§ 2. Численное интегрирование по простейшим формулам
(((2*q-9)*q+11)*q-3)*r4[i]/12;
2: new1 := r2[i]+(q-1)*r3[i]+
((6*q-8)*q+11)*r4[i]/12;
end;
end.
{***Вторая интерполяционная формула Ньютона ***}
function new2 (key:integer;x,y,r1,r2,r3,r4:mas;
x1:real):real;
label 30;
var i,j,k : integer; q : real;
begin
i := 1;
j := 14;
if x1 > x[j] then goto 30;
{ **** Метод двоичного поиска интервала ****}
repeat
k := (i+j) div 2;
if x1< x[k] then j := k;
if x1>=x[k] then i := k;
until j <= i+1;
30:
q := (x1 - x[i]) / (x[2]-x[1]);
{ **** Выбор типа работы процедуры ****}
case key of
0: new2 := y[j]+q*(r1[j-1]+
(q+1)*(r2[j-2]/2.0+(q+2)*r3[j-3]/6.0));
1: new2 := r1[j-1]+0.5*(2*q+
1)*r2[j-2]+((3*q+6)*q+2)*r3[j-3]/6+
(((2*q+9)*q+11)*q+3)*r4[j-4]/12;
2: new2 := r2[j-2]+(q+
1)*r3[j-3]+((6*q+18)*q+11)*r4[j-4]/12;
end;
end.
Для проверки и тестирования процедур вычислялись значения первой и второй производных функции при заданных значениях аргумента с использованием первой или второй интерполяционных формул Ньютона. Расчеты контролировались при составлении таблицы разностей. Вычисления выполнялись с точностью 10-5. Данные взяты из табл. 2.2. Результаты работы процедур приводятся в табл. 3.1.
Таблица 3.1
Первая инт.формула Ньютона |
Вторая инт.формула Ньютона |
||||
xi |
y’ |
y’’ |
xi |
y’ |
y’’ |
1.90800 |
-0.01067 |
-0.00043 |
2.24800 |
-0.01679 |
-0.00041 |
2.13500 |
-0.01449 |
-0.00042 |
2.35900 |
-0.01999 |
0.00025 |
§ 2. ЧИСЛЕННОЕ ИНТЕГРИРОВАНИЕ ПО ПРОСТЕЙШИМ ФОРМУЛАМ
.
Bсякая простая формула, аппроксимирующая отдельный интеграл, назывaется квадратурной. Составная квадратурная формула - это формула, дающая приближение интеграла в виде суммы приближений интегралами по данной квадратурной формуле.
Двумя простейшими квадратурными формулами являются формулы прямоугольников и трапеций, которые в ряде случаев оказываются наиболее эффективными.
Известны три разновидности формул прямоугольников: это формулы левых, правых и средних прямоугольников. Все они основаны на аппроксимации каждого интеграла площадью прямоугольника, одной из сторон которого является hi , а второй - либо значение функции на левом конце отрезка (рис.3.1, а), либо значение функции на правом конце отрезка (рис. 3.1, б), либо значение функции в средней точке отрезка (рис. 3.1, в).
Квадратурные формулы, аппроксимирующие , будут иметь вид:
левых прямоугольников: = hi f (хi);
правых прямоугольников: = hi f (хi+1);
средних прямоугольников: = hi f (хi+1/2).
С учетом представления на элементарном отрезке составные квадpатурные формулы прямоугольников могут быть записаны так:
левых прямоугольников
;
правых прямоугольников
;
средних прямоугольников
.
Тогда площадь фигуры может быть определена из формулы площади прямоугольной трапеции
Si = (fi + fi+1) hi /2.
Если теперь просуммируем последнюю формулу по всем элементарным отрезкам, то получим с учетом выполненных элементарных преобразований следующее выражение:
Заметим, что при бесконечном уменьшении длин элементарных отрезков формулы обоих типов (прямоугольников и трапеций) сходятся к точному значению интеграла. Однако не ясно, как быстро они сходятся ? Попытаемся выяснить данный вопрос, воспользовавшись разложением функции в ряд Тейлора относительно центра элементарного отрезка [хi, хi+1]
f(x) = f(yi) + (x - yi) f ’(yi) + (x - yi)2 f ’’(yi)/2 +
+ (x - yi)3 f’’’(yi)/6 + (x-yi)4 f (IV)(yi)/24 + ...
затем, проинтегрировав полученный ряд по каждому из отрезков в предположении, что оставшиеся члены ряда намного меньше выписанных, с учетом значений коэффициентов ряда Тейлора на элементарном отрезке
получим
Член показывает ошибку формулы прямоугольников без учета членов более высокого порядка. Если теперь подставим в формулу трапеций значения функции в точках х = хi и х = хi+1, то получим
Можно видеть, что ошибки формул средних прямоугольников и трапеций одного порядка, т.е. дают почти одинаковую точность.
Примеров вычисления определенных интегралов по простейшим формулам (прямоугольников и трапеций) можно привести множество. Ни одна БСП не обходится без этих процедур. Одним из примеров такой вычислительной процедуры может служить следующая:
procedure intlps (var s:real; k:integer);
begin
s := 0.0;
x := a;
for i := 0 to n do
begin
if k <>3 then x := a+ h*i
else x:=a+h/2+h*i;
f := func (a,b,x);
if k=4 then f := f/2.0;
if (i=0) and (k<>2) then s := s + f;
if (i=n) and ((k=2) or (k=4)) then s := s + f;
if (i<>0) and (i<>n) then
if k<>4 then s := s+f
else s := s+f*2;
end; s := s * h;
end.
Эта процедура обращается к процедуре-функции func для вычисления подынтегральной функции в заданных