Matrix form for creation and annihilation operator with spins.

This might be little geared towards physics but if anyone can help, it would so much helpful for me. So basically I want the matrix form of the following Hamiltonian ., where creates an electron at site i with spin σ which can be up (↑) or down (↓). For this I need the matrix form of the operators involved. For people who have been working with this type of Hamiltonian are aware that it follows aniti-commutation relations

What I can do now is that I can write the operators without the spin as below but i need help with the above when spin is also involved.

N = 4;
%n = 1:N;
a = diag(sqrt(1:N),1)% Creates matrix with elements of the vector on first super diagonal(i.e 1)
ad = a' %adagger

Is there anyway to do that. How should I go on to create the Hamiltonian for that. Are there any helpful packages to do that. Any help would be appreciated.

Thank you.

Answer 1

[Matt J]

Not sure this will help, but if you can write a function that implements the operator, then his File Exchange offering,

https://www.mathworks.com/matlabcentral/fileexchange/44669-func2mat-convert-linear-function-to-matrix

will convert it to matrix form. For example,

fn=@(x) 7*x(1:3)-3*x(4:6); %A function, doesn't have to be anonymous

A=func2mat(fn, ones(6,1));
full(A)

ans = 3×6

7 0 0 -3 0 0 0 7 0 0 -3 0 0 0 7 0 0 -3

Answer 2

[Marc]

Im sure this comes too late but in case someone else has this question aswell (I know I did) I will try to answer it.

Unlike bosons, the Hilbert space for spins is finite dimensional and if you work with states all you require is a 2x2 matrix: ={0 1; 0 0}. The spin algebra for N spins follows using only the kron product that exists in default MATLAB.

For example, for 2 spins we can define spin lowering (raising) opeartors on each site through . On MATLAB we can define a global spin operator S and then act on it with the kron product to divide the Hilbert space into the 2 distinct spin subspaces:

S = diag(1,1); Sd = S'; %local spin lowering and raising operator

sigma1=kron(S,I); sigma1d=sigma1'; %spin operator on site 1, identity of site 2

sigma2=kron(I,S); sigma2d=sigma2'; %identity of site 1, spin operator on site 2

You can immediately extend this to N spins by iterating this process with the usual tensor product identities, namely for site :

where S is the standard 2x2 spin lowering operator and I the identity matrix.

You can check the spin algebra is satisfied in the 2-spin state above by simply calculating the anticommutators

sigma1d*sigma2d+sigma2d*sigma1d; sigma1*sigma2+sigma2*sigma1; sigma1d*sigma2+sigma2*sigma1d

which should yield zeroes for the diagonals and a kronecker delta for the cross terms.