IPB
ЛогинПароль:

 
 Ответить  Открыть новую тему 
> Численное интегрирование
trminator
сообщение 13.12.2003 12:33
Сообщение #1


Четыре квадратика
****

Группа: Пользователи
Сообщений: 579
Пол: Мужской

Репутация: -  4  +


Численное интегрирование методом прямоугольников:
Код
program integral;
 {        Вычисляет приближенное значение       }
 { интеграла функции F методом прямоугольников  }
 { /b                                           }
 { |  f(x) dx;интервал разбивается на n частей  }
 { /a                                           }

var
 n,i : integer;
 a,b,shag,sum,itog : real;

{======================================}
{Введите сюда нужную функцию}
function F(x:real):real;
 begin
   F:=x*x*x*x+2*x*x+4      
 end;
{======================================}

BEGIN
 write('Начало интегрирования a = '); readln(a);
 write('Конец  интегрирования b = '); readln(b);
 write('Количество разбиений интервала n = '); readln(n);

 shag:=(b-a)/n;
 sum:=0;
 for i:=1 to n-1 do
   sum := sum + F(shag*i+a);
 sum := sum + (F(a)+F(b))/2;

 itog:=(b-a)/n * sum;
 writeLn('Интеграл = ', itog:0:5)
END.


--------------------
Закон добровольного труда Зимерги:
Люди всегда согласны сделать работу, когда необходимость в этом уже отпала
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
hiv
сообщение 18.03.2005 11:59
Сообщение #2


Профи
****

Группа: Пользователи
Сообщений: 660
Пол: Мужской
Реальное имя: Михаил

Репутация: -  11  +


Численное интегрирование методом Симпсона (четное количество разбиений интервала интегрирования):
Код
program Simpson;

{интегрируемая функция}
function F(x:Real):Real;
begin
 F:=2*x;
end;

var  a,b,h,x :real;
    n,i   :integer;
    integ :real;

begin
 write('Введите начало интегрирования a='); readln(a);
 write('Введите конец  интегрирования b='); readln(b);
 write('Введите количество разбиений интервала (четное число) n='); readln(n);
 if (n mod 2)>0 then
 begin
   n:=n+1;
   writeln('Число n было введено нечетное, оно было заменено на n=',n);
 end;

 h:=(b-a)/n;
 integ:=F(a)+F(b)+4*F(a+h);
 for i:=1 to (n div 2)-1 do
 begin
    x:=a+2*h*i;
    integ:=integ+2*F(x)+4*F(x+h);
 end;
 integ:=h*integ/3;

 writeln('Интеграл = ',integ);
end.


--------------------
Никогда не жадничай. Свои проблемы с любовью дари людям!
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
hiv
сообщение 18.03.2005 14:01
Сообщение #3


Профи
****

Группа: Пользователи
Сообщений: 660
Пол: Мужской
Реальное имя: Михаил

Репутация: -  11  +


Самый лучший и быстрый метод интегрирования - десятиточечный метод Гаусса. В программе используется рекурсивная функция для достижения заданной точности вычислений.
Код
program gaussint;
{   Вычисление интегpала десятиточечным методом Гаусса   }

{константы десятиточечного метода Гаусса}
const
 g10c1=0.9739065285/6.2012983932;
 g10c2=0.8650633667/6.2012983932;
 g10c3=0.6794095683/6.2012983932;
 g10c4=0.4333953941/6.2012983932;
 g10c5=0.1488743390/6.2012983932;
 g10x1=0.0666713443/6.2012983932;
 g10x2=0.1494513492/6.2012983932;
 g10x3=0.2190863625/6.2012983932;
 g10x4=0.2692667193/6.2012983932;
 g10x5=0.2955242247/6.2012983932;

function F(x:real):real;  {интегрируемая функция}
 begin
   F:=pi*sin(pi*x);
 end;

function gauss_calc(a,b:real):real; {сам десятиточечный метод Гаусса}
 var  n,m,s,s1,s2,s3,s4,s5  :real;
 begin
   m:=(b+a)/2; n:=(b-a)/2;
   s1:=g10c1*(f(m+n*g10x1)+f(m-n*g10x1));
   s2:=g10c2*(f(m+n*g10x2)+f(m-n*g10x2));
   s3:=g10c3*(f(m+n*g10x3)+f(m-n*g10x3));
   s4:=g10c4*(f(m+n*g10x4)+f(m-n*g10x4));
   s5:=g10c5*(f(m+n*g10x5)+f(m-n*g10x5));
   s:=s1+s2+s3+s4+s5;
   gauss_calc:=s*(b-a);
 end;

{рекурсивная ф-ция подсчета с заданной точностью}
{ gc - ранее посчитаный интеграл на интервале (a,b)}
function gauss(a,b,eps,gc:real):real;  
 var  t,ga,gb :real;                    
 begin
   t:=(a+b)/2;                 {разбиваем интервал на две половинки}
   ga:=gauss_calc(a,t);        {в каждой половинке считаем интеграл}
   gb:=gauss_calc(t,b);
   if abs(ga+gb-gc)>eps then   {проверяем точность вычислений}
     begin
       ga:=gauss(a,t,eps/2,ga);  {рекурсия для первой половинки}
       gb:=gauss(t,b,eps/2,gb);  {рекурсия для второй половинки}
     end;                        {при этом точность повышаем, чтобы }
                                 {при сложении ошибка не накапливалась}
   gauss:=ga+gb;               {интеграл = сумме интегралов половинок}
 end;

var
 Integral :real;
 a,b,eps :real;

begin
 write('Введите начало интервала интегрирования a='); readln(a);
 write('Введите конец  интервала интегрирования b='); readln(b);
 write('Введите точность интегрирования eps='); readln(eps);

 Integral:=gauss(a,b,eps,gauss_calc(a,b));
 writeln('Интеграл = ',Integral);
end.


--------------------
Никогда не жадничай. Свои проблемы с любовью дари людям!
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 
Altair
сообщение 24.05.2005 20:58
Сообщение #4


Ищущий истину
******

Группа: Модераторы
Сообщений: 4 824
Пол: Мужской
Реальное имя: Олег

Репутация: -  45  +


Вычисления определенных интегралов по формуле центральных прямоугольников

UNIT C_Rect;

INTERFACE

TYPE Func=Function(x:real):real;

FUNCTION Rectangles(a,b:Real; f:Func; n:word):Real;
FUNCTION RectanglesRunge(a,b:Real; f:Func; var n: word;eps:Real):Real;

IMPLEMENTATION

FUNCTION Rectangles(a,b:Real; f:Func; n:word):Real;

VAR
h,f1,sum,x,y:Real;
i:Integer;

BEGIN

h := (b-a) / n;
sum := 0;
x := a;
for i := 1 to n do
begin
y := x;
x := x + h;
f1 := f((x + y)/2);
sum := sum + f1;
end;
Rectangles:= sum*h;
END;

FUNCTION RectanglesRunge(a,b:Real; f:Func; var n: word; eps:Real):Real;

VAR
I1, I2: real;
BEGIN
I1 := Rectangles(a,b,f,n);
n := n+n;
I2 := Rectangles(a,b,f,n);
while abs(I1-I2) > 2*eps do
begin
if n >= 16383 then break;
I1 := I2;
n := n+n;
I2 := Rectangles(a,b,f,n);
end;
RectanglesRunge := I2;
END;

BEGIN
END.



ОПИСАНИЕ МОДУЛЯ.
Прикрепленный файл  C_Rect.doc ( 59 килобайт ) Кол-во скачиваний: 2281


--------------------
Помогая друг другу, мы справимся с любыми трудностями!
"Не опускать крылья!" (С)
 Оффлайн  Профиль  PM 
 К началу страницы 
+ Ответить 

 Ответить  Открыть новую тему 
2 чел. читают эту тему (гостей: 2, скрытых пользователей: 0)
Пользователей: 0

 



- Текстовая версия 9.11.2024 21:06
Хостинг предоставлен компанией "Веб Сервис Центр" при поддержке компании "ДокЛаб"