Ohraničená os x, graf integrovateľnej funkcie a úsečky x=a\,\! a x=b\,\!, kde a\,\! a b\,\!- integračné limity (pozri obrázok).

Potreba aplikovať numerickú integráciu môže byť najčastejšie spôsobená absenciou reprezentácie v a teda nemožnosťou analyticky vypočítať hodnotu určitého integrálu nad . Je tiež možné, že tvar primitívnej derivácie je taký zložitý, že je rýchlejšie vypočítať hodnotu integrálu numericky.

Jednorozmerné puzdro

Hlavnou myšlienkou väčšiny metód numerickej integrácie je nahradiť integrand jednoduchším, ktorého integrál sa dá jednoducho analyticky vypočítať. V tomto prípade na odhad hodnoty integrálu použite vzorce formulára

I \približne \sum_(i=1)^(n) w_i\, f(x_i),

kde n\,\! je počet bodov, pri ktorých sa vypočítava hodnota integrandu. bodov x_i\,\! sa nazývajú uzly metódy, čísla w_i\,\!- hmotnosti uzlov. Keď je integrand nahradený polynómom nula, prvého a druhého stupňa, získajú sa metódy , a (Simpson). Vzorce na odhad hodnoty integrálu sa často nazývajú kvadratúrne vzorce.

Metóda obdĺžnika

Metóda obdĺžnika sa získa nahradením integrandu konštantou. Ako konštantu môžete vziať hodnotu funkcie v akomkoľvek bode segmentu \vľavo\,\!. Najčastejšie používané funkčné hodnoty sú v strede segmentu a na jeho koncoch. Zodpovedajúce modifikácie sa nazývajú metódy stredné obdĺžniky, ľavé obdĺžniky a pravé obdĺžniky. Vzorec na približný výpočet hodnoty určitého integrálu metódou obdĺžnikov je

I \približne f(x) (b-a),

kde x=\frac(\left(a+b\right))(2), a\,\! alebo b\,\!, resp.

Lichobežníková metóda

Ak nakreslíme priamku cez konce integračného segmentu, dostaneme lichobežníková metóda. Z geometrických úvah je ľahké ho získať

I \cca \frac(f(a)+f(b))(2) (b-a).

Parabolová metóda

Pomocou troch bodov integračného segmentu môžeme nahradiť integrand parabolou. Zvyčajne sa ako také body používajú konce segmentu a jeho stred. V tomto prípade je vzorec veľmi jednoduchý

I \približne \frac(b-a)(6)\vľavo(f(a)+4f\vľavo(\frac(a+b)(2)\vpravo)+f(b)\vpravo).

Zvýšenie presnosti

Aproximácia funkcie jedným polynómom v celom intervale integrácie spravidla vedie k veľkej chybe v odhade hodnoty integrálu.

Na zníženie chyby je integračný segment rozdelený na časti a na vyhodnotenie integrálu na každej z nich sa používa numerická metóda.

Keďže počet oddielov má tendenciu k nekonečnu, odhad integrálu má tendenciu k jeho skutočnej hodnote pre akúkoľvek numerickú metódu.

Vyššie uvedené metódy umožňujú jednoduchý postup rozdelenia kroku na polovicu, pričom v každom kroku je potrebné počítať hodnoty funkcie len na novo pridaných uzloch. Na odhadnutie chyby výpočtu sa používa.

Gaussova metóda

Vyššie popísané metódy využívajú pevné body segmentov (konce a stredy) a sú nízke (1, 1 a 3). Ak si môžeme vybrať body, v ktorých počítame funkčné hodnoty f(x)\,\!, potom je možné s rovnakým počtom výpočtov integrandu získať metódy vyššieho rádu presnosti. Takže pre dva (ako v lichobežníkovej metóde) výpočty hodnôt integrandu môžete získať metódu už nie 1., ale 3. rádu presnosti:

I \approx \frac(b-a)(2)\left(f\left(\frac(a+b)(2) - \frac(b-a)(2\sqrt(3)) \right)+f\left( \frac(a+b)(2) + \frac(b-a)(2\sqrt(3)) \right) \right).

Vo všeobecnosti pomocou n\,\! bodov, môžete získať metódu s poradím presnosti 2n-1\,\!. Hodnoty uzlov Gaussovej metódy podľa n\,\! body sú koreňmi Legendrovho polynómu stupňa n\,\!.

Hodnoty uzlov Gaussovej metódy a ich váhy sú uvedené v referenčných knihách špeciálnych funkcií. Najznámejšia je Gaussova päťbodová metóda.

Gauss-Kronrodova metóda

Nevýhodou Gaussovej metódy je, že nemá jednoduchý (z výpočtového hľadiska) spôsob odhadu chyby získanej hodnoty integrálu. Použitie Rungeovho pravidla si vyžaduje výpočet integrandu pri približne rovnakom počte bodov, bez toho, aby to prakticky zvýšilo presnosť, na rozdiel od jednoduchých metód, kde sa presnosť niekoľkokrát zvyšuje s každým novým oddielom. Kronrod navrhol nasledujúcu metódu odhadu hodnoty integrálu

I \približne \sum_(i=1)^(n) a_i\, f(x_i) + \sum_(i=1)^(n+1) b_i\, f(y_i),

kde x_i\,\!- Gaussova metóda uzly podľa n\,\! bodov a 3n+2\,\! parametre a_i\,\!, b_i\,\!, y_i\,\! sú zvolené tak, aby sa poradie presnosti metódy rovnalo 3n+1\,\!.

Potom na odhad chyby možno použiť empirický vzorec

\Delta = \left(200 |I - I_G|\right)^(1,5),

kde I_G\,\!- hodnota integrálu, odhadnutá Gaussovou metódou podľa n\,\! bodov. Knižnice [

programovanie numerických integračných vzorcov

Úvod

1. Metódy numerickej integrácie

2. Kvadratúrne vzorce

3. Automatický výber integračného kroku

Záver

Bibliografický zoznam

Úvod

Účelom abstraktu je štúdium a porovnávacia analýza metód numerickej integrácie funkcií; implementácia týchto metód vo forme strojových programov vo vysokoúrovňovom jazyku a praktické riešenie úloh numerickej integrácie na počítači.

Pri riešení inžinierskych problémov je často potrebné vypočítať hodnoty určitého integrálu formulára

. (1)

Ak je funkcia spojitá na intervale [ a , b] a jeho primitívnu vlastnosť možno určiť pomocou známej funkcie, potom sa výpočet takéhoto integrálu vykoná podľa Newton-Leibnizovho vzorca:

.

V inžinierskych problémoch je zriedka možné získať hodnotu integrálu v analytickej forme. Okrem toho funkcia f (X) môže byť dané napríklad tabuľkou experimentálnych údajov. Preto sa v praxi na výpočet určitého integrálu používajú špeciálne metódy, ktoré sú založené na interpolačnom aparáte.

Myšlienka týchto metód je nasledovná. Namiesto výpočtu integrálu pomocou vzorca (1) sa najprv vypočítajú hodnoty funkcie f (x i) = y i v niektorých uzloch x i Î[ a , b]. Potom sa zvolí interpolačný polynóm P (X) prechod cez získané body ( x i , y i), ktorý sa používa pri výpočte približnej hodnoty integrálu (1):

.

Pri implementácii tohto prístupu majú numerické integračné vzorce nasledujúcu všeobecnú formu:

, (2) - uzly interpolácie, Ai sú nejaké koeficienty, R– zvyškový člen charakterizujúci chybu vzorca. Všimnite si, že vzorce vo forme (2) sa nazývajú kvadratúrne vzorce.

Geometrickým významom numerickej integrácie je výpočet plochy krivočiareho lichobežníka ohraničeného grafom funkcie f (X), os x a dve priame čiary x = a a x = b. Približný výpočet plochy vedie k odmietnutiu zvyškového členu v kvadratúrnych vzorcoch R charakterizujúca chybu metódy, ktorá je navyše superponovaná výpočtovou chybou.

1. Metódy numerickej integrácie

V aplikovanom výskume je často potrebné vypočítať hodnotu určitého integrálu

Ako je známe z kurzu matematiky, analytický výpočet integrálu nie je možné vykonať vo všetkých prípadoch. A aj v prípade, že je možné nájsť analytickú formu tohto integrálu, postup výpočtu dáva približný výsledok, takže vzniká problém s približnou hodnotou tohto integrálu.

Podstata približného výpočtu spočíva v dvoch operáciách: 1. vo výbere konečného čísla namiesto n; 2. v bodovom výbere

v príslušnom segmente.

V závislosti od výberu

dostávame rôzne vzorce na výpočet integrálu: Vzorce pre ľavý a pravý obdĺžnik (5), (6) (5) (6)

Trapézový vzorec:

Simpsonov vzorec

b, a - konce uvažovaného segmentu.

Na porovnanie výsledkov výpočtu podľa vyššie uvedených numerických integračných vzorcov vypočítame nasledujúci integrál 3 spôsobmi, pričom segment rozdelíme na 6 rovnakých segmentov: h=

Podľa vzorca ľavých obdĺžnikov:

Podľa lichobežníkového vzorca:

Podľa Simpsonovho vzorca:

A výsledok získaný analyticky sa rovná

=1

Preto môžeme konštatovať, že numerická metóda integrácie podľa Simpsonovho vzorca je presnejšia, ale používa sa vo všeobecnom prípade pri delení segmentu, ktorý sa delí na párny počet intervalov.

2. Kvadratúrne vzorce

Obdĺžnikové vzorce sú najjednoduchšie kvadratúrne vzorce. Rozdeľme integračný segment [ a, b] na P rovnaká dĺžka dielov

. Všimnite si, že hodnota h sa nazýva integračný krok. Na rozdelených bodoch X 0 = a ,X 1 = a + h , ..., x n = b všimnite si súradnice r 0 ,r 1 ,…,y n nepoctivý f (X), t.j. vypočítať i = f (x i)x i = a+ ih = x i -1 + h (i =). Na každom segmente dĺžky h zostavte obdĺžnik so stranami h a y i, kde i =, t.j. hodnotami súradníc vypočítanými na ľavých koncoch segmentov. Potom plochu krivočiareho lichobežníka, ktorá určuje hodnotu integrálu (1), možno približne znázorniť ako súčet plôch obdĺžnikov (obr. 1). Odtiaľ dostaneme vzorec obdĺžnikov:
. (3)

Ak pri výpočte integrálneho súčtu berieme hodnoty funkcie f (X) nie na ľavom, ale na pravom konci segmentov dĺžky h, ktorý je znázornený na obr. 1 s bodkovanou čiarou, potom dostaneme druhú verziu vzorca obdĺžnika:

. (4)

Tretí variant vzorca obdĺžnikov je možné získať pomocou hodnôt funkcie f (X) vypočítané v strede každého segmentu dĺžky h(obr. 2):

. (5)

Vzorce (3), (4) a (4) sa nazývajú vzorce ľavého, pravého a stredného obdĺžnika.

Simpsonov vzorec. Integračný interval delíme na 2 n rovnaká dĺžka dielov

. Na každom segmente [ x i , x i+2] integrand f (X) je nahradená parabolou prechádzajúcou bodmi ( x i , y i), (x i +1 , y i +1), (x i +2 , y i+2). Potom sa približná hodnota integrálu určí podľa Simpsonovho vzorca: . (7)

Pri výpočte na počítači je vhodnejší nasledujúci vzorec:

Simpsonova metóda je jednou z najznámejších a najpoužívanejších metód numerickej integrácie, udáva presné hodnoty integrálu pri integrácii polynómov až do tretieho rádu vrátane.

Newtonov vzorec. Približná hodnota integrálu podľa Newtonovho vzorca sa vypočíta takto:

kde počet segmentov oddielu je násobkom troch, t.j. je 3 n. Pri vývoji počítačových programov je vhodnejšie použiť ekvivalentný vzorec:

Newtonova metóda dáva presné hodnoty integrálu pri integrácii polynómov až do štvrtého rádu vrátane.

3. Automatický výber integračného kroku

Výsledkom výpočtu podľa vzorcov (3) - (8) sa získa približná hodnota integrálu, ktorá sa môže od presnej líšiť o určitú hodnotu, ktorá sa nazýva integračná chyba. Chyba je určená vzorcom zvyšku R rozdielne pre každú z metód integrácie. Ak je potrebné vypočítať hodnotu integrálu s chybou nepresahujúcou e, potom je potrebné zvoliť takýto integračný krok h na uspokojenie nerovnosti R (h) £ e. V praxi sa používa automatický výber hodnoty h, ktorý zabezpečuje dosiahnutie určenej chyby. Najprv vypočítajte hodnotu integrálu ja (n), delením integračného intervalu na P sekcií, potom sa počet sekcií zdvojnásobí a vypočíta sa integrál ja (2n). Proces výpočtu pokračuje, kým sa podmienka nestane pravdivou.

Numerická integrácia

Hlavné otázky diskutované na prednáške:

2. Newton-Cotesove kvadratúrne vzorce

3. Vzorce obdĺžnikov

4. Lichobežníkový vzorec

5. Simpsonov vzorec

6. Kvadratúrne vzorce Gauss

7. Metóda Monte Carlo

1. Vyjadrenie k problému numerickej integrácie

Je potrebné vypočítať určitý integrál tvaru , pričom funkcia môže byť zadaná ako vo forme vzorca, tak aj vo forme tabuľky.

Newton-Cotesove kvadratúrne vzorce

,
kde - Cotesove koeficienty.
Tieto vzorce poskytujú rôzne reprezentácie pre rôzny počet n segmentov oddielu na rovnakom integračnom segmente.

Obdĺžnikové vzorce

Nech je potrebné vypočítať integrál.
Ak je integračný segment dostatočne veľký, musíte ho rozdeliť na menšie segmenty rovnakej dĺžky, kde n je počet segmentov, a nahradením krivočiareho lichobežníka obdĺžnikom na každom zo segmentov, vypočítajte plochy týchto obdĺžnikov. Potom je potrebné výsledné plochy sčítať a toto množstvo bude brané ako približná hodnota požadovaného integrálu.
Čo sa týka konštrukcie obdĺžnikov, môžu byť postavené rôznymi spôsobmi: môžete nakresliť kolmicu na priesečník s krivkou f (x) z pravého konca každého segmentu (obr. 1), môžete - z ľavého konca (obr. 2)

Ryža. jeden

Ryža. 2

V závislosti od toho sú vzorce na výpočet trochu odlišné a nazývajú sa vzorce obdĺžnikov s pravou alebo ľavou súradnicou:

(vzorec „pravých“ obdĺžnikov)

(vzorec „ľavých“ obdĺžnikov)
Existuje aj vzorec pre „stredné“ obdĺžniky: , pre ktoré sa konštrukcia obdĺžnikov vykonáva cez stredy každého zo segmentov oddielu:

· Lichobežníkový vzorec

· Simpsonov vzorec

Nahradením časti krivky y = na každom segmente priečky f(x) na parabolickú krivku, vypočítajúc plochy výsledných čísel a ich súčet, dostaneme Simpsonov vzorec:

·

· Kvadratúrne vzorce Gauss

Tradične sa pri získavaní kvadratúrnych Gaussových vzorcov v pôvodnom integráli vykoná zmena premennej, pričom sa integrál prevedie cez segment na integrál cez segment [-1; jeden]:

.
Potom .
Použijeme lineárnu interpoláciu integrandu.
Ak namiesto segmentu [-1; 1] ak chcete vziať pohyblivé uzly t1, t2 ako interpolačné uzly, musíte tieto hodnoty zvoliť tak, aby plocha lichobežníka ohraničená zhora priamkou prechádzajúcou bodmi A1 (t1, φ(t1) ) a A2 (t2, φ(t2)) sa rovnal integrálu akéhokoľvek polynómu nejakého najvyššieho stupňa.
Za predpokladu, že ide o polynóm tretieho stupňa, vypočítame t1, t2, ktoré sa rovnajú a líšia sa iba očíslovaním hodnôt.
Ďalej, rozdelením integračného segmentu na n častí, aplikovaním vyššie opísanej myšlienky na každú z nich, môžeme získať Gaussov vzorec:

Numerická integrácia

Numerická integrácia(historický názov: (numerický) kvadratúra ) - výpočet hodnoty určitého integrálu (spravidla približný). Numerická integrácia je chápaná ako súbor numerických metód na zistenie hodnoty určitého integrálu.

Numerická integrácia sa používa, keď:

V týchto dvoch prípadoch nie je možné vypočítať integrál pomocou Newtonovho-Leibnizovho vzorca. Je tiež možné, že tvar primitívnej derivácie je taký zložitý, že je rýchlejšie vypočítať hodnotu integrálu numericky.

Jednorozmerné puzdro

Hlavnou myšlienkou väčšiny metód numerickej integrácie je nahradiť integrand jednoduchším, ktorého integrál sa dá jednoducho analyticky vypočítať. V tomto prípade na odhad hodnoty integrálu použite vzorce formulára

kde je počet bodov, v ktorých sa vypočítava hodnota integrandu. Body sa nazývajú uzly metódy, čísla sú váhy uzlov. Keď je integrand nahradený polynómom nula, prvého a druhého stupňa, získajú sa metódy obdĺžnikov, lichobežníkov a parabol (Simpson). Vzorce na odhad hodnoty integrálu sa často nazývajú kvadratúrne vzorce.

Špeciálnym prípadom je metóda na zostavovanie integrálnych kvadratúrnych vzorcov pre rovnomerné siete, tzv Cotesove vzorce. Metóda je pomenovaná po Rogerovi Coatesovi. Hlavnou myšlienkou metódy je nahradiť integrand nejakým druhom interpolačného polynómu. Po prebratí integrálu môžeme písať

kde sa volajú čísla Cotesove koeficienty a sú vypočítané ako integrály zodpovedajúcich polynómov v pôvodnom interpolačnom polynóme pre integrand s hodnotou funkcie v uzle ( je krok mriežky; je počet uzlov mriežky a index uzla je ). Termín je chyba metódy, ktorú možno nájsť rôznymi spôsobmi. V prípade nepárnych možno chybu nájsť integráciou chyby interpolačného polynómu integrandu.

Špeciálne prípady Cotesových vzorcov sú: obdĺžnikové vzorce (n=0), lichobežníkové vzorce (n=1), Simpsonov vzorec (n=2), Newtonov vzorec (n=3) atď.

Metóda obdĺžnika

Nech je potrebné určiť hodnotu integrálu funkcie na intervale . Tento segment je rozdelený bodmi na rovnaké segmenty dĺžky Označme hodnotou funkcie v bodoch Ďalej tvoríme súčty Každý zo súčtov je integrálnym súčtom pre on a teda približne vyjadruje integrál

Ak je daná funkcia kladná a rastúca, potom tento vzorec vyjadruje plochu stupňovitého útvaru zloženého z „prichádzajúcich“ obdĺžnikov, nazývaných aj vzorec ľavých obdĺžnikov, a vzorca

vyjadruje plochu stupňovitého útvaru pozostávajúceho z „vychádzajúcich“ obdĺžnikov, nazývaných aj vzorec pravých obdĺžnikov. Čím kratšia je dĺžka segmentov, na ktoré je segment rozdelený, tým presnejšia je hodnota vypočítaná týmto vzorcom požadovaného integrálu.

Očividne sa oplatí počítať s väčšou presnosťou, ak ako orientačný bod pre zistenie výšky berieme bod v strede medzery. V dôsledku toho dostaneme vzorec pre stredné obdĺžniky:

Vzhľadom na a priori väčšiu presnosť posledného vzorca s rovnakým objemom a povahou výpočtov sa nazýva vzorec obdĺžnikov

Lichobežníková metóda

Ak je funkcia na každom z čiastkových segmentov aproximovaná priamkou prechádzajúcou konečnými hodnotami, získame lichobežníkovú metódu.

Plocha lichobežníka na každom segmente:

Chyba aproximácie na každom segmente:

kde

Úplný vzorec pre lichobežníky v prípade rozdelenia celého integračného intervalu na segmenty rovnakej dĺžky:

kde

Chyba lichobežníkového vzorca:

kde

Parabolová metóda (Simpsonova metóda)

Pomocou troch bodov integračného segmentu môžeme nahradiť integrand parabolou. Zvyčajne sa ako také body používajú konce segmentu a jeho stred. V tomto prípade je vzorec veľmi jednoduchý

.

Ak rozdelíme integračný interval na rovnaké časti, potom máme

Zvýšenie presnosti

Aproximácia funkcie jedným polynómom v celom intervale integrácie spravidla vedie k veľkej chybe v odhade hodnoty integrálu.

Na zníženie chyby je integračný segment rozdelený na časti a na vyhodnotenie integrálu na každej z nich sa používa numerická metóda.

Keďže počet oddielov má tendenciu k nekonečnu, odhad integrálu smeruje k jeho skutočnej hodnote pre analytické funkcie pre akúkoľvek numerickú metódu.

Vyššie uvedené metódy umožňujú jednoduchý postup rozdelenia kroku na polovicu, pričom pri každom kroku je potrebné počítať funkčné hodnoty len na novo pridaných uzloch. Na odhadnutie chyby výpočtu sa používa pravidlo Runge.

Gaussova metóda

Vyššie popísané metódy využívajú pevné segmentové body (konce a stred) a majú nízku presnosť (1 - metóda pravého a ľavého obdĺžnika, 2 - metóda stredného obdĺžnika a lichobežníka, 3 - metóda paraboly (Simpsonova)). Ak si môžeme vybrať body, v ktorých počítame hodnoty funkcie, potom môžeme získať metódy vyššieho rádu presnosti s rovnakým počtom výpočtov integrandu. Takže pre dva (ako v lichobežníkovej metóde) výpočty hodnôt integrandu môžete získať metódu už nie 2., ale 3. rádu presnosti:

.

Vo všeobecnosti pomocou bodov môžete získať metódu s poradím presnosti. Hodnoty uzlov Gaussovej metódy podľa bodov sú koreňmi Legendreho polynómu stupňa.

Hodnoty uzlov Gaussovej metódy a ich váhy sú uvedené v referenčných knihách špeciálnych funkcií. Najznámejšia je Gaussova päťbodová metóda.

Gauss-Kronrodova metóda

Nevýhodou Gaussovej metódy je, že nemá jednoduchý (z výpočtového hľadiska) spôsob odhadu chyby získanej hodnoty integrálu. Použitie Rungeovho pravidla si vyžaduje výpočet integrandu v približne rovnakom počte bodov, bez toho, aby to prakticky zvýšilo presnosť, na rozdiel od jednoduchých metód, kde sa presnosť niekoľkokrát zvyšuje s každým novým oddielom. Kronrod navrhol nasledujúcu metódu odhadu hodnoty integrálu

,

kde sú uzly Gaussovej metódy podľa bodov a parametre , , sú zvolené tak, aby sa poradie presnosti metódy rovnalo .

Potom na odhad chyby môžete použiť empirický vzorec:

,

kde je približná hodnota integrálu získaná Gaussovou metódou nad bodmi. Knižnice gsl a SLATEC na výpočet určitých integrálov obsahujú rutiny využívajúce Gauss-Kronrodovu metódu pre 15, 21, 31, 41, 51 a 61 bodov. Knižnica používa Gaussovu-Kronrodovu metódu pre 15 bodov.

Čebyševova metóda

Integrácia pod nekonečnými hranicami

Ak chcete integrovať cez nekonečné limity, musíte zaviesť nerovnomernú mriežku, ktorej kroky sa zväčšujú, keď idete do nekonečna, alebo môžete urobiť takú zmenu premenných v integráli, po ktorej budú limity konečné. Podobne možno postupovať, ak je funkcia na koncoch integračného intervalu singulárna

Metódy Monte Carlo

Obrázok 3 Numerická integrácia funkcie metódou Monte Carlo

Na určenie oblasti pod funkčným grafom môžete použiť nasledujúci stochastický algoritmus:

Pre malý počet rozmerov integrovateľnej funkcie je výkon integrácie Monte Carlo oveľa nižší ako výkon deterministických metód. Avšak v niektorých prípadoch, keď je funkcia špecifikovaná implicitne, ale je potrebné určiť oblasť špecifikovanú vo forme komplexných nerovností, môže byť vhodnejšia stochastická metóda.

Metódy Runge-Kutta

spline metóda

Viacrozmerný prípad

V malých rozmeroch je možné použiť aj kvadratúrne vzorce založené na interpolačných polynómoch. Vo vyšších dimenziách sa však tieto metódy stávajú neprijateľnými v dôsledku rýchleho nárastu počtu bodov siete a/alebo komplexnej hranice regiónu. V tomto prípade sa používa metóda Monte Carlo. Náhodné body sú generované v našej oblasti a funkčné hodnoty v nich sú spriemerované. Môžete použiť aj zmiešaný prístup - rozdeliť oblasť na niekoľko častí, z ktorých v každej (alebo len v tých, kde sa integrál nedá vypočítať kvôli zložitej hranici) aplikujte metódu Monte Carlo.

Literatúra

1. Kahaner D., Moler K., Nash S. Numerické metódy a softvér (preložené z angličtiny). M.: Mir, 2001, 575 s.