特大地震在瞬时释放巨大能量，产生大面积地表破裂，同时会激发大量的地震波，并导致地球整体的振动，即地球自由振荡(简正模)。简正模的本征频率偏移和分裂是认识地球内部结构、物理参数分布以及主要圈层耦合机制的重要依据。在约束地球大尺度结构层面，如地幔侧向非均匀性、尺度因子变化、密度异常分布、地幔各项异性等方面，简正模有其得天独厚的优势。

地球自转和椭率导致地球自由振荡谱峰分裂形成多谱。根据起始模型，采用Rayleigh变分原理，通过理论模拟，可以非常准确地对其进行预测^[1-3]。因此，自转和椭率导致的谱峰分裂通常称为正常谱峰分裂。除自转和椭率外，地球内部结构横向不均匀性和各向异性等因素同样会导致自由振荡的谱峰分裂^[4-5]，称为异常谱峰分裂。Tanimoto^[6]结合扰动技术和反演方法，利用牛顿方法和变分原理求解非线性方程，分析横向不均匀性和非弹性地球模型的基本模式，研究地球自由振荡谱峰与横向不均匀性的关系。Woodhouse等^[7]通过最小二乘反演法从简正模观测波谱中首次获得谱峰分裂函数。随后，学者们相继给出对地幔、地核敏感的简正模的谱峰分裂函数^[4]，并用于反演、约束地幔横向不均匀性和内核各向异性。Ishii等^[8]联合自由振荡和自由空气重力异常对地球地幔的偶数阶横向不均匀性给出较好的约束。Koelemeijer等^[9]研究对D″区域敏感的简正模，并给出地球核幔边界的介质特性。总之，对不同深度介质敏感的简正模异常谱峰分裂反映了该地区的横向不均匀性、各向异性等介质分布特征，开启了一个新的探测地球内部结构特性的窗口^{[6, 10-15]}。

在利用简正模异常谱峰分裂函数约束地球内部介质参数分布时，首先需要扣除椭率、旋转导致的正常谱峰分裂效应。本文基于MINEOS软件、MATLAB的ode45工具箱和样条插值函数interp1，设计计算简正模正常谱峰分裂的接口。首先根据Dahlen等^[2]的方法计算简正模的一阶、二阶谱峰分裂参数及其单谱频率，采用静力水准重力势能J₂、J₄系数计算地球椭率一阶、二阶项，考虑PREM和1066A两种起始模型，得到₀S₂模态的计算结果，然后与Dahlen^[3]的计算结果进行比较。同时，计算谱峰分裂矩阵，并给出₀S₂与l≤5的球型环型简正模的耦合强度。

1 基本方法

对于任一球对称、非旋转、完全弹性和各项同性(SNREI)的起始模型，都可以得到简正模的本征频率_nω_l(简称ω₀)和本征位移(球型U、V，环型W)。椭率、旋转导致简并频率分裂2l+1单谱，这些单谱频率为_nω_l^m=_nω_l(1+a+mb+m²c)，-l≤m≤l。其中，_nω_l是SNREI简正模简并频率。

$ \begin{array}{l} a = \frac{1}{3}[1 - l(l + 1)\chi ]{\left( {\frac{\mathit{\Omega }}{{{\omega _0}}}} \right)^2} + \\ \;\;\;\;\;\;\frac{1}{2}\left( {v w_0^{ - 2} - \tau } \right) + {\alpha _2}{\left( {\frac{\mathit{\Omega }}{{{\omega _0}}}} \right)^2}\\ b = \chi \left( {\frac{\mathit{\Omega }}{{{\omega _0}}}} \right)\\ c = - \frac{3}{2}l(l + 1)\left[ {vw_0^{ - 2} - \tau } \right] + {\gamma _2}{\left( {\frac{\mathit{\Omega }}{{{\omega _0}}}} \right)^2} \end{array} $

(1)

采用式(2)计算椭率^[16]：

$ \begin{array}{l}{\varepsilon^{1 \mathrm{st}}=-3 f_{2} / 2} \\ {\varepsilon^{2 \mathrm{nd}}=-3 f_{2} / 2-5 f_{4} / 8-3 f_{2}^{2} / 8}\end{array} $

(2)

式中，f₂、f₄与二阶、四阶带状谐波系数J₂、J₄密切相关：

$ \begin{array}{l}{J_{2}=-f_{2}-\overline{m} / 3-11 f_{2}^{2} / 7-\overline{m} f_{2} / 7} \\ {J_{4}=-f_{4}-36 f_{2}^{2} / 35-6 \overline{m} f_{2} / 7}\end{array} $

(3)

值得注意的是，Chambat等^[17]更正了Nakiboglu^[16]关于f₂、f₄方程系数的错误(m=$\frac{{{\mathit{\Omega }^2}{R^3}}}{{GM}}$，G为地球引力常数，R为地球半径，M为地球质量)。式(1)中，第1式右边第1项是离心位的扰动球型部分1/3Ω²r²引起的；第2项是离心位的扰动非球型项Ω²r²P₂(cosθ)/3以及扁率的一级效应(旋转椭球引起的分裂与m²有关，故谱线分裂是l+1谱线分裂，旋转椭球结构引起的谱线分裂是不对称分裂)；第3项是科里奥利力的二阶扰动效应。第2式是科里奥利力的一阶扰动。第3式第1项是离心位和扁率的综合影响，第2项是科里奥利力的二阶扰动效应。如果仅考虑科里奥利力的一阶扰动和椭率的扰动效应，将会得到对角求和法则，即$\frac{1}{{2l + 1}}\sum\limits_{m = - l}^{m = l} {_n\omega _l^m} { = _n}{\omega _l}$。对角求和法则提供了由观测的单谱频率估算SNREI模型频率_nω_l的方法^{[7, 14-15, 18]}。

2 计算结果

首先，计算PREM模型椭率的一阶、二阶项(图 1)以及两者的差值。结果显示，椭率二阶项比一阶项大0.11%。其次，计算科里奥利力对球型简正模₀S₂单谱频率的二阶扰动效应(图 2)。根据Giabdini等^[19]的研究, 定义A=1/2(vωω₀^-2-τ)ω₀/2π，B=bω₀/2π，研究椭率、旋转对简正模_nS_l频率的贡献率(表 1)。地球自转角速度Ω为7.292 115×10⁵，对于低频振型(如₀S₂，角频率ω₀为1.944×10³)，自转效应比扁率大得多，自转效应Ω/ω₀随着振型频率的增加而减小，即简正模与自转的共振效应越来越小，且旋转造成的二阶效应会急剧减小，直至可以忽略；但对于超低频简正模尤其是和自转频率相近的简正模不可忽略。因此，对于高频简正模(如₁₃S₂)，其谱峰分裂的单谱将受椭率计算精度的影响。本文选取₀S₂和₁₃S₂模态，计算1066A和PREM模型一阶、二阶椭率的简正模谱峰分裂参数及其单谱频率，并与Dahlen的结果^[3]进行比较(表 3)。计算发现：1)内核环型振荡对简正模谱峰分裂参数影响极小，其计算结果为1.079 903 577 859 640×10^-5，运算中可以忽略这个影响。2)科里奥利力二阶扰动的主要贡献来自于自耦合，其次是相邻泛音阶球型振荡或环型振荡n=n±1，即在截断无穷耦合链只需考虑相邻阶次。根据Rogister^[20]的研究，位移场耦合链为：

$ s=\sigma_{l-2}^{m}+\tau_{l-2}^{m}+\sigma_{l}^{m}+\tau_{l+1}^{m}+\sigma_{l+2}^{m}, |m| \leqslant l-2 $

(4)

$ s = \tau _{l - 1}^{l - 1} + \sigma _l^{l - 1} + \tau _{l + 1}^{l - 1} + \sigma _{l + 2}^{l - 1}, \;|m| \le l - 1 $

(5)

$ s = \sigma _l^l + \tau _{l + 1}^l + \sigma _{l + 2}^l, \;|m| = l $

(6)

其中σ_l^m、τ_l^m分别为l阶m次的球型/环型位移场。3)PREM模型的谱峰分裂参数要略小于1066A模型。4)通过表 2可知，椭率二阶项得到的单谱频率差异微乎其微，实际计算中可忽略。通过表 3可知，椭率只对谱峰分裂参数a和c有较小影响，这种影响对于单谱频率影响基本可以忽略。

定义简正模耦合强度ε_kk′=‖〈k|δH|k′〉/(ω_k²-ω_k′²)‖，其中ω_k/ω_k′为任意2个简正模的简并频率，δH为谱峰分裂矩阵，左矢〈k和右矢k′〉为相应简正模的本征函数^[21]。得到₀S₂与l≤5的球型环型简正模的耦合强度(表 4)。结果表明，在椭率旋转造成的耦合中，₀S₂与₀T₃的耦合最强，其次是同角序数，不同泛音阶之间的同振型₁S₂耦合，这与选择法则(selection rule)高度吻合。

3 结语

本文计算了基于1066A和PREM模型正常谱峰分裂参数和单谱频率，首先从数值上证明了MATLAB接口的正确性，然后利用重力势能系数J₂、J₄求解椭率一阶项、二阶项，最后对椭率二阶项引起的简正模谱峰分裂参数和单谱频率影响展开讨论。结果表明，椭率二阶项对简正模谱峰分裂参数影响较小。由于谱峰分裂参数很小，所以椭率的二阶项对单谱频率的影响基本可以忽略，为下一步谱峰分裂矩阵可视化、谱峰分裂函数反演和自主开发基于MATLAB的简正模工具箱奠定了基础。