dofile("../definitions.Quanty")

-- we define an Hamiltonian
Hamiltonian = OppU + 0.1 * Oppldots + 0.000001 * (2*OppSz+OppLz)

-- number of states calculated
Npsi=15

-- In order to make sure we have a filling of 2
-- electrons we need to define some restrictions
StartRestrictions = {Nf, Nb,  {"111111",2,2}}

-- We now can create the lowest Npsi eigenstates:
psiList = Eigensystem(Hamiltonian, StartRestrictions, Npsi)

-- We define a list of some operators to look at their expectation values
oppList={Hamiltonian, OppSsqr, OppLsqr, OppJsqr, OppSz, OppLz}

-- after we've created the eigen states we can look 
-- at a list of their expectation values
print(" <E>     <S^2>   <L^2>   <J^2>   <S_z>   <L_z>");
for i=1,#psiList do
  for j=1,#oppList do
    io.write(string.format("%7.4f ",psiList[i]*oppList[j]*psiList[i]))
  end
  io.write("\n")
end