DWA-01 Revisiting MC05
Quantum phase transitions in the Bose-Hubbard model
Attention: this implementation of the directed worm algorithm is deprecated and will be removed in a future version of ALPS. It should currently only be used for Bose-Hubbard models with on-site interactions only. Vistrails input files are not available.
As an example of the directed worm algorithm QMC code, we will study a quantum phase transition in the Bose-Hubbard model.
Superfluid density in the Bose-Hubbard model
Preparing and running the simulation from the command line
We create the parameter file parm1a
to set up Monte Carlo simulations of the Bose-Hubbard model on a square lattice with 4x4 sites for a couple of hopping parameters ($t=0.01, 0.02, …, 0.1$) using the dwa code. It is the same parameter file that was used for the ALPS worm MCO5 tutorial.
LATTICE="square lattice";
L=4;
MODEL="boson Hubbard";
Nmax = 2;
U = 1.0;
mu = 0.5;
T = 0.1;
SWEEPS=20000000;
THERMALIZATION=100000;
SKIP=500;
MEASURE[Winding Number]=1
{ t=0.01; }
{ t=0.02; }
{ t=0.03; }
{ t=0.04; }
{ t=0.05; }
{ t=0.06; }
{ t=0.07; }
{ t=0.08; }
{ t=0.09; }
{ t=0.1; }
Using the standard sequence of commands we can run the simulation using the quantum dwa
code:
parameter2xml parm1a
dwa --Tmin 5 --write-xml parm1a.in.xml
We may also run the code using the usual Python methods.
Evaluating the simulation and preparing plots using Python
We load the results from all output files starting with parm1a
and collect the magnetization density as a function of magnetic field.
import pyalps
import matplotlib.pyplot as plt
import pyalps.plot as aplt
data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm1a'),'Stiffness')
rhos = pyalps.collectXY(data,x='t',y='Stiffness')
plt.figure()
aplt.plot(rhos)
plt.xlabel('Hopping $t/U$')
plt.ylabel('Superfluid density $\\rho _s$')
plt.show()
The transition from the Mott insulator to the superfluid
We next want obtain a better location of the phase transition. For this we simulate a two-dimensional square lattice for various system sizes and look for a crossing of the quantity $\rho_s L$.
Preparing and running the simulation from the command line
In the parameter file parm1b
, we focus on the region around the critical point for three system sizes $L=4,6$, and 8:
LATTICE="square lattice";
MODEL="boson Hubbard";
Nmax =2;
U = 1.0;
mu = 0.5;
T = 0.1;
SWEEPS=2000000;
THERMALIZATION=150000;
SKIP=500;
MEASURE[Winding Number]=1;
{ L=4; t=0.01; }
{ L=4; t=0.02; }
{ L=4; t=0.03; }
{ L=4; t=0.04; }
{ L=4; t=0.05; }
{ L=4; t=0.06; }
{ L=4; t=0.07; }
{ L=4; t=0.08; }
{ L=4; t=0.09; }
{ L=4; t=0.1; }
{ L=6; t=0.01; }
{ L=6; t=0.02; }
{ L=6; t=0.03; }
{ L=6; t=0.04; }
{ L=6; t=0.05; }
{ L=6; t=0.06; }
{ L=6; t=0.07; }
{ L=6; t=0.08; }
{ L=6; t=0.09; }
{ L=6; t=0.1; }
{ L=8; t=0.01; }
{ L=8; t=0.02; }
{ L=8; t=0.03; }
{ L=8; t=0.04; }
{ L=8; t=0.05; }
{ L=8; t=0.06; }
{ L=8; t=0.07; }
{ L=8; t=0.08; }
{ L=8; t=0.09; }
{ L=8; t=0.1; }
Evaluating the simulation and preparing plots using Python
We load the results from all output files starting with parm1b
and collect the scaled stiffness as a function of the hopping parameter.
import pyalps
import matplotlib.pyplot as plt
import pyalps.plot as aplt
data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm1b'),'Stiffness')
rhos = pyalps.collectXY(data,x='t',y='Stiffness',foreach=['L'])
for rho in rhos:
rho.y = rho.y * float(rho.props['L'])
plt.figure()
aplt.plot(rhos)
plt.xlabel('Hopping $t/U$')
plt.ylabel('$\\rho _sL$')
plt.legend()
plt.title('Scaling plot for Bose-Hubbard model')
plt.show()
Note that the legend and labels that are nicely set up.
For the experts: The quantum phase transition between the superfluid and the Mott insulating phase can have two different university classes. First, if we work at constant density, then the transition will have dynamical exponent $z=1$ and can be described by the 3d XY model (more precisely, an emergent (2+1)D CFT). For the simulations, this implies that we must scale the inverse temperature linearly with system size and make sure that we work in the canonical ensemble or enforce that the average density is 1, $\langle n \rangle = 1$, which also works in this case. Second, if we drive the transition by changing the density (as in a drive with chemical potential), then the dynamical exponent $z=2$ and the transition is mean-fieldish describing a few particles (holes) on top of a (rescaled) vacuum. Note that in the tutorial we do neither of these protocols and our results are therefore only approximately locating the phase transition. The scaling form of the stiffness used in the tutorial is the one of the 3D XY model. We refer to this paper by B. Capogrosso-Sansone et al. for a detailed quantum Monte Carlo study of the 2D Bose-Hubbard model.
Contributors
- Matthias Troyer
- Ping Nang Ma
Revisions
- Lode Pollet