Интегрирование методом Симпсона с оценкой точности

Собственное значение определенного интеграла

находится методом Симпсона (парабол). Отрезок [a, b] разбивается на n=2m частей x₀ =a, x₁ =a+h, ..., x_n =b с шагом h=(b-a)/n. Вычисляются значения y_i  = F(x_i ) функции в точках x_i  и находится значение интеграла по формуле Симпсона:

где

Затем количество точек разбиения удваивается и производится оценка точности вычислений:

Если R_n  > e, то количество точек разбиения удваивается. Значение суммы 2(y₁ +y₂ +...+y_2m-1 ) сохраняется, поэтому для вычисления интеграла при удвоении количества точек разбиения требуется вычислять значения y_i  лишь в новых точках разбиения.

Если нашли ошибку в алгоритме - сообщите!

Реализации алгоритма на различных языках:

Блоксхемы:

Реализация алгоритма на AlgoPascal:

unit IntegralSimpsUnit;

declaration
(*
интегрируемая функция
*)
function F(X:Real):Real;

interface
    IntegralSimps;
implementation

(*********************************************************
Интегрирование методом Симпсона с оценкой точности.

Считается интеграл функции F на отрезке [a,b] с погрешностью
порядка Epsilon.

function IntegralSimps(a:Real;b:Real;Epsilon:real):real;
*********************************************************)
function IntegralSimps(
    const   a   :   Real;
    const   b   :   Real;
    const   Epsilon:Real):Real;
var
    i   :   Integer;
    n   :   Integer;
    h   :   Real;
    s   :   Real;
    s1  :   Real;
    s2  :   Real;
    s3  :   Real;
    x   :   Real;
begin
    s2:=1;
    h:=b-a;
    s:=F(a)+F(b);
    repeat
        s3:=s2;
        h:=h/2;
        s1:=0;
        x:=a+h;
        repeat
            s1:=s1+2*F(x);
            x:=x+2*h;
        until not(x<b);
        s:=s+s1;
        s2:=(s+s1)*h/3;
        x:=AbsReal(s3-s2)/15
    until not(x>Epsilon);
    Result:=s2;
end;

end.

Бочканов Сергей, Быстрицкий Владимир
Copyright © 1999-2004
При поддержке проекта MANUAL.RU