Differences

This shows you the differences between two versions of the page.

Link to this comparison view

documentation:tutorials:model_examples_from_physics:hubbard_model_on_a_dimer [2016/10/07 21:12] – created Maurits W. Haverkortdocumentation:tutorials:model_examples_from_physics:hubbard_model_on_a_dimer [2016/10/10 09:41] (current) – external edit 127.0.0.1
Line 1: Line 1:
 +{{indexmenu_n>1}}
 +====== Hubbard model on a dimer ======
  
 +###
 +This example shows the effect of inter-site correlations in a two site Hubbard model, aka a H$_2$ molecule on the smallest basis possible. The solution of this problem is analytically known. The large $U$ limit is known as the Heither-London approximation.
 +###
 +
 +###
 +The input code is:
 +<code Quanty H2.Quanty>
 +Verbosity(0)
 +
 +-- This example calculates the PES and IPES of a H2 molecule in the Hubbard model
 +-- approximation
 +
 +-- we have two s-shells one centered around atom A, one around atom B
 +-- in total there are thus 4 spin-orbitals
 + 
 +NF=4
 +NB=0
 +
 +IndexADn={0}
 +IndexAUp={1}
 +IndexBDn={2}
 +IndexBUp={3}
 +
 +-- The Hamiltonian is given by a hopping from site A to B that conserves spin
 +-- A Coulomb repulsion if the two electrons are on the same site
 +-- A chemical potential that gives an onsite energy of -U/2
 +
 +OppNAdn = NewOperator("Number", NF, 0, 0)
 +OppNAup = NewOperator("Number", NF, 1, 1)
 +OppNBdn = NewOperator("Number", NF, 2, 2)
 +OppNBup = NewOperator("Number", NF, 3, 3)
 +
 +OppNA = OppNAup + OppNAdn
 +OppNB = OppNBup + OppNBdn
 +
 +Oppt = NewOperator("Number", NF, {0,1,2,3}, {2,3,0,1}, {1,1,1,1})
 +OppU = OppNAup * OppNAdn + OppNBup * OppNBdn
 +
 +-- PES and IPES are electron annihilation and creation operators
 +
 +TPes  = NewOperator("An", NF, 0)
 +TIPes = NewOperator("Cr", NF, 0)
 +
 +-- the total number of states is 6.
 +
 +Npsi=6;
 +StartRestrictions = {NF, NB, {"1111",2,2}};
 +psiList = Eigensystem(Oppt, StartRestrictions, Npsi);
 +
 +-- a loop creates U dependent spectra
 +
 +for i=0, 20 do
 +  U = 0.5*i
 +  Hamiltonian=Oppt + U * OppU - U/2 * (OppNA+OppNB)
 +  Eigensystem(Hamiltonian, psiList)
 +  PESIPES= CreateSpectra(Hamiltonian, {TPes, TIPes}, psiList[1], {{"Emin",-1}, {"Emax",9}, {"NE",1000}, {"Gamma",0.25}})
 +  PESIPES.Print({{"file", string.format("PESIPES%i.dat",i)}})
 +end
 +
 +gnuplotInput = [[
 +set autoscale 
 +set xtic auto 
 +set ytic auto 
 +set style line  1 lt 1 lw 1 lc rgb "#000000"
 +set style line  2 lt 1 lw 1 lc rgb "#FF0000"
 +
 +set xlabel "E (eV)" font "Times,12"
 +set ylabel "Intensity (arb. units)" font "Times,12"
 +
 +set out 'H2.ps'
 +set size 1.0, 1.0
 +set terminal postscript portrait enhanced color  "Times" 8
 +
 +
 +plot for [i=0:20] 'PESIPES'.i.'.dat' using ( $1):(-(column(3))+i/2.) notitle with lines ls 1,\
 +     for [i=0:20] 'PESIPES'.i.'.dat' using (-$1):(-(column(5))+i/2.) notitle with lines ls 2
 +]]
 +
 +-- write the gnuplot script to a file
 +file = io.open("H2.gnuplot", "w")
 +file:write(gnuplotInput)
 +file:close()
 +
 +-- call gnuplot to execute the script
 +os.execute("gnuplot H2.gnuplot ; ps2pdf H2.ps")
 +</code>
 +###
 +
 +###
 +The example produces a picture of the spectral function as output:
 +
 +{{:documentation:tutorials:model_examples_from_physics:h2.png?nolink |}}
 +###
 +===== Table of contents =====
 +{{indexmenu>.#1|msort}}
Print/export