Dynamics of a quantum phase transition in a ferromagnetic Bose-Einstein condensate

TL;DR

This paper investigates the dynamics of a slow quantum phase transition in a spin-1 Bose-Einstein condensate, analyzing how the system's magnetization scales during the transition through analytical and numerical methods.

Contribution

It provides an analytical derivation of the scaling laws for magnetization during a quantum phase transition in a spinor BEC, supported by numerical simulations.

Findings

Analytical scaling laws for magnetization dynamics

Numerical verification in a 1D model

Insights into slow quantum phase transition behavior

Abstract

We discuss dynamics of a slow quantum phase transition in a spin-1 Bose-Einstein condensate. We determine analytically the scaling properties of the system magnetization and verify them with numerical simulations in a one dimensional model.

Figures

Equations

$$E[\Psi ]=\int$$

$$q=Q/(n|c_1|),n={\Psi }^{†}\Psi$$

$${\Psi }_B^T\sim (\sqrt{4-2q}{e}^{i{\chi }_1},\sqrt{8+4q}{e}^{i({\chi }_1+{\chi }_{-1})/2},\sqrt{4-2q}{e}^{i{\chi }_{-1}})$$

$$q(t)=2-t/{\tau }_Q,$$

$$\hat{t}\sim {\tau }_Q^{1/3}.$$

$$f_{\alpha }=⟨\Psi |F_{\alpha }|\Psi ⟩,\alpha =x,y,z.$$

$$M_T(t)=\int dz[f_x^2(z,t)+f_y^2(z,t)]=\int dz{m}_T,$$

$$\hat{\xi }\sim {\tau }_Q^{1/3}.$$

$$i\frac{d}{dt}{G}_{\chi }=-\frac{1}{2}\frac{d^2}{d{z}^2}{G}_{\chi }+\frac{\alpha }{2}q{G}_{\chi }-\frac{\alpha }{2}(G_{\chi }+G_{\chi }^{\ast }),$$

$$\frac{d}{dt}\left[\begin{array}[]{c}a_{\chi}\\ b_{\chi}\end{array}\right]=\frac{1}{2}\left(\begin{array}[]{cc}0&k^{2}+\alpha q\\ 2\alpha-k^{2}-\alpha q&0\end{array}\right)\left[\begin{array}[]{c}a_{\chi}\\ b_{\chi}\end{array}\right].$$

$$a_{\chi }(k,t)={\alpha }_{k\chi }\mathrm{Ai}(s)+{\beta }_{k\chi }\mathrm{Bi}(s),\frac{s}{\kappa }=\frac{t}{{\tau }_Q^{1/3}}-\frac{k^2{\tau }_Q^{2/3}}{\alpha },$$

Full text

Studies of phase transitions have traditionally focused on equilibrium scalings of various properties near the critical point. Dynamics of the phase transition presents new challenges and there is a strong motivation for analyzing it. Nonequilibrium phase transitions may play a role in the evolution of the early Universe kibble . Their analogues can be studied in the condensed matter experiments. The latter observation led to development of the theory based on the universality of critical behavior zurek , which in turn resulted in a series of beautiful experiments eksperymenty . The recent progress in the cold atom experiments allows for time dependent realizations of different models undergoing a quantum phase transition (QPT) lewenstein ; nature_kurn . These experimental developments are only a proverbial tip of the iceberg, but they call for an in-depth theoretical understanding of the QPT dynamics.

A QPT is a fundamental change in ground state (GS) of the system as a result of small variations of an external parameter, e.g., a magnetic field sachdev . It takes place ideally at zero absolute temperature, which is in striking contrast to thermodynamical phase transitions. The most complete description of the QPT dynamics has been obtained so far in spin models dorner ; spiny that are exactly solvable. In these systems the gap in the excitation spectrum goes to zero at the critical point, which precludes the adiabatic evolution across the phase boundary. It leads to creation of excitations whose density and scaling with a quench rate follow from a quantum version dorner ; bodzio of the Kibble-Zurek (KZ) theory kibble ; zurek .

We study dynamics of a ferromagnetic condensate of spin-1 particles fer . For simplicity, we consider 1D homogeneous (untrapped) clouds: atoms in a box as in the experiment raizen with spinless bosons. We adopt the parameters for our 1D model such that the length and time scales are comparable to experimental ones units . Assuming that the system is placed in a magnetic field $B$ aligned in the $z$ direction, one gets the following dimensionless mean-field energy functional units

[TABLE]

where $\Psi^{T}=(\psi_{1},\psi_{0},\psi_{-1})$ describes the $m=0,\pm 1$ condensate components, $\int dz\sum_{m}|\psi_{m}|^{2}=1$, and $F_{x,y,z}$ are spin-1 matrices ueda_broken . The first term in (1) is the kinetic energy, the second and the fourth term describe spin-independent and spin-dependent atom interactions respectively, the third term is a quadratic Zeeman shift coming from atom interactions with a magnetic field. For 87Rb atoms considered here $c_{1}<0$, which results in an interesting phase diagram due to the competition between the last two terms in (1). Restricting analysis to zero longitudinal magnetization case, and introducing

[TABLE]

one finds a polar phase for $q>2$, described by $\Psi_{P}^{T}\sim(0,1,0)$, and the broken-symmetry phase where

[TABLE]

for $0\leq q<2$. The freedom of choosing the $\chi_{\pm 1}$ results in rotational symmetry of the transverse magnetization on the $(x,y)$ plane. The transition between these phases can be driven by the change of the magnetic field $B$ imposed on the atom cloud, $q\sim Q\sim B^{2}$ ueda_pra2007 , which was experimentally done in nature_kurn .

The dynamics of a QPT depends on the rate of quench driving the system across the phase boundary. For very fast “impulse” transition, the system has no time to adjust to the changes of the Hamiltonian and arrives in a region where a new phase is expected with the “old” wave-function untouched during the evolution. Slow transitions are different: the system has time to “probe” various broken symmetry “vacua” in the neighborhood of the critical point where it gets excited. We are interested in evolutions that are fast enough to produce macroscopic excitations of the system, but slow enough to reflect scalings of the critical region. By comparing analytical findings to numerical simulations for experimentally relevant parameters we provide the first complete description of QPT dynamics in a ferromagnetic condensate.

Fast transitions were realized in the Berkeley experiment nature_kurn . The 3D numerical simulations closely following this experiment were reported in ueda_pra2007 . Analytical studies of the evolution after “impulse” quench were presented in austen ; uwe . The paper of Lamacraft austen also discusses dynamics of non instantaneous transitions in 2D spinor condensates focusing on analytical predictions on the growth of the transverse magnetization correlation functions.

We start with a qualitative discussion. Considering small perturbations around the GS of the broken-symmetry phase one finds three Bogolubov modes as in ueda_broken where quantum fluctuations are studied. In the long wavelength limit (important for slow transitions) there is only one nonzero eigenvalue mode: the gapped mode having eigenenergy $\Delta\sim\sqrt{4-q^{2}}$. Suppose now that we drive the system from polar to broken-symmetry phase. The system reaction time to Hamiltonian changes in the broken-symmetry phase is given by the inverse of the gap: $\tau\sim\frac{1}{\Delta}$ dorner ; bodzio . For example, when $\tau$ is small enough the evolution becomes adiabatic so the system adjusts fast to parameter changes. Right after entering the broken-symmetry phase, the reaction time is large with respect to the transition time, $\Delta/\frac{d\Delta}{dt}$, and so the system undergoes the “impulse” evolution where its state is “frozen”. The gapped mode starts to be populated around the instant $\hat{t}$ after entering the broken-symmetry phase: the system leaves the “impulse” regime to catch up with instantaneous GS solution. This happens when the two time scales become comparable: $1/\Delta(\hat{t})\sim\Delta/\frac{d\Delta}{dt}|_{t=\hat{t}}$. We consider here transitions driven by

[TABLE]

where $\tau_{Q}$ is the quench time inversely proportional to the speed of driving the system through the phase transition. For slow transitions of interest here, $\tau_{Q}\gg 1$, we obtain

[TABLE]

In the following we analyze dynamics induced by a linear decrease of $q(t)$ (2). The evolution starts from $t<0$, i.e., in the polar phase, and ends at $t=2\tau_{Q}$ ($q=0$). Such $q(t)$ dependence is achieved by ramping down the magnetic field as $\sim\sqrt{2-t/\tau_{Q}}$. The initial state is chosen as a slightly perturbed GS in the polar phase, $\Psi^{T}\sim(\delta\psi_{1},1/\sqrt{L}+\delta\psi_{0},\delta\psi_{-1}),$ where $|\delta\psi_{m}|\ll 1/\sqrt{L}$ are random. We generate the real and imaginary part of $\delta\Psi_{m}$ at different grid points with the probability distribution $p(x)=\exp(-x^{2}/2\sigma^{2})/\sqrt{2\pi}\sigma$. We take $\sigma=10^{-4}$ to start evolution closely to the polar phase GS.

To find the full numerical solution within the mean-field approximation, we integrate three coupled nonlinear Schrödinger equations for the $\psi_{m}$ condensates that can be easily obtained by the variation of (1). During evolution we look at the magnetization of the sample

[TABLE]

The transverse magnetization. A total transverse (to the magnetic field in the $z$ direction) magnetization reads

[TABLE]

and is experimentally measurable. It disappears in the GS of the polar phase and equals $(1-q^{2}/4)/L$ in the broken-symmetry GS. Its typical evolution is depicted in Fig. 1. We see there that nothing happens in the polar phase. The system starts nontrivial evolution in the broken-symmetry phase at a distance $\hat{t}/\tau_{Q}$ after the critical point was passed. The magnetization grows fast from that point until it exceeds the static prediction and starts oscillations with the amplitude decreasing in time. We consider slow transitions. Therefore, by the end of time evolution, when $q=0$, the system is in the slightly perturbed ferromagnetic GS: globally $M_{T}L\approx 1$ (Fig. 1) and locally $L^{2}m_{T}(z)\approx 1$ (Fig. 3). We can now ask: Does the scaling (3) hold? To find out we define arbitrarily $\hat{t}$ as the instant when $M_{T}L$ intersects $1\%$. A fit to numerics for $\tau_{Q}\geq 10$ yields $\ln\hat{t}=(0.056\pm 0.01)+(0.332\pm 0.002)\ln\tau_{Q}$ which confirms prediction (3). This fit is presented in Fig. 1a, where the gradual departure of the numerical data for $\tau_{Q}<10$ from $\hat{t}\sim\tau_{Q}^{1/3}$ indicates that $\tau_{Q}\gg 1$ or $37ms$ has to be taken for the observation of $1/3$ exponent: quench has to be slow enough to reflect the critical dynamics.

In the GS configuration of the broken-symmetry phase the vector $(f_{x},f_{y})$ can have arbitrary orientation, so in the dynamical problem considered here it is interesting to find out how is this rotational symmetry broken. When unstable evolution starts, spatial correlations in magnetization appear (Fig. 2a). In the subsequent evolution these correlations evolve such that the correlation length increases: see Fig. 2b obtained by the end of time evolution. This is a generic picture though the details depend on the quench time $\tau_{Q}$ and initial state of the system. This behavior suggests creation of spin textures stephens ; vilenkin . In our case, topological textures are spin configurations where the magnetization direction varies in space so that the kinetic energy term in (1) is not minimized, but magnetization magnitude follows closely a GS result. Such structures appear in 1D when the first homotopy group of the vacuum manifold $\cal M$ is nontrivial, which happens here: $\pi_{1}({\cal M})={\mathbb{Z}}$ makela . These textures are characterized by the winding number, $\frac{1}{2\pi}\int dz\frac{d}{dz}{\rm Arg}(f_{x}+if_{y}),$ which is not conserved. Indeed, it reads $+1$ in Fig. 2a, while by the end of that evolution (Fig. 2b) it equals [math].

Are different stages of this evolution experimentally observable? Let’s look at $\tau_{Q}=10$ case presented in Figs. 1-3. The evolution from the phase boundary to the first peak in magnetization $M_{T}$ (the $q=0$ point) takes $2.81\times 37ms\cong 104ms$ ($2\tau_{Q}=740ms$). Both these time scales are well within the reach of the experiment nature_kurn .

The longitudinal magnetization. Initially, $f_{z}(z)\approx 0$ so that $\int dzf_{z}\approx 0$. The conservation of the latter allows only for creation of a network of magnetic domains (nontopological structures with fixed $f_{z}$ sign) having opposite polarizations. The domains appear by the time when the system enters unstable evolution and the maxima of $|f_{z}|$ tend to move towards the minima of $m_{T}$ (Fig. 3). More quantitatively, we performed $N_{r}$ evolutions starting from different initial conditions, but fixed $\sigma$. As in the experiment nature_kurn , we average over these runs to wash out shot-to-shot fluctuations. In Fig. 4 we plot the mean domain size: $\xi=\sum_{i}\xi_{z}(i)/N_{r}$, where $i=1,...,N_{r}$ and $\xi_{z}(i)$ is the mean domain size in the $i$-th run. As shown in Fig. 4a, for $t\lesssim\hat{t}$ we observe $\xi\approx f(t/\tau_{Q}^{1/3})$ as for $M_{T}(t)$. The domains are formed on a time scale of $\sim\hat{t}$. A simple analysis based on KZ theory kibble ; zurek suggests that their characteristic post-transition size, $\hat{\xi}$, should be roughly given by $\int_{0}^{\sim\hat{t}}dtv_{s}(t)$, where $v_{s}(t)$ is a sound velocity. There are two sound modes in the broken-symmetry phase that propagate both spin and density fluctuations ueda_broken : the faster (slower) one has velocity $\sim\sqrt{c_{0}}$ ($\sim\sqrt{q|c_{1}|}$). Putting any of these as $v_{s}$ into the integral, and assuming $\tau_{Q}\gg 1$ for the slower mode, we get

[TABLE]

This result correctly predicts the scaling property of the size of post-transition “defects” as is evident from the overlap of different curves in Fig. 4, which shows up for $\tau_{Q}\geq 25$ or $0.9s$. Quantitatively, we define $\hat{\xi}$ as the value of $\xi$ averaged over $q\in[1/2,1]$ to wash out post-transition fluctuations. A fit got us $\ln\hat{\xi}=(-0.38\pm 0.03)+(0.30\pm 0.01)\ln\tau_{Q}$, in good agreement with (5). The fit was done to $\tau_{Q}\geq 30$ data and is presented in Fig. 4b which illustrates that smaller $\tau_{Q}$ data gradually departs from $1/3$ scaling law.

Now we focus on the analytical calculations providing predictions about early stages of time-evolution. We assume that the wave-function stays close to the polar phase GS, $\Psi^{T}=(\delta\psi_{1}(t),1/\sqrt{L}+\delta\psi_{0}(t),\delta\psi_{-1}(t))\exp(-i\mu t)$, where the chemical potential is $\mu=c_{0}/L$, $|\delta\psi_{m}|\ll 1/\sqrt{L}$, and $\int dz(\delta\Psi_{0}+\delta\Psi_{0}^{*})\equiv 0$ to keep $\int dz\Psi^{\dagger}\Psi=1+O(\delta\Psi^{2})$. Linearizing the coupled nonlinear-Schrödinger equations that describe the system we get $f_{\chi}={\rm Re}G_{\chi}$, where $\chi=x,y$, $G_{x}=\sqrt{2}(\delta\Psi_{1}+\delta\Psi_{-1})/\sqrt{L}$, $G_{y}=i\sqrt{2}(\delta\Psi_{1}-\delta\Psi_{-1})/\sqrt{L}$, and

[TABLE]

where $\alpha=2|c_{1}|/L$. To solve this equation we go to momentum space, $a_{\chi}(k)=\int dzf_{\chi}\exp(ikz)$ and $b_{\chi}(k)=\int dz{\rm Im}G_{\chi}\exp(ikz)$, getting

[TABLE]

Diagonalizing the matrix from Eq. (6) we see that there is instability for $k^{2}/\alpha<2-q$ as in the Bogolubov spectrum of this model ueda_broken . Thus, the system is stable in the polar phase, and so small initial perturbations do not grow during the evolution towards broken-symmetry phase. The instability for $q<2$ is responsible for the magnetization jump in Fig. 1 and the subsequent breakdown of the linear approach (Fig. 1b).

To solve Eq. (6) with $q(t)$ given by (2) we derive the equation for $d^{2}a_{\chi}(t)/dt^{2}$, keep leading order terms in the slow transition ($\tau_{Q}\gg 1$) and long-wavelength ($k^{2}/\alpha\ll 2$) limits, and get that

[TABLE]

where $\kappa=(\alpha^{2}/2)^{1/3}$, $\alpha_{k\chi}$ and $\beta_{k\chi}$ are constants given by initial conditions, while Ai and Bi are Airy functions. From (7) we see that the instability arises from unbounded increase of the Bi$(s)$ function happening for $s>0$, i.e., $k^{2}/\alpha<2-q(t)$, which is a dynamical manifestation of the static result for unstable modes. This solution works till $t\sim\hat{t}\sim\tau_{Q}^{1/3}$ when a significant increase of $f_{\chi}$ invalidates the linearized theory: this calculation rigorously derives scaling (3). Additionally, the solution (7) can be reliably used as long as $\tau_{Q}\gg 1$ or $37ms$, which is also supported by numerics (Fig. 1a). The quench time scale in the experiment nature_kurn is much smaller than this bound. Finally, these results hold for any initial state spread over the $k$ modes.

The (re)scalings $t/\tau_{Q}^{1/3}$ and $\hat{\xi}\sim\hat{t}\sim\tau_{Q}^{1/3}$ derived above in a 1D system were also found by different means in a 2D spinor condensate austen . A trivial extension of our mean-field analytical calculations to 2D and 3D systems shows that they hold for any number of spatial dimensions.

To summarize, we have developed a theory of the dynamics of symmetry-breaking in the quantum phase transition inspired by the experiment nature_kurn , but for the range of quench rates that are sufficiently slow so that the critical scalings can determine phase transition dynamics. This regime should be accessible by a “slower” version of the quench nature_kurn . Our analysis points to a Kibble-Zurek-like scenario, where the state of the system departs from the old symmetric vacuum with a delay $\sim\hat{t}$ after the critical point was crossed. This sets up an initial post-transition state with a characteristic length scale $\hat{\xi}$ (5). This scale should determine the initial density of topological features. In our 1D simulations textures appear, but we predict that in real 3D experiments other topological defects are created (as they were in nature_kurn ), and the distance between them should be initially $\sim\hat{\xi}$. Such topological defects are more stable than textures so measurement of their density should be possible and would be a good test of the theory we have presented.