$W_{N}^{A}=W_{N}^{e1}+Y_{N}^{e1}+W_{N}^{e}+Y_{N}^{e}$ and, consequently,
$w_{A}=w_{e1}+y_{e1}+w_{e}+y_{e}$. Finding all eight functions resolves the
problem of interest.

The functions of interest satisfy the equations describing their evolution
with increasing the number of spins $N$ as

$\displaystyle
Y_{N}^{o}=\sum_{k=2}^{\infty}Y_{N-k}^{e1},~{}Y_{N}^{e}=\sum_{k=2}^{\infty}Y_{N-k}^{o1},~{}Y_{N}^{o1}=Y_{N-1}^{e}+Y_{N-1}^{e1}+2\delta_{N1},~{}Y_{N}^{e1}=Y_{N-1}^{o}+Y_{N-1}^{o1},$
$\displaystyle
W_{N}^{o}=\sum_{k=2}^{\infty}W_{N-k}^{e1}+2(1-\delta_{N,1}),~{}W_{N}^{e}=\sum_{k=2}^{\infty}W_{N-k}^{o1},~{}W_{N}^{o1}=W_{N-1}^{e}+W_{N-1}^{e1},~{}W_{N}^{e1}=W_{N-1}^{o}+W_{N-1}^{o1}.$
(64)

These equations can be reformulated in the algebraic form for the generating
functions

$\displaystyle
y_{o}=\frac{x^{2}}{1-x}y_{e1},~{}y_{e}=\frac{x^{2}}{1-x}y_{o1},~{}y_{o1}=xy_{e}+xy_{e1}+2x,~{}y_{e1}=xy_{o}+xy_{o1},$
$\displaystyle
w_{o}=\frac{x^{2}}{1-x}w_{e1}+\frac{2x^{2}}{1-x},~{}w_{e}=\frac{x^{2}}{1-x}w_{o1},~{}w_{o1}=xw_{e}+xw_{e1},~{}w_{e1}=xw_{o}+xw_{o1}.$
(65)

Solving these equations we obtained the generating function in Eq. (63) in the
form

$\displaystyle
w_{im}=\frac{2x^{2}(x^{5}-2x^{4}+x^{3}-2x^{2}+1)}{(1-x)(x^{6}+x^{4}-2x+1)}.$
(66)

The exponential asymptotic $W_{N}\propto\lambda_{min}^{-N}$ is determined by
the root $\lambda_{min}$ of the equation
$x^{7}+x^{6}+x^{5}+x^{4}-2x^{2}-x+1=0$ possessing a minimum absolute value.
Solving this equation we get $\lambda_{min}=0.5698$. Consequently, the
asymptotic behavior of interest can be expressed as $W_{N}^{im}=1.134\cdot
1.7549^{N}$. This result is quoted in the main text.

We evaluate the number of all frozen states that should be extracted from the
group II in a similar manner. All odd-numbered sequences within these states
contain a single spin and all even-numbered sequences contain even numbers of
spins or vice versa. Similarly to the previous consideration one can introduce
numbers of spin states $Z_{N}^{oo}$, $Z_{N}^{eo}$, $Z_{N}^{oe}$, $Z_{N}^{ee}$
for states consisting of odd or even numbers of sequences (first superscript)
and having one spin in odd- or even-numbered sequences (second superscript).
One can introduce corresponding generating functions satisfying the equations
below

$\displaystyle z_{oo}=xz_{eo}+2x,~{}z_{eo}=\frac{x^{2}}{1-x^{2}}z_{oo},$
$\displaystyle
z_{oe}=\frac{x^{2}}{1-x^{2}}z_{ee}+\frac{2x^{2}}{1-x^{2}},~{}z_{ee}=xz_{oe}.$
(67)

The generating function for the numbers of states of interest can be expressed
as (cf. Eq. (60))

$\displaystyle
w_{II,im}=z_{eo}+z_{ee}+z_{oo}\frac{x^{2}+x^{4}}{(1-x^{2})^{2}}=\frac{2x^{3}(2+x)}{1-2x^{2}}.$
(68)

Consequently, the numbers of states is defined by the series expansion of Eq.
(68) as

$\displaystyle
W_{N}^{II,im}=2^{N/2}\left[\frac{1+2\sqrt{2}}{4}+(-1)^{N}\frac{1-2\sqrt{2}}{4}\right]-\delta_{2,N}.$
(69)

###  B.3 Counting group III states

Here we evaluate the total number of spin states having no frozen sets with
the structure $\\{odd,1,odd\\}$, $\\{odd,1,(even,1)_{k},odd\\}$ ($k=1,2,...$,
see the main text). We consider the problem for open boundary conditions
(OBC). Problem with periodic boundary conditions (PBC) can be also solved
using generating function method, but it is overcomplicated so we leave it to
readers. The result for OBC gives the right asymptotic exponential dependence
on the system size for the number of states in the case of PBC and we estimate
the preexponential factor by comparing the analytical asymptotic behavior with
the numerical results.

To evaluate the desirable number of states we use iteration procedure similar
to that in the previous sections for the numbers $W_{N}^{e}$ and $W_{N}^{o}$
of spin states with the last sequence containing odd or even numbers of spins,
respectively, with all mobile spins (except for the very last one in the last
odd sequence) and numbers $Y_{N}^{e}$, $Y_{N}^{o}$ for states ended by even or
odd spin sequences, respectively, with two or more spins in the end being
frozen. Those states must have at the end the sequences $odd,(1,even)_{k}$
($k=1,2,...$) for $Y_{N}^{e}$ or the sequence $odd,(1,even)_{k},1$
($k=0,1,2,...$) for $Y_{N}^{o}$ so the sequence in the end must contain one
spin only. The first and last sequences in the state containing all mobile
spins with OBC must have even numbers of spins because the edge spin in the
last odd sequence is always frozen.

Then the iteration equations connecting the numbers of states with $N$ spins
with the numbers of states with smaller numbers of spins can be written as
(cf. Eq. (51))

$\displaystyle
W_{N}^{e}=\sum_{k=1}^{\infty}(W_{N-2k}^{e}+W_{N-2k}^{o}+Y_{N-2k}^{e}),~{}W_{N}^{o}=\sum_{k=0}^{\infty}W_{N-2k-1}^{e}+\sum_{k=1}^{\infty}(W_{N-2k-1}^{o}+Y_{N-2k-1}^{e}),$
$\displaystyle
Y_{N}^{e}=\sum_{k=1}^{\infty}Y_{N-2k}^{o}+2\sum_{k=1}^{\infty}\delta_{N,2k},~{}Y_{N}^{o}=Y_{N-1}^{e}+W_{N-1}^{o}.$
(70)

Our target is the main exponential asymptotic for the total number of mobile
spins $W_{m}\propto c\cdot A^{N}$. We wish to determine the factor $A$
analytically, while the prefactor will be estimated numerically.

Similarly to the previous section, Eq. (52), we introduce generating functions
using low case letter names with the same subscripts as superscripts in
capital letter notations. Eq. (70) can be rewritten in terms of generating
functions as

$\displaystyle\left(1-\frac{x^{2}}{1-x^{2}}\right)w_{e}-\frac{x^{2}}{1-x^{2}}w_{o}-\frac{x^{2}}{1-x^{2}}y_{e}=0,~{}-w_{e}\frac{x}{1-x^{2}}+w_{o}\left(1-\frac{x^{3}}{1-x^{2}}\right)-\frac{x^{3}}{1-x^{2}}y_{e}=0,$
$\displaystyle
y_{e}-\frac{x^{2}}{1-x^{2}}y_{o}=\frac{2x^{2}}{1-x^{2}},~{}-xw_{o}-xy_{e}+y_{o}=0.$
(71)

Eq. (71) has the standard form of the linear equation
$\widehat{M}\mathbf{x}=\mathbf{v}$ with the matrix $\widehat{M}$ defined as

$\displaystyle\widehat{M}=\begin{pmatrix}1-\frac{x^{2}}{1-x^{2}}&-\frac{x^{2}}{1-x^{2}}&-\frac{x^{2}}{1-x^{2}}&0\\\
-\frac{x}{1-x^{2}}&1-\frac{x^{3}}{1-x^{2}}&-\frac{x^{3}}{1-x^{2}}&0\\\
0&0&1&-\frac{x^{2}}{1-x^{2}}\\\ 0&-x&-x&1\end{pmatrix}$ (72)

The solutions of Eq. (71) are inversely proportional to the determinant of
this matrix. They have poles $1/(x-\lambda)$ at zeros $\lambda$ of that
determinant. The large $N$ asymptotic of the solution $W_{No,e}$ can be
expressed using the root $\lambda_{min}$ having the smallest absolute value in
the form $W_{No,e}\propto\lambda^{-N}$. This can be proved using the partial
fraction expansion similarly to previous sections.

The determinant of the matrix $\widehat{M}$ can be expressed as

$\displaystyle{\rm
det}\widehat{M}=\frac{4x^{5}-2x^{4}-3x^{3}-3x^{2}+1}{(1-x^{2})^{2}}.$ (73)

We found the solutions of the equation det $\widehat{M}=0$ numerically. The
solution with the minimum absolute value is $\lambda_{min}\approx 0.5423$.
Consequently, one can expect the number of states without frozen spins to
scale as $W_{N}\propto\lambda_{min}^{-N}=1.8442^{N}$. This number increases
faster then the numbers of states in groups I and II and, therefore,
subtracting their numbers will not affect the main exponential asymptotic.

The solution for the generating function of interest $w_{e}(x)$ representing
the states with all mobile spins takes the form

$\displaystyle w_{e}(x)=\frac{x^{4}}{(4x^{5}+2x^{4}-3x^{3}-3x^{2}+1)}.$ (74)

Expanding this result into simple fractions one can define the preexponential
factor in the dominating exponential asymptotic for the OBC problem that is
$W^{e}_{N}=0.0551\cdot 1.8442^{N}$

To derive the exponential asymptotic for the number of group III states in the
periodic system we compare analytical results vs. numerical calculations of
the numbers of states having no frozen spins as reported in Fig. 11. It turns
out that the best exponential fit can be obtained for the total number of
states belonging to groups I, II and III and this fit is $W^{e}_{N}=1.0\cdot
1.8442^{N}$. The number of group III states approaches this asymptotic with
increasing the number of spins as it supposes to be because the relative
weight of states of groups I and II decreases exponentially with increasing
the number of spins compared to that of the group III.

Figure 11: Numbers of states with all mobile spins vs the number of spins
compared to the analytical fit $1.0\cdot 1.8442^{N}$

##  Appendix C Imbalance as a measure of system dynamic properties

Here we consider the definition of an infinite time imbalance in greater
detail with special attention to its experimental determination. The imbalance
is defined following the conditions of a typical experiment similarly to Ref.
[2] (Eq. (39) in the main text). It is determined by the system evolution from
some initial product state $|a>$, characterized by fixed projections of spins
$\pm 1/2$ to the $z$ axis ($t=0$). The system evolution can be expressed in
terms of average spin projections $<S_{k}^{z}(t)>$ at time $t$. The extraction
of the infinite time limit for $<S_{k}^{z}(t)>$ is not straightforward using
the experimental data available at a finite time only. Yet it can be estimated
reasonably well as described below.

According to our observations imbalance converges to its infinite time limit
quite fast for the initial state belonging to the groups I or II with
delocalized eigenstates. For the family of initial states belonging to the
group II, shown in Fig. 12.a, imbalances approach its infinite time limit
during the time of order of few inverse spin hopping amplitudes $\Delta^{-1}$
Eq. (1). The convergence gets better with increasing the number of spins as
illustrated in Fig. 12.a,.

For the localized state considered in Fig. 12.b, originated from the group III
state, shown in the graph, the oscillations of imbalance around its infinite
time limit are very strong. However, one can estimate the infinite time limit
of imbalance using time averaged imbalance defined as

$\displaystyle I_{av}(t)=\frac{\int_{0}^{t}I(\tau)d\tau}{t}.$ (75)

According to Fig. 12.b the time averaged imbalance $I_{av}(t)$ converges to
its infinite time limit after the time around few inverse hopping amplitudes
$\Delta^{-1}$ similarly to Fig. 12.a. Thus the infinite time limit of
imbalance can be determined experimentally and compared with the predictions
of theory. One should notice that this definition can be not applicable in
marginal situations of very slow dynamics emerging for instance for unlikely
realized states of anomalously high or low energy [76].

(a)  (b)

Figure 12: (a) Convergence of the imbalances for the family of states
belonging to the group II to the infinite time limit. (b) Convergence of the
time averaged imbalance to its infinite time limit for the localized state
belonging to the group III. The initial product states are shown within the
graphs in the inverted representation.

Eq. (39) determines the imbalance measured in Ref. [2] and we use it here for
a characterization of different groups of states. The Hamming distances
reported in Ref. [1] can be expressed as

$\displaystyle HD(t)=1-\frac{4}{N}\sum_{k=1}^{N}<S_{k}^{z}(t)><S_{k}^{z}(0)>.$
(76)

In a thermodynamic limit of $N\rightarrow\infty$ one can set $<S^{z}>=0$ in
the definition of imbalance (Eq. (39) in the main text). Then two definitions
becomes equivalent, since $HD(t)+I(t)=1$. However, we prefer to use Eq. (39)
or its generalizations above since they are less sensitive to finite size
effects due to the subtraction of the finite size expectation values for
average spin projections.


# Anisotropic Coherence properties in a trapped quasi2D dipolar gas

Christopher Ticknor  Theoretical Division, Los Alamos National Laboratory, Los
Alamos, New Mexico 87545, USA

###### Abstract

We consider a trapped quasi2D dipolar Bose Einstein condensate (q2D DBEC) with
a polarization tilted into the plane of motion. We show that by tilting the
polarization axis, the coherence properties are anisotropic. Such a system
will have density fluctuations whose amplitude depend on their trap location.
Additionally, interference contrast will also be anisotropic despite an
isotropic density profile. The anisotropy is related to a roton like mode that
becomes unstable and supports local collapse along the polarization axis.

###### pacs:

67.85.-d, 03.75.Kk, 03.75.Lm, 47.37.+q

Coherence properties underlie collective quantum behavior. An example of a
need to further understand coherence is in high $T_{c}$ superconductors
phaseTC , where the origin of the phase coherence leading to superconducting
behavior is still mysterious. Such solid state systems where this occurs are
complex, making the study of coherence properties a challenge. One avenue to
study coherence of quantum systems more clearly is in ultracold gases, which
have offered direct insight and unprecedented control to the study of quantum
collective behavior. For example, two macroscopic Bose Einstein condensates
(BECs) were interfered displaying phase coherence int . More refined
experiments have directly detected higher order correlations correlations
which cannot be easily accessed by electronic systems. Additionally, reduced
dimensional ultracold gases have been the focus of intense interest because
fluctuations are enhanced in such geometries. Phase coherence has been studied
at the Berezinskii Kosterlitz Thouless (BKT) transition BKT ; BKTrev  where
two independent q2D ultracold gases were interfered, revealing thermally
activated phase defects are the origin of decoherence zoran . The BKT
transition is a transition where the phase coherence changes character, and
insights into controlling such coherence properties are of general importance.
In this paper we show the coherence properties of a trapped dipolar gas can be
controlled.

The progress of experiments to probe ultracold gases is remarkable. Recent
advances in imaging q2D ultracold gases have been used to detect density
fluctuations chicago ; jila ; ENS . Additionally, ultracold atomic systems in
one dimensional systems have also amazingly been used to study quantum
fluctuations directly armijo2  and other coherence properties armijo1 ; manz ;
imam ; kinoshita . For these examples of ultracold gases, the interaction
between the constituents are short range. There are ultracold dipolar systems
which have anisotropic, long range interactions, which offer an entirely new
avenue to study quantum correlations and collective behavior. Amazing
experimental progress has been made with the strongly magnetic atoms: Chromium
(Cr), Dysprosium (Dy), and Erbium (Er). All three have been Bose condensed and
displayed strong dipolar effects Cr ; Dy ; Er . The study of dipoles in
reduced geometries have begun with Cr, where it has been put into a one
dimensional lattice, forming a series of coupled q2D dipolar BECs layer .
Additionally, polar molecules have been loaded into q2D geometries to control
the rate of their chemical reaction rbk .

Correlations of ultracold dipolar gases have been studied with a perpendicular
polarization for weakly interacting gases near the roton instability klaw ;
sykes ; blair  and strongly interacting gases, e.g. buchler . Using the
direction of the polarization axis as a means to control the collective
behavior and coherence properties, a q2D dipolar gas has been studied in the
many body dipolar systems at the mean field level nath ; aniso  and for 2D
dipolar scattering tilt . Additional theories have investigated a tilted
polarization axis in a q2D dipolar fermi gas tilt2 .

In this work, we go beyond mean field theory and examine the correlations of a
trapped q2D DBEC at finite temperature with a partially tilted polarization
into the plane of motion. This situation leads to anisotropic coherence
properties, in which the gas maintains coherence along the direction of
polarization, but more rapidly loses coherence perpendicular to the
polarization direction. This anisotropy is evident in both the phase and
density-density coherence properties. This emerges because of a roton-like
mode in the excitation spectrum and is strongly anisotropic in character. We
consider the implications of this theory through interference experiments and
density fluctuations, both of which will exhibit anisotropic correlated
behavior. We also calculate the compressibility of the dipolar gas.

To study the finite temperature coherence properties, we employ the Hartree
Fock Bogoliubov method within the Popov approximation (HFBP) with nonlocal
interactions griffin ; CThfb . The HFB breaks the wavefunction into a
condensate and thermal component:
$\hat{\Psi}=[\sqrt{N_{0}}\phi_{0}(\vec{\rho})+\hat{\theta}(\vec{\rho})]$ where
we have replaced $\hat{a}_{0}\rightarrow\sqrt{N_{0}}$. We use the Bogoliubov
transformation:
$\hat{\theta}(\vec{\rho})=\sum_{\gamma}[u_{\gamma}(\vec{\rho})\hat{a}_{\gamma}e^{-i\omega_{\gamma}t}-v_{\gamma}^{*}(\vec{\rho})\hat{a}^{*}_{\gamma}e^{i\omega_{\gamma}t}]$
where $\hat{a}_{\gamma}$ ($\hat{a}^{*}_{\gamma}$) is the bosonic annihilation
(creation) operator for the $\gamma^{th}$ quasiparticle. The HFBP solves a
non-local, generalized Gross-Pitaevskii equation and a nonlocal Bogoliubov de
Gennes set of equations from which the eigenvalues and vectors
(quasiparticles) are self-consistently obtained CThfb .

To evaluate the interaction, we assume a dipole moment of the form:
$\vec{d}=d[\hat{x}\sin(\alpha)+\hat{z}\cos(\alpha)]$ where $\alpha$ is the
angle between $\hat{z}$ and $\hat{d}$, and we assume a single transverse
Gaussian wavefunction of width $l_{z}$ is occupied in the $z$ direction, where
$l_{i}=\sqrt{\hbar/\omega_{i}m}$ and $\omega_{i}$ is the trapping frequency.
This requires $\mu,T\ll\hbar\omega_{z}$, and this leads to a simple form of
the interaction. We use an isotropic trap ($x,y$) with $l_{z}/l_{\rho}=0.15$
($\omega_{z}/\omega_{\rho}$=44.4), and an interaction strength of
$g_{d}=0.025$ and $g=0$, where $g_{d}=d^{2}/\sqrt{2\pi}l_{z}$ and
$g=\sqrt{8\pi}a/l_{z}$ with $a$ being the 3D s-wave scattering length ($a\ll
l_{z}$). The polarization angle is set to $\alpha/\pi=0.25$, this is beyond
the critical angle. The critical angle is when the interaction is zero along
the $x$ axis, so $\sin(\alpha_{c})=1/\sqrt{3}$. For $\alpha>\alpha_{c}\sim
0.2\pi$ there is an attractive region to the interaction. This value of
$\alpha$ is picked to optimize the collapse with respect to particle number.
If it were steeper, the critical number would be less, and the depletion would
not be large enough to observe the desired behavior. In addition, if $\alpha$
were smaller then roton the mode would require many more particles to
collapse, and the collapses would be more 3D in nature. For $\alpha/\pi=0.25$,
the system is well described by the q2D formalism. For these parameters, the
trapping potential for Dy [Er] is $(\omega_{\rho},\omega_{z})/2\pi$= (12,533)
[(48,2133)] Hz. In a trapped ideal 2D gas, Bose condensation can occur trap ,
and the the critical temperature is at $T_{C}/\hbar\omega=\sqrt{6N}/\pi\sim
0.78\sqrt{N}$. This is for the pure 2D trapped gas, we only use it as a
reference temperature. Throughout this paper we work in trap units. It is
important to mention that the examples given here are weakly interacting; the
HFBP method cannot accurately handle strongly interacting dipolar gases
($\mu<10\hbar\omega_{\rho}$) CThfb . This criteria is from the deviation of
the Kohn mode, or slosh mode, from 1$\hbar\omega_{\rho}$. If the interaction
is too strong then the static thermal cloud alters the slosh mode of the
condensate; this is a numerical signature that the HFBP has broken down. This
restricted the work from entering the Thomas Fermi regime where
$l_{\rho}>\xi$, where $\xi$ is the healing length ($\sqrt{m\mu/\hbar^{2}}$).
Because of this, trapping effects are always important in this study.

Once $g_{d}$, $g$, $l_{z}$, $T$, $\alpha$ and the condensate number, $N_{0}$,
are selected, the HFBP calculation finds: the number of thermal atoms,
$\tilde{N}$; the chemical potential, $\mu$; the condensate wavefunction,
$\phi_{0}$; and the quasiparticle wavefunctions, $u_{\gamma}$ and
$v_{\gamma}$. From these wavefunctions, we construct the correlation
functions: the nonlocal thermal correlation function,
$\tilde{n}(\vec{\rho},\vec{\rho}^{\prime})=\sum_{\gamma}[u^{*}_{\gamma}(\vec{\rho})u_{\gamma}(\vec{\rho}^{\prime})+v^{*}_{\gamma}(\vec{\rho})v_{\gamma}(\vec{\rho}^{\prime})]N_{\gamma}+v^{*}_{\gamma}(\vec{\rho})v_{\gamma}(\vec{\rho}^{\prime})$
with $N_{\gamma}$ being the Bose-Einstein occupation, and the condensate
correlation function, $n_{0}(\vec{\rho},\vec{\rho}^{\prime})$, given by
$N_{0}\phi^{*}_{0}(x)\phi_{0}(x^{\prime})$. The total local density is
$n(\vec{\rho})=n_{0}(\vec{\rho})+\tilde{n}(\vec{\rho})$ where e.g.
$n_{0}(\vec{\rho})=n_{0}(\vec{\rho},\vec{\rho})$. The HFBP method uses these
correlation functions in the self-consistent calculation, and they can be
related to the more common $g_{1}$ and $g_{2}$ correlation functions phase :

$\displaystyle g_{1}(\vec{\rho};\vec{\rho}^{\prime})$
$\displaystyle={\langle\Psi^{*}(\vec{\rho})\Psi(\vec{\rho}^{\prime})\rangle\over\sqrt{n(\vec{\rho})n(\vec{\rho}^{\prime})}}={n_{0}(\vec{\rho},\vec{\rho}^{\prime})+\tilde{n}(\vec{\rho},\vec{\rho}^{\prime})\over\sqrt{n(\vec{\rho})n(\vec{\rho}^{\prime})}},$
(1) $\displaystyle g_{2}(\vec{\rho};\vec{\rho}^{\prime})$
$\displaystyle={\langle\Psi^{*}(\vec{\rho})\Psi^{*}(\vec{\rho}^{\prime})\Psi(\vec{\rho})\Psi(\vec{\rho}^{\prime})\rangle\over{n(\vec{\rho})n(\vec{\rho}^{\prime})}}$
$\displaystyle=1+g_{1}(\vec{\rho},\vec{\rho}^{\prime})^{2}-{n_{0}(\vec{\rho},\vec{\rho}^{\prime})^{2}\over
n(\vec{\rho})n(\vec{\rho}^{\prime})}.$

$g_{1}$ details phase correlations, for example it determines the fringe
contrast if a BEC were split and interfered. $g_{2}$ is the density-density
correlation function, and is related to density fluctuations and the static
structure factor zam .

We study a system with a temperature of $T/T_{C}=0.5$. At this temperature,
the condensate fraction is $N_{0}/N>0.57$ for all $N$ considered. This is well
within the region of validity for the HFBP. This temperature thermally
populates the quasiparticles, so that the anisotropy in the coherence
properties is significant. Furthermore, this temperature makes the density
more isotropic. Values of $T/T_{C}$ between $0.4$ and $0.6$ produces similar
results. Raising the temperature would lead to populating more thermal modes,
and would make correlations isotropic, reflecting the nature of the trap
($T/T_{C}\sim 1$). For lower temperatures, the correlations are determined by
the condensate’s properties alone, as the depletion and thermal population
become negligible. A pure condensate has total phase coherence, it is the
thermal fraction that leads to reduced coherence.

We now look at the correlation functions, $g_{1}(\vec{\rho};0)$ and
$g_{2}(\vec{\rho};0)$. These quantify the coherence of the system from the
origin to $\vec{\rho}$. To illustrate the anisotropic coherence of the gas, we
plot them along the polarization axis and in the perpendicular direction in
figure 1.

Figure 1: (a) The total density, (b) $g_{1}(x,0;0,0)$ (solid line) and
$g_{1}(0,y;0,0)$ (dashed), and (c) $g_{2}(x,0;0,0)$ (solid line) and
$g_{2}(0,y;0,0)$ (dashed) at $T/T_{C}$=0.50 and for $N$ =2000 (blue), 3000
(red), and 4000 (black) from top to bottom. The total density is isotropic in
the center and yet the correlation function is significantly anisotropic at
short range.

Fig. 1 (a) shows the total density and Fig. 1 (b) shows $g_{1}(\vec{\rho};0)$
along the $x$ (solid line) and $y$ (dashed) axes. The important feature is
that $g_{1}$ quickly decreases in the $y$ direction, but only steadily
decreases its value in the $x$ direction or polarization direction. We have
shown the system for three different number of total particles 2000 (blue),
3000 (red), and 4000 (black), from top to bottom at the temperature
$T/T_{C}=0.5$. Varying the number is like varying the interaction strength.
The density is isotropic near the center of the trap, even when the $g_{1}$ is
anisotropic. Near the edge of the trap the total density becomes anisotropic
and is more extended in the $x$ axis.

A recent Bragg spectroscopy experiment observed an anisotropic speed of sound
for a 3D dipolar gas bragg . A similar experiment for a q2D DBEC could be used
to probe these predicted anisotropic coherence properties. If a DBEC were
split into two equal clouds such that one is moved by $\vec{a}$ with a
relative phase of $S(\vec{\rho})$, then the total density would be
$n(\vec{\rho})=0.5(n(\vec{\rho})+n(\vec{\rho}+a))+\sqrt{n(\vec{\rho})n(\vec{\rho}+\vec{a})}g_{1}(\vec{\rho},\vec{\rho}+\vec{a})\cos[S(\vec{\rho})]$.
The fringe contrast will depend on the angle between $\vec{a}$ and $\vec{d}$
because of the anisotropy in $g_{1}(\vec{\rho},\vec{\rho}+\vec{a})$. To
maximize the difference in fringe contrast, the phase difference should be
$\pi/2$ when $\vec{a}$ is near 3/4 $l_{\rho}$. This signal is only a few
percent and would be challenging to measure.

We now move on to the density-density correlation function,
$g_{2}(\vec{\rho};0)$. In Fig. 1 (c), $g_{2}$ this is plotted along the $x$
(solid line) and $y$ (dashed) axis. The results are similar to those of
$g_{1}$, the extent of the coherence is further in the $x$ direction than the
$y$. The anisotropy is more significant in this correlation function. The fact
that the $g_{2}$ goes below one is a trap effect and is most significant when
correlating with the trap center manz ; imam . This correlation function is
important for density fluctuations, this will be discussed below.

To understand the nature of the anisotropy we look at the Bogoliubov de Gennes
excitation spectrum as a function of particle number, this is shown in Fig. 2
(a). The spectrum goes soft or an eigenvalue becomes complex implying the
system is unstable when there are more than 4850 particles. Additionally, we
show contour plots of the wavefunctions of the two modes which go soft
($u_{\gamma}$) in Fig. 2 (b,c), and the dashed black lines are $\sqrt{n_{0}}$.
The contours at 0.25, 0.5, and 0.75 of the maximum values of the individual
wavefunctions, and for $u_{\gamma}$ the shaded regions signify a region less
than zero. These quasiparticle modes start with an energy at or above
$4\hbar\omega_{\rho}$. Then as the particle number is increased, their energy
decreases. In a homogeneous q2D DBEC, the roton is a local minimum in the
excitation spectrum, which occurs near the $k_{\rho}l_{z}\sim 1$ roton1 . If
the spectrum goes soft, the system locally collapses into plane waves which
form dense stripes along the polarization axis.

The Bogoliubov de Gennes excitation spectrum has been analyzed in detail in
Ref. CThfb  for a dipolar gas with $\alpha=0$ or a cylindrically symmetric
system. In that case there are degenerate modes with azimuthal symmetry $\pm
m$. In the present case this degeneracy has been broken by the tilted
polarization axis, leading to slightly different energies depending on the
excitation shape relative to the $x$ axis. This leads to band like structures
forming in Fig. 2 (a).

The mode in 2 (b) makes a significant contribution to the correlation function
for two reasons: first, it is a low energy mode, and is more thermally
populated. Second, the correlation function, $g_{i}(\vec{\rho}:0)$, has a
factor of $u_{\gamma}(0)$ and $v_{\gamma}(0)$ in it. Only a few low energy
modes have a significant contribution at the origin. The breathing mode just
below 2$\hbar\omega$ and the roton mode in 2 (b) are two low energy modes with
a maximum in the center of the trap. Other low energy modes have small or zero
amplitude in the center of the cloud because symmetry.

In the trapped case, the finite extent of the gas alters the nature of the
“roton” mode roton ; Wilson ; Wilson2 . For $\alpha=0$, the gas become 3D in
nature, $\mu\sim\hbar\omega_{z}$ blair . When the gas becomes unstable, it
locally collapses and forms regions of high density. A good demonstration of
this is given in Ref. parker . For $\alpha/\pi\sim 0.25$, the stability of the
system displays an interplay of the interaction and the trap geometry, while
maintaining $\mu\ll\hbar\omega_{z}$. The collapse is local and forms a series
density regions parallel to the polarization. In Fig. 2 (b,c), we see this
character in the quasiparticles, which have an extended amplitude along the
polarization axis and a series of nodes in the perpendicular direction. It is
this character of the quasiparticles that leads to anisotropic character to
the correlations. The nodes in the roton modes lead to destructive
interference and reduce correlation in that direction.

Figure 2: (a) The Bogoliubov de Gennes spectrum is shown as a function of
particle number for $g_{d}=0.025$ for $\alpha/\pi=0.25$. (b,c) A contour plot
of the roton quasiparticle modes ($u_{\gamma}$) which go soft for a dipolar
system with tilted polarization axis are also shown. The dashed black line is
$\sqrt{n_{0}}$ and for $u_{\gamma}$ the shaded regions are less than zero. The
contours at 0.25, 0.5, and 0.75 the maximum values of the individual
wavefunctions.

To further study the implications of the anisotropic correlations, we consider
the density fluctuations of the gas by calculating the compressibility:
$kT\kappa(\vec{\rho})=$ $\langle\int d\vec{\rho}^{\prime}\delta
n(\vec{\rho})\delta n(\vec{\rho}^{\prime})\rangle$ =$\int
d\vec{\rho}^{\prime}n(\vec{\rho})n(\vec{\rho}^{\prime})g_{2}(\vec{\rho};\vec{\rho}^{\prime})-n(\vec{\rho})(N-1)$
where $\delta
n(\vec{\rho})=\hat{\Psi}^{*}(\vec{\rho})\hat{\Psi}(\vec{\rho})-n(\vec{\rho})$
and this is the density fluctuations of the system. The compressibility and
number fluctuations have been measured in BEC experiments chicago ; jila . We
present the compressibility for a trapped dipolar gas in Fig. 3 (b) at
$T/T_{0}=0.5$ for 2000 (blue), 3000 (red), and 4000 (black) in the $x$ (solid)
and $y$ (dashed). These results show that the density fluctuations roughly
follow the anisotropic density profile. Another way to state this is: that
depending on the local value of the density, the compressibility is determined
almost uniquely.

Figure 3: (a) The total density (b) the compressibility,
$kT\kappa(\vec{\rho})$, for the dipolar gas along the $x$ (solid) and $y$
(dashed) axis for particles number 2000 (blue), 3000 (red), and 4000 (black).

In this work, we found that the correlations are anisotropic in a q2D DBEC
with a polarization axis tilted into the plane of motion. We studied the DBEC
with the HFBP at at temperature of $T/T_{C}\sim 0.5$. At this temperature,
quasiparticles are thermally populated and total density is nearly isotropic
in the central region of the trap. We computed the $g_{1}(\vec{\rho};0)$ and
$g_{2}(\vec{\rho};0)$ correlations functions, which have anisotropic character
in the center of the gas even though the density is isotropic in the central
region of the DBEC. This has measurable implications for interference
properties and density fluctuations of the gas. We found that the anisotropic
correlations occur due to roton like quasiparticle modes that can lead to
local collapse along the polarization axis.

An immediate extension of this work is to study the static structure factor
such as was recently done in Ref. blair  at zero temperature. An important
difference is that the static structure will become anisotropic, i.e.
$S(k_{x},k_{y})$, and with present methods, it could be done at finite
temperature stemp . This also presents a promising method to observe the
anisotropic coherences predicted here.

Future work will be to investigate the impact of correlations at higher
temperatures and strongly interactions where they might impact the BKT
transition BKT ; BKTrev ; dBKT . This will require an improvement to the
method mora  or new approach altogether. Either way, it will be intriguing to
study the impact of the anisotropic dipolar interaction on vortex correlations
near the BKT transition and the time dependence of such a gas. The present
work was not able to measure the decay of the correlation functions due to the
relatively weak interaction strength or small number the numerical method
could handle. If the next method can probe stronger interactions, we would
hope to observe a different character in the decay of the correlation
functions along and perpendicular to the polarization axis.

Acknowledgments, the author is grateful for support from the Advanced
Simulation and Computing Program (ASC) and LANL which is operated by LANS, LLC
for the NNSA of the U.S. DOE under Contract No. DE-AC52-06NA25396. The author
is grateful for discussion with A. Sykes, R. Behunin, and L. A. Collins.

## References

  * (1) T. Timusk and B. Statt, Rep. Prog. Phys. 62, 61 (1999); V. M. Loktev, R. M. Quick, S. G. Sharapov, Physics Reports, 349, 1 (2001).
  * (2) M. R. Andrews et al., Science 275, 553 (1997).
  * (3) For example: A. Ottl et al., Phys. Rev. Lett. 95 090404 (2005); S. S. Hodgman et al., Science 331, 1046 (2011); S. Foelling et al., Nature 434, 481 (2005).
  * (4) Berezinskii, V. L. Sov. Phys. JETP 34, 610 (1972). L. M. Kosterlitz, D. J. Thouless, J. Phys. C 6, 1181 (1973);
  * (5) Z. Hadzibabic and J. Dalibard, Riv. Nuo. Cime., 34, 389 (2011).
  * (6) Z. Hadzibabic et al., Nature 441, 1118 (2006).
  * (7) T. Yefsah et al., Phys. Rev. Lett. 107, 130401 (2011).
  * (8) C.L. Hung et al., Nature 470, 236 (2011); N. Gemelke et al., Nature 460, 995 (2009).
  * (9) S. Tung et al., Phys. Rev. Lett. 105, 230408 (2010).
  * (10) J. Armijo, Phys. Rev. Lett. 108, 225306 (2012).
  * (11) S. Manz et al., Phys. Rev. A 81, 031610 (2010).
  * (12) A. Imambekov et al., Phys. Rev. A 80, 033612 (2009).
  * (13) T. Kinoshita, T. Wenger, and D. S. Weiss, Phys. Rev. Lett. 95, 190406 (2005).
  * (14) J. Armijo et al., Phys. Rev. Lett. 105, 230402 (2010).
  * (15) Th. Lahaye et al., Nature 448, 672 (2007); Th. Lahaye et al., Nature Physics 4, 218 (2008); Th. Lahaye et al., Phys. Rev. Lett 101, 080401 (2008).
  * (16) K. Aikawa et al., Phys. Rev. Lett. 108, 210401 (2012).
  * (17) M. Lu et al.,, Phys. Rev. Lett., 107, 109401 (2011); M. Lu et al.,, Phys. Rev. Lett., 108, 215301 (2012).
  * (18) S. Meuller et al., Phys. Rev. A 84, 053601 (2011).
  * (19) M. H. G. de Miranda et al., Nat. Phys. 7, 502 (2010).
  * (20) M. Klawunn et al. Phys. Rev. A 84, 033612 (2012).
  * (21) A. Sykes and C. Ticknor, arXiv1205.1350.
  * (22) P. B. Blakie, D. Baillie, and R. N. Bisset, Phys Rev. A, 86 021604(R) (2012).
  * (23) H. P. Buechler et al., Phys Rev. Lett 98, 060404 (2007);G. E. Astrakharchick et al., Phys Rev. Lett 98, 060405 (2007).
  * (24) R. Nath, P. Pedri, and L. Santos, Phys. Rev. Lett. 102 050401 (2009).
  * (25) C. Ticknor, R. M. Wilson, J. L. Bohn, Phys. Rev. Lett. 106, 065301 (2011).
  * (26) C. Ticknor, Phys. Rev A 84, 032702 (2011).
  * (27) M. M. Parish and F. M. Marchetti, Phys. Rev. Lett. 108, 145304 (2012); L. M. Sieberer and M. A. Baranov, Phys. Rev. A 84, 063633 (2012).
  * (28) A. Griffin, Phys. Rev. B, 53, 9341 (1996).
  * (29) C. Ticknor, Phys. Rev. A 85, 033629 (2012).
  * (30) V. Bagnato and D. Kleppner, Phys. Rev. A 44, 7439 (1991).
  * (31) M. Naraschewski and R. J. Glauber, Phys. Rev. A 59, 4595 (1999).
  * (32) F. Zambelli et al., Phys. Rev. A 61 063608 (2000).
  * (33) G. Bismut et al., arXiv 1205.6305.
  * (34) L. Santos, G. V. Shlyapnikov, and M. Lewenstein, Phys. Rev. Lett. 90, 250403 (2003).
  * (35) S. Ronen, D. C. E. Bortolotti, and J. L. Bohn, Phys. Rev. Lett. 98, 030406 (2007).
  * (36) R. Wilson, S. Ronen, and J.L. Bohn Phys. Rev. Lett. 104, 094501 (2010).
  * (37) R. M. Wilson, S. Ronen, and J. L. Bohn, Phys. Rev. Lett. 100, 245302 (2008).
  * (38) N. G. Parker et al. Phys. Rev. A 79, 013617 (2009).
  * (39) A. Filinov et al. Phys. Rev. Lett. 105, 070401 (2010).
  * (40) C. Mora and Y. Castin, Phys. Rev. A 67, 053615 (2003).
  * (41) A. Csordas, R. Graham, and P. Szepfalusy, Phys. Rev. A 54, R2543 (1996).