在不同的地幔储库中Li同位素的组成存在较大的差异，同时Li在不同矿物中的储存能力也不相同，锂元素及其同位素(⁶Li和⁷Li)因其特殊的地球化学性质，近二十年来作为地球化学示踪剂被广泛的用于研究地幔物质演化和俯冲带壳幔循环的相关问题中(Marschall et al., 2007; Ionov and Seitz, 2008; Magna et al., 2015)。由于Li的轻重同位素⁶Li和⁷Li拥有较大的相对质量差(~17%)(Seitz and Woodland, 2000; Woodland et al., 2002; Seitz et al., 2003, 2004; Teng et al., 2006; Trail et al., 2016)，而且具有不同的扩散速率，所以在高温条件下，会发生较大的与质量相关的同位素分馏效应(Gallagher and Elliott, 2009; Pogge von Strandmann et al., 2011)。斜顽辉石(Mg₂Si₂O₆)、镁铝石榴石(Mg₃Al₂(SiO₄)₃)以及镁铝尖晶石(MgAlO₂)被认为是地球上地幔圈层的重要组成矿物，从成因上，三者主要是岩浆作用的产物，是基性、超基性岩的主要造岩矿物，常与橄榄石、单斜辉石等矿物共生(Magna et al., 2006; Jeffcoate et al., 2007)。在上地幔和转换带的化学元素组成中Li元素是不可分割的一部分，研究发现Li可以以取代矿物晶格中Mg空位的方式进入辉石、石榴石和尖晶石等地幔矿物中(Lodders, 2003; Grant and Wood, 2010)，并可通过填隙位与取代空缺的机制在其中迁移(Richter et al., 2014)。但是对于Li在原子水平上是如何在地幔矿物中迁移的，在上地幔共生矿物相中Li同位素组成存在差异的原因，这些问题的研究都存在着争议。本文通过经典力场经验势的方法计算模拟了Li同位素在上地幔主要矿物，斜顽辉石、镁铝榴石和镁铝尖晶石中原子水平上的扩散机制和矿物内Li同位素的分馏情况，从而探讨Li同位素在斜顽辉石、镁铝榴石和镁铝尖晶石中由于扩散导致的同位素组成特征以及Li同位素在这些上地幔主要矿物中的分布等相关问题。

1 计算方法

在Li同位素迁移扩散以及分馏的计算模拟中, 我们采用了由Gale J. D.教授提供的General Utility Lattice Program (GULP)软件包(Gale, 1997)。由于矿物晶体的周期性边界条件，为了避免原胞内缺陷-缺陷间的相互作，对于所有的Li同位素迁移计算和Li同位素分馏计算，我们均使用了：斜顽辉石使用了2×2×4(19.08Å×17.28Å×21.25Å)，含有640个原子；镁铝石榴石，2×2×2(22.65Å×22.65Å×22.65Å)，含有1280个原子；镁铝尖晶石，2×2×2(16.16Å×16.16Å×16.16Å)，含有448个原子的Supercell来进行模拟计算，其中对于Li同位素的分馏计算时，Supercell中的Li/Mg值为1/128，在建立好的Supercell模型的基础上我们模拟计算了Li在各矿物中不同位置的点缺陷，包括取代空位和填隙位的点缺陷，通过计算不同点缺陷体系的能量来确定其稳定性。

本文中完整晶体和缺陷晶体的能量都是在Born近似条件下计算得到的。在这种近似下，用下式描述相距为r的两离子间的相互作用：

(1)

(1) 式中的第一项为长程库仑相互作用，后两项是用Buckingham势描述的短程原子间相互作用能。其中Z_i、Z_j分别为原子i和j的核电荷数，A_ij、ρ_ij和C_ij是需要程序拟合的势参数，具有不唯一性和可转移性。由于晶格中存在的缺陷会极化周围的其他离子，因此考虑了离子极化率α。离子极化率采用壳模型(shell model)来处理：把原子或离子分成核(core)和壳层(shell)两部分，其中，核具有质量和正电荷Y，壳层则没有质量，但有负电荷Z-Y(Z为表观电荷)。核和壳通过一个弹性系数K相关联。通过壳模型处理方法，由于电子极化而引起的壳层电子相对于核的位移就能通过以下公式得到很好的处理：

(2)

本文选取了Catlow (1977)、Lewis and Catlow (1985)对MgO、Al₂O₃、SiO₂和Li₂O拟合的Buckingham短程势参数来描述O原子与各阳离子间的非共价键相互作用，并用三体项(Harmonic three-body term)来描述Si-O四面体中O-Si-O的三体相互作用(Sanders et al., 1984)，可表示为：1/2k (θ-θ₀)，其中θ代表O-Si-O的夹角。这组势参数对氧化物的模拟已十分成熟，并成功用来模拟存在多种阳离子的材料的相关性质(Palin et al., 2003, 2005; Zhang and Wright, 2012a)。所选用势参数的详细相关参数值如表 1所示。

我们对短程互作用势的最大作用范围(cutoff)分别选取10Å和20Å进行收敛性测试，得到相应的Supercell的晶格能的相对误差小于10^-5，因此可以确认选取的势作用截断半径cutoff为10Å时已经达到收敛。另外所建立的Supercell最小边长为16.16Å，大于相互作用势的截断半径(cutoff)10Å；故认为所模拟计算的Supercell已经足够大，有效避免了晶格中可能存在的缺陷-缺陷间的相互作用。同时，为了考察这组经验势参数是否能够较好的模拟地幔矿物的晶体结构，我们对计算得到的晶格常数与实验值进行了比较，二者之间的误差小于3%。

2 计算结果 2.1 斜顽辉石，镁铝榴石和镁铝尖晶石中Li的点缺陷模拟计算

斜顽辉石(Mg₂Si₂O₆)晶格中共有两种Mg的位置，Mg1和Mg2，所以Li以空缺取代机制存在的可能性有两种——Li′_Mg1和Li′_Mg2(用Kroger-Vink缺陷表示法表示)，而Li在填隙位的点缺陷可表示为：Li_i^·。在Mg1和Mg2位的空缺表示为：V″_Mg1和V″_Mg2。计算2×2×4的斜顽辉石Supercell中各种Li的点缺陷的能量得到表 2。我们发现Li取代Mg2位上的Mg所需要的能量比取代Mg1位上的Mg低0.78eV，证明Li位于Mg2位上时比位于Mg1位更为稳定。

镁铝榴石(Mg₃Al₂(SiO₄)₃)和镁铝尖晶石(MgAlO₂)的晶格中只有一种镁，Mg1，所以Li以空位取代机制存在的可能性只有一种—Li′_Mg1，Li在填隙位的点缺陷可表示为：Li_i^·。Li在Mg1位的空缺表示为：V″_Mg1。我们计算得到了镁铝榴石和镁铝尖晶石Supercell中各种Li的点缺陷的能量(表 2)。

在所有的计算中，通过经典力场方法计算得到的能量均为矿物体系的吉布斯自由能G：

其中G为吉布斯自由能，H为体系的焓，U为体系的内能，T为开尔文温度，P、V分别为压力和体积。

Li在地幔矿物中可以在填隙位和Mg空位存在，两者均需要进行电荷平衡。在纯净的矿物晶格中Li的点缺陷可能的电荷平衡方案有两种

① 2Li_i^·+V″_Mg，两个填隙位Li+一个Mg空位；

② Li_i^·+Li′_Mg，一个填隙位Li+一个取代Mg位的Li

为比较这两种方案哪种更为稳定，针对斜顽辉石、镁铝榴石和镁铝尖晶石体系我们构造下面三个反应：

(ⅰ)

(ⅱ)

(ⅲ)

其中(ⅰ)、(ⅱ)和(ⅲ)分别代表斜顽辉石、镁铝榴石和镁铝尖晶石体系，方程左边代表方案①，方程右边代表方案②。对于斜顽辉石，考虑到Li在Mg2位更加稳定，所以Li取代Mg空位的点缺陷能量值选用的Li′_Mg2数据。根据表 2中得到的能量值，我们发现，在三种矿物体系中，都是反应式右边(即方案②)的能量值相对较低，证明在不考虑除Li以外其他杂质元素的影响时，方案②为Li在斜顽辉石、镁铝榴石和镁铝尖晶石中最可能的电价平衡方式。这与Zhang and Wright (2012b)计算的Li在镁橄榄石中的电价平衡方式是一致的，同时也与Dohmen and Milke (2010)及Grant and Wood (2010)在实验中的假设条件相符。在本次研究中，我们仅考虑了Li在纯净矿物中的存在机制，但是在自然界的地幔矿物中还有很多其他离子元素及微量元素，因此有许多不同的电荷平衡方式。前人工作表明三价阳离子，如Al³⁺、Fe³⁺、Sc³⁺及REE(Purton et al., 1997; Wood and Blundy, 2001; Grant and Wood, 2010; Zhang and Wright, 2010, 2012a)，或者单价阳离子，如Li⁺及H⁺(Grant and Wood, 2010; Zhang and Wright, 2010)等都可参与地幔矿物中的电荷平衡过程。在自然界的斜顽辉石晶格中，由于多种阳离子(Fe³⁺，Al³⁺等)的存在，会有多种电价平衡的方式，这些离子的存在能够平衡更多的Li离子，从而增加辉石对Li的吸收(Grant and Wood, 2010)。

2.2 锂离子在斜顽辉石、镁铝榴石及尖晶石中的微观扩散机制

扩散是引起地幔物质中同位素和化学组成不均一性的一个重要原因，讨论离子在矿物中的扩散有助于了解元素分配的原因和离子扩散的微观机制等方面的性质。因为Li会以点缺陷的形式在矿物晶格中的填隙位和取代空缺位存在，所以对于斜顽辉石，我们考虑了Li在斜顽辉石晶格中Mg空位间的迁移路径；对于镁铝榴石和尖晶石，我们考虑了Li在矿物晶格中Mg空位以及填隙位间的迁移路径。并计算出Li沿这些可能的路径迁移时的活化能，以此来衡量沿着该路径Li的迁移能力的大小，从而确定在能量趋向上，Li更倾向以哪条路径进行迁移。

我们用Li取代Mg的位置代表取代空位机制的点缺陷，将Li放置在晶格间隙中代表填隙位机制的点缺陷。对于填隙位点缺陷的模拟，我们将Li放在多个填隙位位置上分别进行晶格结构优化，计算得到相应的能量值，以能量最小值的方法确定Li在晶格中最稳定的填隙位位置，并将其作为Li在矿物晶格中以填隙位点缺陷存在的形式。我们所模拟的Li在矿物中沿扩散路径的迁移是从一个稳定的低能量位置迁移若干步后到另一个等价的低能量位置，这里的低能量位置可以是填隙位，也可以是空缺取代位。Li沿迁移路径迁移的每一步，我们都固定Li原子平行于迁移路径的空间方向，同时对所有原子进行弛豫。我们将Li迁移的每一步的能量作图得到Li沿该路径迁移时能量的变化，并将最大和最小能量值的差值作为Li沿该路径迁移时所需的迁移活化能。

对于斜顽辉石，Li以空位机制迁移的路径有3条：A(Mg1-Mg2)，B(Mg1-Mg1)，C(Mg2-Mg2)。这三条路径大致沿c轴方向分布，如图 1。计算发现，Li以空缺取代机制迁移，当沿路径A(Mg1-Mg2)迁移时，可得到最低的迁移活化能(1.05eV)。Li沿这三条路径迁移的活化能值及晶格能变化如表 3，图 2所示。通过计算我们发现Li在斜顽辉石中的迁移具有方向性，均沿c轴方向。对于空位机制，路径A(Mg1-Mg2)的迁移路径具有最低的迁移活化能，路径B(Mg1-Mg1)路径次之，而路径C(Mg2-Mg2)的迁移活化能则最高(图 2)，因此Li难以沿该路径进行迁移。

对于镁铝榴石及镁铝尖晶石，由于两种矿物中的Mg位仅有一种，再加上镁铝尖晶石结构对称性较高，Li在其中的迁移路径较为简单。我们模拟计算了镁铝榴石和镁铝尖晶石中的可能的空位机制(A)和填隙机制(B)迁移路径(如图 3所示)。Li在各矿物中沿不同路径迁移的活化能值(表 4)表明，在镁铝尖晶石中，空位机制(A)和填隙机制(B)的迁移活化能分别为0.493eV和1.445eV，后者迁移方向沿矿物几何对角线方向[111]方向。而在镁铝石榴石中空位机制迁移路径(A)迁移活化能为0.896eV，迁移方向为沿着[100]方向，填隙机制迁移路径(B)迁移活化能为1.655eV，迁移方向为沿着[010]方向。Li在镁铝榴石和镁铝尖晶石中沿各路径迁移活化能的变化如图 4所示。

我们发现Li在镁铝尖晶石中的最稳定的填隙位位于两个Mg位之间，而在镁铝榴石中则位于三个Mg位所组成的等边三角形中心。且以空位机制进行迁移时具有较低的迁移活化能。

2.3 锂同位素在斜顽辉石和镁铝榴石以及镁铝尖晶石中不同晶格位上的分馏作用

Li同位素在斜顽辉石和镁铝榴石矿物晶格内填隙位和空缺取代位间的分馏系数将用公式(3)进行计算：

(3)

其中，α_{interstitial-vacancy}代表Li同位素在不同晶格位间的分馏系数，β代表与质量相关的同位素配分函数比，可表示为：

(4)

(5)

那么，

(6)

公式(6)中，G代表Li的两个同位素于特定温度压力条件下，在不同晶格位上存在时的吉布斯自由能，Q代表两个同位素在某一种晶格位上的配分函数。kB代表玻尔兹曼常数，T代表开尔文温度。通过以上公式，我们就可以通过得到相应的吉布斯自由能的值来计算Li同位素的分馏系数α_{interstitial-vacancy}。

对于斜顽辉石，我们选取300K、1000K、1500K三个温度和4GPa、10GPa两个压力条件进行了相应的分馏系数的计算，计算结果见表 5。我们的结果发现，在10GPa条件下，当温度从300K升高到1500K，分馏系数从1.0882 (88.2‰)迅速降低到1.0040(4.0‰)，该趋势在4GPa压力条件下也符合，这表明温度是Li同位素在斜顽辉石晶格中不同晶格位上的分馏作用的一个重要的影响因素。而且比较同一温度不同压力的分馏程度我们发现，例如，在300K条件下，斜顽辉石中Li同位素在4GPa和10GPa两个压力下的分馏程度相差达到了~39‰，这表明Li同位素在高压环境下的斜顽辉石中分馏程度更高。而且从表 5的数据我们同样可以看出，Li的重同位素更倾向于以填隙位的方式存在于斜顽辉石晶格中，相反的轻同位素则更易以Mg空位的形式存在。

对于镁铝榴石，我们选取1000K、1400K、1700K三个温度点和4GPa、10GPa、13GPa三个压力点，并在这些条件下计算了Li同位素在镁铝榴石中相应的分馏系数值(表 6)，通过计算发现，当压力为4GPa，温度从1000K升高至1700K时，Li同位素的分馏系数由1.0087(8.7‰)降低为1.0028(2.8‰)，分馏程度随着温度的升高而减弱，减弱程度与斜顽辉石相比较小，这种趋势在压力为10GPa和13GPa时同样存在。比较同一温度下不同压力对于Li同位素分馏程度的影响，我们发现，压力作用对Li同位素的分馏影响较为明显，以1000K为例，当压力从4GPa升高到13GPa时，分馏程度从8.7‰升高到9.7‰(表 6)，这说明在镁铝榴石中，压力升高使得Li同位素的分馏程度升高。同时我们还发现，与斜顽辉石类似，由于分馏系数均大于1，可以得到⁷Li更倾向进入镁铝石榴石中的填隙位，⁶Li则更容易进入Mg空位这个结论。

对于镁铝尖晶石，我们选择了300K、1000K、1500K三个温度点在0GPa、4GPa、10GPa的压力条件下计算了Li同位素在镁铝尖晶石中相关的分馏情况，见表 7。

我们的计算结果表明，在压力为0GPa时，当温度从300K升高到1500K，Li同位素的分馏系数从1.0358(35.8‰)迅速降低至1.0019(1.9‰)，随着温度的升高，分馏系数减小，分馏程度减弱，在4GPa和10GPa的压力条件下，同样也出现了这样的趋势，这种趋势和斜顽辉石以及镁铝榴石中分馏系数变化趋势相同。在镁铝尖晶石中，与温度对分馏程度的影响相比，压力的影响较弱，例如在温度同为300K时，压力从0GPa升高到10GPa，Li同位素的分馏系数仅从1.0358(35.8‰)降低到了1.0327(32.7‰)(表 7)。此外我们还发现，在镁铝尖晶石中，Li同位素的分馏系数均大于1(表 7)，这说明⁷Li更倾向于以填隙位的方式存在于镁铝尖晶石中，而⁶Li则倾向于以Mg空位的方式存在于镁铝尖晶石。

3 讨论与结论

斜顽辉石中Li同位素的分馏效应不仅对温度敏感，对压力也十分敏感，从计算数据中可看出，随着温度的增加，Li在辉石晶格中的分馏程度降低，温度对Li同位素分馏程度的这种影响，可能是由于温度升高后粒子运动速率加快，化学不均一性减弱所造成。当体系的压力升高时，与体系温度升高对矿物晶格具有相反的影响，故当压力增加时，Li同位素在辉石晶格中的分馏程度与温度效应相反，即随着压力升高，Li同位素的分馏系数增加。这与McDade et al. (2003)讨论的压力对Li元素的分配及同位素的分馏效应的影响相一致。在低温条件下Li的重同位素⁷Li在斜顽辉石晶格中会被固定在填隙位上，从而使得⁶Li更易扩散出矿物。Xu et al. (2013)发现在辉石晶体的边缘δ⁷Li值呈极负值，这可能是由于⁶Li的迁出矿物而⁷Li则滞留在矿物内，使得观察到矿物边缘比矿物内部具有更低的⁷Li值。

对比同一压力条件下温度对Li同位素分馏效应的影响程度，我们发现：在高压条件下，与镁铝石榴石相比，斜顽辉石的这种影响程度更大。如压力为10GPa时，对于镁铝石榴石，温度由1000K升高到1400K，分馏效应由9.4‰减小到4.7‰，变化了4.7‰，而对于斜顽辉石，温度由1000K升高到1500K：分馏效应由9.2‰减小到4.0‰，变化了5.2‰；而在低压条件下，镁铝石榴石的这种影响则变得较大，斜顽辉石相对较小，如压力为4GPa时，对于镁铝榴石，温度由1000K升高到1400K，分馏效应由8.7‰减小到4.3‰，变化了4.4‰，对于斜顽辉石，温度由1000K升高到1500K：分馏效应由5.2‰减小到2.3‰，变化了2.9‰。

同样的，当温度一致时，压力升高对Li同位素在地幔矿物中的分馏效应的影响程度表现为：斜顽辉石中最强，镁铝石榴石次之。如1000K温度下，压力从4GPa升至10GPa，斜顽辉石中分馏效应由5.2‰升高到9.2‰，变化了4.0‰。而在镁铝石榴石中，Li同位素的分馏效应由8.7‰升高到9.4‰，变化了0.7‰。与斜顽辉石和镁铝榴石不同，在镁铝尖晶石中，压力与Li同位素的分馏程度为反相关关系，但是压力对于尖晶石中Li同位素分馏的影响非常小。1000K时，压力从4GPa升高到10GPa，Li同位素的分馏效应由4.2‰减小到4.0‰，变化了0.2‰，1500K时，压力由4GPa升高到10GPa，同位素的分馏效应由1.91‰减小到1.87‰，变化了0.04‰

在相同的温度压力条件下，斜顽辉石、镁铝石榴石以及镁铝尖晶石中Li同位素的分馏效应存在明显的差异，镁铝石榴子石具有最强的分馏效应，斜顽辉石次之，镁铝尖晶石则最弱，斜顽辉石中的Li同位素分馏效应对压力很敏感，若在高压条件下，分馏程度可与镁铝石榴石近似，而在低压条件时，分馏程度则小于镁铝榴石，反观镁铝尖晶石，Li同位素的分馏对于压力的响应并不明显，而温度成为影响尖晶石中Li同位素分馏的重要因素。

对于同一温度下，压力升高对不同矿物中Li同位素分馏效应的影响程度是不同的，在斜顽辉石中最强，镁铝石榴石次之，即压力升高，斜顽辉石中的Li同位素分馏效应升高的幅度最大，镁铝石榴石中则较低。但是，镁铝尖晶石中则出现相反的趋势，压力的升高，使得尖晶石中Li同位素的分馏程度减弱，我们认为，这是由于镁铝尖晶石，镁铝榴石，斜顽辉石不同矿物之间不同的化学成分所造成。对于同一压力下，温度升高对不同矿物中Li同位素分馏效应的影响程度表现为：高压下，斜顽辉石中Li同位素分馏效应降低的幅度较大，镁铝尖晶石中次之，镁铝石榴石中的Li同位素分馏效应降低的幅度最小；低压下同样也是斜顽辉石减小程度最大，镁铝尖晶石次之，镁铝榴石最小。

我们的计算发现Li在斜顽辉石中取代空位机制的路径进行迁移时，更倾向于以路径A(Mg1-Mg2)进行迁移，迁移方向为c轴方向。在镁铝榴石和镁铝尖晶石中与填隙位机制相比，Li更倾向于以空位机制来进行迁移。Li以取代空位机制迁移时一定程度上与Li-Mg交换速度，Mg空位浓度有关。Li同位素在不同晶格位上的分馏作用计算表明，⁷Li会更多的存在于晶格的填隙位上，⁶Li则会相对更多的进入Mg空位中。Li同位素这样的分馏行为，是共生矿物相间存在Li同位素的组成差异的微观机理。

在低温条件下，⁷Li会被固定在填隙位上，当温度升高时，其移动性逐渐变强，而⁶Li的迁移速率一定程度上会受控于Mg空位的浓度。我们认为这种Li在地幔矿物晶体内部不同晶格位上分馏的行为是绝大多数幔源物质中Li同位素分馏的主要原因，可用于解释浅部幔源岩石中低δ⁷Li的特征(Nishio et al., 2004)及冷却驱动的同位素分馏作用模型等科学问题。

致谢      感谢Julian Gale教授提供最新版本的GULP软件包。计算资源由中国科学院地质与地球物理所计算机模拟实验室(Computer Simulation Lab, IGGCAS)及国家超算深圳中心提供。感谢二位匿名审稿人对本文提出的宝贵修改意见。