I would like to ask how to properly calculate the off-diagonal terms in covariance matrix for the entangled Gaussian state?
E.g. from https://arxiv.org/abs/0810.0534v1 we have a coherent Gaussian state in the following form $$ |\psi\rangle = \sum_{n=0}^{\infty}\sqrt{\frac{N^n}{(N+1)^{n+1}}}|n\rangle_A|n\rangle_B $$ and the covariance matrix $$ V=\frac{1}{4} \left( \begin{array}{cccc} S & 0 & C & 0 \\ 0 & S & 0 &-C \\ C & 0 & S & 0 \\ 0 &-C & 0 & S \\ \end{array} \right) $$ where $S=2N+1$ and $C=2\sqrt{N(N+1)}$.
Using the definition of the covariance matrix $$ V_{ij}=\frac{1}{2}Tr[\hat \rho\{\hat q_i;\hat q_j\}] $$ (assuming zero displacement) where $\hat \rho$ is the appropriate density operator and vector $\hat q=(\hat X_A, \hat P_A, \hat X_B, \hat P_B)$ can be expressed using kvadrature, i.e. $\hat X = \hat a + \hat a^\dagger$, $\hat P = \hat a^\dagger - \hat a$. It is clear to me how I get the diagonal terms, but not the off-diagonal.
For example $$ V_{13}=\frac{1}{2}Tr[\hat \rho(\hat X_A\hat X_B + \hat X_B\hat X_A ) ] $$ since operators $A$ and $B$ commute $$ V_{13}=Tr[\hat \rho(\hat X_A\hat X_B ) ]\\ =Tr[\hat \rho(\hat a_A + \hat a_A^\dagger )(\hat a_B + \hat a_B^\dagger ) ] $$ but this combination of creation and annihilation operators change the states but never "return" and, therefore, trace will be zero.
Probably I do something trivially wrong, but I'm blind. Thanks.
Edit 1: in more details: My density matrix reads $$ \hat \rho = \sum_{n=0}^{\infty}\frac{N^n}{(N+1)^{n+1}} |n\rangle_A\langle n| |n\rangle_B \langle n| $$ Then the above described term 13 of the covariance matrix is $$ V_{13}=Tr\left[\sum_{n=0}^{\infty}\frac{N^n}{(N+1)^{n+1}} _B\langle n|_A\langle n| \hat a_A\hat a_B + \hat a_A^\dagger \hat a_B^\dagger + \hat a_A^\dagger\hat a_B + \hat a_A^\dagger\hat a_B |n\rangle_A|n\rangle_B \right]\\ Tr\left[\sum_{n=0}^{\infty}\frac{N^n}{(N+1)^{n+1}} _B\langle n|_A\langle n| \left( n|n-1\rangle_A|n-1\rangle_B + (n+1)|n+1\rangle_A|n+1\rangle_B + \sqrt{n(n+1)}|n+1\rangle_A|n-1\rangle_B + \sqrt{n(n+1)}|n-1\rangle_A|n+1\rangle_B \right) \right] $$ then the trace gives zero. Where I do a mistake? What is wrong?
