一种加速预测梯度流方法，用于计算具有Lee-Hang-Yang修正的旋转二维玻色-爱因斯坦凝聚体的基态

《Communications in Nonlinear Science and Numerical Simulation》：An accelerated projected gradient flow method for computing the ground states of rotating two-dimensional Bose-Einstein condensates with Lee-Hang-Yang correction

　　唐凯智|冯志旭|舒清洲|胡冰
四川大学数学学院，成都，610064，中国

**摘要**
在本文中，我们引入了一种对Gross-Pitaevskii方程（GPE）的正则化方法，并加入了Lee-Huang-Yang（LHY）修正项，以消除由其对数非线性引起的奇异性。这种正则化的GPE包含一个小的正则化参数0 ≤ ε ? 1，且其正则化能量随着ε → 0而线性收敛到原始GPE的能量。此外，我们还严格证明了正则化模型基态的存在性。基于这种正则化公式，我们提出了一种加速投影梯度流（APGF）方法，用于计算带有LHY修正的旋转二维玻色-爱因斯坦凝聚体（BECs）的基态。APGF中的拉格朗日乘子具有一个显式表达式，确保了该方法在连续水平上保持质量守恒性。然后，我们提出了两种时间半离散化方案来离散化APGF。通过结合傅里叶伪谱方法进行空间离散化，这两种方案都可以高效求解。最后，我们提供了大量的数值示例，以比较两种方案的性能，并研究了在不同设置下带有LHY修正的旋转二维BECs的基态。

**引言**
玻色-爱因斯坦凝聚体（BECs）是一种在超冷原子气体中实验实现的宏观量子现象[1]，[2]，已成为探索基本量子多体物理和模拟复杂量子系统的关键平台。它们对参数（如原子间相互作用、外部势和维度）的高度可控性，使我们能够深入了解超流性、量子涡旋和量子相变。在各种现象中，受到外部旋转驱动的旋转BECs特别引人注目，因为它们可以模拟旋转的超流体和中子星，产生复杂的涡旋 lattice 和拓扑激发[3]，[4]，[5]。近年来，对量子液滴的研究（由量子涨落稳定的自束缚态）重新激发了这一领域[6]，[7]，[8]，[9]，[10]。这些液滴通常在包含Lee-Huang-Yang（LHY）修正的扩展平均场框架内进行描述，该修正考虑了超出平均场的量子涨落。LHY修正从根本上改变了BECs的稳定性和形态，使得在二元或偶极凝聚体中可能存在稳定的液滴相，并修改了旋转系统中的涡旋态结构[6]，[8]，[10]，[11]。因此，理解带有LHY修正的旋转BECs的基态对于探索新型量子流体相及其在量子信息和精密计量学中的应用至关重要。

在二维情况下，Lee-Huang-Yang（LHY）修正源自平均场之外的量子涨落，表现为对数非线性。这种有效的描述对于二维量子液滴尤为重要，其中平均场相互作用和LHY修正相互平衡形成自束缚态[12]，[13]。需要指出的是，LHY修正在三维中的形式有所不同。我们的研究聚焦于二维情况，二维LHY泛函的原始推导可以在[13]中找到：
$$
i\partial_t\psi = \left(-\frac{1}{2}\Delta + V(x) + \beta|>\psi|^2 + \delta|>\psi|^2\ln|>\psi|^2 - \omega L_z\right)\psi, \quad x \in \mathbb{R}^2, \quad t \geq 0,
$$
其中 $\psi \equiv \psi(x, t)$。角动量定义为 $L_z = -i(x\partial_y - y\partial_x)$。该方程包含一个表示LHY修正的非线性项。系统的特征参数包括：旋转速度 $\omega$、两体相互作用和LHY相互作用的强度（即 $\beta$ 和 $\delta$），以及实值外部陷阱势 $V(x)$。参数 $\beta$ 表示两体s波相互作用的强度，与散射长度成正比；$\beta > 0$ 对应排斥相互作用，而 $\beta < 0$ 对应吸引相互作用。参数 $\delta$ 表示LHY量子涨落修正的强度，在二维中通常是正的。$\beta$ 和 $\delta$ 之间的相互作用决定了凝聚体的稳定性和形态。动力学守恒总质量（或粒子数）：
$$
N(\psi(\cdot, t)) = \|\psi(\cdot, t)\|^2 := \int_{\mathbb{R}^2|>\psi(\cdot, t)|^2 \, dx \equiv 1, \quad t \geq 0,
$$
以及总能量：
$$
E(\psi(\cdot, t)) := \int_{\mathbb{R}^2\left(\frac{1}{2}|\nabla\psi|^2 + V(x)|>\psi|^2 + \beta^2|>\psi|^4 - \omega\psi'\nabla_z\psi + \delta^2|>\psi|^4\ln|>\psi|^2 - 1/2\right) \, dx \equiv E(\psi(\cdot, 0)), \quad \psi' \表示 \psi$ 的复共轭。
基态被定义为在质量约束下的能量泛函的最小化值：
$$
\phi_g(x) := \argmin_{\phi \in \mathcal{E}(\phi), \quad \mathcal{E} := \{ \phi \in \H_1(\mathbb{R}^2) \mid \|\phi|\|^2 = 1, \quad E(\phi) < \infty \}。
$$
对于旋转系统，它对应于旋转框架中能量最低的稳定平衡配置，通常具有以特定模式排列的量子化涡旋。与（1.4）相关的欧拉-拉格朗日方程为：
$$
\mu\psi = -\frac{1}{2}\Delta\psi + V(x)\psi + \beta|>\psi|^2\psi + \delta|>\psi|^2\ln(|>\psi|^2)\psi - \omega L_z\psi, \quad \psi \in \mathcal{S}, \quad x \in \mathbb{R}^2,
$$
这是一个对于（$\mu$, $\psi$）的$L^2$-归一化非线性特征值问题。相应的特征值 $\mu$ 可以通过下式计算：
$$
\mu(\psi) = \int_{\mathbb{R}^2\left(\frac{1}{2}|\nabla\psi|^2 + V(x)|>\psi|^2 + \beta|>\psi|^4 + \delta|>\psi|^4\ln|>\psi|^2 - \omega\psi'\nabla_z\psi\right) \, dx = E(\psi) + \int_{\mathbb{R}^2\left(\beta^2|>\psi|^4 + \delta^2|>\psi|^4\ln|>\psi|^2 + 1/2\right) \, dx。
$$
在本文中，我们始终假设 $(1 + K)\omega^2|x|^2/2 \leq V(x) \in \mathcal{L}_{\text{loc}}(\mathbb{R}^2)$，其中 $K$ 是某个常数且 $K > 0$，这也是存在基态的常用充分条件[14]。当 $\delta > 0$ 时，问题（1.4）的基态存在。基态存在性的证明细节可以在附录A中找到。值得注意的是，即使当 $\beta < 0$（吸引相互作用）时，对数项在低密度下也表现出超二次增长，这防止了体系坍缩并确保了在适当参数范围内存在稳定的基态。我们的分析和数值方法涵盖了 $\beta > 0$ 和 $\beta < 0$ 的情况。

由于能量泛函的非凸性、高度非线性（包括LHY项）、由对数非线性引起的奇异性以及质量守恒和旋转对称性的约束，数值计算这样的基态面临着重大挑战。此外，旋转项使得结构变得复杂，需要仔细处理以准确解析涡旋 lattice，而奇异性要求适当的正则化。已经进行了大量的数学和数值工作来研究旋转BEC的基态；有关不完整的列表，请参考[14]，[15]，[16]，[17]，[18]，[19]，[20]，[21]，[22]，[23]，[24]及其中的参考文献。BEC的基态可以通过优化问题（1.4）或$L^2$-归一化非线性特征值问题（1.5）来数值计算。基于（1.5）的非线性特征值求解器可以在[25]，[26]，[27]，[28]及其中的参考文献中找到。在最小化方法中，有两种主要策略：一种是基于相关偏微分方程（PDEs）离散化的方法，如$L^2$-投影梯度流[16]，[18]，[19]，[29]，[30]，[31]，[32]，[33]，[34]和Sobolev投影梯度流[20]，[28]，[35]，[36]；另一种是直接通过特定优化方法最小化能量泛函（1.4）[17]，[21]，[22]，[24]，[28]，[37]，[38]。然而，关于带有LHY修正的旋转BEC的基态的研究相对较少。

本文的目的是开发和分析专门用于计算带有LHY修正的旋转二维BECs基态的健壮、高效且保持质量的数值方法。为此，我们首先开发了一个（1.1）的正则化模型，以避免对数非线性项的奇异性。然后基于正则化模型，我们提出了一种结合拉格朗日乘子的新型加速梯度流公式，以实现精确的质量守恒。我们将通过全面的数值实验来证明我们方法的效率，包括比较所提出的方法、正则化模型的收敛性测试以及探索各种物理状态。这项工作旨在为研究带有LHY修正的旋转二维BECs的丰富基态相图提供一个可靠的计算工具。

**本文的其余部分组织如下。**
第2节中，我们首先对GPE和能量进行正则化，以消除对数非线性项的奇异性，然后分析正则化能量如何逼近原始能量。此外，我们证明了正则化模型基态的存在性。第3节中，我们提出了一种加速投影梯度流（APGF）来计算基态，该方法在连续水平上保持质量并降低修正后的能量。然后，构建了两种半离散化方案来离散化APGF的时间变量。最后，我们提出了一种傅里叶伪谱方法来进一步离散化APGF的空间变量[39]。第4节中，进行了大量的数值示例，以展示这些方案的性能。第5节得出了几个结论。

**正则化模型及其性质**
众所周知，对数非线性的奇异行为（即 $\ln \rho \rightarrow \infty$ 随着 $\rho \rightarrow 0$）可能导致数值爆炸或显著的四舍五入误差。为了解决对数奇异性带来的数值挑战，我们直接用一个小的正则化参数 $0 < \epsilon \leq 1$ 对（1.1）（以及（1.3）中的对数非线性 $f(\rho)$）进行正则化[40]：
$$
i\partial_t\psi_{\epsilon} = \left(-\frac{1}{2}\Delta + V(x) + \beta|>\psi_{\epsilon}|^2 + \delta|>\psi_{\epsilon}|^2\ln(|>\psi_{\epsilon}|^2 + \epsilon) - \omega L_z\right)\psi_{\epsilon}, \quad x \in \mathbb{R}^2, \quad t \geq 0,
$$
其中 $\psi_{\epsilon}$ 表示正则化后的波函数。

**数值方法及其性质**
为了解决快速旋转和LHY修正带来的数值挑战，本节我们提出了一种加速投影梯度流（APGF）方法，用于计算带有LHY修正的旋转二维BEC的基态。我们开发了两种时间半离散化方案来离散化APGF，同时使用傅里叶伪谱方法进行空间离散化。

**数值结果**
在本节中，我们首先比较了两种提出的方法。然后我们研究了在不同模型参数下的基态模式。本文采用的停止准则为：
$$
\| \psi_{n+1} - \psi_n \|_{\infty} \leq \epsilon_r.
$$
在实践中，我们选择停止准则（4.21），其中 $\epsilon_r = 10^{-8}$。这些算法在Ubuntu GNU/Linux系统中使用MATLAB（2022b）实现，并在Intel(R) Xeon(R) Gold 6136R CPU（3.00GHz，25 MB缓存）上运行。

**结论**
在本文中，我们开发了用于计算带有Lee-Huang-Yang修正的旋转二维玻色-爱因斯坦凝聚体基态的健壮且高效的数值方法。为了克服对数非线性的奇异性，我们引入了对Gross-Pitaevskii方程的正则化，并证明了正则化能量随着 $\epsilon \rightarrow 0$ 而线性收敛到原始能量。正则化模型的基态存在性得到了证明。

**未引用的参考文献**
[43]

**作者贡献声明**
唐凯智：撰写 - 审查与编辑、撰写 - 原始草稿、可视化、软件、方法论、调查、形式分析、数据管理、概念化。
冯志旭：撰写 - 审查与编辑、验证、方法论、形式分析。
舒清洲：撰写 - 审查与编辑、撰写 - 原始草稿、验证、软件、方法论、调查、形式分析、数据管理、概念化。
胡冰：撰写 - 审查与编辑、监督、资源、项目。

**利益冲突声明**
作者声明他们没有已知的竞争性财务利益或个人关系，这些利益或关系可能会影响本文报告的工作。