Thermodynamics of spin-orbit coupled bosons in two dimensions from complex Langevin

TL;DR

This paper uses the complex Langevin method to study the thermodynamics of two-dimensional spin-orbit coupled bosons, revealing limitations of mean-field theory and effects of spin-orbit interactions on pseudo-condensates.

Contribution

It introduces a non-perturbative computational approach to analyze the thermodynamics of spin-orbit coupled bosons in two dimensions, addressing the sign problem.

Findings

Mean-field solutions underestimate the density at strong spin-orbit coupling.

Finite volume effects lead to pseudo-condensate formation.

Spin-orbit interactions destroy pseudo-condensates.

Abstract

We investigate the thermal properties of interacting spin-orbit coupled bosons with contact interactions in two spatial dimensions. To that end, we implement the complex Langevin method, motivated by the appearance of a sign problem, on a square lattice with periodic boundary conditions. We calculate the density equation of state non-perturbatively in a range of spin-orbit couplings and chemical potentials. Our results show that mean-field solutions tend to underestimate the average density, especially for stronger values of the spin-orbit coupling. Additionally, the finite nature of the simulation volume induces the formation of pseudo-condensates. These have been observed to be destroyed by the spin-orbit interactions.

Equations

$$S=\int d^{2}xd\tau \mathcal{L},$$

$${\mathcal{L}}_0=$$

$$-{\Phi }^{†}\frac{(-i\text{openone}{\partial }_x-{\kappa }_x{\sigma }_x{)}^2+(-i\text{openone}{\partial }_y-{\kappa }_y{\sigma }_y{)}^2}{2m}\Phi ,$$

$${\mathcal{L}}_{\mathrm{int}}=\frac{\lambda }{4}\left[{\left({\varphi }_{\uparrow }^{†}{\varphi }_{\uparrow }\right)}^2+{\left({\varphi }_{\downarrow }^{†}{\varphi }_{\downarrow }\right)}^2\right]+\frac{g}{2}\left({\varphi }_{\uparrow }^{†}{\varphi }_{\uparrow }{\varphi }_{\downarrow }^{†}{\varphi }_{\downarrow }\right),$$

$$S=$$

$$-\frac{1}{2m}\sum _j{\Phi }_{(\vec{x},\tau )}^{†}\left[v_j{\Phi }_{(\vec{x}+a\hat{j},\tau )}+v_j^{†}{\Phi }_{(\vec{x}-a\hat{j},\tau )}-2{\Phi }_{(\vec{x},\tau )}\right]$$

$$+\frac{\lambda +g}{8}{\left({\Phi }_{(\vec{x},\tau )}^{†}{\Phi }_{(\vec{x},\tau -a_{\tau })}\right)}^2$$

$$+\frac{\lambda -g}{8}{\left({\Phi }_{(\vec{x},\tau )}^{†}{\sigma }_z{\Phi }_{(\vec{x},\tau -a_{\tau })}\right)}^2\},$$

$$S_{\mathrm{free}}=\sum _{\vec{p},\omega }{\tilde{\Phi }}_{(\vec{p},\omega )}^{†}M(\vec{p},\omega ,\mu ,\vec{\kappa }){\tilde{\Phi }}_{(\vec{p},\omega )},$$

$$M=\left[\begin{array}{cc}T+v_x^{+}+v_x^{-}+v_y^{+}+v_y^{-}& v_x^{+}-v_x^{-}-i(v_y^{+}-v_y^{-})\\ v_x^{+}-v_x^{-}+i(v_y^{+}-v_y^{-})& T+v_x^{+}+v_x^{-}+v_y^{+}+v_y^{-}\end{array}\right],$$

$$\left[\begin{array}{c}{\phi }_{+}\\ {\phi }_{-}\end{array}\right]=\left[\begin{array}{cc}\frac{(v_x^{+}-v_x^{-})+i(v_y^{+}-v_y^{-})}{\sqrt{2[(v_x^{+}-v_x^{-}{)}^2+(v_y^{+}-v_y^{-}{)}^2]}}& \frac{1}{\sqrt{2}}\\ -\frac{(v_x^{+}-v_x^{-})+i(v_y^{+}-v_y^{-})}{\sqrt{2[(v_x^{+}-v_x^{-}{)}^2+(v_y^{+}-v_y^{-}{)}^2]}}& \frac{1}{\sqrt{2}}\end{array}\right]\left[\begin{array}{c}{\varphi }_{\uparrow }\\ {\varphi }_{\downarrow }\end{array}\right],$$

$${\lambda }_{\pm }(\vec{p},\omega ,\mu ,\vec{\kappa })$$

$$\pm \sqrt{(v_x^{+}-v_x^{-}{)}^2+(v_y^{+}-v_y^{-}{)}^2}.$$

$$\beta {\Omega }_{\mathrm{free}}(\mu ,\vec{\kappa })=-\ln {\mathcal{Z}}_{\mathrm{free}}$$

$$=\sum _{\vec{p},\omega }[\ln {\lambda }_{+}(\vec{p},\omega ,\mu ,\vec{\kappa })+\ln {\lambda }_{-}(\vec{p},\omega ,\mu ,\vec{\kappa })].$$

$$⟨n⟩=\frac{1}{V}\frac{\partial \ln \mathcal{Z}}{\partial (\beta \mu )}=\sum _{s=\pm }\sum _{\vec{p},q}\frac{1}{1-X_s(\vec{p},\mu )e^{-i{\omega }_q}},$$

$$X_s(\vec{p},\mu )=e^{-\xi \mu }(1+{\lambda }_s-T).$$

$$⟨n⟩=\frac{1}{V}\sum _{s=\pm }\sum _{\vec{p}}\frac{1}{[X_s(\vec{p},\mu ){]}^{N_{\tau }}-1}.$$

$$[X_s(\vec{p},\mu ){]}^{N_{\tau }}\to e^{\beta ({ϵ}_s(\vec{p},\vec{\kappa })-\mu )}$$

$${ϵ}_{\pm }(\vec{p},\vec{\kappa })$$

$$\pm \frac{|\vec{p}||\vec{\kappa }|}{m}\sqrt{{\cos }^2({\theta }_{\kappa }){\cos }^2({\theta }_p)+{\sin }^2({\theta }_{\kappa }){\sin }^2({\theta }_p)},$$

$${\lambda }_{\pm }$$

$$\frac{\partial {\varphi }_s(x,\tau )}{\partial \theta }=-\frac{\delta S}{\delta {\varphi }_s(x,\tau )}+{\eta }_s(x,\tau ),$$

$$⟨{\eta }_s(x,\tau ){⟩}_{\eta }$$

$$⟨{\eta }_s(x,\tau ){\eta }_{s^{\prime }}(x^{\prime },{\tau }^{\prime }){⟩}_{\eta }$$

$$⟨\hat{O}⟩=\underset{\theta \to \infty }{\lim }⟨\hat{O}({\varphi }_{\uparrow }(\theta ),{\varphi }_{\downarrow }(\theta )){⟩}_{\eta },$$

$$⟨n⟩$$

$$=\frac{e^{\xi \mu }}{V}\sum _{x,\tau }⟨{\Phi }_{(x,\tau )}^{†}{\Phi }_{(x,\tau -a_{\tau })}⟩,$$

$$⟨n_s⟩=\frac{2e^{\xi \mu }(e^{\xi \mu }-1-\xi /m(\cos ({\kappa }_x)+\cos ({\kappa }_y)-2))}{\xi \lambda },$$

$$P(\mu )-P({\mu }_0)={\int }_{{\mu }_0}^{\mu }⟨n({\mu }^{\prime })⟩d{\mu }^{\prime }.$$

$$G_{s{s}^{\prime }}\left(|\vec{x}-\vec{y}|\right)=⟨{\varphi }_{s(\vec{x},\tau )}^{†}{\varphi }_{s^{\prime }(\vec{y},\tau )}⟩.$$

$$R_{s{s}^{\prime }}=\frac{G_{s{s}^{\prime }}(a{N}_x/2)}{G_{s{s}^{\prime }}(0)},$$

Full text

NT@UW-19-10, INT-PUB-19-035

I Introduction

The experimental realization of ultracold atomic systems with spin-orbit coupling (SOC), nearly a decade ago lin_spinorbit-coupled_2011 ; wang_spin-orbit_2012 ; cheuk_spin-injection_2012 , opened an exciting new set of directions for the exploration of the properties of matter in extreme (yet highly controllable) conditions. The SOC, which couples the atomic pseudo-spin (which itself denotes the particle species or ‘flavor’) to momentum, is realized as the coupling of conventional nonrelativistic neutral particles to a synthetic non-abelian background gauge field lin_synthetic_2009 ; lin_synthetic_2011 ; Dalibard:2010ph . Such a construction has potential applications for the exploration of a variety of physical situations (including Rashba- ISI:A1960WT33500031 and Dresselhaus-type PhysRev.100.580 couplings). Notably, these systems have been under intense scrutiny in recent years as they may be used experimentally to realize exotic phases of matter such as supersolids wu_realization_2016 ; leonard_supersolid_2017 ; 2017Natur.543…91L , superfluids with a crystalline structure, and study exotic topological properties 2016NatPh..12..540H ; PhysRevLett.117.235304 ; 2019arXiv190708637V .

More generally, ultracold bosons subject to SOC are known, at mean field level, to exhibit stripe or plane wave phases PhysRevLett.105.160403 ; Zhai:2014gna : the ground state wave function in the former is composed of two plane waves propagating in opposite directions, leading to an interference pattern, while on the latter it has only one plane wave. Many recent studies have investigated mean field properties of different types of SOC 2017PhRvA..95e1601K , and also first order (one loop) quantum corrections 2018PhRvL.120n0403L ; 2019arXiv190308182L . Recent theoretical studies of spin-orbit coupling include its effect combined with finite angular momentum 2017PhRvL.118n5302S and a harmonic trap 2019arXiv190713355M .

In this work we characterize, in a non-perturbative fashion, some of the basic thermodynamic observables of spin-orbit coupled bosons with quartic interactions. Using the complex Langevin method, we perform a study of such a non-relativistic Bose gas on a spacetime lattice, determining the density equation of state, pressure and pseudo-condensate fraction. Notably, the lattice formulation compactifies the spin-orbit interaction, which appears as a constant background non-abelian gauge field. This allows for an easy study of the limit of very large SOC compared to the momentum.

II Model and lattice formulation

We consider a system of $(2+1)$-dimensional non-relativistic bosons with two hyperfine (pseudo-spin) states, denoted by $\uparrow$ and $\downarrow$, in Euclidean spacetime. They are subject to a Rashba-Dresselhaus spin-orbit coupling (SOC) and density-density contact interactions. More specifically, we will use the Euclidean action

[TABLE]

where $\mathcal{L}=\mathcal{L}_{0}+\mathcal{L}_{\mathrm{int}}$,

[TABLE]

is the free part, and

[TABLE]

is the interacting part. In the above equations, $\Phi=[\phi_{\uparrow}\quad\phi_{\downarrow}]^{T}$, $\mu$ is the chemical potential for both species, $m$ is the mass, $\lambda$ and $g$ are the intra- and inter-species couplings, and $\kappa_{x}$ and $\kappa_{y}$ characterize the spin-orbit interaction coupled to the $\sigma_{x}$ and $\sigma_{y}$ Pauli matrices. Note that, by virtue of the SOC, neither particle number is individually conserved, but the total particle number is, such that $\mu$ is a sensible chemical potential.

In order to perform our lattice studies, we discretise the system in a hypercubic lattice of size $N_{x}^{2}\times N_{\tau}$ and spacings $a$ and $a_{\tau}$ in the spatial and temporal directions, respectively. The temperature is given by $T=\beta^{-1}=1/a_{\tau}N_{\tau}$. We assume periodic boundary conditions in all directions. The SOC enters the action in the same way as a background SU($2$) gauge field, similar to a minimal coupling. On the lattice, the SOC is treated as a background non-abelian gauge field and it is discretized in the same way. The discretized version of the action, with the fields and couplings being rescaled by appropriate powers of the lattice spacing, is given by

[TABLE]

where $\hat{j}$ represents a unit vector in the $j$-th direction, $v_{j}=e^{-i\kappa_{j}\sigma_{j}}$, and $\xi=a_{\tau}/a$ is the lattice spacing anisotropy factor. The contact interactions have been regularized in the same way as the number density operator. This formulation displays explicitly the global SU($2$) flavor symmetry of the contact interactions when $\lambda=g$, and the conservation of the total particle number by all interactions due to the global U($1$) symmetry.

In this work, we study the interplay between the self-coupling $\lambda$ and the spin-orbit couplings $\kappa_{x}$ and $\kappa_{y}$. The coupling $g$ between different pseudo-spins is left for a future publication.

II.1 Exact solution in the quadratic case

When the quartic terms are not present, the action (hereafter referred to as $S_{\mathrm{free}}$) can be written in momentum space as

[TABLE]

with

[TABLE]

where $v_{j}^{\pm}=\xi\sin^{2}((p_{j}\pm\kappa_{j})/2)/m$, and $T=(1-e^{\xi\mu}e^{i\omega_{q}})$. The lattice momenta and Matsubara frequencies are given by $p_{j}=2\pi k_{j}/N_{x}$, with $k_{j}=0,\ldots,N_{x}-1$, and $\omega_{q}=2\pi q/N_{\tau}$ with $q=0,\ldots,N_{\tau}-1$, respectively. The tildes represent Fourier transforms. The matrix $M$ can be diagonalized via a change to the helicity basis:

[TABLE]

with eigenvalues

[TABLE]

leading to the grand thermodynamic potential $\Omega$ via

[TABLE]

It is clear from the above equation that, since $p_{x}$ and $p_{y}$ run over the same interval, the free energy is symmetric under $\kappa_{x}\leftrightarrow\kappa_{y}$. Moreover, the partition function on the lattice is an even, periodic function of the spin-orbit couplings, with period $\pi$ (in lattice units). This periodicity is a lattice artifact that disappears in the continuum limit.

The average density of the noninteracting case can be obtained by differentiation with respect to $\beta\mu$:

[TABLE]

where

[TABLE]

It can be shown that the sum over the Matsubara frequencies can be carried out, which yields

[TABLE]

In the $\xi\to 0$ limit,

[TABLE]

where, in the continuum limit,

[TABLE]

with $\theta_{p}=\tan^{-1}(p_{y}/p_{x})$ and $\theta_{\kappa}=\tan^{-1}(\kappa_{y}/\kappa_{x})$.

The continuum eigenvalues are given by

[TABLE]

Notice that, for $|\vec{\kappa}|=0$, $\lambda_{\pm}$ vanishes at zero momentum and $\mu=0$, signaling the well known instability associated with Bose-Einstein condensation. For $|\vec{\kappa}|\neq 0$, the instability shifts to $\mu\neq 0$ and exists for both $|\vec{p}|=0$ and $|\vec{p}|\neq 0$. The continuum eigenvalues have been studied in three dimensions in the Hamiltonian formulation for isotropic spin-orbit coupling in Ref. Ozawa:2011kj . Additionally, one can see that $N_{\pm}=\varphi_{\pm}^{\dagger}\varphi_{\pm}$ are conserved, and that $N_{+}+N_{-}=N_{\uparrow}+N_{\downarrow}$.

Using trigonometric identities it is possible to show that $v^{+}_{i}+v^{-}_{i}=\xi(1-\cos(p_{i})\cos(\kappa_{i}))/m$ and $v^{+}_{i}-v^{-}_{i}=\xi\sin(p_{i})\sin(\kappa_{i})/m$, from which it is clear that $\kappa_{i}=\pm\pi/2$ corresponds to the limit of $|\kappa_{i}|\gg|p_{i}|$ in the continuum. When $|\vec{\kappa}|\gg|\vec{p}|$, such that the $|\vec{p}|^{2}$ term can be ignored, the single particle Hamiltonian becomes an anisotropic Weyl Hamiltonian with $\vec{v}_{0}=\vec{\kappa}/m$ playing the role of the (anisotropic) speed of light and an effective chemical potential of $-m|\vec{v}_{0}|^{2}/2$. A similar case has been discussed, in the context of ultracold fermionic atoms, in Ref. 2008PhRvA..77a1802J .

III Many-body method

The first order time derivative in the action is a non-Hermitian operator, making $e^{-S}$ a complex weight for the path integral; this is known as phase (or sign) problem. This prevents the use of traditional Monte Carlo methods, since they use $e^{-S}$ as a probability weight. One alternative in this scenario is the complex Langevin technique, which has been used to study theories with sign problems such as those with repulsive interactions Loheac:2017yar as well as polarized Loheac:2018yjh ; Rammelmuller:2018hnk and mass-imbalanced fermions Rammelmuller:2017vqn (see Ref. Berger:2019odf for a review), finite density QCD with staggered quarks Sexty:2013ica ; Attanasio:2018rtq ; Kogut:2019qmi ; Sexty:2019vqx , random matrix models Nagata:2016alq ; Bloch:2017sex , rotating bosons Hayata:2014kra ; Berger:2018xwy , superstring-inspired matrix models Nishimura:2019qal , among others.

The complex Langevin method is an extension of stochastic quantisation Parisi:1980ys . The latter method consists of evolving the fields along a fictitious time dimension, $\theta$, according to the Langevin equation

[TABLE]

where $\eta_{s}(x,\tau)$ is a Gaussian white noise field satisfying

[TABLE]

with $\langle\cdot\rangle_{\eta}$ indicating an ensemble average over the noise field. Quantum expectation values are obtained as

[TABLE]

where $O$ is some observable. In practice, the Langevin equations are solved numerically with a step size $\varepsilon>0$, chosen adaptively Aarts:2009dg . We follow an Euler-like discretization scheme in this work. This generates a sequence of field configurations. Ensemble averages are performed as simple averages of the observables calculated using the configurations generated after the system reaches its steady state.

In order to deal with theories that have a complex action, each of the fields has to be complexified Parisi:1984cs ; Klauder:1983nn ; Klauder:1983zm ; Klauder:1983sp ; Aarts:2008rr ; Aarts:2008wh ; Aarts:2009uq . For complex fields, both the real and imaginary parts become complex and obey the Langevin equation (16). We choose the noise to remain real Aarts:2011ax . Expectation values are calculated in the same way as in the case with real action. Note that real observables do become complex in this method, but their imaginary parts should be statistically compatible with zero.

The average density was calculated as

[TABLE]

where the angular brackets on the right-hand side indicate an average over configurations generated by the Langevin process. Mean-field results can be obtained by solving the Langevin equation without noise, which finds the minimum of the action.

We have performed our simulations on a periodic lattice of volume $20^{2}\times 64$, a spacing anisotropy of $\xi=1/8$, and mass $m=1$. These parameters lead to a thermal wavelength of $\sim 7$ in lattice units, which is consistent with the continuum limit window $1\ll\lambda_{T}/a\ll N_{x}$. The Langevin step size was chosen adaptively, with average of $\mathcal{O}(10^{-4})$. We have estimated the auto-correlation time via the method proposed in Wolff:2003sm .

IV Results

IV.1 Quadratic case

We have studied the case of non-interacting ($\lambda=0$) exactly on the lattice. The average particle number density can be seen in figure 1.

One can see that the density has its minimum when the system becomes Weyl-like, i.e., when the $p^{2}$ term in the Hamiltonian is much smaller than $\vec{\sigma}\cdot\vec{p}$ and can be neglected.

At $\beta\mu\geq 0$ and $\kappa_{x}=\kappa_{y}=0$ the condensation of bosons in the ground state makes the simulation unstable. This condensation is sharper at $T=0$ and softer at $T>0$, which is our case. For non-zero spin-orbit couplings, however, the chemical potential where such condensation happens is pushed to larger values due to the $\kappa^{2}$ term behaving as an effective chemical potential. To better visualise this, we have calculated the average density as a function of $\kappa_{x}$ and $\kappa_{y}$ by including only their effect on the chemical potential, i.e., the terms proportional to $\cos(\kappa_{j})$, and ignoring the $\sin(\kappa_{j})$ terms ($\vec{\sigma}\cdot\vec{p}$ in the continuum). The result is shown in Fig. 2.

The figure further shows that the $\vec{\sigma}\cdot\vec{p}$ term in the continuum action has a non-trivial effect on the average density.

IV.2 Interacting case I – isotropic SOC

As an initial check on the ability of the complex Langevin method to give correct results, we have looked at the imaginary part of the density. Despite the sign problem in the Euclidean formulation, the density is expected to be (statistically compatible with) zero. This can be verified in figure 3. Throughout this section we have used $\kappa_{x}=\kappa_{y}\equiv\kappa$. Complex-valued densities would indicate a failure of the simulations.

The repulsive quartic couplings, similarly to the spin-orbit couplings, have the effect of keeping the system stable at small but positive chemical potentials. This can be seen in Fig. 4, where we show the average density as function of $\beta\mu$ for different values of (isotropic) SOC at $\lambda/a=0.5$. As in the case of $\lambda=0$, indicated by the long-dashed, short-dashed, dotted and dash-dotted lines for each $\kappa$, we observe smaller densities as the spin-orbit coupling increases.

The mean field average density is also shown in the figure. For each flavour, it is given by

[TABLE]

and is zero when the right-hand side is negative.

For the larger values of $\kappa$ shown in the figure, the mean field density becomes positive at much higher values of $\beta\mu$. The mean field average density for $a\kappa=0.4\pi/2$ and $0.6\pi/2$ is very small for the chemical potentials considered and not shown on the figure. It is clear that the average density is not well described by the mean field result for $\mu\leq 0$ and/or $\kappa>0$.

In order to have a better look at the effect of the SOC over the bosonic system, we show in Fig. 5 how the average density changes as we vary $\kappa$. For comparison, we also show the mean field density at $\beta\mu=0.4$ and non-interacting results at $\beta\mu=-0.4$ and $\beta\mu=0.0$.

In all cases, an increase in the spin-orbit coupling has led to decreasing average densities as the system gets closer to the Weyl-like state.

We have also calculated the pressure difference from a reference value via the Gibbs-Duhem relation

[TABLE]

The numerical integration has been carried out using the trapezoid rule, with $\beta\mu_{0}=-1.6$ as reference point. We have estimated the errors via bootstrapping with $1000$ samples. The results are shown in figure 6. The pressure difference in the plot is shown in units of $1/\beta\lambda_{T}^{2}$.

The confined nature of the simulation volume can induce the formation of (pseudo-)condensates at finite temperature. These condensates manifest themselves as off-diagonal long range order in the correlation function between similar spins. The correlation function is, in general, given by

[TABLE]

We show, in figure 7, the condensate fraction,

[TABLE]

between two pseudo-spin “up” fields, as a function of $\kappa$ and different values of $\beta\mu$. A similar result is obtained for two “down” spins.

At $\kappa=0$ we observe a finite condensate fraction, which tends to zero as $\kappa$ increases, and is larger for higher values of $\beta\mu$. Beyond $a\kappa>0.1\pi$, our results for $R$ are compatible with zero and thus excluded from the plot. The correlations between fields of different spins have been measured to be statistically compatible with zero for all values of $\kappa$ and $\mu$ considered.

The above results indicate a destructive interplay between the spin-orbit coupling and condensation in finite systems in two spatial dimensions. A similar phenomenon has been observed in three dimensions with an inter-species coupling in 2018PhRvL.120n0403L , and isotropic $s$-wave coupling PhysRevLett.109.025301 .

IV.3 Interacting case II – anisotropic SOC

We have investigated the effects of anisotropic spin-orbit coupling, i.e., $\kappa_{x}\neq\kappa_{y}$, on the density equation of state by using $\kappa_{y}=\eta_{\mathrm{soc}}\kappa_{x}$, with $0\leq\eta_{\mathrm{soc}}\leq 1$. Because of the symmetry $\kappa_{x}\leftrightarrow\kappa_{y}$ in the partition function there is no need to consider $\eta_{\mathrm{soc}}>1$.

In figure 8 we plot the average density as function of $|\vec{\kappa}|$. The $x$-axis has been normalized so that the maximum value of $|\vec{\kappa}|$ is the same for all anisotropies. We observe a slower decay of $\langle n\rangle$ as a function of $|\vec{\kappa}|$ for $\eta_{\mathrm{soc}}<1$. We remind the reader that the mean field densities for negative chemical potentials are zero.

At $\eta_{\mathrm{soc}}=0.0$ there is a periodic behaviour of period $2\pi/N_{x}$ on the average density for both non-interacting and interacting cases. In the former, we have verified it to be a property of the partition function as a whole, and not of the eigenvalues. Some remnant of this periodicity can be seen at $\eta_{\mathrm{soc}}=0.5$ for small values of $|\vec{\kappa}|$.

The density equation of state for different SOC anisotropies is shown in figure 9 for both interacting (points) and non-interacting (lines) cases. We have used $a\kappa_{x}=0.25\pi/2$. Mean field results are very small in comparison to those in the figure and therefore omitted.

The figure shows that the distincion between the interacting and non-interacting equation of states becomes more pronounced as $\eta_{\mathrm{soc}}$ decreases, and is stronger for larger chemical potentials.

As with the isotropic spin-orbit case, we have investigated the condensate fractions for $\eta_{\mathrm{soc}}<1$. Figure 10 shows the results for $R_{\uparrow\uparrow}$ as a function of $|\vec{\kappa}|$. Similar to the $\eta_{\mathrm{soc}}=1$ case, $R_{\downarrow\downarrow}$ displayed a similar behaviour. We observe a slower decay of both of the (pseudo)-condensate fraction for larger SOC anisotropies. This effect has been seen in three dimensions in PhysRevA.87.031604 .

V Summary and outlook

We have investigated the effects of an artificial spin-orbit coupling on the density equation of state for a bosonic system of two pseudo-spins, as well as contact interactions between similar boson species using non-perturbative numerical methods. The action is complex in Euclidean spacetime due to the first order time derivative, and therefore our method of choice for the simulations was the complex Langevin technique.

Due to how the spin-orbit term enters the lattice formulation, it was possible to investigate how the equation of state changes in the range $0\leq\kappa_{x},\kappa_{y}<\infty$. In particular, the Hamiltonian becomes Weyl-like when the spin-orbit coupling is much larger than the momentum.

We have obtained the density and pressure equations of state for different values of the spin-orbit coupling in the isotropic case, where $\kappa_{x}=\kappa_{y}\equiv\kappa$. The average number density has been seen to be a decreasing function of $\kappa$, having its minimum value when $a\kappa=\pi/2$ (which corresponds to $\kappa\to\infty$ in the continuum). A comparison with mean field results has shown that quantum effects play a bigger role for larger values of $\kappa$ and positive chemical potentials. For $\mu\leq 0$ the mean field average density is zero, in clear contrast with the simulations, which include all quantum effects.

We have also investigated the case of anisotropic spin-orbit coupling, with $\kappa_{y}=\eta_{\mathrm{soc}}\kappa_{x}$. A periodic behaviour of the number density as a function of the spin-orbit coupling, induced by the finite volume, has been observed in both interacting and non-interacting cases when the SOC anisotropy $\eta_{\mathrm{soc}}=0$. As $\eta_{\mathrm{soc}}$ is increased from [math] towards $1$, the aforementioned decaying behaviour of the density is recovered.

For both isotropic and anisotropic spin-orbit couplings, we have carried out comparisons between the interacting and non-interacting equations of state. The distinction between them is smallest for large SOC, and becomes more apparent as $|\vec{\kappa}|\to 0$. Moreover, the finite volume of the simulations allows the formation of a pseudo-condensate, which is depleted by the spin-orbit coupling.

Possible future studies include non-perturbative investigations of the interplay between the spin-orbit coupling and rotation, as well as different types of SOC PhysRevLett.107.150403 . Moreover, the determination of physical quantities such as the scattering length or binding energy, which help connecting with experimental results, can be done via Lüscher’s method luscher1986 or the second virial coefficient Shill:2018tan .

VI Acknowledgements

The authors would like to thank T. Ozawa for useful discussions. This work was facilitated though the use of advanced computational, storage, and networking infrastructure provided by the Hyak supercomputer system at the University of Washington. The work of FA was supported by US DOE Grant No. DE-FG02-97ER-41014. This material is based upon work supported by the National Science Foundation under Grant No. PHY1452635 (Computational Physics Program). F.A and acknowledges the hospitality of the University of North Carolina at Chapel Hill where part of this work was performed.