Finite Temperature Phases of Two Dimensional Spin-Orbit Coupled Bosons

TL;DR

This paper explores the finite temperature phase diagram of two-dimensional spin-orbit coupled bosons, revealing various phases including superfluid, stripe, and fractionalized quasi-condensates, influenced by anisotropy and isotropy in interactions and coupling.

Contribution

It provides the first detailed Monte Carlo analysis of finite temperature phases of 2D spin-orbit coupled bosons, highlighting the effects of anisotropy and isotropy on phase stability and order.

Findings

Identification of BKT superfluid transition with anisotropic spin-orbit coupling.

Discovery of spin-order driven by interaction anisotropy, leading to homogeneous or stripe phases.

Observation of fractionalized quasi-condensate in isotropic interaction conditions.

Abstract

We determine the finite temperature phase diagram of two dimensional bosons with two hyperfine (pseudo-spin) states coupled via Rashba-Dresselhaus spin-orbit interaction using classical field Monte Carlo calculations. For anisotropic spin-orbit coupling, we find a transition to a Berenzinskii-Kosterlitz-Thouless superfluid phase with quasi-long range order. We show that the spin-order of the quasi-condensate is driven by the anisotropy of interparticle interaction, favoring either a homogeneous plane wave state or stripe phase with broken translational symmetry. Both phases show characteristic behavior in the algebraically decaying spin density correlation function. For fully isotropic interparticle interaction, our calculations indicate a fractionalized quasi-condensate where the mean-field degeneracy of plane wave and stripe phase remains robust against critical fluctuations. In the case of fully isotropic spin-orbit coupling, the circular degeneracy of the single particle ground state destroys the algebraic ordered phase in the thermodynamic limit, but a cross-over remains for finite size systems.

Equations

$$H_0=\int d^2\mathbf{r}{\hat{\Psi }}^{†}(\mathbf{r})\left[-\frac{{\hslash }^2{\nabla }^2}{2m}-i\frac{{\hslash }^2\kappa }{m}\left({\sigma }_x{\partial }_x+{\eta }_{soc}{\sigma }_y{\partial }_y\right)\right]\hat{\Psi }(\mathbf{r})$$

$$V=\frac{1}{2}\sum _{\sigma ,{\sigma }^{\prime }=\uparrow \downarrow }{g}_{\sigma {\sigma }^{\prime }}{\hat{\Psi }}_{\sigma }^{†}(\mathbf{r}){\hat{\Psi }}_{{\sigma }^{\prime }}^{†}(\mathbf{r}){\hat{\Psi }}_{{\sigma }^{\prime }}(\mathbf{r}){\hat{\Psi }}_{\sigma }(\mathbf{r})$$

$$\displaystyle\psi^{MF}_{{{\boldsymbol{\kappa}}}}=\frac{1}{\sqrt{2}}\left[\cos(\phi)e^{i{{\boldsymbol{\kappa}}}\cdot{{\bf r}}}\left(\begin{array}[]{c}1\\ -e^{i\theta_{{\boldsymbol{\kappa}}}}\end{array}\right)+\sin(\phi)e^{-i{{\boldsymbol{\kappa}}}\cdot{{\bf r}}}\left(\begin{array}[]{c}1\\ e^{i\theta_{{\boldsymbol{\kappa}}}}\end{array}\right)\right]$$

$$S[\Psi (\mathbf{r})]$$

$$\Delta n=\frac{1}{L^2}\sum _{\mathbf{k},\alpha =\pm }\left[n_B({\epsilon }_{\mathbf{k}\alpha }^{mf}-\mu )-n_{cf}({\epsilon }_{\mathbf{k}\alpha }^{mf}-\mu )\right]$$

$$n_0^{\kappa }=\sum _{\mathbf{k}=(\pm \kappa ,0)}\sum _{\sigma {\sigma }^{\prime }}\int \frac{d\mathbf{r}d{\mathbf{r}}^{\prime }}{L^2}{\psi }_{\mathbf{k}\sigma }^{MF}(\mathbf{r})G_{\sigma ,{\sigma }^{\prime }}(\mathbf{r},{\mathbf{r}}^{\prime }){\psi }_{\mathbf{k}{\sigma }^{\prime }}^{MF\ast }({\mathbf{r}}^{\prime })$$

$${\rho }_n=\frac{1}{k_{B}T{L}^2}⟨{\left[{\mathbf{P}}_{\alpha }^{\text{tot}}+\hslash \kappa {\mathbf{S}}_{\alpha }^{\text{tot}}\right]}^2⟩$$

Full text

Introduction. The coupling of artificial gauge fields to ultracold atomic gases Dalibard et al. (2011) has opened the possibility of studying spin-orbit coupled Bose gases Lin et al. (2011); Jiménez-García et al. (2015); Galitski and Spielman (2013); Li et al. (2015) where translational symmetry may be broken spontaneously in the superfluid ground state Li et al. (2016). At the mean-field level, spin-orbit coupling (SOC) introduces degenerate ground states expected to enhance fluctuation effects and giving rise to new, exotic quantum phases. The occurrence and nature of finite temperature transitions in these systems have not yet been fully established Su et al. (2016); Ozawa and Baym (2012a, 2013); Zhai (2015); Jian and Zhai (2011); Liao et al. (2014).

In the following we consider a two-dimensional homogeneous gas of Rashba-Dresselhaus spin-orbit coupled bosons. Mean-field calculations Zhai (2015); Jian and Zhai (2011); Liao et al. (2014); Ho and Zhang (2011); Wang et al. (2010) indicate a Bose condensed ground state of a single plane wave with non-vanishing momentum or a linear superposition of two plane waves with opposite momenta, called plane wave state (PW) and stripe phase (SP), respectively. For spin-independent particle interaction, PW and SP remain degenerate at the mean-field level. In addition, in the case of isotropic SOC, ground states with momenta lying on a circle are connected by symmetry. These degeneracies may be broken by classical or quantum fluctuations.

In this work we explore the phase diagram using classical field Monte Carlo calculations. We show that for anisotropic SOC, the systems undergoes a Kosterlitz-Thouless phase transition from a normal to superfluid state. In the superfluid state, the single particle density matrix decays algebraically and directly reflects the PW/SP character of the mean-field ground state. In the limit of isotropic interparticle interaction, the PW/SP degeneracy is unaffected by the transition. Thus, at large but finite system sizes, fragmentation Mueller et al. (2006) of the condensate occurs. In the case of isotropic SOC, we show that the transition temperature decreases with increasing system size due to the increasing number of degenerate mean-field ground states and eventually vanishes in the thermodynamic limit.

Since continuous phase transitions at finite temperatures are driven by classical long wave length fluctuations, our classical field calculations allows us to establish quantitatively the finite temperature phase diagram in a parameter regime relevant for ultracold atom experiments.

Model. We consider bosons in two hyperfine states, $\sigma$, labeled $\uparrow$ and $\downarrow$, described by the Hamiltonian $H=H_{0}+V$, where $H_{0}$ describes ideal spin-orbit coupled bosons and $V$ the interparticle interaction. We have

[TABLE]

where $\hat{\Psi}^{\dagger}\equiv(\hat{\Psi}_{\uparrow}^{\dagger},\hat{\Psi}_{\downarrow}^{\dagger})$ is the field creation operator, $\sigma_{\alpha}$ are the Pauli matrices, and $m$ is the atomic mass. The SOC is characterized by its strength $\kappa$ and the anisotropy $0\leq\eta_{soc}\leq 1$. For anisotropic SOC, $\eta_{soc}<1$, the minimum of the single particle energies is reached at two wave vectors $(\pm\kappa,0)$, whereas all wave vectors on a circle of absolute value $\kappa$ are degenerate for isotropic SOC $\eta_{soc}=1$.

The interparticle interaction is described by

[TABLE]

where the coupling strengths $g_{\sigma\sigma^{\prime}}$ are determined by the scattering amplitudes of different hyperfine states. For simplification, we only consider $g_{\uparrow\uparrow}=g_{\downarrow\downarrow}$ in the following and use $g=g_{\uparrow\downarrow}-g_{\uparrow\uparrow}$ to characterize the interaction imbalance.

The mean field ground state wave function can be written in terms of

[TABLE]

with $|{{\boldsymbol{\kappa}}}|=\kappa$ and $e^{i\theta_{{\bf k}}}\equiv(k_{x}+i\eta_{soc}k_{y})/\sqrt{k_{x}^{2}+\eta_{soc}^{2}k_{y}^{2}}$. The circular degeneracy of the wave vector is broken for anisotropic SOC where ${{\boldsymbol{\kappa}}}$ points along the $x$ direction. For $g<0$, the mean field energy is minimized by a single plane wave state, $\phi=0$ (PW), whereas stripes corresponding to $\phi=\pi/4$ (SP) are formed in the real space density for $g>0$. The mean-field energy is independent of $\phi$ for isotropic interactions, $g=0$; PW and SP are degenerate at the mean-field level in this case.

Expanding the effective action at low temperatures around the mean-field state at $g\neq 0$ to include thermal fluctuations Jian and Zhai (2011); Liao et al. (2014), a Kosterlitz-Thouless transition is expected for anisotropic SOC, but the occurrence of a superfluid phase for $\eta_{soc}\to 1$ is controversial. Further, in the limit of isotropic interaction, $g=0$, the effective action approach breaks down, since PW and SP become degenerate at the mean-field level and the character of the low temperature phase is still an open problem. To overcome these limitations, we have performed classical field Monte Carlo calculations which are expected to correctly describe the finite temperature behavior close to a possible phase transition where quantum fluctuations should not play an essential role and can therefore be neglected Baym et al. (1999); Holzmann and Baym (2003); Prokof’ev and Svistunov (2002); Holzmann et al. (2007).

We have calculated the reduced single particle density matrix, $G_{\sigma,\sigma^{\prime}}({{\bf r}},{{\bf r}}^{\prime})=\langle\hat{\Psi}_{\sigma^{\prime}}^{\dagger}({{\bf r}}^{\prime})\hat{\Psi}_{\sigma}({{\bf r}})\rangle$ expected to reveal algebraic order for a quasi-Bose condensed state Berezinskiǐ (1972). Superfluid mass density, $\rho_{s}$, can be directly related to the phase stiffness, $\rho_{s}=\partial^{2}F(\theta)/\partial\theta^{2}$, where $F(\theta)$ is the free energy density where the momentum operator, $\hat{p}$, is replaced by $\hat{p}-\theta$ in the Hamiltonian Baym (1969); Pollock and Ceperley (1987); Holzmann and Baym (2007).

The spin-order of the low temperature phase is characterized by the spin-density correlation functions, $M_{\alpha}({{\bf r}})=\langle\hat{{{\bf S}}}_{\alpha}({{\bf r}})\hat{{{\bf S}}}_{\alpha}(0)\rangle$, where $\hat{{{\bf S}}}_{\alpha}({{\bf r}})=\hat{\Psi}^{\dagger}({{\bf r}})\sigma_{\alpha}\hat{\Psi}({{\bf r}})$ is the local magnetization operator. From the mean-field solution, we expect $M_{\alpha}({{\bf r}})$ to monotonically reach a constant value at large distances for PW states when $\alpha$ points in the direction of ${{\boldsymbol{\kappa}}}$. In the stripe phase oscillating behavior occurs, since spins rotate in the plane orthogonal to the vector ${{\boldsymbol{\kappa}}}$ with a spatial frequency equal to $2\kappa$ Ozawa and Baym (2012b).

Method. In the following we are interested to establish the phase diagram in the limit of small interaction strength, $mg\ll 1$ and small spin-orbit coupling, $\kappa\lambda_{T}\ll 1$ where $\lambda_{T}=\sqrt{2\pi\hbar^{2}/mk_{B}T}$ is the thermal wave length at temperature $T$. In this limit, the leading order corrections to mean-field are captured within classical field theory Baym et al. (1999); Holzmann and Baym (2003); Giorgetti et al. (2007); Prokof’ev and Svistunov (2002); Holzmann et al. (2007) where the occupation of low energy modes is high such that commutators like $[\hat{\Psi}^{\dagger}({{\bf r}}),\hat{\Psi}({{\bf r}}^{\prime})]$ can be neglected. In this approximation, the field operator, $\hat{\Psi}({{\bf r}})$, can be replaced by two complex fields, $\Psi({{\bf r}})\equiv(\Psi_{\uparrow}({{\bf r}}),\Psi_{\downarrow}({{\bf r}}))$, one for each spin, and the theory is regularized by discretizing space on a lattice of linear extension $L$.

The probability of a given field configuration is then proportional to $\exp(-S[\Psi({{\bf r}})])$ and the action writes

[TABLE]

where $a$ is the lattice spacing, $\mu$ is the chemical potential, and $\nabla_{D}$ and $\partial_{\alpha}^{D}$ are finite difference expressions approximating the derivatives. It can be shown that the action is real for all field configurations so that the partition function can be sampled by Monte Carlo methods. We have implemented local moves in real and Fourier space to efficiently sample around the mean-field energy minima.

Within classical field theory, the occupation of eigenmodes of energy $\varepsilon$ is by the equipartition theorem, $n_{cf}(\varepsilon)=k_{B}T/\varepsilon$ instead of the full Bose distribution, $n_{B}(\varepsilon)=[\exp(\varepsilon/k_{B}T)-1]^{-1}$. Since we expect mean-field theory to quantitatively describe high energy modes, we have corrected the densities of our classical field calculations to account for the correct ultraviolet behavior Holzmann et al. (2007) adding the difference

[TABLE]

Here, the single particle mean-field energies are given by $\varepsilon_{{{\bf k}}\alpha}^{mf}=\varepsilon_{{{\bf k}}\alpha}+2\sum_{\alpha^{\prime}}g_{\alpha\alpha^{\prime}}n_{\alpha^{\prime}}^{mf}$, where $\varepsilon_{{{\bf k}}\alpha}$ are the eigen energies of the ideal SOC gas. The corresponding mean-field densities, $n_{\alpha}^{mf}=L^{-2}\sum_{{{\bf k}}}n_{cf/B}(\varepsilon_{{{\bf k}}\alpha}^{mf}-\mu)$, have to be determined self-consistently

Results. In order to study the competition between SOC and interparticle interaction, we have fixed $\sum_{\sigma\sigma^{\prime}=\uparrow\downarrow}mg_{\sigma\sigma^{\prime}}/4=\kappa/\sqrt{mk_{B}T}=\pi/20$ with $mg=0$ to address isotropic interaction and $mg=\pm\pi/100$ to slightly break the spin isotropy of scattering particles. In the following we will study the phase diagram as a function of the SOC anisotropy $\eta_{soc}$ and observe signatures of the different phases in the condensate and superfluid fraction and in the spin-resolved density correlation function.

In two dimensions, without SOC, the single particle density matrix, $G_{\sigma,\sigma^{\prime}}({{\bf r}},{{\bf r}}^{\prime})$, is expected to decay algebraically in the low temperature superfluid phase Berezinskiǐ (1972); Kosterlitz and Thouless (1973); Kosterlitz (1974). For SOC bosons, quasi-long range order for SOC bosons occurs when the single particle density matrix is dominated by one or few highly occupied modes. We therefore project $G_{\sigma,\sigma^{\prime}}({{\bf r}},{{\bf r}}^{\prime})$ over all degenerate PW mean-field ground states

[TABLE]

to estimate the condensate fraction $n_{0}^{\kappa}/n$ where $n=\sum_{\sigma}G_{\sigma\sigma}({{\bf r}},{{\bf r}})$ is the total particle density. Although $n_{0}^{\kappa}$ is a direct indicator for a phase transition it does not distinguish PW or SP character.

In Fig. 1 we show $n_{0}^{\kappa}$ as a function of density for a finite system of extension $L/a=80$ with $a=\hbar/\sqrt{mk_{B}T}$. The condensate fraction grows rapidly around a cross-over density which decreases from $\eta_{soc}=0$ to $\eta_{soc}=0.9$ 111Notice that for $\eta_{soc}=0$ the SOC can be eliminated via spin-rotation combined with a gauge transformation whenever ${{\boldsymbol{\kappa}}}$ is commensurate with the boundary conditions. The resulting $\psi^{4}$ theory with $N=4$ field components does not show a Kosterlitz-Thouless phase transition. Here, we address the limit $\eta_{soc}\to 0$ continuously connected to non vanishing $\eta_{soc}>0$, which corresponds to non commensurate values of ${{\boldsymbol{\kappa}}}$.. However, no differences are visible changing the sign of our small anisotropic interaction from negative to positive $g$. We therefore expect that the cross-over/transition temperature is a smooth, continuous function of $g$ around $g=0$.

In order to determine a possible sharp phase transition in the thermodynamic limit, we have determined the behavior of the condensate fraction increasing the system size, $n_{0}^{\kappa}/n\sim L^{-(2-\eta(T))}$. Whereas in the high temperature, normal phase the condensate fraction decreases with the volume, $\eta(T)=0$, the exponent changes in the low temperature phase indicating a Berenzinskii-Kosterlitz-Thouless transition Berezinskiǐ (1972); Kosterlitz and Thouless (1973); Kosterlitz (1974). Assuming the transition to be within the Kosterlitz-Thouless class, the critical temperature can be estimated to occur when $\eta(T_{KT})=1/4$. Our calculations indicate that for anisotropic SOC, $\eta_{soc}<1$ the Berenzinskii-Kosterlitz-Thouless phase occurs at finite temperature, independent of the sign of $g$ (see Fig. 2 for $\eta_{soc}=0$). Further, the limit of isotropic interaction, $g=0$, is indeed approached smoothly from both sides, $g>0$ and $g<0$.

For isotropic SOC, $\eta_{soc}=1$, however, we do not observe the onset of quasi-long range order for $g<0$ in the considered density regime and system sizes. For $g>0$, a cross-over occurs, but the onset of algebraic order strongly depends on the number of degenerate mean-field ground states. Although these energy minima form a circle in the thermodynamic limit, in a finite geometry only a certain number of single particle states are strictly degenerate. As shown in Fig. 3, the behavior of the condensate fraction is qualitatively and quantitatively affected by the number of degenerate states. In particular, the onset of algebraic order is shifted towards considerable higher densities (lower temperatures) increasing the degeneracy from $4$ to $8$ degenerate modes. For $\eta_{soc}=1$ and infinity system sizes, the transition will therefore be shifted to zero temperature and no finite temperature transition with algebraic order in the single particle channel should occur due to the circular degeneracy.

We have further calculated the superfluid and normal mass density, $\rho_{n}=mn-\rho_{s}$ from the phase stiffness

[TABLE]

where ${{\bf P}}_{\alpha}^{\text{tot}}$ is the total momentum and ${{\bf S}}_{\alpha}^{\text{tot}}=\int d{{\bf r}}S_{\alpha}({{\bf r}})$ the total magnetization of the system in the $\alpha=x/y$ direction. Deviations from a Boltzmann distribution of $[{{\bf P}}^{\text{tot}}_{\alpha}+\hbar\kappa{{\bf S}}^{\text{tot}}_{\alpha}]^{2}/(2mnL^{2})$ are directly connected to the quantization of the center of mass motion. Consistent with the prediction of Berenzinskii, Kosterlitz and Thouless, the low temperature, algebraically ordered phase is superfluid for $\eta_{soc}<1$ whereas it vanishes for isotropic SOC with increasing degeneracy. Josephson’s scaling relates the algebraic decay of the single particle density matrix to the condensate and superfluid density Josephson (1966); Holzmann and Baym (2007), and our values of $\eta(T)$ from $n_{0}^{\kappa}$ are consistent with $\eta(T)=m/\rho_{s}\lambda^{2}$ within the numerical uncertainty.

The absence of a transition for isotropic SOC is consistent with rather general considerations particular to two dimensional systems. From the analysis of the non-linear $\sigma$-model Polyakov (1975); Brézin and Zinn-Justin (1976), a phase transition of the Kosterlitz-Thouless type is only expected for $N=2$ component fields, e.g. real and imaginary part of a spinless complex field. For $\eta_{soc}=1$, the PW (SP) order parameter of the infinite system is characterized its phase and by its spin-direction, $-{{\boldsymbol{\kappa}}}/\kappa$, and can be mapped to the two-dimensional (half) unit-sphere. Therefore, vortices are not any more topologically stable objects and destroy algebraic order and superfluidity at any finite temperature. For $\eta_{soc}<1$, the two possible spin-directions are disconnected, so that algebraic order and superfluidity is protected against vortex excitations at low temperatures.

We now want to characterize the order – PW or SP – of the superfluid phase for $\eta_{soc}<1$. Therefore, we diagonalize the single particle density matrix, $G_{\sigma\sigma^{\prime}}({{\bf r}},{{\bf r}}^{\prime})$, calculated by a single realization, not ensemble averaged over different initial conditions. The spin structure of the eigenmodes of $G$ is resolved by taking into account the degenerate Fourier modes with ${{\bf k}}=(\pm\kappa,0)$. As shown in Fig 2, for $g\neq 0$, above the critical density, the single particle density matrix is dominated by a single, highly occupied mode, yielding PW and SP order for $g<0$ and $g>0$, respectively.

The spin structure of this condensate mode is directly reflected in the spin correlation function, $M_{\alpha}({{\bf r}})$, shown in Fig. 4. For $g<0$, $M_{x}(x,0)$ develops quasi-long range order, whereas it remains short ranged in the SP for $g>0$. Quasi-long range stripe order is reflected in slowly decaying oscillations of period $2\kappa$ in $M_{y}(x,0)$. In both cases, the exponent of the algebraic decay is given by the scaling exponent of $n_{0}^{\kappa}$ and compatible with $\eta(T)$ obtained from the superfluid density. Therefore, the quasi-long range spin order results from the spin structure of the underlying quasi-condensate.

In the limit of isotropic interaction, $g=0$, we always obtained two highly occupied modes of the single particle density matrix, degenerate within our numerical precision. Therefore, PW and SP remain degenerate and robust against thermal, critical fluctuations and we obtain a fractionalized quasi-condensate. Both spin correlation functions, $M_{x}(x,0)$ and $M_{y}(x,0)$ become quasi-long ranged and indicate simultaneous PW and SP character. From the Bogoliubov approximation around the $T=0$ mean-field ground states, we expect that quantum fluctuations lift the degeneracy and favor the PW character decreasing the temperature without further phase transition Barnett et al. (2012).

Conclusion. Using classical field simulations we have studied the phase diagram of two-dimensional, spin-orbit coupled Bosons as a function of the spin-orbit anisotropy $\eta_{soc}$ with vanishing or small spin-anisotropy, $g$, in the interparticle interaction. For $\eta_{soc}<1$, our calculations are consistent with a Berenzinskii-Kosterlitz-Thouless transition. The low temperature superfluid phase is characterized by the PW or SP character of the underlying quasi-condensate determined by the sign of $g$. For isotropic interactions, $g=0$, we obtained a fractionalized quasi-condensate with two degenerate modes at the transition showing both, PW and SP character. For fully isotropic SOC, $\eta_{soc}=1$, a cross-over occurs for finite systems, but no superfluid transition is expected in the thermodynamic limit.

Acknowledgments. We thank Gordon Baym for seminal discussions, Frank Hekking and Tomoki Ozawa for useful comments and a critical reading of the manuscript, and Francesco Calcavecchia for continuous support and discussions.