function Tz(indexup, indexdn, NF) if #indexup ~= #indexdn then error("Length of index up must be equal to length of index dn and equal to 2l+1 in Tz") end local l=(#indexup-1)/2 if not IntegerQ(l) then error("Length of index must be equal to 2l+1 in Tz") end local function Sx(m1,m2) return NewOperator("Number",NF,{indexup[m1+l+1],indexdn[m1+l+1]}, {indexdn[m2+l+1],indexup[m2+l+1]},{1/2,1/2}) end local function Sy(m1,m2) return NewOperator("Number",NF,{indexup[m1+l+1],indexdn[m1+l+1]}, {indexdn[m2+l+1],indexup[m2+l+1]},{-I/2,I/2}) end local function Sz(m1,m2) return NewOperator("Number",NF,{indexup[m1+l+1],indexdn[m1+l+1]}, {indexup[m2+l+1],indexdn[m2+l+1]},{1/2,-1/2}) end local opp = 0 * NewOperator("Number",NF,0,0) for m1 = -l,l do for m2 = -l,l do opp = opp -2 * ( SlaterCoefficientC({l,m1},{2, 0},{l,m2}) )*Sz(m1,m2) - ( sqrt(3/2)*SlaterCoefficientC({l,m1},{2,-1},{l,m2}) -sqrt(3/2)*SlaterCoefficientC({l,m1},{2, 1},{l,m2}) )*Sx(m1,m2) -I * ( sqrt(3/2)*SlaterCoefficientC({l,m1},{2,-1},{l,m2}) +sqrt(3/2)*SlaterCoefficientC({l,m1},{2, 1},{l,m2}) )*Sy(m1,m2) end end opp.Chop() opp.Name="Tz" return opp end