该研究由哈佛大学地球与行星科学系的Dimitri Komatitsch和Jeroen Tromp合作完成，发表于《Geophysical Journal International》期刊，1999年7月接受发表（正式发表于1999年第139卷第3期）。
这项研究属于计算地震学领域，旨在解决三维地球模型中合成地震图计算的数值方法问题。传统方法如有限差分法（finite difference method）在存在地形起伏和各向异性介质时面临显著困难，而伪谱方法（pseudospectral method）则受限于模型的平滑变化要求。谱元法（spectral element method，SEM）结合了有限元法的灵活性和谱方法的高精度，为解决这些问题提供了创新方案。
研究基于运动方程的弱形式（weak formulation）求解，采用适应自由表面和主要内部不连续面的六面体单元网格（hexahedral elements mesh）。单元上的波场离散采用高阶拉格朗日插值（high-degree Lagrange interpolants），而单元积分则基于Gauss-Lobatto-Legendre积分法则。这种方法结合离散化和积分技术，产生了对角质量矩阵（diagonal mass matrix），大大简化了算法。
研究首先建立了弹性或非弹性固体中地震波传播的运动方程，可以采用强形式（strong form）或弱形式（weak form）求解。在强形式中直接处理微分形式的运动方程和边界条件；而在弱形式中使用运动方程的积分形式，类似于有限元法（finite element method，FEM）和直接解法。
动量方程表达为：ρ∂²s/∂t²=∇·T+f。应力张量T通过胡克定律与位移梯度∇s线性相关：T=C:∇s（对于弹性各向异性固体）。对于衰减介质，胡克定律需修改为考虑整个应变历史的形式。
自由表面边界条件为：T·n̂=0。在人工模型边界Γ上，使用近似吸收边界条件将牵引力与速度相关。
通过对动量方程进行点积测试向量w，在模型体积Ω上积分并施加边界条件，得到弱形式方程。这种方法自然地满足自由表面边界条件，不需要显式施加，因此能更准确地模拟表面波。
模型体积被细分为不重叠的六面体单元Ωₑ，吸收边界Γ由四边形表面单元Γ_b表示。每个单元通过控制点和形函数（shape functions）定义，形函数是1次或2次拉格朗日多项式的乘积。
在边界和体积单元上，函数通过拉格朗日多项式插值表示。在边界单元Γ_b上，函数f插值为拉格朗日多项式乘积的和；在体积单元Ωₑ上，函数f表示为三重拉格朗日多项式乘积的和。
全局系统的位移向量u服从常微分方程：MÜ+CṺ+KU=F。谱元法的一个显著特点是质量矩阵M在构造上是对角的，这大大降低了算法的复杂度。研究使用显式二阶有限差分（finite difference，FD）方案进行时间离散，将刚度和吸收项移到右侧。
对于衰减介质，采用吸收带本构关系（absorption-band constitutive relation）确定应力张量，并单独处理记忆变量方程。这些一阶时间方程使用已知高效的Runge-Kutta四阶方案积分。
研究首先测试了一个覆盖半空间的单层模型，比较了谱元法与离散波数/反射率方法（discrete wavenumber/reflectivity method）的结果。源为垂直力，位于半空间中块体中间25.05km深度处。结果表明，对于直接P波和S波以及众多多次波，两种方法吻合良好。
在相同模型中，将源置于靠近地表（536.1m深度）时，响应包含显著的表面波贡献。谱元法解与离散波数/反射率参考解再次显示出极好的一致性，验证了该方法对表面波的准确表示能力。
第二个模拟测试了三层覆盖半空间的更复杂模型。研究表明，谱元法能精确模拟垂直入射平面P波或S波从模型底部来的情况下地表记录的位移。与高精度一维有限差分参考解相比，谱元法结果显示出极好的一致性。
在均匀半空间中测试了纯倾滑源（dip-slip source）情况，比较了谱元法、频率-波数（frequency-wavenumber，FK）法和有限差分法的结果。三种技术获得的速度分量总体上吻合良好，验证了谱元法中矩张量源（moment tensor source）的实现。
通过半球形陨石坑（hemispherical crater）嵌入均匀半空间的模型，验证了自由表面边界条件在非常陡峭地形情况下的准确性。比较了谱元法与近似边界法的结果，两种方法显示出极好的一致性，特别是陨石坑边缘附近的强烈放大效应被很好地再现。
在二维均匀强衰减介质（Qₚ≈30和Qₛ≈20）中，比较了谱元法与Carcione等（1988）推导的点源力解析解。结果显示很好的一致性，验证了在谱元法中引入衰减的方法。特别展示了与相同松弛材料性质的弹性介质解相比，S波振幅减少了超过两倍的强烈衰减效应。
该研究详细介绍了谱元法在地震波传播领域的应用，展示了该方法在三维地球模型中计算合成地震图的强大能力。谱元法能够自然地结合自由表面地形、衰减和各向异性效应，并准确表示表面波传播。
研究的科学价值在于： 1. 提供了处理复杂三维地震波传播问题的高精度数值工具 2. 验证了谱元法对各类源（点力、矩张量）和介质（层状、各向异性、衰减）的适用性 3. 展示了该方法在陡峭地形和强衰减等挑战性问题中的可靠性
应用价值体现在： 1. 为区域和全球地震学中真实三维地球模型的精确地震图计算提供了必要工具 2. 为地震风险评估等实际应用问题奠定了基础 3. 展示了该方法在并行计算环境中的良好适应性
研究还讨论了网格设计的重要考量，指出六面体单元需要谨慎设计以适应速度变化和地形特征。提出了基于”立方球”（cubed sphere）的创新型网格划分方法，用于处理半球形陨石坑等复杂几何形状。
此外，论文探讨了多项式阶数选择（通常5-10为最佳）和网格分辨率（每个最小波长约5个点）的优化问题，为实际应用提供了指导。研究还指出，虽然空间精度高，但使用简单二阶时间方案可能最终限制整体精度，建议考虑更高阶时间推进方案的可能性。