Numerička integracija

U numerici, numerička integracija se sastoji od velike porodice algoritama za računanje numeričke vrijednosti određenog integrala, a pojam se također koristi za opis numeričkog rješenja diferencijalne jednačine. Termin numerička kvadratura je manje-više sinonim za numeričku integraciju, pogotovo za jednodimenzionalne integrale. Numerička integracija u više dimenzija opisuje se kao kubatura,^[1] iako se izrazom kvadratura podrazumijeva i višedimenziona integracija.

Osnovni problem u numeričkoj integraciji je izračunati približnu vrijednost određenog integrala:

I=\int _{a}^{b}\!f(x)\,dx

sa zadatim stepenom tačnosti. Ako je $f(x)$ glatka funkcija integrisana preko malog broja dimenzija i ako je domen integracije zatvoren, tada postoji više metoda za aproksimaciju integrala na željenu preciznost.

Historija

Kvadratura je historijski matematički pojam, a znači računanje površine. Kvadraturni problemi su bili glavni zadaci koji su zadavani kao izvor za matematičku analizu. Starogrčki matematičari su prema Pitagorinoj doktrini shvatili računanje površine kao proces konstrukcije geometrijskog kvadrata koji ima jednaku površinu (kvadriranje). To je razlog zašto je proces nazvan kvadratura. Primjeri: kvadratura kruga, hipokratov mjesec i kvadratura parabole. Ova konstrukcija može biti izvedena jedino upotrebom šestara i lenjira.

Za kvadraturu pravougaonika sa stranicama a i b potrebno je konstruisati kvadrat sa stranicom: $x={\sqrt {ab}}$

Za ovu svrhu moguće je koristiti sljedeću činjenicu: ako se nacrta krug sa zbirom stranica a i b kao prečnik, onda je visina BH (sa tačke njihovog dodira do presijecanja sa krugom) jednaka njihovoj geometrijskoj sredini. Jednaka geometrijska konstrukcija rješava problem kvadrature paralelograma i trougla.

Problemi kvadrature za krivolinijske figure su mnogo složeniji. Kvadratura jednog kruga sa šestarom i lenjirom je bila dokazana u 19. vijeku kao nemoguća. Ipak, za neke oblike (npr. hipokratov mjesec) kvadratura se može obaviti. Kvadrature sferične površine i dijela parabole koje je uradio Arhimed postao je najveći uspjeh u antičkoj analizi.

Površina lopte je jednaka četverostrukoj površini velikog kruga ove lopte.
Površina dijela parabole odsječena pravcem je 4/3 površine od trougla koji je upisan u ovaj segment.

Za dokaz rezultata, Arhimed je koristio metodu Eudoksa.

U srednjovjekovnoj Evropi, kvadratura je značila proračun površine koristeći bilo koju metodu. Češće je metoda nedjeljivosti korištena; manje je rigorozna, a jednostavnija i dobra. Pomoću nje, Galileo Galilei i Gilles de Roberval su našli površinu svoda cikloide, Grégoire de Saint-Vincent je istražio površinu ispod hiperbole (Opus Geometricum, 1647), a Alphonse Antonio de Sarasa, de Saint-Vincentov učenik i komentator, uspostavio je relaciju ove površine sa logaritmima.

John Wallis je "algebrizirao" ovaj metod: napisao je u svojim Arithmetica Infinitorum (1656.) serijama ono što se danas zove određenim integralom, te izračunao njihove vrijednosti. Isaac Barrow i James Gregory su napravili dalji napredak: kvadrature za neke algebarske krivulje i spirale. Christiaan Huygens je uspješno izveo kvadrature nekih "materijala revolucije".

Kvadratura hiperbole od Saint-Vincenta i de Sarasae je ustupila novu funkciju, prirodni logaritam, kritičnog značaja.

S otkrićem integralnog računa došla je univerzalna metoda za računanje površine. Kao odgovor, termin kvadratura je postao tradicionalan, a umjesto toga moderna fraza "računanje određenog integrala" je češće upotrebljavana.

Razlozi za numeričku integraciju

Postoji više razloga zašto se koristi numerička integracija. Funkcija koja se integriše f(x) može biti poznata samo na pojedinim mjestima, što se radi uzimanjem uzorka. Neki super računari i ostale računarske aplikacije nekad trebaju numeričku integraciju baš radi ovog razloga.

Formula za funkciju koja se integriše može biti poznata, ali može biti teško ili nemoguće naći antiderivaciju koja je elementarna funkcija. Jedan primjer je funkcija f(x) = exp(−x²), antiderivacija koja se ne može napisati u elementarnoj formi.

Moguće je naći antiderivaciju simbolično, ali je mnogo jednostavnije naći numeričku aproskimaciju nego računati antiderivaciju (anti-izvod). Ovo može biti korišteno ako je antiderivacija data kao neograničeni niz proizvoda ili ako bi proračun zahtijevao specijalne funkcije koje nisu dostupne računarima.

Metode za jednodimenzione integrale

Metode numeričke integracije se općenito mogu opisati kao kombinacije procjena integranda (funkcije) da se dobije aproksimacija integrala. Integrand je definisan kao konačan skup tačaka nazivanih integracione tačke i zbir ovih vrijednosti se koristi za aproksimaciju integrala. Integracione tačke i zbir zavise od posebnih metoda koje se koriste i preciznosti koja se traži aproksimacijom.

Važan dio analize bilo koje numeričke integracione metode je proučavanje ponašanja greške aproksimacije kao funkcija broja procjene integranda. Metoda koja ima malu grešku za mali broj procjena se često smatra superiornom. Smanjenjem broja procjena integranda, smanjuje se i broj korištenih aritmetičkih operacija te se smanjuje ukupna greška proračuna. Također, svaka procjena uzima vrijeme, a integrand može biti svojevoljno komplikovan.

Brute force metoda ("nasilna metoda" koja koristi sve kombinacije) numeričkog integrisanja može biti obavljena ako integrand ima "dobro ponašanje" (tj. ako je funkcija neprekidna i iz ograničene varijacije), sa procjenom integranda sa veoma malim korakom - inkrementom.

Kvadraturna pravila zasnovana na interpolaciji funkcija

Veliki broj kvadraturnih pravila može se naći izvođenjem funkcija koje su jednostavne za integriranje. Najčešće funkcije koje se interpoliraju su polinomi.

Najjednostavnija metoda ovog tipa je dopuštanje da interpolirana funkcija bude konstantna funkcija (polinom nultog stepena) koja prolazi kroz tačku sa koordinatama ((a+b)/2, f((a+b)/2)). Ovo se zove pravilo sredine (srednje tačke) ili kvadratno pravilo.

\int _{a}^{b}f(x)\,dx\approx (b-a)\,f\left({\frac {a+b}{2}}\right)

Interpolacijska funkcija može biti polinom prvog stepena koja prolazi kroz tačke: (a, f(a)) i (b, f(b)). Ovo se zove trapezno pravilo.

\int _{a}^{b}f(x)\,dx\approx (b-a)\,{\frac {f(a)+f(b)}{2}}.

Za bilo koje od ovih pravila može se napraviti preciznija aproksimacija prekidanjem intervala [a, b] na veći broj podintervala, računajući aproksimaciju za svaki podinterval te sabirajući sve u jedan rezultat. Ovo se zove kompozitno pravilo, prošireno pravilo ili iterativno pravilo. Naprimjer, kompozitno trapezno pravilo ima slijedeći oblik:

\int _{a}^{b}f(x)\,dx\approx {\frac {b-a}{n}}\left({f(a) \over 2}+\sum _{k=1}^{n-1}\left(f\left(a+k{\frac {b-a}{n}}\right)\right)+{f(b) \over 2}\right)

gdje podintervali imaju oblik [k h, (k+1) h], sa h = (b−a)/n i k = 0, 1, 2, ..., n−1.

Interpolacija sa polinomima procijenjena tačkama sa jednakim rastojanjima u zatvorenom sektoru [a, b] se radi Newton-Cotesovim formulama (njutn-kotes), od koje su primjeri trapezno i kvadratno pravilo. Simpsonovo pravilo, koje se bazira na polinomima drugog reda, također je Newton-Cotesova formula.

Kvadraturna pravila sa jednako razmaknutim tačkama imaju veoma pogodnu osobinu gniježđenja. Odgovarajuće pravilo sa svakim potpodijeljenim intervalom sadrži sve trenutne tačke, pa ove vrijednosti integranda mogu biti ponovo iskorištene.

Ako se dopusti da intervali između interpolacionih tačaka variraju, onda se nalazi nova grupa kvadraturnih formula, kao što je Gaussova kvadraturna formula. Gaussovo kvadraturno pravilo je još preciznije od Newton-Cotesovog pravila koje zahtjeva jednak broj iteracija ako je funkcija koja se integriše (integrand) glatka krivulja (npr. ako je dovoljno diferencijabilna). Ostale kvadraturne metode sa varijacijom intervala uključuju Clenshaw–Curtis kvadraturne metode (također zvane Fejér kvarature), koje se gnijezde.

Gaussova kvadraturna pravila nemaju osobinu gniježđenja, ali Gauss–Kronrod kvadraturne formule imaju.

Prilagodljivi algoritmi

Ako f(x) nema puno izvoda na svim tačkama, ili ako izvodi postanu ogromni, onda je Gaussova kvadratura često nedovoljna. U tom slučaju, algoritam ispod će bolje uraditi zadatak:

def racunanje_odredjenog_integrala(f, initial_step_size):
    '''
    Ovaj algoritam računa određeni integral (onaj koji ima granice)
    funkcije od 0 do 1, birajući manje korake kod problematičnih
    tačaka.
    '''
    x = 0.0
    h = initial_step_size
    akumulator= 0.0
    while x < 1.0:
        if x + h > 1.0:
            h = 1.0 - x
        quad_this_step =
        if error_too_big_in_quadrature_of_over_range(f, [x,x+h]):
            h = make_h_smaller(h)
        else:
            akumulator+= quadrature_of_f_over_range(f, [x,x+h])
            x += h
            if error_too_small_in_quadrature_of_over_range(f, [x,x+h]):
                h = make_h_larger(h) # Izbjegavanje dangubljenja na malim koracima.
    return akumulator

Neki detalji algoritma zahtijevaju pažljiv pristup. Za više slučajeva, pronalazak greške kvadrature na intervalu za funkciju f(x) nije očit. Jedno poznato rješenje je koristiti dva pravila kvadratura i korištenjem razlike se može procijeniti greška kvadrature. Čest problem je odlučiti šta označava "premalo" ili "preveliko". Lokalni kriterij za "preveliko" je da kvadraturna greška ne smije biti veća od t · h, gdje je t, realni broj, tolerancija koja se želi podesiti za globalnu grešku. Opet, ako je h već malo, onda ga može biti bespotrebno praviti čak manjim iako je kvadraturna greška velika. Globalni kriterij je da zbir grešaka na svim intervalima bude manji od t. Ovaj tip analize greške se često naziva "posteriori" zato što se nakon svake procjene računa i greška.

Metoda ekstrapolacije

Preciznost kvadraturnog pravila kod Newton-Cotesovog tipa je uglavnom funkcija broja na tačkama procjene. Rezultat je često više pouzdan kad se broj tačaka poveća ili, ekvivalentno, kad se razmaci između tačaka smanje. Nameće se pitanje kakav bi rezultat bio kada bi razmaci među tačkama bili jednaki nula. Ovo se može odgovoriti ekstrapolacijom rezultata od 2 do n razmaka koristeći ubrzavajuće metode kao npr. Richardsonove ekstrapolacije. Ekstrapolaciona funkcija može biti polinom ili racionalna funkcija. Ekstrapolacione metode su opisanje detaljnije od strane Stoera i Bulirscha (Sekcija 3.4) i implementirane su u više metoda u QUADPACK biblioteci.

Konzervativna (a priori) procjena greške

Ako se pusti da f ima ograničen prvi izvod na zatvorenom intervalu [a, b], teorema srednje vrijednosti za f, gdje je x < b, daje:

(x-a)f'(v_{x})=f(x)-f(a)\,

za neke v_x iz intervala [a,x] zavisne od x. Ako se integriše po x od a do b na obje strane i uzmu apsolutne vrijednosti, dobije se:

\left|\int _{a}^{b}f(x)\,dx-(b-a)f(a)\right|=\left|\int _{a}^{b}(x-a)f'(v_{x})\,dx\right|

(*)

Može se dalje aproksimirati integral na desnoj strani nejednakosti ubacujući apsolutnu vrijednost u integrand i smjenom slova u f' prema gornjoj granici:

\left|\int _{a}^{b}f(x)\,dx-(b-a)f(a)\right|\leq {(b-a)^{2} \over 2}\sup _{a\leq x\leq b}\left|f'(x)\right|

(**)

Ako se aproksimira integral ∫_a^b f(x) dx sa kvadraturnim pravilom (b − a)f(a), tada greška nije veća od one na desnoj strani kod (**). Ovo se može pretvoriti u analizu greške za Riemannovu sumu (*), dajući gornju granicu za izražavanje greške od određene aproksimacije:

{n^{-1} \over 2}\sup _{0\leq x\leq 1}\left|f'(x)\right|

(Komentar: U ovom primjeru je ovo greška koja je izračunata za primjer $f(x)=x$ .)

Koristeći više izvoda i popravljanja kvadrature može se dobiti ista greška analize koristeći se Taylorovim redovima (koristeći parcijalnu sumu sa ostatkom) za f. Ova analiza greške daje striktnu gornju granicu greške ako postoje izvodi za f.

Ova metoda integracije može biti kombinovana sa aritmetikom intervala da kreira računarski dokaz i verificirane proračune.

Integrali sa beskonačnim intervalima

Postoji nekoliko metoda za aproksimaciju integracije kod neograničenih intervala. Uobičajena tehnika uključuje posebno izvedena kvadraturna pravila, poput Gauss-Hermite kvadratura za integrale na čitavoj liniji realnih brojeva i Gauss-Laguerre kvadrature za integrale na pozitivnim realnim brojevima.^[2] Monte Carlo metode se također mogu koristiti ili zamjenom promjenljivih na konačan interval. To se može predstaviti za neograničene intervale s obje strane:

\int _{-\infty }^{+\infty }f(x)\,dx=\int _{-1}^{+1}f\left({\frac {t}{1-t^{2}}}\right){\frac {1+t^{2}}{(1-t^{2})^{2}}}\,dt,

ili polu-beskonačne intervale:

\int _{a}^{+\infty }f(x)\,dx=\int _{0}^{1}f\left(a+{\frac {t}{1-t}}\right){\frac {dt}{(1-t)^{2}}}

\int _{-\infty }^{a}f(x)\,dx=\int _{0}^{1}f\left(a-{\frac {1-t}{t}}\right){\frac {dt}{t^{2}}}

kao moguće transformacije.

Višedimenzioni integrali

Kvadraturna pravila koja su razmatrana dosad pravljena su s ciljem računanja jednodimenzionih integrala. Za računanje istih u više dimenzija, jedan pristup je izraziti višestruki integral kao ponovljeni jednodimenzioni integral koristeći se Fubinijevom teoremom. Ovaj pristup traži proračune funkcije da rastu eksponencijalno kako se povećava broj dimenzija. Dvije metode su poznate za postizanje ovoga tzv. prokletstva dimenzionalnosti.

Monte Carlo

Metode Monte Carla i kvazi metode istog su jednostavne za koristiti na višedimenzionim integralima i mogu dobivati veću preciznost za isti broj od proračuna funkcija nego metodom ponovljenih integracija koristeći se jednodimenzionim metodama.

Veliki broj korisnih Monte Carlo metoda su tzv. Markov lanac za Monte Carlo algoritme, koji sadrže Metropolis-Hastings algoritam i Gibbsovo uzimanje uzorka.

Sparsova rešetka

Sparsove rešetke je originalno razvio Smolyak za kvadrature višedimenzionih funkcija. Ova metoda je uvijek zasnovana na jednodimenzionom pravilu kvadrature, ali obavlja sofisticiraniju kombinaciju jednovarijantnih rezultata.

Veza sa diferencijalnim jednačinama

Problem traženja integrala:

F=\int _{a}^{b}f(x)\,dx

može biti reduciran na problem početne vrijednost za uobičajenu diferencijalnu jednačinu. Ako je gornji integral obilježen sa F(x), onda funkcija F zadovoljava:

{\frac {dF(x)}{dx}}=f(x),\quad F(a)=0.

Metode razvijene za uobičajene diferencijalne jednačine, kao Runge-Kutta, mogu biti iskorištene na ponovljeni problem i za proračun integrala. Naprimjer, standardna Runge-Kutta metoda iskorištena na diferencijalnu jednačinu kreira Simpsonovo pravilo odozgo.

Diferencijalna jednačina F ' (x) = ƒ(x) ima specifičan oblik: desna strana sadrži samo zavisnu varijablu (ovdje x) i nezavisnu (ovdje F). Ovo znatno pojednostavljuje teoriju i algoritme.

Također pogledajte

Reference

^ Eric W. Weisstein, Numerička integracija na MathWorld-u.
^ Leader, Jeffery J. (2004). Numerical Analysis and Scientific Computation. Addison Wesley. ISBN 0-201-73499-0. CS1 održavanje: nepreporučeni parametar (link)

Philip J. Davis and Philip Rabinowitz, Methods of Numerical Integration
George E. Forsythe, Michael A. Malcolm, and Cleve B. Moler, Computer Methods for Mathematical Computations. Englewood Cliffs, NJ: Prentice-Hall, 1977. (v. Poglavlje 5.)
Press, W.H.; Teukolsky, S.A.; Vetterling, W.T.; Flannery, B.P. (2007), "Chapter 4. Integration of Functions", Numerical Recipes: The Art of Scientific Computing (3rd izd.), New York: Cambridge University Press, ISBN 978-0-521-88068-8, arhivirano s originala, 4. 3. 2016, pristupljeno 20. 6. 2014
Josef Stoer and Roland Bulirsch, Introduction to Numerical Analysis. New York: Springer-Verlag, 1980. (v. Poglavlje 3.)
Boyer, C. B., A History of Mathematics, 2nd ed. rev. by Uta C. Merzbach, New York: Wiley, 1989 ISBN 0-471-09763-2 (1991 pbk ed. ISBN 0-471-54397-7).
Eves, Howard, An Introduction to the History of Mathematics, Saunders, 1990, ISBN 0-03-029558-0,

Vanjski linkovi

Integracija: Pozadina, simulacije, itd. (en) na Holistic Numerical Methods Institute

Besplatan softver za numeričku integraciju

Numerička integracija je jedan od najviše istraživanih problema u numerici. Od velikog broja softverskih rješenja, izdvojeno je nekoliko besplatnih paketa otvorenog koda:

QUADPACK (dio SLATEC-a): opis [1], izvorni kod [2]. QUADPACK je kolekcija algoritama, u Fortranu, za numeričku integraciju baziranu na Gausovim kvadraturama.
interalg: softver od OpenOpt/FuncDesigner frameworksa, baziran na analizi intervala, garantovanoj preciznosti, licenci: BSD (besplatan za svaku svrhu)
GSL: The GNU Scientific Library (GSL) jeste numerička biblioteka napisana u programskom jeziku C koja pruža veliki obim matematičkih rutina, poput Monte Carlo integracije.
Algoritmi numeričke integracije se mogu naći u GAMS razredu H2.
ALGLIB je kolekcija algoritama napisanih u C# / C++ / Delphi / Visual Basic / itd. za numeričku integraciju (uključuje Bulirsch-Stoer i Runge-Kutta integratore).
Cuba Arhivirano 7. 4. 2014. na Wayback Machine je biblioteka besplatnog softvera od nekoliko višedimenzionih integracionih algoritama.
Cubature kod za višedimenzionu integraciju.
Scilab softver otvorenog koda pod CeCILL licencom (kompatibilna sa GPL licencom), pruža moćne mogućnosti uključujući numeričku integraciju.

[1] Eric W. Weisstein, Numerička integracija na MathWorld-u.

[2] Leader, Jeffery J. (2004). Numerical Analysis and Scientific Computation. Addison Wesley. ISBN 0-201-73499-0. CS1 održavanje: nepreporučeni parametar (link)

[1]

[2]