Інтеграли зі змінними границями

2



Міністерство освіти  і науки України

Дніпропетровський національний університет



Механіко-математичний факультет

Кафедра обчислювальної механіки і міцності конструкцій





























Курсова робота 



з чисельних методiв

на тему:



Інтеграли зі змінними границями











                                                                                                            



    Виконавець      студент групи  МД-01-1                 Ромащук Р. В.

    Керівник          старший викладач                            Гарт Е.Л.                                                                                                   























Дніпропетровськ

2003 р.

     Ця курсова робота мiстить в собi такi теоретичнi питання, як « Визначений інтеграл зі змінною верхньої межею. Властивості визначеного інтегралу зі змінною верхньої межею. Чисельні методи знаходження визначеного інтегралу зі змінною верхньої межею», розв’язок за допомогою обчислювальної машини задачi для знаходження визначеного інтеграла зі змінними границями інтегрування, а також наведенi висновки, на основi отриманих результатiв.



























































































З М I С Т

Постановка задачi ………………………………………………………………………........4

		 Вступ…………………………………………………………………………………………....5

  1.  Постановка задачі чисельного інтегрування.............……............................6

  2.  Квадратурні формули………...........................................................................6

         2.1.  Формула прямокутників.......................................................................6

         2.2.  Формула трапецій..................................................................................7

         2.3.  Формула парабол (Сімпсона)...............................................................9  

 3.  Чисельні методи знаходження визначеного 

       інтеграла зі змінною верхньою межею.........................................................10

 4. Опис обчислювального алгоритму………………………………………….10

 5. Обговорювання результатів…………………………………………………11

Висновки…………………………………………………………………………12

Список посилань………………………………………………………………...13

Додатки:………………………………………………………………………….14

               А  Опис вихiдних даних та результатiв розрахунку………………...14

               В  Схема обчислювального алгоритму……………………………….15

               С  Лiстiнг програми…………………………………………………....18

































































Постановка задачі



     За допомогою квадратурних формул обчислити визначений інтеграл зі мінною границею 

                                                                      (1)



Побудувати сітку, і скласти таблицю значень інтеграла на цій сітці fn=f(x)

За квадратурною формулою високої точності. Тоді

                                            xnxxn+1





































































В С Т У П



     В практичних розрахунках, у т.ч. в задачах механіки, нерідко виникає необхідність обчислення визначених інтегралів

    

 де під інтегральна функція f(x) неперервна відрізку [a,b], а вагова функція (x) неперервна на інтервалі (a,b).

     До чисельного знаходження інтеграла звертаються тоді, коли його або неможливо виразити через елементарні функції, або підінтегральна функція задана таблично, а також коли внаслідок інтегрування приходять до незручного для використання виразу. Формули чисельного знаходження визначених інтегралів називаються квадратурними формулами. Побудова квадратурних формул ґрунтується на заміні складної підінтегральної функції деякою більш простою, інтеграл від якої легше обчислити. Виникаюча при цьому похибка називається похибкою квадратурної формули. Най простіші квадратурні формули можуть бути отримані із простих геометричних міркувань. 



























Постановка задачі чисельного інтегрування

     Нехай потрібно знайти визначений інтеграл

                                                                                                    (1.1)

де функція f(x) неперервна відрізку [a,b], а вагова функція (x) неперервна на інтервалі (a,b). Тоді f(x) наближають такою функцією (x;C) від якої інтеграл легко взяти в елементарних функціях. Завдяки лінійності такої апроксимації відносно параметрів ci функцію можна записати так:

                                                                                        (1.2)

де r(x) – залишковий член апроксимації. Підставляємо (1.2) в (1.1), отримаємо загальну формулу чисельного інтегрування – квадратурну формулу:



;          

де хi  - вузли, сi – ваги, R – залишковий член. Інтеграл приблизно заміняється сумою, схожою на інтегральну суму, причому вузли та коефіцієнти цієї суми не залежать від f(x).  



2.  Квадратурні формули.

2.1.  Формула прямокутників.

Припустимо, що fC2[-h/2,h/2], h>0 . 

                                                 		           	(2.1.1)

де f0=f(0), тобто площа криволінійної трапеції, обмеженої зверху графіком функції f(x) , апроксимується площею прямокутника, висота якого дорівнює значенню f(x)  в середній точці трапеції (мал. 2.1.1).



мал. 2.1.1. Формула прямокутників



	Знайдемо залишковий член , тобто похибку формули (2.1.1) . 	Нехай 

                                                       				(2.1.2)

  Тому що F(0)=0, F/(0)=f0, F//(0)=f/0, F///(x)=f//0,   

то відповідно до формули Тейлора з залишковим членом у формі Лагранжа маємо 

	                (2.1.3)

де - , +  - деякі точки , -h/-<+<h/2.  

	Функція F(x) є первісної для f(x). Тому для інтеграла, що стоїть в лівій частині наближеної рівності (2.1.1), з формули Ньютона-Лейбница з розрахунком (2.1.3) випливає наступне співвідношення 



     Ззвідси одержуємо формулу прямокутників із залишковим членом:

	                                                  (2.1.4)



2.2.  Формула трапецій.

Нехай fC2[0,h], h>0   

	                                              				(2.2.1)

де f0=f(0), f1=f(h)  тобто інтеграл  приблизно заміняється площею заштрихованої трапеції,  показаної  на малюнку (мал. 2.2.1).

  

мал. 2.2.1. Формула трапецій.



	Знайдемо залишковий член, тобто похибку формули (2.2.1). Виразимо f1 та F1=F(h) де F - функція (2.1.2), по формулі Тейлора з залишковим членом в інтегральній формі (*):

    (*)

	                                    			 (2.2.2)

(2.2.3)

Згідно (2.2.1) маємо

	           		(2.2.4)

Відокремивши в правій частині (2.2.3) доданок hf0/2  і замінивши його вираженням (2.2.4), з урахуванням того, що 



 знаходимо 



	Перетворимо тепер другий доданок у правій частині, використовуючи узагальнену  теорему про середнє. 	Тому що (h-t)t0,  t[0,t]  то за  теоремою 



де [a,b]  - деяка точка . Підставляючи отримане в (*), приходимо  до  формули трапецій із залишковим членом :

	                             (2.2.5)

	



2.3.  Формула Сімпсона .

Припустимо, що fC4[-h,h]. Тоді інтеграл



наближеного заміняємо площею заштрихованої криволінійної трапеції, обмеженою зверху  параболою, що проходить через точки (-h,f-1), (0,f0), (h,f1),   де  fi=f(ih) (мал. 2.3.1) 



мал. 2.3.1 Формула парабол (Сімпсона)



	Зазначена парабола задається рівнянням 



у цьому неважко переконатися, поклавши x=-h, x=0, x=h  (її можна також одержати, побудувавши інтерполяційний багаточлен другого ступеня і приводячи подібні ).	Звідси знаходь 



     Таким чином, формула Сімпсона , називають також формулою парабол, має вид 

                                 		           	 (2.3.1)

     Покладемо F1=F(h), де F функція (2.1.2). Оскільки F(0)=0, F(k)(x)=

f(k-1 ) (x), 1k5 то згідно формули Тейлора з залишковим членом в інтегральній формі маємо



Звідси одержуємо 

		(2.3.2)

тому що інші члени взаємно знищуються.

     Оскільки   , t[0,h] то застосовуючи до інтеграла (2.3.2) узагальнену теорему про середнє, знаходимо 



                                                      (2.3.3)

де [0,h], [-h,h] - деякі точки. Приймаючи до уваги, що 



з (2.3.2), (2.3.3) приходимо до формули

	                                 		           (2.3.4) 

тобто  до  формули Сімпсона з залишковим членом.	





3.  Чисельні методи знаходження визначеного

Інтеграла зі змінною верхньою межею



     У деяких випадках необхідно обчислити такі інтеграли



Можна, звичайно, розглядати його для кожного значення верхньої границі х як інтеграл зі сталими границями і обчислювати за однією з квадратурних формул, що невигідно у випадку великої кількості значень x. Краще вибрати деяку сітку і скласти таблицю значень інтеграла на цій сітці Fn=F(x) за квадратурної формули високої точності. Тоді 

                                                    (3.1)

причому останній інтеграл можна одчислювати за простими квадратурними формулами. 

     Окрім того, маючи таблицю F(xn), можна знаходити  F(x) інтерполяцією за цією таблицею. Природно, маючи і похідну інтеграла F (x)=(x)f(x). Краще скористатись інтерполяційним поліномом Ерміта.  







4. Опис обчислювального алгоритму



     При реалізаціі алгоритму обчислення визначеного інтеграла зі змінними границями інтегрування використовуються процедури та функцiї, для того щоб скоротити витрати машинного часу при обчислюваннi, та для компактностi программи. Программа для знаходження написана на мовi Delphi5, стан пограмми – вiдлажена.  



Обговорювання результатів

                                                                                                                   Таблиця 1

     



Формула (3.1)

Формула Сімпсона

Формула трапецій

Дійсне значення інтеграла

a=0;  b=1;

-0.7974398040

Різниця 0.0000012883

		-0.7974386790

		Різниця 0.0000001633

-0.7993252434

Різниця 0.00188672780

-0.7974385156

a=0;  b=2;

3.9190337956      Різниця 0.0000062805

3.9190353338     Різниця 0.0000047422

3.90875628130     Різниця  0.01028379486

3.9190400761

a=0;  b=3;

10.5498688094      Різниця 0.00002744251

10.5498688094     Різниця 0.00002744251

10.5247085565     Різниця 0.02518769537

10.5498962519

a=0;  b=4;

17.8842287345      Різниця 0.0000804723

17.8842201707     Різниця 0.00008903613

17.8382724576     Різниця 0.0460367491

17.8843092068

a=0;  b=5;

25.5043003647      Різниця 0.0001835185

25.5042688642     Різниця 0.00021501907

25.4318420115     Різниця 0.0726418717

25.5044838833

a=0;  b=6;

33.2576007639      Різниця

0.00035637138

33.2575244054    Різниця 0.00043272988

33.1529684530     Різниця 0.1049886822

33.2579571352



     Таблиця 1 була отримана при наступних вхідних даних:

Кількість вузлів при побудові таблиці значень інтегралу (1) =20

Кількість вузлів при застосуванні формули трапецій =20

Кількість вузлів при застосуванні формули Сімпсона =20





























Висновки



Таким чином з таблиці 1 видно, що чим більший проміжок ми беремо тим кращу точність отримаємо, навіть краще за формулу Сімпсона, але загальна похибка (відносно дійсного значення) також збільшується. Формулу (1) доцільно використовувати, якщо потрібно обчислити інтеграл на відносно великому проміжку та якщо треба обчислити відразу декілька інтегралів.  















































































Список посилань



Каліткін Н.Н. ‘Чисельні методи’ – М.: Наука, 1978. – 512 с.



Балашова С.Д. ‘Тексты лекций по курсу “ Численные методы”’. – Днепропетровск: Из – во ДГУ, 1989. – 206 с



Мусiяка В.Г. Основи чисельних методiв механiки. – Днiпропетровськ: Вид – во ДДУ, 1993. – 156 с.



 Методические рекомендации по курсу “ Методы вычислений в инженерных расчётах”/ Составитель В.Г. Мусияка. – Днепропетровск: Из – во ДГУ, 1992. – 40 с.



Фіхтегольц Г.М. ‘Основи математичного аналізу’– М.: Наука, 1968. – 440 с.  





























































Д О Д А Т К И



А Опис вихiдних даних та результатiв розрахунку



Вихiднi данi  

Кількість вузлів при побудові таблиці значень інтегралу (1)  nGrid - integer;

Кількість вузлів при застосуванні формули трапецій nTrap - integer;

Кількість вузлів при застосуванні формули Сімпсона nSim - integer;

Границі інтегрування a і b – real;



Наслiдки виконання програми друкуються у виглядi:



Вихідні дані це функції типа real.

FullIntegral(L,R);        

integralSimpsona(L,R);

integralTrapeciay(L,R);

first(L,R);

































B Схемаобчислювального алгоритму



Функція     y(x) тип real                             Функція     first(x1,x2) тип real



Початок







Result:=2*arctan(x)-11*Power(3,(-x))+5;







Кінець





Result:=5*x2+2*x2*ArcTan(x2)+11*Power(3,(-x2))/Ln(3)-Ln(1+x2*x2)-(5*x1+2*x1*ArcTan(x1)+11*Power(3,(-x1))/Ln(3)-Ln(1+x1*x1));

Початок

Кінець









Процедура createGrid                                 Процедура createGridOfInt

















Початок

  n:=High(Xgrid);

  SetLength(GridOfInt,n+1);

  for i:=1 to n do

    GridOfInt[i]:=integralSimpsona(a,Xgrid[i]);

Кінець



h:=(b-a)/nGrid;     SetLength(Xgrid,nGrid+1);

  for i:=0 to nGrid do

                Хgrid[i]:=a+h*i;

Початок

Кінець



							             			







Процедура setN(nG, nT, nS )        





















Початок

nG, nT, nS :integer

  nGrid:=nG;  nSim:=nS;

  nTrap:=nT;

Кінець







Функція integralSimpsona(aSim,bSim) тип real













aSim<>bSim

Початок

aSim, bSim : real 

           так                                            ні         







      h:=(bSim-aSim)/nSim;//2

      SetLength(X,nSim+1);

      for i:=0 to nSim do

        X[i]:=aSim+h*i;      sum1:=0;

      sum2:=0;

      for i:=1 to (nSim div 2) do

        sum1:=sum1+y(X[2*i-1]);  

      for i:=1 to ((nSim div 2)-1) do

        sum2:=sum2+y(X[2*i]);

      Result:=(bSim-aSim)*(y(X[0])+y(X[nSim])+4*sum1+2*sum2)

    

Result:=0;

























Кінець







Функція     integralTrapeciay(aTrap,bTrap) тип real



Кінець

aTrap,bTrap : real

Початок













  h:=(bTrap-aTrap)/nTrap;

  SetLength(X,nTrap+1);

  for i:=0 to nTrap do

    X[i]:=aTrap+h*i;

  sum:=(y(X[0])+y(X[nTrap]))/2;

  for i:=1 to (nTrap-1) do

      sum:=sum+y(X[i]);

  Result:=sum*h;



















Кінець





Функція     FullIntegral(aFull,bFull:real) тип real





raznost:=GridOfInt[ai]+integralTrapeciay(Xgrid[ai],aFull);

Result:=GridOfInt[bi]+integralTrapeciay(Xgrid[bi],bFull)-raznost;

bi:=i;

(Xgrid[i]<=bFull)

and

(bFull<=Xgrid[i+1])

ai:=i;

(Xgrid[i]<=aFull) and

(aFull<=Xgrid[i+1])

Кінець

i<(High(Xgrid)-1)

i:=i+1

i:=-1

Result:=0;

aFull<>bFull

Кінець

aFull, bFull : real

Початок









                                            так                                        ні









                      так                                                             ні









    так                                                           ні









                                          так                                                            ні













































С   Лiстiнг програми



unit funct;



interface

uses

  Windows, Messages, SysUtils, Classes, Graphics, Controls, ExtCtrls, Math;



const

  a = 0;

  b = 6;



type

     one_array = array of real;

     Parametrs = record

     nSimp, nTrap :integer;

     end;



function  y(x:real):real;

function  first(x1,x2:real):real;

procedure setN(nG, nT, nS :integer);

procedure createGrid;

procedure createGridOfInt;

function  integralSimpsona(aSim,bSim:real):real;

function  integralTrapeciay(aTrap,bTrap:real):real;

function  FullIntegral(aFull,bFull:real):real;



var

  Xgrid:one_array;     GridOfInt:one_array;      nGrid, nSim, nTrap :integer;



implementation

//--------------------------------------------------------------

function y(x:real):real;

begin

  Result:=2*arctan(x)-11*Power(3,(-x))+5;

end;

//--------------------------------------------------------------

function  first(x1,x2:real):real;

begin

  result:=5*x2+2*x2*ArcTan(x2)+11*Power(3,(-x2))/Ln(3)-Ln(1+x2*x2)-(5*x1+2*x1*ArcTan(x1)+11*Power(3,(-x1))/Ln(3)-Ln(1+x1*x1));

end;

//--------------------------------------------------------------

procedure setN(nG, nT, nS :integer);

begin

  nGrid:=nG;  nSim:=nS;  nTrap:=nT;

end;

//--------------------------------------------------------------

procedure createGrid;

var i:integer;    h:real;

begin

  h:=(b-a)/nGrid;

  SetLength(Xgrid,nGrid+1);

  for i:=0 to nGrid do

    Xgrid[i]:=a+h*i;

end;

//--------------------------------------------------------------

procedure createGridOfInt;

var i, n :integer;

begin

  n:=High(Xgrid);

  SetLength(GridOfInt,n+1);

  for i:=1 to n do

    GridOfInt[i]:=integralSimpsona(a,Xgrid[i]);

end;

//--------------------------------------------------------------

function  integralSimpsona(aSim,bSim:real):real;

var X:one_array;      i:integer;      sum1, sum2, h:real;

begin

  if(aSim<>bSim)then

    begin

      h:=(bSim-aSim)/nSim;

      SetLength(X,nSim+1);

      for i:=0 to nSim do

        X[i]:=aSim+h*i; 

      sum1:=0;

      sum2:=0;

      for i:=1 to (nSim div 2) do

        sum1:=sum1+y(X[2*i-1]);  

      for i:=1 to ((nSim div 2)-1) do

        sum2:=sum2+y(X[2*i]);

      Result:=(bSim-aSim)*(y(X[0])+y(X[nSim])+4*sum1+2*sum2)/(3*nSim);

    end

  else

    Result:=0;

end;

//--------------------------------------------------------------

function  integralTrapeciay(aTrap,bTrap:real):real;

var        i:integer;      sum, h :real;       X:one_array;

begin

  h:=(bTrap-aTrap)/nTrap;

  SetLength(X,nTrap+1);

  for i:=0 to nTrap do

    X[i]:=aTrap+h*i;

  sum:=(y(X[0])+y(X[nTrap]))/2;

  for i:=1 to (nTrap-1) do

      sum:=sum+y(X[i]);

  Result:=sum*h;

end;

//--------------------------------------------------------------

function  FullIntegral(aFull,bFull:real):real;

var        z1, z2, z3, raznost :real;          i, ai, bi :integer;

begin

  if(aFull<>bFull)then

    begin

      for i:=0 to (High(Xgrid)-1) do

        begin

          if((Xgrid[i]<=aFull)and(aFull<=Xgrid[i+1]))then ai:=i;

          if((Xgrid[i]<=bFull)and(bFull<=Xgrid[i+1]))then bi:=i;

        end;

      raznost:=GridOfInt[ai]+integralTrapeciay(Xgrid[ai],aFull);

      Result:=GridOfInt[bi]+integralTrapeciay(Xgrid[bi],bFull)-raznost;

    end

  else

    Result:=0;

end;

//--------------------------------------------------------------

end.





  unit UnitMAIN;

interface

uses

  Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms, Dialogs,  jpeg, ExtCtrls, StdCtrls,  funct;



type

  TForm1 = class(TForm)

    Image1: TImage;    Image2: TImage;    Button1: TButton;    Image3: TImage;    Edit1: TEdit;    Edit2: TEdit;

    Label1: TLabel;    Label2: TLabel;    Label3: TLabel;    Memo1: TMemo;    Panel1: TPanel;    Button2: TButton;

    Button3: TButton;    Label4: TLabel;    Label5: TLabel;    Label6: TLabel;    Label7: TLabel;    Edit3: TEdit;

    Edit4: TEdit;    Edit5: TEdit;    Label8: TLabel;



    procedure Button1Click(Sender: TObject);

    procedure Button2Click(Sender: TObject);

    procedure Button3Click(Sender: TObject);

    procedure FormCreate(Sender: TObject);

  end;



var

  Form1: TForm1;



implementation



{$R *.DFM}



procedure TForm1.Button1Click(Sender: TObject);

var          i:integer;            intF, intS, intT, intR, L,R:real;          file1:TextFile;

begin

  L:=StrToFloat(Edit1.Text);            R:=StrToFloat(Edit2.Text);

  if(((0<=L)and(L<=6))and(((0<=R)and(R<=6))))

    then

      begin

        createGrid;

        createGridOfInt;

        intF:=FullIntegral(L,R);

        Label1.Caption:=FloatToStr(intF);

        intS:=integralSimpsona(L,R);

        intT:=integralTrapeciay(L,R);

        intR:=first(L,R);

        AssignFile(file1,'Result.txt');

        Rewrite(file1);

        writeln(file1,'Дйствительное значение интеграла            ',intR:17:15);

        writeln(file1,'Значение интеграла по формуле (1)             ',intF:17:15,'      Расность по модулю ',abs(intF-intR):17:15);

        writeln(file1,'Значение интеграла по формуле трапеций  ',intT:17:15,'     Расность по модулю ',abs(intT-intR):17:15);

        writeln(file1,'Значение интеграла по формуле Симпсона ',intS:17:15,'     Расность по модулю ',abs(intS-intR):17:15);

        CloseFile(file1);

        Memo1.Lines.LoadFromFile('Result.txt');

      end

  else ShowMessage('граници должна быть в пределах от 0 до 6');



end;



procedure TForm1.Button2Click(Sender: TObject);

begin

  Panel1.Visible:=true;    Button2.Visible:=false;

end;



procedure TForm1.Button3Click(Sender: TObject);

begin

   Panel1.Visible:=false;     Button2.Visible:=true;

   setN(StrToInt(Edit5.Text),StrToInt(Edit3.Text),StrToInt(Edit4.Text));

end;



procedure TForm1.FormCreate(Sender: TObject);

begin

  setN(StrToInt(Edit5.Text),StrToInt(Edit3.Text),StrToInt(Edit4.Text));

end;

end.