1.   引　言

一直以来, 少电子原子分子体系(氢原子、氦原子、氢分子离子等)的精密光谱研究为检验基本物理理论、测定基本物理常数(里德伯常数${R_\infty }$、精细结构常数$\alpha$等)、确定原子核的电荷分布半径和磁矩分布半径等提供了有效的方法^[1-4]. 锂原子以及类锂的铍离子(Be⁺)作为三电子原子体系, 受到了理论和实验研究者的广泛关注. 锂原子精密光谱的一个重要应用是确定原子核的性质, 如核电核分布半径和磁矩分布半径(Zemach半径)等. 例如, Yan和Drake^[5]发展了锂原子能级的精密计算理论和方法. 根据他们的方法, 锂原子最低的三个能级(2S, 3S和2P)间的跃迁频率的同位素移动的高精度理论计算值和实验值的对比被用于确定锂原子核同位素的电荷分布半径, 特别是具有“晕核”结构的¹¹Li的电荷分布半径^[4] (这里的电荷分布半径指的是相对于某个电荷半径为已知的同位素的相对电荷分布半径, 目前的理论计算值的精度还不足以确定锂原子核的绝对电荷分布半径). 对于铍离子(Be⁺), 也有类似的应用^[4]. Puchalski和Pachucki^[6]则发展了锂原子超精细结构的高精度计算理论和方法, 根据他们计算的锂原子基态的超精细结构劈裂理论值与实验值的比较发现, ⁷Li原子核的Zemach半径比⁶Li的大40%, 为原子核结构的研究提出了一个重要的问题. 最近, 基于锂离子(Li⁺)的超精细结构的研究也得到了相似的结论^[7].

超精细结构是原子精密谱研究的一个重要内容. 对于原子体系的磁偶极超精细结构来说, 其最主要的来源是费米接触相互作用. 因此, 要想得到高精度的超精细结构劈裂, 必须要把费米接触项计算到足够高的精度^[6,7]. 对于类氢原子体系, 费米接触项(Fermi contact term, ${f_{\text{c}}}$)可以由解析计算获得; 而对于含有两个电子或更多电子的原子体系, ${f_{\text{c}}}$只能通过数值计算得到, 其数值精度严重依赖波函数在原点处的质量. 除此之外, 电子间的关联效应也会对${f_{\text{c}}}$产生重要影响. 通过使用全局化方法^[8,9], 可以降低${f_{\text{c}}}$对波函数在原点附近行为的依赖性; 而通过使用Hylleraas型变分波函数, 则可以充分照顾到电子关联效应对${f_{\text{c}}}$的影响^[10]. 目前, 利用全局化方法和Hylleraas变分波函数, 锂原子基态(${\text{1}}{{\text{s}}^{\text{2}}}{\text{2s}} {\;^{\text{2}}}{\text{S}}$)非相对论${f_{\text{c}}}$的理论计算值的精度已经达到了10^–8[6,11]; 对于锂原子的其他几个较低的自旋双重态, ${f_{\text{c}}}$的值也可以达到很高的精度^[11]. 对于锂原子的自旋四重态, Yan^[12]用Hylleraas变分波函数计算了${\text{1s2s3s}}{\;^{\text{4}}}{\text{S}}$, ${\text{1s2s4s}}{\;^{\text{4}}}{\text{S}}$和${\text{1s2s5s}}{\;^{\text{4}}}{\text{S}}$这三个态的${f_{\text{c}}}$, 精度高达10^–8; King^[13]也做了类似的计算(包括$1{\text{s}}2{\text{s}}6{\text{s}}{\;^4}{\text{S}}$), 但是与Yan的结果相比, 精度稍差. Zhuo及其合作者^[14,15]用组态相互作用方法计算了类锂离子体系(Z = 7—10)几个较低的自旋四重P态(${\text{1s}}nl nl'{\;^{\text{4}}}{{\text{P}}^{{o}}}$和${\text{1s2p}}n{\text{p}}{\;^{\text{4}}}{\text{P}}$)的费米接触项以及其他超精细结构常数. 对于${f_{\text{c}}}$, 他们的结果的精度在10^–5(通过和本研究小组的结果对比而发现). 本文对锂原子及类锂离子体系(Z = 4—10)的${\text{1s2s3s}} {\;^{\text{4}}}{\text{S}}$, ${\text{1s2s4s}}{\;^{\text{4}}}{\text{S}}$和${\text{1s2s2p}}{\;^{\text{4}}}{\text{P}}$的费米接触项做了计算, 一方面弥补了在Z = 4—6这个区域中相关数据的缺失; 另一方面, 对于文献上已有的结果, 在精度上有明显的提高(1—3个量级). 除此之外, 还研究了原子核的有限质量对${f_{\text{c}}}$的影响. 本文所有${f_{\text{c}}}$的精度都在10^–10以上, 可以作为其他计算方法的基准, 同时也为相关的实验研究提供了有价值的参考数据.

2.   理论方法

2.1   非相对论哈密顿量

在质心坐标系中, 三电子原子体系的非相对论哈密顿量为^[16] (本文使用原子单位, 即$\hbar = e = {m_{\text{e}}} = 1$)

$$H = \sum\limits_{i = 1}^3 {\left( { - \frac{1}{2}\nabla _i^2 - \frac{Z}{{{r_i}}}} \right) + \sum\limits_{i < j}^3 {\frac{1}{{{r_{ij}}}}} } - \frac{\mu }{M}\sum\limits_{i < j}^3 {{\nabla _i} \cdot {\nabla _j}} ,$$

(1)

其中, $M$代表原子核的质量, $\mu = {{{m_{\text{e}}}M} \mathord{\left/ {\vphantom {{{m_{\text{e}}}M} {({m_{\text{e}}} + M)}}} \right. } {({m_{\text{e}}} + M)}}$代表电子的折合质量. (1)式中能量的单位是$2{R_M}$, 而${R_M} = (1 - {\mu \mathord{\left/ {\vphantom {\mu M}} \right. } M}){R_\infty }$. 原子核的有限质量对能量的影响通过能量单位$2{R_M}$(mormal mass shift)和(1)式中的最后一项(specific mass shift)体现出来.

2.2   三电子体系Hylleraas型变分波函数

本文使用瑞利-李兹变分法求解体系的薛定谔方程, 选择三电子体系的Hylleraas型基函数构造变分波函数. Hylleraas型基函数是目前求解两电子和三电子原子体系薛定谔方程的精度最高的基函数^[17]. 把体系的变分波函数做如下展开:

$$\psi = \sum\limits_{i = 1}^N {{c_i}} {\varPhi _i}({{\boldsymbol{r}}_1},{{\boldsymbol{r}}_2},{{\boldsymbol{r}}_3}) ,$$

(2)

而

$${\varPhi _i}({{\boldsymbol{r}}_1},{{\boldsymbol{r}}_2},{{\boldsymbol{r}}_3}) = \mathcal{A}\left[ {{\phi _i}({{\boldsymbol{r}}_1},{{\boldsymbol{r}}_2},{{\boldsymbol{r}}_3})\chi (1,2,3)} \right] ,$$

(3)

其中, $\mathcal{A}$表示三电子体系的反对称化算符, ${\phi _i} ({\boldsymbol r}_1, {\boldsymbol r}_2, {\boldsymbol r}_3)$表示基函数的轨道部分:

$${\phi _i}({{\boldsymbol{r}}_1},{{\boldsymbol{r}}_2},{{\boldsymbol{r}}_3})= r_1^{{j_1}}{r_2}^{{j_2}}{r_3}^{{j_3}}r_{12}^{{j_{12}}}r_{23}^{{j_{23}}}r_{31}^{{j_{31}}}{{\rm{e}}^{ - \alpha {r_1} - \beta {r_2} - \gamma {r_3}}}\mathcal{Y}_{({l_1}{l_2}){l_{12}},{l_3}}^{LM}({{\boldsymbol{r}}_1},{{\boldsymbol{r}}_2},{{\boldsymbol{r}}_3}) ,$$

(4)

而

$$\begin{split} \;&\mathcal{Y}_{({l_1}{l_2}){l_{12}},{l_3}}^{LM}({{\boldsymbol{r}}_1},{{\boldsymbol{r}}_2},{{\boldsymbol{r}}_3}) \\ =\;& {r_1}^{{l_1}}{r_2}^{{l_2}}{r_3}^{{l_3}} \sum\limits_{all\,{m_i}} {\left\langle {{{l_1}{m_1};{l_2}{m_2}}} \mathrel{\left | {\vphantom {{{l_1}{m_1};{l_2}{m_2}} {{l_1}{l_2};{l_{12}}{m_{12}}}}} \right. } {{{l_1}{l_2};{l_{12}}{m_{12}}}} \right\rangle } \left\langle {{{l_{12}}{m_{12}};{l_3}{m_3}}} \mathrel{\left | {\vphantom {{{l_{12}}{m_{12}};{l_3}{m_3}} {{l_{12}}{l_3};LM}}} \right. } {{{l_{12}}{l_3};LM}} \right\rangle {Y_{{l_1}{m_1}}}({{\boldsymbol{r}}_1}){Y_{{l_2}{m_2}}}({{\boldsymbol{r}}_2}){Y_{{l_3}{m_3}}}({{\boldsymbol{r}}_3}), \end{split}$$

(5)

表示三个角动量函数耦合出的基函数的角向部分; $\chi (1, 2, 3)$表示自旋函数, 对于自旋四重态($S = {3 \mathord{\left/ {\vphantom {3 2}} \right. } 2}, \;{M_S} = {3 \mathord{\left/ {\vphantom {3 2}} \right. } 2}$), 有

$$\chi (1,2,3) = \alpha \left( 1 \right)\alpha \left( 2 \right)\alpha \left( 3 \right) .$$

(6)

使用(2)式的变分波函数, 根据瑞利-李兹变分法^[17], 可以得到广义本征值方程:

$${\boldsymbol{Hc}} = E{\boldsymbol{Oc}} ,$$

(7)

其中${\boldsymbol{H }}$和${\boldsymbol{O}}$分别表示哈密顿矩阵和交叠矩阵. 求解广义本征值方程(7), 即可得到体系的能量本征值和本征波函数的最佳近似值. 为了提高变分能量的收敛率, 参考Yan的方法^[12], 根据以下条件把变分基底分为5组:

$$\begin{split}& \text{sector}\text{1}: \quad\text{all}\;{j}_{12},\quad{j}_{23}=0,\quad{j}_{31}=0;\\ & \text{sector}2: \quad\text{all}\;{j}_{12},\quad{j}_{23}=0,\quad{j}_{31}\ne 0;\\ & \text{sector}3: \quad\text{all}\;{j}_{12},\quad{j}_{23}\ne 0,\quad{j}_{31}=0;\\ & \text{sector}4:\, {j}_{12}=0,\quad{j}_{23}\ne 0,\quad{j}_{31}\ne 0;\\ &\text{sector}5:\, {j}_{12}\ne 0,\quad{j}_{23}\ne 0,\quad{j}_{31}\ne 0.\end{split}$$

(8)

每组基底的规模用下面的条件控制:

$${j_1} + {j_2} + {j_3} + {j_{12}} + {j_{23}} + {j_{31}} \leqslant {\varOmega _i},$$

(9)

其中, ${j_\mu }$, ${j_{\mu \nu }}$和${\varOmega _i}$均取非负整数.

2.3   费米接触项与全局化方法

费米接触项的定义为

$${f_{\text{c}}} = 4{\text{π}}\left\langle {\psi \left| {\sum\limits_{i = 1}^3 {\delta \left( {{{\boldsymbol{r}}_i}} \right){{\boldsymbol{\sigma}} _{i z}}} } \right|\psi } \right\rangle ,$$

(10)

其中${{\boldsymbol{\sigma}} _z}$是泡利矩阵. 从(10)式可知, ${f_{\text{c}}}$的数值精度对体系波函数在坐标原点处的取值有比较强的依赖关系, 这在数值计算中是不利的. 为了克服这一缺点, Hiller及其合作者以及Drachman发展了“全局化方法”^[8,9]. 本文采用Drachman的全局化方法^[9], 其主要思想是用一个等价的表达式代替$\delta$函数的期望值, 即

$$\begin{split} \left\langle \psi \right.\left| {\delta \left( {{{\boldsymbol{r}}_i}} \right)} \right|\left. \psi \right\rangle = \;& - \frac{1}{{\text{π }}}\left\langle {\psi \left| {\frac{1}{{{r_i}}}\left( {{V_{\text{d}}} - E} \right)} \right|\psi } \right\rangle \\& - \frac{1}{{{\text{2π}}}}\sum\limits_{s = 1}^3 {\left\langle {{\nabla _s}\psi \left| {\frac{1}{{{r_i}}}} \right|{\nabla _s}\psi } \right\rangle } , \end{split}$$

(11)

其中

$${V_{\text{d}}} = - \sum\limits_{i = 1}^3 {\frac{Z}{{{r_i}}}} + \sum\limits_{i < j}^3 {\frac{1}{{{r_{ij}}}}} - \frac{\mu }{M}\sum\limits_{i < j}^3 {{\nabla _i} \cdot {\nabla _j}} .$$

(12)

不难发现, 利用全局化方法计算得到的${f_{\text{c}}}$比直接利用定义(10)式计算在精度上可以提高1—2个量级.

3.   结果与讨论

使用2.2节中描述的变分基底, 首先用瑞利-李兹变分法求解锂原子和类锂离子体系(Z = 4—10)的${\text{1s2s3s}}{\;^{\text{4}}}{\text{S}}$, ${\text{1s2s4s}}{\;^{\text{4}}}{\text{S}}$和${\text{1s2s2p}}{\;^{\text{4}}}{\text{P}}$的薛定谔方程, 得到了这些体系的变分能量和变分波函数. 把变分基底的规模从$N{\text{ = }}210$($N$表示基函数的数目)逐步扩大到$N = 27130$, 在每一组基底下, 用牛顿法对5组非线性参数进行充分优化, 以得到尽可能低的变分能量. 在$N < 10000$的情况下, 优化参数能够带来比较明显的效果, 但是, 随着变分基底的进一步增大, 优化参数带来的能量降低很有限. 这时, 主要靠增大基底的规模带来能量的进一步降低.

表1列出了锂原子${\text{1s2s3s}}{\;^{\text{4}}}{\text{S}}$态的变分能量(原子核质量取为无穷大)在不同基底下的值, 从收敛序列可以确认变分能量的收敛精度达到了10^–13. 对于其他两个态, 变分能量也达到了相似的精度. 表2列出了从锂原子到Ne⁷⁺离子的${\text{1s2s3s}}{\;^{\text{4}}}{\text{S}}$, ${\text{1s2s4s}}{\;^{\text{4}}}{\text{S}}$和${\text{1s2s2p}}{\;^{\text{4}}}{\text{P}}$的变分能量, 并与文献上已有的结果进行了对比. 相较于Yan^[12]的结果, 本研究组的结果在精度上提高了1—3个量级; 相较于Barrois等^[18]的结果, 则提高了5—6个量级. 精度提高的主要原因是使用了更多的变分基函数.

表1和表2中的变分能量都是在原子核质量取无穷大的假设下得到的, 即在方程(1)中取${\mu \mathord{\left/ {\vphantom {\mu M}} \right. } M} = 0$. 对于原子核质量取有限值的情况, 可以通过直接对角化方程(1)中的哈密顿量来得到体系的能量, 也可以用微扰法来处理有限核质量带来的影响, 即把方程(1)中的最后一项当作微扰算符来处理. 为了得到各种同位素的能量, 使用计算公式

$${E_M} = {\varepsilon _0} + \left( {{\varepsilon _0} + {\varepsilon _1}} \right)\lambda + \left( {{\varepsilon _1} + {\varepsilon _2}} \right){\lambda ^2} + \left( {{\varepsilon _2} + {\varepsilon _3}} \right){\lambda ^3} + \cdots,$$

(13)

其中, ${E_M}$表示核质量为$M$的同位素的变分能量; ${\varepsilon _0}$表示核质量取无穷大时体系系的能量; ${\varepsilon _1}$, ${\varepsilon _2}$和${\varepsilon _3}$分别表示一阶、二阶和三阶质量极化系数; $\lambda = - {\mu \mathord{\left/ {\vphantom {\mu M}} \right. } M}$. 质量极化系数可以通过微扰论计算得到. 用微扰论计算了一阶、二阶和三阶质量极化系数, 数据列于表3. 根据表2和表3中的数据, 只要知道了核质量, 就可以根据方程(13)计算出该同位素的非相对论能量.

根据变分法得到的高精度波函数以及方程(10)和方程(11), 即可得到${f_{\text{c}}}$的数值. 表4列出了锂原子的${f_{\text{c}}}$随基底逐渐增大而收敛的情况. 根据收敛序列, 可以确定${f_{\text{c}}}$的收敛精度达到了10^–10. 表5列出的是从锂原子到Ne⁷⁺三个最低的自旋四重态的${f_{\text{c}}}$. 这些结果的收敛精度普遍达到了10^–10. 与早期使用Hylleraas变分波函数计算的锂原子费米接触项^[12,21]相比, 本研究组的结果在精度上提高了2—5个量级. 精度提高的原因有两点: 一是在变分波函数中包括了更多的基函数; 二是使用了Drachman全局化方法^[9], 这种方法有利于加速结果的收敛. 通过比较Zhuo等^[14,15]使用组态相互作用方法计算的N⁴⁺, O⁵⁺, F⁶⁺和Ne⁷⁺的${\text{1s2s2p}}{\;^{\text{4}}}{\text{P}}$态的费米接触项, 发现他们的结果的精度在10^–5左右.

表4和表5中的${f_{\text{c}}}$是在核质量取无穷大的情况下得到的, 与核的有限质量对变分能量的影响类似, 核的有限质量对${f_{\rm{c}}}$的影响也可以通过质量极化系数表示为

$$f_{\text{c}}^M = f_{\text{c}}^0 + \left( {3f_{\text{c}}^0 + f_{\text{c}}^1} \right)\lambda ,$$

(14)

其中, $f_{\text{c}}^0$表示无穷大核质量情况下计算得到的费米接触项, $f_{\text{c}}^M$表示有限核质量下的费米接触项, $f_{\text{c}}^1$表示费米接触项的一阶质量极化系数. 方程(14)精确到$O(\lambda )$. 计算得到的费米接触项的一阶质量极化系数列于表6. 根据表5和表6所列的数据可以看出, 原子核的有限质量对${f_{\text{c}}}$的影响主要来自于$3 f_{\text{c}}^0\lambda$, 即normal mass shift, 质量极化效应$f_{\text{c}}^1\lambda$的影响相对来说要小很多.

4.   结　论

本文利用Hylleraas坐标下的瑞利-李兹变分法对锂原子及类锂离子(Z = 4—10)的自旋四重态${\text{1s2s3s}}{\;^{\text{4}}}{\text{S}}$, ${\text{1s2s4s}}{\;^{\text{4}}}{\text{S}}$和${\text{1s2s2p}}{\;^{\text{4}}}{\text{P}}$的非相对论薛定谔方程进行了高精度求解, 得到的变分能量的精度达到了10^–13, 并给出了变分能量的一阶、二阶和三阶质量极化系数. 与文献上已有的结果相比, 本文的结果在精度上提高了1—6个数量级. 使用高精度的变分波函数, 计算了这些态的费米接触项, ${f_{\text{c}}}$的精度达到了10^–10, 是目前已有结果中最精确的. 本文的结果可以作为其他计算方法参照的基准, 也为以后要开展的相对论修正和QED修正计算工作打下了基础, 同时也为相关的实验研究提供了有价值的参考数据.