Difference between revisions of "Documentation:Bosons in an optical lattice"

From ALPS
Jump to: navigation, search
(→‎Momentum distribution and time-of-flight images)
(→‎Bosons in an optical lattice trap)
Line 159: Line 159:
 
Here, the units will be cleverly normalized later on.
 
Here, the units will be cleverly normalized later on.
  
−
=== Momentum distribution ===
 
  
−
The momentum distribution <math>\langle n(\vec{k})\rangle</math> of bosons in an optical lattice at equilbrium can be expressed as:
+
(c)
−
<math>
 
−
\langle \hat{\psi}^+ (\vec{k}) \,  \hat{\psi} (\vec{k}) \rangle  = \frac{1}{N} \int d\vec{r} \, d\vec{r}' \langle \hat{\psi}^+ (\vec{r}) \,  \hat{\psi} (\vec{r}') \rangle \, e^{i \vec{k} \cdot (\vec{r} - \vec{r}')} = \frac{1}{N} \sum_{i,j} \int d\vec{r} \, d\vec{r}' \, w(\vec{r} - \vec{r}_i) e^{i\vec{k}\cdot(\vec{r} - \vec{r}_i)} \, w(\vec{r}' - \vec{r}_j) e^{-i\vec{k}\cdot(\vec{r}' - \vec{r}_j)}    \langle \hat{b}^+_i \hat{b}_j \rangle \, e^{i \vec{k} \cdot (\vec{r}_i - \vec{r}_j)}
 
−
</math>
 
−
 
 
−
Expressing the fourier transform of the wannier function <math>w(\vec{r}) </math> as
 
−
<blockquote>
 
−
<math>
 
−
\tilde{w}(\vec{k}) = \frac{1}{\sqrt{N}} \int d\vec{r} \, w(\vec{r}) e^{-i\vec{k}\cdot\vec{r}}  \,\, ,
 
−
</math>
 
−
</blockquote>
 
−
and the interference term:
 
−
<blockquote>
 
−
<math>
 
−
S(\vec{k}) = \sum_{i,j} \langle \hat{b}^+_i \hat{b}_j \rangle \, e^{i \vec{k} \cdot (\vec{r}_i - \vec{r}_j)} \,\, ,
 
−
</math>
 
−
</blockquote>
 
−
the momentum distribution reduces to the following form:
 
−
<blockquote>
 
−
<math>
 
−
\langle n(\vec{k})\rangle = \left| \tilde{w}(\vec{k})\right|^2 S(\vec{k})  \,\, .
 
−
</math>
 
−
</blockquote>
 
−
 
 
−
=== Time of flight images (in experiment) ===
 
−
 
 
−
 
 
−
To measure it in experiments, the optical lattice is momentarily switched off and the bosons expand freely by assumption with momentum gained from the lattice momentum previously.
 
−
Measurements are performed in our classical world, and therefore semiclassical treatment is already sufficient, i.e.
 
−
<blockquote>
 
−
<math>
 
−
\hbar \vec{k} = m \left( \frac{\vec{r} }{t_f} \right)
 
−
</math>
 
−
</blockquote>
 
−
where <math>t_f</math> is the time-of-flight taken by the bosons to move from the origin (experiment) to the detector probe at position <math>\vec{r} </math>.
 
−
Here, we assume the simplest picture, ie. 1) no interaction, 2) no collision, and 3) uniform motion.
 
−
 
 
−
 
 
−
Of course, at equilibrium, i.e. "after the bosons have travelled for a long time", the time-of-flight image will capture a density of
 
−
<blockquote>
 
−
<math>
 
−
\langle n_\mathrm{tof} (\vec{r} ) \rangle = \left( \frac{m}{\hbar t_f} \right)^3  \left| \tilde{w}\left( \frac{m}{\hbar t} \vec{r} \right) \right|^2 S(\vec{k})
 
−
</math>
 
−
</blockquote>
 
−
which indirectly probes the original momentum distribution of the lattice bosons.
 
−
This form is only correct in the far-field limit or in the limit of infinite time-of-flight (tof).
 
−
 
 
−
 
 
−
To improve upon the accuracy, we shall correct it with semiclassical dynamics during the time of flight.
 
−
For the bosons that originate from lattice site at <math>\vec{r}_j</math> to the detector probe at <math>\vec{r}</math>, their time-dependent wavefunction must be:
 
−
<blockquote>
 
−
<math>
 
−
\Psi(\vec{r},t) = \left( \frac{m}{\hbar t_f} \right)^{3/2} \tilde{w}\left( \frac{m}{\hbar t_f} (\vec{r} - \vec{r}_j) \right)  \, \exp \left( \frac{i \epsilon_k t_f}{\hbar} \right)
 
−
</math>
 
−
</blockquote>
 
−
where the kinetic energy of the free bosons must be:
 
−
<blockquote>
 
−
<math>
 
−
\epsilon_k = \frac{\hbar^2 k^2}{2m} = \frac{\hbar^2}{2m}  \left( \frac{m}{\hbar t_f} \right)^2 (\vec{r} - \vec{r}_j )^2
 
−
</math>
 
−
</blockquote>
 
−
Therefore, the corrected time-of-flight image must be:
 
−
<blockquote>
 
−
<math>
 
−
\langle n_\mathrm{tof} (\vec{r} ) \rangle \approx \left( \frac{m}{\hbar t_f} \right)^3  \left| \tilde{w}\left( \frac{m}{\hbar t} \vec{r} \right) \right|^2  \sum_{i,j} \, \langle \hat{b}^+_i \hat{b}_j \rangle \exp \left( - \frac{i}{\hbar} \frac{\hbar^2}{2m}  \left( \frac{m}{\hbar t_f} \right)^2 (\vec{r} - \vec{r}_i )^2 t_f \right) \, \exp \left(  \frac{i}{\hbar} \frac{\hbar^2}{2m}  \left( \frac{m}{\hbar t_f} \right)^2 (\vec{r} - \vec{r}_j )^2 t_f \right)
 
−
</math>
 
−
</blockquote>
 
−
where the dependence of the initial site on the wannier enveloped function is neglected.
 
−
Finally, we arrive at:
 
−
<blockquote>
 
−
<math>
 
−
\langle n_\mathrm{tof} (\vec{k} ) \rangle \approx  \left| \tilde{w} ( \vec{k} ) \right|^2  \sum_{i,j} \, \langle \hat{b}^+_i \hat{b}_j \rangle \, e^{i \vec{k} \cdot (\vec{r}_i - \vec{r}_j) - i \gamma_f (\vec{r}_i^2 - \vec{r}_j^2) }  \,\, .
 
−
</math>
 
−
</blockquote>
 
−
with time-of-flight phase <math> \gamma_f = \frac{m}{2\hbar t_f} </math>.
 
−
 
 
−
 
 
−
For example, the time-of-flight phase <math>\gamma</math> for Rb-87 bosons (mass 86.99 a.m.u.) with a time-of-flight of 15.5 ms
 
−
 
 
−
a. expanded from an isotropic optical lattice of wavelength 765 nm:
 
−
 
 
−
>>> pyalps.dwa.tofPhase(time_of_flight=15.5, wavelength=765, mass=86.99)
 
−
0.006464602556863682
 
−
 
 
−
b. expanded from an optical lattice of wavelength 843 nm along x-direction and 765 nm along y- and z- directions:
 
−
 
 
−
>>> pyalps.dwa.tofPhase(time_of_flight=15.5, wavelength=[843,765,765], mass=86.99)
 
−
[0.007850080468935228, 0.006464602556863682, 0.006464602556863682]
 
−
 
 
−
 
 
−
The following table summarizes the basic option to turn on time-of-flight imaging options.
 
−
 
 
−
<div id="links_option3"></div>
 
−
{| border="1"
 
−
|- style="background:gray; color:white"
 
−
! &nbsp; Option &nbsp;
 
−
! &nbsp; Default &nbsp;
 
−
! &nbsp; Remark &nbsp;
 
−
|-
 
−
| &nbsp; tof_phase &nbsp;
 
−
| align="center" | 0.
 
−
| &nbsp; time-of-flight phase <math>\gamma</math> &nbsp;
 
−
|}
 
−
 
 
−
=== Symmetry leads to further simplification ===
 
−
 
 
−
Symmetries:
 
−
<blockquote>
 
−
<math>
 
−
\langle n_\mathrm{tof} (\vec{k} ) \rangle  = \langle n_\mathrm{tof} (-\vec{k} ) \rangle
 
−
</math>
 
−
</blockquote>
 
−
<blockquote>
 
−
<math>
 
−
\langle \hat{b}^+_i \hat{b}_j \rangle  = \langle \hat{b}^+_j \hat{b}_i \rangle
 
−
</math>
 
−
</blockquote>
 
−
 
 
−
Therefore:
 
−
<blockquote>
 
−
<math>
 
−
\langle n_\mathrm{tof} (\vec{k} ) \rangle \approx  \left| \tilde{w} ( \vec{k} ) \right|^2  \sum_{i,j} \, \langle \hat{b}^+_i \hat{b}_j \rangle \, \cos( \vec{k} \cdot (\vec{r}_i - \vec{r}_j) ) \, \cos(\gamma_f (r_i^2 - r_j^2))  \,\, .
 
−
</math>
 
−
</blockquote>
 
−
 
 
−
Define green function:
 
−
<blockquote>
 
−
<math>
 
−
g_f(i,j) = \langle \hat{b}^+_i \hat{b}_j \rangle \, \cos(\gamma_f (r_i^2 - r_j^2)) 
 
−
</math>
 
−
</blockquote>
 
−
 
 
−
In reduced coordinates, <math> \vec{r}_\alpha = \vec{r}_i - \vec{r}_j </math> , we redefine green function:
 
−
<blockquote>
 
−
<math>
 
−
g_f(\alpha) = \sum_{i,j : \, \vec{r}_i - \vec{r}_j \, = \, \pm \vec{r}_\alpha} \langle \hat{b}^+_i \hat{b}_j \rangle \, \cos(\gamma_f (r_i^2 - r_j^2)) 
 
−
</math>
 
−
</blockquote>
 
−
and therefore the momentum distribution/ time-of-flight image:
 
−
<blockquote>
 
−
<math>
 
−
\langle n_\mathrm{tof} (\vec{k} ) \rangle \approx  \left| \tilde{w} ( \vec{k} ) \right|^2  \sum_{\alpha : \, |r_\alpha| \geq 0} \, g_f (\alpha) \cos( \vec{k} \cdot (\vec{r}_\alpha) )  \,\, .
 
−
</math>
 
−
</blockquote>
 
−
 
 
−
Here, we distinguish the momentum distribution from the interference
 
−
<blockquote>
 
−
<math>
 
−
\langle n_k (\vec{k} ) \rangle =  \sum_{\alpha : \, |r_\alpha| \geq 0} \, g_f (\alpha) \cos( \vec{k} \cdot (\vec{r}_\alpha) )  \,\, .
 
−
</math>
 
−
</blockquote>
 
−
 
 
−
Please refer to the list of [[#link_measurements|measurement observables]] provided by the DWA code.
 

Revision as of 18:34, 13 September 2013

Bandstructure of an homogeneous optical lattice

Theory

At this first moment, we shall look at the simplest case, i.e. a single particle of mass m which experiences a periodic potential V(\vec{r}), where

V(\vec{r}) = \sum_{x_\alpha = x,y,z} V_0^{x_\alpha} \sin^2 (\pi x_\alpha)

in the units of recoil energy E_r^\alpha = \frac{\hbar^2}{2m} \left( \frac{2\pi}{\lambda_\alpha} \right)^2 and lattice spacing \frac{\lambda_\alpha}{2}.

The quantum mechanical behaviour of the single particle follows

\left[\frac{1}{\pi^2} \left( -i \nabla + 2\pi \vec{k} \right)^2 + \sum_{x_\alpha = x,y,z} V_0^{x_\alpha} \sin^2 (\pi x_\alpha) \right] u_k (\vec{r}) = \epsilon_k u_k(\vec{r})

which is clearly separable to say the x-component:

\left[\frac{1}{\pi^2} \left( -i \partial_x + 2\pi k_x \right)^2 + V_0^{x} \sin^2 (\pi x) \right] u_{k_x} (x) = \epsilon_{k_x} u_{k_x}(x) \,\,\,\,\, ,\,\,\,\,\, k_x = 0, \frac{1}{L_x} ,\cdots \frac{L_x-1}{L_x}

In the plane wave basis,

u_{k_x} (x) = \frac{1}{\sqrt{L_x}} \sum_{m \in \mathbf{Z}} c_m^{(k_x)} e^{i2m\pi x}

we arrive at a tridiagonal diagonalization problem:

\left[ 4(m + k_x)^2 + \frac{V_0^x}{2} \right] c_m^{(k_x)} - \frac{V_0^x}{4} c_{m-1}^{(k_x)} - \frac{V_0^x}{4} c_{m+1}^{(k_x)} = \epsilon_{k_x} c_m^{(k_x)}

The wannier function is defined as:

w(x) = \frac{1}{\sqrt{L_x}} \sum_{k_x} u_{k_x} (x) e^{i 2\pi k_x x} = \frac{1}{L_x} \sum_{k_x} \sum_{m \in \mathbf{Z}} c_m^{(k_x)} e^{i 2\pi (m+k_x) x}

and from there, one can calculate the onsite interaction:

U = g \int | w(x) |^4 dx = \frac{4 \pi a_s \hbar^2}{m} \int | w(x) |^4 dx

After a little bit of algebra, we arrive at the hopping strength:

t = -\frac{1}{L_x} \sum_{k_x} \epsilon_{k_x} e^{-i2\pi k_x}

Finally, the Fourier transform of the wannier function is:

\tilde{w}(q_x) = \frac{1}{\sqrt{L_x}} \int w(x) e^{-i2\pi q_x x} dx = \frac{1}{\sqrt{L_x}} \sum_{k_x} \sum_{m \in \mathbf{Z}} c_m^{(k_x)} \delta_{q_x, k_x+m}

Implementation in Python

An example

For instance: 

import numpy;
import pyalps.dwa;

V0   = numpy.array([8. , 8. , 8.]);      # in recoil energies
wlen = numpy.array([843., 843., 843.]);  # in nanometer
a    = 114.8;                            # s-wave scattering length in bohr radius
m    = 86.99;                            # mass in atomic mass unit
L    = 200;                              # lattice size (along 1 direction)

band = pyalps.dwa.bandstructure(V0, wlen, a, m, L);

A first glance of the band structure: 

>>> band

Optical lattice: 
================
V0    [Er] = 8	8	8	
lamda [nm] = 843	843	843	
Er2nK      = 154.89	154.89	154.89	
L          = 200 
g          = 5.68473

Band structure:
===============
t [nK] : 4.77051	4.77051	4.77051	
U [nK] : 38.7018
U/t    : 8.11272	8.11272	8.11272	

wk2[0 ,0 ,0 ] : 5.81884e-08
wk2[pi,pi,pi] : 1.39558e-08

Well, the values of t (nK), U (nK), and U/t can be obtained via: 

>>> numpy.array(band.t())
array([ 4.77050984,  4.77050984,  4.77050984])
>>>
>>> numpy.array(band.U())
array(38.7018197381118)
>>>
>>> numpy.array(band.Ut())
array([ 8.11272192,  8.11272192,  8.11272192])

In momentum (\vec{q} ) space, the (squared) wannier function |\tilde{w}(\vec{q})|^2 can be obtained in the x-direction from: 

>>> numpy.array(band.q(0))
array([-5.   , -4.995, -4.99 , ...,  5.985,  5.99 ,  5.995])
>>> 
>>> numpy.array(band.wk2(0))
array([  7.57249518e-15,   7.88189086e-15,   8.20434507e-15, ...,
         1.62988573e-18,   1.56057426e-18,   1.49429285e-18])

and the y- or z- direction by replacing the index 0 to 1 and 2 respectively.

Bosons in an optical lattice trap

Boson Hubbard model

Hamiltonian

Bosons in an optical lattice trap can be effectively described by the single band boson Hubbard model

\hat{H} = -t \sum_{\langle i,j \rangle} \hat{b}_i^+ \hat{b}_j + \frac{U}{2} \sum_i \hat{n}_i (\hat{n}_i - 1) - \sum_i ( \mu - V_T ( \vec{r}_i) ) \hat{n}_i

with hopping strength t, onsite interaction strength U, and chemical potential \mu at finite temperature T via Quantum Monte Carlo implemented in the directed worm algorithm. Here, \hat{b} (\hat{b}^+) is the annihilation (creation) operator, and \hat{n}_i being the number operator at site i. Bosons in an optical lattice are confined, say in a 3D parabolic trapping potential, i.e.

V_T (\vec{r}_i) = K_x x_i^2 + K_y y_i^2 + K_z z_i^2 .

due to the gaussian beam waists as well as other sources of trapping.

Finite temperature

At finite temperature T, the physics is essentially captured by the partition function

Z = \mathrm{Tr} \, \exp \left(-\beta \hat{H} \right)

and physical quantities such as the local density

\langle n_i \rangle = \frac{1}{Z} \, \mathrm{Tr} \, \hat{n}_i \exp \left(-\beta \hat{H} \right) = \frac{1}{Z} \sum_{\mathcal{C}} n_i (\mathcal{C}) Z(\mathcal{C})

for some configuration \mathcal{C} in the complete configuration space, with inverse temperature \beta = 1/T . Here, the units will be cleverly normalized later on.


(c)

Retrieved from "https://alps.comp-phys.org/mediawiki/index.php?title=Documentation:Bosons_in_an_optical_lattice&oldid=7066"