水镁石与方镁石的热力学及热导特征的第一性原理研究

於造林， 宋亚兰， 马超杰， 何开华*， 姬广富

（1.中国地质大学 数学与物理学院，湖北 武汉430074；3.中国工程物理研究院流体物理研究所冲击波物理与爆轰物理国防科技重点实验室，四川绵阳621900）

水镁石由富镁的硅酸盐以及富镁碳酸盐矿物经低温热液蚀变而成，它的结构和化学成份都比较简单，尽管在地幔岩石中含量较少，但它作为含水层状矿物基本代表类型，可以为复杂的含水矿物研究提供依据.水镁石分子式Mg（OH）2，结构分层简单，其结构为三方晶系，空间群为，由镁配位组成的八面体，八面体里面含有2个羟基，每层八面体间间距0.48 nm，羟基键0.1 nm，并垂直于八面体层［1-6］.近年来，对于水镁石的研究越来越广泛，不管是在实验方面还是理论计算方面都取得了重大进展.文献［3］通过拉曼光谱的研究给出了水镁石高温高压作用下的结构相变；文献［2］通过布里渊区散射的方法得到了单晶水镁石在15 GPa时的弹性常数以及其它弹性性质；文献［1］通过实验方法得到了高压下水镁石的热膨胀系数以及其压强与熵的变化关系等.文献［7］通过中子衍射的方法研究了同位素氢和氘对Mg（OH）2性质的影响.含结构水（OH）对矿物物理性质的影响很大，近年来对此问题的研究吸引了人们的兴趣，比如会引起地震波、电导以及弹性性质较大的变化等［8-10］.Mg（OH）2中水的质量分数达31.03%，在高温高压的作用下会脱水形成方镁石，这一过程必将影响矿物的物理性质.揭示脱水前后矿物的物理性质有利于理解地球内部结构，解释地球物理观测的现象.关于2种矿物的热力学性质以及热传导性质之间的异同目前还未见报道.本文采用第一性原理计算，分别对水镁石和方镁石的热力学性质以及热传导性质进行了研究，通过对比分析发现2种矿物的热力学性质及热导性质有着很大的差异.

1 研究方法

本文中结构优化和力常数的计算（二阶和三阶）采用基于密度泛函理论的第一性原理计算，该方法在众多研究领域中被证实是行之有效的［11-12］.本文采用计算精度和效率较好的VASP软件［13］.粒子间的相互作用采用缀加投影平面波来描述［14］，交换关联能选取局域密度近似处理.布里渊区积分采用Monkhorst-Pack方法［15］，K 点网格大小为10×10×10.截断能选取500 eV，能量收敛精度设定为10-6eV，力的收敛精度为0.5×10-3eV／nm.优化过程中考虑晶格常数和原子位置均充分弛豫.

2种矿物声子谱及热力学性质的计算采用基于准谐近似的 PHONOPY软件［16］，此软件采用有限差分法并且基于原子力可以算出二阶力常数.求解高对称点的声子频率的本征值可以构建动力学矩阵，而其他点的声子模式和频率通过对动力学矩阵在实空间中进行傅里叶变换可以得到.对于原子力的计算，需要在2种矿物的超胞（水镁石与方镁石均为3×3×3超胞）中建多个独立的原子位移结构，每个结构需要利用VASP计算单点能进而得到原子相互作用力.

三阶力常数的计算采用基于实空间有限差分法的ShengBTE软件［17］，晶格热导的计算通过求解线性玻尔兹曼传输方程.晶格热导的表达式

CV、vg、τ分别为等体热容、声子群速度、声子寿命.其中声子的寿命与散射率成反比，因此可以用散射率的大小来反映声子寿命的长短.本文三阶力常数的计算选取了与二阶力常数相同的超胞大小，并且考虑到了四次近邻相互作用.为得到可靠的晶格热导率计算中选取了较大的q网格点（41×41×41）.

2 结果与分析

2.1 热力学性质 水镁石在高温高压条件下脱水为方镁石，为了更好的理解水对矿物性质，以及压力和温度等外界环境对其热力学性质的影响，计算给出了2种矿物在不同条件下的热膨胀系数（见图1）.方镁石的热膨胀系数在压强一定时会随温度增大而增大，低温时增长较快，然后逐渐变缓.当温度给定时，热膨胀系数会随压力的增大而减小.在压力为0 GPa，温度为300 K时，MgO的热膨胀系数为3.26×10-5K-1，文献［18］通过实验得到 MgO 的热膨胀系数为3.15×10-5K-1，文献［19］通过实验得到MgO的热膨胀系数为3.16×10-5K-1，文献［20］通过计算得到其热膨胀系数为3.11×10-5K-1.本文计算得到的MgO的热膨胀系数与前面的实验和计算值都符合得很好，也保证了此次的计算的可靠性.水镁石的热膨胀系数在低温时随温度增长很快，在大概600 K时达到最大值，但是温度继续升高后，其热膨胀系数反而减小，这个趋势和方镁石有区别.另外还有一个值得注意的是，当给定温度时，Mg（OH）2的热膨胀系数随压力的增大而增加，这与MgO的变化趋势相反.因此以上热膨胀的结果揭示了2种矿物间的差异，体现了结构水对物理性质的调控作用.为了研究清楚2种矿物热膨胀系数的区别，近一步对影响热膨胀的其他几个参数进行了探讨.

图1 水镁石与方镁石热膨胀系数随压强和温度的变化关系Fig.1 Thermal expansivity of brucite and periclase

对比图2与图3，可以发现2种矿物的CV、k0及V的变化趋势基本相同.CV随着温度升高逐渐增大，在高温时（800 K）与温度的依赖关系减弱，最终会达到饱和，符合杜隆-伯特定律.在压力的作用下，2种矿物的热容会略微减小.水镁石脱水的压力不是很高，在本文考虑的压力范围内，热容变化不大.在压力作用下水镁石的体变模量会从147 GPa增加到213 GPa，增幅为45%；而方镁石的体变模量会从170 GPa增加到220 GPa，增幅为29%.对比2种矿物的体积来看，他们在压力和温度作用下的变化非常类似.现在来关注格林艾森常数，从图2中可以看到水镁石的格林艾森常数随着压力的增加而增加.在低温时，水镁石的格林艾森常数随温度增加而增大，在较高温度时与温度依赖关系减弱，进一步升高温度，发现此参数会出现略微减小.需要注意的是在低温端（～40 K）水镁石的格林艾森常数出现了一个峰，方镁石也有类似的特点，存在尖峰的问题应该是少数低频率的声子能提供较大的格林艾森常数（微观格林艾森常数）.因给定温度下的格林艾森常数可表示为［21］：

其中，ωi为第i个声子频率，其对应的格林艾森常数为 γ（ωi）.不同频率的声子对应的 γ（ωi）大小不同.在压力的作用下，2种矿物的结构有较大的变化（特别是水镁石属三方晶系），可能出现一些低频声子所对应的γ（ωi）会比较大，进而导致在低温时出现尖峰.方镁石的格林艾森常数在压力的作用下会逐渐减小，这个趋势恰好和与水镁石相反.另外在温度的作用下，特别是高温部分，方镁石的格林艾森常数会达到饱和，没有呈现出减小的趋势，这也和水镁石有差别.

图2 水镁石的热力学参数，C V，k 0，V以及γ随温度与压力的变化关系Fig.2 Thermodynamic parameters（CV，k 0，V，γ）of brucite

因此从上述的4个参数分析来看，格林艾森常数与体变模量是主要决定热膨胀系数的物理参数，可以用来解释2种矿物的热膨胀的不同特征.在低温部分，水镁石的格林艾森常数随压力增大而增大，与方镁石的变化恰好相反，因此决定了水镁石的热膨胀系数随压力增加而增加，方镁石随压力增加而减小.需要说明的是在高温部分，水镁石的格林艾森常数会随着温度的降低而出现了小幅小降，并且相对方镁石而言，水镁石的体变模量变化较大，2个参数共同决定了水镁石在高温时热膨胀系数减小的趋势.另外水镁石的格林艾森常数低温端的峰（12 GPa，～40 K）在热膨胀系数中并没有体现出来，是因为在此温度附近的定体热容很小接近于0，此时热膨胀系数主要由定体热容决定，而格林艾森常数对热膨胀系数的影响有限，进而导致热膨胀系数在低温时非常小甚至等于0（见热膨胀系数图3）.随着温度升高，格林艾森常数对热膨胀系数的影响增强，其变化趋势与热膨胀变化趋势相似.

2.2 热导性质 图4分别给出了方镁石与水镁石的晶格热导率.可以看出在环境条件下，方镁石的热导率为80 W／mK，而水镁石的热导率仅为0.97 W／mK，因此可以推断含水对矿物的热导率产生了极大的影响.随着温度的升高，2种矿物的热导率都会降低，符合热导率与1／T近似正比的关系.在压力的作用下，方镁石的热导率会增加，而水镁石却表现出了相反的特征，其热导率会随着压力的增加而减小.2种矿物热导率的相异特征与前述的热膨胀特征相似.

图3 方镁石热力学参数，CV，k 0，V以及γ随温度与压力的变化关系.Fig.3 Thermodynamic parameters （CV，k 0，V，γ）of periclase.

图4 水镁石与方镁石的晶格热导率与温度压力的关系Fig.4 Lattice thermal conductivity of brucite and periclase

图5 水镁石与方镁石的非谐散射率Fig.5 Anharmonic scattering of brucite and periclase

3 结束语

1）2种矿物的热膨胀系数随温度压力的变化存在较大的差异.低温时，2种矿物的热膨胀系数都会随着温度的升高而升高，但水镁石在较高温度时会呈现出微弱减小的趋势.另在压力的作用下，水镁石的热膨胀系数会增大，而方镁石的热膨胀系数会减小.2种矿物的格林艾森常数和体变模量的变化决定了他们热膨胀系数的变化.

2）通过求解玻尔兹曼方程得到了2种矿物的晶格热导率，计算显示在压力作用下的2种矿物的晶格热导率的变化与热膨胀系数的变化相似，但是两者间的幅度相差很大，方镁石的晶格热导率比水镁石的要高2个数量级，说明含水或脱水对热导性质的影响很大.脱水后形成的方镁石具有较小的非谐散射率，因此声子寿命增大，导致方镁石的热导率增大.