Strona główna » Algorytmy » Artykuły » Całkowanie Numeryczne
 

Całkowanie Numeryczne

https://mattomatti.com/a54

Wstęp

Całkę można intuicyjnie wytłumaczyć jako obszar ograniczony przez osie odcięte oraz wykres funkcji f(x) na pewnym przedziale [a, b]. Istnieje wiele rodzajów całek. Jednak nie zawsze można wyznaczyć całkę dla funkcji f(x). Szukane pole można wyznaczyć przy pomocy tzw. oznaczonej całki Riemanna. W przypadku obliczeń komputerowych nie zawsze istnieje potrzeba wyznaczania wzoru całki. W aplikacjach możemy otrzymać przybliżone wyniki dzięki Całkowaniu numerycznemu.

Całkowanie Numeryczne

Rozpatrujemy całkę nieujemną funkcji f(x) na przedziale [a, b]. Całkę traktujemy tutaj jako pole powierzchni ograniczone przez proste x = a, x = b, oś odciętą oraz od góry poprzez wykres funkcji. Podaną interpretację można przedstawić graficznie jako:

W Całkowaniu numerycznym, aby obliczyć pole powierzchni zastępujemy funkcję podcałkową prostszą funkcją, którą potrafimy policzyć. Jednak takie uproszczenie powoduje, że wyniki są tylko przybliżone. Istnieją jednak różne warianty pozwalające obliczyć pole tą metodą. Dzięki temu przybliżenia są na tyle małe, że można je zaniedbać.

Wstępna implementacja

Ze względu na fakt, że całkowanie wymaga podania funkcji f(x) do programu dołączę bibliotekę wielomiany.h. Jest to biblioteka napisana przeze mnie, a proces jej powstawania można prześledzić tutaj. Dzięki niej program będzie umiał wczytać wielomian dowolnego stopnia.

Oczywiście z tej biblioteki nie trzeba korzystać. Jednak pozwoli ona znacząco usprawnić pisanie programu poprzez gotowe funkcje. Istnieje możliwość napisania oddzielnej funkcji np. f_wartosc(double x), która będzie obliczała wartość funkcji dla podanego argumentu x.

Warianty Całkowania Numerycznego

Metoda prostokątów

Przedział całkowania dzielimy na n równych odcinków. Następnie obliczamy wartości funkcji f(x) w punktach a i b oraz w punktach, które są końcami dwóch odcinków. Wtedy całkę obliczamy następująco:

Funkcja, która będzie wyliczać pole musi dla każdego odcinka wyliczyć średnią arytmetyczną dla wartości funkcji na początku i na końcu każdego odcinka. W ten sposób uzyskujemy dla i-tego odcinka wysokość h. Potem pozostaje tylko pomnożyć przez podstawę - w tym przypadku jest to długość odcinka czyli .

Metoda trapezów

Przedział całkowania dzielimy na n równych odcinków. Następnie obliczamy wartości funkcji f(x) w punktach a i b oraz w punktach, które są końcami dwóch odcinków. Wtedy całkę obliczamy następująco:

W tym przypadku funkcja, aby wyliczyć pole i-tego trapezu musi przyjąć za jedną podstawę wartość f(xi), a za drugą f(xi + 1). Wysokością jest tu długość i-tego odcinka (). Oczywiście samo pole obliczymy wtedy z podstawowego wzoru trapezu co w tym przypadku przybierze formę:

Implementacja

Chociaż oba warianty wyliczania całki różnią się to w rzeczywistości otrzymujemy ten sam wynik dla podanej dokładności. Jako dowód wystarczy porównać pole i-tego prostokąta oraz trapezu. Obie strony są sobie równe, dlatego napiszemy funkcję calka, która będzie obliczać pole powierzchni zgodnie z algorytmem całkowania numerycznego.

Funkcja calka() będzie przyjmowała cztery argumenty: wskaźnik na wielomian w, a i b oznaczające kolejno początek i koniec zakresu całkowania oraz n - określający długość pojedynczego przedziału. n to w rzeczywistości dokładność. Im argument będzie mniejszy tym dokładniejszy będzie wynik:

  1. double calka(wielomian* w, double a, double b, double n = 0.01){
  2.   double suma = 0;
  3.   for(double x = a; x < b; x += n)
  4.     suma += abs(((wartosc_wielomianu(w, x) + wartosc_wielomianu(w, x + n)) / 2) * n);
  5.   return suma;
  6. }

(2.) Deklarujemy zmienną sumę, która będzie przechowywała dotychczas policzone pole. (3.) Dla każdego odcinka długość n w przedziale [a, b]: (4.) wyliczamy pole i dodajemy do suma. (4.) Zwracamy wartość wyliczonego pola.

Testowanie funkcji

Napisane wcześniej funkcję można przetestować zamieszczoną poniżej funkcją main().

  1. int main () {
  2.   wielomian* w = wczytajWieloman();
  3.   wypiszWieloman(w);
  4.   double a, b, n;
  5.   cin >> a >> b >> n;
  6.   cout << calka(w, a, b, n) << endl;
  7.   system("pause");
  8.   return 0;
  9. }

Przypuśćmy, że mamy funkcję postaci: f(x) = x2 + 1 i chcemy obliczyć całkę na przedziale [0, 1]. Wtedy program wywołamy następująco:

  1. 2
  2. 1 0 1
  3. 0 1
  4. 0.001

(1.) Określamy stopień wielomianu. (2.) Współczynniki kolejnych stopni wielomianu. (3.) Zakres a, b. (4.) Dokładność obliczeń (długość odcinka xi). Jako wynik otrzymamy obliczone pole dla podanej funkcji:

  1. 1.33333

Zadania

Zadanie 1

Napisz funkcję calka(), gdzie argument n nie będzie długością odcinka, a będzie określał na ile odcinków został podzielony dany zakres [a, b]

Dla przykładu podanego przy testowaniu funkcji tym razem będzie trzeba wpisać:

  1. 2
  2. 1 0 1
  3. 0 1
  4. 1000

Kolejne linijki oznaczają: (1.) stopień wielomianu, (2.) współczynniki wielomianu, (3.) zakres całkowania i (4.) na ile odcinków ma zostać podzielony zakres.