Вычисление определённых интегралов

Пусть задан интеграл S.

картинка

Значение интеграла численно равно значению площади криволинейной траппеции, ограниченной кривой y=f(x) (непрерывная функция) и прямыми x=a и x=b.

картинка

Разделим отрезок [a,b] на n равных частей с точками X0=a, X1, X2 … Xn=b. Шаг интегрирования (приращение) dx = (b-a)/n. Вычислим значение функции в каждой точке деления Y0=F(X0), Y1=F(X1), .. Yn=F(Xn).

картинка

Заменим каждую полоску криволинейной трапеции прямоугольником с высотой, равной значению функции на левом конце промежутка. (На рисунке красным цветом показаны дополнительные построения). Площадь криволинейной фигуры заменится площадью сумм прямоугольников и будет равна

S = Y0*h+Y1*h+..+Yn-1*h = ((b-a)/n)*(Y0+Y1+...Yn-1)

картинка

Если же каждую полоску криволинейной трапеции заменим прямоугольником с высотой, равной значению функции на правом конце промежутка, то формула вычисления площади получит вид

S = Y1*h+Y2*h+..+Yn*h = ((b-a)/n)*(Y1+Y2+...Yn)

{ Вычисление интеграла методом средних прямоугольников }
var
  n : integer;
  s: real;
 
function f(x:real):real;
begin
  f:=1-0.25*sin(x)*sin(x);
end;
 
procedure rect (a,b:real; n:integer; var s:real);
var
  i : integer;
  h,x: real;
begin
  h:=(b-a)/n;
  x:=a+h/2;
  s:=0.0;
  for i:=1 to n do
    begin
      s:=s+f(x);
      x:=x+h;
    end;
  s:=s*h;
end;
 
begin
  rect (0,pi/2,40, s);
  write('Интеграл = ',s:16:4);
end.
{ Ответ:            }
{ Интеграл = 1.3744 }

картинка

Если каждую полоску криволинейной трапеции заменить обычной трапецией, то площадь фигуры будет равна сумме площадей трапеций:

S = h*((Y0+Yn)/2 +Y1+Y2+...+Yn)

картинка

Для вычисления этим методом необходимо, чтобы число частей было четным. Через каждые три точки (Xi,Yi), (Xi+1,Yi+1), (Xi+2,Yi+2) проведем пораболу и заменим соответствующий кусок криволинейной трапеции на параболическую трапецию. Получим формулу:

S = ((b-a)/(3*n))*(Y0+Yn+4*(Y1+Y3+...Yn-1)+2*(Y2+Y4+..Yn))

uses crt;
 
{Функция, площадь которой нужно вычислить}
function f(x: real): real;
begin
 f := x*x+3;
end;
 
var a,b,dx,x,y,y1,y2: real;
    n: integer;
begin
 clrscr;
 
 {Вводим промежуток}
 write ('Введите a:'); readln (a);
 write ('Введите b:'); readln (b);
 
 {Вводим количество частей}
 write ('Введите n:'); readln (n);
 
 {Определяем шаг интегрирования}
 dx := (b-a)/n;
 
 {Вычисление по формуле левых прямоугольников}
 y := 0; x := a;      {x = a = X0}
 while x < b do begin {x < Xn}
  y := y + f(x);
  x := x + dx;        {x = X0, X1, X2 ... Xn-1}
 end;
 y := y * dx;
 writeln ('Площадь по формуле левых прямоугольников: ', y: 10: 3);
 
 {Вычисления по формуле правых многоугольников}
 y := 0; x := a + dx;  {x = X1}
 while x <= b do begin {x <= Xn}
  y := y + f(x);
  x := x + dx;         {x = X1, X2, X3 ... Xn}
 end;
 y := y * dx;
 writeln ('Площадь по формуле правых прямоугольников: ', y:10:3);
 
 {Вычисление по формуле трапеций}
 y := 0; x := a + dx;  {x = X1}
 while x < b do begin  {x < Xn}
  y := y + f(x);
  x := x + dx;         {x = X1, X2, X3 ... Xn-1}
 end;
 y := (y + (f(a)+f(b))/2) * dx;
 writeln ('Площадь по формуле трапеций: ', y:10:3);
 
 {Вычисление по методу Симпсона}
 y1 := 0; x := a + dx;       {x = X1}
 while x < b - dx do begin   {x < Xn}
  y1 := y1 + f(x);
  x := x + 2*dx;             {x = X1, X3, X5 ... Xn-1}
 end;
 y2 := 0; x := a + 2*dx;     {x = X2}
 while x < b - 2*dx do begin {x < Xn-1}
  y2 := y2 + f(x);
  x := x + 2*dx;             {x = X2, X4, X6 ... Xn-2}
 end;
 y := ((b-a)/(3*n)) * (f(a)+f(b) + 4*y1 + 2*y2);
 writeln ('Площадь по формуле Симпсона: ', y:10:3);
 
 readkey;
end.