Coding of Spin Spin correlation function

[woodydewar]

Homework Statement

Code the spin spin correlation function (see below) for a Spin chain. The states are represented by integers.

Relevant Equations

$\langle \vec{S}_i \vec{S}_j /rangle$

I tried to code spinoperators who act like $S_x^iS_x^j$ (y and z too) and to apply them to the states, which works fine. I am not sure about how to code the expectation value in the product Space. Has anyone pseudo Code to demonstrate that?

 

[yezia]

Hello woodydewar,

Are we talking about quantum averages or statistical averages ?
Quantum average is pretty straightforward, you just give as input a state and apply a really big (but sparse) matrix on it. I believe your problem is on how to build such a matrix. For instance, on the $i$-th site of a $N$ $1/2$ spin-chain, in the product space, the spin operator $\mathrm{S}_x$ takes the form :
$$\mathrm{S}_x \leftarrow \mathrm{1}_2^{\otimes i}\otimes\mathrm{S}_x\otimes\mathrm{1}_2^{\otimes N-i-1}.$$
You can easilly generalize for any spin value, any direction and subsequently build the vector of matrices $\mathbf{S}_i$. Computing $\bra{\Psi}\mathbf{S}_i\cdot\mathbf{S}_j\ket{\Psi}$ finally consists of matrix multiplications.

On the other hand, statistical average requires the use of Monte Carlo methods (for very large $N$). If you are instead working at finite temperature, I can add details on how to sample the configuration space.

 

[woodydewar]

Hi yezia,

thanks! I already built the matrix and have now the groundstate, where I want to measure these correlations. I want to measure $$ \braket{\psi \vert \mathbf{S}_i \mathbf{S}_j \vert \psi} $$
but am not sure, how to do it and to get all contributions. Could you write some pseudocode and explain why it works?

What I did until now:
I have a state $$ \ket{\psi} =\alpha \ket{\phi} $$ where phi is representing the basis spinstate (up down up ...) in the z-basis and alpha is the coefficient in the groundstate of the matrix weighting each basis state.

I defined two functions: one is simulating the action of the pauli matrices on phi (called state), the other gives the factor (called factor). I did this for x,y, and z and want to add them at the end.

I tried something like (here an example for y)

for i in spins:
for j in spins:
for alpha in range(length(Groundstate)):
SpincorrelationYY += GS[alpha] * GS[state(alpha,i,j)] * Factor(alpha, i,j)

 

[yezia]

I did not perfectly understand what you did, but for now I'll just try to explain you one way to compute this expectation value. Don't worry : I will try to go more in depth with your approach later.

You managed to build the matrix representation of $\mathbf{S}_i\cdot\mathbf{S}_j$, I will call it $\mathrm{M}$.
You already know the state you are working with - the ground-state $\ket{\Psi}$ - , you probably already have its vectorial representation : I will call it $u$ (if you numerically computed it, what your program furnishes is already its vectorial representation).
$\mathrm{M}$ is a matrix, $u$ is a column vector. All what is left is basic matrix multiplication :
$$\langle\mathbf{S}_i\cdot\mathbf{S}_j\rangle_{\Psi} = u^\dagger \mathrm{M} u.$$
In python, it will look like something like this (everything is a $\mathtt{numpy}$ matrix):
correl = numpy.dot(u.H, np.dot(M, u))

 

[yezia]

I still don't really understand what you did, but I guess you're basically trying to code matrix multiplication. Sure, why not; but there a lot of libraries already doing that for you.

 

[woodydewar]

Sorry for my late reply and thanks you. I did it as you proposed, but all the results are negative. I am pretty sure that my matrix is correct, as I build my hamiltonian exactly the same way using the same structures and the groundstate energy matches the right value. Since I want to calculate this for a Heisenberg Model, the matrix structure of H and the "corelaton matrix" $$ \mathbf{S}_i \mathbf{S}_j$$ is similar.

Thats why I am so confused.

 

[yezia]

I will assume your hamiltonian to be of the simplest form:
$$
\mathrm{H} = -J\sum_i \mathbf{S}_i\cdot\mathbf{S}_{i+1}.
$$
The correlations depend not only on the sign of $J$ but also the spin magnitude. Can you at least test the influence of $J$ ? Please also provide the way you build the matrix representation of $\mathbf{S}_i\cdot\mathbf{S}_j$.

 

[woodydewar]

Hi,
I tested the influence of J regarding to the Groundstate energy of the ferromagnetic Heisenberg chain, it is as it should be.

Of course you can rewrite the Heisenberg Interaction as
$$S_xS_x +S_yS_y +S_zS_z$$

I built the matrix using the following c++ functions

   int sxsq(int state, int i, int j){
    state ^= 1UL << i;
    state ^= 1UL << j;
    return state;    };

This function simulates the acting of the $$S^x_iS^x_j$$ operators. The "weight" is always 0.25.
As the y-Pauli matrix has the same shape, the resulting state is the same, but with another weight. Therefore I defined:

    std::complex<double> sysqvalue(int state, int i, int j){
    if (((state >> i) & 1U) == ((state >> j) & 1U))
    {return std::complex<double>(-0.25,0.0);} else {return std::complex<double>(0.25,0.0);}
    };

For z-part I got:

std::complex<double> szsqvalue(int state, int i, int j){
    if (((state >> i) & 1U) == ((state >> j) & 1U))
    {return std::complex<double>(0.25,0.0);} else {return std::complex<double>(-0.25,0.0);}
};

All this put together is multiplied with the groundstate vector, as you said before, into a vector "vec" with the spins i,j

            for (int state = 0; state < dim; state++)
            {  
                vec[state] = NeutralElement;
                vec[state] += GS[sxsq(state, i,j)] * 0.25;
                vec[state] += GS[sxsq(state, i,j)] * sysqvalue(state, i,j);
                vec[state] += GS[state] * szsqvalue(state, i,j);

            };

which is then multiplied with the conjugated groundstate:

            for (int state = 0; state < dim; state++)
            {  
                SC += std::conj(GS[state])*vec[state];
            };

A far as I know this corresponds to the static structure factor of a spin chain, which should not be negative. But I get negative values for the sum of all terms.

 

[woodydewar]

Hi yezia, do you see any error/mistake in my concept?

 

[yezia]

Dear woodydewar,

I'm sorry for the late reply. I forgot this thread was not solved and did not really pay attention to it/my mails.
I am currently not home, but I will take a closer look to your method before next monday.

 

[yezia]

I'm starting to understand what you did. I believe the general outline of your answer (mainly the reduction from $2^{2N}$ to $2^N$ operations) is correct, I don't see any blatant errors. Now I'm not sure about the way you "select" the non-zero elements (for $\mathrm{S}_x$ and $\mathrm{S}_y$).
Could you please explain this to me ?

And finally: can you tell me how does the sign of $J$ affect your correlations ?

 

[woodydewar]

Thank you!

Sorry, I thought you were talking about how J affects my eigenvalues of the Heisenberg Model. I expect a large positive summed-over correlation/Structurefactor in the ferromagnetic case and a positive but smaller Structurefactor in the antiferromagnetic case. But I get both times small negative values, which has to be wrong.

For the scalarproduct to be nonzero between two states, every spin has to be the same oriented in the two states (they have to be the same component in the vector rep). Therefore I am selecting the "position"/integer representation (equivalent) of the state which corresponds to the state which results from application of the sxsx operators to the state I am currently interested in.

 

[yezia]

Yes, in the FM case, the correlations should be positive. In the AFM case, the sign should alternate. The fact that your answer does not depend on the sign of $J$ suggests a problem. The GS could be wrong, or even the functions you coded. But if you say you built the Hamiltonian with these same functions, and you obtained the correct GS; I believe they must be efficient. You can still test for $N=3$ for instance, and compare the acting of your functions with the formal expressions of the matrices (with $\hbar=2$):
$$
\mathrm{S}_1^x\mathrm{S}_2^x =
\begin{pmatrix}
0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \\
1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 1 & 0 & 0 & 0 & 0
\end{pmatrix}
\begin{pmatrix}
0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 \\
1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \\
0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 1 & 0 & 0
\end{pmatrix}.
$$

I have real trouble understanding this sentence:

Therefore I am selecting the "position"/integer representation (equivalent) of the state which corresponds to the state which results from application of the sxsx operators to the state I am currently interested in.

And how does that explain the selection of the non-zero elements for $\mathrm{S}_{x/y}$ (the $\mathtt{sxsq}$ function) ?

I usually solved the spin chain by the use of sparse matrices in Python, and your explanations are at the moment outside the scope of my understanding.
Someone more competent than me is welcome to help our friend here.

 

[woodydewar]

Oh sorry, I thought you were talking about something differently. How the operator works is simply explained:
A state is represented as an integer in bit representation. The "shape" of Sx or Sy is something like $$ \begin{pmatrix} 0 & a\\ b & 0\end{pmatrix}$$ which means applying Sx/Sy to up results in down and the way around with different factors. This can be simply implemented in this bit representation by a bitwise xor.

I'll test for the two spin case and tell you.

 

[woodydewar]

Hi Yezia, many thanks! I solved the issue. The error was in another loop. Thank you for your help and patience with me!

 

[yezia]

I'm glad our long discussion managed to help you!