磁阱中超冷玻色氣體的 BEC 相變與 Kibble - Zurek 動力學(下)

  • 物理專文
  • 撰文者:劉翼綱、郭西川
  • 發文日期:2022-10-23
  • 點閱次數:180

壹. 古典場模型與數值模擬

GP 方程式是平均場近似下的理論,它僅能描述凝聚體波函數在絕對零度下的同調演化。我們採用古典場方法來處理 BEC 在非零溫環境下的演化,包含了玻色氣體中同調與非同調的交互作用過程。古典場方法的原理簡述如下[3]

根據給定的能量閥值Ecut,如圖三(a)所示,系統可分成兩部分;(一)古典場區(C-field region): 由所有能量低於Ecut 的模態組成,包含了 BEC 及其他高佔有數的低能量模態。(二)非同調區(incoherent region):由能量高於 Ecut 之低佔有數激發態所組成(可視為一處於熱平衡溫度為 T 的熱庫)。系統的演化主要發生在古典場區,與 GP 方程有類似的形式。這也是此方法最大的優點,不僅在數值計算上容易處理並且可以適用在大部分的實驗參數範圍。

picture3.png

如圖三(b)、(c)所示,非同調區與古典場區粒子間的碰撞模式可分為:(i)非同調區中兩粒子碰撞後,其中之一進入古典場區;(ii)一古典場與一非同調區粒子碰撞,古典場區粒子數不變但佔有數重新分配。由於效應(ii)在古典場佔有數極大時可忽略,故我們僅考慮(i),此時古典場可以由一隨機投影GP 方程(stochastic projected GP equation,SPGPE)描述[3,4]

\[d\psi=\hat{\mathcal{P}}\left\{\frac{i+\gamma}{\hbar}\left[\mu-\hat{\mathcal{L}}\right]\psi\right\}dt+dW,
\hspace{1cm}(1)\]

其中投影算符 \( \hat{\mathcal{P}}\) 排除非古典場模態進行演化,而

\[\hat{\mathcal{L}}=-\frac{\hbar\nabla^2}{2M}+V_{\mathrm{trap}}\left(\mathbf{r}\right)+g\left|\psi\right|^2。 \hspace{1cm}(2)\]

\(dW\left(\mathbf{r},t\right)\)為經散射進入古典區的粒子所帶進之複數高斯白噪音,它滿足以下相干關係

\[\left\langle d W^\ast\left(\mathbf{r},t\right)dW\left(\mathbf{r}\prime,t\right)\right\rangle=\frac{2\gamma k_BT}{\hbar}\delta\left(\mathbf{r}-\mathbf{r}^\prime\right)dt,  \hspace{1cm}(3)\]

其中 \(\gamma \) 為古典區粒子的增加率,\(\mu\) 為化學勢,\(k_B\) 為波茲曼常數。依據“ 各態歷經假說”(ergodic assumption),可利用短時間平均來計算密度矩陣\(\rho\left(\mathbf{r},\mathbf{r}^\prime,t\right)=\left\langle\psi^\ast\left(\mathbf{r}^\prime,t\right)\psi\left(\mathbf{r},t\right)\right\rangle \)。再依據 Penrose - Onsager 定理, 亦即 BEC 模態(或是 PO 模態)對應到密度矩陣之最大本徵值的解,我們便可以計算 BEC 在隨機噪訊環境下的演化。

根據 Trento 團隊的實驗參數[1,2],我們考慮於 \( -\tau_Q < t < \tau_Q\) 期間化學勢與溫度呈線性變化:

picture4.png

取\(\mu_f=\mathrm{22}\hbar \omega_{⊥},T_0=500 \hspace{1mm}\rm nK\) 與\(\Delta T = 290 \hspace{1mm}\rm nK\) ;位能勢\(V_{\mathrm{trap}}\left(\mathbf{r}\right)=\left({\mathrm{M}}/{\mathrm{2}}\right)\left[\omega_x^2x^2+\omega_\bot^2\left(y^2+z^2\right)\right] \)的振盪頻率,\(\omega_x=\mathrm{2}\pi×13 \) Hz,\(\omega_\bot=\mathrm{2}\pi×131.4 \) Hz。被捕捉在此位能阱的原子團呈現狹長(\(\omega_\bot\approx10\omega_x \))的雪茄型,屬三維系統,故可觀測到的缺陷為渦漩線。我們透過 SPGPE 模擬來探索該系統相變點之動力學行為[4,5]。系統的初態設定為一熱平衡的玻色氣體,末態為一內含 BEC 的超冷玻色氣(凝聚體的粒子數約佔原有氣體粒子數的75%)。圖四中,我們考慮\(\tau_Q= 9 \sim 720 \hspace{1mm}\rm ms\),\(\gamma=0.005\) 。根據不同的 \( \tau Q\),我們歸納出以下三種降溫模式:(a)緩慢降溫—系統在通過相變點後,逐一生成分散的 BEC,但無大量的缺陷生成。因為熱擾動影響,缺陷將會離開系統,故系統有很大的機會觀測不到缺陷;(b)典型降溫—系統生成數個拓樸缺陷,隨後因為熱漲落與缺陷間的交互作用,僅有少數能夠存活長時間。(c)快速降溫—系統出現大量的缺陷,即便缺陷的平均生命期有限,但依然能夠有數個能夠存活的非常久。


貳. 線性近似與臨界暴漲

系統粒子數由化學勢 \(\mu\) 決定。對於均勻系統,當 \(\mu > 0\) 時,系統將發生凝聚態現象,故可以藉由μ 來做為相變調控參數。當系統尚未發生凝聚態現象時,即 \(\mu > 0\)  ( \(t < 0\) ) ,古典場波函數僅由無序的相位構成且 \( \psi \sim 0\)(及其期望值),故我們能夠忽略式(2)中的非線性項。對SPGPE 式(1)做傅立葉變換得到\[\frac{d\widetilde{\psi}}{dt}=\left[\frac{\left(\gamma+i\ \right)}{\hbar}\frac{\hbar^2k^2}{2M}+\frac{\gamma}{\hbar}\mu\left(t\right)\right]\widetilde{\psi}+\widetilde{\zeta},\hspace{1cm} (5)\] 其中 \(\widetilde{\psi} \) 為古典場波函數的傅立葉轉換,而噪訊項 \(\widetilde{\xi}\) 則遵守相干條件,

\[\left\langle{\widetilde{\zeta}}^\ast\left(\mathbf{k},t\right)\widetilde{\zeta}\left(\mathbf{k}^\prime,t\right)\right\rangle\ =\frac{2\gamma k_BT}{\hbar}\delta\left(\mathbf{k}-\mathbf{k}^\prime\right)\delta\left(t-t\prime\right)。
\hspace{1cm}(6)\]

隨機方程式(5)於動量空間中有精確解且僅躁訊項有主要貢獻,其精確解形式為
\[\widetilde{\psi}\left(\mathbf{k},t\right)=\int_{-\infty}^{t}{dt^\prime\widetilde{\zeta}\left(\mathbf{k},t^\prime\right)e^{-\frac{\gamma\hbar k^2}{2M}\left(t-t^\prime\right)+\frac{\gamma\mu_f}{2\hbar\tau_Q}\left(t^2-t^{\prime2}\right)}} \hspace{1cm}(7)\] 

古典場波函數之動量分布(譜函數)為可量測之物理量,其定義為
 \[f\left(\mathbf{k},t\right)=\left\langle{\widetilde{\psi}}^\ast\left(\mathbf{k},t\right)\widetilde{\psi}\left(\mathbf{k},t\right)\right\rangle。 \hspace{1cm}(8)\] 


picture5.png
我們可以假定當 \(t > 0\),系統發生BEC 相變,故上述的解僅適用於接近相變點。利用高斯積分將其近似為
\[f\left(\mathbf{k},t\right)\approx\frac{2\sqrt\pi\gamma k_BT}{\hbar}\hat{t}e^{\left({t}/{\hat{t}}-\hat{\xi}k^2\right)^2}, \hspace{1cm}(9)\] 
其中
 \[\hat{t}=\left(\frac{\hbar\tau_Q}{\gamma\mu_f}\right)^{-\frac{1}{2}}, \hat{ξ}=τ_Q^{1/4}  \left(\frac{4M^2μ_f}{γℏ^3} \right)^{-\frac{1}{4}}。 \hspace{1cm}(10)\] 

由式(9),我們發現譜函數滿足標度不變性(scaling invariant)
 \[f\left(\mathbf{k},t\right)=\hat{t}F\left(\frac{t}{\hat{t}},\hat{\xi}\mathbf{k}\right), \hspace{1cm}(11)\] 
另外,\({\hat{t}}/{{\hat{\xi}}^3}\propto{\hat{\xi}}^{-1} \)符合式(5)的標度行為假設。


上述標度行為適用於\({\ \tau}_Q\gg\hat{t}\  \)之模擬,因其不受式(4)之降溫參數停止變化之影響。此外要明確地知道相變點的位置,在考慮位能阱及粒子交互作用的情況下,我們知道相變點不會發生在 \(\mu = 0\)。我們透過研究在不同時刻下給定 \(\mu\) 與 \(T\) 之古典場平衡態性質發現相變發生在 \(t_c ≈ (0.188 ± 0.013)τ_Q \) [4],於此我們將對於數值計算中之 \(t\) 平移至 \(t_c\)。在引入相變點位移後,我們發現模擬數值計算滿足式(11)。圖五中以 \(k = 0\) 為例,將譜函數經根據式(11)進行標度轉換後,譜函數至\(\left(t-t_c\right)/\hat{t}\approx\mathrm{1.25}\  \)左右皆遵守標度律。

由式(9),譜函數在通過相變點後會隨時間增加呈指數發散,其特徵時間為 \( \hat{t}\),但不同波數的發散時刻會隨著 \(k\) 增加而延後。由譜函數我們回推對應的系統密度為
\[\left\langle\left|\psi\left(\mathbf{r},t\right)\right|^2\right\rangle\approx\frac{\gamma k_BT}{16\pi^{3/2}\hbar}\frac{\hat{t}}{{\hat{\xi}}^3}\left(\frac{\hat{t}}{t}\right)^\frac{3}{2}{e^{\left(t/\hat{t}\right)}}^2 \hspace{2cm},(12)\]

本近似預期可精確至 \(t\sim\hat{t}\),因隨著時間演化密度成指數發散,非線性項將不可忽略,系統也將不再臨界發散,故可以被激發的最高動量波數\(~\ {\hat{\xi}}^{-1}\)將決定系統的最小尺度,即相干長度之尺度。我們發現與KZ 圖像不同是,缺陷之相干長度並非由相變點前之凍結時域決定,而是由過相變點後,系統所能激發的最大波數決定。故此處的 \( \hat{t}\) 與KZ圖像之凍結時域大小相符但僅發生在相變後,且恰巧對應到 \(zν = 1\) 其中 \(z = 2\) 之平均場預測。

參. 降溫相圖

以線性近似分析SPGPE,我們可以更進一步利用局域密度近似來分析Trento 團隊的實驗。為利於討論與計算,我們依然假定相變點發生在 \(\mu = 0 \)。在陷阱中,空間中各區間發生相變的條件為
\[\mu\left(t\right)-V_{\mathrm{trap}}\left(\mathbf{r}\right)>0,\hspace{2cm} (13)\]
可推知系統中相變波前尺寸為
\[a_{x/\bot}\left(t\right)=\sqrt{\frac{2\mu_f}{M\omega_{x/\bot}^2}\frac{t}{\tau_Q}},\hspace{2cm} (14)\]
與其對應的波前速率
\[v_{x/\bot}=\frac{da_{x/\bot}}{dt}=\sqrt{\frac{\mu_f}{2M\omega_{x/\bot}^2\tau_Qt}},\hspace{2cm} (15)\]
比較系統發生相變後,聲速地平線波前(相干長度增長)增長的速率
\[\hat{v}\approx\frac{\hat{\xi}}{\hat{t}}=\left(\frac{\gamma^3\hbar\mu_f}{4M^2}\right)^\frac{1}{4}\tau_Q^{-\frac{1}{4}},\hspace{2cm} (16)\]
當聲速地平線擴展速度小於相變波前時,我們可得勻相軸長為
\[{\hat{a}}_{x/\bot}=\left(\frac{4\mu_f^3}{M^2\omega_{x / \bot}^8\gamma\hbar\tau_Q^3}\right)^{\frac{1}{4}\ }\hspace{2cm} (17)\]
反之當聲速地平線擴展速度大於相變波前時,即當 \( a_x (t)/ \hat{a}_x\)或 \(a_{\perp}(t)/\hat{a}_{\perp}\) 大於1 時,勻相近似將會失準。我們整理出以下幾個降溫模式

\( \tau_Q\lesssim\hat{t} \),          準瞬間降溫
 \( \hat{t}<\tau_Q<\frac{\mu_f}{\gamma\hbar\omega_\bot^2} \),    勻向降溫
\(\frac{\mu_f}{\gamma\hbar\omega_\bot^2}<\tau_Q<\frac{\mu_f}{\gamma\hbar\omega_x^2}\),長軸勻向降溫
\(\frac{\mu_f}{\gamma\hbar\omega_x^2}<\tau_Q\)。  非勻向降溫

圖六為上述不同降溫模式所對應的相圖。我們的降溫範圍剛好處於準瞬間降溫與勻向降溫區間,其中前節所述之” 快速降溫”對應到準瞬間降溫,而典型與緩慢降溫屬於勻向降溫,可以解釋即便我們在線性近似之預測所中忽略的外加陷阱,卻可以依然適用於我們的數值計算中。想要探討長軸勻向與非勻向降溫動力學,根據本參數需要極緩慢的降溫,可預期會有更加複雜有趣的結果。
picture6.png


肆. 後相變時期

最後,我們比較Trento 實驗進行KZ 缺陷數目量測的結果[1]。該實驗使用飛行時間(Time of flight)量測法, 觀測在特定時刻凝聚體由位能阱釋放後膨脹的影像。他們僅能在凝聚體生成後約250 ms 左右進行有效的缺陷數測量。由上述理論與模擬,在臨界暴漲期間,即便已過了相變點,PO 模態尚未能夠顯示完整的凝聚體[4],故可以預期實驗上可觀測到凝聚體的時間是已過了臨界暴漲期。為了便於與實驗觀測比較,我們模擬實驗的量測方式,以PO 模態粒子數達最終平衡態粒子數的5% 時為實驗可能量測到凝聚體生成的時間點,標記為 \(t_bec\)。在圖七我們比較 \(t − t_bec\) 在50 ms 與200 ms 所能觀測到缺陷數目與實驗量測之結果。因實驗僅能控制溫度隨時間之變化,故我們將可觀測之渦旋數 \(N_{\mu}\) 對 \(dT/dt \) 作圖比較,發現數值計算結果整體行為與實驗吻合,實驗中所量測到的標度律為 \(N_{\mu}\propto \tau Q^{−7/3}\),符合受侷限系統之預測[6]。即便因為熱漲落,數值模擬在後期所能觀測到的粒子數較實驗結果為少,但整體趨勢一致。


快速降溫出現極值之平台,可以理解為因降溫過程生成大量之缺陷,在SPGPE 模擬中我們發現因受交互作用影響、熱漲落而將缺陷排出於系統外,導致其數目快速地減少,故僅存數個可觀測的缺陷。這些缺陷因距離遙遠,猶如相互隔離的缺陷,即便在熱漲落下,依舊可以存活相當長的時間。而較緩的降溫則因本身能夠生成的缺陷數變不多,且如上述相互隔離的缺陷可以相當地長壽,故能夠吻合KZ 圖像預測。總而言之,想要簡單地驗證KZ 機制,需要具有穩定不受熱漲落影響的缺陷系統[7]。

然而即便這些隔離的長壽缺陷,最終還是會受到熱漲落而離開系統。根據模擬與實驗顯示,單一的渦漩最長可存活超過 1 秒,而接著系統將演化至一有限溫度、處於平衡態的凝聚體,如圖三(a)所示。


伍. 總結


本文利用古典場近似探索一受侷限之超冷玻色氣體經過 BEC 相變點之動力學行為。我們發現缺陷數並非如 KZ 圖像所描述,由被凍結的相位決定,而是在通過相變點後,系統將發生臨界暴漲,此時所能激發之最大動量模態,或是最小系統尺度決定。其伴隨著密度暴增,將離開暴漲時域,而此時域恰與 KZ 圖像之凍結時域相符。透過在 SPGPE 線性分析中引入局域密度近似,得出勻向-非勻向降溫相圖,而本文探索的參數剛好處於勻向區。

透過 SPGPE 數值模擬(這是首宗引入化學勢與溫度參量的數值模擬)與分析古典場中的 PO 凝聚體模態,我們得以更進一步探究凝聚態在相變後之動力學行為,得知了生成的拓樸缺陷會因為彼此的交互作用而排出系統,導致有大量的自發缺陷無法順利在實驗中被觀測到,數值模擬結果與實驗有相當好地吻合,並了解在系統如何逐漸演化至平衡態。

本研究由英國 Newcastle 大學、義大利 Trento 大學、波蘭 Jagiellonian 大學與台灣國立彰化師範大學四方合作而成,論文詳見[4,5]。

參考文獻:

[1] S. Donadello, S. Serafini, T. Bienaimé, F. Dalfovo, G. Lamporesi and G. Ferrari, “Creation and counting of defects in a temperature-quenched Bose-Einstein conden-sate”, Phys. Rev. A 94, 023528(2016).

[2] S. Donadello, S. Serafini, M. Tylutki, L. P. Pitaevskii, F. Dalfovo, G. Lamporesi and G. Ferrari, “Observation of Solitonic Vortices in Bose-Einstein Condensates” Phys. Rev. Lett. 113, 065302(2014).

[3] 蘇士煒、劉翼綱、郭西川, “ 冷原子不冷-非零溫玻色-愛因斯坦凝聚體的動力學,” 物理雙月刊,第34 卷

3 期,184(2012)。

[4] I.-K. Liu, S. Donadello, G. Lamporesi, G. Ferrari, S.-C. Gou, F. Dalfovo and N. P. Proukakis, “Dynamical equilibration across a quenched phase transition in a trapped quantum gas”, Comm. Phys. 1, 24(2018).

[5] I.-K. Liu, J. Dziamaga, S.-C. Gou, F. Dalfovo and N. P. Proukakis, “Kibble-Zurek dynamics in a trapped ultracold Bose gas” Phys. Rev. Research 2, 033183(2020).

[6] A. del Campo, A. Retzker, and M. B. Plenio, “The inhomogeneous Kibble-Zurek mechanism: Vortex nucleation during Bose-Einstein condensation”, New Journal of Physics 13, 083022(2010).

[7] S.-W. Su, S.-C. Gou, A. Bradley, O. Fialko, J. Brand, “Kibble-Zurek scaling and its breakdown for spontaneous generation of Josephson vortices in Bose-Einstein condensates”, Phys. Rev. Lett. 110, 215302(2013).