Hartree-Fock for Helium excited state II

October 29, 2017

GoLuckyRyan Math , Leave a comment

This times, we will show the energy level diagram for helium atom. We already calculated the ground state, the 1s2s singlet and triplet excited states. We will calculate higher excited states, for example, 1s3s, 1s2p or 1s3p, etc, and included in the diagram.

The most difficult part is the evaluation of the mutual interaction matrix element

\displaystyle G_{ij}^{hk} = \langle b_i(x) b_h(y) | \frac{1}{r_{xy}} | b_j(x) b_k(y) \rangle

The angular integral was evaluated in the last post. And the radial integral is done using Mathematica.

I will use hydrogen 1s, 2s, 2p, 3s, 3p, 4s, and 4p states for basis. Thus, the diagram will contain some of the excited states from some possible combination. The 1s4s and 1s4p state will not be calculated, because there is no room for mixing with higher excited states but only for lower excited states. This will make the eigen state be unbound.

During the calculation, one of the tricky point is the identify of the n-th s-state. It is because the mixing is always there, and the mixing can be very large. In order to determine the principle quantum number, we have to check the wave function can see how many peaks in there.

Here is the energy level diagram obtained by Hartree-Fock method using limited hydrogen wave functions as basis set.

He_excite.PNG

 

Hartree-Fock method for Helium excited state

October 22, 2017

GoLuckyRyan Basic 6 Comments

There are two excited wavefunctions, singlet and triplet state. In the singlet state, the total spin is zero, while in the triplet state, the total spin is 1.

In any case, the Slater determinant is

\displaystyle \Psi(x,y) = \frac{1}{\sqrt{2}} \begin{vmatrix} \phi_1(x) & \phi_2(x) \\ \phi_1(y) & \phi_2(y) \end{vmatrix}

For spin-singlet state, we set,

\phi_1(q) = \uparrow \phi_{1s}(r),  \phi_2(q) = \downarrow \phi_{2s}(r)

For triplet state,

\phi_1(q) = \uparrow \phi_{1s}(r),  \phi_2(q) = \uparrow \phi_{2s}(r)

Since the Hamitonian is spin independent, we only take care of the spatial part.

For the triplet state, if we put in the wavefucntion in the Hamiltonian. Using the basis set approximation, The Fock matrix is

\displaystyle F_{ij} = H_{ij} + \sum_{hk}^M a_{\nu h}^*a_{\nu k} (G_{ij}^{hk} - \delta(\sigma_h, \sigma_k)G_{ik}^{hj})

The factor $latex \delta(\sigma_h, \sigma_k) $ is due to the spin component.

From the beginning, the Fock matrix is

\displaystyle \sum_{j}^M F_{ij} a_{\mu j} = \sum_{j}^M H_{ij} a_{\mu j} + \sum_{j}^M a_{\mu j} \sum_{hk}^M a_{\nu h}^*a_{\nu k} G_{ij}^{hk} - \sum_{j}^M a_{\nu j} \sum_{hk}^M a_{\nu h}^*a_{\mu k} G_{ij}^{hk}

We can see the exchange term is difference. We can interchange k \leftrightarrow j

\displaystyle \sum_{j}^M a_{\nu j} \sum_{hk}^M a_{\nu h}^*a_{\mu k} G_{ij}^{hk} = \sum_{k}^M a_{\nu k} \sum_{hj}^M a_{\nu h}^*a_{\mu j} G_{ik}^{hj}

and pull the a_{\mu j} out.

The calculation is straight forward. First, we calculate the 1s state due to the mean field created by the 2s orbital, then calculate the 2s state due to the mean field by the 1s orbital. And then repeat until converged.

Using the Hydrogen 1s, 2s, 3s, and 4s orbtial, we calculated

1s KE+PE = -2 Hartree = -54.42 eV
2s KE+PE = -0.4349 Hartree = -11.83 eV
The direct term = 0.2800 Hartree = 7.62 eV
The exchange term = 0.016 Hartree = 0.44 eV

The total energy = -2.171 Hartree = -59.07 eV.

Compare to the calculated ground state energy ( -77.11 eV), the excitation energy is 18.04 eV. And the ionization energy is 4.65 eV.

The experimental value for the triplet state excitation energy is 19.77 eV, and the ionization energy is 4.70 eV.

Consider the ground state energy is different from the experimental value by 2 eV. the calculation is fairly good.

 

For the singlet state, because the spin state of the two electron are opposite, the Fock matrix becomes Hartree matrix.

\displaystyle F_{ij} = H_{ij} + \sum_{hk}^M a_{\nu h}^*a_{\nu k} G_{ij}^{hk}

The calculated result is,

1s KE+PE = -2.000 Hartree = -54.42 eV
2s KE+PE = -0.400 Hartree = -10.89 eV
The direct term = 0.2531 Hartree = 6.89 eV
The exchange term = 0.040 Hartree = 1.09 eV

The total energy = -2.147 Hartree = -58.43 eV.

Compare to the calculated ground state energy ( -77.11 eV), the excitation energy is 18.68 eV. And the ionization energy is 4.00 eV.

The experimental value for the triplet state excitation energy is 20.55 eV, and the ionization energy is 3.92 eV.

It is interesting to compare with other formalism. Many other textbook use singlet wavefunction,

\displaystyle \Psi_S(x,y) = \frac{1}{\sqrt{2}}( \phi_1(x) \phi_2(y) + \phi_2(x) \phi_1(y)) \frac{1}{\sqrt{2}}(\uparrow \downarrow - \downarrow \uparrow) 

The triplet state is the same,

\displaystyle \Psi_T(x,y) = \frac{1}{\sqrt{2}}( \phi_1(x) \phi_2(y) - \phi_2(x) \phi_1(y)) \frac{1}{\sqrt{2}}(\uparrow \uparrow) 

Since the Hamiltonian is spin independent, we can see, when we calculate the total energy of the singlet state and triplet state are

E_S = E_0 + J + K

E_T = E_0 + J - K

Thus, the energy difference between singlet and triplet state is 2 exchange term.

In our calculation, the wave function of the singlet state is difference. We directly calculate from the Slater determinant. And the Fock matrix is lack of the exchange term due to the opposite spin. However, the triplet state is the same, the exchange term is 0.44 eV. If we do not do Hartree-Fock on the singlet state, but calculate the singlet state using the exchange term, the singlet state energy would be -58.19 eV. Comparing with the Hartree-Fock result (-58.43 eV), I would say they agree. The small difference is due to the difference of the mean fields.

Hartree method for Helium ground state

October 17, 2017

GoLuckyRyan atomic, Basic, Math , , , , , , 7 Comments

After long preparation, I am ready to do this problem.

The two electron in the helium ground state occupy same spacial orbital but difference spin. Thus, the total wavefunction is

\displaystyle \Psi(x,y) = \frac{1}{\sqrt{2}}(\uparrow \downarrow - \downarrow \uparrow) \psi(x) \psi(y)

Since the Coulomb potential is spin-independent, the Hartree-Fock method reduce to Hartree method. The Hartree operator is

F(x) = H(x) + \langle \psi(y)|G(x,y) |\psi(y) \rangle

where the single-particle Hamiltonian and mutual interaction are

\displaystyle H(x) = -\frac{\hbar^2}{2m} \nabla^2 - \frac{Ze^2}{4\pi\epsilon_0 x} = -\frac{1}{2}\nabla^2 - \frac{Z}{x}

\displaystyle G(x,y) = \frac{e^2}{4\pi\epsilon_0|x-y|} = \frac{1}{|x-y|}

In the last step, we use atomic unit, such that \hbar = 1, m=1, e^2 = 4\pi\epsilon_0. And the energy is in unit of Hartree, 1 \textrm{H} = 27.2114 \textrm{eV}.

We are going to use Hydrogen-like orbital as a basis set.

\displaystyle b_i(r) = R_{nl}(r)Y_{lm}(\Omega) \\= \sqrt{\frac{(n-l-1)!Z}{n^2(n+l)!}}e^{-\frac{Z}{n}r} \left( \frac{2Z}{n}r \right)^{l+1} L_{n-l-1}^{2l+1}\left( \frac{2Z}{n} r \right) \frac{1}{r} Y_{lm}(\Omega)

I like the left the 1/r, because in the integration \int b^2 r^2 dr, the r^2 can be cancelled. Also, the i = nlm is a compact index of the orbital.

Using basis set expansion, we need to calculate the matrix elements of

\displaystyle H_{ij}=\langle b_i(x) |H(x)|b_j(x)\rangle = -\delta \frac{Z^2}{2n^2}

\displaystyle G_{ij}^{hk} = \langle b_i(x) b_h(y) |G(x,y) |b_j(x) b_k(y) \rangle

Now, we will concentrate on evaluate the mutual interaction integral.

Using the well-known expansion,

\displaystyle G(x,y) = \frac{1}{|x-y|}=\frac{1}{r_{12}} = \sum_{l=0}^{\infty} \sum_{m=-l}^{l} \frac{4\pi}{2l+1} \frac{r_<^l}{r_>^{l+1}} Y_{lm}^{*}(\Omega_1)Y_{lm}(\Omega_2)

The angular integral

\displaystyle \langle Y_i(x) Y_h(y)| Y_{lm}^{*}(x) Y_{lm}(y) | Y_j(x) Y_k(y) \rangle \\ = \big( \int Y_i^{*}(x) Y_{lm}^{*}(x) Y_j(x) dx \big) \big( \int Y_h^{*}(y) Y_{lm}(y) Y_k(y) dy \big)

where the integral \int dx = \int_{0}^{\pi} \int_{0}^{2\pi} \sin(\theta_x) d\theta_x d\phi_x .

From this post, the triplet integral of spherical harmonic is easy to compute.

\displaystyle \int Y_h^{*}(y) Y_{lm}(y) Y_k(y) dy = \sqrt{\frac{(2l+1)(2l_k+1)}{4\pi (2l_h+1)}} C_{l0l_k0}^{l_h0} C_{lm l_km_k}^{l_hm_h}

The Clebsch-Gordon coefficient imposed a restriction on l,m.

The radial part,

\displaystyle \langle R_i(x) R_h(y)| \frac{r_<^l}{r_>^{l+1}} | R_j(x) R_k(y) \rangle \\ = \int_0^{\infty} \int_{0}^{\infty} R_i(x) R_h(y) \frac{r_<^l}{r_>^{l+1}} R_j(x) R_k(y) y^2 x^2 dy dx \\ = \int_0^{\infty} R_i(x) R_j(x) \\ \left( \int_{0}^{x} R_h(y) R_k(y) \frac{y^l}{x^{l+1}} y^2dy  + \int_{x}^{\infty} R_h(x)R_k(x) \frac{x^l}{y^{l+1}}  y^2 dy   \right) x^2 dx

The algebraic calculation of the integral is complicated, but after the restriction of l from the Clebsch-Gordon coefficient, only a few terms need to be calculated.

The general consideration is done. now, we use the first 2 even states as a basis set.

\displaystyle b_{1s}(r) = R_{10}(r)Y_{00}(\Omega) = 2Z^{3/2}e^{-Zr}Y_{00}(\Omega)

\displaystyle b_{2s}(r) = R_{20}(r)Y_{00}(\Omega) = \frac{1}{\sqrt{8}}Z^{3/2}(2-Zr)e^{-Zr/2}Y_{00}(\Omega)

These are both s-state orbital. Thus, the Clebsch-Gordon coefficient

\displaystyle C_{lm l_k m_k}^{l_h m_h} = C_{lm00}^{00}

The radial sum only has 1 term. And the mutual interaction becomes

\displaystyle G(x,y) = \frac{1}{|x-y|}=\frac{1}{r_{12}} = 4\pi \frac{1}{r_>} Y_{00}^{*}(\Omega_1)Y_{00}(\Omega_2)

The angular part

\displaystyle \langle Y_i(x) Y_h(y)| Y_{lm}^{*}(x) Y_{lm}(y) | Y_j(x) Y_k(y) \rangle = \frac{1}{4\pi}

Thus, the mutual interaction energy is

G_{ij}^{hk} = \displaystyle \langle b_i(x) b_h(y) |G(x,y) |b_j(x) b_k(y) \rangle = \langle R_i(x) R_h(y)| \frac{1}{r_>} | R_j(x) R_k(y) \rangle

The radial part

G_{ij}^{hk} = \displaystyle \langle R_i(x) R_h(y)| \frac{1}{r_>} | R_j(x) R_k(y) \rangle \\ \begin{pmatrix} G_{11}^{hk} & G_{12}^{hk} \\ G_{21}^{hk} & G_{22}^{hk} \end{pmatrix} = \begin{pmatrix} \begin{pmatrix} G_{11}^{11} & G_{11}^{12} \\ G_{11}^{21} & G_{11}^{22} \end{pmatrix} & \begin{pmatrix} G_{12}^{11} & G_{12}^{12} \\ G_{12}^{21} & G_{12}^{22} \end{pmatrix} \\ \begin{pmatrix} G_{21}^{11} & G_{21}^{12} \\ G_{21}^{21} & G_{21}^{22} \end{pmatrix} & \begin{pmatrix} G_{22}^{11} & G_{22}^{12} \\ G_{22}^{21} & G_{22}^{22} \end{pmatrix} \end{pmatrix} \\= \begin{pmatrix} \begin{pmatrix} 1.25 & 0.17871 \\ 0.17871 & 0.419753 \end{pmatrix} & \begin{pmatrix} 0.17871 & 0.0438957 \\ 0.0439857 & 0.0171633 \end{pmatrix} \\ \begin{pmatrix} 0.17871 & 0.0438957 \\ 0.0438957 & 0.0171633 \end{pmatrix} & \begin{pmatrix} 0.419753 & 0.0171633 \\ 0.0171633 & 0.300781 \end{pmatrix} \end{pmatrix}

We can easy to see that G_{ij}^{hk} = G_{ji}^{hk} = G_{ij}^{kh} = G_{hk}^{ij} = G_{ji}^{kh} . Thus, if we flatten the matrix of matrixes, it is Hermitian, or symmetric.

We also notice that G_{ij}^{hk} is a matrix of matrixes. The inner matrix runs for hk indexes, while the bigger matrix runs for ij index. This will simplify our minds.

Now, we can start doing the Hartree method.

The general solution of the wave function is

\psi(x) = a_1 b_{1s}(x) + a_2 b_{2s}(x)

The Hartree matrix is

F_{ij} = H_{ij} + \sum_{h,k} a_h a_k G_{ij}^{hk}

In this simple 2-states example, index h, k runs from 1 to 2. And remember the G_{ij}^{hk} is a matrix of matrixes, using the identify G_{ij}^{hk} = G_{hk}^{ij} the sum becomes:

\sum_{h,k} a_h a_k G_{ij}^{hk} = a_1^2 G_{11} + a1 a_2 (G_{12} + G_{21}) + a_2^2 G_{22},

where the G_{ij} represent the bigger matrix in G_{ij}^{hk}.

The first trial wave function is the Hydrogen-like orbital,

\psi^{(0)}(x) = b_{1s}(r)

F_{ij}^{(0)} = \begin{pmatrix} -2 & 0 \\ 0 & -0.5 \end{pmatrix}  + \begin{pmatrix} 1.25 & 0.17871 \\ 0.17817 & 0.419753 \end{pmatrix}

Solve for eigen system, we have the energy after 1st trial,

\epsilon^{(1)} = -0.794702 , (a_1^{(1)}, a_2^{(1)}) = (-0.970112, 0.242659)

After 13th trial,

\epsilon^{(13)} = -0.880049 , (a_1^{(13)}, a_2^{(13)}) = (-0.981015, 0.193931)

F_{ij}^{(13)} = \begin{pmatrix} -2 & 0 \\ 0 & -0.5 \end{pmatrix}  + \begin{pmatrix} 1.15078 & 0.155932 \\ 0.155932 & 0.408748 \end{pmatrix}

Thus, the mixing of the 2s state is only 3.7%.

Since the eigen energy contains the 1-body energy and 2-body energy. So, the total energy for 2 electrons is

E_2 = 2 * \epsilon^{(13)} - G = -2.82364 \textrm{H} = -76.835 \textrm{eV}

In which ,

G = \langle \psi(x) \psi(y) |G(x,y) |\psi(x) \psi(y) \rangle = 1.06354 \textrm{H} = 28.9403 \textrm{eV}

So the energies for

From He to He++.  E_2 = -2.82364 \textrm{H} = -76.835 \textrm{eV}
From He+ to He++, E_1^+ = -Z^2/2 = 2 \textrm{H} = -54.422 \textrm{eV} .
From He to He+, is E_1 = E_2 - E_1^+ = -0.823635 \textrm{H} =  -22.4123 \textrm{eV}

The experimental 1 electron ionization energy for Helium atom is

E_1(exp) = -0.90357 \textrm{H} = -24.587 \textrm{eV}
E_1^+(exp) = -1.99982 \textrm{H} = -54.418 \textrm{eV}
E_2(exp) = -2.90339 \textrm{H} = -79.005 \textrm{eV}

The difference with experimental value is 2.175 eV. The following plot shows the Coulomb potential, the screening due to the existence of the other electron, the resultant mean field, the energy, and r \psi(x)

Helium.PNG

Usually, the Hartree method will under estimate the energy, because it neglected the correlation, for example, pairing and spin dependence. In our calculation, the E_2 energy is under estimated.

From the (a_1^{(13)}, a_2^{(13)}) = (-0.981015, 0.193931) , we can see, the mutual interaction between 1s and 2s state is attractive. While the interaction between 1s-1s and 2s-2s states are repulsive. The repulsive can be easily understood. But I am not sure how to explain the attractive between 1s-2s state.

Since the mass correction and the fine structure correction is in order of 10^{-3} \textrm{eV} , so the missing 0.2 eV should be due to something else, for example, the incomplete basis set.

If the basis set only contain the 1s orbit, the mutual interaction is 1.25 Hartree = 34.014 eV. Thus, the mixing reduce the interaction by 5.07 eV, just for 3.7% mixing

I included the 3s state,

\epsilon^{(13)} = -0.888475 , (a_1^{(13)}, a_2^{(13)}, a_3^{(13)}) = (0.981096, -0.181995, -0.06579)

The mutual energy is further reduced to 1.05415 Hartree = 28.6848 eV. The E_2 = -77.038 \textrm{eV} . If 4s orbital included, the E_2 = -77.1058 \textrm{eV} . We can expect, if more orbital in included, the E_2 will approach to E_2(exp).

Hartree-Fock method for 1D infinite potential well

October 12, 2017

GoLuckyRyan Basic , Leave a comment

Following the previous post, I tested my understanding on the Hartree method. Now, I move to the Hartree-Fock method. The “mean field” energy of the Hartree-Fock is

\displaystyle G_{\alpha \beta}= \sum_{j=i+1}^{N} \langle \psi_{\alpha}(i) \psi_{\nu} (j) | G(i,j) \left(| \phi_{\beta}(i) \phi_{\nu}(j) \rangle  - |\phi_{\nu}(i) \phi_{\beta}(j) \rangle \right) \\ = \langle \alpha \nu | \beta \nu \rangle - \langle \alpha \nu | \nu \beta \rangle

I also use the same method, in which, the trial wave function is replaced every iteration and the integration is calculated when needed.

Since the total wave function must be anti-symmetry under permutation, therefore one state can be occupied by only one particle. Thus, if we use the ground state in the mean field, the “meaningful” wave functions are the other states.

It is interesting that the mean field energy is zero when \mu = \nu , the consequence is no mean field for the same state. Suppose the mean field energy is constructed using the ground state, and we only use 3 states, the direct term is

G_D = \begin{pmatrix} \langle 11|11 \rangle & \langle 11|21 \rangle & \langle 11|31 \rangle \\  \langle 21|11 \rangle & \langle 21|21 \rangle & \langle 21|31 \rangle \\ \langle 31|11 \rangle & \langle 31|21 \rangle & \langle 31|31 \rangle \end{pmatrix}

The exchange term is

G_E = \begin{pmatrix} \langle 11|11 \rangle & \langle 11|12 \rangle & \langle 11|13 \rangle \\  \langle 21|11 \rangle & \langle 21|12 \rangle & \langle 21|13 \rangle \\ \langle 31|11 \rangle & \langle 31|12 \rangle & \langle 31|13 \rangle \end{pmatrix}

Due to the symmetry of the mutual interaction. We can see that some off-diagonal terms are cancelled. For example,

\displaystyle \langle 1 1 | 3 1 \rangle = \int \psi_1^*(x) \psi_1^*(y) \cos(x-y) \psi_3(x) \psi_1(y) dy dx

\displaystyle \langle 1 1 | 1 3 \rangle = \int \psi_1^*(x) \psi_1^*(y) \cos(x-y) \psi_1(x) \psi_3(y) dy dx

These two integrals are the same. In fact,

$latex \langle \alpha \nu | \beta \nu \rangle – \langle \alpha \nu | \nu \beta \rangle = $

whenever \alpha = \nu

\displaystyle \langle \nu \nu | \beta \nu \rangle = \int \psi_\nu^*(x) \psi_\nu^*(y) \cos(x-y) \psi_\beta(x) \psi_\nu(y) dy dx

\displaystyle \langle \nu \nu | \nu \beta \rangle = \int \psi_\nu^*(x) \psi_\nu^*(y) \cos(x-y) \psi_\nu(x) \psi_\beta(y) dy dx

We can see, when interchange x \leftrightarrow y , the direct term and the exchange term are identical, and then the mean field energy is zero. Also, when \beta = \nu the mean field energy is also zero.

Due to the zero mean field, the off-diagonal terms of the Hamiltonian H_{\alpha \beta} with \alpha = \nu or \beta = \nu are zero. Then, the eigen energy is the same as the diagonal term and the eigen vector is un-change.

Back to the case, the direct matrix at the 1st trial is,

G_D = \begin{pmatrix} 0.720506 & 0 & -0.144101 \\  0 & 0.576405 & 0 \\ -0.144101 & 0 & 0.555819 \end{pmatrix}

The exchange matrix is

G_E = \begin{pmatrix} 0.720506 & 0 & -0.144101 \\  0 & 0.25 & 0 \\ -0.144101 & 0 & 0.0288202 \end{pmatrix}

Thus, the Fock matrix is

F = \begin{pmatrix} 1 & 0 & 0  \\  0 & 4.3264 & 0 \\ 0 & 0 & 9.527 \end{pmatrix}

Therefore, the eigen states are the basis, unchanged

\displaystyle \psi_1(x) = \sqrt{\frac{2}{\pi}} \sin(x)

\displaystyle \psi_2(x) = \sqrt{\frac{2}{\pi}} \sin(2x)

\displaystyle \psi_3(x) = \sqrt{\frac{2}{\pi}} \sin(3x)

Only the eigen energies are changed, as \epsilon_1 = 1 , \epsilon_2 = 4.3264 , \epsilon_3 = 9.527

The total wave function for 2 particles at state 1 and state-\mu is

\Psi(x,y) = \frac{1}{\sqrt{2}} ( \psi_1(x) \psi_\mu(y) - \psi_\mu(x) \psi_1(y) )

I found that the “mean field” function is not as trivial as in the Hartree case, because of the exchange term. In principle, the mean field for particle-i at state-\mu is,

G(i) = \int \phi_\nu^*(j) G(i,j) \phi_\nu(j) dj -  \int \phi_\nu^*(j) G(i,j) \phi_\mu(j) dj

However, the direct term is multiply with \psi_\mu(i) , but the exchange term is multiply with \psi_\nu(i) , which are two different functions. i.e., the “mean field” is affecting two functions, or the “mean field” is shared by two states.

Although the mean field for single state can be defined using exchange operator symbolically, but I don’t know how to really implement into a calculation. Thus, I don’t know how to cross-check the result.

Some comments on Hartree-Fock method

October 10, 2017

GoLuckyRyan Basic , Leave a comment

Most confusing thing for the method is that, weather the index is represent particle or state.

Look back the method. The total wave function is product of the wave function of each particle.

\Psi = \phi_1(1) \phi_2(2) .. \phi_N(N)

Here, we fixed our notation that the subscript is state, the parentheses is particle ID.

The single particle energy is,

\displaystyle \langle \Psi |\sum_{i}^{N} H(i)|\Psi \rangle = \sum_{i}^{N} \langle \phi_\mu(i) | H(i)|\phi_\mu(i) \rangle

Here, the i-particle in state \mu

The mutual interaction energy is

\displaystyle \langle \Psi |\sum_{i}^{N-1}\sum_{j=i+1}^{N} G(i,j)|\Psi \rangle = \sum_{i}^{N-1}\sum_{j=i+1}^{N} \langle \phi_\mu(i) \phi_\nu(j) | G(i,j)|\phi_\mu(i) \phi_\nu(j) \rangle

In both energy terms, we can factor out the sum of particle, we get the Hartree equation

\displaystyle \left( H(i) + \sum_{j=i+1}^{N} \langle \phi_\nu(j) | G(i,j)| \phi_\nu(j) \rangle  \right) |\phi_\mu(i) \rangle = \epsilon_\mu | \phi_\mu(i) \rangle

Now, we can see, the sum is sum over the number of particle, not state.

When solving the Hartree equation using a given basis b_k(x) , we set the first trial of the solution to be

| \phi_\mu(i) \rangle = | b_\mu(i) \rangle

Then, in the first trial, the Hartree matrix is

\displaystyle \hat{F}_{\alpha \beta} = \langle b_\alpha(i) | \left( H(i) + \sum_{j=i+1}^{N} \langle b_\nu(j) |G(i,j) | b_\nu(j) \rangle \right) | b_\beta(i) \rangle

In here, we keep the particle ID. Inside the mean field integral, the trial state b_\nu(j) is used for the particle-j.

After a trial, we will get the eigen vector for each state v_\mu, then we construct a new trial function and iterate until converge.

We can see, for multi-particle wave function, the sum is sum over the particle, even thought there are particles share the same state. In fact, since each particle should be at a state, the state should be a function of particle ID. So, the general Hartree matrix is

\displaystyle \hat{F}_{\alpha \beta} = \langle \psi_\alpha(i)|\left( H(i) + \sum_{j=i+1}^{N} \langle \psi_{\nu(j)}(j) | G(i,j)| \psi_{\nu(j)}(j) \rangle  \right) |\psi_\beta(i) \rangle

where \alpha, \beta = 1,2,3,..., k, k is number state to use.

In Hartree-Fock method, the sum is also over all other particle, not state.

 

Hartree method for 1D infinite potential well

October 7, 2017

GoLuckyRyan theory , , , Leave a comment

In order to test my understanding on the Hartree method. I played a simple case. The potential is

V(x) = \left\{ \begin{matrix} 0 , 0 \leq x \leq \pi  \\  \infty , else \end{matrix} \right\}

The mutual interaction between 2 particles is, for simple calculation,

G(x_1, x_2) = g_0 \cos(x_1 - x_2)

This mutual interaction is repulsive, we can expect the eigen wave function will be more spread out.

We use the eigen state of the V(x) potential to be the basis,

\displaystyle \phi_n(x) = \sqrt{\frac{2}{\pi}} \sin(nx)

The Hartree method is basically solving the Hartree equation

\displaystyle \left( -\frac{\hbar^2}{2m} \frac{d^2}{dx_1^2} + V(x_1) + \int G(x_1, x_2) (\psi_k(x_2))^2 dx_2 \right) \psi_n(x_1) = \epsilon_n \psi_n(x_1)

We assume the two particles are in the same ground state, thus, k = n = 1. We assume,

\psi_1(x) = \sum a_i \phi_i(x)

The algorithm is that, for the 1 particles matrix

\displaystyle F_{ij} = \int \psi_i(x) \left( -\frac{\hbar^2}{2m} \frac{d^2}{dx^2} + V(x)  \right) \psi_j(x) dx

Define

\displaystyle G_0(x) = \int G(x, x_2) (\psi_1(x_2))^2 dx_2

The 2-particle matrix

\displaystyle G_{ij} = \int \psi_i(x) G_0(x) \psi_j(x) dx

For the 1st trial, we set

\psi_i^{(0)}(x) = \phi_{(2i+1)}(x)

We only use “even” function with respect to x = \pi/2, because mixing parity will result non-reasonable result. ( why ? )

Then we get the matrix of the Hartree Hamiltonian H_{ij} = F_{ij} + G_{ij}, and solve for the eigen value and eigen vector,

\vec{v}^{(0)} = (v_1^{(0)}, v_2^{(0)}, .... v_n^{(0)})

The next trial is

\displaystyle \psi_i^{(1)}(x) = \sum_i v_i^{(0)} \psi_i^{(0)}

The process repeated until converge. The convergent can be exam using the Hartree Hamiltonian that the non-diagonal terms are vanished, i.e.

F_{ij} + G_{ij} = 0, i\neq j

I tried with 2 basis and 6 basis, I set \hbar^2/2m = 1 for simplicity, the Hartree ground state is 1.7125 for both basis. The eigen wavefunction for the 2-basis is

\psi(x) = 0.797752 \sin(x) + 0.0145564 \sin(3x)

with the mean field

V_m(x) = 0.842569 \sin(x)

For 6 basis, the wave function is

\psi(x) = 0.797751 \sin(x) + 0.0145729 \sin(3 x) + 0.000789143 \sin(5 x) + 0.000125871 \sin(7 x) + 0.0000336832 \sin(9 x) + 0.0000119719 \sin(11 x)

we can see, the contribution of the higher energy states are getting smaller and smaller. The higher energy states only contributes 0.03%. Following is the plot compare \phi_1(x) (orange) and \psi(x) (blue). We can see that the wave function becomes little spread out.

psi.PNG

If we increase the strength of the mutual interaction by 10 times, the wave function becomes more spread out.

psi3.PNG

To cross check the result,  we can use the mean field and solve the Schrodinger equation using variational method. The method gives the ground state energy be 1.7125 and converge.

vm.PNG

This is interesting that the variational method needs 4 or 5 basis to converge, while Hartree method only need 2.

At last, we substitute the wave functon back to the Hartree equation, it turns out the result is not satisfied 100%. The total wave function is

\Psi(x_1, x_2) = \psi(x_1) \psi(x_2)

The Schrodinger equation is

\displaystyle  H \Psi = -\frac{\hbar^2}{2m}\left(\frac{d^2}{dx_1^2} + \frac{d^2}{dx_2^2}\right) \Psi(x_1,x_2) + G(x_1,x_2)\Psi(x_1,x_2) = E \Psi(x_1,x_2)

However, the left side and the right side are not equal, but the integration

\langle \Psi |H|\Psi\rangle  = \langle |\Psi|\Psi \rangle

is equal.

I plot the \Psi H\Psi (orange) and $ E \Psi^2 $ (blue), and there difference is the following

diff.PNG

We can imagine, the integration of the bottom plot is zero. I don’t know the reason, but I guess it is not because of the number of basis, because the higher energy states only contributes for 0.03%. I guess, it is because of the energy is the quantity that is optimized rather then the wave function. And I am not sure that the reverse of the following is true.

H|\Psi\rangle = E |\Psi\rangle  \rightarrow \langle \Psi|H|\Psi \rangle = E \langle \Psi|\Psi \rangle

For the case of many particles, the G_0 has to sum over all other particle. When all particles are in the same states, the sum is simple. For particles distribute among different states, in general, the sum is

\displaystyle G_0 =\sum_{k\neq i} \int G(x_1,x_2) \psi_k(x_2)^2 dx

Here, the sum is excluded the particle x_1, which is in state i. For example, if the system consists of 2 particles on state-1, and 2 particles on state-2. The mean field for the particle in state-1 is constructed using 1 particles on state-1 and 2 particles on state-2. Thus the energy for the state-2 may be not correct.

The mean field for the particle on state-2, the mean field should be constructed by 2 particles on state-1 and 1 particle on state-2.

For higher states that without particle, the mean field is built by 2 particles on state-1 and 2 particles on state-2.

In general, the mean field is constructed using other particles.

This fixed mean field method is simple but we have to define the mean field for each state and calculate again and again.

Hartree-Fock and Mean-field

July 12, 2017

GoLuckyRyan Basic , , , , Leave a comment

For long time, i know the Hartree-Fock method is a way to find the mean-field, and the method is a mean-field theory. But how exactly the Hartree-Fock connects to the mean field, I have no idea. And occasionally, I mentioned the Hartree-Fock in its simplest form. Here, I will makes the connection crystal clear.

Hartree = self-consistence field

Fock = anti-symmetry wave function

Mean-field approximation, or in chemistry, self-consistence field approximation.

The full Hamiltonian is

\displaystyle H = H_1 + H_2

\displaystyle H_1 = \sum_{i}^{N} \left( -\frac{1}{2}\nabla^2_i - \frac{Z}{r_i} \right) = \sum_{i} H_i

\displaystyle H_2 =\sum_{i<j} \frac{1}{r_{ij}}= \sum_{i<j} G_{ij}

Here, H_1 is one-body operator, and H_2 is two-body operator. I used Coulomb potential in the one-body operator, but it can be generalized as H_i, and the mutual interaction can also be generalized as G_{ij}.

The idea of the mean-field approximation is that, what if, we can find a one-body potential V(r) , the mean-field, such that

\displaystyle H = H_0 + H_R

\displaystyle H_0 = \sum_{i}^N \left( -\frac{1}{2}\nabla^2_i + V(r_i) \right)

\displaystyle H_R = \sum_{i<j}\frac{1}{r_{ij}} - \sum_{i}\left( \frac{Z}{r_i} - V(r_i) \right)

Here, H_0 is the mean-field Hamiltonian, which represent most of the effective interaction to a particle, such that H_R is the residual interaction, which is very small, and can be later treated as perturbation.

So, the problem is, how to find this mean field V(r)?

Back to the form H = H_1 + H_2 , we first construct a trial wave function,

\displaystyle \Phi = \frac{1}{\sqrt{N!}} \sum_{P} (-1)^P P \Phi_H = \sqrt{N!} A \Phi_H

\displaystyle \Phi_H = \phi_\alpha(1)\phi_\beta(2) ... \phi_{\nu}(N)

\displaystyle A = \frac{1}{N!} \sum_{P}(-1)^P P = \frac{1}{N!} (1 - P_{ij} + P_{ijk} + .. )

Here, P is the permutation operator, it can be 1-body exchange, 2-body exchange, and so on, but we will see that, only 1-body exchange (which is no change at all) and 2-body exchange are needed. \Phi_H is the Hartree wave function, it is a simple product of wave function of difference particles of difference states, or the diagonal product of the Slater determinant. In $\phi_\lambda(i) $, $\latex \lambda$ represents the state and i is the “id” of the particle. Notice that

\langle \phi_\mu | \phi_\nu \rangle = \delta_{\mu \nu}

A is the anti-symmetrization operator, it commute with Hamiltonian and a kind of projector operator,

[A,H] = 0,  A^2 = A.

Now, we evaluate the energy using this trial wave function.

E_{\Phi} = \langle \Phi | H_1 | \Phi \rangle + \langle \Phi | H_2 | \Phi \rangle

\displaystyle \langle \Phi | H_1 | \Phi \rangle  = N! \langle \Phi_H | A H_1 A | \Phi_H \rangle = N! \langle \Phi_H | H_1 A | \Phi_H \rangle

\displaystyle N! \langle \Phi | H_1 | \Phi \rangle  = \sum_{i} \langle \Phi_H| H_i (1 - P_{ij} + P_{ijk} + .. ) | \phi_\alpha(1)\phi_\beta(2) ... \phi_{\nu}(N) \rangle

Since the one-body operator H_i only acts on the i particle, any exchange will make the operator do nothing on the j particle, then the orthogonality of the wave function makes the integration zero. Lets demonstrate on 2 particles case.

\displaystyle \langle \phi_a(1) \phi_b(2) | H_1 | \phi_a(2) \phi_b(1) \rangle = \langle \phi_a(1)| H_1 | \phi_b(1) \rangle \langle \phi_b(2) | \phi_a(2) \rangle = 0

Notice that this H_1 is the one-body operator act on particle 1.

OK,

\displaystyle \langle \Phi | H_1 | \Phi \rangle  = \sum_{i} \langle \phi_\mu(i) | H_i | \phi_\mu(i) \rangle = \sum_\mu e_\mu

\displaystyle \langle \Phi | H_2 | \Phi \rangle  =\frac{1}{2}\sum_{ij} \langle \Phi_H | H_2 (1 - P_ij) | \Phi_H \rangle

\displaystyle = \frac{1}{2} \sum_{ij} \begin{matrix} \langle \phi_\mu(i) \phi_\nu(j)|G_{ij}|\phi_\mu(i) \phi_\nu(j) \rangle \\ - \langle \phi_\mu(i) \phi_\nu(j)|G_{ij}|\phi_\nu(i) \phi_\mu(j) \rangle \end{matrix}  = \sum_{\mu<\nu} \left( \langle \mu\nu | \mu\nu \rangle - \langle \mu\nu | \nu\mu \rangle \right)

The \langle \mu\nu | \mu\nu \rangle is direct term, and the \langle \mu\nu | \nu\mu \rangle is the exchange term. This is a simplified notation, the first position is for the i-th particle, and the second position is always for the j-th particle.

Thus, the total energy is

\displaystyle E_{\Phi} = \sum_{\mu} e_{\mu} + \sum_{\mu<\nu} \left(\langle \mu\nu | \mu\nu \rangle - \langle \mu\nu | \nu\mu \rangle\right)

We can factor out the ket |\phi_\mu(i) \rangle in the above equation and get the Fock-operator, with a notation for the exchange term operator

\displaystyle F_i = H_i + \sum_{i<j} \left(J_{ij} - K_{ij}\right)

J_{ij} |\phi_\mu(i) \rangle = \langle \phi_\nu(j)|G_{ij} | \phi_\nu(j) \rangle |\phi_\mu(i) \rangle

K_{ij} |\phi_\mu(i) \rangle = \langle \phi_\nu(j)|G_{ij} | \phi_\mu(j) \rangle |\phi_\nu(i) \rangle

Becareful on the \mu , \nu ! I know the notation is messy, I know….

The Fock-operator is an effective one-body operator.

First, we put the trial basis wave function \phi_\lambda(i) can get the Fock-matrix \langle \phi_\nu(i) | F_i | \phi_\mu(i) \rangle , then diagonalize it, get the new eigen-states. Use this set of new eigen-states to calculate agian and again until converged!

So, where is the mean-field? Lets expand the Fock-operator into the Hartree-Fock equation.

F_i |\phi_\mu(i) \rangle = \epsilon_\mu |\phi_\mu(i) \rangle

\displaystyle \left( -\frac{1}{2} \nabla^2_i -  \frac{Z}{r_i} +  \frac{1}{2}\sum_{j} \left( J_{ij} - K_{ij} \right) \right) \phi_\mu(i)  = \left( -\frac{1}{2} \nabla^2_i + V(r_i)  \right) \phi_\mu(i) = E_\mu \phi_\mu(i)

\displaystyle V(r_i) = -  \frac{Z}{r_i} +  \sum_{i<j} \left( J_{ij} - K_{ij} \right)

This is the mean-field!

Using the trial basis, we evaluate the Fock-matrix, that is equivalent to evaluate the mean-field.

Capture.PNG

In this post, the Hartree-fock for 2-body ground state is discussed. Unfortunately, that method is not the same in here. I would said, that method is only Hartree but not Fock. Since the method in that post can find a consistence field, but the ground state spatial is not anti-symmetric.

Since the spin-state is factored out, the spatial wave function of the case is identical. Thus, the exchange term is gone. The Slater determinant is

\displaystyle \Phi = \frac{1}{\sqrt{2}}\begin{pmatrix} \phi(1)\alpha & \phi(1) \beta \\ \phi(2) \alpha & \phi(2) \beta \end{pmatrix} = \phi(1) \phi(2) \frac{1}{\sqrt{2}}(\alpha(1) \beta(2) - \beta(1) \alpha(2))

Here, \alpha, \beta are spin-state.

In general, the 2 particle system, the energy is

\langle \Phi | H_1 |\phi \rangle = \langle \phi_a(1) | H_1|\phi_a(1) \rangle + \langle \phi_b(2) | H_1|\phi_b(2) \rangle

\langle \Phi | H_2 |\phi \rangle = \langle \phi_a(1) \phi_b(2) | G_{12}|\phi_a(1)  \phi_b(2) \rangle - \langle \phi_a(1) \phi_b(2) | G_{12}|\phi_b(1) \phi_a(2) \rangle

When a = b , the direct term and exchange term cancelled. Thus, the “mean-field” is simply the Coulomb potential. Therefore, the method in that post is kind of getting around.

Some people may found that the Fock operator is

\displaystyle F_i = H_i + \frac{1}{2} \sum_{j} \left(2 J_{ij} - K_{ij} \right)

In this way, the direct term for i = j will not be cancelled.  In the case of 2 particle, the mean field is the Coulomb potential plus the average mutual interaction from the other particle.

 

Hartree-Fock for 2-body ground state

June 8, 2017

GoLuckyRyan Basic , Leave a comment

Long time ago, we talked about the mean-field calculation, an touched Hartree-Fock method. In that time, we explained excited-state approach. Now, we explain another approach by variation of the wave functions. This approach is inevitable in atomic physics, because the potential is fixed.

The Hamiltonian is

H = H_1 + H_2 + V_{12}

Since the spin component is anti-parallel, the space part of the total wave function is

\Psi = \phi_1(r_1) \phi_2(r_2)

The Schrodinger equation is

H\Psi = E\Psi

Integrate both side with \phi_1(r_1)

\int dr_1 \phi_1(r_1) H \Psi = \int dr_1 \phi_1(r_1) E \Psi

using the normalization and define \left<x| 1\right> = \phi_1(r_1)

\left<1|H|1\right> \phi_2(r_2) = E \phi_2(r_2)

similarly

\left<2|H|2\right> \phi_1(r_1) = E \phi_1 (r_1)

expand H = H_1 + H_2 + V_{12}

(H_2 + \left<1|V_{12}|1\right>) \phi_2(r_2) = (E - \left<1|H_1|1\right>) \phi_2(r_2)

(H_1 + \left<2|V_{12}|2\right>) \phi_1(r_1) = (E - \left<2|H_2|2\right>) \phi_1(r_1)

Now, we have 2 equations,  an initial guess of \Psi = \phi_1(r_1) \phi_2(r_2) ,

The difficulty is that, the \left<1|V_{12}|1\right> contains 2 variables.

Capture.PNG

From Mean field calculation to Independent particle model

February 28, 2016

GoLuckyRyan Basic, Math, nuclear, theory , , Leave a comment

The independent particle model (IPM) plays a fundamental role in nuclear model. The total Hamiltonian of a nucleus is:

H = \sum{\frac{P_i^2}{2m_i}} + \sum_{i,j}{V_{ij}}

The mean field is constructed by Hartree-Fork method, so that, the total Hamiltonian,

H = \sum{h_i} + \sum_{i,j}{V_{ij} - U_i}, h_i = \frac{P_i^2}{2m_i} + U_i

The second term is called residual interaction. The residual should be as small as possible. Under the mean filed, each nucleon can be treated as independent particle model that experience by h_i, so that,

h_i \left|\phi_i\right> = \left|\phi_i\right> \epsilon_i,

where \epsilon_i is the single particle energy, and \left|\phi_i\right> is the single particle wave function.

Hartree-Fock method

Start with a trial single particle function \phi_i, construct a first trial total wavefuction \Psi_0

\Psi_0=\frac{1}{\sqrt{A!}} \left| \begin{array}{ccc} \phi_1(r_1) & \phi_1(r_2) & ... \\ \phi_2(r_1) & \phi_2(r_2) & ... \\ ... & ... & \phi_A(r_A) \end{array} \right|

where A is the number of particle. Using variation method,

\delta \left<\Psi_0|H|\Psi_o\right> = 0 \Rightarrow \left<\delta \Psi_0|H|\Psi_0\right>

H = \sum\limits_{i} h_i + \sum\limits_{i\neq j} V_{ij}

we don’t need to variate the ket and bar, since they are related.

The variation can be made on

  1. \phi_i or
  2. the particle-hole excitation.

For the particle-hole excitation,

\left|\delta \Psi_0\right> = \sum \eta_{kt} \left|\Psi_{kt}\right>

\left|\Psi_{kt}\right>=\left|1,2,...,t-1, t+1,...,A,k\right>

Then, the variation of the energy becomes,

\sum \eta_{kt} \left< \Psi_{kt}|H| \Psi_0 \right> =0

\left<\Psi_{kt}|H|\Psi_0\right> = \left<k|\sum h |t\right>+ \sum\limits_{r} \left< kr | V_{ij}| rt \right>=0

we used the one-body expression for the two-body interaction,

\sum \limits_{ij} V_{ij} = \sum \limits_{\alpha \beta \gamma \delta} \left|\alpha \beta\right> V_{\alpha \beta \gamma \delta} \left<\gamma \delta\right|

V_{\alpha \beta \gamma \delta} = \left<\alpha \beta |V_ij |\gamma \delta\right>

Thus, we take out \left<k\right| and \left|t\right>, we get the Hartree-Fock single particle Hamiltonian,

h_{HF} = h + \sum \limits_{r} \left< r | V_{ij} | r \right>

using this new single particle Hamiltonian, we have a better single particle wavefunction

h_{HF} \phi_i^1 = \epsilon_i^1 \phi_1^1

U = \sum \limits_{i} h_{HF}

with this new wavefunction, the process start again until convergence. After that, we will get a consistence mean field (in the sense that the wavefunction and the potential are consistence), the single particle energy and the total binding energy.

I am not sure how and why this process can minimized the mean field U