ALPS 2 Tutorials:MC-06 QWL

From ALPS
Jump to: navigation, search


In this tutorial we will introduce QMC simulations using the quantum version of the Wang-Landau algorithm.

Thermodynamics of quantum Heisenberg spin chains

The ferromagnetic Heisenberg chain

We will start with a quick tutorial on using the ALPS qwl code for a spin chain.

Using the command line

The parameter file parm6a sets up a Monte Carlo simulation of the quantum mechanical Heisenberg ferromagnet on a one-dimensional chain with 40 sites, using the quantum Wang-Landau (QWL) method.

LATTICE="chain lattice" 
MODEL="spin"
local_S   = 1/2
L       = 40
CUTOFF  = 500
{J = -1}

After preparing the input files and running the qwl code using the standard commands:

parameter2xml parm6a
qwl  parm6a.in.xml

You can produce XML plot files for the thermodynamic and magnetic observables using qwl_evaluate:

qwl_evaluate --T_MIN 0.1 --T_MAX 10 --DELTA_T  0.1 parm6a.task1.out.xml

This will generate the following XML plot files:

parm6a.task1.plot.energy.xml
parm6a.task1.plot.free_energy.xml
parm6a.task1.plot.entropy.xml
parm6a.task1.plot.specific_heat.xml
parm6a.task1.plot.uniform_structure_factor.xml
parm6a.task1.plot.staggered_structure_factor.xml
parm6a.task1.plot.uniform_susceptibility.xml

To extract the calculated results from the XML plot files generated by qwl_evaluate, you can use the plot2text tool, and then view this data with your favorite plotting tool. For example, to extact the data of the energy density vs. temperature, use

plot2text parm6a.task1.plot.energy.xml

In a similar way, you can extract the data from the other XML plot files.

When Grace is your favorite plotting tool, you can also directly generate a Grace project file from the XML plot file using the plot2xmgr tool. For example, to generate a Grace project file of the energy vs. temperature, use

plot2xmgr parm6a.task1.plot.energy.xml > energy.agr

Similarly the tool plot2gp produces Gnuplot scripts and plot2text converts the file to plain text. However, the preferred method for data evaluation and plotting is using Python.

Using Python

To set up and run the simulation in Python we use the script tutorial6a.py. The first parts of this script imports the required modules and then prepares the input files as a list of Python dictionaries, and then runs the simulation

import pyalps
import matplotlib.pyplot as plt
import pyalps.pyplot
 parms = [{ 
         'LATTICE'        : "chain lattice", 
         'MODEL'          : "spin",
         'local_S'        : 0.5,
         'L'              : 40,
         'J'              : -1 ,
         'CUTOFF'         : 1000
       }]
input_file = pyalps.writeInputFiles('parm6a',parms)
res = pyalps.runApplication('qwl',input_file)

We next run the evaluation program on all output files

data = pyalps.evaluateQWL(pyalps.getResultFiles(prefix='parm6a'),DELTA_T=0.1, T_MIN=0.1, T_MAX=10.0)

and finally show all the plots:

for s in pyalps.flatten(data):
  plt.figure()
  plt.title("Ferromagnetic Heisenberg chain")
  pyalps.pyplot.plot(s)
plt.show()

Using Vistrails

To run the simulation in Vistrails open the file mc-06-qwl.vt and look at the workflow labeled "Ferromagnetic chain". Click on "Execute" to prepare the input file, run the simulation and create the output figures similar to the ones below:

vt_id:12 version:138vt_id:12 version:138vt_id:12 version:138vt_id:12 version:138vt_id:12 version:138

The antiferromagnetic Heisenberg chain

To simulate the antiferromagnetic chain prepare new simulations setting J=1 instead of J=-1. The parameters are in parm6b, the Python script in tutorial6b.py and the Vistrails workflow is in mc-06-qwl.vt.

Questions

  • Where are differences between the two cases most pronounced?
  • Why are there only minor differences at high temperatures?
  • What is the value of the entropy at zero and for infinite temperature in both cases (if not sure, perform a simulation for a 8 sites chain to obtain further data)?
  • How is this compatible with the third law of thermodynamics?
  • Why does the uniform susceptibility behave so differently in the two cases?

The three-dimensional Heisenberg antiferromagnet

Simulating the 3D quantum Heisenberg antiferromegnet

The parameter file parm6c sets up a Monte Carlo simulation of the quantum mechanical Heisenberg antiferromagnet on a three-dimensional simple cubic lattice with 43 sites, using the QWL method. The Python script is tutorial6c.py and the Vistrails workflow is in mc-06-qwl.vt.

The simulations are set up and run as above.

Questions

  • Why does the staggered structure factor start to increase near T≈1?
  • What are further indications of this phenomena in the thermodynamic properties?

Finite size scaling analysis to determine the critical point

Finite size scaling theory predics the staggered structure factor S(L) for this transition to scale at the critical point as L2-η, where η≈0.034. A scaling plot of S(L)/L2-η vs. temperature is expected to show a crossing of curves for different linear system sizes L at the critical temperature Tc. In order to produce such a scaling plot, we set up a further simulation of the cubic antiferromagnet, for a larger system with L=4 and a cutoff at 1000, given in the parameter file parm6d

Evaluation now requires multiplication of the results with L2-η which is easiest done in Python or Vistrails.After running the simulation we first load the results:

results = pyalps.evaluateQWL(pyalps.getResultFiles(prefix='parm6d'),DELTA_T=0.05, T_MIN=0.5, T_MAX=1.5)

Next we extract just the staggered structure factor S(Q) for any system size L, rescale it by L^{-2+\eta}, and set a label according to the system size:

data = []
for s in pyalps.flatten(results):
  if s.props['ylabel']=='Staggered Structure Factor per Site':
    print 'yes'
    d = copy.deepcopy(s) # make a deep copy to not change the original
    l = s.props['L']
    d.props['label']='L='+str(l)
    d.y = d.y * pow(float(l),-1.97)
    data.append(d)

And finally we make the plot

plt.figure()
plt.title("Scaling plot for cubic lattice Heisenberg antiferromagnet")
pyalps.pyplot.plot(data)
plt.legend()
plt.xlabel('Temperature $T/J$')
plt.ylabel('$S(\pi,\pi,\pi) L^{-2+\eta}$')
plt.show()

The Vistrails workflow is in mc-06-qwl.vt and will produce a figure like:

vt_id:12 version:143


Questions

  • Do the curves indeed cross?
  • What is your estimated value of the critical temperature? Compare your estimate to Tc=0.946.
  • How could you improve your estimated value?
  • Would you expect the critical temerature for the quantum ferromagnet to be the same?
  • How would you proceed to obtain a guess for its value? (Give it a try!)


© 2004-2010 by Stefan Wessel and Matthias Troyer