Il metodo di Numerov

Vogliamo ora considerare il problema della risoluzione dell'equazione di Schrödinger indipendente dal tempo e unidimensionale in maniera numerica. Questo ci permetterà di apprendere la metodologia generale da applicare per casi specifici (ad esempio quello dell'oscillatore armonico), e di comprendere la potenza e le limitazioni del metodo numerico.

Il metodo di Numerov è utile per integrare equazioni differenziali del secondo ordine della forma generale

$$\frac{d^2 y}{dx^2} = -g(x) y(x) + s(x)$$ (2.32)

dove $g(x)$ e $s(x)$ sono funzioni date, e condizioni iniziali della forma

$$y(x_0)=y_0 , y'(x_0)=y'_0$$ (2.33)

L'equazione di Schrödinger (2.2) ha questa forma, con $g(x)\equiv \frac{2m}{\hbar^2}[ E - V(x) ]$ e $s(x) = 0$. Si vedrà in seguito che anche le parti radiali di equazioni di Schrödinger in tre dimensioni a simmetria sferica appartengono a questa classe. Un'altra importante equazione che ricade in questa categoria è l'equazione di Poisson dell'elettromagnetismo,

$$\frac{d^2 \phi}{dx^2} = -4\pi \rho(x)$$ (2.34)

dove $\rho(x)$ è una densità di carica. In questo caso $g(x)=0$ e $s(x)=-4\pi \rho(x)$.

La metodologia è simile a quella dell'algoritmo di Størmer-Verlet per seguire l'evoluzione temporale di un punto materiale, sostituendo la coordinata spaziale a quella temporale.

Dividiamo dunque l'intervallo spaziale di interesse in $N$ intervallini di ampiezza $\Delta x$, siano $x_i$ i nodi della griglia così ottenuta e $y_i=y(x_i)$ i valori della funzione incognita $y(x)$ in corrispondenza di tali punti. Analogamente indichiamo con $g_i$ e $s_i$ i valori delle funzioni (date) $g(x)$ e $s(x)$ negli stessi punti. Al fine di ottenere una equazione alle differenze finite espandiamo in serie di Taylor attorno ad un punto $x_n$, spingendoci fino al quinto ordine:

$$\begin{array}{ll} y_{n-1} = & y_n - y'_n \Delta x + \frac{1}... ...20} y'''''_n (\Delta x)^5 \\ & + O[(\Delta x)^6] \end{array}$$ (2.35)

Sommiamo le due equazioni:

$$y_{n+1} + y_{n-1} = 2y_n + y''_n (\Delta x)^2 + \frac{1}{12} y''''_n (\Delta x)^4 + O[(\Delta x)^6]$$ (2.36)

La (2.32) ci dice che

$$y''_n = -g_n y_n + s_n$$ (2.37)

Inoltre, indicando temporaneamente con $z_n$ questa quantità, sarà anche vero

$$z_{n+1} + z_{n-1} = 2z_n + z''_n (\Delta x)^2 + O[(\Delta x)^4]$$ (2.38)

e quindi

$$y''''_n \equiv z''_n = \frac{z_{n+1} + z_{n-1} - 2z_n}{(\Delta x)^2} + O[(\Delta x)^2]$$ (2.39)

Inserendo questi risultati nella (2.36)

$$\begin{array}{ll} y_{n+1} & = 2y_n - y_{n-1} + ( -g_n y_n +... ... y_n - 2 s_n ) (\Delta x)^2 \\ & + O[(\Delta x)^6] \end{array}$$ (2.40)

da cui la formula di Numerov

$$\begin{array}{ll} y_{n+1} \left[ 1 + g_{n+1} \frac{(\Delta x... ..._{n-1} ) \frac{(\Delta x)^2}{12} + O[(\Delta x)^6] \end{array}$$ (2.41)

che permette di ottenere $y_{n+1}$ a partire da $y_n$ e $y_{n-1}$ e quindi ricorsivamente--dalle condizioni iniziali date--tutta la funzione sull'intervallo di interesse.

Dalle condizioni iniziali (2.33) è ovviamente possibile integrare muovendosi sia nella direzione degli $x$ positivi che in quella degli $x$ negativi, e in presenza di simmetria rispetto ad un punto di inversione basterà integrare in una direzione sola.

Nel caso di nostro interesse--l'equazione di Schrödinger--tutti i termini $s_n$ sono assenti, e in pratica è conveniente porre

$$f_n \equiv 1 + g_n \frac{(\Delta x)^2}{12}$$ (2.42)

dove

$$g_n = \frac{2m}{\hbar^2} [ E - V(x_n) ]$$ (2.43)

Con questa posizione la formula di Numerov si riduce a

$$y_{n+1} = \frac{(12-10f_n) y_n - f_{n-1}y_{n-1}}{f_{n+1}}$$ (2.44)