Programma: hydrogen_gauss

Il programma hydrogen_gauss risolve il problema secolare per l'atomo di idrogeno utilizzando due diverse basi non ortonormali:

1. una base gaussiana ``onda S'':
   

   $$b_i({\bf r}) = {\rm e}^{-\alpha_i r^2}$$ (5.91)

   
   

2. una base gaussiana ``onda P'':
   

   $$b_i({\bf r}) = x {\rm e}^{-\alpha_i r^2}$$ (5.92)

   $$b_i({\bf r}) = y {\rm e}^{-\alpha_i r^2}$$ (5.93)

   $$b_i({\bf r}) = z {\rm e}^{-\alpha_i r^2}$$ (5.94)

   
   

L'operatore hamiltoniano del problema è ovviamente (in unità atomiche con energie date in Rydberg)

$$H = -\nabla^2 -\frac{2Z}{r}$$ (5.95)

e per l'atomo di idrogeno $Z=1$.

Il programma permette di assegnare gli esponenti $\alpha_i$, sotto forma di larghezze $\lambda_i$ delle gaussiane:

$$\alpha_i = \left( \frac{Z}{\lambda_i} \right)^2$$ (5.96)

Gli esponenti vengono fatti corrispondere a una griglia di $N$ valori di $\lambda$ equispaziati da un minimo ad un massimo specificati. I calcoli per l'onda S e l'onda P sono del tutto indipendenti: si tratta di due basi distinte. La base in onda P è chiaramente inadatta a descrivere lo stato fondamentale, perchè non è possibile riprodurre il giusto andamento per piccoli $r$, e viene inclusa a scopo didattico.

Il codice procede quindi a valutare tutti gli elementi delle matrici $H_{ij}$ e $S_{ij}$. Il calcolo è basato sulle espressioni degli integrali calcolati analiticamente. In particolare si trova per l'onda S

$$S_{ij} = \int {\rm e}^{-(\alpha_i+\alpha_j)r^2} d^3{\bf r} = \left(\frac{\pi}{\alpha_i+\alpha_j}\right)^{3/2}$$ (5.97)

e inoltre i termini cinetico e coulombiano di $H_{ij}$ sono rispettivamente

$$H^K_{ij} = \int {\rm e}^{-\alpha_i r^2} [-\nabla^2] {\rm e}^... ..._i+\alpha_j} \left(\frac{\pi}{\alpha_i+\alpha_j}\right)^{3/2}$$ (5.98)

$$H^V_{ij} = \int {\rm e}^{-\alpha_i r^2} \left[-\frac{2Z}{r}\... ...{-\alpha_j r^2} d^3{\bf r} = -\frac{2\pi Z}{\alpha_i+\alpha_j}$$ (5.99)

Per l'onda P si procede analogamente, utilizzando le corrispondenti espressioni per gli integrali.

Il codice procede dunque a chiamare la subroutine diaghg che risolve il problema secolare generalizzato (ossia applica il principio variazionale), ritornando il vettore e contenente gli autovalori in ordine crescente di energia, e la matrice v contenente gli autovettori, ossia i coefficienti dell'espansione.

Al suo interno, diaghg effettua il calcolo in due stadi descritto nella sezione precedente, delegando la risoluzione del problema agli autovalori semplice alla subroutine dsyev della libreria di algebra lineare LAPACK, ed utilizzando alcune chiamate addizionali a subroutines della libreria BLAS. Informazioni addizionali su queste librerie sono disponibili a partire dalla pagina del corso lapack.html.

Da notare che dopo la prima diagonalizzazione (relativa a $S$) diaghg getta via tutti gli autovettori corrispondenti ad autovalori nulli (entro una soglia di tolleranza numerica). Questi corrispondono a combinazioni linearmente dipendenti dalle altre. Il numero di vettori di base linearmente indipendenti trovati viene riportato nell'output.

Il programma procede quindi a scrivere i coefficienti dell'espansione nel file s-coeff.out (p-coeff.out): per ogni funzione $j$, e la funzione d'onda dello stato fondamentale nel file s-wfc.out (p-wfc.out).