# Generating Wick's theorem contractions

## Preliminaries

Let's start by importing wick&d and defining a Slater determinant reference


```python
import wicked as w
from IPython.display import display, Math, Latex

def latex(expr):
    """Function to render any object that has a member latex() function"""
    display(Math(expr.latex()))
    
w.reset_space()
w.add_space("o", "fermion", "occupied", ['i','j','k','l','m'])
w.add_space("v", "fermion", "unoccupied", ['a','b','c','d','e','f'])    
```

## Wick's theorem contractions

Consider the product of operators
\begin{equation}
\hat{F} \hat{T}_1 = 
\left( \sum_{pq}^\mathbb{O \cup V} f_{p}^{q} \{ \hat{a}^{p}_{q} \} \right)
\left( \sum_{i}^\mathbb{O} \sum_{a}^\mathbb{V} t_{a}^{i} \{ \hat{a}^{a}_{i} \} \right)
\end{equation}
which we can create with the following code:


```python
F = w.utils.gen_op('f',1,'ov','ov')
T1 = w.op('t',['v+ o'])
FT1 = F @ T1
print(FT1)
```

    + f { v+ v } t { v+ o }
    + f { v+ o } t { v+ o }
    + f { o+ v } t { v+ o }
    + f { o+ o } t { v+ o }


To apply Wick's theorem to this product, you need to create a WickTheorem object. You can compute the contraction with the `contraction` function.
This function takes as arguments the operator expression to contract and the minimum and maximum number of second quantized operators in the final expression (`minrank/maxrank`). Here we request all terms with zero to four second quantized operators:


```python
wt = w.WickTheorem()
expr = wt.contract(FT1,minrank=0,maxrank=4)
expr
```




    f^{o1}_{o0} t^{o2}_{v0} a+(o0) a+(v0) a-(o2) a-(o1)
    +f^{v0}_{o0} t^{o0}_{v0}
    +f^{v0}_{o0} t^{o1}_{v0} a+(o0) a-(o1)
    -f^{v0}_{o0} t^{o0}_{v1} a+(v1) a-(v0)
    -f^{v0}_{o0} t^{o1}_{v1} a+(o0) a+(v1) a-(v0) a-(o1)
    -f^{o0}_{o1} t^{o1}_{v0} a+(v0) a-(o0)
    +f^{o0}_{v0} t^{o1}_{v1} a+(v0) a+(v1) a-(o1) a-(o0)
    +f^{v1}_{v0} t^{o0}_{v1} a+(v0) a-(o0)
    -f^{v1}_{v0} t^{o0}_{v2} a+(v0) a+(v2) a-(v1) a-(o0)



You can get a LaTeX expression with the `latex()` function


```python
print(expr.latex())
```

    + {f}^{j}_{i} {t}^{k}_{a} \hat{a}^{i} \hat{a}^{a} \hat{a}_{k} \hat{a}_{j} \\ 
    + {f}^{a}_{i} {t}^{i}_{a} \\ 
    + {f}^{a}_{i} {t}^{j}_{a} \hat{a}^{i} \hat{a}_{j} \\ 
    - {f}^{a}_{i} {t}^{i}_{b} \hat{a}^{b} \hat{a}_{a} \\ 
    - {f}^{a}_{i} {t}^{j}_{b} \hat{a}^{i} \hat{a}^{b} \hat{a}_{a} \hat{a}_{j} \\ 
    - {f}^{i}_{j} {t}^{j}_{a} \hat{a}^{a} \hat{a}_{i} \\ 
    + {f}^{i}_{a} {t}^{j}_{b} \hat{a}^{a} \hat{a}^{b} \hat{a}_{j} \hat{a}_{i} \\ 
    + {f}^{b}_{a} {t}^{i}_{b} \hat{a}^{a} \hat{a}_{i} \\ 
    - {f}^{b}_{a} {t}^{i}_{c} \hat{a}^{a} \hat{a}^{c} \hat{a}_{b} \hat{a}_{i}


Or render the expression directly in a notebook


```python
latex(expr)
```


$\displaystyle + {f}^{j}_{i} {t}^{k}_{a} \hat{a}^{i} \hat{a}^{a} \hat{a}_{k} \hat{a}_{j} \\ 
+ {f}^{a}_{i} {t}^{i}_{a} \\ 
+ {f}^{a}_{i} {t}^{j}_{a} \hat{a}^{i} \hat{a}_{j} \\ 
- {f}^{a}_{i} {t}^{i}_{b} \hat{a}^{b} \hat{a}_{a} \\ 
- {f}^{a}_{i} {t}^{j}_{b} \hat{a}^{i} \hat{a}^{b} \hat{a}_{a} \hat{a}_{j} \\ 
- {f}^{i}_{j} {t}^{j}_{a} \hat{a}^{a} \hat{a}_{i} \\ 
+ {f}^{i}_{a} {t}^{j}_{b} \hat{a}^{a} \hat{a}^{b} \hat{a}_{j} \hat{a}_{i} \\ 
+ {f}^{b}_{a} {t}^{i}_{b} \hat{a}^{a} \hat{a}_{i} \\ 
- {f}^{b}_{a} {t}^{i}_{c} \hat{a}^{a} \hat{a}^{c} \hat{a}_{b} \hat{a}_{i}$


Let's evaluate a more complex expression
\begin{equation}
\langle \Phi | \hat{\Lambda}_1 [\hat{F}, \hat{T}_1] | \Phi \rangle
\end{equation}
where 
\begin{equation}
\hat{\Lambda}_1 =
\sum_{i}^\mathbb{O} \sum_{a}^\mathbb{V} \lambda_{i}^{a} \{ \hat{a}^{i}_{a} \} 
\end{equation}


```python
L1 = w.op('lambda',['o+ v'])
opexpr = L1 @ w.commutator(F,T1)
expr = wt.contract(opexpr,0,0)
latex(expr)
```


$\displaystyle - {f}^{j}_{i} {\lambda}^{a}_{j} {t}^{i}_{a} \\ 
+ {f}^{b}_{a} {\lambda}^{a}_{i} {t}^{i}_{b}$


## Example: Coupled cluster energy expression

Here we compute the coupled cluster energy expression
\begin{equation}
E = \langle \Phi | e^{-\hat{T}} \hat{H} e^{\hat{T}} | \Phi \rangle
= E_0 + \sum_{i}^\mathbb{O} \sum_{a}^\mathbb{V} f^{a}_{i} t^{i}_{a} + 
\frac{1}{4} \sum_{ij}^\mathbb{O} \sum_{ab}^\mathbb{V}
 (t^{i j}_{a b} + 2 t^{i}_{a} t^{j}_{b}) v^{a b}_{i j} 
\end{equation}
We first define the full Hamiltonian operator and the cluster operator


```python
E0 = w.op("E_0",[""])
F = w.utils.gen_op('f',1,'ov','ov')
V = w.utils.gen_op('v',2,'ov','ov')
H = E0 + F + V

T = w.op("t",["v+ o", "v+ v+ o o", "v+ v+ v+ o o o"])
```

Next, we compute the operator $e^{-\hat{T}} \hat{H} e^{\hat{T}}$ up to the term quadratic in $\hat{T}$ and compute the fully contracted term


```python
Hbar = w.bch_series(H,T,2)
expr = wt.contract(Hbar,0,0)
latex(expr)
```


$\displaystyle + {E_0}^{}_{} \\ 
+ {f}^{a}_{i} {t}^{i}_{a} \\ 
+\frac{1}{2} {t}^{i}_{a} {t}^{j}_{b} {v}^{a b}_{i j} \\ 
+\frac{1}{4} {t}^{i j}_{a b} {v}^{a b}_{i j}$


ipynb
c1e2dc3b607c0c64580d3bd932ec5b05301b37b8
Jupyter Notebook
9bc0e13f6e45c86222ea95fdadf1cb66eb59862f
MIT
fevangelista/pyWicked
tutorials/03-WickTheorem.ipynb
9bc0e13f6e45c86222ea95fdadf1cb66eb59862f
MIT
fevangelista/pyWicked
tutorials/03-WickTheorem.ipynb
9bc0e13f6e45c86222ea95fdadf1cb66eb59862f
MIT
fevangelista/pyWicked
tutorials/03-WickTheorem.ipynb