基于脉冲泵浦的固体激光制冷研究

Research on Solid-State Laser Refrigeration Based on Pulsed Pumping

摘要: 固体激光制冷因其具有静音、无污染和微型化等优点，在光电探测器制冷、传感器制冷和空间飞行器冷却等领域中有着广阔的应用前景。本文将讨论以脉冲光作为泵浦源来进行固体激光制冷。我们以掺Yb³⁺氟化物晶体(Yb³⁺:YLiF₄)为例，基于密度矩阵，进一步建立更精确的四能级模型，系统分析脉冲泵浦下的动力学，给出制冷系统的优化激发方案。通过简化二能级模型解析求解出激发态粒子数极值对应的最优泵浦持续时间，并与四能级模型的数值解作误差比较。结合共振吸收截面和背景吸收截面与温度的变化关系，并带入制冷效率公式，最后模拟了制冷功率密度和制冷效率随脉冲持续时间和温度的变化关系。

Abstract: Solid-state laser cooling, with its advantages of being silent, pollution-free, and miniaturized, holds broad application prospects in fields such as photodetector cooling, sensor cooling, and spacecraft cooling. This paper discusses solid-state laser cooling using pulsed light as the pump source. Taking Yb³⁺-doped fluoride crystal (Yb³⁺:YLiF₄) as an example, we further establish a more precise four-level model based on the density matrix, systematically analyze the dynamics under pulsed pumping, and propose an optimized excitation scheme for the cooling system. By simplifying the two-level model, we analytically determine the optimal pump duration corresponding to the extremum of excited-state particle numbers and compare it with the numerical solution of the four-level model. Combining the variation relationships of resonance absorption cross-section and background absorption cross-section with temperature, we incorporate them into the cooling efficiency formula. Finally, we simulate the variation relationships of cooling power density and cooling efficiency with pulse duration and temperature.

1. 引言

激光作为20世纪的一项重大发明，其在军事、通信、医学、工业等等领域都起重要作用，例如激光焊接、激光切割、光纤通信、激光测距、激光雷达等。此外，激光照射在高品质光学样品上，可以诱导产生反斯托克斯荧光，从而产生制冷效应。

固体激光制冷，又被称为光学制冷，基本概念最早在1929年由德国物理学家Pringsheim提出，他当时预言光与物质相互作用可能产生冷却效应，并给出了反斯托克斯荧光制冷的基本原理[1]。所谓反斯托克斯荧光，指能量高于入射光能量的荧光。直到1995年，美国洛斯阿拉莫斯国家实验室的R.I. Epstein团队使用波长为1010纳米的激光照射掺杂稀土镱离子的氟化物玻璃ZBLANP，实现了0.3摄氏度的降温[2]。2011年，新墨西哥大学的Mansoor Sheik-Bahae教授带领的研究组用波长为1020纳米的激光成功将掺镱氟化钇锂晶体的温度从室温降到零下160摄氏度(约113 K) [3] [4]。2020年来自墨西哥大学的Azzurra Volpi等人将Yb³⁺和Tm³⁺共掺杂的LiYF₄冷却到87 K，并首次证实在室温下光学制冷效率超过4% [5]。

对基于脉冲泵浦的固体激光制冷领域，2015年，我国哈尔滨工程大学张新陆研究团队针对铥掺杂氟化钇锂(Tm:YLF)晶体开发了脉冲激光制冷理论与实验方案[6]。并且该团队提出并论证了一种基于脉冲泵浦的固体激光制冷新方案，解决传统连续波泵浦在极低温下制冷效率低的问题[7]。2015年Andrei等人提出新型光学制冷模型比直接4f-4f泵浦的“传统”模型的冷却效率高14% [8]。随后在2016年，他们又提出了Yb³⁺离子掺杂晶体的受激拉曼绝热通道激光冷却振动模型，激光冷却的振动机制提供了更高的冷却效率和功率[9]。

本文讨论基于脉冲光泵浦的固体激光制冷，首先我们以掺Yb³⁺氟化物晶体(如Yb³⁺:YLiF₄)为例，基于密度矩阵，建立更精确的四能级模型，以处理脉冲泵浦下的相干效应和退相位过程。然后，通过系统分析脉冲泵浦下的动力学，给出制冷系统的优化激发方案。通过简化二能级模型解析求解出激发态粒子数极值对应的最优泵浦持续时间，并与四能级模型的数值解作误差比较。再结合共振吸收截面和背景吸收截面均随温度变化，模拟制冷功率密度和制冷效率随脉冲持续时间和温度的变化关系。

2. 基本原理

光学制冷的基本原理是利用通过光泵浦晶体物质，进行能量交换，吸收晶体声子能，辐射出能量更高的荧光，从而达到制冷的目的。制冷功率P_c和制冷效率η_c分别表示为：

$P_c=P_f-P_r-P_b$ (1)

${\eta }_c={\eta }_{ext}{\eta }_{abs}\frac{\lambda }{{\lambda }_f}-1$ (2)

其中，P_f为单位时间内荧光从材料中带走的热量，P_r、P_b分别为共振吸收功率密度和背景吸收功率密度。λ_p和λ_f分别为泵浦光和辐射光的波长。η_abs为对泵浦光的吸收效率，表示为：

${\eta }_{abs}=\frac{a_r}{a_r+a_b}=\frac{1}{1+\frac{a_r}{a_b}}$ (3)

其中a_r为共振吸收系数，大小与温度和波长有关。a_b为背景吸收系数。η_ext为外量子效率，表示能通过辐射跃迁将热量从材料中带出去的激发粒子的比例：

${\eta }_{ext}=\frac{{\eta }_e{w}_r}{{\eta }_e{w}_r+w_{nr}}$ (4)

η_e为荧光光子从材料中能辐射出去的比例，w_r和w_nr分别为辐射跃迁速率和非辐射跃迁速率。

从上式中利用看出要想使η_c > 0，必须让η_ext和η_abs尽可能的大，稀土离子掺杂到低声子能晶体中，其电子激发态具有较高的外量子效率η_ext (>99%)。

由于光学制冷是泵浦光与晶体物质相互作用的结果，我们以掺Yb³⁺氟化物晶体(如Yb³⁺:YLiF₄)为例，可以建立一个简易的四能级模型来描述该过程。如图1，离子能级结构中每一层都有多条斯塔克子能级组成，其中能级$|0\rangle$ 和能级$|1\rangle$ 为基态能级的两条不同的斯塔克子能级，$|2\rangle$ 和$|3\rangle$ 为激发态能级中的两条不同的斯塔克子能级，两条斯塔克子能级上的粒子数分布符合玻尔兹曼分布。制冷过程如下：(1) 用频率为ν_p的泵浦光将基态能级$|1\rangle$ 上的粒子抽运到激发态能级$|2\rangle$ 上。(2) 由于打破了激发态能级上的玻尔兹曼统计分布，导致能级$|2\rangle$ 上的粒子会吸收晶体的声子能，向能级$|3\rangle$ 进行热弛豫，恢复玻尔兹曼统计的平衡，(3) 能级$|3\rangle$ 上的激发态粒子是不稳定的，会再辐射到基态能级$|0\rangle$ ，辐射跃迁释放出频率更高为ν_f的荧光或者非辐射跃迁产生热量。(4) 能级$|0\rangle$ 上的玻尔兹曼分布被打破，电子吸收声子能量跃迁到能级$|1\rangle$ 上，完成一个制冷循环。

Figure 1. Schematic diagram of the four-level refrigeration model

图1. 四能级的制冷模型示意图

在制冷循环的过程中，主要有两条能级参与了共振吸收和跃迁过程，对于一个二能级系统，其演化满足线性主方程：

$\dot{\rho }=-\frac{i}{\hslash }\left[H,\rho \right]+{\displaystyle \sum _i^{}\left(L_i\rho L_i^{+}-\frac{1}{2}{L}_i^{+}{L}_i\rho -\frac{1}{2}\rho L_i^{+}{L}_i\right)}$ (5)

我们可以写出基于密度矩阵的动力学方程并对掺Yb³⁺和Ho³⁺的晶体的制冷系统进行建模，可以得到：

${\dot{\rho }}_{33}={\Gamma }_2{n}_2{\rho }_{22}-{\Gamma }_2\left(n_2+1\right){\rho }_{33}-\left({\beta }_{31}+{\beta }_{30}\right)\gamma {\rho }_{33}$ (6a)

${\dot{\rho }}_{22}=\text{i}\Omega \left({\rho }_{21}-{\rho }_{12}\right)+{\Gamma }_2\left(n_2+1\right){\rho }_{33}-{\Gamma }_2{n}_2{\rho }_{22}-\left({\beta }_{21}+{\beta }_{20}\right)\gamma {\rho }_{22}$ (6b)

${\dot{\rho }}_{11}=\text{i}\Omega \left({\rho }_{12}-{\rho }_{21}\right)-{\Gamma }_1\left(n_1+1\right){\rho }_{11}+{\Gamma }_1{n}_1{\rho }_{00}+{\beta }_{31}\gamma {\rho }_{33}+{\beta }_{21}\gamma {\rho }_{22}$ (6c)

${\dot{\rho }}_{00}={\Gamma }_1\left(n_1+1\right){\rho }_{11}-{\Gamma }_1{n}_1{\rho }_{00}+{\beta }_{30}\gamma {\rho }_{33}+{\beta }_{20}\gamma {\rho }_{22}$ (6d)

${\dot{\rho }}_{21}=\text{i}\Omega \left({\rho }_{22}-{\rho }_{11}\right)-\left({\gamma }_1+{\gamma }_2\right){\rho }_{21}-0.5\left[\left({\beta }_{21}+{\beta }_{20}\right)\gamma +\left(n_1+1\right){\Gamma }_1+n_2{\Gamma }_2\right]{\rho }_{21}$ (6e)

${\dot{\rho }}_{01}={\dot{\rho }}_{10}^{*},{\dot{\rho }}_{12}={\dot{\rho }}_{21}=0$ (6f)

其中${\beta }_{ij}$ 为粒子从能级$|i\rangle$ 到$|j\rangle$ 的跃迁分支比，$\gamma ={\eta }_{ext}{\gamma }_r+{\gamma }_{nr}$ 为衰减速率，${\gamma }_1$ 和${\gamma }_2$ 为退相干速率。对于脉冲泵浦，利用密度矩阵来描述四能级模型中粒子的动力学过程，可以更准确的表示出两条激发态能级之间的相干作用。

3. 泵浦持续时间的判定

Figure 2. The number of particles with different light intensities varies over time, and the light intensities are as follows: (a) 5 × 10¹⁷ W/m²; (b) 2 × 10¹⁸ W/m²; (c) 5 × 10¹⁸ W/m²; (d) 1 × 10¹⁹ W/m²

图2. 不同光强下粒子数随时间变化，光强分别为：(a) 5 × 10¹⁷ W/m²；(b) 2 × 10¹⁸ W/m²；(c) 5 × 10¹⁸ W/m²；(d) 1 × 10¹⁹ W/m²

在脉冲激发的条件下，两个参与相互作用的能级上的粒子布局数会产生拉比振荡效应，表现为该两能级的粒子数周期性反转要想实现尽可能多的制冷，要求参与制冷的粒子数尽可能的多，即在脉冲抽运结束时，激发态能级粒子数最大。所以我们需要研究：(1) 判断何时激发态粒子数最大，即确定脉冲持续时间τ_p；(2) 分析脉冲持续时间τ_p与光强的变化关系；(3) 分析激发态能级粒子数与光强最大值的变化关系。

根据给出的四能级动力学模型，为了更清晰的了解各个能级的粒子数变化，需要特别注意在泵浦过程中，两个主要能级产生的拉比振荡现象，所以我们对能级$|2\rangle$ 和能级$|1\rangle$ 上的粒子数进行模拟。设温度为300 K情况下，用不同大小的光强对掺杂浓度为10%的Yb³⁺:YLiF₄晶体进行泵浦，并对能级粒子数进行模拟，下面是光强分别为5 × 10¹⁷ W/m²、2 × 10¹⁸ W/m²、5 × 10¹⁸ W/m²和1 × 10¹⁹ W/m²条件下能级$|2\rangle$ 和能级$|1\rangle$ 上的粒子数随时间的变化。

如图2，当光强为5 × 10¹⁷ W/m²时，能级$|2\rangle$ 上的粒子数随时间先增加后趋于稳定且没有出现粒子数反转，当光强的大小在2 × 10¹⁸ W/m²附近时，能级$|2\rangle$ 和能级$|1\rangle$ 上的粒子数会有明显的拉比振荡，并伴随有粒子数反转的现象。继续提高光强到5 × 10¹⁸ W/m²和1 × 10¹⁹ W/m²时，拉比振荡的次数会更多，振荡会更加明显，并伴有多次粒子数反转，然后慢慢趋于稳定。

当能级间出现显著的拉比振荡时，我们可以尝试先将四能级模型简化为二能级模型，即主要考虑参与泵浦的两条能级，如下图3。

Figure 3. Simple two-level model

图3. 简易二能级模型

由于脉冲泵浦的时间(皮秒级)远快于激发态粒子的寿命(毫秒级)，即脉冲抽运的速度远快于辐射跃迁的速度，所以当我们只讨论泵浦阶段时，我们可以暂时将激发态粒子的辐射跃迁忽略。基于密度矩阵的动力学方程如下所示。

${\dot{\rho }}_{22}=\text{i}\Omega \left({\rho }_{21}-{\rho }_{12}\right)$ (7a)

${\dot{\rho }}_{11}=\text{i}\Omega \left({\rho }_{12}-{\rho }_{21}\right)$ (7b)

${\dot{\rho }}_{21}=\text{i}\Omega \left({\rho }_{22}-{\rho }_{11}\right)-{\gamma }_1{\rho }_{21}$ (7c)

${\dot{\rho }}_{12}=\text{i}\Omega \left({\rho }_{11}-{\rho }_{22}\right)-{\gamma }_1{\rho }_{12}$ (7d)

${\dot{\rho }}_{12}={\dot{\rho }}_{21}=0$ (7e)

其中γ₁为退相干速率，我们设最初基态粒子数为A，则激发态粒子数即为1-A，并设：

$\sqrt{{\gamma }_1^2-16{\Omega }^2}=2i{\Omega }_t$ (8)

其中的等效拉比振荡Ω_t为实数，显而易见的要求：

${\gamma }_1^2-16{\Omega }^2<0$ (9)

其中Ω为拉比振荡频率，表示为：

$\Omega =\frac{\mu E}{2\hslash }$ (10)

其中μ为电偶极矩，E为电场振幅，光强I可由电场振幅来表示：

$I=\frac{nc{\epsilon }_0{E}^2}{2}$ (11)

其中，n为相互作用折射率，ε₀为真空中的电导率，再结合式(8)、式(9)和式(10)可以得到要求：并由此可以给出激发态粒子数的解析解为：

$I>\frac{{\gamma }_1^2{\hslash }^{2}nc{\epsilon }_0}{8{\mu }^2}$ (12)

并由此可以给出激发态粒子数的解析解为：

${\rho }_{22}=\frac{1}{2}-\frac{\left(-1+2A\right){\text{e}}^{-\frac{t{\gamma }_1}{2}}\left[2{\Omega }_t\cos \left(t{\Omega }_t\right)+{\gamma }_1\sin \left(t{\Omega }_t\right)\right]}{4{\Omega }_t}$ (13)

对上式求导可得：

${\dot{\rho }}_{22}=\frac{\left(1-2A\right){\text{e}}^{-\frac{t{\gamma }_1}{2}}\left({\gamma }_1^2+4{\Omega }_t^2\right)\sin \left(t{\Omega }_t\right)}{8{\Omega }_t}$ (14)

由上式可以看出，当

$t=\frac{\left(2n+1\right)\pi }{{\Omega }_t}$ (15)

时激发态粒子数取极值，且由于存在衰减机制，拉比振荡的峰值会随时间演化会越来越小，粒子数趋于常数，所以在第一次极值时，即n = 0时，激发态能级粒子数取最大值。这个结论非常重要，它指出，在存在退相干的情况下，拉比振荡的幅度是衰减的，因此第一个峰值通常是整个演化过程中的最大值。为了在脉冲结束时获得最多的激发态粒子，最优的脉冲持续时间应该对应于第一个拉比振荡峰值出现的时间。这为脉冲激光器的参数设定提供了直接的理论依据。

根据以上分析，我们已经知道在τ_p1 = π/Ω_t时，激发态粒子数取最大值。为了确定持续时间τ_p1的具体数值，且τ_p1与光强I有关，所以我们需要讨论在不同温度下，激发态能级粒子数最大值和脉冲持续时间τ_p1与光强的变化关系。把n = 0和式(12)带入式(8)中可得激发态粒子数与等效拉比振荡Ω_t的关系：

${\rho }_{22}=\frac{1}{2}+\frac{\left(-1+2A\right){\text{e}}^{-\frac{\pi {\gamma }_1}{{\Omega }_{t}2}}}{2}$ (16)

因为等效拉比振荡Ω_t与光强I_p呈正比关系，当光强不断增大的情况下，等效拉比振荡Ω_t不断增大，激发态粒子数无限接近静态基态能级粒子数A。

下面为当温度分别为150 K、200 K、250 K和300 K时下激发态能级粒子数第一次达到峰值时，脉冲持续时间τ_p1随光强I_p的变化关系。

如图4所示，脉冲持续时间τ_p1随泵浦光强的变大先快速变小，最后趋于稳定。不同温度下，脉冲持续时间τ_p1的变化趋势保持一致，温度越高，脉冲持续时间τ_p整体越小。

Figure 4. Relationship between pulse duration τ_p and light intensity I_p

图4. 脉冲持续时间τ_p随光强I_p的变化关系

下图为温度分别为150 K、200 K、250 K和300 K时激发态能级粒子数第一次达到峰值时，能级粒子数极值与光强关系：

Figure 5. Relationship between extreme particle number and light intensity I_p: (a) Excited state, (b) Ground state

图5. 粒子数极值随光强I_p的变化关系：(a) 激发态，(b) 基态

从图5中可以看到，激发态能级粒子数的极值会先随泵浦光强I_p的变大而变大，最后趋于稳定，在τ_p1时基态能级粒子随泵浦光强I_p的变大先减小趋于稳定。

综上我们通过解析解，就可以先确定泵浦光的脉冲持续时间，合理的去调整光强的大小，达到一个合理的泵浦条件，使激发态粒子数达到一个极值，得到最佳的制冷性能。

但是，由于上述模拟结果本质是简化了能级结果，并忽略了能级辐射跃迁再通过数值推导计算出的解析解，与实际结果存在着误差。为了严谨性，需要讨论解析解与数值解之间误差的大小。我们利用基于密度矩阵的四能级动力学模型，对Yb³⁺:YLiF₄晶体的粒子数随时间变化进行模拟再取其数值解的极值，并于数值解进行比较。

Figure 6. Numerical and analytical solutions of particle number at different temperatures: (a) 300 K; (b) 250 K; (c) 200 K; (d) 150 K

图6. 不同温度下粒子数的数值解与解析解：(a) 300 K；(b) 250 K；(c) 200 K；(d) 150 K

Figure 7. Difference between numerical and analytical solutions of extreme particle numbers in excited states at different temperatures

图7. 不同温度激发态粒子数极值的数值解与解析解差值

如图6所示，其中的实线为取得数值解，虚线为解析解的数值，可以看到两条曲线差距很小，并随着光强越大，数值解与解析解之间的差值越来越小。为了更清晰展示差值的大小，继续用两者之间的差值占数值解的占比来随光强的变化来表示。

如图7所示，不同温度下粒子数极值的解析解与数值解的误差占比不到1%，所以这微小的误差我们可以忽略不计。

4. 制冷模拟分析

本节我们将探讨共振吸收截面和背景吸收截面与温度的变化关系，以及对吸收效率以及制冷效率的影响。

对于背景吸收，在近红外波段，声子吸收占主导地位，背景吸收截面可表示为：

${\sigma }_i={\sigma }_0\exp \left(-\frac{E_{ph}}{k_{B}T}\right)$ (17)

其清晰的变化规律如图8所示。

Figure 8. Relationship between background absorption cross-section σ_i and temperature

图8. 背景吸收截面σ_i随温度的变化关系

为了更清晰的讨论共振吸收截面与温度的变化关系，可以选择用吸收光谱法测量吸收截面与温度的变化关系，并根据测量出来的曲线拟定一个近似函数[10]：

${\sigma }_a=b_0+b_{1}T+b_2{T}^2+b_3{T}^3+b_4{T}^4+b_5{T}^5+b_6{T}^6$ (18)

其中$b_i\left(i=1-6\right)$ 的值与波长相关，再通过查表，我们可以确定在泵浦光波长为1020 nm时$b_i\left(i=1-6\right)$ 的值分别为−5.676、2.379、−3.940、3.265、−1.420、3.166和−2.860。所选用的六阶多项式提供了足够的灵活性来精确拟合实验曲线。通过此函数模拟两者的变化曲线为：

Figure 9. Relationship between resonance absorption cross section σ_a and temperature

图9. 共振吸收截面σ_a随温度的变化关系

由图9可知无论是背景吸收截面σ_i还是背景吸收截面σ_a都是跟温度呈正相关的。分别代入可得不饱和吸收情况下吸收效率η_abs跟温度相关的图像，如图10：

Figure 10. Relationship between absorption efficiency η_abs and temperature

图10. 吸收效率η_abs随温度的变化关系

在用连续进行抽运时，稀土离子的寿命τ_r为毫秒级，杂质离子的寿命τ_i为微秒级。设连续光的光强为5 × 10⁷ W/m²，τ_r为1.5 ms，τ_i为25 μs；在超快脉冲泵浦中，τ_r和τ_i均为τ_p，设光强为5 × 10¹⁸ W/m²，脉冲持续时间为1ps情况下，可分别得到连续抽运和脉冲抽运下吸收效率跟温度相关的图像：

Figure 11. Relationship between absorption efficiency η_abs and temperature under continuous pumping and pulse pumping conditions

图11. 连续抽运和脉冲抽运条件下吸收效率η_abs随温度的变化关系

从图11中可以看到连续抽运条件下吸收效率≥96.5%，并随温度增大先减小然后趋于稳定。而脉冲抽运条件下，由于对背景吸收的降低，是得吸收效率得到提高(≥98.8%)，并随温度增大而减小。

为了更直观的表现分析脉冲持续时间和温度共同对制冷的影响，我们可以分别以两者为x、y坐标综合得到制冷效率和制冷功率密度的三维模拟图。

Figure 12. Relationship between extreme cooling power density and pulse duration and temperature

图12. 制冷功率密度极值随脉冲持续时间和温度的变化关系

如图12所示，脉冲时间较小且温度较高时，得到制冷功率密度的极大值。

Figure 13. Relationship between extreme refrigeration efficiency density and pulse duration and temperature

图13. 制冷效率极值随脉冲持续时间和温度的变化关系

如图13所示，在饱和吸收的情况下，制冷功率的极大值与脉冲持续时间无关，随温度变大而变大。

5. 总结

本文首先介绍了光学制冷基本原理，建立简易的四能级模型，并基于密度矩阵建立的动力学方程。为了确定不同制冷条件下的泵浦持续时间，使参与制冷的粒子数尽可能多，粒子数由于拉比振荡会导致有波动性，我们先对不同光强下，激发态粒子数进行模拟，得到不同振荡幅度的粒子数情况。为了方便分析，我们通过简化的二能级模型，得到激发态粒子数的解析解，并求极值的得到了脉冲持续时间与光强的变化关系。并与四能级方程模拟得到的数值解进行误差分析，发现误差极小以证明方案的可行性。其次，我们引入共振吸收界面和背景吸收截面与温度的变化关系的函数，模拟发现两者均随温度变化而升高，并对吸收效率进行分析。最后，对Yb³⁺:YLiF₄晶体基于脉冲泵浦下进行制冷模拟，得到制冷功率密度极值和制冷效率极值随泵浦持续时间的变化关系并进行总结：随着脉冲持续时间不断提高，不同温度下制冷功率密度极值大致为先降低后升高的趋势，随着脉冲持续时间不断提高，不同温度下制冷效率极值均保持不变；随着温度不断提高，不同脉冲持续时间下制冷效率和制冷功率密度均随着温度的提高而提高。