Difference between revisions of "ALPS 2 Tutorials:MC-01 Equilibration"

From ALPS
Jump to: navigation, search
(Using Vistrails)
(Example: Quantum Monte Carlo (directed worm algorithm) simulations)
Line 3: Line 3:
 
'''Rule of thumb: All Monte Carlo simulations have to be equilibrated before taking measurements.'''
 
'''Rule of thumb: All Monte Carlo simulations have to be equilibrated before taking measurements.'''
  
== Example: Quantum Monte Carlo (directed worm algorithm) simulations ==
+
== Example: Classical Monte Carlo (local updates) simulations ==
  
 
As an example, we consider a Quantum Monte Carlo simulation implemented in the directed worm algorithm for boson Hubbard model in square lattice geometry of size 20<sup>2</sup>.
 
As an example, we consider a Quantum Monte Carlo simulation implemented in the directed worm algorithm for boson Hubbard model in square lattice geometry of size 20<sup>2</sup>.
Line 11: Line 11:
 
The parameter file [http://alps.comp-phys.org/static/tutorials2.1.0/mc-01b-equilibration/parm1a parm1a]:
 
The parameter file [http://alps.comp-phys.org/static/tutorials2.1.0/mc-01b-equilibration/parm1a parm1a]:
  
  LATTICE="square lattice"
+
  T=2.269186
MODEL="boson Hubbard"
+
  J=1
 
L=20 
 
Nmax=20
 
 
  t=1.
 
U=16.
 
mu=32.
 
 
 
  THERMALIZATION=10000
 
  THERMALIZATION=10000
  SWEEPS=100000
+
  SWEEPS=50000 
  SKIP=400
+
UPDATE="local"
   
+
MODEL="Ising"
  {T=1.0}
+
{L=2;}
 +
{L=4;}
 +
{L=8;}
 +
  {L=16;}
 +
  {L=32;}
 +
  {L=48;}
 +
 
  
 
We first convert the input parameters to XML and then run the application '''dwa''':
 
We first convert the input parameters to XML and then run the application '''dwa''':

Revision as of 12:11, 5 September 2013

Equilibration

Rule of thumb: All Monte Carlo simulations have to be equilibrated before taking measurements.

Example: Classical Monte Carlo (local updates) simulations

As an example, we consider a Quantum Monte Carlo simulation implemented in the directed worm algorithm for boson Hubbard model in square lattice geometry of size 202.

Using command line

The parameter file parm1a:

T=2.269186
J=1
THERMALIZATION=10000
SWEEPS=50000  
UPDATE="local"
MODEL="Ising"
{L=2;}
{L=4;}
{L=8;}
{L=16;}
{L=32;}
{L=48;}


We first convert the input parameters to XML and then run the application dwa:

parameter2xml parm1a
dwa parm1a.in.xml

Detailed information regarding the Density measurements, for example, can be extracted from:

h5dump -g /simulation/results/Density parm1a.task1.out.h5

or its (binned) timeseries from:

h5dump -g /simulation/results/Density/timeseries parm1a.task1.out.h5

We can extract the timeseries data into a CSV file:

h5dump -d /simulation/results/Density/timeseries/data -w 1 -y -o parm1a.task1.out.density.timeseries.csv parm1a.task1.out.h5

and plot it using your favourite graphical plotter, say xmgrace:

xmgrace parm1a.task1.out.density.timeseries.csv

or, say gnuplot:

gnuplot <<@@
set datafile separator ','
plot "parm1a.task1.out.density.timeseries.csv"
@@

Based on the timeseries, the user will then judge for himself/herself whether the simulation has reached equilibration.

Using Python

The following describes what is going on within the script file tutorial1a.py.

The headers:

import pyalps

Set up a python list of parameters (python) dictionaries:

parms = [{
  'LATTICE'         : "square lattice",          
  'MODEL'           : "boson Hubbard",
  'L'               : 20,
  'Nmax'            : 20,
  't'               : 1.,
  'U'               : 16.,
  'mu'              : 32.,
  'T'               : 1.,
  'THERMALIZATION'  : 10000,
  'SWEEPS'          : 100000,
  'SKIP'            : 400
}]

Write into XML input file:

input_file = pyalps.writeInputFiles('parm1a',parms)

and run the application dwa:

pyalps.runApplication('dwa', input_file, Tmin=10, writexml=True)

We first get the list of all hdf5 result files via:

files = pyalps.getResultFiles(prefix='parm1a', format='hdf5')

and then extract, say the timeseries of the Density measurements:

ar = pyalps.hdf5.h5ar(files[0])
density_timeseries = ar['/simulation/results']['Density']['timeseries']['data']

We can then visualize graphically:

import matplotlib.pyplot as plt
plt.plot(density_timeseries)
plt.show()

Based on the timeseries, the user will then judge for himself/herself whether the simulation has reached equilibration.

A convenient tool: steady_state_check

ALPS Python provides a convenient tool to check whether a measurement observable(s) has (have) reached steady state equilibrium.

Here is one example:

pyalps.steady_state_check(files[0], 'Density')

and another one:

pyalps.steady_state_check(files[0], ['Density', 'Energy Density'])

Description
1. steady_state_check first performs a linear fit on the timeseries, and decides whether the measurement observable has reached steady state equilibrium based on the gradient/slope of the fitted line.

2. The optional arguments of steady_state_check are:

argument default remark
tolerance 0.01  \mathrm{tolerance} = \frac{X^\mathrm{(fit)} (t_\mathrm{final}) - X^\mathrm{(fit)} (t_\mathrm{initial})}{\bar{X}}
simplified False shall we combine the checks of all observables as 1 final boolean answer?
includeLog False shall we print the detailed log?

3. To see the complete log for instance:

pyalps.steady_state_check(files[0], ['Density', 'Energy Density'], includeLog=True)

Using Vistrails

To run the simulation in Vistrails open the file mc-01b-equilibration.vt.