ZCE 111 Computational Physics Semester Project
Ng Wei Chun, 105525

Solving the Schrödinger equation for Helium Atom with Hatree-Fock Method.

Abstract

The Hartree-Fock method is a method of approximation for solving the Hamiltonian of a many-body system. It is based on a simple approximation to the many-body wave function: that the wave function as given by Slater determinant. The wave function captures much of the physics required for accurate solution of the Hamiltonian, while being anti-symmetric with respect to an exchange correlation. This is as well required by the Pauli exclusion principle. When solving many-body system with computation, seek ahead a self-consistent eigenvalue problem, which can be solve by an some self-consistency iterative procedure. In the case of Hartree approximation, this involved solving a generalised eigenvalue problem, where the normalisation is govern by an overlap matrix of the chosen basis. The choice of basis functions ultimately alter the outcome and included with it are complication when dealing with pre-factors (Subscript[C, p]) and exponents parameters, (Subscript[α, p]).

Helium ground state, self-consistent approximation

Helium ground state with Gaussian s-basis function

(*initialise*)
α[1]=0.298073;
α[2]=1.242567;
α[3]=5.782948;
α[4]=38.474970;

(*in-loop parameters*)
tol=0.00000001;n=0;Enew=100*tol;ΔE=100*tol;

(*fixed matrix components*)
(*χ[p_,vecr_]:=Exp[-α[p]*vecr^2];*)
For[p=1,p≤4,p++,
For[q=1,q≤4,q++,
h[p,q]=Integrate[4π*Exp[-vecr^2*(α[p]+α[q])]*(3*α[q]*vecr^2-2*α[q]^2*vecr^4-2*vecr),{vecr,0,999}];
S[p,q]=Integrate[4π*vecr^2*Exp[-vecr^2*(α[p]+α[q])],{vecr,0,999}];
For[r=1,r≤4,r++,
For[s=1,s≤4,s++,
Q[p,r,q,s]=((2(π^(5/2)))/((α[p]+α[q])(α[r]+α[s])(Sqrt[α[p]+α[q]+α[r]+α[s]])));
];];];];


(*generate initial C_{p}*)
Cnew={1,1,1,1};

Print["loop counts = ", Dynamic[n]];(*tracer*)

(*relaxation looping*)
While[ΔE>tol,
n++;
Eold=Enew;

(*replace matrix C_{p}*)
For[i=1,i≤4,i++,
Cc[i]=Cnew[[i]];
];

(*normalise C_{p}*)
norm=Sum[Cc[p]*S[p,q]*Cc[q],{p,1,4},{q,1,4}];
For[i=1,i≤4,i++,
Cc[i]=Cc[i]/Sqrt[norm];
];(*end For i*)

(*constructing F_{pq}*)
For[p=1,p≤4,p++,
For[q=1,q≤4,q++,
F[p,q]=h[p,q]+Sum[Q[p,r,q,s]*Cc[r]*Cc[s],{r,1,4},{s,1,4}];
];];

(*matrix*)
mF=Table[F[p,q],{p,1,4},{q,1,4}];
mS=Table[S[p,q],{p,1,4},{q,1,4}];
(*Eigensystem[{mF,mS}];*)(*solving the generalised eigenvalue problem*)
Enew=Min[Eigenvalues[{mF,mS}]];
For[i=1,i≤4,i++,
If[Enew==Eigenvalues[{mF,mS}][[i]],Cnew=Eigenvectors[{mF,mS}][[i]];];
];(*end For i*)
ΔE=Abs[Eold-Enew];

];(*end While*)

EG=2*Sum[Cc[p]*Cc[q]*h[p,q],{p,1,4},{q,1,4}]+Sum[Q[p,r,q,s]*Cc[p]*Cc[q]*Cc[r]*Cc[s],{p,1,4},{q,1,4},{r,1,4},{s,1,4}];
Eth=-2.85516038;(*given Eth from text book*)
Print["Resulting ground state energy for Helium atom is E_{G} = ", EG, " a.u."];
Print["in atomic energy unit of Hartree; % of deviation is ", 100*(EG-Eth)/Eth, "%; Time used = ", TimeUsed[], "s"];

General Hatree-Fock method

Helium ground state with Gaussian s-basis function with exchange correlation

(*initialise*)
α[1]=0.298073;
α[2]=1.242567;
α[3]=5.782948;
α[4]=38.474970;

(*in-loop parameters*)
tol=0.00000001;n=0;Enew=100*tol;ΔE=100*tol;

(*fixed matrix components*)
(*χ[p_,vecr_]:=Exp[-α[p]*vecr^2];*)
For[p=1,p≤4,p++,
For[q=1,q≤4,q++,
h[p,q]=Integrate[4π*Exp[-vecr^2*(α[p]+α[q])]*(3*α[q]*vecr^2-2*α[q]^2*vecr^4-2*vecr),{vecr,0,999}];
S[p,q]=Integrate[4π*vecr^2*Exp[-vecr^2*(α[p]+α[q])],{vecr,0,999}];
For[r=1,r≤4,r++,
For[s=1,s≤4,s++,
(*exchange G_{pq} matrix element*)
GJ[p,r,q,s]=((2(π^(5/2)))/((α[p]+α[q])(α[r]+α[s])(Sqrt[α[p]+α[q]+α[r]+α[s]])));
GK[p,r,s,q]=((2(π^(5/2)))/((α[p]+α[s])(α[r]+α[q])(Sqrt[α[p]+α[q]+α[r]+α[s]])));
];];];];


(*generate initial C_{p}*)
Cnew={1,1,1,1};

Print["loop counts = ", Dynamic[n]];(*tracer*)

(*relaxation looping*)
While[ΔE>tol,
n++;
Eold=Enew;

(*replace matrix C_{p}*)
For[i=1,i≤4,i++,
Cc[i]=Cnew[[i]];
];

(*normalise C_{p}*)
norm=Sum[Cc[p]*S[p,q]*Cc[q],{p,1,4},{q,1,4}];
For[i=1,i≤4,i++,
Cc[i]=Cc[i]/Sqrt[norm];
];(*end For i*)

(*constructing F_{pq}*)
For[p=1,p≤4,p++,
For[q=1,q≤4,q++,
F[p,q]=h[p,q]+Sum[Cc[r]*Cc[s]*(2GJ[p,r,q,s]-GK[p,r,s,q]),{r,1,4},{s,1,4}];
];];

(*matrix*)
mF=Table[F[p,q],{p,1,4},{q,1,4}];
mS=Table[S[p,q],{p,1,4},{q,1,4}];
(*Eigensystem[{mF,mS}];*)(*solving the generalised eigenvalue problem*)
Enew=Min[Eigenvalues[{mF,mS}]];
For[i=1,i≤4,i++,
If[Enew==Eigenvalues[{mF,mS}][[i]],Cnew=Eigenvectors[{mF,mS}][[i]];];
];(*end For i*)
ΔE=Abs[Eold-Enew];

];(*end While*)

EG=2*Sum[Cc[p]*Cc[q]*h[p,q],{p,1,4},{q,1,4}]+Sum[Cc[p]*Cc[q]*Cc[r]*Cc[s]*(2GJ[p,r,q,s]-GK[p,r,s,q]),{p,1,4},{q,1,4},{r,1,4},{s,1,4}];
Eth=-2.85516038;(*given Eth from text book*)
Print["Resulting ground state energy for Helium atom is E_{G} = ", EG, " a.u."];
Print["in atomic energy unit of Hartree; % of deviation is ", 100*(EG-Eth)/Eth, "%; Time used = ", TimeUsed[], "s"];

Conclusion

In general, Hartree-Fock method predict in it’s own accord the ground state of many body system. In this project, we are dealing with single nucleus particle Helium, whereby the spatial coordinate of nucleus can be fixed at origin. Hence, the Hamiltonian is much simplified, leaving the electrons-electrons and electrons-nucleus (not electrons-nuclei) potential. Since electrons are fermions, they obey Pauli exclusion principle and any particle-exchange leads to a (-1) to the eigenvalue. In accord to exchange corellation, Helium have a permutation that does not change the ground state energy level.