====== MeanFieldOperator ======

###
//MeanFieldOperator($O$, $\rho$)// creates the mean-field version of operator $O$ with the corresponding density matrix $\rho$.
$rho$ stores the expectation values of $a^{\dagger}_{\tau}a^{\phantom{\dagger}}_{\tau'}$, a table of dimensions $NFermion$ by $NFermion$.

Any two particle parts of the operator will be replaced in mean-field, using the Hartree-Fock approximation by:
\begin{eqnarray}
a^{\dagger}_{i}a^{\dagger}_{j}a^{\phantom{\dagger}}_{k}a^{\phantom{\dagger}}_{l} &\to&\\
\nonumber &-& a^{\dagger}_{i}a^{\phantom{\dagger}}_{k} \langle a^{\dagger}_{j}a^{\phantom{\dagger}}_{l} \rangle \\
\nonumber &+& a^{\dagger}_{i}a^{\phantom{\dagger}}_{l} \langle a^{\dagger}_{j}a^{\phantom{\dagger}}_{k} \rangle \\
\nonumber &+& a^{\dagger}_{j}a^{\phantom{\dagger}}_{k} \langle a^{\dagger}_{i}a^{\phantom{\dagger}}_{l} \rangle \\
\nonumber &-& a^{\dagger}_{j}a^{\phantom{\dagger}}_{l} \langle a^{\dagger}_{i}a^{\phantom{\dagger}}_{k} \rangle \\
\nonumber &-& \langle a^{\dagger}_{i}a^{\phantom{\dagger}}_{l} \rangle \langle a^{\dagger}_{j}a^{\phantom{\dagger}}_{k} \rangle \\
\nonumber &+& \langle a^{\dagger}_{i}a^{\phantom{\dagger}}_{k} \rangle \langle a^{\dagger}_{j}a^{\phantom{\dagger}}_{l} \rangle 
\end{eqnarray}

If the option AddDFTSelfInteraction was set to true more terms are added to the Mean-Field Operator, namely
\begin{equation}
\sum_{m} U \langle a^\dagger_m a^{\phantom{\dagger}}_m \rangle a^\dagger_m a^{\phantom{\dagger}}_m
\end{equation}
where
\begin{equation}
U
=
\left(
\frac{N_{Fermion} (N_{Fermion}-1)}{2}
\right)^{-1}
\sum_{m,n}
\left(
U_{m\,n\,n\,m}
-
U_{m\,n\,m\,n}
\right)
\end{equation}
is the average interaction energy electrons have with one another.
###

===== Input =====

  * $O$ : Operator
  * $rho$ : Matrix (Table of Table of length $O.NF$) of doubles
  * Possible options are:
    * "AddDFTSelfInteraction" bool defining if the electron self-interaction is to be included. (Standard false)

===== Output =====

  * $O_{MF}$ The mean-field approximated operator

===== Example =====



==== Input ====
<code Quanty Example.Quanty>
NF = 4
op = NewOperator("Number",NF,{1},{1},{0.1+I}) + NewOperator("U",NF,{0},{1},{5}) + 3
rho = {{0.7,0.3+I,0,0},{0.3-I,0.4,0,0},{0,0,0,0},{0,0,0,0}}

print("Full Operator:")
print(op)
print("\nDensity:")
print(rho)
print("\nMeanFieldOperator:")
print( MeanFieldOperator(op, rho) )
print("\nMeanFieldOperator with electron self-interaction:")
print( MeanFieldOperator(op, rho, {{"AddDFTSelfInteraction",true}}) )
</code>

==== Result ====
<file Quanty_Output>
Full Operator:

Operator: CrAn
QComplex         =          2 (Real==0 or Complex==1 or Mixed==2)
MaxLength        =          4 (largest number of product of lader operators)
NFermionic modes =          4 (Number of fermionic modes (site, spin, orbital, ...) in the one particle basis)
NBosonic modes   =          0 (Number of bosonic modes (phonon modes, ...) in the one particle basis)

Operator of Length   0
QComplex      =          0 (Real==0 or Complex==1)
N             =          1 (number of operators of length   0)
|  3.000000000000000E+00

Operator of Length   2
QComplex      =          1 (Real==0 or Complex==1)
N             =          1 (number of operators of length   2)
C  1 A  1 |  1.000000000000000E-01  1.000000000000000E+00

Operator of Length   4
QComplex      =          0 (Real==0 or Complex==1)
N             =          1 (number of operators of length   4)
C  1 C  0 A  1 A  0 | -5.000000000000000E+00



Density:
{ { 0.7 , (0.3 + 1 I) , 0 , 0 } , 
  { (0.3 - 1 I) , 0.4 , 0 , 0 } , 
  { 0 , 0 , 0 , 0 } , 
  { 0 , 0 , 0 , 0 } }

MeanFieldOperator:

Operator: 
QComplex         =          0 (Real==0 or Complex==1 or Mixed==2)
MaxLength        =          4 (largest number of product of lader operators)
NFermionic modes =          4 (Number of fermionic modes (site, spin, orbital, ...) in the one particle basis)
NBosonic modes   =          0 (Number of bosonic modes (phonon modes, ...) in the one particle basis)

Operator of Length   0
QComplex      =          0 (Real==0 or Complex==1)
N             =          1 (number of operators of length   0)
|  1.255000000000000E+01

Operator of Length   2
QComplex      =          1 (Real==0 or Complex==1)
N             =          4 (number of operators of length   2)
C  1 A  1 | -3.400000000000000E+00  1.000000000000000E+00
C  1 A  0 |  1.500000000000000E+00  5.000000000000000E+00
C  0 A  1 |  1.500000000000000E+00 -5.000000000000000E+00
C  0 A  0 | -2.000000000000000E+00  0.000000000000000E+00



MeanFieldOperator with electron self-interaction:

Operator: 
QComplex         =          0 (Real==0 or Complex==1 or Mixed==2)
MaxLength        =          4 (largest number of product of lader operators)
NFermionic modes =          4 (Number of fermionic modes (site, spin, orbital, ...) in the one particle basis)
NBosonic modes   =          0 (Number of bosonic modes (phonon modes, ...) in the one particle basis)

Operator of Length   0
QComplex      =          0 (Real==0 or Complex==1)
N             =          1 (number of operators of length   0)
|  1.255000000000000E+01

Operator of Length   2
QComplex      =          1 (Real==0 or Complex==1)
N             =          4 (number of operators of length   2)
C  1 A  1 | -3.066666666666666E+00  1.000000000000000E+00
C  1 A  0 |  1.500000000000000E+00  5.000000000000000E+00
C  0 A  1 |  1.500000000000000E+00 -5.000000000000000E+00
C  0 A  0 | -1.416666666666667E+00  0.000000000000000E+00
</file>

===== Table of contents =====
{{indexmenu>.#1}}