\documentclass{report}
\usepackage[utf8x]{inputenc}
\usepackage[a4paper]{geometry}

\usepackage{hyperref}

\usepackage{amsmath,amsthm}

\usepackage{tikz}
\usepackage{verbatim}

\usepackage[lf]{MinionPro}
\usepackage[scr=rsfso,calscaled=.96]{mathalfa}


\newcommand{\One}{\hat{\mathbf{1}}} \newcommand{\eff}{\text{eff}}
\newcommand{\Heff}{\hat{H}_\text{eff}}
\newcommand{\Veff}{\hat{V}_\text{eff}}
\newcommand{\braket}[1]{\langle#1\rangle}
\newcommand{\Span}{\operatorname{sp}}
\newcommand{\tr}{\operatorname{trace}}
\newcommand{\diag}{\operatorname{diag}}
\newcommand{\bra}[1]{\left\langle #1 \right|}
\newcommand{\ket}[1]{\left| #1 \right\rangle} \newcommand{\element}[3]
{\bra{#1}#2\ket{#3}}

\usepackage{amsmath}

\begin{document}


\section*{First project fall 2015 FYS-KJM4480/9480}

{\bf Date given: Thursday October 8, 2015.\\
Date revised: Monday October 19, 2015.
Deadline: Monday November 2 at 12pm (noon). It counts 30\% of the final mark.}

\subsection*{Introduction and background}


This project contains 3 exercises. 

A \emph{quantum dot} is a system of electrons
confined in semiconducting heterostructures. Such systems exhibit, due
to their small size, discrete quantum levels. Small confined systems
such as quantum dots have become very popular for experimental study.
Beyond their possible relevance for nanotechnology, they are highly
tunable in experiments and introduce level quantization and quantum
interference in a controlled way.  The possibility to manufacture
systems with a tailored electronic structure may improve electrical or
optical properties of materials, a reason why quantum dots are good
candidates as components in such diverse applications as quantum
computers, optimized solar cells, laser technology and medical
imaging, to name a few. On the other hand, strongly confined electrons
offer a wide variety of complex and subtle phenomena which pose severe
challenges to existing many-body methods.

We study a model called \emph{the parabolic quantum dot}: A system of
electrons confined in two dimensions by a harmonic oscillator
potential.  Our Hamiltonian for $N$ electrons takes the following form
\begin{equation}
  \hat{H} = \sum_{i=1}^{N} \left(\frac{-\hbar^2}{2m_e}\nabla_i^2 \right)
  + \sum_{i=1}^{N} \frac{1}{2}m_e\omega^2r^2_i + \sum_{i<j}^{N}
  \frac{ke^2}{r_{ij}} = \hat{H}_0 + \hat{W},
  \label{eq:ham}
\end{equation}
with
\begin{equation}
  \hat{H}_0 = \sum_{i=1}^N \hat{h}(\vec{r}_i), \quad \hat{h}(\vec{r})
  = -\frac{\hbar^2}{2m_e} \nabla^2 + \frac{1}{2}m_e \omega^2 r^2.
\end{equation}
Here, $r\equiv |\vec{r}|$, $r_i \equiv |\vec{r}_i|$, and $r_{ij} \equiv |\vec{r}_i -
\vec{r}_j|$. 

We will consider a system of $N=2$ and $N=6$ electrons. These
represent so-called closed shell systems for the paraboic quantum
dot. The next ``magic numbers'' for the parabolic quantum dot is
$N=12$, $N=20$, $N=30$, $N=42$, and these electron numbers are
analogous to noble gas electron numbers in atoms and closed-shell
structures in atomic nuclei.

As single-particle basis functions for our calculations we will use
spin-orbitals made from the harmonic oscillator eigenfunctions
$\psi_{nm}(r,\theta)$ in polar coordinates,
\begin{equation}
  \hat{h}(\vec{r})\psi_{nm}(r,\theta) = e_{nm} \psi_{nm}(r,\theta)
\end{equation}
with
\begin{equation}
  e_{nm} = \frac{1}{2}\hbar\omega (2n + |m| + 1).
\end{equation}
and (with $a = \sqrt{m_e\omega/\hbar}$), 
\begin{equation}
  \psi_{nm}(r,\theta)=a e^{im\theta}\sqrt{\frac{n!}{\pi(n+|m|)!}}(a
  r)^{|m|}L_{n}^{|m|}(a^2 r^2))\exp{(-a^2 r^2/2)}. 
\end{equation}
Here, $L_n^{|m|}$ are the so-called associated Laguerre polynomials, the
first few of which are
\begin{align}
  L_0^{|m|}(x) &= 1, \\
  L_1^{|m|}(x) & = -x + 1 + |m|.
\end{align}
The quantum number $n$ is the number of nodes in the radial part,
while $\hbar m$ is the orbital angular momentum, the eigenvalue of the
operator
\begin{equation}
  \hat{l}_z(\vec{r}) = -i\hbar \frac{\partial}{\partial \theta}.
\end{equation}
You will not need to use all the details of the form of $\psi_{nm}$ in this
project.

The spin-orbitals become
\begin{equation}
  \phi_p(x) = \psi_{nm}(r,\theta)\chi_\alpha(\sigma), \quad p =
  p(nm\alpha),
\end{equation}
where $\alpha = \pm 1/2$ is the quantum number for the $z$-projection
of the electron spin, and $\chi_\alpha$ is the corresponding spinor basis
function, $\braket{\chi_\alpha|\chi_\beta} =
\delta_{\alpha\beta}$. Thus $p$ counts the spin-orbitals, and can be
considered as in invertible function of the quantum numbers $nm\alpha$.

(We sometimes write $n_p$, $m_p$, and $\alpha_p$ for the quantum numbers
belonging to a given spin-orbital index $p$.)

The spin-projection operator can be written in terms of its action on
the spin-orbital basis,
\begin{equation}
  \hat{s}_z \ket{\phi_p} = \hbar\alpha_p \ket{\phi_p}.
\end{equation}


The orbital angular momentum operator is, on first-quantized form,
\begin{equation}
  \hat{L}_z = \sum_{i=1}^N \hat{l}_z(i),
\end{equation}
while the spin-projection operator is
\begin{equation}
  \hat{S}_z = \sum_{i=1}^N \hat{s}_z(i).
\end{equation}

Note that each eigenvalue $e_{nm}$
of the harmonic-oscillator $\hat{h}$ is now doubly degenerate due to spin.
The one-particle operator $\hat{h}$ is now diagonal in the
spin-orbital basis. Writing $p = p(\alpha n m)$ and $p' = p'(\alpha'
n' m')$, we get
\begin{equation}
  \braket{\phi_{p}|\hat{h}(\vec{r})|\phi_{p'}} =
  \braket{\chi_{\alpha}|\chi_{\alpha'}} (\psi_{n
    m}|\hat{h}|\psi_{n' m'}) = \delta_{\alpha,\alpha'}\delta_{n,n'} \delta_{m,m'}
  e_{n m} = \delta_{p,p'} e_{n m}. \label{eq:h-diagonal}
\end{equation}
Note that
\begin{equation}
  \delta_{p,p'} = \delta_{\alpha,\alpha'}\delta_{n,n'}\delta_{m,m'}
\end{equation}
since $p=p'$ if and only if $\alpha=\alpha'$, $n=n'$, and $m=m'$.


The matrix elements of the Coulomb interaction $\hat{W}$ in the
spin-orbital basis can be written in terms of the integrals
\begin{equation}
  \begin{split}
(\psi_{n_1m_1}\psi_{n_2m_2}|&\hat{w}|\psi_{n_3m_3}\psi_{n_4m_4}) = \\
&\iint r_1dr_1 d\theta_1 r_2dr_2
d\theta_2\psi_{n_{1}m_{1}}^*(r_1,\theta_1)\psi_{n_{2}m_{2}}^*(r_2,\theta_2) 
  \frac{ke^2}{r_{12}}\psi_{n_{3}m_{3}}(r_1,\theta_1)\psi_{n_{4}m_{4}}(r_2,\theta_2).
\end{split}
\label{eq:w-integrals}
\end{equation}
The explicit values for these matrix elements for $ke^2=1$ are found in the
appendix of this project.  Note well that these matrix elements
do not include spin. When setting up, say, the matrix elements
$\braket{\phi_{p_1}\phi_{p_2}|\hat{w}|\phi_{p_3}\phi_{p_4}}$,
you need to consider the spin degrees of freedom as well.


Our computational Hilbert space will be defined by the three lowest harmonic
oscillator single-particle energy levels. In total, when we add spin, we
have $L=12$ single-particle states. Thus, we have 6 orbitals and 12
spin-orbitals. These are illustrated in Fig.~\ref{fig:fig1}.


\subsection*{Exercise 1 (10 points)}

\begin{enumerate}
\item[1a)]
  Let $c^\dag_p$ and
  $c_p$ be creation and annihilation operators for the spin-orbitals $\phi_p$.

  Define the second-quantized forms of $\hat{L}_z$ and
  $\hat{S}_z$. Using that the spin-orbitals are eigenfunctions,
  simplify the expressions as much as possible.
\item[1b)]
  Show that any Slater determinant $\ket{p_1\cdots p_N}$ (built from
  the spin-orbitals $\phi_p$) is an eigenfunction for $\hat{L}_z$
  with eigenvalue $\hbar M = \hbar \sum_i m_i$, where $m_i =
  m_{p_i}$. Show that the Slater determinant is also an
  eigenfunction for $\hat{S}_z$ with eigenvalue $\hbar S_z = \hbar\sum_i
  \alpha_i$, where $\alpha_i = \alpha_{p_i}$.
\item[1c)]
  Show that
  \begin{equation}
    [\hat{S}_z, \hat{L}_z] = 0,
  \end{equation}
  by computing the action of each commutator on an arbitrary Slater
  determinant.
\end{enumerate}


From 1c) and the fact that $[\hat{H},\hat{L}_z]= [\hat{H},\hat{S}_z] = 0$ (you do not have to
show this) we obtain that the ground-state wavefunction of $\hat{H}$ can
be chosen to be an eigenfunction of both $\hat{S}_z$ and
$\hat{L}_z$. We will in the rest of the exercise consider only the
case $M = 0$ and $S_z = 0$. From 1c) we then obtain that the
ground-state wavefunction can be constructed as a linear combination
of those Slater determinants that have $M = S_z = 0$ only: all the
$m_i$'s sum to zero and all the $\alpha_i$'s sum to zero for these determinants. All other
Slater determinants are eliminated \emph{a priori}.

\subsection*{Exercise 2 (20 points)}

We can transform the Hamiltonian to dimensionless form by introducing
a suitable set of units. We will compute in a set of units where
\begin{equation}
  \hbar = m_e = ke^2 = 1.
\end{equation}
leaving $\omega$ as the last independent
parameter, controlling the strength of the confinement potential. In
the rest of the project, we use these units for simplicity.

We also limit ourselves to the case
\begin{equation}
  \omega = 1.
\end{equation}

\begin{enumerate}
\item[2a)] We start with the two-electron system and define our
  single-particle Hilbert space to consist of the orbitals $(n=0,
  m=0)$, $(n=1, m=0)$, $(n=0, m=\pm 1)$ and $(n=0,m=\pm 2)$, with
  their corresponding spin degeneracies, see Fig.~\ref{fig:fig1} for a
  schematic picture. Without spin we have thus six single-particle
  states.  Adding the two spin degeneracies, we end up with 12
  states. In the figure, the indices $P = P(n,m,\alpha)$ are
  listed below each spin-orbital.
\begin{figure}
\begin{center}
\begin{tikzpicture}
  \draw[->] (-3.5,-1) -- (3.5,-1) node[anchor=west] {$m$};
  \draw[->] (-3.5,-1) -- (-3.5,3) node[anchor=south] {$e$};

  \draw (-3.6,0) node[anchor=east] {$1\omega$} -- (-3.4,0);
  \draw (-3.6,1) node[anchor=east] {$2\omega$} -- (-3.4,1);
  \draw (-3.6,2) node[anchor=east] {$3\omega$} -- (-3.4,2);

  \draw (-2,-1.2) node[anchor=north] {$-2$} -- (-2,-.8);
  \draw (-1,-1.2) node[anchor=north] {$-1$} -- (-1,-.8);
  \draw (0,-1.2) node[anchor=north] {$0$} -- (0,-.8);
  \draw (0,-1.2) node[anchor=north] {$0$} -- (0,-.8);
  \draw (1,-1.2) node[anchor=north] {$+1$} -- (1,-.8);
  \draw (2,-1.2) node[anchor=north] {$+2$} -- (2,-.8);

  \begin{scope}
    \draw (-.5,0) -- (.5,0);
    \draw[->] (-.25,-.25) node[anchor=north]{$1$} -- (-.25,.25);
    \draw[<-] (.25,-.25) node[anchor=north]{$2$} -- (.25,.25);
  \end{scope}

  \begin{scope}[xshift=-1cm,yshift=1cm]
    \draw (-.5,0) -- (.5,0);
    \draw[->] (-.25,-.25) node[anchor=north]{$3$} -- (-.25,.25);
    \draw[<-] (.25,-.25) node[anchor=north]{$4$} -- (.25,.25);
  \end{scope}


  \begin{scope}[xshift=1cm,yshift=1cm]
    \draw (-.5,0) -- (.5,0);
    \draw[->] (-.25,-.25) node[anchor=north]{$5$} -- (-.25,.25);
    \draw[<-] (.25,-.25) node[anchor=north]{$6$} -- (.25,.25);
  \end{scope}


  \begin{scope}[xshift=-2cm,yshift=2cm]
    \draw (-.5,0) -- (.5,0);
    \draw[->] (-.25,-.25) node[anchor=north]{$7$} -- (-.25,.25);
    \draw[<-] (.25,-.25) node[anchor=north]{$8$} -- (.25,.25);
  \end{scope}


  \begin{scope}[xshift=0cm,yshift=2cm]
    \draw (-.5,0) -- (.5,0);
    \draw[->] (-.25,-.25) node[anchor=north]{$9$} -- (-.25,.25);
    \draw[<-] (.25,-.25) node[anchor=north]{$10$} -- (.25,.25);
  \end{scope}


  \begin{scope}[xshift=2cm,yshift=2cm]
    \draw (-.5,0) -- (.5,0);
    \draw[->] (-.25,-.25) node[anchor=north]{$11$} -- (-.25,.25);
    \draw[<-] (.25,-.25) node[anchor=north]{$12$} -- (.25,.25);
  \end{scope}

\end{tikzpicture}
\caption{Spin-orbitals for an electron in a two-dimensional
  oscillator well using a so-called Fock--Darwin representation.  Here
  we show the states and their oscillator energies for three
  harmonic-oscillator levels. The oscillator energies are given by
  $e_{nm}=\hbar\omega(2n+|m|+1)$. }
\label{fig:fig1}
\end{center}
\end{figure}

Set up the $M=S_z = 0$ ground-state of the noninteracting problem
($\hat{H}_0$) as a reference Slater determinant $\ket{\Phi}$. Use second
quantization. Define the Fermi level/energy. Draw a diagram of the
reference state, in the same fashion as Fig.~\ref{fig:fig1}. (If you
use \LaTeX, such diagrams can be easily made, see the source code for
this project.)

Define quasiparticle creation- and annihilation operators for this problem.

Construct thereafter all possible one-particle-one-hole excitations
$\ket{\Phi_i^a}$, where $i$ is an index below the Fermi level, and $a$
is an index above the Fermi level. Write the Slater determinants in
terms of quasiparticle operators. Ensure that you only write down
those states that have $M = m_1 + m_2 = 0$ and $S_z = \alpha_1 +
\alpha_2 = 0$. Draw diagrams.

Construct thereafter all possible two-particle-two-hole excitations
$\ket{\Phi_{ij}^{ab}}$ in a similar manner. Draw diagrams.


\item[2b)] Write down the Hamiltonian $\hat{H}$ on second-quantized
  form. Use
  Eqns.~\eqref{eq:h-diagonal} and \eqref{eq:w-integrals}. Include only
  nonzero terms in the sums. (Do not separate the sums over
  spin-orbitals into particles and holes; keep the sums over general
  spin-orbitals.)

  Compute the expectation value
  \begin{equation}
    \braket{\Phi|\hat{H}|\Phi}
  \end{equation}
  and obtain a numerical value by using the integral table in the
  Appendix. Explain why this expectation value can be considered an
  approximation to the exact ground-state energy.



\item[2c)] We now limit ourselves to the configuration-interaction
  with singles (CIS) basis, consisting of $\ket{\Phi}$ and the
  one-particle-one-hole functions $\ket{\Phi_i^a}$ you found in 2a).

  Compute the matrix elements (without inserting the numerical
  values of the $\hat{w}$ matrix elements) of  
  \[
  \langle \Phi | \hat{H}| \Phi_j^b \rangle, \quad \braket{\Phi_i^a|\hat{H}|\Phi},
  \] 
  and 
  \[
  \langle \Phi_i^a | \hat{H}| \Phi_j^b \rangle.
  \]
  Use Wick's Theorem for quasiparticles.

  % Represent these expectation values in a diagrammatic form, both 
%for the onebody part and the two-body part of the Hamiltonian. 
 
Insert then the explicit values for the various matrix elements and 
set up the final CIS Hamiltonian matrix and diagonalize it using your
favorite tool, for example Octave, Matlab, Python, C++ or Fortran.

Compare your results from those of exercise 2b) and comment your results. 
The exact ground-state energy with our Hamiltonian is $E_0=3$ for the
two-electron quantum dot. ``Exact'' means in the limit of a complete infinite-dimensional
harmonic-oscillator basis.


\item[2d)] We repeat exercises b) and c) but now for the six electron quantum dot.
Define the reference wavefunction $|\Phi\rangle$ and limit yourself again to
one-particle-one-hole excitations.    

Compute the reference energy $\langle \Phi | \hat{H}| \Phi \rangle $,
inserting the numerical values for the $\hat{w}$ matrix
elements. Thereafter you will need to set up the appropriate
Hamiltonian matrix which involves also one-particle-one-hole
excitations, also using numerical values for matrix
elements. Diagonalize this matrix and compare your eigenvalues with
$\langle \Phi | \hat{H}| \Phi \rangle$ and comment the result.  The
best numerical result for our Hamiltonian is obtained with Diffusion
Monte Carlo calculations, and for $\omega=1$ the energy for the ground
state is $E_0 = 20.1597$ atomic units for the six electron quantum
dot, see for example Lohne Pedersen {\em et al}, Physical Review B
{\bf 84}, 115302 (2011). 
\end{enumerate}


\subsection*{Exercise 3 (35 points)}

Having treated the parabolic quantum dot with
configuration-interaction, we now turn to a restricted Hartree--Fock (RHF)
treatment. 

Recall that we have 12 spin-orbitals, but only 6 orbitals. Define an index
$P = P(n,m)$ to count the orbitals without spin. Thus,
$P=1,2,\cdots,6$. 
\[ P(0,0) = 1, \quad P(0,-1) = 2, \quad P(0,1) = 3, \quad P(0,-2) = 4,
\quad P(1,0) = 5, \quad P(0,4) = 6. \]
This is the numbering used in the Appendix. Note well, that this is
for the \emph{spatial orbitals} 
only, hence the capital ``$P$'' instead of the lowercase ``$p$''.


\begin{enumerate}
\item[3a)] Define the RHF wavefunction $\ket{\Phi_\text{RHF}}$ for
  $N=2$ and $N=6$. Use
  the notation $\tilde{\psi}_I(\vec{r})$, $I=1,2,\cdots,N/2$, for the
  unknown spatial orbitals to be found.

  Compute, in terms of the spatial matrix elements of $\hat{h}$ and $\hat{w}$,
  \begin{equation}
    E_\text{RHF} = \braket{\Phi_\text{RHF}|\hat{H}|\Phi_\text{RHF}},
  \end{equation}
  
  In what manner does the RHF wavefunction
  ansatz and energy expression differ from the general Hartree--Fock wavefunction for fermions?
\item[3b)] In the remainder, we will expand the RHF spatial
  orbitals in the harmonic oscillator functions $\psi_{P}(\vec{r})$,
  $P=1,\cdots,6$. Define therefore an $6\times 6$ matrix $U$ and write
  \begin{equation}
    \tilde{\psi}_Q(\vec{r}) = \sum_{P=1}^6 \psi_P(\vec{r}) U_{PQ}.\label{eq:tildepsi}
  \end{equation}
  Explain why the matrix $U$ must be a unitary matrix. What part of
  $U$ corresponds to the occupied RHF orbitals? What part corresponds
  to the virtual RHF orbitals?
\item[3c)]
  
  Using your result in 3a) and the definition~\eqref{eq:tildepsi},
  show that the RHF energy of can be written
  \begin{equation}
    E_\text{RHF} =
    \braket{\Phi_\text{RHF}|\hat{H}|\Phi_{\text{RHF}}} = 2 \sum_{PQI}
    U_{QI}^* (\psi_Q|\hat{h}|\psi_P) U_{PI} + \sum_{PQI}
    U_{QI}^* \left(\sum_{RS} D_{SR}
     [QR|PS]\right) U_{PI} .
    \label{eq:rhf-energy}
  \end{equation}
  The sums over $I$ and $J$ go from 1 to $N/2$, while the sums over
  $P$, $Q$, $R$, and $S$ go from $1$ to $L/2$, and we have defined
  \begin{equation}
    [QR|PS] \equiv (\psi_Q\psi_R|\hat{w}|\psi_P\psi_S) -\frac{1}{2}
    (\psi_Q\psi_R|\hat{w}|\psi_S\psi_P).
  \end{equation}
  We have also defined the reduced density matrix,
  \begin{equation}
    D_{SR} = 2 \sum_I U_{SJ} U_{RJ}^*.
  \end{equation}
  
  Split the $\hat{W}$-part of $E_\text{RHF}$ into a direct and an
  exchange energy.  Explain their origin, and compare with standard HF
  theory formulated in an orthonormal spin-orbital basis.
\item[3d)]
  The RHF equations in the given basis can be written
  \begin{equation}
    F(D) U = U\epsilon, \qquad \text{(Roothan--Hall equation)}
  \end{equation}
  (you do not have to show this),  where $F(D)$ is the Fock matrix, 
  \begin{equation}
    F_{QP}(D) = (\psi_Q|\hat{h}|\psi_P) + \sum_{RS} D_{SR} [QR|PS],
  \end{equation}
  and where $\epsilon =
  \operatorname{diag}(\epsilon_1,\cdots,\epsilon_{L/2})$ is a diagonal
  matrix of eigenvalues.
  
  Explain in what way the Roothan--Hall equation differs from an ordinary
  matrix eigenvalue problem.

\item[3e)]
  Show that the RHF energy can be written
  \begin{equation}
    E_\text{RHF} = 2\sum_I \epsilon_I - \sum_{IPQ} U_{QI}^* \left(\sum_{SR}
    D_{SR} [QR|PS] \right) U_{PI} .
  \end{equation}
  Hint: Multiply The Roothan--Hall equation with $U^H$ from the left,
  and compere with Eq.~\eqref{eq:rhf-energy}.
\item[3f)]
  We are now going to solve the Roothan--Hall equation numerically.
  Use your favorite programming language and environment. Use comments
  to document your code and make it readable.
  
  You will need to read the integrals
  $(\psi_P\psi_Q|\hat{w}|\psi_R\psi_S)$ from a text file
  \texttt{coulomb.dat} downloadable
  from the course web page. This file is listed in the Appendix. You
  can then tabulate $[PQ|RS]$.

  Write first a subroutine/function/module that, given a reduced
  density matrix $D_0$, computes the corresponding Fock matrix
  $F(D_0)$, and then diagonalizes it.

  Thus, your subroutine should solve the standard eigenvalue problem
  $F(D_0)U = U\epsilon$ and return $U$ and the diagonal elements of
  $\epsilon$.

  (A useful check on the sanity of your implementation is that $F$ is
  Hermitian, $F_{PQ} = F_{QP}^*$.)

  Make sure the eigenvalues are sorted in increasing order (the
  eigenvectors also need to be in corresponding order). Hint: if you
  use LAPACK library calls, it can be useful to eliminate some
  numerical noise that can cause LAPACK to call a non-symmetric
  eigenvalue solver. For example, you can explicitly symmetrize
  the matrix by diagonalizing $(F+F^H)/2$.

\item[3g)]
  Write a program that performs SCF iterations,
  \begin{equation}
    F(U^{(k-1)}) U^{(k)} = U^{(k)}\epsilon^{(k)},
  \end{equation}
  with $U^{(0)}=I$, the identity matrix.

  You will need approximately 15 iterations for convergence. Measure
  the convergence via
  \begin{equation}
    \delta_k = \max_P\{|\epsilon^{(k)}_P - \epsilon^{(k-1)}_P|\}.
  \end{equation}
  Plot $\delta_k$ versus $k$, with $\delta_k$ on a log-axis.
  
  At each iteration $k=1,2,\cdots$, compute $E_\text{RHF}$. Discuss the
  convergence properties.

  Compare the converged RHF energy with the results obtained for CIS,
  and with the numerically exact values listed elsewhere in this project.

  Make sure your program handles both $N=2$ and $N=6$, and that you
  discuss both cases.
\end{enumerate}


\section*{Appendix:  Table of matrix elements and expressions for the expectation value of the two-body matrix elements}


\subsection*{Explicit values for the two-body matrix elements} 

The expectation values for the two-body matrix elements are computed
using the closed-form expression in the article of
E.~Anisimovas and A.~Matulis, J.~Phys.: Condens.~Matter {\bf 10}, 601
(1998).

Each line in the table below is on the form
\begin{verbatim}
P Q R S ( psi_P psi_Q|w|psi_R psi_S ).
\end{verbatim}
Only the non-zero matrix elements are listed. Note that the integrals
below do not include spin. The indices $P$, $Q$, etc, are 
related to the quantum numbers $(n,m)$ as $1=(n=0,m=0)$, $2=(n=0,
m=-1)$, $3=(n=0, m=+1)$, $4=(n=0,m=-2)$, $5=(n=1,m=0)$, and
$6=(n=0,m=+2)$ with no spin degrees of freedom. This table can also be
found in text format at the webpage of the course, see the file
\texttt{coulomb.dat} available from the course page.  
\begin{verbatim}
1 1 1 1     1.253314137
1 1 1 5    0.3133285343
1 1 2 3    0.3133285343
1 1 3 2    0.3133285343
1 1 4 6    0.1174982004
1 1 5 1    0.3133285343
1 1 5 5    0.2349964007
1 1 6 4    0.1174982004
1 2 1 2     0.939985603
1 2 2 1    0.3133285343
1 2 2 5    -0.07833213358
1 2 3 4    0.2769459142
1 2 4 3    0.1661675485
1 2 5 2    0.07833213358
1 3 1 3     0.939985603
1 3 2 6    0.2769459142
1 3 3 1    0.3133285343
1 3 3 5    -0.07833213358
1 3 5 3    0.07833213358
1 3 6 2    0.1661675485
1 4 1 4     0.744155269
1 4 2 2    0.2769459142
1 4 4 1    0.1174982004
1 4 4 5    -0.08812365028
1 4 5 4    -0.009791516698
1 5 1 1    0.3133285343
1 5 1 5    0.8616534694
1 5 2 3    -0.07833213358
1 5 3 2    -0.07833213358
1 5 4 6    -0.08812365028
1 5 5 1    0.2349964007
1 5 5 5    0.1370812338
1 5 6 4    -0.08812365028
1 6 1 6     0.744155269
1 6 3 3    0.2769459142
1 6 5 6    -0.009791516698
1 6 6 1    0.1174982004
1 6 6 5    -0.08812365028
2 1 1 2    0.3133285343
2 1 2 1     0.939985603
2 1 2 5    0.07833213358
2 1 3 4    0.1661675485
2 1 4 3    0.2769459142
2 1 5 2    -0.07833213358
2 2 1 4    0.2769459142
2 2 2 2    0.8616534694
2 2 4 1    0.2769459142
2 2 4 5    -0.1523202528
2 2 5 4    -0.1523202528
2 3 1 1    0.3133285343
2 3 1 5    -0.07833213358
2 3 2 3    0.8616534694
2 3 3 2    0.2349964007
2 3 4 6    0.3035370176
2 3 5 1    -0.07833213358
2 3 5 5    0.1370812338
2 3 6 4    0.1468727505
2 4 2 4    0.7539467857
2 4 4 2    0.3035370176
2 5 1 2    -0.07833213358
2 5 2 1    0.07833213358
2 5 2 5    0.7245722356
2 5 3 4    -0.1246256614
2 5 4 3    -0.1523202528
2 5 5 2    0.1370812338
2 6 1 3    0.2769459142
2 6 2 6    0.7539467857
2 6 3 1    0.1661675485
2 6 3 5    -0.1246256614
2 6 5 3    -0.1523202528
2 6 6 2    0.1468727505
3 1 1 3    0.3133285343
3 1 2 6    0.1661675485
3 1 3 1     0.939985603
3 1 3 5    0.07833213358
3 1 5 3    -0.07833213358
3 1 6 2    0.2769459142
3 2 1 1    0.3133285343
3 2 1 5    -0.07833213358
3 2 2 3    0.2349964007
3 2 3 2    0.8616534694
3 2 4 6    0.1468727505
3 2 5 1    -0.07833213358
3 2 5 5    0.1370812338
3 2 6 4    0.3035370176
3 3 1 6    0.2769459142
3 3 3 3    0.8616534694
3 3 5 6    -0.1523202528
3 3 6 1    0.2769459142
3 3 6 5    -0.1523202528
3 4 1 2    0.2769459142
3 4 2 1    0.1661675485
3 4 2 5    -0.1246256614
3 4 3 4    0.7539467857
3 4 4 3    0.1468727505
3 4 5 2    -0.1523202528
3 5 1 3    -0.07833213358
3 5 2 6    -0.1246256614
3 5 3 1    0.07833213358
3 5 3 5    0.7245722356
3 5 5 3    0.1370812338
3 5 6 2    -0.1523202528
3 6 3 6    0.7539467857
3 6 6 3    0.3035370176
4 1 1 4    0.1174982004
4 1 2 2    0.2769459142
4 1 4 1     0.744155269
4 1 4 5    -0.009791516698
4 1 5 4    -0.08812365028
4 2 2 4    0.3035370176
4 2 4 2    0.7539467857
4 3 1 2    0.1661675485
4 3 2 1    0.2769459142
4 3 2 5    -0.1523202528
4 3 3 4    0.1468727505
4 3 4 3    0.7539467857
4 3 5 2    -0.1246256614
4 4 4 4    0.7160046585
4 5 1 4    -0.08812365028
4 5 2 2    -0.1523202528
4 5 4 1    -0.009791516698
4 5 4 5    0.6682710146
4 5 5 4    0.1395291129
4 6 1 1    0.1174982004
4 6 1 5    -0.08812365028
4 6 2 3    0.3035370176
4 6 3 2    0.1468727505
4 6 4 6    0.7160046585
4 6 5 1    -0.08812365028
4 6 5 5    0.1395291129
4 6 6 4    0.1285136567
5 1 1 1    0.3133285343
5 1 1 5    0.2349964007
5 1 2 3    -0.07833213358
5 1 3 2    -0.07833213358
5 1 4 6    -0.08812365028
5 1 5 1    0.8616534694
5 1 5 5    0.1370812338
5 1 6 4    -0.08812365028
5 2 1 2    0.07833213358
5 2 2 1    -0.07833213358
5 2 2 5    0.1370812338
5 2 3 4    -0.1523202528
5 2 4 3    -0.1246256614
5 2 5 2    0.7245722356
5 3 1 3    0.07833213358
5 3 2 6    -0.1523202528
5 3 3 1    -0.07833213358
5 3 3 5    0.1370812338
5 3 5 3    0.7245722356
5 3 6 2    -0.1246256614
5 4 1 4    -0.009791516698
5 4 2 2    -0.1523202528
5 4 4 1    -0.08812365028
5 4 4 5    0.1395291129
5 4 5 4    0.6682710146
5 5 1 1    0.2349964007
5 5 1 5    0.1370812338
5 5 2 3    0.1370812338
5 5 3 2    0.1370812338
5 5 4 6    0.1395291129
5 5 5 1    0.1370812338
5 5 5 5    0.7490510274
5 5 6 4    0.1395291129
5 6 1 6    -0.009791516698
5 6 3 3    -0.1523202528
5 6 5 6    0.6682710146
5 6 6 1    -0.08812365028
5 6 6 5    0.1395291129
6 1 1 6    0.1174982004
6 1 3 3    0.2769459142
6 1 5 6    -0.08812365028
6 1 6 1     0.744155269
6 1 6 5    -0.009791516698
6 2 1 3    0.1661675485
6 2 2 6    0.1468727505
6 2 3 1    0.2769459142
6 2 3 5    -0.1523202528
6 2 5 3    -0.1246256614
6 2 6 2    0.7539467857
6 3 3 6    0.3035370176
6 3 6 3    0.7539467857
6 4 1 1    0.1174982004
6 4 1 5    -0.08812365028
6 4 2 3    0.1468727505
6 4 3 2    0.3035370176
6 4 4 6    0.1285136567
6 4 5 1    -0.08812365028
6 4 5 5    0.1395291129
6 4 6 4    0.7160046585
6 5 1 6    -0.08812365028
6 5 3 3    -0.1523202528
6 5 5 6    0.1395291129
6 5 6 1    -0.009791516698
6 5 6 5    0.6682710146
6 6 6 6    0.7160046585
\end{verbatim}
\end{document}







