Форум «Всё о Паскале» _ FAQ _ Численное интегрирование
Автор: trminator 13.12.2003 12:33
Численное интегрирование методом прямоугольников:
Код
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;
Численное интегрирование методом Симпсона (четное количество разбиений интервала интегрирования):
Код
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.
Автор: hiv 18.03.2005 14:01
Самый лучший и быстрый метод интегрирования - десятиточечный метод Гаусса. В программе используется рекурсивная функция для достижения заданной точности вычислений.
Код
program gaussint; { Вычисление интегpала десятиточечным методом Гаусса }
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);
Вычисления определенных интегралов по формуле центральных прямоугольников
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 килобайт )
Кол-во скачиваний: 2249