Numerical simulation of dynamic process of earthquake
基于Fluent内燃机压缩过程数值模拟-孙宝庆
(a) smoothing 参数设置
(b) remeshing 参数设置
图 4 动网格参数设置
图 5 活塞运动参数设置
(a)5oCA 网格
(b)60oCA 网格
(c)120oCA 网格 图 6 不同曲轴转角下动网格模型
(d)180oCA 网格
4
模拟结果分析
4.1 速度场分析
如图 7 所示，本文给出了，在不同曲轴转角下，速度场云图。由图 7 可知，在不同曲 轴转角下速度场云图变化较大，在 0-180oCA 进气过程中，活塞由上止点到下止点运动 140mm,速度最大值变化较大，在 60oCA 时，最大速度为 237ms/s,在 120oCA 时最大速度为 148m/s，而到 180oCA 最大速度为 55.8m/s。由此可得在进气过程中，速度由零逐渐增大， 而后逐渐降低，在进气末尾速度逐渐降低。且由云图可知，在运动过程中存在着漩涡。这 更加有利有在发动机后期燃烧过程中，燃气混合。
[参考文献] [1] 杨忠敏.车用动力柴油化的技术研究与开发方向[J].柴油机设计与制造，2004,(01) [2] 王东,黄震,张有等.螺旋进气道三维流场数值模拟[J].计算机辅助工程， 2007,(09) [3] 曹永悌. 基于 Fluent 对康明斯 ISM370 柴油掺氢发动机燃烧过程的数值模拟[C].郑州： 华北水利水电大学， 2012 [4] 张国庆. 内燃机进排气过程数值模拟[C].太原：中北大学,2013
1
前言
气缸是发动机压缩过程中最主要的工作部件，对缸内气流流动数值模拟也显的格外重 要，同时由于目前随着柴油机转速逐渐提高，柴油机的动力、经济性指标和环保指标的提 高也越来越受到人们的重视。 而这些性能指标的提高在很大程度上受限于燃料的燃烧程度。 气流动状态则对柴油机燃烧室内的混合气生成的质量有着很大的影响。进气道的形状对柴 油机的进气性能有着很大影响，进而也影响着整机的性能 随着有限元技术的成熟，在包括汽车发动机在内的几乎所有工程领域得到愈来愈广泛 的运用 [1] 。有限元技术在汽车发动机零部件设计上的应用也越发成熟。计算流体动力学的 应用可以提高直喷柴油机进气道设计的整体性能，探求气道内每一点气体流动的变化，使 我们了解相应流场的微观结构。同时可以缩短研发周期，大大推动柴油机工业的发 [2]。
人工冻土冻结过程中热_力耦合的数值模拟方法研究
；c 为比热（J/kg·℃） ；t 为温度（℃） 。在非 （kg/m3） 相变区，式（3）中的 qv / λ = 0 。 （2）本构方程 当土体中作用一个温度场时，除考虑力场作用下 的应力和应变外，还需要考虑温度场对土体力场的影 响。土中的应力会随着温度的变化而变化，这种由于 温度变化产生的应力被称为热应力。 按照弹性力学方法，建立热应力方程。其中平衡 方程只与物体受力有关，而与产生力的原因无关；几 何方程中应变只与位移有关，而与引起位移的原因无 关， 所以平衡方程和几何方程的形式不变。 在变温时， 弹性体的应变由两部分组成。 （1）自由膨胀引起的正应变分量 ε 0 = αΔt ，其中 α 为膨胀系数，Δt 为温度变化值。自由膨胀引起的剪 应变分量为 0。 （2） 热膨胀时， 土体的应力–应变关系符合广义 胡克定律： 1 ε x = E (σ x − μ (σ y + σ z )) + αΔt ， 1 (4) ε y = (σ y − μ (σ x + σ z )) + αΔt ， E 1 ε z = E (σ z − μ (σ y + σ x )) + αΔt 。 联合温度场微分方程和热力本构关系，可以得到 热力耦合控制方程。 在进行冻胀位移计算时，地铁隧道衬砌管片设为 线弹性材料；土体亦假设为线弹性材料，但其弹性模 量的取值随温度变化而变化；土体的冻胀率通过原状 土的冻胀试验确定。 2.2 计算参数的取值方法 土样的物理力学性质指标、热物理参数随冻结过 程存在较大的变化，若采用常量指标进行运算不符合 实际情况，本文根据室内试验得出导热系数、比热以 及其他指标随温度场的变化规律， 在进行数值计算时， 通过建立比热、导热系数随温度场变化的场变量赋予 冻土的热物理参数，同时也考虑了模量随温度场的变 化。从而建立了计算参数随温度场变化的热力耦合的 数值计算方法。 本文结合天津地铁某联络通道的人工冻土工程， 利用上述方法代入随温度场变化的参数指标进行计 算，具体取值见下文第 3.2 节。
注水诱导裂缝动态特征及数值模拟研究
嵌入式离散裂缝模型( EDFM) 由 Lee 等 [14] 提
出,其主要优点为裂缝系统和基质系统相对独立,
改变裂缝的分布时无需重新对基质网格进行划分。
Yang 等 [15] 基于 EDFM 提出了一种适用于注水诱
导裂缝的高效数值模拟方法,即动态嵌入式离散裂
缝方法( dEDFM) ,但仅研究了张性裂缝的直线扩
展过程,而受古构造运动的影响,砂岩油藏中一般
发育天然裂缝,在注水诱导裂缝的扩展过程中发生
主裂缝与天然裂缝相遇形成交叉裂缝扩展的情况,
导致 dEDFM 的准确性降低。 因此, 对 dEDFM 进
行前处理算法的动态改进,使其能够模拟复杂天然
裂缝分布下注水诱导裂缝的扩展过程。
模型中考虑油水两相的渗流问题。 以油相为
力,Pa;p of 为与井筒相连裂缝单元的流体压力,Pa;
r em 为与基质网格相连井筒的等效半径,m;r ef 为与
1 高效注水诱导裂缝数值模拟方法
▽(
q ofw 为油相由井筒向裂缝的注入量,m 3 / s;p wf 为井
(1)
(2)
式中:K m 、K f 分别为基质渗透率和裂缝渗透率,mD;
的距离,m;θ 为目标方向与 x 轴的夹角,°;K Ⅰ 为Ⅰ
Key words: ultra-low permeability reservoir; waterflood-induced fracture; cross fracture; numerical simulation
0 引 言
鄂尔多斯盆地 Y2 超低渗透油藏存在初期产
量低、产能递减迅速等问题,通常采用面积注水的
模拟方法能够高效模拟油藏注水开发过程中裂缝的动态变化,反映天然裂缝对于主裂缝扩展
莺琼盆地特殊岩性识别方法研究及应用
莺-琼盆地位于印支半岛与南海北部大陆架交接区，在板块构造位置上位于欧亚板块、印支板块和太平洋板块的交汇处。
莺歌海盆地形态呈 NNW 走向的长条形，由东南部的莺歌海凹陷和西北部的河内凹陷组成，两凹陷间被临高凸起所分隔[1]。
莺琼盆地是南海西部最重要的天然气勘探战场。
近年来随着二维、三维资料的覆盖及连续的钻探，发现该盆地存在较多的特殊岩性。
在莺歌海盆地，钻探的主要特殊岩性为低速泥岩，诸多的低速泥岩在地震特征上均有“亮点”强振幅特征，因此给勘探目标搜寻及评价带来了巨大的挑战。
从莺歌海及琼东南盆地已钻井出发，主要研究低速泥岩识别。
首先，理清了各区岩石物理规律，总结了低速泥岩形成机理，形成了一套识别低速泥岩的方法[2-3]，排除了特殊岩性对目标的影响，确定了莺琼盆地有利目标特征，为南海西部天然气钻探奠定了基础。
传统 “低频强振幅亮点”反射特征的地震异常体通常认为是含气砂岩的响应，而莺歌海盆地钻探结果表明，该种特征泥岩广泛发育。
经过统计发现盆地中此类泥岩多呈现低速度特征[4]。
其中乐东区钻遇低速泥岩井十余口，而东方区则达到二十余口。
因此了解低速泥岩形成机理，准确识别泥岩，确定砂岩储层，降低勘探风险十分关键。
1 泥岩速度影响因素地质、地球物理、测井等资料综合研究发现地层压力、泥质含量、有机质丰度均会对泥岩的纵波速度产生影响。
压力一定时，泥质含量超过25%时，泥质含量越高，纵波速度越低。
这也代表泥岩越纯，纵波速度越低。
当泥质含量相差不大时，有效压力（上覆压力减去孔隙压力）的减小时，泥岩速度也随之变小。
在泥质含量、孔隙压力等差异较小时，泥岩的有机质增加，即干酪根含量增大，泥岩的纵波速度减小。
带着上述认识，从本区已钻井出发，来剖析三种影响泥岩速度的因素。
由于地层孔隙压力大导致了泥岩中的微裂缝增加了孔隙流体间的相互联系，造成泥岩低电阻率和含水导致的地震波声传播异常低速[5]。
而实际中，莺歌海盆地在黄流组二段普遍发育异常超压地层，甚至部分在黄流组一段亦存在高压。
用随机方法模拟唐山7.1级余震地震动场
用随机方法模拟唐山7.1级余震地震动场李启成;宋志勇【摘要】在系统介绍用随机方法模拟地震动的基础上,用该方法模拟了唐山7.1级余震地震动场.模拟结果表明,唐山地区发生Ms7.1级余震时,该地区产生了强烈的地震动,很多区域加速度峰值超过400 cm/s2,地震断层面长65 km,宽13.5 km,地震过程中断层面的平均滑动为1.08 m.模拟与记录的地震动持时符合得不好,其原因是随机方法中的经验关系大都是在美国大陆得到的,所以建议加强中国大陆相应的地震动经验关系的研究,使随机方法能更好地适用于中国大陆.尽管如此,由于用随机方法模拟的地震动峰值和反应谱与记录符合得比较好,这些模拟结果可以作为建筑物抗震设防的地震动输入.【期刊名称】《国际地震动态》【年(卷),期】2012(000)007【总页数】7页(P13-19)【关键词】地震动;迁安;随机方法【作者】李启成;宋志勇【作者单位】辽宁工程技术大学,阜新123000;黑龙江科技学院,哈尔滨150027;黑龙江科技学院,哈尔滨150027【正文语种】中文【中图分类】P315.9引言地震动的模拟是重要的基础研究项目，对于理解地震断层的破裂过程、地震波的传播过程和地震动的场地反应有十分重要的意义。
在地震记录不十分丰富的今天，地震动模拟研究显得十分重要，尤其是对建筑物有针对性的抗震设防提供地震动输入所起到的作用更是少量的地震记录无法替代的。
长期以来，日本更多地采用经验格林函数方法模拟地震动[1]。
上世纪末，美国和加拿大等国家开始更多地采用随机方法模拟地震动[2，3]。
经验格林函数方法在我国研究比较早[4，5]，随机方法在本世纪初才引进到中国[6]。
经验格林函数方法模拟地震动的效果比较好，但存在着缺乏小震记录的问题。
中国科学家在本世纪初已经开始用随机方法模拟中国大陆地震动，但有待于在更大范围内推广。
本文通过用随机方法模拟唐山MS 7.1级余震的地震动场来进一步证明用该方法模拟中国大陆地震动的有效性。
一种近场地震结构振动半主动控制方法★
第 55 卷 第 2 期 2024 年 3 月Vol.55 No.2Mar. 2024太原理工大学学报J OURNAL OF TAIYUAN UNIVERSITY OF TECHNOLOGY一种近场地震结构振动半主动控制方法樊长林1，项雪琦2（1.山西工程科技职业大学 建筑工程学院，山西 晋中 030619；2.太原理工大学 土木工程学院，太原 030024）摘 要：【目的】发展一种近场地震结构振动半主动控制方法。
采用神经模糊模型来描述SMA 应变率、温度相关的非线性力学行为。
【方法】基于非支配遗传算法设计了自组织模糊控制器实现对结构的半主动控制，结合近场地震特点，以最大层间位移和各层绝对加速度峰值为性能控制指标，调整模糊管理规则实现多目标启发式优化。
最后采用非线性时程分析程序，比对分析了受控和无控三层框架，研究了该方法对近场地震结构振动控制性能。
【结果】研究结果表明，近场地震作用下，在不增加结构加速度响应的基础上，有效地减小主体结构位移和损伤，实现了对结构振动的有效控制。
关键词：半主动控制；变摩擦阻尼器；形状记忆合金；近场地震；非支配遗传算法中图分类号：TU352.1 文献标识码：ADOI ：10.16355/j.tyut.1007-9432.20220841 文章编号：1007-9432（2024）02-0348-08A Semi -active Control Method for Structural Vibration Subject toNear -field EarthquakesFAN Changlin 1, XIANG Xueqi 2（1.College of Architectural Engineering ， Shanxi Vocational University of Engineer Science and Technology ，Jinzhong 030619， China ； 2.College of Civil Engineering ， Taiyuan University of Technology ， Taiyuan030024， China ）Abstract ：【Purposes 】 In this paper ， a semi -active control method for structural vibration subject to near -field earthquakes is proposed. A neuro -fuzzy model is used to capture rate and temperature de⁃pendent nonlinear behavior of the SMA components of the hybrid device. 【Methods 】According to characteristics of near -field earthquakes ， a non -dominated genetic algorithm is employed to design a self -organizing fuzzy controller. To obtain optimal controllers for concurrent mitigation of displace⁃ment and acceleration responses ， fuzzy governing rules is tuned with the multi -objective heuristic opti⁃mization. Then ， numerical simulation of a multi -story building is conducted to evaluate the perfor⁃mance of the hybrid device. 【Findings 】 Results show that a self -centering variable friction device modulated with a fuzzy logic control strategy can effectively reduce structural deformations without in⁃creasing acceleration response during near -field earthquakes.Keywords ：semi -active control ； variable friction damper ； shape memory alloy ； near -field引文格式：樊长林，项雪琦.一种近场地震结构振动半主动控制方法［J ］.太原理工大学学报，2024，55（2）：348-355.FAN Changlin ，XIANG Xueqi.A semi -active control method for structural vibration subject to near -field earthquakes ［J ］.Journal of Taiyuan University of Technology ，2024，55（2）：348-355.收稿日期：2022-10-29；修回日期：2023-03-01基金项目：国家自然科学基金资助项目（12172244）；山西省教育教学规划项目（GH -19266）通信作者：樊长林（1975-），博士，教授，国家一级注册结构师，主要从事结构工程抗震及振动控制研究，（E -mail ）fanc17538@第 2 期樊长林，等：一种近场地震结构振动半主动控制方法earthquake； non-dominated genetic algorithm随着经济的发展和城市人口密度的增加，地震造成的灾难越来越严重，人们对结构抗震性能的要求随之提高。
冲击波和气泡作用下舰船结构动态响应的数值模拟
冲击波和气泡作用下舰船结构动态响应的数值模拟王诗平;孙士丽;张阿漫;陈海龙【摘要】针对水下爆炸载荷、瞬态流固耦合效应在舰船水下爆炸数值模拟中的困难,在现有水下爆炸载荷计算模型(Geers and Hunter)的基础上,结合边界元法,修正水下爆炸气泡载荷计算方法.针对用二阶双渐近法(the second-order doubly-asymptotic approximation,DAA2)在处理低频大幅运动流固耦合问题时的局限性,提出非线性双渐近法(nonlinear doubly-asymptotic approximation,NDAA),计入了舰船大幅低频运动、自由面效应等非线性因素,并与结构有限元程序相结合,形成一套适合于工程应用研究的舰船水下爆炸数值模拟程序.数值模拟结果与实验结果的平均误差在16.8％左右,说明NDAA方法可为冲击波和气泡对舰船结构毁伤数值模拟提供参考.%Aimed at the difficulties of underwater explosion load and instantaneous fluid-structure interaction in the numerical simulation of dynamic response of warship structures, based on the Geers-Hupter model, a modified model for calculating the underwater explosion bubble load was proposed by combining with the boundary element method. To overcome the shortcoming of the second-order doubly-asymptotic approximation (DAA2) method in treating the fluid-structure interaction with low frequency and large amplitude motion, a nonlinear double asymptotic approximation (NDAA) method was presented by considering the nonlinear factors such as the violent low-frequency response of the warship and its free surface effect. And a code was developed by combining the presented NDAA method with the finite element code to numerically investigate the dynamic response of warship structuressubjected to underwater explosion Shockwaves and bubbles. The mean error between numerical and experimental vertical vibration velocities at the different positions is about 16. 8% , and it displays that the presented NDAA method is feasible.【期刊名称】《爆炸与冲击》【年(卷),期】2011(031)004【总页数】6页(P367-372)【关键词】爆炸力学;非线性双渐近法;边界元法;舰船;水下爆炸;冲击波;气泡【作者】王诗平;孙士丽;张阿漫;陈海龙【作者单位】哈尔滨工程大学船舶工程学院,黑龙江哈尔滨150001;哈尔滨工程大学船舶工程学院,黑龙江哈尔滨150001;哈尔滨工程大学船舶工程学院,黑龙江哈尔滨150001;哈尔滨工程大学船舶工程学院,黑龙江哈尔滨150001【正文语种】中文【中图分类】O382.1水下爆炸主要包括冲击波和气泡［1-13］，通常冲击波对舰船结构造成严重的局部毁伤，气泡对舰船造成总体和局部双重毁伤。
【16】Numerical simulation of blasting-induced rock fractures(1)
School of Civil and Environmental Engineering, Nanyang Technological University, Singapore 639798, Singapore Received 3 January 2007; received in revised form 3 December 2007; accepted 4 December 2007 Available online 20 February 2008
Abstract In the present study, the Johnson–Holmquist (J–H) material model is implemented into the commercial software LS-DYNA through user-subroutines to simulate the blasting-induced rock fractures. The J–H model consists of strength models for both intact and fully fractured materials, a polynomial equation of state, and a damage model that represents the material from an intact state to a fully fractured state. Influences of the key parameters in smooth blasting, viz., loading rate, distance from a free face, earth stress, and preexisting joint planes, etc., on fracture patterns are explored. According to the simulation results, the rock fracture pattern is significantly influenced by the loading rate. Fracture control techniques, namely, notched borehole and charge holder with slits are also simulated. Effectiveness of the fracture control techniques is demonstrated. The numerical simulation in the present study reproduces some of the well-known phenomena observed by other researchers. It has the potential to be applied in practical blast control and gas and hydraulic fracturing engineering. r 2007 Elsevier Ltd. All rights reserved.
2017年9月8日墨西哥沿岸Mw8.2级地震海啸观测数据分析与模拟
第41卷 第5期海 洋 学 报V o l .41,N o .52019年5月H a i y a n g Xu e b a o M a y 2019孙立宁,于福江,王培涛.2017年9月8日墨西哥沿岸M w 8.2级地震海啸观测数据分析与模拟[J ].海洋学报,2019,41(5):35-46,d o i :10.3969/j.i s s n .0253-4193.2019.05.004S u nL i n i n g ,Y uF u j i a n g ,W a n g P e i t a o .N u m e r i c a l s i m u l a t i o n a n d o b s e r v a t i o n a n a l y s i s o fM e x i c o M w 8.2e a r t h q u a k e t s u n a m i o n 8S e p-t e m b e r ,2017[J ].H a i y a n g X u e b a o ,2019,41(5):35-46,d o i :10.3969/j.i s s n .0253-4193.2019.05.0042017年9月8日墨西哥沿岸M w 8.2级地震海啸观测数据分析与模拟孙立宁1,于福江1,2*,王培涛1,2(1.国家海洋环境预报中心,北京100081;2.海洋灾害预报技术研究国家海洋局重点实验室,北京100081)收稿日期:2018-05-07;修订日期:2018-07-05㊂基金项目:国家重点研发计划项目(2018Y F C 1407000,2016Y F C 1401500)㊂作者简介:孙立宁(1993 ),男,河北省沧州市人,主要从事海洋灾害预警技术研究㊂E -m a i l :s 396647854@y e a h .n e t*通信作者:于福江(1965 ),男,研究员,主要从事海洋灾害应急管理及预警技术研究㊂E -m a i l :y v f u j i a n g@163.c o m 摘要:2017年9月8日4时49分(U T C ),墨西哥瓦哈卡州沿岸海域(15.21ʎN ,93.64ʎW )发生M w 8.2级地震,震源深度30k m ㊂强震在该海域引发海啸,海啸对震源附近数百千米范围内造成了严重影响㊂位于太平洋上的多个海啸监测网络捕捉到了海啸信号并详细记录了此次海啸的传播过程㊂本文选用了近场2个D A R T 浮标和6个验潮站的水位数据,通过潮汐调和分析和滤波分离出海啸信号,对近场海啸特征值进行了统计分析,并采用小波变换分析方法进一步分析了海啸的波频特征㊂基于O k a d a弹性位错理论断层模型计算得到了强震引发的海底形变分布,并采用M O S T 海啸模式对本次海啸事件近场传播特征进行了模拟,模拟结果与观测吻合较好㊂最后,基于实测和模拟结果,详细分析了此次地震海啸的近场分布特征,发现除受海啸源的强度和几何分布特征影响外,近岸海啸波还主要受地形特征控制,在与特定地形相互作用后波幅产生放大效应,会进一步加剧海啸造成的灾害㊂关键词:地震海啸;数值模拟;海啸波中图分类号:P 731.25文献标志码:A文章编号:0253-4193(2019)05-0035-121 引言海啸一般是由海底地震㊁火山喷发㊁海底或海岸滑坡及小行星撞击等因素产生的破坏性波浪㊂地震海啸主要因为由海底形变导致海水自海底至海表的巨大扰动所形成的长周期重力波㊂海啸波抵达近岸时,由于水深变浅,波速变缓,波长变短导致波高迅速增加,对海岸带产生剧烈的破坏作用[1]㊂海啸作为一种突发性且破坏力极强的海洋灾害,近年来受到广泛关注㊂据统计,全球范围内平均每年均有破坏性海啸事件发生,远高于20世纪每6年1次的频率,海啸灾害已造成245622人死亡,129658栋房屋受损,直接经济损失超过2500亿美元㊂可见海啸已成为世界上最严重的的海洋灾害之一,严重威胁着全球沿海居民生命财产安全[1-3]㊂在经历了2004年印度洋大海啸㊁2010年智利大海啸㊁2011年日本东北大海啸等重大海啸事件后,各沿海国家逐步提高了应对海啸灾害的关注度,相继建立了适合本区域的海啸预警系统,初步具备了应对越洋海啸和区域海啸的能力,但有效应对近场海啸的能力仍存在较大不足[4]㊂主要因为近场海啸留给预警系统响应的时间有限,另一方面的原因在于近场海啸与海岸地形间存在着复杂的相互作用关系,增加了预报难度㊂近年来,随着现代观测技术和海啸数值模式模拟能力的提高,对海啸在近岸传播的特征也不断有新的发现,于福江等[5]通过谱分析和小波分析方法发现海啸初波激发的海湾或港口内的局地振荡叠加在原始波形上,往往会使水位大幅升高;O m i r a等[6]通过对智利沿岸的3起地震海啸事件的对比得出海啸在近岸的传播明显受到智利大陆架的影响;张金凤[7]通过对不同地震源产生的海啸波在近岸产生的共振响应分析,认为海啸波与近岸的共振响应不仅与近岸复杂地形有关,而且与海啸波传播到近岸时的波浪入射方向有关;潘文亮等[8]通过分析2010年智利和2011年日本海啸在中国近岸传播的海啸波形,发现近岸海啸波相比于深海,能量衰减的更慢,持续时间更长,能量分布更集中;D a o和T k a l i c h[9]基于数值模拟,以2004年印度洋海啸事件为案例进行了关于天文潮㊁底摩擦和频散效应等因素的敏感性实验,研究中发现,海啸在浅水传播过程中底摩擦起到了非常重要的作用,海啸波与海岸的相互作用需要得到更多的关注;Y a m a z a k i等[10]利用有限断层模型很好地模拟了2011年日本东北地震海啸近场传播特征及陆架共振效应㊂以上研究对于认识海啸在近岸的传播特征和能量分布规律起到了推动作用㊂2017年9月8日墨西哥M w8.2级地震的震源位于墨西哥近海区域,地震发生后引发的海啸对墨西哥西海岸区域造成了灾害性的影响㊂通过对近岸验潮站水位数据的初步分析,发现此次海啸事件中部分验潮站的海啸波幅明显受到近岸地形的影响,需要进行进一步的分析和探讨㊂因此,研究分析此次墨西哥地震海啸事件,有助于进一步揭示海啸在近岸的传播规律以及近岸地区海啸波与地形的相互作用规律,具有一定的科学意义㊂本文基于近场深水D A R T(D e e p-o c e a nA s s e s s m e n t a n dR e p o r t i n g o fT s u n a m i s)浮标与近岸验潮站的观测数据,统计分析了海啸波在近岸的波浪特征,并采用小波变换进一步分析了海啸波的时频特征㊂同时采用M O S T海啸模型模拟重现了海啸在近场的传播演变过程㊂基于实测与模拟结果,详细分析了此次地震海啸的近场分布特征,初步解释了空间上邻近的区域对同一海啸事件的响应差异的原因㊂基于港湾共振理论,重点研究分析了港湾内海啸长期振荡特征及海啸波的放大效应,对可能的共振周期进行了估算㊂22017年9月8日墨西哥M w8.2级地震海啸事件概述2017年9月8日4时49分(U T C),墨西哥瓦哈卡州沿岸近海海域(15.21ʎN,93.64ʎW)发生M w8.2级地震,震源深度为30k m,震中位置距离墨西哥西海岸50k m,震后监测到余震700余次,最强余震达6级以上㊂根据历史资料统计,本次墨西哥地震是自1995年10月9日强震以及1985年9月21日地震后,墨西哥沿岸遭受的最严重的的地震,此次地震诱发了海啸㊂墨西哥M w8.2级地震发生在太平洋东部的可可斯板块(C o c o s)与北美板块交界的附近,可可斯板块以大约72~81m m/a的速度向东北方向俯冲,俯冲速度由南至北逐渐减小并俯冲至北美板块之下㊂俯冲的结果导致了该区域强烈的地震活动以及一系列活跃的火山链,该区域的中源地震大多发生在可可斯板块300k m深的位置㊂自1900年以来,该区域内发生过多次中源深度的中强震事件,包括1915年9月7日发生于萨尔瓦多的7.4级地震以及1950年10月5日发生于哥斯达黎加的7.5级地震㊂区域内历史最大地震事件为1932年6月3日发生在墨西哥的8.1级地震,震后当月在该区域又分别发生了1次7.7级和1次7.8级强震㊂据监测资料显示,墨西哥地震引发的海啸在震中附近产生了较为严重的影响,其中墨西哥H u a t u l c o 站在地震发生12分钟后最早监测到海啸波,墨西哥南部C h j a p a s站监测到最大海啸波高为1.75m㊂地震发生后的10余小时里,墨西哥㊁东南太平洋马克萨斯群岛以及萨摩亚群岛等地先后监测到海啸波,最大波高超过1m,最小仅有几厘米㊂从R a mír e z-H e r r e r a 等[11]对震后受灾最严重的恰帕斯州实地考察结果来看,本次海啸在恰帕斯州最大爬高为3m,最大淹没距离为190m㊂根据全球地震震源机制解目录(G C M T)和美国地质调查局(U S G S)公布的震源机制解(表1)可知,本次地震机制为典型的高角度正断层性质,两机构反演的震源机制解大致相同㊂反演结果如图1所示,地震破裂方向为北西向,破裂长度大约为200k m,宽约50k m㊂3数据处理及分析在本次海啸事件中,位于太平洋上的多个海啸监测网络捕捉到了海啸信号并详细记录了此次海啸的传播过程,为研究此次海啸的传播特征提供了珍贵的数据㊂D A R T浮标布放在深海,水位记录仪放置在海底平坦处,用于实时监测水位变化㊂本文选取了来自63海洋学报41卷美国国家海洋大气局国家浮标数据中心(N O A A ,N D B C )(h t t p ://w w w .n d b c .n o a a .go v /)的D A R T 浮标数据和联合国教科文组织政府间海洋学委员会(I O C -U N E S C O )网站(h t t p ://w w w .i o c -s e a l e v e l m o n i t o r i n g .o r g/)提供的近岸海洋台站的水位数据㊂常规状态下,D A R T 回传的观测数据采样间隔为15m i n ,发生海啸事件时回传频率提高为1m i n 或15s,本文中对数据进行差值处理,统一调整至1m i n 采样率;近岸海洋台站的采样周期为1m i n ㊂其中D A R T -43413站点记录的海啸最大波幅为7c m ,为D A R T 浮标记录的本次海啸最大波幅;墨西哥南部的C h j a p a s 站记录的最大海啸波幅为1.75m ,为本次海啸在近岸被记录到的最大波幅㊂图2为选取站点的地理信息㊂表1 地震断层节面参数T a b .1 N o d a l pl a n e p a r a m e t e r s o f s e i s m i c f a u l t s 序号机构节面Ⅰ节面Ⅱ走向/(ʎ)倾角/(ʎ)滑动角/(ʎ)走向/(ʎ)倾角/(ʎ)滑动角/(ʎ)M w矩心位置纬度经度深度/k m 1G C M T 14813-8332077-928.215.34ʎN 94.62ʎW 502U S G S16420-6131473-1008.215.02ʎN93.90ʎW47注:G C M T 为全球地震震源机制解目录,U S G S 为美国地质调查局㊂图1 震源位置及震源机制(五角星位置为震源区域)F i g .1 H y po c e n t e r l o c a t i o n a n d f o c a lm e c h a n i s m (ʀr e pr e n s t t h e s e i s m i c s o u r c e )3.1 海啸波提取在近岸站点观测到的水位时间序列中,除了海啸引起的水位变化,还有其他的波动成分㊂其中潮汐是主要的组成部分,除此之外还包括气象因素㊁过往船只等引起的水位扰动㊂为了在水位信号中更准确地分离出海啸波,本文首先根据对应站点9月份水位数据进行调和分析,分离出潮汐,然后对余下的信号采用巴特沃斯滤波器进行带通滤波,从而得到较为纯净的海啸波信号㊂图3为各站点的海啸波信号,表2为提取的海啸波信号的特征值㊂735期 孙立宁等:2017年9月8日墨西哥沿岸M w 8.2级地震海啸观测数据分析与模拟图2 水深分布及站点地理信息图F i g .2 B a t h y m e t r y an d l o c a t i o n o f s t a t i o ns 图3 海啸波信号F i g.3 T s u n a m iw a v e f o r m 83海洋学报 41卷从海啸波信号及其特征值中可以看出,在深海, D A R T-43413站点最早观测到了海啸的传播,海啸传播时间为1h8m i n,且该站点7c m的最大波幅为海啸在深海观测到的最大波幅;相比于D A R T-43413站,D A R T-32411站距离海啸源更远,在海啸传播1h 48m i n后观测到海啸波,最大波幅为3c m㊂在近岸区域,H u a t u l c o站最早观测到海啸波,海啸传播时间为24m i n;近岸观测到的最大波幅站点为C h i a p a s站,最大波幅为1.75m㊂从图3中可以看到,S a l i n a C r u z和H u a t u l c o两站海啸波震荡的持续时间要明显长于其他近岸站点,其原因在于这两个站点处于半封闭海湾地形中,海啸波传播到海湾中会受到地形的影响,产生反射或共振使得海啸波的持续时间变长㊂海湾地形与海啸波的相互作用同样对海啸波周期产生影响,H u a t u l c o站海啸波显著周期为12m i n,为近岸站点中的最小值㊂在近岸与海啸源距离较近的4个站点(L a z a r o C a r d e n a s㊁S a l i n aC r u z㊁H u a t u l c o㊁P u e r t oA n g e l),都是最早监测到负的水位信号,这与海啸初始场增水方向远离海岸线㊁减水方向面向海岸线的特征相一致,表明近岸海啸首波受到海啸源初始场的控制㊂在近岸站点中,S a l i n aC r u z与H u a t u l c o两站位于海湾中,这种特殊的海岸地形使得这两站具有相对较大的海啸最大波幅和较小的海啸波动周期㊂从地理空间分布来看,H u a t u l c o与P u e r t oA n g e l两站距离十分近(相距仅41k m),两站也几乎在同一时间观测到了海啸波的信号(前后相差3m i n),但是两站的观测信息在最大波幅和海啸波信号显著周期上都存在明显差别, H u a t u l c o站的最大波幅要明显大于P u e r t oA n g e l站,并且海啸波具有更短的周期㊂造成这种差异的原因与海啸本身关系不大,主要是由于两站所处位置的地形差异㊂H u a t u l c o站处于明显的海湾地形之中,这种特殊的地形使得海啸波在传播至海湾内时发生了反射㊁形成驻波,产生共振,使得海啸波幅变大㊂两站的差异也说明了海啸波在近岸的传播过程中,除了受到海啸波初始场的影响外,近岸地形特征起到了更为重要的作用㊂表2各站点海啸波主要信息T a b.2T h e i n f o r m a t i o n o f s t a t i o n s a n d t s u n a m i c h a r a c t e r i s t i c s站点位置到达时间(U T C)传播时间最大波幅/m显著周期/m i n D A R T-4341310.842ʎN,100.137ʎW5:581h8m i n0.0730,60 D A R T-324114.995ʎN,90.850ʎW6:381h48m i n0.0326L a z a r oC a r d e n a s17.940ʎN,102.178ʎW6:461h56m i n0.2060S a l i n aC r u z16.168ʎN,95.197ʎW6:001h10m i n0.9120,60 H u a t u l c o15.750ʎN,96.117ʎW5:1424m i n0.7110,20 P u e r t oA n g e l15.667ʎN,96.492ʎW5:1727m i n0.2220,60 A c a j u t l a13.574ʎN,89.838ʎW15:3710h47m i n0.5040L aL i b e r t a d13.485ʎN,89.319ʎW11:326h42m i n0.3774C h i a p a s14.697ʎN,92.411ʎW9:494h59m i n1.7532注:C h i a p a s站海啸波特征值为P T W C报告给出,I O C网站中无该站点水位数据㊂3.2小波分析为了进一步探究近岸地形特征在海啸波近岸传播过程中的影响,本文采用小波分析方法来分析近岸站点海啸波的时频特征㊂小波变换具有多分辨率分析的特点,并且在时频两域都具有表征信号局部特征的能力㊂小波变换通过将时间序列分解到时间频率域内,从而得到时间序列的显著的波动模式,即周期变化动态,以及周期变化动态的时间格局[12]㊂小波分析将信号分解成一系列小波函数的叠加,而这些小波都是由一个母小波函数经过平移与尺度伸缩得来的,用这种不规则的小波函数可以逼近那些非稳态信号中尖锐变化的部分,也可以去逼近离散不连续具有局部特性的信号,从而更为真实地反应原信号在某一时间尺度上的变化㊂本文采用935期孙立宁等:2017年9月8日墨西哥沿岸M w8.2级地震海啸观测数据分析与模拟M o r l e t小波变换对各站点的海啸波信号进行分析,从而进一步分析海啸波在近岸传播过程中的时频特征㊂小波变换分析结果如图4所示㊂可以看到, D A R T-43413站在6:00前后能量最强,对应海啸波首波的到达时间,在4:00 14:00之间能量主要集中在10~60m i n周期域内,14:00后高频部分能量减少㊂对应图3中该站点的海啸波信号发现,4:00 14:00为海啸波经过该站点的时间,说明此次海啸事件中海啸波在传播过程中能量集中在20~30m i n之间㊂与远海相比,近岸站点的海啸波主频波谱分布更广,高频部分能量持续时间更长㊂其中H u a t u l c o站在能量分布上与其他站点存在明显差别,近岸站点的能量峰值所在周期为60m i n,而H u a t u c l o站其能量峰值所在周期在20~30m i n之间,60m i n周期能量很少㊂图4海啸波小波功率谱和小波全谱图F i g.4 W a v e l e t a n a l y s i s o f t s u n a m iw a v e s根据M e r i a n提出的1/4波长共振理论[13]中的共振周期计算公式T=2L/n(g h)1/2(式中,L㊁n㊁g㊁H 分别表示海湾长度㊁模态数㊁重力加速度㊁平均水深),可以算出H u a t u l c o站的共振周期为21.04m i n(n=1时),可见海湾的共振周期与海啸波的主要能量周期相一致㊂综合前文对H u a t u l c o站的地形特点和海啸波特征值的分析可知,由于H u a t u l c o站所处的海湾地形特征,使得海啸波受地形影响激发了海湾固有周04海洋学报41卷期的震荡[14],而后相互作用,形成共振,从而大幅度提高了海啸的最大波幅,并延长了海啸波持续时间,对该地区造成了更为严重的灾害影响㊂3.3海啸波分窗口滤波分析在海啸波的小波分析中,我们发现海啸波的能量主要集中在10~60m i n之间,并且高频部分的能量衰减的更快㊂为了探究海啸波传播过程中能量和频率的关系,本文对各站点海啸波信号进行了不同窗口的滤波,分别为10~30m i n频段和60~80m i n频段,结果如图5所示㊂可以发现,相比于60~80m i n频段,10~30m i n频段海啸波的最大波幅更大,并且最大波幅与海啸波信号的最大波幅对应较为一致,由此可见海啸在传播过程中能量主要集中在10~30m i n 频段部分㊂10~30m i n频段携带的海啸波能量更多,但相比于60~80m i n频段其衰减速度更快㊂在深海区域(D A R T-43413)10~30m i n频段海啸波能量只维持几个小时,之后迅速衰减,相比而言60~80m i n频段海啸波的能量可以长时间维持㊂在近岸10~30 m i n频段海啸波的能量相比于远海可以维持更长的时间,但是衰减速度仍大于60~80m i n频段海啸波㊂在海湾地形内(H u a t u l c o和S a l i n aC r u z站),10~30 m i n频段海啸波能量可以维持更长时间㊂由此说明海啸波在传播初期主要能量集中在高频部分,之后迅速衰减,低频部分海啸波能量相对较小但是可以维持更长时间,近岸的特殊地形作用于海啸波可以使高频部分海啸波在特殊地形区域内持续震荡㊂图5分窗口带通滤波F i g.5 B a n d-p a s s f i l t e r i nd i f f e r e n t b a n d s 145期孙立宁等:2017年9月8日墨西哥沿岸M w8.2级地震海啸观测数据分析与模拟4地震海啸数值模拟海啸数值预报模型是探索海啸生成机制㊁评估海啸灾害和建立实时海啸预警系统的重要工具和技术手段㊂根据海啸波的物理特性,海啸数值模型通常有:计算海啸初始位移场的断层模型(G e n e r a t i o n);基于非线性浅水方程或者类B o u s s i n e s q方程计算海啸波传播的传播模型(P r o p a g a t i o n);考虑沿岸地形㊁地物㊁采用动态移动边界技术来实现海啸爬高及与海岸相互作用的淹没模型(I n u n d a t i o n)[4]㊂为了进一步研究海啸波在近岸的传播特征,本文采用数值模拟方法,通过将观测值与模拟值进行对比来探究近岸海啸波的传播特点㊂4.1海啸源模型海啸源模型是海啸数值模型研究的基础,传统上采用的方法是瞬时响应方法,即将断裂带滑移引起的海底地形变化直接作为海表面的初始扰动[15]㊂针对地震产生的海底地形变化,可采用弹性半空间位错理论计算㊂该理论首先由S t e k e t e e提出,并经过许多学者的深入研究,最后经过O k a d a的发展,得到在弹性半空间内断层滑动引起的位错及应变表达式,并在其后得到不断完善[1,16-17]㊂基于有限断层破裂模型,利用O k a d a弹性位错理论计算了海啸传播所需的初始扰动场,具体的参数采用表3中给出U S G S参数值(采用G C M T参数计算出的初始场与之相似,不再另行给出)㊂海啸的初始位移场如图6所示㊂从海啸源的形态特征我们可以看出,海面形变场能量均成具有西北向带状分布特征,海啸初始场具有显著的极性特征,最大抬升约为0.6m,最大水位下沉约为1.1m ㊂图6海啸初始位移场F i g.6S e a s u r f a c e d e f o r m a t i o n表3海啸初始场模型参数T a b.3T s u n a m i s o u r c e p a r a m e t e r s参数符号参数值(U S G S)参数值(G C M T)断层长度/k m L189.09189.09断层宽度/k m W70.2570.25倾角/(ʎ)δ7377滑动角/(ʎ)λ-100-92走向角/(ʎ)φ314320滑移量/m S l i p4.824.82深度/k m d45.550.0震源位置(x0,y0)(15.07ʎN,93.90ʎW)(15.34ʎN,94.62ʎW)4.2海啸传播模型M O S T(M e t h o do f S p l i t t i n g T s u n a m i s)模型采用非线性浅水方程(N o n l i n e a r S h a l l o w-W a t e rE q u a t i o n, N S W)来研究海啸波的演变过程㊂非线性方程包括两个空间变量和一个时间变量,即向岸方向空间变量x,沿岸方向空间变量y以及时间变量t,T i t o v与S y n o l a k i s采用分裂法将2+1维问题转化成两个1+1维问题,分别包含向岸方向变量x和时间变量t㊁沿岸方向变量y和时间变量t㊂描述海啸演变过程的2 +1维浅水方程具体如下:h t+(u h)x+(v h)y=0,u t+u u x+v u y+g h x=g d x,v t+u v x+v v y+g h y=g d y,ìîíïïïï(1)式中,h=η(x,y,t)+d(x,y,t);η(x,y,t)为自由水24海洋学报41卷面水位高度;d (x ,y ,t )为静水深;u (x ,y ,t )和v (x ,y ,t )分别为向岸方向和沿岸方向的深度平均速度;g 为重力加速度㊂M O S T 模型采用分裂法解2+1维浅水方程,将原本包含两个空间变量的方程分解成两个水动力系统,各自包含一个空间变量和时间变量t ㊂h t +(u h )x =0u t +u u x +g h x =gd x ,v t +u v x =0ìîíïïïï(2)h t +(v h )x =0v t +v v y +g h y =gd y .u t +v u y =0ìîíïïïï(3) M O S T 模型采用多重嵌套网格来进行海啸波传播模拟的计算,从而做到在深远海采用较粗的分辨率提高计算效率,在近岸区域采用更精细的网格保证模拟的准确性㊂在本次海啸的模拟中,本文采用了3层嵌套网格㊂第一层网格分辨率为1a r c m i n ,计算区域为2ʎ~22ʎN ,85ʎ~108ʎW (计算区域及水深分布如图2所示);第二层网格分辨率为0.5a r c m i n ,计算区域为13ʎ~19ʎN ,88ʎ~103ʎW ;第三层网格分辨率为3a r c s e c ,有3个计算区域,分别为C 1(15.5012ʎ~15.7112ʎN ,96.1052ʎ~96.6011ʎW )㊁C 2(14.4579ʎ~14.7804ʎN ,92.2503ʎ~92.5728ʎW )和C 3(16.0435ʎ~16.3485ʎN ,94.9221ʎ~95.2746ʎW )㊂前两层网格的地理信息数据来自G E B C O (h t t p s ://w w w .ge b c o .n e t/),第三层网格的地理信息数据来自美国国家海洋和大气管理局海啸研究中心(N O A A C e n t e r f o rT s u n a m iR e a s e a r c h)㊂4.3 模拟结果与观测数据的对比从最大海啸波幅分布情况来看(图7),本次海啸波的主要能量传播方向为远海向西南方向传播,近岸向东北方向传播㊂这种海啸能量分布特征主要受断层走向的影响㊂在近岸,最大海啸波幅主要分布在恰帕斯州沿岸,最大可达2m 左右,在瓦哈卡州的东南部也有0.6m 左右的最大波幅分布㊂恰帕斯州也是距离海啸源最近的区域,距离海啸源较远的区域最大海啸波幅相对较小,说明此次海啸事件影响区域主要为近海啸源区域,对远离海啸源区域的地方影响不大㊂图7为采用U S G S 参数模拟的最大海啸波幅分布场,采用G C M T 参数的模拟结果几乎一致,不再另行给出㊂图7 最大海啸波幅分布F i g .7 M a x i m u mt s u n a m iw a v e a m pl i t u d e 345期 孙立宁等:2017年9月8日墨西哥沿岸M w 8.2级地震海啸观测数据分析与模拟从观测值与模拟值的对比(图8)情况来看,海啸传播模型较好地模拟了海啸波在深海中的传播情况, D A R T-43413站的首波到达时间与波幅观测值与模拟值几乎完全一致,D A R T-32411站首波到达时间模拟值略早于观测值,并且模拟值的最大波幅略大于观测值,总体来看模拟结果在深海区域与观测值较为一致㊂相比而言,近岸的模拟结果较为复杂,L a z a r o C a r d e n a s站首波最大增水模拟值和观测值相同,模拟值首波到达时间稍早于观测值㊂S a l i n aC u r z站模拟值首波到达时间要比观测值晚20~30m i n,并且最大波幅明显大于观测值;H u a t u l c o站最大减水出现时间模拟值稍晚于观测值,并且浅水较小㊂造成S a l i n a C r u z和H u a t u l c o两站观测值与模拟值差距较大的原因在于两站处于海湾的地理位置,数值模拟采用的地形数据精度有限,无法刻画其真实的海湾地形,导致模拟值与观测值出现了一定的差距㊂从数值模拟的情况看,H u a t u l c o与P u e r t oA n g e l站的模拟值波形特征相差不大,没有实测值之间那种较大的差异,这也从侧面印证了数值模拟所采用的地形数据无法精确刻画实际地形特征㊂P u e r t oA n g e l站海啸最大波幅观测值为减水,模拟值的最大减水与观测相比出现稍晚,最大减水值相差不大㊂A c a j u t l a和L aL i b e r t a d由于距离海啸源较远,海啸波首波到达时模式已运行8h以上,导致模拟的海啸波信号与观测值有相对较大的差别,说明海啸传播模型在长时间运行后其误差也会累加,导致模拟结果变差㊂图8海啸波观测数据与模拟数据对比(C h i a p a s站无观测值)F i g.8 C o m p a r i s o n s o f o b s e r v e d a n d s i m u l a t e d t i m e s e r i e s(n o o b s e r v e dd a t a i nC h i a p a s)5总结本文采用D A R T和近岸站点水位数据,从原始水位信息中通过潮汐调和分析和滤波的手段分离出海啸波,并采用小波分析和分窗口滤波方法进一步分析海啸波的波频和能量特征㊂同时基于O k a d a弹性44海洋学报41卷位错理论断层模型计算得到了地震引发的海底形变分布,并采用M O S T 海啸传播模型对本次海啸事件进行模拟㊂通过对本次海啸事件的分析,得出以下结论:(1)本次海啸主要对墨西哥沿岸造成灾害性影响,对远洋区域影响不大㊂由于震中位置距离墨西哥海岸较近,导致海啸波在地震发生后短时间内就到达了墨西哥沿岸,沿岸监测到的最大波高为1.75m ㊂海啸的初始水位为减水面向海岸,因此在距离海啸源较近的近岸站点都是首先观测到负的水位变化㊂(2)海啸除受海啸源初始场控制外,在近岸的传播过程中也会受到地形的影响,并且近岸地形在海啸波能量分布上往往起到更大的作用㊂在本次地震海啸事件中,S a l i n a C r u z 和H u a t u l c o 两站处于海湾地形中,导致两站产生了较大的最大海啸波幅,尤其是H u a t u l c o 站海湾地形的共振周期与海啸占优周期相一致,导致其最大海啸波幅明显大于与之地理位置十分接近的P u e r t oA n g e l 站㊂海湾地形与海啸波相互作用的结果也使得海啸波的持续时间变长㊂这种特殊地形的作用,可能会加重该区域海啸的灾害性,造成更严重的损失,在之后的研究和工作过程中要警惕㊂(3)海啸波主要能量集中在10~30m i n 频段部分,在该频段衰减迅速,60~80m i n 频段海啸波能量虽然相对较小但是可以维持更长时间㊂在本次海啸事件中,H u a t u l c o 站与当地的海湾地形发生共振效应,海啸波地形的相互作用使得高频海啸波能量在该区域内长期维持㊂(4)海啸传播模型在深海地区可以得到较好的模拟效果㊂而在近岸地区,一方面由于近岸区域地形更为精细和复杂,地形数据无法准确刻画;另一方面由于近岸区域地形与海啸波存在着更为复杂的相互作用,导致模拟结果相对较差㊂对于拥有较复杂地形的区域,需要有更高分辨率的地形数据作为基础,同时需要进一步考虑近岸区域海啸波与地形间的复杂相互作用,才能取得较好的模拟效果㊂参考文献:[1] 任智源.南海海啸数值模拟研究[D ].上海:上海交通大学,2015.R e nZ h i y u a n .N u m e r i c a l s i m u l a t i o no f t s u n a m i i nS o u t hC h i n a S e a [D ].S h a n g h a i :S h a n g h a i J i a oT o n g U n i v e r s i t y,2015.[2] 任智源,原野,赵联大,等.2016年全球地震海啸监测预警与数值模拟研究[J ].海洋科学,2017,41(6):98-110.R e nZ h i y u a n ,Y u a nY e ,Z h a oL i a n d a ,e t a l .M o n i t o r i n g ,e a r l y w a r n i n g a n d n u m e r i c a l s t u d y o f gl o b a l t s u n a m i s i n 2016[J ].M a r i n e S c i e n c e s ,2017,41(6):98-110.[3] 闪迪,王培涛,任智源,等.有限断层模型在2015年9月16日智利M w 8.3级地震海啸数值模拟中的应用与评估[J ].海洋学报,2017,39(11):49-60.S h a nD i ,W a n g P e i t a o ,R e nZ h i y u a n ,e t a l .A p p l i c a t i o n a n d e v a l u a t i o n o f t h e 16S e p t e m b e r 2015I l l a p e l ,C h i l e M w 8.3e a r t h q u a k e f i n i t e f a u l t r u p-t u r em o d e l f r o mn u m e r i c a l s i m u l a t i o n [J ].H a i y a n g Xu e b a o ,2017,39(11):49-60.[4] 王培涛,于福江,赵联大,等.2011年3月11日日本地震海啸越洋传播及对中国影响的数值分析[J ].地球物理学报,2012,55(9):3088-3096.W a n g P e i t a o ,Y uF u j i a n g ,Z h a oL i a n d a ,e t a l .N u m e r i c a l a n a l y s i s o f t s u n a m i p r o p a g a t i n g g e n e r a t e d b y t h e j a p a n M W 9.0e a r t h qu a k e o nM a r .11i n 2011a n d i t s i m p a c t o nC h i n a c o a s t s [J ].C h i n e s e J o u r n a l o fG e o p h ys i c s ,2012,55(9):3088-3096.[5] 于福江,原野,赵联大,等.2010年2月27日智利8.8级地震海啸对我国影响分析[J ].科学通报,2011,56(3):239-246.Y uF u j i a n g ,Y u a nY e ,Z h a o L i a n d a ,e t a l .E v a l u a t i o n o f p o t e n t i a l h a z a r d s f r o mt e l e t s u n a m i i nC h i n a :T i d a l o b s e r v a t i o n s o f a t e l e t s u n a m i ge n e r a t e d b y t h eC h i l e 8.8M w e a r t h qu a k e [J ].C h i n e s e S c i e n c eB u l l e t i n ,2011,56(11):1108-1116.[6] O m i r aR ,B a p t i s t aM A ,L i s b o aF .T s u n a m ic h a r a c t e r i s t i c s a l o n g t h eP e r u -C h i l e t r e n c h :a n a l y s i s o f t h e 2015M W 8.3I l l a p e l ,t h e 2014M W 8.2I q u i q u e a n d t h e 2010M W 8.8M a u l eT s u n a m i s i n t h eN e a r -f i e l d [J ].P u r e a n dA p p l i e dG e o p h y s i c s ,2016,173(4):1063-1077.[7] 张金凤.海啸波近岸共振响应的数值模拟及分析[J ].防灾减灾工程学报,2010,30(2):147-151.Z h a n g J i n f e n g .N u m e r i c a l s i m u l a t i o n a n d a n a l y s i s o f o f f s h o r e r e s o n a n c e o f t s u n a m iw a v e s [J ].J o u r n a l o fD i s a s t e rP r e v e n t i o na n d M i t i g a t i o nE n gi -n e e r i n g,2010,30(2):147-151.[8] 潘文亮,王盛安,孙璐,等.2010年智利和2011年日本海啸在华南沿岸的实测海啸波形和特征[J ].热带海洋学报,2014,33(6):17-23.P a n W e n l i a n g ,W a n g S h e n g a n ,S u nL u ,e t a l .O b s e r v e dw a v e f o r ma n d c h a r a c t e r i s t i c s o f t h e 2010C h i l e a n d 2011J a pa n t s u n a m i s n e a r t h e c o a s t o f S o u t hC h i n a [J ].J o u r n a l o fT r o p i c a lO c e a n o g r a p h y,2014,33(6):17-23.[9] D a oM H ,T k a l i c hP .T s u n a m i p r o p a g a t i o nm o d e l l i n g -a s e n s i t i v i t y s t u d y [J ].N a t u r a lH a z a r d s a n dE a r t hS y s t e mS c i e n c e s ,2007,7(6):741-754.[10] Y a m a z a k i Y ,C h e u n g K ,L a y T .M o d e l i n g o f t h e 2011t o h o k un e a r -f i e l d t s u n a m i f r o mf i n i t e -f a u l t i n v e r s i o no f s e i s m i cw a v e s [J ].B u l l e t i no f t h e S e i s m o l o g i c a l S o c i e t y o fA m e r i c a ,2013,103(2B ):1444-1455.[11] R a m ír e z -H e r r e r aM T ,C o r o n aN ,R u i z -A n g u l oA ,e t a l .T h e 8S e p t e m b e r 2017T s u n a m iT r i g g e r e db y t h e M W 8.2I n t r a p l a t e E a r t h qu a k e ,C h i a -545期 孙立宁等:2017年9月8日墨西哥沿岸M w 8.2级地震海啸观测数据分析与模拟。
广义塑性Cosserat连续体模型及其应用
摘要岩土材料在自然条件下普遍是非均匀的，这种非均匀性对土体的应力应变关系、承载力与变形特性、剪切带的发展有着显著影响。
当前的主流理论仍然假设土体在空间中是均匀的，并且很少有本构理论可以同时考虑静动力问题。
另外颗粒破碎现象在土体受荷加载过程中也是十分常见的，对于土体应力应变关系的影响不可忽视。
另外，在土石坝等实际工程中，大颗粒材料如堆石料等也经常由于高围压下的颗粒破碎而导致应力应变关系变化，进而影响大坝整体的沉降变形。
同时有研究表明，颗粒破碎是地震荷载影响坝体变形的重要因素。
为了发展一种能有效进行静动力分析并考虑颗粒破碎的本构理论，本文做了如下工作：（1）将Pastor，Zienkiewicz等人提出的广义塑性理论以及后人对这一理论的完善进行ABAQUS软件中的本构程序二次开发，这一程序被用于常规三轴压缩试验、动三轴试验以及土石坝的填筑与地震分析模拟中。
通过与试验结果、工程实例结果对比，得到了较好的数值模拟结果，验证了广义塑性模型的正确性，证明了颗粒破碎对应力应变关系以及承载力位移关系有着重要影响，也验证了广义塑性理论可以在同时考虑静动力荷载的条件下保持模型参数的一致性。
通过土石坝静动力条件下的分析，进一步揭示了颗粒破碎对沉降变形的影响，同时从填筑层数的角度说明了沉降与分层数量的反相关关系。
（2）通过模拟Borja等人的平面应变试验，将考虑了空间非均匀性分布的土体，与常规的均匀性土体同时用广义塑性模型进行平面应变模拟，进一步验证广义塑性模型的正确性。
通过对比揭示了非均匀性对于剪切带发展规律以及承载力软化的影响，指出通常认知中的X型剪切带在非均匀性作用下的发展规律；通过使用不同类型的单元进行平面应变模拟，得到单元类型对于非均匀性土体剪切带发展的影响；经过多重比较，确定了C3D20R单元在模拟非均匀性土体剪切带时的优势；通过对模型参数的校正，揭示了颗粒破碎对剪切带数量的影响，也揭示了关键参数对剪切带的控制作用。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
261Numerical simulation of dynamic process of anearthquake
Yongen Cai and Ren WangGeodynamic Research Center of Peking University Department of Geology, PekingUniversity, Beijing,100871, China (e-mail: yonen@pku.edu.cn; phone: 86 10 62755441).
AbstractIn this investigation, we use a new numerical method LDDA ( Lagrangian DiscontinuousDeformation Analysis)[1-2] to simulate a dynamic process of an earthquake occurring onstrike slip fault, The method can be used to solve the dynamic problem with discontinuousdeformation block system, it grows up on the base of the finite element method with do-main decomposition algorithm [3-4] and DDA (Discontinuous Deformation Analysis) [5]and overcomes the weakness of the classical Lagrangian multiplier method by using do-main decomposition algorithm and obviously prevails over the penalty method in the satis-faction of constrain conditions. The process of solving the problem by LDDA is as fol-lows:1) find total contact points (Lagrangian multiplier points) between blocks according tothe contact criteria of the DDA,2) solve the forces on the contact points between blocks by domain decompositionmethod.3) calculate the displacement and stresses caused by the contact forces in term of theFEM for each block.In mechanics, the earthquake occurrence can be considered as the results when theshear stress on the fault overcomes the frictional resistance of the fault. The LDDAmethod will be used to simulate the dynamic process by changing the frictional state of thefault[6] in the following example. The geometry of the model used in the study is shown inFig. 1. The fault is simulated by contact element across which the displacement is discon-tinuous when slip occurs. The initial stress field is induced by applying displacements onthe boundaries of the model. Under the state the frictional coefficient of the middle section(thick line) is changed from 1.0 to 0.1, it starts the dynamic process of fault sliting. Thecomputed results of the dynamic process are shown in Figs. 2 and 3. Fig. 2 shows the“Seismographs” of the x-, y- components of the displacement at “station” A, B, C and D. Itshows that the “elastic rebound” is about 7cm at A at time 0.03 sec. after the start of slid-ing, but the residual displacement found at time 5.0 sec.(quasi-static state) is less than thismagnitude. At other locations the displacements are also not recovered to their initial dis-placement. Fig. 3 shows the contours of the change in displacement and maximum shearstress caused by earthquake at different time.262
0.5m0.5m1m1m700m700m70m25m25mµ=1.0µ=0.1µ=0.1D+A+B++C
Fig. 1 Model geometry and boundary conditionsWe see from the results that 1) the relative displacement of the fault obtained after thissimulated earthquake is less than the maximum reached at the initial moment of the dy-namic process, 2) The maximum shear stress drop obtained is about 38Mpa at the fault anddecreases away from the fault, 3) P and S wave can not be identified at the earthquakesource region.Although the material model and geometry shapes in the examples are very simple theresults still reveal some very interesting phenomena and are helpful for understanding thedynamic process of earthquake events, especially for the discontinuous deformation dy-namic problems.The method will be used to simulate the dynamic process of earthquake migration in afault network263
(((((((((7LPHVHF
'LVSODFHPHQW;P1RGH$(((((((
(
(
7LPHVHF7LPHVHF
'LVSODFHPHQW<P1RGH$
(((((((7LPH7LPHVHFVHF
'LVSODFHPHQW;P1RGH%(((((
(
7LPHVVHHFF
'LVSODFHPHQW<P1RGH%
(((((((7LPHPHVHVHFF
'LVSODFHPHQW;P1RGH&((((((