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

(→An example) |
(→Bandstructure of an homogeneous optical lattice) |
||

Line 119: | Line 119: | ||

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

+ | |||

+ | === Model hamiltonian === | ||

+ | |||

+ | The DWA code simulates the single band boson Hubbard model<br/> | ||

+ | <blockquote> | ||

+ | <math> | ||

+ | \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 | ||

+ | </math> | ||

+ | </blockquote> | ||

+ | with hopping strength <math>t</math>, onsite interaction strength <math>U</math>, and chemical potential <math>\mu</math> at finite temperature <math>T</math> via Quantum Monte Carlo implemented in the directed worm algorithm. | ||

+ | Here, <math>\hat{b}</math> (<math>\hat{b}^+</math>) is the annihilation (creation) operator, and <math>\hat{n}_i</math> being the number operator at site i. | ||

+ | Bosons in an optical lattice are confined, say in a 3D parabolic trapping potential, i.e. | ||

+ | <blockquote> | ||

+ | <math> | ||

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

+ | </math>. | ||

+ | </blockquote> | ||

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

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

+ | <blockquote> | ||

+ | <math> | ||

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

+ | </math> | ||

+ | </blockquote> | ||

+ | and physical quantities such as the local density | ||

+ | <blockquote> | ||

+ | <math> | ||

+ | \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}) | ||

+ | </math> | ||

+ | </blockquote> | ||

+ | for some configuration <math>\mathcal{C}</math> in the complete configuration space, with inverse temperature <math>\beta = 1/T </math>. | ||

+ | Here, the units will be cleverly normalized later on. | ||

+ | |||

+ | |||

+ | === Momentum distribution and time-of-flight images === | ||

+ | |||

+ | The momentum distribution <math>\langle n(\vec{k})\rangle</math> of bosons in an optical lattice at equilbrium can be expressed as: | ||

+ | <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> | ||

+ | |||

+ | |||

+ | 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" | ||

+ | ! Option | ||

+ | ! Default | ||

+ | ! Remark | ||

+ | |- | ||

+ | | tof_phase | ||

+ | | align="center" | 0. | ||

+ | | time-of-flight phase <math>\gamma</math> | ||

+ | |} | ||

+ | |||

+ | |||

+ | === 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:27, 13 September 2013

## Contents

# 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 , where

in the units of recoil energy and lattice spacing .

The quantum mechanical behaviour of the single particle follows

which is clearly separable to say the x-component:

In the plane wave basis,

we arrive at a tridiagonal diagonalization problem:

The wannier function is defined as:

and from there, one can calculate the onsite interaction:

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

Finally, the Fourier transform of the wannier function is:

## 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 () space, the (squared) wannier function 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.

### Model hamiltonian

The DWA code simulates the single band boson Hubbard model

with hopping strength , onsite interaction strength , and chemical potential at finite temperature via Quantum Monte Carlo implemented in the directed worm algorithm. Here, () is the annihilation (creation) operator, and being the number operator at site i. Bosons in an optical lattice are confined, say in a 3D parabolic trapping potential, i.e.

.

due to the gaussian beam waists as well as other sources of trapping. At finite temperature T, the physics is essentially captured by the partition function

and physical quantities such as the local density

for some configuration in the complete configuration space, with inverse temperature . Here, the units will be cleverly normalized later on.

### Momentum distribution and time-of-flight images

The momentum distribution of bosons in an optical lattice at equilbrium can be expressed as:

Expressing the fourier transform of the wannier function as

and the interference term:

the momentum distribution reduces to the following form:

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.

where is the time-of-flight taken by the bosons to move from the origin (experiment) to the detector probe at position . 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

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 to the detector probe at , their time-dependent wavefunction must be:

where the kinetic energy of the free bosons must be:

Therefore, the corrected time-of-flight image must be:

where the dependence of the initial site on the wannier enveloped function is neglected. Finally, we arrive at:

with time-of-flight phase .

For example, the time-of-flight phase 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.

Option | Default | Remark |
---|---|---|

tof_phase | 0. | time-of-flight phase |

### Symmetry leads to further simplification

Symmetries:

Therefore:

Define green function:

In reduced coordinates, , we redefine green function:

and therefore the momentum distribution/ time-of-flight image:

Here, we distinguish the momentum distribution from the interference

Please refer to the list of measurement observables provided by the DWA code.