Skip to content

Commit ca10307

Browse files
OoolabOoolab
authored
Merge pull request #57 from ALPSim/fix/mc05-rewrite
Rewrite MC-05 bosons tutorial
2 parents 1c64727 + 9d19aee commit ca10307

1 file changed

Lines changed: 167 additions & 93 deletions

File tree

‎content/en/tutorials/mcs/mc05.md‎

Lines changed: 167 additions & 93 deletions
Original file line numberDiff line numberDiff line change
@@ -6,140 +6,214 @@ toc: true
66
weight: 7
77
---
88

9-
## Quantum phase transitions in the Bose-Hubbard model
9+
The Bose-Hubbard model describes interacting bosons on a lattice:
1010

11-
As an example of the worm QMC code, we will study a quantum phase transition in the Bose-Hubbard mode.
11+
$$H = -t \sum_{\langle i,j \rangle} (a_i^\dagger a_j + \text{h.c.}) + \frac{U}{2} \sum_i n_i(n_i-1) - \mu \sum_i n_i,$$
1212

13-
### Superfluid density in the Bose Hubbard model
13+
where $t$ is the hopping amplitude, $U$ the on-site repulsion, and $\mu$ the chemical potential.
14+
At integer filling and large $U/t$ the system is a **Mott insulator**: bosons are localized by interactions and the superfluid density $\rho_s = 0$.
15+
As $t/U$ increases, quantum fluctuations eventually drive a transition to a **superfluid** phase with $\rho_s > 0$.
16+
This tutorial uses the ALPS worm QMC code to locate this quantum phase transition on a two-dimensional square lattice at filling $\langle n \rangle = 1$ (set by $\mu = U/2 = 0.5$).
1417

15-
#### Preparing and running the simulation from the command line
18+
## Superfluid density across the transition
1619

17-
The parameter file <a href="https://github.com/ALPSim/ALPS/blob/master/tutorials/mc-05-bosons/parm5a" download>`parm5a`</a> with the following contents sets up Monte Carlo simulations of the quantum 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 worm code.
20+
We first scan a wide range of hopping values on a $4 \times 4$ lattice to observe how the superfluid density $\rho_s$ (called "Stiffness" in ALPS) evolves across the transition.
21+
The Hilbert space is truncated at `Nmax=2` bosons per site, which is a good approximation near the Mott lobe at unit filling.
1822

19-
```
20-
LATTICE="square lattice";
21-
L=4;
22-
MODEL="boson Hubbard";
23-
NONLOCAL=0;
24-
U = 1.0;
25-
mu = 0.5;
26-
Nmax = 2;
27-
T = 0.1;
28-
SWEEPS=500000;
29-
THERMALIZATION=10000;
30-
{ t=0.01; }
31-
{ t=0.02; }
32-
{ t=0.03; }
33-
{ t=0.04; }
34-
{ t=0.05; }
35-
{ t=0.06; }
36-
{ t=0.07; }
37-
{ t=0.08; }
38-
{ t=0.09; }
39-
{ t=0.1; }
40-
```
23+
#### Setting up and running on the command line
4124

42-
The corresponding Python script is found at <a href="https://github.com/ALPSim/ALPS/blob/master/tutorials/mc-05-bosons/tutorial5a.py" download>`tutorial5a.py`</a>.
25+
The parameter file <a href="https://github.com/ALPSim/ALPS/blob/master/tutorials/mc-05-bosons/parm5a" download>`parm5a`</a>:
4326

44-
#### Evaluating the simulation and preparing plots using Python
27+
```
28+
LATTICE="square lattice"
29+
L=4
30+
MODEL="boson Hubbard"
31+
NONLOCAL=0
32+
U=1.0
33+
mu=0.5
34+
Nmax=2
35+
T=0.1
36+
SWEEPS=500000
37+
THERMALIZATION=10000
38+
{t=0.01;}
39+
{t=0.02;}
40+
{t=0.03;}
41+
{t=0.04;}
42+
{t=0.05;}
43+
{t=0.06;}
44+
{t=0.07;}
45+
{t=0.08;}
46+
{t=0.09;}
47+
{t=0.1;}
48+
```
4549

46-
To load the results and prepare plots we load the results from the output files and collect the magntization density as a function of magnetic field from all output files starting with `parm5a`.
50+
`NONLOCAL=0` disables non-local measurements to keep the output compact.
4751

4852
```
49-
data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm5a'),'Stiffness')
50-
magnetization = pyalps.collectXY(data,x='h',y='Stiffness')
53+
parameter2xml parm5a
54+
worm --Tmin 10 --write-xml parm5a.in.xml
5155
```
5256

53-
To make plots we call the pyalps.plot.plot and then set some nice labels, a title, and a range of y-values:
54-
57+
#### Setting up and running in Python
58+
59+
The script <a href="https://github.com/ALPSim/ALPS/blob/master/tutorials/mc-05-bosons/tutorial5a.py" download>`tutorial5a.py`</a>:
60+
61+
```Python
62+
import pyalps
63+
import matplotlib.pyplot as plt
64+
import pyalps.plot
65+
66+
parms = []
67+
for t in [0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1]:
68+
parms.append(
69+
{
70+
'LATTICE' : "square lattice",
71+
'L' : 4,
72+
'MODEL' : "boson Hubbard",
73+
'NONLOCAL' : 0,
74+
'U' : 1.0,
75+
'mu' : 0.5,
76+
'Nmax' : 2,
77+
'T' : 0.1,
78+
'SWEEPS' : 500000,
79+
'THERMALIZATION' : 10000,
80+
't' : t
81+
}
82+
)
83+
84+
input_file = pyalps.writeInputFiles('parm5a', parms)
85+
pyalps.runApplication('worm', input_file, Tmin=5)
5586
```
87+
88+
#### Evaluating and plotting
89+
90+
Load the superfluid stiffness and plot it as a function of hopping:
91+
92+
```Python
93+
data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm5a'), 'Stiffness')
94+
rhos = pyalps.collectXY(data, x='t', y='Stiffness')
95+
5696
plt.figure()
5797
pyalps.plot.plot(rhos)
5898
plt.xlabel('Hopping $t/U$')
59-
plt.ylabel('Superfluid density $\\rho _s$')
99+
plt.ylabel('Superfluid density $\\rho_s$')
100+
plt.title('Bose-Hubbard model on a $4\\times 4$ lattice')
60101
plt.show()
61102
```
62103

63-
#### Questions
104+
$\rho_s$ should be small (consistent with zero) for small $t/U$ and grow to a finite value for large $t/U$, with a crossover near the critical hopping $(t/U)_c \approx 0.060$.
64105

65-
What is the signature of the phase transition?
106+
## Locating the critical point
66107

67-
### The transition from the Mott insulator to the superfluid
108+
To pin down $(t/U)_c$ more precisely we exploit finite-size scaling.
109+
At the quantum critical point, $\rho_s \sim L^{-(d+z-2)}$ where $d=2$ is the dimension and $z=1$ is the dynamical exponent for this universality class (3D XY).
110+
For $d=2$, $z=1$ this gives $\rho_s \sim L^{-1}$, so the combination $\rho_s L$ is dimensionless at criticality and curves for different $L$ cross at $t_c$.
68111

69-
We next want to pin down the location of the phase transition more accurately. For this we simulate a two-dimensional square lattice for various system sizes and look for a crossing of the quantity $\rho_s L$.
112+
We simulate three system sizes $L = 4, 6, 8$ on a fine grid of hopping values around the expected critical point.
70113

71-
#### Preparing and running the simulation from the command line
114+
#### Setting up and running on the command line
72115

73-
In the parameter file <a href="https://github.com/ALPSim/ALPS/blob/master/tutorials/mc-05-bosons/parm5b" download>`parm5b`</a> we focus on the region around the critical point for three system sizes L=4, 6, and 8:
116+
The parameter file <a href="https://github.com/ALPSim/ALPS/blob/master/tutorials/mc-05-bosons/parm5b" download>`parm5b`</a>:
74117

75118
```
76-
LATTICE="square lattice";
77-
MODEL="boson Hubbard";
78-
NONLOCAL=0;
79-
U = 1.0;
80-
mu = 0.5;
81-
Nmax = 2;
82-
T = 0.05;
83-
SWEEPS=600000;
84-
THERMALIZATION=150000;
85-
{ L=4; t=0.045; }
86-
{ L=4; t=0.05; }
87-
{ L=4; t=0.0525; }
88-
{ L=4; t=0.055; }
89-
{ L=4; t=0.0575; }
90-
{ L=4; t=0.06; }
91-
{ L=4; t=0.065; }
92-
{ L=6; t=0.045; }
93-
{ L=6; t=0.05; }
94-
{ L=6; t=0.0525; }
95-
{ L=6; t=0.055; }
96-
{ L=6; t=0.0575; }
97-
{ L=6; t=0.06; }
98-
{ L=6; t=0.065; }
99-
{ L=8; t=0.045; }
100-
{ L=8; t=0.05; }
101-
{ L=8; t=0.0525; }
102-
{ L=8; t=0.055; }
103-
{ L=8; t=0.0575; }
104-
{ L=8; t=0.06; }
105-
{ L=8; t=0.065; }
119+
LATTICE="square lattice"
120+
MODEL="boson Hubbard"
121+
NONLOCAL=0
122+
U=1.0
123+
mu=0.5
124+
Nmax=2
125+
T=0.05
126+
SWEEPS=600000
127+
THERMALIZATION=150000
128+
{L=4; t=0.045;}
129+
{L=4; t=0.05;}
130+
{L=4; t=0.0525;}
131+
{L=4; t=0.055;}
132+
{L=4; t=0.0575;}
133+
{L=4; t=0.06;}
134+
{L=4; t=0.065;}
135+
{L=6; t=0.045;}
136+
{L=6; t=0.05;}
137+
{L=6; t=0.0525;}
138+
{L=6; t=0.055;}
139+
{L=6; t=0.0575;}
140+
{L=6; t=0.06;}
141+
{L=6; t=0.065;}
142+
{L=8; t=0.045;}
143+
{L=8; t=0.05;}
144+
{L=8; t=0.0525;}
145+
{L=8; t=0.055;}
146+
{L=8; t=0.0575;}
147+
{L=8; t=0.06;}
148+
{L=8; t=0.065;}
106149
```
107-
108-
The corresponding Python script is found at <a href="https://github.com/ALPSim/ALPS/blob/master/tutorials/mc-05-bosons/tutorial5b.py" download>`tutorial5b.py`</a>.
109-
110-
#### Evaluating the simulation using Python
111150

112-
We first load the superfluid density (stiffness) into three different data sets, one for each system size L:
151+
The lower temperature ($T = 0.05$) and longer runs compared to `parm5a` are needed to resolve the crossing clearly.
113152

114153
```
115-
data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm5b'),'Stiffness')
116-
rhos = pyalps.collectXY(data,x='t',y='Stiffness',foreach=['L'])
154+
parameter2xml parm5b
155+
worm --Tmin 10 --write-xml parm5b.in.xml
117156
```
118157

119-
Next we multiply each data set by the size L:
120-
158+
#### Setting up and running in Python
159+
160+
The script <a href="https://github.com/ALPSim/ALPS/blob/master/tutorials/mc-05-bosons/tutorial5b.py" download>`tutorial5b.py`</a> adapts `tutorial5a.py`: rename the prefix to `parm5b`, lower `T` to 0.05, increase `THERMALIZATION` to 150000 and `SWEEPS` to 600000, and loop over both `L` and `t`:
161+
162+
```Python
163+
import pyalps
164+
import matplotlib.pyplot as plt
165+
import pyalps.plot
166+
167+
parms = []
168+
for L in [4, 6, 8]:
169+
for t in [0.045, 0.05, 0.0525, 0.055, 0.0575, 0.06, 0.065]:
170+
parms.append(
171+
{
172+
'LATTICE' : "square lattice",
173+
'L' : L,
174+
'MODEL' : "boson Hubbard",
175+
'NONLOCAL' : 0,
176+
'U' : 1.0,
177+
'mu' : 0.5,
178+
'Nmax' : 2,
179+
'T' : 0.05,
180+
'SWEEPS' : 600000,
181+
'THERMALIZATION' : 150000,
182+
't' : t
183+
}
184+
)
185+
186+
input_file = pyalps.writeInputFiles('parm5b', parms)
187+
pyalps.runApplication('worm', input_file, Tmin=5)
121188
```
189+
190+
#### Evaluating and plotting
191+
192+
Load the stiffness for each system size, multiply by $L$, and plot:
193+
194+
```Python
195+
data = pyalps.loadMeasurements(pyalps.getResultFiles(prefix='parm5b'), 'Stiffness')
196+
rhos = pyalps.collectXY(data, x='t', y='Stiffness', foreach=['L'])
197+
122198
for s in rhos:
123199
s.y = s.y * float(s.props['L'])
124-
```
125-
126-
And finally we make a plot in the usual way:
127200

128-
```
129201
plt.figure()
130-
pyalps.pyplot.plot(rhos)
202+
pyalps.plot.plot(rhos)
131203
plt.xlabel('Hopping $t/U$')
132-
plt.ylabel('$\\rho _sL$')
204+
plt.ylabel('$\\rho_s L$')
133205
plt.legend()
134-
plt.title('Scaling plot for Bose-Hubbard model')
206+
plt.title('Finite-size scaling: Bose-Hubbard model')
135207
plt.show()
136208
```
137209

138-
Note the legend and labels that are nicely set up.
210+
The curves for $L = 4, 6, 8$ should cross near $(t/U)_c$.
211+
The exact result for the 2D Bose-Hubbard model at unit filling is $(t/U)_c = 0.05974\ldots$
139212

140-
#### Questions
213+
## Questions
141214

142-
- How can you determine the location of the quantum phase transition in the thermodynamic limit?
143-
- Tip: Multiply your results for the superfluid stiffness by the respective linear system size L.
144-
- Compare your result to the exact result (t/U)c = 0.05974...
145-
- Why does the Monte Carlo simulation overestimate the critical point of the transition?
215+
- At what hopping value does $\rho_s$ first become clearly nonzero in the coarse scan? Does this agree with the crossing in the $\rho_s L$ plot?
216+
- Where do the $\rho_s L$ curves cross? How does your estimate compare with the exact value $(t/U)_c = 0.05974$?
217+
- The finite-temperature simulations systematically overestimate $(t/U)_c$. Why? What would happen to the crossing point if you lowered $T$ further?
218+
- Try increasing `Nmax` to 3. How much do the results change? At what filling or interaction strength would `Nmax=2` become a poor approximation?
219+
- (Bonus) Repeat the finite-size scaling analysis with larger system sizes. How does the quality of the crossing improve?

0 commit comments

Comments
 (0)