原子核处电子密度、Isomer shifts、Mössbauer 四极分裂与 Mössbauer 谱

关键词：EFG-Tensor；Q-Tensor；Isomer shifts；Mössbauer quadrupole splittings；穆斯堡尔谱

材料中的某些原子有一定的几率能够发射 γ 射线，而该原子不受到反冲，处于基态的其他原子对该 γ 射线也有一定的几率产生共振吸收，并且也不受到反冲。原子的这种无反冲发射或共振吸收 γ 射线的现象，叫做 Mössbauer 效应。

当 γ 射线照射某种材料时，如果 γ 光子的能量与材料中某些原子的能级跃迁能量相等，γ 光子就有一定几率被该原子共振吸收。能量差异较大的 γ 光子，一般不会产生共振吸收现象。吸收后测得的 γ 光子的数量与能量的对应关系，就是 Mössbauer 谱。由于 Mössbauer 效与原子核与内层电子性质相关，包括原子核与芯电子能级结构，及其所处化学环境，因此可以应用 Mössbauer 谱研究原子的价态、化学键的离子性与配位数、晶体结构、电子密度和磁性质等。

默认情况下，计算原子核处的电子密度，不需要额外的输入，任何计算都会默认输出这些数值，因此用户只需要考虑合适的高精度参数，进行结构优化、单点计算即可。如下文所示，ADF 输出的原子核处电子密度，并不是在原子核的精确中心点的电子密度值（精确中心的值受基组影响非常剧烈，不宜采用），而是在原子核中心一个微小半径上的小球表面上的点的电子密度求均值得到的。

原子核处的电子密度，可用于解释 Mössbauer 谱中的Isomer shifts，通常需要对实验测量的Isomer shifts和计算的电子密度进行拟合，具体可以参考文献中所述：

• W.-G. Han, T. Liu, T. Lovell and L. Noodleman, DFT calculations of 57Fe Mössbauer isomer shifts and quadrupole splittings for iron complexes in polar dielectric media: Applications to methane monooxygenase and ribonucleotide reductase, Journal of Computational Chemistry 27, 1292 (2006)

对于 Mössbauer 活性的元素，应该使用全电子基组。如果比较不同分子中原子核的电子密度，则应使用相同的参数，包括基组、交换相关函数、积分精度和核模型（AMSinput → Details → Relativity → Nuclear Model：点电荷模型/核电荷Gaussian分布模型）。原子核处的电子密度在很大程度上取决于该原子核核心区域中基组的准确性，尤其是当考虑了相对论效应的时候。

如果相对论方法使用 ZORA 方法（AMSinput → Details → Relativity → Formalism：ZORA），则程序调用的是 ZORA-4 计算原子核处的电子密度，该方法包括“小分量”密度，而老版本中使用的是 ZORA 密度，不包括小分量密度。2020版及其以上使用的 ZORA-4 密度。

在ADF2017中，ZORA 默认使用 MAPA，而不是 SAPA，与SAPA相比，MAPA的优势在于减少了ZORA的规范依赖性。不过除重原子核的电子密度外，几乎所有性质的 ZORA 规范依赖性都很小。

与相对论 Dirac 方法相比，ZORA 计算给出的原子核密度太大，过大的来源主要是来自 1s 轨道。不过，如果对实验测量的 Isomer shifts 与计算的电子密度进行拟合，那么 ZORA 电子密度过大的问题就不重要了。

相对论 X2C 方法（AMSinput → Details → Relativity → Formalism：X2C）得到的原子核处电子密度的绝对值，则显著更佳。此时需要注意使用具有更紧密的 STO 1s 函数的大基组，来获得更准确的结果（TZ2P-J或QZ4P-J）。

本文采用AMS2020.101完成计算。

ADF采用Slater基组，能够很好地描述电子在原子核附近的渐进行为。而Gaussian基组则不能很好的描述。至少，原子核处电子波函数一阶导数是不连续的，Gaussian基组得到的电子波函数一阶导数是连续的。

模型创建好之后，点击窗口底部的★符号，程序计算的时候，才能自动识别点群，同时计算效率也会大大提高。参数的设置，可以参考ADF参数设置详解。

ADF计算原子核处的电子密度，设置非常简单，因为默认就会计算。因此做一个Single Point计算就会得到原子核处电子密度。但需要设置非常高的精度，例如下图所示，计算HgCl₂的原子核处电子密度：

该功能也支持开壳层计算，即Charge、Spin polarization可以根据需要设置，如果后者不等于0，则需要勾选Unrestricted。需要注意的是，基组的Frozen Core一定要设置为None。

另外，如果要非常精确的计算，还可以设置原子核的模型为Gaussian型，而不是点电荷模型，相对论方法也可以选择X2C。

计算 Mössbauer quadrupole splittings，则需要勾选EFG Q-tensor。如果只需要分析部分原子，则选中这些原子点击 EFG atoms 后面的 + 即可。

保存并运行作业。

SCM - Output - Properties - Electron Density at Nuclei显示下面一段：

 Electron Density at Nuclei
 
 
 The electron density is calculated at points on a small sphere around the center of a nucleus.
 The printed electron density is the average electron density on these points.
 The radius of the sphere is the printed approximate finite nuclear radius.
 Note: density is calculated in the Dirac picture
 
         Atom        Nuclear Radius (Angstrom)             Electron Density (a.u.)
         ----        -------------------------             -----------------------
     1)   Cl                0.0000426624                             3299.78224
     2)   Hg                0.0000706845                          1545276.89486
     3)   Cl                0.0000426624                             3299.78224

即是。

Nuclear Electric Quadrupole Hyperfine 相互作用即 Q-tensor，而 Q-tensor 与 EFG 密切相关。 Q-tensor 张量矩阵元（单位MHz）

在原子核处电子密度网上翻，就可以看到所有原子的EFG-tensor 和 Q-tensor（由于该分子具有对称性，两个Cl原子等价，因此此处只列出了一个Cl原子），也可以在 SCM → Output → Other properties → EFG Q-Tensor 看到：

Cl原子

HYPERFINE INTERACTIONS  ***  EFG and ESR Q-TENSOR  ***

  
 Small component electron density is included in the calculation of EFG and ESR Q-tensor
 Note: density is calculated in the Dirac picture


 ------------------------------------------------
  isotope             =  35Cl 
  nuclear spin I      =  1.5
  magnetic moment mu  =  0.8218743
          g_n = mu/I  =  0.5479162
  quadrupole moment Q = -0.08165 e 10-24 cm2
 ------------------------------------------------

 Atom    1        Cl Input number    1  xyznuc(Angstrom) =     0.0000    0.0000   -2.2599

 ==== EFG-tensor (a.u.) Cl
            X             Y             Z
     X   -0.141751E+01  0.000000E+00  0.000000E+00
     Y    0.000000E+00 -0.141751E+01  0.000000E+00
     Z    0.000000E+00  0.000000E+00  0.283503E+01

 ==== principal axes Q-tensor (EFG) Cl
           11            22            33
     X    1.000000      0.000000      0.000000
     Y    0.000000      1.000000      0.000000
     Z    0.000000      0.000000      1.000000
 principal values EFG Cl
          V11           V22           V33
         -0.141751E+01 -0.141751E+01  0.283503E+01  a.u.
         -0.137745E+02 -0.137745E+02  0.275490E+02  10^21 V/m^2
 principal values  35Cl  Q-tensor (includes factor Q/(2I(2I-1)))
          q11           q22           q33
          0.453247E+01  0.453247E+01 -0.906494E+01  MHz
          0.151187E+01  0.151187E+01 -0.302374E+01  10-4 cm-1
  35Cl Vzz =   2.83503     a.u.,   NQCC =  -54.3896     MHz,  asym. parameter eta =  0.00000

Hg原子：

 ------------------------------------------------
  isotope             = 199Hg 
  nuclear spin I      =  0.5
  magnetic moment mu  =  0.5058855
          g_n = mu/I  =  1.0117710
  quadrupole moment Q =  0.00000 e 10-24 cm2
 ------------------------------------------------

 Atom    3        Hg Input number    2  xyznuc(Angstrom) =     0.0000    0.0000    0.0000

 ==== EFG-tensor (a.u.) Hg
            X             Y             Z
     X    0.268584E+01  0.000000E+00  0.000000E+00
     Y    0.000000E+00  0.268584E+01  0.000000E+00
     Z    0.000000E+00  0.000000E+00 -0.537169E+01

 ==== principal axes Q-tensor (EFG) Hg
           11            22            33
     X    0.000000      0.000000      1.000000
     Y    0.000000      1.000000      0.000000
     Z    1.000000      0.000000      0.000000
 principal values EFG Hg
          V11           V22           V33
         -0.537169E+01  0.268584E+01  0.268584E+01  a.u.
         -0.521986E+02  0.260993E+02  0.260993E+02  10^21 V/m^2
 principal values 199Hg  Q-tensor (includes factor Q/(2I(2I-1)))
          q11           q22           q33
          0.000000E+00  0.000000E+00  0.000000E+00  MHz
          0.000000E+00  0.000000E+00  0.000000E+00  10-4 cm-1
 199Hg Vzz =  -5.37169     a.u.,   NQCC =   0.00000     MHz,  asym. parameter eta =  0.00000

其中 Cl 的 Q-tensor 张量（注意，0.453247E+01 表示 0.453247×10$^{1}$ = 4.53247）：

• q$_{11}$=0.453247E+01 MHz
• q$_{22}$=0.453247E+01 MHz
• q$_{33}$=-0.906494E+01 MHz

上面还用另一种单位10$^{-4}$cm$^{-1}$列出：

• q$_{11}$=0.151187E01 × 10$^{-4}$cm$^{-1}$
• q$_{22}$=0.151187E01 × 10$^{-4}$cm$^{-1}$
• q$_{33}$=-0.302374E01 × 10$^{-4}$cm$^{-1}$

这二者之间的关系：Q张量元 = EFG张量元 × 234.9647 × 核四极矩quadrupole moment Q(以 barn为单位，1 barn = 10-24cm$^2$) ÷ [2 × nuclear spin I × (2 × nuclear spin I -1)]

Cl 的核四极耦合常数 NQCC = -54.3896 MHz，与 EFG的关系：NQCC = 以 a.u.为单位的 EFG 张量中最大的元(本例中是0.268584E+01 a.u.) × 234.9647 × 核四极矩quadrupole moment Q(以 barn为单位，1 barn = 10-24cm$^2$)

计算 $^{57}$Fe、$^{119}$Sn、$^{125}$Te、$^{193}$Ir和$^{197}$Au 的 Q-Tensor 时，quadrupole splittings 的输出中，单位改为 mm/s 用于 Mössbauer 谱。

如果要进行 EFG 主成分的 Mulliken 类型分析，以及根据正则 MO 的分析，参考：使用NBO、NLMO分析NMR化学位移、EFG贡献来源

1. 该数值对基组很敏感，因此需要使用尽量大的基组；
2. 原子核模型有点电荷模型和Gaussian型两种，对结果有影响，但一般不影响定性的大小关系，在Details > Relativity > Nuclear model里面可以设置；
3. 是否采用相对论方法Scalar或者Spin-Orbit也有影响，后者更精确；
4. Details > Relativity > Formallism > ZORA也可以修改为X2C，X2C效率和精度均高于ZORA，使用X2C时，建议使用TZ2P-J或QZ4P-J基组。

AMS软件提供免费试用，试用申请方式参见：AMS免费试用