地形效应影响下地震动参数与斜坡稳定性的相关性研究*

刘甲美 高孟潭 陈 鲲

（中国北京 100081 中国地震局地球物理研究所）

引言

震害调查表明，局部地形对震害分布具有重要影响．在1976年唐山大地震中，位于迁西县景忠山顶部的庙宇式建筑被严重破坏，烈度达Ⅸ度，而与其高度相差约300m的山脚周围村庄的破坏程度则明显较轻，烈度仅为Ⅵ度（杨宇等，2011）．2008年汶川MS8．0地震诱发了数以万计的滑坡（吴树仁等，2009），其所在的龙门山地区地形地貌比较复杂，地形高差悬殊，这种复杂的地形对地震动的影响是地震滑坡多发的一个重要原因（罗永红，王运生，2013）．与风暴形成的降雨型滑坡沿山体均匀分布不同，地震滑坡多分布在脊峰及山顶部位，具有明显的地形分布特征（Meunieretal，2008）．强震观测和数值模拟研究进一步表明，不规则地形对地震动的幅值、频谱等特性具有显著的影响（刘晶波，1996；Bouckovalas，2005；姜慧等，2007；Leeetal，2009；丁志华等，2014），这是地形影响震害分布的根本原因．

地震动是影响地震作用下斜坡稳定性的重要因素（王秀英等，2011，2012），地形对地震动的影响也会影响斜坡的稳定性（Pengetal，2009）．实际工作中，由于台站分布稀疏不一，往往无法获取不同位置的地震动记录而采用地震动参数表示地震动强度，其中峰值加速度和阿里亚斯（Arias）烈度为两个较为常用的地震动参数（Jibson，2007；王秀英等，2012；徐光兴等，2012）．山区具有复杂的地形地貌，地形对地震动的持时、振幅和频率等特性均会产生显著影响且作用较为复杂（刘洪兵，朱晞，1999；石玉成等，1999；卢育霞等，2011；郭明珠等，2013），而地震动峰值加速度和阿里亚斯烈度均不能完全涵盖地震动持时、振幅和频率等特性，因此考虑地形效应影响的地震动作用下的斜坡稳定性，需要选取与斜坡稳定性相关性较好的地震动参数．在评估地震滑坡稳定性的方法中，Newmark（1965）提出的累积位移模型，在斜坡稳定性评价中取得了较好的效果（Jibson，1993；Jibsonetal，2000；陈晓利等，2013；葛华等，2013；王涛等，2013）．该方法将滑坡体视为在斜坡上滑动的刚滑块，该滑块具有临界加速度，当地震动加速度记录值大于斜坡临界加速度时，滑块产生位移，该位移值的大小则反映了地震动作用下斜坡的稳定性（Jibson，2007）．

本文采用Newmark（1965）滑坡稳定评价模型，利用自贡地形影响台阵获取的汶川地震主震加速度记录，分析地形对地震动峰值加速度、阿里亚斯烈度以及90%地震动能量持时的影响，以期对比在地形效应影响下地震动峰值加速度和阿里亚斯烈度与地震作用下斜坡稳定性的相关性．

1 数据

图1 自贡地形影响台阵的台站位置剖面示意图（引自王海云和谢礼立，2010）Fig．1 Schematic illustration of the stations location of Zigong topography array（after Wang and Xie，2010）

自贡地形影响台阵位于四川省自贡市西山公园．该公园地处自流井凹陷区，总体构造形态为一个复式向斜，由一系列北东向背斜和向斜组成．该台阵场址位于自流井背斜附近，基本沿着山脊地形的轴线，观测对象为山脊地形．该台阵由8个台站（Z0，Z1，…，Z7）组成，均配置内置ES-T型加速度计的ETNA型数字强震动仪，其中台站Z0位于山底土层上，台站Z1位于同一高度的山底基岩上，其它6个台站位于山脊不同高度的侏罗系基岩上（杨宇等，2011）．图1为自贡地形影响台阵的台站分布剖面示意图．

汶川地震中自贡台阵的8个台站均获得了良好的三分量主震加速度记录．该台阵距离汶川地震震中约226．4km，而距离最远的Z0与Z7两个台站的水平距离约为0．4km（王海云，谢礼立，2010；凌代俭等，2013），因此可以认为地震波从震源传播到各个台站的传播路径是相同的，并且具有相同的震源效应，不同台站间地震动的差异受到局部场地效应的影响．台站Z0位于土层之上，其地震记录受到上覆土层和地形的共同影响，而台站Z1—Z7由于位于不同高度的基岩上，可以认为地形是影响其地震动记录差异的主要因素（杨宇等，2011）．

本文使用台站Z1—Z7的三分量加速度记录，以台站Z1为参考台站，研究地形对地震动峰值加速度、阿里亚斯烈度及90%地震动能量持时的影响，分析在地形效应影响下地震动峰值加速度和阿里亚斯烈度与斜坡稳定性的相关性．阿里亚斯烈度 （Arias，1970）定义为一系列单位质量的单自由度体系（自振周期0—∞）所吸收的地震动能量总和，综合了地震动多持时及强度信息，具体表示为

式中：IAx为x方向的阿里亚斯烈度，单位为m／s；ax（t）为x方向的地震加速度时程，单位为m／s2；g为重力加速度；Td为地震加速度记录总的持续时间，单位为s．90%地震动能量持时Ds（Dobryetal，1978）定义为阿里亚斯烈度达到5%—95%的时间间隔．各台站计算所得到的三分量地震动参数值如表1所示．

表1 各台站地震动峰值加速度、阿里亚斯烈度IA及90%地震动能量持时Ds一览表Table 1 Measured PGA，Arias intensity（IA），Dobry duration（Ds）for the seven stations

2 峰值加速度、阿里亚斯烈度及持时随地形高度的变化

以台站Z1为地形参考台站，将地形放大系数K定义为

式中，Yj为台站Z2—Z7的三分量地震动参数值（地震动峰值加速度、阿里亚斯烈度或地震动持时），Y1为台站Z1所对应分量的地震动参数值．

从图2可以看出，各台站的峰值加速度及其放大系数整体上随台站高度的增加而增大，并且地形效应对峰值加速度水平向的影响大于竖直向．水平向峰值加速度EW分量由台站Z1（高度0m）的0．225m／s2增大到台站Z6（高度72m）的0．412m／s2，放大系数为1．18—1．77；NS分量由台站Z1的0．265m／s2增大到台站Z7（高度52m）的0．451m／s2，放大系数为1．14—1．73．竖直向峰值加速度由台站Z1的0．147m／s2增大到台站Z6的0．196m／s2，放大系数为1．08—1．34．

图2 自贡地形影响台阵各台站峰值加速度PGA（a）及其放大系数K（b）随台站高度变化Fig．2 Variation of PGA（a）and its amplification factor K （b）with height of the stations of Zigong topography array

图3给出了阿里亚斯烈度及其放大系数随台站高度的变化．可以看出，地形场地对阿里亚斯烈度也有显著的放大效应，整体上阿里亚斯烈度及其放大系数随台站高度的增加而增大，地形效应对阿里亚斯烈度水平向的影响大于竖直向．水平向阿里亚斯烈度EW分量由台站Z1（高度0m）的0．046m／s增大到台站Z6（高度72m）的0．150m／s，放大系数为1．48—3．28；NS分量由台站Z1的0．055m／s增大到台站Z6的0．163m／s，放大系数为1．27—2．99．竖直向阿里亚斯烈度由台站Z1的0．022m／s增大到台站Z6的0．039m／s，放大系数为1．26—1．77．

图3 自贡地形影响台阵各台站阿里亚斯烈度IA（a）及其放大系数K（b）随台站高度的变化Fig．3 Variation of Arias intensity IA（a）and its amplification factor K （b）with height of the stations of Zigong topography array

不同台站的峰值加速度和阿里亚斯烈度并非完全随台站高度的增加而增大．例如：台站Z4的高度为47m，大于台站Z5的高度（45m），从表1可以看出，台站Z4NS分量和EW分量的峰值加速度和阿里亚斯烈度均小于台站Z5的相应数值；台站Z6（高度72m）NS分量的峰值加速度也小于台站Z7（高度52m）的相应数值．这表明地形场地对地震动的影响机制较为复杂，高度并非影响地形放大效应的唯一因素．此外，地形效应对地震动影响是多方面的，影响到地震动幅值和频谱等特性（刘晶波，1996；王海云，谢礼立，2010），因此单一的地震动参数（如峰值加速度或阿里亚斯烈度）无法全面地刻画地形对地震动的影响．对比图2与图3可以看出，峰值加速度和阿里亚斯烈度随台站的变化具有不同的特征．例如：台站Z7NS分量的峰值加速度放大系数为1．7，是所有台站NS分量峰值加速度放大系数的最大值；而其阿里亚斯烈度放大系数为2．3，并不是所有台站NS分量阿里亚斯烈度放大系数的最大值．因此，采用不同的地震动参数表征地形效应，其分析结果会有所差异．在地震滑坡的研究中，考虑地形放大效应的地震滑坡稳定性分析需要选择合适的地震动参数，以合理地反映地形效应对斜坡稳定性的影响．

图4给出了地震动持时及其放大系数随台站高度的变化．可以看出，地形场地对地震动持时也有一定的影响．山脊地形增大了地震动持时，不同地形高度的地震动记录持时并没有显著差异，不同分量的地震动持时也没有显著差异．参考台站Z1的三分量地震动持时约为90s，EW分量持时为90．46s，NS分量持时为89．89s，UD分量持时为91．46s；台站Z2—Z7的三分量地震动持时均接近120s，其放大系数约为1．33．

图4 各台站地震动持时Ds（a）及其放大系数K（b）随台站高度的变化Fig．4 Variation of Dobry duration Ds（a）and its amplification factor K （b）with height of the stations of Zigong topography array

3 地震动参数与斜坡稳定性的相关性

实际上，诸如台站Z1—Z7所记录到的地震动强度很难诱发斜坡失稳．王秀英等（2012）的研究表明，汶川地震中整个龙门山震区诱发地震滑坡的地震动峰值加速度最小值平均约为0．98m／s2，诱发地震滑坡的阿里亚斯烈度最小值约为0．2m／s．但现实中可能存在一些斜坡，在某些因素（如斜坡已遭受一次强震破坏）影响下已处于临界滑动状态，随之而来的微弱震动（余震等）有可能触发其失稳．鉴于此，考虑到台站Z1—Z7台站水平分量地震动峰值加速度最小值约为0．20m／s2，采用Newmark方法，本文假定存在一个临界加速度值为0．20m／s2而处于临界滑动状态的斜坡，使各台站的地震动强度均能对斜坡稳定性产生影响，计算其在台站Z1—Z7水平分量地震动作用下所产生的Newmark累积位移（表2），以对比考虑地形效应时峰值加速度和阿里亚斯烈度与斜坡稳定性的相关性．

由于台站布设密度不均匀等因素，实际工作中很难获取不同位置的地震动记录．故考虑地形效应进行区域上的地震滑坡评估，通过地震动参数乘以放大系数进行分析较为适宜．对比地形效应对阿里亚斯烈度和峰值加速度的影响（图2，3）可以看出，总体上二者随台站的高度变化趋势一致．例如台站Z3的阿里亚斯烈度NS和EW分量均高于台站Z2，其峰值加速度NS和EW分量也同样高于台站Z2．但是阿里亚斯烈度与峰值加速度的变化趋势也不尽相同．例如：台站Z6的阿里亚斯烈度NS分量为0．163m／s，大于台站Z7（0．125m／s），而台站Z6的峰值加速度 NS分量为0．421m／s2，小于台站Z7（0．451m／s2）．这是由于地形对地震动作用较为复杂，对地震动不同频段具有不同的影响所致（王海云，谢礼立，2010；杨宇等，2011；唐晖等，2012；凌代俭等，2013）．而峰值加速度和阿里亚斯烈度仅表征地震动的部分信息，地形效应反映在不同的地震动参数上则有所差异，这种差异会导致选用不同的地震动参数分析叠加地形效应影响的地震滑坡稳定性时，不能反映真实的地形影响．

表2 各台站EW和NS分量的Newmark累积位移（DN）Table 2 Measured Newmark’s cumulative displacements DNof EW and NS components at the stations

图5给出了利用各台站EW分量加速度记录计算的Newmark累积位移随峰值加速度和阿里亚斯烈度的变化．可以看出：Newmark位移值按台站Z1，Z2，Z4，Z3，Z7，Z5，Z6的顺序依次增大，峰值加速度值按台站Z1，Z2，Z4，Z7，Z5，Z3，Z6的顺序增大，除台站Z3具有明显偏离外，Newmark位移值随峰值加速度的增大而增大；阿里亚斯烈度值按台站Z1，Z2，Z4，Z3，Z5，Z7，Z6的顺序增大，除台站Z5（相对台站Z7）稍有偏离外，Newmark位移值也随阿里亚斯烈度的增大而增大．

图5 自贡地形影响台阵各台站EW分量Newmark位移DN随峰值加速度PGA（a）和阿里亚斯烈度IA（b）的变化Fig．5 Variation of EW component Newmark displacement DNwith PGA （a）and Arias intensity IA（b）for the stations of Zigong topography array

图6 自贡地形影响台阵各台站NS分量Newmark位移DN随峰值加速度（a）和阿里亚斯烈度IA（b）的变化Fig．6 Variation of NS component Newmark displacements DNwith PGA （a）and Arias intensity IA（b）for the stations of Zigong topography array

图6给出了利用各台站NS分量加速度记录计算的Newmark累积位移随峰值加速度和阿里亚斯烈度的变化．可以看出：Newmark位移总体随着峰值加速度增大而增大，但二者的变化趋势在一些台站有所不同，台站Z6NS分量峰值加速度为0．421m／s2，与台站Z5NS分量峰值加速度相同而小于台站Z7（0．451m／s2），而斜坡在台站Z6地震动加速度NS分量作用下所产生的Newmark位移为1．07cm，大于斜坡在台站Z5的0．88cm和Z7台站的0．59cm；阿里亚斯烈度值和Newmark位移值以相同的台站顺序同时增大，即Newmark位移完全随阿里亚斯烈度的增大而增大．

总体上看：Newmark位移随峰值加速度和阿里亚斯烈度的增大而增大，表明以峰值加速度或阿里亚斯烈度为地震动参数都可以很好地表征地形效应对斜坡稳定性的影响；相对于Newmark位移随峰值加速度曲折型的变化趋势，Newmark位移随阿里亚斯烈度增大而增大的趋势近乎呈直线型，表明阿里亚斯烈度可以更好地表征地形对地震动的影响，与在地震动作用下斜坡稳定性的相关性更强．

4 讨论与结论

使用自贡地形影响台阵记录到的2008年汶川地震主震加速度记录，分析了地震动峰值加速度、90%能量持时及阿里亚斯烈度随地形高度的变化，探讨了地形效应作用下峰值加速度和阿里亚斯烈度与地震作用下斜坡稳定性的相关性，得到以下结论：

1）地形对地震动的影响较为复杂，反映在不同地震动参数如峰值加速度和阿里亚斯烈度上具有不同的特征．地形场地对峰值加速度和阿里亚斯烈度均有显著的放大效应，整体上其放大效应随台站高度的增加而增大，水平向的放大效应大于竖直向，即水平向峰值加速度放大系数为1．1—1．8，阿里亚斯烈度放大系数为1．2—3．3；竖直向峰值加速度放大系数为1．1—1．3，阿里亚斯烈度放大系数为1．2—1．7．

2）地形场地对地震动持时也有一定的放大效应，不同高度及不同分量的放大效应没有显著差异，其放大系数均为1．3左右，这表明持时的地形放大效应或许不是地形加剧滑坡震害的重要影响因素．

3）地形对地震动强度、持时具有放大效应，且影响较为复杂，因此考虑地震滑坡稳定性分析时需要选择合适的地震动参数来表征地形放大效应对斜坡破坏的影响．本文结果表明，阿里亚斯烈度和峰值加速度均能很好地表征地形对地震动的影响，与地震动对斜坡稳定性的影响具有很强的相关性；与峰值加速度相比，阿里亚斯烈度可以更好地表征地形对地震动的影响，与地震动作用下斜坡稳定性的相关性更强，是考虑地形效应影响的地震滑坡危险性的优选地震动参数．

国家强震动台网中心（CSMNC）为本文提供数据，中国地震局地球物理研究所王玉石副研究员提出建议，审稿专家对本文提出修改意见，作者在此一并表示感谢．

陈晓利，袁仁茂，庾露．2013．Newmark方法在芦山地震诱发滑坡分布预测研究中的应用［J］．地震地质，35（3）：661-670．

Chen X L，Yuan R M，Yu L．2013．Applying the Newmark’s model to the assessment of earthquake-triggered landslides during the Lushan earthquake［J］．SeismologyandGeology，35（3）：661-670（in Chinese）．

丁志华，周红，蒋涵．2014．三维台阶地形地震动效应研究［J］．地震学报，36（2）：184-199．

Ding Z H，Zhou H，Jiang H．2014．Effect of 3-D step topography on ground motion［J］．ActaSeismologicaSinica，36（2）：184-199（in Chinese）．

葛华，陈启国，王德伟．2013．地震滑坡危险性评价及编图：以映秀震中区为例［J］．中国地质，40（2）：644-652．

Ge H，Chen G Q，Wang D W．2013．The assessment and mapping of seismic landslide hazards：A case study of Yingxiu area，Sichuan Province［J］．GeologyinChina，40（2）：644-652（in Chinese）．

郭明珠，赵芳，赵凤仙．2013．场地地震动局部地形效应研究进展［J］．震灾防御技术，8（3）：311-318．

Guo M Z，Zhao F，Zhao F X．2013．A review of the effect of small-scale surface topography on ground motions［J］．TechnologyforEarthquakeDisasterPrevention，8（3）：311-318（in Chinese）．

姜慧，高孟潭，俞言祥，唐丽华．2007．北天山地区S波非弹性衰减和场地效应研究［J］．地震学报，29（2）：181-186．

Jiang H，Gao M T，Yu Y X，Tang L H．2007．S wave inelastic attenuation and site effect in the northern Tianshan area［J］．ActaSeismologicaSinica，29（2）：181-186（in Chinese）．

凌代俭，陈伟，袁建力．2013．用双向反应谱比法估计地形对结构反应的影响［J］．中国地震，29（2）：256-264．

Ling D J，Chen W，Yuan J L．2013．Evaluation of the seismic amplification of topography on structure with bidirectional response spectra ratio method［J］．EarthquakeResearchinChina，29（2）：256-264（in Chinese）．

刘洪兵，朱晞．1999．地震中地形放大效应的观测和研究进展［J］．世界地震工程，15（3）：20-25．

Liu H B，Zhu X．1999．Advance on topographic amplification effects of seismic response［J］．WorldInformationon EarthquakeEngineering，15（3）：20-25（in Chinese）．

刘晶波．1996．局部不规则地形对地震地面运动的影响［J］．地震学报，18（2）：239-245．

Liu J B．1996．Local effects of irregular topography on earthquake ground motion［J］．ActaSeismologicaSinica，18（2）：239-245（in Chinese）．

卢育霞，刘琨，姚凯，石玉成，李倩．2011．甘肃文县上城台地的地震记录分析［J］．西北地震学报，33（S1）：393-397．

Lu Y X，Liu K，Yao K，Shi Y C，Li Q．2011．Analysis on seismic records observed on Shangcheng platform in Wenxian，Gansu Province［J］．NorthwesternSeismologicalJournal，33（S1）：393-397（in Chinese）．

罗永红，王运生．2013．汶川地震诱发山地斜坡震动的地形放大效应［J］．山地学报，31（2）：200-210．

Luo Y H，Wang Y S．2013．Mountain slope ground motion topography amplification effect induced by Wenchuan earthquake［J］．JournalofMountainScience，31（2）：200-210（in Chinese）．

石玉成，王兰民，张颖．1999．黄土场地覆盖层厚度和地形条件对地震动放大效应的影响［J］．西北地震学报，21（2）：203-208．

Shi Y C，Wang L M，Zhang Y．1999．Effects of overburden thickness and topography of loessial sites on earthquake ground motion［J］．NorthwesternSeismologicalJournal，21（2）：203-208（in Chinese）．

唐晖，李小军，李亚琦．2012．自贡西山公园山脊地形场地效应分析［J］．振动与冲击，31（8）：74-79．

Tang H，Li X J，Li Y Q．2012．Site effect of topograghy on ground motions of Xishan park of Zigong city［J］．Journalof VibrationandShock，31（8）：74-79（in Chinese）．

王海云，谢礼立．2010．自贡市西山公园地形对地震动的影响［J］．地球物理学报，53（7）：1631-1638．

Wang H Y，Xie L L．2010．Effects of topography on ground motion in the Xishan park，Zigong city［J］．ChineseJournal ofGeophysics，53（7）：1631-1638（in Chinese）．

王涛，吴树仁，石菊松，辛鹏．2013．基于简化Newmark位移模型的区域地震滑坡危险性快速评估：以汶川MS8．0级地震为例［J］．工程地质学报，21（1）：16-24．

Wang T，Wu S R，Shi J S，Xin P．2013．Case study on rapid assessment of regional seismic landslide hazard based on simplified Newmark displacement model：WenchuanMS8．0earthquake［J］．JournalofEngineeringGeology，21（1）：16-24（in Chinese）．

王秀英，聂高众，王松．2011．汶川地震诱发滑坡的地震动加速度评判标准［J］．地震学报，33（1）：82-90．

Wang X Y，Nie G Z，Wang S．2011．Ground motion acceleration criterion for judging landslide induced by the 2008Wenchuan earthquake［J］．ActaSeismologicaSinica，33（1）：82-90（in Chinese）．

王秀英，聂高众，马牧军．2012．地震滑坡灾害评估中地震影响因素的联合应用［J］．地震学报，34（1）：76-84．

Wang X Y，Nie G Z，Ma M J．2012．Application of multiple ground motion factors in earthquake-induced landslide hazard evaluation［J］．ActaSeismologicaSinica，34（1）：76-84（in Chinese）．

吴树仁，石菊松，姚鑫，王涛，汪华斌．2009．四川汶川地震地质灾害活动强度分析评价［J］．地质通报，27（11）：1900-1906．

Wu S R，Shi J S，Yao X，Wang T，Wang H B．2009．Analysis and evaluation of geohazard intensity of the Wenchuan earthquake，Sichuan，China［J］．GeologicalBulletinofChina，27（11）：1900-1906（in Chinese）．

徐光兴，姚令侃，李朝红，王晓放．2012．基于汶川地震强震动记录的边坡永久位移预测模型［J］．岩土工程学报，34（6）：1131-1136．

Xu G X，Yao L K，Li C H，Wang X F．2012．Predictive models for permanent displacement of slopes based on recorded strong-motion data of Wenchuan earthquake［J］．ChineseJournalofGeotechnicalEngineering，34（6）：1131-1136（in Chinese）．

杨宇，李小军，贺秋梅．2011．自贡西山公园山脊场地地形和土层效应数值模拟［J］．震灾防御技术，6（4）：436-447．

Yang Y，Li X J，He Q M．2011．Numerical simulation for site effect of ridge terrain and overlaying soil in Zigong Xishan park［J］．TechnologyforEarthquakeDisasterPrevention，6（4）：436-447（in Chinese）．

Arias A．1970．A measure of earthquake intensity［G］∥SeismicDesignforNuclearPowerPlants．Cambridge，Mass：The MIT Press：438-483．

Bouckovalas G D，Papadimitriou A G．2005．Numerical evaluation of slope topography effects on seismic ground motion［J］．SoilDynEarthqEng，25（7／8／9／10）：547-558．

Dobry R，Idriss I M，Ng E．1978．Duration characteristics of horizontal components of strong-motion earthquake records［J］．BullSeismolSocAm，68（5）：1487-1520．

Jibson R W．1993．Predicting earthquake-induced landslide displacements using Newmark’s sliding block analysis［G］∥TranspResRec1441．Washington：Nat Acad Press：9-17．

Jibson R W．2007．Regression models for estimating coseismic landslide displacement［J］．EngGeol，91（2）：209-218．

Jibson R W，Harp E L，Michael J A．2000．A method for producing digital probabilistic seismic landslide hazard maps［J］．EngGeol，58（3）：271-289．

Lee S J，Komatitsch D，Huang B S，Tromp J．2009．Effects of topography on seismic-wave propagation：An example from northern Taiwan［J］．BullSeismolSocAm，99（1）：314-325．

Meunier P，Hovius N，Haines J A．2008．Topographic site effects and the location of earthquake induced landslides［J］．EarthPlanetSciLett，275（3／4）：221-232．

Newmark N M．1965．Effects of earthquakes on dams and embankments［J］．Géotechnique，15（2）：139-160．

Peng W F，Wang C L，Chen S T，Lee S T．2009．Incorporating the effects of topographic amplification and sliding areas in the modeling of earthquake-induced landslide hazards，using the cumulative displacement method［J］．CompGeosci，35（5）：946-966．