deletions | additions
diff --git a/RDMFT1.tex b/RDMFT1.tex
index 6bcaf33..3c453cd 100644
--- a/RDMFT1.tex
+++ b/RDMFT1.tex
...
Within Reduced Density Matrix Functional Theory (RDMFT) the total ground state energy is given as a functional of the one body reduced density matrix (1-RDM)
\begin{eqnarray}
\gamma(\mathbf{r},\mathbf{r'})=N\int\cdots\int d\mathbf{r_2}...d\mathbf{r_N} \Psi^*(\mathbf{r'},\mathbf{r_2}...\mathbf{r_N})\Psi(\mathbf{r},\mathbf{r_2}...\mathbf{r_N})
...
\end{eqnarray}
The natural orbitals $\phi_i(\mathbf{r})$ and their occupation numbers $n_i$ are the eigenfunctions and eigenvalues of the 1-RDM, respectively.
As the exact functional is unknown different approximate functionals are employed and minimized with respect to the 1-RDM in order to find the ground state energy. However, from the total energy
\begin{align} \begin{eqnarray}
E=\sum_{i=1}^\infty\int d\mathbf{r} n_{i}\phi^{*}_{i}(\mathbf{r})\left(-\frac{\nabla^2}{2}\right) \phi_{i}(\mathbf{r})+\sum_{i=1}^\infty \int d\mathbf{r}
V_{\mathrm{ext}}(\mathbf{r})n_{i}|\phi_{i}(\mathbf{r})|^{2}\\ V_{\mathrm{ext}}(\mathbf{r})n_{i}|\phi_{i}(\mathbf{r})|^{2}\nonumber\\
+\frac{1}{2}\sum_{i,j=1}^\infty n_{i} n_{j}\int d\mathbf{r} d\mathbf{r'} \frac{|\phi_{i}(\mathbf{r})|^{2} |\phi_{j}(\mathbf{r})|^{2}}{|\mathbf{r}-\mathbf{r'}|} + E_{xc}\left[\{n_{j}\},\{\phi_{j}\}\right]\label{eqn:energy}
\end{align} \end{eqnarray}
the part that needs to be approximated $E_{xc}[\gamma]$ comes only from the interaction term (contrary to KS-DFT), as the interacting kinetic energy can be explicitely expressed in terms of $\gamma$.
For closed-shell systems the necessary and sufficient conditions for the 1-RDM to be $N$-representable, i.e.\ to correspond to a $N$-electron wavefunction is that $ 0 \leq n_{i} \leq 2$ and
\begin{align} \begin{eqnarray}
\sum_{i=1}^{\infty}n_{i}=N.
\end{align} \end{eqnarray}
Note that within the RDMFT implementation in octopus only closed-shell systems are treated at the momment. Minimization of the energy functional of Eq.
\eqref{eqn:energy} ref{eqn:energy} is performed under the $N$-representability constraints and the orthonormality requierements of the natural orbitals,
\begin{eqnarray}
\langle \phi_{i} | \phi_{j}\rangle = \delta_{ij}
\end{eqnarray}
In practice, the minimization of the energy is not performed with respect to the 1-RDM but with respect to $n_{i}$ and $\phi_{i}$, separately. This does not lead to eigenvalue equation..The The bounds on the occupation numbers are automatically satisfied by setting $n_{i}$=2sin$^2
2\pi\theta_i$ 2\pi\vartheta_i$ and varying
$\theta_{i}$ $\vartheta_{i}$ without constraints. The occupation numbers summing up to the number of electrons is taken into account by using a Lagrange multiplyer
$\mu$. The $\mu$ and Lagrange multiplyers $\lambda_{ik}$ are used for the orthonormality constraints. Then one can define the following functional
\begin{eqnarray}
\Omega(N,\{\vartheta\} ,\{\phi_i(\mathbf{r})\})= E - \mu (\sum_i 2sin^2( 2\pi\vartheta_i)-N)-\sum_{ik} \lambda_{ik}(\langle\phi_k|\phi_i\rangle-\delta_{ki})
\end{eqnarray}
which has to be stationary with respect to variations in $\{\vartheta_i\}$, $\{\phi_i(\mathbf{r})\}$ and $\{\phi_i^{*}(\mathbf{r})\}$.
%\begin{eqnarray}
%\delta \Omega = 4sin(4\pi\vartheta_i)\big[\mu-\frac{\partial E}{\partial n_{i}}\big] d\vartheta+
%\sum_i \int d\mathbf{r} \delta \phi_i^{*}(\mathbf{r})\Big[\frac{\delta E}{\delta \phi_{i}^{*}(\mathbf{r})}-\sum_{k}\lambda_{ki}\phi_{k}(\mathbf{r})\Big ]+\nonumber\\
% \sum_i \int d\mathbf{r} \delta \phi_i(\mathbf{r})\Big[\frac{\delta E}{\delta \phi_{i}(\mathbf{r})}-\sum_{k}\lambda_{ki}\phi_{k}^{*}(\mathbf{r})\Big] = 0
%\end{eqnarray}
In practice, the minimization
of the energy is not performed with respect to the 1-RDM but with respect to $n_{i}$ and $\phi_{i}$, separately. Thus, for a fixed set of orbitals the energy functional is minimized with respect to occupation numbers
corresponds to and accordingly for a fixed set of occupations the
minimization energy functional is minimized with respect to variations of the
objectivefunctional $E(\{\theta_i\})-\mu (\sum_{i}(2sin(2\pi\theta_i)-N)$. orbitals until overall convergence is achieved. For the first step of occupation numbers minimization HF orbitals are used. As the correct $\mu$ is not known bisection is used. For every $\mu$, the objective functional is minimized with respect to
$\theta_i$ $\vartheta_i$ until
$\sum_{i}(2sin(2\pi\theta_i)-N=0$ $\sum_{i}(2sin(2\pi\vartheta_i)-N=0$ is
satisfied . satisfied.\par
The implementation of the natural orbital minimization follows a method by Piris and Ugalde (\cite{Piris}).
As one can show for fixed occupation numbers
\begin{eqnarray}
\lambda_{ki}=h_{ki}n_i+\int d\mathbf{r} \frac{\delta V_{ee}}{\delta \phi_i^{*}(\mathbf{r})}\phi_k^{*}(\mathbf{r}).
\end{eqnarray}
At the extremum, the matrix of the Lagrange multiplyers must be Hermitian, i.e. $\lambda_{ki}=\lambda_{ik}^{*}$
\begin{eqnarray}
\lambda_{ki}=h_{ki}n_i+\int d\mathbf{r} \frac{\delta V_{ee}}{\delta \phi_i^{*}(\mathbf{r})}\phi_k^{*}(\mathbf{r}).
\end{eqnarray}
At the extremum, the matrix of the Lagrange multiplyers must be Hermitian, i.e. $\lambda_{ki}=\lambda_{ik}^{*}$
\begin{eqnarray}
F_{ki}=\theta(i-k)(\lambda_{ki}-\lambda^{*}_{ik})+\theta(k-i)(\lambda^{*}_{ik}-\lambda_{ki})
\end{eqnarray}
