Приближенное вычисление определенного интеграла

В настоящем пункте рассматриваются способы приближенного вычисления определенных интегралов

Введем на [а, в] равномерную сетку с шагом h, т.е. множества точек

и представим интеграл в виде суммы интегралов по частным отрезкам:

Для построения формулы численного интегрирования на всем отрезке [а, в] достаточно построить квадратичную формулу для на частном отрезке [х_i-1, х_i].

Формула прямоугольников. Заменим интеграл S_i выражением Геометрический такая замена означает, что площадь криволинейной трапеции АВСД заменяется площадью прямоугольника АВС₁Д₁ (см. рис. 1).

Рис. 1

Тогда получим формулу

(7)

которая называется формулой прямоугольников на частичном отрезке [х_i-1, х_i].

Погрешность метода (7) определяется величиной

которую легко оценить с помощью формулы Тейлора. Действительно, запишем ψ _i в виде

и воспользуемся разложением

Обозначая оценим ψ _i следующим образом:

Таким образом, для погрешности формулы прямоугольников на частичном отрезке справедлива формула

т.е. формула имеет погрешность О(h³) при h→ 0.

Суммируя равенства (26) по I от 1 до N, получим составную формулу прямоугольников

Погрешность этой формулы

Отсюда, обозначая получим

т.е. погрешность формулы прямоугольников на всем отрезке есть величина О(h²). В этом случае говорят, что квадратурная формула имеет второй порядок точности.

Определение. Приближенное равенство

.

Называется квадратурной формулой.

Формула трапеции. На частичном отрезке (х_i-1, x_i) площадь криволинейной трапеции АВСД заменяется площадью прямоугольной трапеции АВСД (рис. 41).

Рис. 41

Тогда

Для оценки погрешности

Представим его в виде

Отсюда получим

Составная формула трапеции имеет вид

Погрешность этой формулы оценивается следующим образом:

Таким образом, формула трапеции имеет вид, так же как и формула прямоугольников, второй порядок точности, но ее погрешность оценивается величиной в два раза большей.

Применение формулы трапеции или прямоугольников требует оценки второй производной на отрезке [а, в]. Если такая оценка затруднительна (или вообще невозможно, например, в случае функции определяемых опытным путем), то в предположений малого изменения (или монотонности) второй производной можно во всех полученных оценках заменить множителя М₂h² наибольшей величиной

Отсюда видно, что формула прямоугольников и трапеции дает достаточную точность только при достаточно малых разностях второго порядка (а именно, когда произведения не превосходят допустимой погрешности расчета).

Для уточнения величины интеграла можно использовать, то обстоятельство, что с уменьшением шага h в два раза погрешность формулы трапеций уменьшается примерно в четыре раза. Отсюда следует, что совпадающие знаки в значениях интеграла, вычисленных с шагом h и можно считать верным. Действительно, если погрешность значения интеграла, вычисленного с шагом обозначить через ε, то погрешность значения интеграла, вычисленного с шагом h, будет приближенно равна 4ε, и значить, разность указанных значений интеграла будет не менее чем 3ε. Поэтому из совпадения m десятичных знаков у рассматриваемых значений интеграла можно заключить, что погрешность , а это означает, что в значений интеграла вычисленном с шагом , все m десятичных знаков верны (здесь предполагается, что погрешность исходных данных пренебрежимо мало).

Формула Симпсона. При аппроксимации интеграла заменяем функцию параболой, проходящей через точки , , т.е. представим приближенно в виде

Тогда

(8)

Вычислим

Из (8) получим, что

Таким образом, приходим к приближенному равенству

которое называется формулой Симпсона.

Погрешность этой формулы ψ _i оценивается так [1]:

На всем отрезке [a, в] формула Симпсона имеет вид

Погрешность этой формулы оценивается неравенством:

Из этой оценки видно, что с уменьшением шага h в два раза погрешность формулы Симпсона уменьшается примерно в 16 раз; поэтому значение интеграла, вычисленное с шагом содержащий на один верный знак больше, чем значение интеграла, вычисленное с шагом h. Это правило на практике очень удобно при оценке точности интеграла.

Программа

program arai5;

uses crt;

const T0=307; Tk=414; T1_=407; n=60; Dt=(T1_-T0)/n;

E=-0.93; AA=96.5; a0=2.7;

Var xx, xp, tt, tx, d: array [0..n] of real;

i: integer; f1, f3: text;

Ti, Ma, Si, k1, k2, k3, k4, y_: real;

function F(T: real): real;

begin F: =exp(AA/T-E)/(Tk-T); end;

function k(T: real): real;

begin k: =exp(E-AA/T); end;

function Fx(T, x: real): real;

begin fx: =exp(E-AA/T)*(a0-x); end;

{Vy4islenie integrala, formula simpsona}

begin TT[0]: =0; tx[0]: =0;

for i: =0 to n-1 do begin Ti: =T0+I*Dt;

TT[i+1]: =TT[i]+(F(ti)+4*F(ti+dt/2)+F(ti+dt))*dt/6; end;

{Formula pryamougolnikov}

for i: =0 to n-1 do

Tx[i+1]: =Tx[i]+F(T0+I*Dt+Dt/2)*dt;

Ma: =0; Si: =0;

for I: =1 to n do D[i]: =Abs(TT[i]-TX[i]);

for I: =1 to n do Ma: =Ma+D[i]/n;

for i: =1 to n do Si: = Si+sqr(D[i]-Ma)/n;

Si: =sqrt(si);

writeln('n=', n, ' Mt=', Ma: 6: 6, ' tsi=', si: 6: 6);

assign(f1, 'out1.dat');

rewrite(f1);

for i: =0 to n do

writeln(f1, Tt[i]: 6: 4, ' ', Tx[i]: 6: 6);

close(f1);

{reshenie zada4i koshi-kinetika}

xx[0]: =0; xp[0]: =0;

{metodom R-K 4-go p}

for i: =0 to n-1 do begin

k1: =Fx(T0+i*dt, xx[i]);

k2: =Fx(T0+i*dt+dt/2, xx[i]+k1*dt/2);

k3: =Fx(T0+i*dt+dt/2, xx[i]+k2*dt/2);

k4: =Fx(T0+i*dt+dt, xx[i]+k3*dt);

xx[i+1]: =xx[i]+(k1+2*k2+2*k3+k4)*(TT[i+1]-TT[i])/6; end;

{R-K 2go p}

For i: =0 to n-1 do begin

y_: =xp[i]+0.5*(TT[i+1]-TT[i])*Fx(T0+i*dt, xp[i]);

xp[i+1]: =xp[i]+(TT[i+1]-TT[i])*Fx(T0+i*dt+dt/2, y_); end;

Ma: =0; Si: =0;

for I: =1 to n do D[i]: =Abs(xx[i]-xp[i]);

for I: =1 to n do Ma: =Ma+D[i]/n;

for i: =1 to n do Si: = Si+sqr(D[i]-Ma)/n;

Si: =sqrt(si);

writeln('n=', n, ' Mx=', Ma: 6: 6, ' xsi=', si: 6: 6);

assign(f3, 'out3.dat');

rewrite(f3);

for i: =0 to n do

writeln(f3, xx[i]: 6: 4, ' ', xp[i]: 6: 6);

close(f3);

end.

Вывод: