Regula Cavalerilor-Simpson

Prin regula Knights-Simpson sau regula Knights sau regula lui Simpson se înțelege o metodă pentru calculul numeric aproximativ al integralelor definite ale formei:

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

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

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

Ca toate procedurile pentru calculul aproximativ al integralelor definite și pentru alte calcule aproximative pornind de la funcțiile unei variabile reale, această metodă este utilizată pentru funcții $\,f\left(x\right)$ ${\ displaystyle \, f \ left (x \ right)}$ ${\ displaystyle \, f \ left (x \ right)}$ dintre care nu se cunoaște funcția primitivă sau ale cărei primitive sunt cunoscute numai caracteristici din care nu este posibil să se obțină o expresie prin funcții elementare care pot fi utilizate în mod rezonabil pentru calculele solicitate. Aceste metode aproximative sunt folosite și în cazurile în care nu se cunoaște o expresie analitică a funcției care urmează să fie integrată, dar doar unele dintre valorile sale sunt cunoscute (obținute experimental sau obținute din alte surse), sau când este cunoscută doar diagrama sa (trasat cu ajutorul unor instrumente speciale sau obținut din literatură).

Formula de cvadratură sau metoda parabolelor

Regula Knights-Simpson aproximează integrala funcției necesare (în albastru) cu cea a parabolei care o interpola în noduri (în roșu)

Regula Cavalieri-Simpson prevede subdivizarea intervalului de integrare în subintervale și înlocuirea în aceste subintervale a funcției integrandului prin intermediul arcurilor parabolei , adică prin intermediul polinoamelor pătratice .

Prin urmare, să luăm în considerare $\,\int _{a}^{b}f(x)\,dx.$ ${\ displaystyle \, \ int _ {a} ^ {b} f (x) \, dx.}$ ${\ displaystyle \, \ int _ {a} ^ {b} f (x) \, dx.}$ Pentru simplitate presupunem $f(x)\geq 0$ ${\ displaystyle f (x) \ geq 0}$ ${\ displaystyle f (x) \ geq 0}$ pe tot parcursul intervalului de integrare $\left[a,b\right].$ ${\ displaystyle \ left [a, b \ right].}$ ${\ displaystyle \ left [a, b \ right].}$

Să ne descompunem $\left[a,b\right]$ ${\ displaystyle \ left [a, b \ right]}$ ${\ displaystyle \ left [a, b \ right]}$ într-un număr par $n=2m$ ${\ displaystyle n = 2m}$ ${\ displaystyle n = 2m}$ de sub-intervale fiecare de amplitudine $h:={\frac {b-a}{n}}.$ ${\ displaystyle h: = {\ frac {ba} {n}}.}$ ${\ displaystyle h: = {\ frac {b-a} {n}}.}$ Sa spunem:

x_{i}:=a+i\cdot h,\qquad y_{i}:=f\left(x_{i}\right),\qquad {\text{per }}i=0,1,\ldots ,n,

{\ displaystyle x_ {i}: = a + i \ cdot h, \ qquad y_ {i}: = f \ left (x_ {i} \ right), \ qquad {\ text {per}} i = 0.1, \ ldots, n,}

{\ displaystyle x_ {i}: = a + i \ cdot h, \ qquad y_ {i}: = f \ left (x_ {i} \ right), \ qquad {\ text {per}} i = 0.1, \ ldots, n,}

pentru extremele următoarelor subintervale și pentru valorile pe care funcția le asumă în corespondența lor.

De asemenea, considerăm intervalul parțial format din două subintervaluri consecutive având ca extreme $x_{0}$ ${\ displaystyle x_ {0}}$ $x_0$ Și $x_{2}.$ ${\ displaystyle x_ {2}.}$ ${\ displaystyle x_ {2}.}$ În plus, luăm în considerare și următoarele $m-1$ ${\ displaystyle m-1}$ $m-1$ intervale parțiale având ca extreme respectiv $x_{2}$ ${\ displaystyle x_ {2}}$ $x_2$ Și $x_{4},\ldots ,x_{n-2}$ ${\ displaystyle x_ {4}, \ ldots, x_ {n-2}}$ ${\ displaystyle x_ {4}, \ ldots, x_ {n-2}}$ Și $x_{n}.$ ${\ displaystyle x_ {n}.}$ ${\ displaystyle x_ {n}.}$

În fiecare dintre aceste intervale parțiale ne propunem să le substituim $f\left(x\right)$ ${\ displaystyle f \ left (x \ right)}$ ${\ displaystyle f \ left (x \ right)}$ cu o funcție rațională de număr întreg de gradul doi. Începem de la primul interval parțial și alegem un polinom al formei

y\left(x\right):=A\left(x-x_{1}\right)^{2}+B\left(x-x_{1}\right)+C

{\ displaystyle y \ left (x \ right): = A \ left (x-x_ {1} \ right) ^ {2} + B \ left (x-x_ {1} \ right) + C}

{\ displaystyle y \ left (x \ right): = A \ left (x-x_ {1} \ right) ^ {2} + B \ left (x-x_ {1} \ right) + C}

astfel încât integrala sa între $x_{0}$ ${\ displaystyle x_ {0}}$ $x_0$ Și $x_{2}$ ${\ displaystyle x_ {2}}$ $x_2$ diferă de cea a funcției inițiale printr-o cantitate care poate fi neglijabilă.

Expresia polinomului surogat reprezintă o parabolă generică cu axă de simetrie verticală. Pentru a determina valoarea constantelor $LA, B.$ ${\ displaystyle A, B}$ $A, B$ Și $C.$ ${\ displaystyle C}$ $C.$ trecerea parabolei este impusă pentru punctele coordonate:

\left(x_{0},\,y_{0}:=f(x_{0})\right),~(x_{1},\,y_{1}:=f(x_{1})),~(x_{2},\,y_{2}:=f(x_{2})).

{\ displaystyle \ left (x_ {0}, \, y_ {0}: = f (x_ {0}) \ right), ~ (x_ {1}, \, y_ {1}: = f (x_ {1 })), ~ (x_ {2}, \, y_ {2}: = f (x_ {2})).}

{\ displaystyle \ left (x_ {0}, \, y_ {0}: = f (x_ {0}) \ right), ~ (x_ {1}, \, y_ {1}: = f (x_ {1 })), ~ (x_ {2}, \, y_ {2}: = f (x_ {2})).}

În acest fel, parabola este determinată în mod unic de următorul sistem de ecuații liniare:

{\begin{cases}y_{0}=A(x_{0}-x_{1})^{2}+B(x_{0}-x_{1})+C\\y_{1}=C\\y_{2}=A(x_{2}-x_{1})^{2}+B(x_{2}-x_{1})+C,\end{cases}}

{\ displaystyle {\ begin {cases} y_ {0} = A (x_ {0} -x_ {1}) ^ {2} + B (x_ {0} -x_ {1}) + C \\ y_ {1 } = C \\ y_ {2} = A (x_ {2} -x_ {1}) ^ {2} + B (x_ {2} -x_ {1}) + C, \ end {cases}}}

{\ displaystyle {\ begin {cases} y_ {0} = A (x_ {0} -x_ {1}) ^ {2} + B (x_ {0} -x_ {1}) + C \\ y_ {1 } = C \\ y_ {2} = A (x_ {2} -x_ {1}) ^ {2} + B (x_ {2} -x_ {1}) + C, \ end {cases}}}

din care rezultă:

A={\frac {y_{0}+y_{2}-2y_{1}}{2h^{2}}},\quad B={\frac {y_{2}-y_{0}}{2h}},\quad C=y_{1}.

{\ displaystyle A = {\ frac {y_ {0} + y_ {2} -2y_ {1}} {2h ^ {2}}}, \ quad B = {\ frac {y_ {2} -y_ {0} } {2h}}, \ quad C = y_ {1}.}

{\ displaystyle A = {\ frac {y_ {0} + y_ {2} -2y_ {1}} {2h ^ {2}}}, \ quad B = {\ frac {y_ {2} -y_ {0} } {2h}}, \ quad C = y_ {1}.}

Prin valoarea integralei acestui polinom $J_{1}^{'}$ ${\ displaystyle J_ {1} ^ {'}}$ ${\ displaystyle J_ {1} ^ {'}}$ avem:

\int _{x_{0}}^{x_{2}}\left[A\left(x-x_{1}\right)^{2}+B\left(x-x_{1}\right)+C\right]\,dx

{\ displaystyle \ int _ {x_ {0}} ^ {x_ {2}} \ left [A \ left (x-x_ {1} \ right) ^ {2} + B \ left (x-x_ {1} \ right) + C \ right] \, dx}

{\ displaystyle \ int _ {x_ {0}} ^ {x_ {2}} \ left [A \ left (x-x_ {1} \ right) ^ {2} + B \ left (x-x_ {1} \ right) + C \ right] \, dx}

=

\left[{\frac {A\left(x-x_{1}\right)^{3}}{3}}+{\frac {B\left(x^{2}-2x_{1}x+x_{1}^{2}\right)}{2}}+Cx\right]_{x_{0}}^{x_{2}}={\frac {2Ah^{3}}{3}}+2Ch.

{\ displaystyle \ left [{\ frac {A \ left (x-x_ {1} \ right) ^ {3}} {3}} + {\ frac {B \ left (x ^ {2} -2x_ {1 } x + x_ {1} ^ {2} \ right)} {2}} + Cx \ right] _ {x_ {0}} ^ {x_ {2}} = {\ frac {2Ah ^ {3}} { 3}} + 2Ch.}

{\ displaystyle \ left [{\ frac {A \ left (x-x_ {1} \ right) ^ {3}} {3}} + {\ frac {B \ left (x ^ {2} -2x_ {1 } x + x_ {1} ^ {2} \ right)} {2}} + Cx \ right] _ {x_ {0}} ^ {x_ {2}} = {\ frac {2Ah ^ {3}} { 3}} + 2Ch.}

Înlocuind valorile $LA, B.$ ${\ displaystyle A, B}$ $A, B$ Și $C.$ ${\ displaystyle C}$ $C.$ obținută din sistem, se obține valoarea aproximativă

J_{1}^{'}=\int _{x_{0}}^{x_{2}}f(x)\,dx={\frac {h}{3}}\left(y_{0}+4y_{1}+y_{2}\right).

{\ displaystyle J_ {1} ^ {'} = \ int _ {x_ {0}} ^ {x_ {2}} f (x) \, dx = {\ frac {h} {3}} \ left (y_ {0} + 4y_ {1} + y_ {2} \ dreapta).}

{\ displaystyle J_ {1} ^ {'} = \ int _ {x_ {0}} ^ {x_ {2}} f (x) \, dx = {\ frac {h} {3}} \ left (y_ {0} + 4y_ {1} + y_ {2} \ dreapta).}

Operăm într-un mod analog pentru calcularea integralelor polinoamelor în următoarele $m-1$ ${\ displaystyle m-1}$ $m-1$ intervale parțiale. Valorile obținute pe $m$ ${\ displaystyle m}$ $m$ intervale parțiale și pentru întregul interval de integrare și se obține o valoare aproximativă $J$ ${\ displaystyle J}$ $J$ pentru a se calcula integralul:

\int _{a}^{b}f(x)\,dx\simeq {\frac {h}{3}}\left[\left(y_{0}+4y_{1}+y_{2}\right)+\left(y_{2}+4y_{3}+y_{4}\right)+\ldots +\left(y_{n-2}+4y_{n-1}+y_{n}\right)\right].

{\ displaystyle \ int _ {a} ^ {b} f (x) \, dx \ simeq {\ frac {h} {3}} \ left [\ left (y_ {0} + 4y_ {1} + y_ { 2} \ right) + \ left (y_ {2} + 4y_ {3} + y_ {4} \ right) + \ ldots + \ left (y_ {n-2} + 4y_ {n-1} + y_ {n } \ corect corect].}

{\ displaystyle \ int _ {a} ^ {b} f (x) \, dx \ simeq {\ frac {h} {3}} \ left [\ left (y_ {0} + 4y_ {1} + y_ { 2} \ right) + \ left (y_ {2} + 4y_ {3} + y_ {4} \ right) + \ ldots + \ left (y_ {n-2} + 4y_ {n-1} + y_ {n } \ corect corect].}

Asa de:

J={\frac {h}{3}}\left(y_{0}+4y_{1}+2y_{2}+4y_{3}+2y_{4}+\ldots +2y_{n-2}+4y_{n-1}+y_{n}\right).

{\ displaystyle J = {\ frac {h} {3}} \ left (y_ {0} + 4y_ {1} + 2y_ {2} + 4y_ {3} + 2y_ {4} + \ ldots + 2y_ {n- 2} + 4y_ {n-1} + y_ {n} \ dreapta).}

{\ displaystyle J = {\ frac {h} {3}} \ left (y_ {0} + 4y_ {1} + 2y_ {2} + 4y_ {3} + 2y_ {4} + \ ldots + 2y_ {n- 2} + 4y_ {n-1} + y_ {n} \ dreapta).}

Această formulă se numește formula Knights-Simpson sau formula de cvadratură a parabolelor .

Greșeli

Metoda de cuadratură Knights-Simpson, ca orice metodă de aproximare numerică, este predispusă la erori . În plus față de eroarea datorată înlocuirii funcției integrand cu o succesiune de funcții aproximative, polinoamele pătratice, intrinseci metodei utilizate, întâmpină și erori datorate rotunjirii valorilor $y_{i}$ ${\ displaystyle y_ {i}}$ $y_ {i}$ care sunt calculate în mod concret cu instrumente care operează inevitabil cu o precizie limitată.

Pentru a reduce la minim acest lucru, este recomandabil:

alegeți un pas de integrare $h$ ${\ displaystyle h}$ $h$ cu un număr finit de zecimale;
efectuați calcule cu cel puțin dublul numărului de zecimale pe care le doriți exact în rezultat.

Indicarea erorii intrinseci a metodei $e:=J-I$ ${\ displaystyle e: = JI}$ ${\ displaystyle e: = J-I}$ , se poate arăta că:

e\simeq kh^{4},

{\ displaystyle și \ simeq kh ^ {4},}

{\ displaystyle și \ simeq kh ^ {4},}

unde este $k$ ${\ displaystyle k}$ $k$ este o constantă care depinde de funcția integrand și de intervalul de integrare. Prin urmare, regula Knights-Simpson este o metodă de ordinul patru.

Poate fi foarte util să cunoașteți o suprataxă pentru această eroare. Evaluarea exactă a acestei creșteri nu este ușoară, deoarece necesită calcularea celei de-a patra derivate a funcției integrand. Pentru multe funcții date analitic, calculul celei de-a patra derivate este foarte oneros; pentru funcții cunoscute empiric, aceeași evaluare a derivatei a patra constituie în sine o problemă de calcul numeric care tinde să fie oneroasă. În consecință, în general, metodele empirice sunt preferate pentru a evalua eroarea: cea mai cunoscută și cea mai utilizată este metoda de înjumătățire a pasului .

Din cele observate anterior rezultă că prin aplicarea metodei Cavalieri-Simpson cu o etapă de integrare $h,$ ${\ displaystyle h,}$ ${\ displaystyle h,}$ obținem aproximarea cu care indicăm acum $J_{(1)}$ ${\ displaystyle J _ {(1)}}$ ${\ displaystyle J _ {(1)}}$ afectate de o eroare pe care o scriem $e_{(1)}=J_{(1)}-I\simeq kh^{4}.$ ${\ displaystyle e _ {(1)} = J _ {(1)} - I \ simeq kh ^ {4}.}$ ${\ displaystyle e _ {(1)} = J _ {(1)} - I \ simeq kh ^ {4}.}$

Folosind pasul de integrare ${\frac {h}{2}}$ ${\ displaystyle {\ frac {h} {2}}}$ ${\ frac {h} {2}}$ , se va obține valoarea aproximativă $J_{(2)}$ ${\ displaystyle J _ {(2)}}$ ${\ displaystyle J _ {(2)}}$ cu eroare: $e_{(2)}=J_{(2)}-I\simeq k{\frac {h^{4}}{16}}.$ ${\ displaystyle e _ {(2)} = J _ {(2)} - I \ simeq k {\ frac {h ^ {4}} {16}}.}$ ${\ displaystyle e _ {(2)} = J _ {(2)} - I \ simeq k {\ frac {h ^ {4}} {16}}.}$

Din aceste rapoarte rezultă

e_{(1)}-e_{(2)}=J_{(1)}-J_{(2)}\simeq k\left(h^{4}-{\frac {h^{4}}{16}}\right),

{\ displaystyle e _ {(1)} - e _ {(2)} = J _ {(1)} - J _ {(2)} \ simeq k \ left (h ^ {4} - {\ frac { h ^ {4}} {16}} \ dreapta),}

{\ displaystyle e _ {(1)} - e _ {(2)} = J _ {(1)} - J _ {(2)} \ simeq k \ left (h ^ {4} - {\ frac { h ^ {4}} {16}} \ dreapta),}

din care rezultă:

k\simeq {\frac {16\left(J_{(1)}-J_{(2)}\right)}{15h^{4}}}.

{\ displaystyle k \ simeq {\ frac {16 \ left (J _ {(1)} - J _ {(2)} \ right)} {15h ^ {4}}}.}

{\ displaystyle k \ simeq {\ frac {16 \ left (J _ {(1)} - J _ {(2)} \ right)} {15h ^ {4}}}.}

Deoarece, neglijând aproximările de rotunjire, cea mai bună aproximare este dată de $J_{(2)},$ ${\ displaystyle J _ {(2)},}$ ${\ displaystyle J _ {(2)},}$ substituind valoarea lui $k$ ${\ displaystyle k}$ $k$ în $e_{(2)}\simeq k{\frac {h^{4}}{16}},$ ${\ displaystyle e _ {(2)} \ simeq k {\ frac {h ^ {4}} {16}},}$ ${\ displaystyle e _ {(2)} \ simeq k {\ frac {h ^ {4}} {16}},}$ primesti:

e_{(2)}\simeq {\frac {16\left(J_{(1)}-J_{(2)}\right)}{15h^{4}}}\cdot {\frac {h^{4}}{16}}

{\ displaystyle e _ {(2)} \ simeq {\ frac {16 \ left (J _ {(1)} - J _ {(2)} \ right)} {15h ^ {4}}} \ cdot { \ frac {h ^ {4}} {16}}}

{\ displaystyle e _ {(2)} \ simeq {\ frac {16 \ left (J _ {(1)} - J _ {(2)} \ right)} {15h ^ {4}}} \ cdot { \ frac {h ^ {4}} {16}}}

de la care

e_{(2)}\simeq {\frac {J_{(1)}-J_{(2)}}{15}}.

{\ displaystyle e _ {(2)} \ simeq {\ frac {J _ {(1)} - J _ {(2)}} {15}}.}

{\ displaystyle e _ {(2)} \ simeq {\ frac {J _ {(1)} - J _ {(2)}} {15}}.}

Prin urmare, putem asuma valoarea absolută a $e_{(2)}$ ${\ displaystyle e _ {(2)}}$ ${\ displaystyle e _ {(2)}}$

{\frac {\left|J_{(1)}-J_{(2)}\right|}{15}}

{\ displaystyle {\ frac {\ left | J _ {(1)} - J _ {(2)} \ right |} {15}}}

{\ displaystyle {\ frac {\ left | J _ {(1)} - J _ {(2)} \ right |} {15}}}

ca o creștere a erorii absolute

e_{a}=\left|J_{(2)}-I\right|.

{\ displaystyle e_ {a} = \ left | J _ {(2)} - I \ right |.}

{\ displaystyle e_ {a} = \ left | J _ {(2)} - I \ right |.}

Este interesant de observat că dacă aproximările $J_{(1)}$ ${\ displaystyle J _ {(1)}}$ ${\ displaystyle J _ {(1)}}$ Și $J_{(2)}$ ${\ displaystyle J _ {(2)}}$ ${\ displaystyle J _ {(2)}}$ coincid pentru prima $r$ ${\ displaystyle r}$ $r$ cifre zecimale, rezultă:

\left|J_{(1)}-J_{(2)}\right|<10^{-r},

{\ displaystyle \ left | J _ {(1)} - J _ {(2)} \ right | <10 ^ {- r},}

{\ displaystyle \ left | J _ {(1)} - J _ {(2)} \ right | <10 ^ {- r},}

din care rezultă

{\frac {\left|J_{(1)}-J_{(2)}\right|}{15}}<{\frac {10^{-r}}{15}}=0,666\ldots \cdot 10^{-(r+1)},

{\ displaystyle {\ frac {\ left | J _ {(1)} - J _ {(2)} \ right |} {15}} <{\ frac {10 ^ {- r}} {15}} = 0,666 \ ldots \ cdot 10 ^ {- (r + 1)},}

{\ displaystyle {\ frac {\ left | J _ {(1)} - J _ {(2)} \ right |} {15}} <{\ frac {10 ^ {- r}} {15}} = 0,666 \ ldots \ cdot 10 ^ {- (r + 1)},}

ceea ce înseamnă că primul $r$ ${\ displaystyle r}$ $r$ zecimale nu sunt afectate de eroare.

Prin urmare, se poate concluziona că dacă două aproximări ale unei integrale, dintre care a doua se obține prin înjumătățirea etapei de integrare folosite pentru a calcula prima, coincid pentru prima $r$ ${\ displaystyle r}$ $r$ cifre zecimale, aceste cifre pot fi considerate exacte.

Mai general dacă doriți să cunoașteți o aproximare a unei integrale cu garanția acurateței pentru un număr dat $s$ ${\ displaystyle s}$ $s$ de cifre zecimale, trebuie calculat un anumit număr de aproximări succesive, înjumătățind pasul de fiecare dată, până se obțin două care coincid pentru $s$ ${\ displaystyle s}$ $s$ cifre.

Observăm că se poate întâmpla să ajungă la prima pereche de valori aproximative satisfăcătoare care au mai mult de s cifre coincidente. De asemenea, observăm că, procedând la reducerea amplitudinii subintervalelor, pe lângă timpul de calcul mai mare necesar, erorile de rotunjire pot fi departe de a fi neglijabile din cauza creșterii numărului de operații necesare; poate duce, de asemenea, la situații în care eroarea generală crește odată cu reducerea etapei de integrare. Din acest motiv, au fost concepute metode de integrare numerică „adaptive” (sau adaptive) care măresc numărul de sub-intervale numai în zonele indicate de un test de eroare specific.

Notă istorică

Regula este unul dintre cele mai bune exemple ale legii eponimiei lui Stigler : se pare că această regulă era deja cunoscută de Torricelli , în timp ce Cavalieri și- a demonstrat deja formularea geometrică în 1635 . Kepler îl tratase și el: de fapt, multe texte germane îl numesc Keplersche Fassregel cu 100 de ani înainte de Simpson. Cu toate acestea, la nivel internațional astăzi este raportat conform desemnării anglo-saxone ca regulă a lui Simpson : chiar și în acest mediu fusese deja folosit, de exemplu, de Gregory . Pe de altă parte, spre deosebire de predecesorii săi, Thomas Simpson astăzi este în esență cunoscut pentru această regulă, în ciuda faptului că are doar meritul de a formaliza o metodă deja cunoscută.