Modern numerikus módszerek
Bene Gyula
Eötvös Loránd Tudományegyetem, Elméleti Fizikai Tanszék
1117 Budapest, Pázmány Péter sétány 1/A
9. hét csütörtöki előadás
Speciális függvények értékeinek kiszámítása
Differenciálegyenletek, hatványsorok, aszimptotikus sorok
Rekurziós formulák
Pl.:
Legendre-polinomok:
$$(n+1)P_{n+1}(x)=(2n+1)xP_n(x)-nP_{n-1}(x)$$
Bessel-függvények:
$$J_{n+1}(x)=\frac{2n}{x}J_n(x)-J_{n-1}(x)$$
Exponenciális integrálok:
$$nE_{n+1}(x)={\rm e}^{-x}-xE_n(x)$$
Koszinusz- és szinuszfüggvények:
$$\cos nx=2\cos x \cos (n-1)x-\cos(n-2)x$$
$$\sin nx=2\cos x \sin (n-1)x-\sin(n-2)x$$
Általában egy
$$y_{n+1}+a_ny_{n}+b_ny_{n-1}=0$$
típusú rekurziós összefüggésnek két lineárisan független megoldása van,
jelöljük ezeket $f_n$-nel
és $g_n$-nel. Ha
$$\lim_{n\rightarrow \infty}\frac{f_n}{g_n}=0\;,$$
akkor $f_n$ a minimális, $g_n$ pedig a domináns megoldás. Növekvő $n$-ek
irányában alkalmazva a rekurziót, a domináns megoldást kapjuk, csökkenő $n$-ek
esetén pedig a minimális megoldást.
Lánctörtek
$$b_0+\frac{a_1}{b_1+\frac{a_2}{b_2+\frac{a_3}{b_3+\dots}}}$$
Kényelmesebb jelölés:
$$b_0+\frac{a_1}{b_1+}\;\frac{a_2}{b_2+}\;\frac{a_3}{b_3+}\dots$$
Pl.:
$$\frac{\sqrt{5}-1}{2}=1+\frac{1}{1+}\;\frac{1}{1+}\;\frac{1}{1+}\dots$$
$$\tan x =\frac{x}{1-}\;\frac{x^2}{3-}\;\frac{x^2}{5-}\;\frac{x^2}{7-}\dots$$
Lánctörtek kiszámítása:
$$t_n\equiv b_0+\frac{a_1}{b_1+}\;\frac{a_2}{b_2+}\dots\frac{a_n}{b_n}
=\frac{A_n}{B_n}$$
Ekkor az $A_n$, $B_n$ mennyiségek a következő rekurzióval kaphatók meg
(J.Wallis, 1655):
$$A_{-1}=1\;,\quad B_{-1}=0 \;,\quad A_{0}=b_0\;,\quad B_{0}=1$$
$$A_j=b_jA_{j-1}+a_jA_{j-2}\;,\quad B_j=b_jB_{j-1}+a_jB_{j-2}\quad (j\ge 1)$$
Ui.
$$\left(\begin{array}{c}A_{j-2}\\B_{j-2}\end{array}\right)=\left(\begin{array}{cc}\alpha&\beta\\\gamma&\delta\end{array}\right)\left(\begin{array}{c}a_{j-2}\\b_{j-2}\end{array}\right)$$
$$\left(\begin{array}{c}A_{j-1}\\B_{j-1}\end{array}\right)=\left(\begin{array}{cc}\alpha&\beta\\\gamma&\delta\end{array}\right)\left(\begin{array}{c}a_{j-2}b_{j-1}\\b_{j-2}b_{j-1}+a_{j-1}\end{array}\right)=\left(\begin{array}{cc}\alpha&\beta\\\gamma&\delta\end{array}\right)\left(\begin{array}{cc}0&
a_{j-2}\\1&
b_{j-2}\end{array}\right)\left(\begin{array}{c}a_{j-1}\\b_{j-1}\end{array}\right)$$
$$\left(\begin{array}{c}A_{j}\\B_{j}\end{array}\right)=\left(\begin{array}{cc}\alpha&\beta\\\gamma&\delta\end{array}\right)\left(\begin{array}{cc}0&
a_{j-2}\\1& b_{j-2}\end{array}\right)\left(\begin{array}{cc}0&
a_{j-1}\\1& b_{j-1}\end{array}\right)\left(\begin{array}{c}a_{j}\\b_{j}\end{array}\right)=\left(\begin{array}{cc}A_{j-2}&
A_{j-1}\\ B_{j-2}& B_{j-1}\end{array}\right)\left(\begin{array}{c}a_{j}\\b_{j}\end{array}\right)$$
Lánctörtek numerikus meghatározása (Lentz módszere):
Legyen
$$C_j=\frac{A_j}{A_{j-1}}\;,\quad D_j=\frac{B_{j-1}}{B_j}\;.$$
Ezzel
$$t_j=C_jD_jt_{j-1}\;.$$
A $C_j$, $D_j$ mennyiségekre
$$C_j=\frac{b_jA_{j-1}+a_jA_{j-2}}{A_{j-1}}=b_j+\frac{a_j}{C_{j-1}}$$
és
$$D_j=\frac{B_{j-1}}{b_jB_{j-1}+a_jB_{j-2}}=\frac{1}{b_j+a_jD_{j-1}}$$
teljesül
$$C_0=b_0\;,\quad D_0=0\;,\quad t_0=b_0$$
kezdeti feltételekkel. Ha valamelyik nevező nullává válna, akkor egy kis
számmal (pl. $10^{-30}$-cal) helyettesítjük, a beavatkozás a nullával osztást
kivédi, egy újabb ciklus múlva pedig a hatása elenyészik.
Lánctört kifejtés a rekurziós formula alapján:
$$y_{n+1}+a_ny_{n}+b_ny_{n-1}=0$$
$$\frac{y_{n+1}}{y_n}+a_n+b_n\frac{y_{n-1}}{y_n}=0$$
$$\frac{y_n}{y_{n-1}}=-\frac{b_n}{a_n+\frac{y_{n+1}}{y_n}}=-\frac{b_n}{a_n-}\;\frac{b_{n+1}}{a_{n+1}-}\dots$$
Pincherle tétele: a fenti lánctört akkor és csak akkor konvergál, ha létezik a
rekurziónak az $f_n$ minimális megoldása. Ebben az esetben a lánctört értéke
$f_n/f_{n-1}$.
Gamma-függvény, béta-függvény, faktoriális, binomiális
együtthatók
$$\Gamma(z)=\int_0^\infty {\rm e}^{-t}\;t^{z-1}dt$$
Néhány fontos összefüggés:
$$\Gamma(n+1)=n!$$
$$\Gamma(z+1)=z\Gamma(z)$$
Bizonyítás: parciális integrálással.
Reflexiós összefüggés:
$$\Gamma(1-z)=\frac{\pi}{\Gamma(z)\sin(\pi z)}$$
Közelítő formula a numerikus kiértékeléshez (az $Re(z)>0$ komplex félsíkon):
$\begin{eqnarray}\Gamma(z)=\sqrt{2\pi}\left(z+\gamma+\frac{1}{2}\right)^{z+\frac{1}{2}}{\rm
e}^{-\left(z+\gamma+\frac{1}{2}\right)}\left[c_0+\frac{c_1}{z+1}+\frac{c_2}{z+2}+\cdots +\frac{c_N}{z+N}+\epsilon\right]
\end{eqnarray}$
Itt $\gamma$ egész szám. $6$ tag már elegendő $|\epsilon|<2\times 10^{-10}$
pontossághoz
az $Re(z)>0$
komplex
félsíkon
mindenütt. Gyakorlatban
a függvény
gyors
növekedése
miatt
célszerűbb a
logaritmusát
kiszámítani.
Euler-féle béta-függvény:
$$B(z,w)=B(w,z)=\int_0^1t^{z-1}(1-t)^{w-1}dt$$
$$B(z,w)=\frac{\Gamma(z)\Gamma(w)}{\Gamma(z+w)}$$
Bizonyítás:
$$\Gamma(z)\Gamma(w)=\int_0^\infty\int_0^\infty {\rm
e}^{-(t+u)}\;t^{z-1}u^{w-1}dtdu$$
Változócsere:
$$t=xy\;,\quad u=x(1-y)$$
Belátható, hogy $x\in (0,\infty)$ és $y\in (0,1)$, és
$$dt du = x dx dy$$
Ezzel
$$\Gamma(z)\Gamma(w)=\int_0^\infty dx\int_0^1 dy\; x\; {\rm
e}^{-x} x^{z-1}y^{z-1}x^{w-1}(1-y)^{w-1}=\underbrace{\int_0^\infty dx\;{\rm
e}^{-x} \;x^{(z+w)-1}}_{\Gamma(z+w)}\underbrace{\int_0^1 dy\;y^{z-1}(1-y)^{w-1}}_{B(z,w)}$$QED.
Nem teljes gamma-függvény és alkalmazásai: hibafüggvény,
$\chi^2$-eloszlás, Poisson-eloszlás, exponenciális
integrálok
$$\gamma(a,x)=\int_0^x {\rm e}^{-t}\;t^{a-1}dt$$
$$\Gamma(a,x)=\int_x^\infty {\rm e}^{-t}\;t^{a-1}dt$$
$$P(a,x)=\frac{\gamma(x,a)}{\Gamma(a)}$$
$$Q(a,x)=\frac{\Gamma(x,a)}{\Gamma(a)}=1-P(a,x)$$
Sorfejtés ($x< a+1$ esetén gyorsan konvergál):
$$\gamma(a,x)={\rm e}^{-x}x^a\sum_{n_0}^\infty
\frac{\Gamma(a)}{\Gamma(a+1+n)}x^n$$
Lánctört-kifejtés ($x> a+1$ esetén gyorsan konvergál):
$$\Gamma(a,x)={\rm
e}^{-x}x^a\left(\frac{1}{x+}\;\frac{1-a}{1+}\;\frac{1}{x+}\;\frac{2-a}{1+}\;\cdots\right)$$
A gyorsabb konvergencia érdekében célszerű a lánctört páros részét (melynek
közelítő törtjei $t_{2n}$-ek, ha az eredetié $t_n$-ek) venni.
Hibafüggvény:
$${\rm erf}(x)=\frac{2}{\sqrt{\pi}}\int_0^x{\rm e}^{-t^2}dt$$
$${\rm erfc}(x)=1-{\rm erf}(x)=\frac{2}{\sqrt{\pi}}\int_x^\infty{\rm
e}^{-t^2}dt$$
Kapcsolat a nem teljes gamma-függvénnyel:
$${\rm erf}(x)=P\left(\frac{1}{2},x^2\right)$$
$${\rm erfc}(x)=Q\left(\frac{1}{2},x^2\right)$$
Poisson-eloszlás és a nem teljes gamma-függvény:
$$\sum_{n=0}^{k-1}{\rm e}^{-\lambda}\; \frac{\lambda^n}{n!}=Q(k,\lambda)$$
$\chi^2$-eloszlás és a nem teljes gamma-függvény:
$$P(\chi^2|\nu)=P\left(\frac{\nu}{2},\frac{\chi^2}{2}\right)$$
Exponenciális integrálok:
$$E_n(x)=\int_1^\infty \frac{{\rm e}^{-xt}}{t^n}dt=x^{n-1}\Gamma(1-n,x)$$
Nem teljes béta-függvény és alkalmazásai: Student-eloszlás,
F-eloszlás, binomiális eloszlás
$$I_x(a,b)=\frac{B_x(a,b)}{B(a,b)}=\frac{\int_0^xt^{a-1}(1-t)^{b-1}dt}{\int_0^1t^{a-1}(1-t)^{b-1}dt}$$
Nyilván
$$I_x(a,b)=1-I_{1-x}(b,a)\;.$$
Lánctört-kifejtés:
$\begin{eqnarray}I_x(a,b)&=&\frac{x^a(1-x)^b}{aB(a,b)}\left[\frac{1}{1+}\;\frac{d_1}{1+}\;\frac{d_2}{1+}\;\cdots\right]\\
d_{2m+1}&=&-\frac{(a+m)(a+b+m)x}{(a+2m)(a+2m+1)}\\
d_{2m}&=&\frac{m(b-m)x}{(a+2m)(a+2m-1)}\end{eqnarray}$
Gyorsan konvergál, ha $x<(a+1)/(a+b+2)$. Ha ez nem teljesül, előbb a fenti reflexiós
összefüggést alkalmazzuk.
Student-eloszlás:
$$A(t|\nu)=\frac{1}{\nu^{1/2}B\left(\frac{1}{2},\frac{\nu}{2}\right)}\int_{-t}^t\left(1+\frac{x^2}{\nu}\right)^{-\frac{\nu+1}{2}}dx=1-I_{\frac{\nu}{\nu+t^2}}\left(\frac{\nu}{2},\frac{1}{2}\right)$$
F-eloszlás:
$$Q(F|\nu_1,\nu_2)=I_{\frac{\nu_2}{\nu_2+\nu_1
F}}\left(\frac{\nu_2}{2},\frac{\nu_1}{2}\right)$$
Binomiális eloszlás:
$$\sum_{j=k}^n\left(\begin{array}{c}n\\j\end{array}\right)p^j(1-p)^{n-j}=I_p(k,n-k+1)$$
Bessel-függvények
Gömbharmonikusok
Elliptikus integrálok és elliptikus függvények
Hipergeometrikus függvény
$$\phantom{a}_{2}F_1(a,b,c,z)=1+\frac{ab}{c}z+\frac{a(a+1)b(b+1)}{c(c+1)}\frac{z^2}{2!}+\cdots
+\frac{a(a+1)\cdots (a+j-1)b(b+1)\cdots (b+j-1)}{c(c+1)\cdots
(c+j-1)}\frac{z^j}{j!}+\cdots$$
Differenciálegyenlete:
$$z(z-1)F''=abF-\left[c-(a+b+1)z\right]F'$$
bene@arpad.elte.hu