四川地区地壳上地幔速度结构的初步研究

赵珠, 张润生

赵珠, 张润生. 1987: 四川地区地壳上地幔速度结构的初步研究. 地震学报, 9(2): 154-166.
引用本文: 赵珠, 张润生. 1987: 四川地区地壳上地幔速度结构的初步研究. 地震学报, 9(2): 154-166.
ZHAO ZHU, ZHANG RUNSHENG h. 1987: PRIMARY STUDY OF CRUSTAL AND UPPER MANTLE VELOCITY STRUCTURE OF SICHUAN PROVINCE. Acta Seismologica Sinica, 9(2): 154-166.
Citation: ZHAO ZHU, ZHANG RUNSHENG h. 1987: PRIMARY STUDY OF CRUSTAL AND UPPER MANTLE VELOCITY STRUCTURE OF SICHUAN PROVINCE. Acta Seismologica Sinica, 9(2): 154-166.

四川地区地壳上地幔速度结构的初步研究

PRIMARY STUDY OF CRUSTAL AND UPPER MANTLE VELOCITY STRUCTURE OF SICHUAN PROVINCE

  • 摘要: 本文使用10个工业爆破和154个天然地震,以及四川台网50个台站记录的 P 波组到时等资料,以龙门山断裂和二次大地构造单元分界线为界(图1),得出四川东部盆地和西部高原不同的地壳上地幔平均速度模型.若简单地采用双层地壳模型,则东部地壳厚40-41km,壳下 P 波速度为8.15-8.2km/s,壳内上层平均速度为5.82-5.9km/s,厚18km,下层平均速度为6.47-6.54km/s,厚22-23km;西部地壳厚61-64km,壳下 P 波速度为7.8-7.84km/s,壳内上层平均速度为5.82-5.98km/s,厚27-28.5km,下层平均速度为6.94-7.0km/s,厚34-35.5km.此模型为四川地区走时表提供了依据,也为研究地壳上地幔结构与地震的关系,研究我国大陆地块的构造演化及形成等,提供了有用的约束.
    Abstract: Arrival-time data from 10 industrical explosions and 154 earthquakes recorded by 50 stations of the Sichuan Seismic Network are analyezd to formulate the different average curstal and mantle velocity models for the east basin and west tableland of Sichuan divided by the Longmenshan fault and secondary boundaries of land structure blocks. If a simple two-layer crust is assumed, we have obtained: (1) the velocity model in the eastern part has a crustal thickness of 40-41km, an upper mantle velocity (P wave) of 8.15-8.2 km/s, an 18 km thick upper crust with an average velocity of 5.82-5.9 km/s, and a 22-23 km thick lower crust with an average velocity of 6.47-6.54 km/s; (2) the velocity model in the western part has a crustal thickness of 61-64km, an upper mantle velocity (P wave) of 7.8-7.84 km/s, a 27-28.5 km thick upper crust with an average velocity of 5.82-5.96 km/s and a 34-35.5 km thick lower crust with an average velocity of 6.94-7.0 km/s. This model provides the hasis for the regional travel-time tables of Sichuan and useful constraints for researches the relationship of the crustal and mantle structure with earthquakes, the structural evolution and formation of the China land blocks and so on.
  • 岩土及地下结构的静动力耦合模拟直接影响着其地震反应的合理性。在考虑土体重力引起的初始应力和地震动力作用的静动力耦合地震动力反应分析中,静力计算采用的是静力边界,如固定边界、滚轴边界等;有别于静力分析,动力计算使用的是动力边界,如透射边界(廖振鹏等,1984)、无限元边界(Ungless,1973)、黏弹性边界(Deeks,Randolph,1994)等。一般的动力边界不能用于静力分析,于是数值模拟计算过程中出现了边界不统一的矛盾。以黏弹性边界为例,若在静力分析中使用常规黏弹性动力边界,会影响后续动力反应(郭亚然等,2018)。为解决这一问题,目前发展了三类方法。第一类是实现静动力耦合分析时边界条件的统一性,将动力边界用于静力分析中。此类方法利用了黏弹性边界良好的鲁棒性,通过修正弹性系数Kb来发展静动力统一边界。刘晶波和李彬(2005)提出了一种三维黏弹性静动力统一人工边界以解决静动力耦合动力分析中静动力边界不统一的问题,但此种边界在静力计算时模拟结果存在一定误差(高峰,赵冯兵,2011)。廖维张等(2006)提出了一种开放系统非线性动力反应分析算法,在进行地下结构非线性动力地震反应分析时,将无零频限制的黏弹性边界用于静力分析中,以此实现静动力边界条件的统一。第二类方法是采用动力学解法计算静力反应(杜修力等,2000)。此方法将静力荷载视为零荷载时刻所施加的阶跃荷载,以阶跃函数的形式施加在结构上,用动力方法进行计算,由于模型阻尼的作用,阶跃荷载引起的动力反应会很快衰减并处于平衡状态,然后再进行动力计算分析。第三类处理静动力边界矛盾的传统方法是静动力边界转换方法。此方法在静力计算完成后,将静力边界转换成动力边界,再进行动力计算。具体来讲,先用固定边界计算重力引起的单元应力和边界节点反力,以此作为动力分析的初始条件,最后进行地震动力分析。

    使用这种传统的静动力边界转换方法来进行地震动力反应数值模拟比较普遍,但是静动力边界转换的处理细节不尽相同(高峰,赵冯兵,2011王苏等,2012付培帅,2015郭亚然等,2018),且转换过程中静力边界转换成动力边界时模型会出现振荡反应(Su et al,2021)。Su等(2021)提出利用ABAQUS中生死单元功能来杀死和激活黏弹性边界(弹簧和阻尼单元)以实现静动力边界的切换,但未解释振荡的原因和机理。王苏等(2012)给出了静动力边界转换的分析思路。高峰和赵冯兵(2011)针对静动力边界转换方法开展了相应的工作,提出了一种验证静动力人工边界转换合理性的方法,但不能直观地反应出静动力边界转换结果对后续动力计算的影响。

    静动力边界转换过程较为繁琐,涉及到应力的提取和施加,静动力边界的切换以及地震动的输入,人为操作过程中容易出错,出现振荡反应,因此有必要在进行动力分析前对静动力边界转换的合理性进行检验。为此,本文将围绕静动力边界转换中出现的振荡问题及转换合理性验证方法开展相关研究工作:根据岩土材料的特性提出参考解法,并通过算例对该方法的适应性进行探讨,以验证静动力边界转换的合理性;通过对算例的分析揭示振荡的原因和机理,以加深对边界转换过程中振荡的理解。

    黏弹性边界是众多动力边界中的一种,其本质是为了保证计算区域内反射波和散射波在边界上不发生反射,方法是在截断的模型边界上施加连续的并联弹簧和阻尼器。弹簧刚度系数和阻尼系数由下式获得(刘晶波等,2007):

    $$ \left\{ \begin{array}{l} {K_{{\rm{BT}}}} = {\alpha _{\rm{T}}}\dfrac{G}{R},{C_{{\rm{BT}}}} = \rho {C_{\rm{S}}} \text{,} \\ {K_{{\rm{BN}}}} = {\alpha _{\rm{N}}}\dfrac{G}{R},{C_{{\rm{BN}}}} = \rho {C_{\rm{P}}} \text{,} \end{array} \right. $$ (1)

    式中:$ K_{{\rm{BT}}} $和$ K_{{\rm{BN}}} $分别为切向弹簧和法向弹簧的刚度系数;$ C_{{\rm{BT}}} $和$ C_{{\rm{BN}}} $分别为切向阻尼器和法向阻尼器的阻尼系数;$ \alpha_{{\rm{T}}} $和$ {\alpha}_{{\rm{N}}} $为黏弹性边界的修正系数,二维情况下,$ \alpha_{{\rm{T}}} $和$ {\alpha}_{{\rm{N}}} $的建议取值为$ \alpha_{{\rm{T}}}=0.5 $,$\alpha_{{\rm{N}}}=1$ (刘晶波等,2007);R为波源至人工边界点的距离;$\; \rho $和G分别为介质的密度和剪切模量;$ C_{{\rm{S}}} $和$ C_{{\rm{P}}} $分别为介质的剪切波速和纵波波速。各边界节点地震动等效节点力计算公式如式(2)—(5)所示。本文在马笙杰等(2020)和马笙杰(2021)文献基础上,进一步编制了适用于成层介质的等效节点力计算及黏弹性边界施加的程序VBEA2.0(program for viscous-spring boundary and equivalent nodal load application 2.0),以实现地震动在黏弹性边界上的垂直入射。等效节点力Fb的计算方法为

    $$ {F_{{\rm{b}}}}= ( {K_{{\rm{}}b}} {{\boldsymbol{u}_{{\rm{b}}}}}+{{C_{{\rm{b}}}} {{\dot{\boldsymbol{u}}}}_{{\rm{b}}}}+{{\sigma_{{\rm{b}}}} {\boldsymbol{n}}} ) { A_{{\rm{b}}} }\text{,} $$ (2)

    底边界,

    $$ \left\{\begin{aligned} F_{{\rm{b}} x}^{-y}= & A_{{\rm{b}}}\left[K_{\mathrm{BT}} {\boldsymbol{u}} ( t ) +C_{\mathrm{BT}} \dot{{\boldsymbol{u}}} ( t ) - \frac{G}{{\rm{d}} y} ( {\boldsymbol{u}}_{x}^{i}-{\boldsymbol{u}}_{x}^{j} ) -\frac{1}{2} \rho {\rm{d}} y \ddot{{\boldsymbol{u}}}_{x}^{j}\right] \\ F_{{\rm{b}} y}^{-y}= & 0\text{,} \end{aligned}\right. \text{,} $$ (3)

    左边界,

    $$ \left\{\begin{array}{l} F_{{{\rm{b}}x}}^{-x}=A_{{\rm{b}}} [ K_{\mathrm{BN}} {\boldsymbol{u}} ( t ) +C_{\mathrm{BN}} \dot{{\boldsymbol{u}}} ( t ) ] \text{,} \\ F_{{\rm{b}} y}^{-x}=-A_{{\rm{b}}}\left[\dfrac{G}{{\rm{d}} y} ( {\boldsymbol{u}}_{x}^{i}-{\boldsymbol{u}}_{x}^{j} ) +\dfrac{1}{2} \rho {\rm{d}} y \ddot{{\boldsymbol{u}}}_{x}^{j}\right]\text{,} \end{array}\right. $$ (4)

    右边界,

    $$ \left\{\begin{array}{l} F_{{{\rm{b}}x}}^{x}=A_{{\rm{b}}} [ K_{\mathrm{BN}} {\boldsymbol{u}} ( t ) +C_{\mathrm{BN}} \dot{{\boldsymbol{u}}} ( t ) ] \text{,} \\ F_{{\rm{b}} y}^{x}=A_{{\rm{b}}}\left[\dfrac{G}{{\rm{d}} y} ( {\boldsymbol{u}}_{x}^{i}-{\boldsymbol{u}}_{x}^{j} ) +\dfrac{1}{2} \rho {\rm{d}} y \ddot {\boldsymbol{u}}_{x}^{j}\right]\text{,} \end{array}\right. $$ (5)

    式中:$K_{\rm{b}} {\boldsymbol{u}}_{\rm{b}}$表示为了克服位移引起弹簧单元所产生的附加应力;$C_{\rm{b}} \dot{{\boldsymbol{u}}}_{\rm{b}}$表示为了克服速度引起阻尼器产生的附加应力;$\sigma_{\rm{b}} {\boldsymbol{n}}$为自由场地震动在边界产生的应力张量;${\boldsymbol{u}}_{\rm{b}}= [ \begin{array}{lll}u & v & w\end{array} ] ^{\mathrm{T}}$为入射波位移向量;$\dot {\boldsymbol{u}}_{\rm{b}}= [ \begin{array}{lll}\dot u & \dot v & \dot w\end{array} ] ^{\mathrm{T}}$为入射波速度向量;Ab为边界节点控制面积,二维模型中为节点控制边长;n为边界外法线方向余弦向量;二维情况下,边界节点入射波位移向量为${\boldsymbol{u}}_{\rm{b}}= [ u \;\;v ] ^{\mathrm{T}}$,入射波速度向量为$\dot {\boldsymbol{u}}_{\rm{b}}= [ \dot u \;\;\dot v ] ^{\mathrm{T}}$;$K_{\rm{b}}$,$ C_{\rm{b}} $分别为黏弹性边界上弹簧刚度系数和阻尼系数。

    本文采用如下静动力边界转换方法来进行静动力耦合计算:先用固定边界计算由重力引起的单元应力和边界节点反力,以此作为动力分析的初始条件;然后进行地震动力分析,即先进行地应力平衡计算,再进行静动力边界转换,最后进行地震动力计算。具体过程如图1所示:

    图  1  静动力边界转换方法示意图
    (a) 地应力平衡;(b) 静动力边界转化;(c) 施加地震荷载
    Figure  1.  Schematic diagram of static-dynamic boundary switch
    (a) In-situ stress balance;(b) Switching between static and dynamic boundary;(c) Application of earthquake loading

    1) 首先约束模型侧边界水平方向、底部水平及竖向位移,施加重力,通过静力分析进行地应力平衡,使模型的初始位移为零,一般10−5即满足要求。提取侧边界节点的水平反力$ R F_{x}^{1} $、底部节点的水平反力$ R F_{x}^{2} $和竖向反力$ R F_{y}^{2} $以及模型所有单元的应力$ \sigma $。其中反力RF (reaction force)的下标表示力的方向,上标1表示侧边界,上标2表示底边界。对于二维问题,$ \sigma $为在ABAQUS中提取的积分点(centroid)应力分量S11,S22,S33,S12,对于三维问题,还应提取S13和S23。

    2) 对地下结构地震反应分析模型施加黏弹性边界,取消1)中静力边界条件并施加地应力平衡后获得的边界节点反力、所有单元应力和重力,以此实现静动力边界的转换。

    3) 在前面静力分析步的基础上创建隐式动力分析步,并在该分析步施加地震等效节点力和重力荷载。

    由2.1节可知,静动力转换过程较为繁琐,需要人为操作的步骤较多,因此在进行动力分析前有必要进行边界转换成功与否的验证,以保证后续动力分析计算的正确性。高峰和赵冯兵(2011)提出施加零动力荷载的方法来检验静动力人工边界转换误差。零动力荷载法能观察到地震动力反应时程曲线(位移、应力等)的初期是否有波动,但不能直观反映此波动的相对大小以及对后续地震动力反应的影响程度。本文中的参考解法能够直接得到静动力边界转换后与参考解的误差数值,为合理判断静动力边界合理转换提供一种有效的方法。

    根据弹性叠加原理,现给出一种判断静动力边界转换合理性方法:参考解法。由该原理可知,在线弹性无阻尼条件下,考虑静动力耦合与不考虑静动力耦合的自由场场地地震动力计算模型在各深度处的地震动力反应(加速度时程、位移时程等)是一致的。因此可以把不考虑静动力耦合的模型计算结果作为动力反应分析的参考解来判断静动力边界转换的正确性,这样就能直接得出静动力耦合计算模型在静边界转换后计算结果与理论值之间的误差,可以确保后续非线性静动力耦合模拟的合理性和正确性。具体方法如下:先建立不施加重力引起的初始应力的自由场弹性计算模型(不考虑静动力耦合作用的计算模型),施加地震动力荷载,进行动力分析;另外建立一个施加重力引起的初始应力并进行静动力边界转换的弹性计算模型(考虑静动力耦合作用的计算模型),然后施加同样的地震动力荷载,进行动力分析。通过比较两模型在各深度处地震动力反应(位移、加速度时程),即可判断考虑静动力耦合模型的静动力边界转换是否合理,若两者同一模型位置处的位移、加速度时程等动力反应一致,则静动力边界转换合理;反之,不合理。在具体操作过程中,建议选取侧边界底部、顶部节点以及模型中部、顶部和底部的节点的加速度动力响应来比较验证,具体原因将在第4节讨论。

    建立两类不同土层厚度的成层自由场模型(100 m×50 m和200 m×100 m),模型都均分为5层土层且两个模型都使用图2中的材料参数。不考虑土体的阻尼,采用弹性本构模型来描述土层的地震动力反应。有限元网格尺寸只有足够小,高频的地震波才能在离散体系中传播,进而保证计算精度。要获得波长为$ \lambda_{T} $波动的可靠数值结果,网格尺寸$ \Delta x $应当满足式(廖振鹏,刘晶波,1986):$\Delta x {\text{≤}}$(1/10—1/8)λT,本研究中采用0.5 m的四节点应变单元(CPE4)对模型进行离散。两种不同尺寸的模型分别考虑施加重力引起的初始应力(考虑静动力耦合作用)和不考虑施加初始应力(不考虑静动力耦合作用),一共4种不同计算工况。使用适用于成层地基的黏弹性边界施加和地震输入程序VBEA2.0,在模型两侧、底部施加黏弹性边界和等效节点力,将计算模型考虑为柔性地基,即:在地震波输入后,计算区域内的反射波和散射波不会在边界上发生反射。输入的地震动为神户海洋气象台处获得的地震南北分量前20 s地震动记录作为模型的输入,如图3所示。

    图  2  二维成层半空间模型
    Figure  2.  Two-dimensional layered half-space model
    图  3  日本神户地震南北分量加速度时程
    Figure  3.  Acceleration time history of the north and south components of the Kobe earthquake

    在不同工况下两种不同尺寸模型的中部沿深度方向的地震反应计算结果如图4所示。由图4a4b可知,50 m高的多土层自由场模型在地震作用下,考虑静动力耦合作用和不考虑静动力耦合作用两种工况下,随深度变化的峰值加速度和最大相对位移曲线吻合得很好。对于100 m高的多土层自由场模型,土层各深度节点峰值加速度以及最大相对位移的趋势与50 m高的自由场模型相似,考虑静动力耦合作用和不考虑静动力耦合工况的加速度、位移动力计算结果是一致的,也吻合得很好,如图4c4d所示。

    图  4  自由场地震反应计算结果
    (a) 100 m×50 m模型峰值加速度;(b) 100 m×50 m模型最大相对位移;(c) 200 m×100 m模型峰值加速度;(d) 200 m×100 m模型最大相对位移;(e) 随深度变化的加速度、位移时程误差
    Figure  4.  Seismic dynamic response of free field
    (a) Maximum acceleration for 100 m×50 m model;(b) Maximum relative displacement for 100 m×50 m model;(c) Maximum acceleration for 200 m×100 m model;(d) Maximum relative displacement for 200 m×100 m model;(e) Acceleration and displacement time-history errors with varying depths

    为进一步确定静动力耦合作用对自由场地震动力反应的影响,采用二阶欧几里得范数来量化和分析自由场模型中部沿深度方向节点的加速度时程误差值,如式(6)所示。

    $$ R=\frac{\|X-Y\|^{2}}{\|Y\|^{2}}=\frac{\sqrt{\displaystyle\sum_{i=1}^{n}\left|x_{i}-y_{i}\right|^{2}}}{\sqrt{\displaystyle\sum_{i=1}^{n}\left|y_{i}\right|^{2}}} {\text{×}} 100 \text{%} \text{,} $$ (6)

    式中,X为未考虑静动力耦合作用的加速度时程计算结果,Y表示考虑静动力耦合作用的加速度时程计算结果,将其作为参考解。图4e给出了两种不同尺寸模型不考虑静动力耦合作用时,误差随深度变化的曲线,可以看出,四个算例模型中部沿深度方向节点的加速度和位移时程误差都接近于零。

    综合上述结果,在黏弹性边界下,弹性模型考虑静动力耦合作用与否,不影响模型各深度的峰值加速度、相对最大位移的动力计算结果。因此,在弹性自由场静动力耦合计算模型中,把不考虑静动力耦合作用的计算模型的动力反应结果作为参考解来判断静动力边界转换的合理性和正确性是可行的。

    本节以日本大开地铁车站及其场地条件(曹炳政等,2002)作为研究对象,并以其为基准,通过调整土层剪切波速来获得另外较软地层和较硬地层场地,以此建立三种不同坚硬程度土层的弹性土-地下结构静动力耦合计算模型,场地条件如表1所示。车站主体为C30混凝土,将其考虑为线弹性,弹性模量为$ 3 {\text{×}} 10^{4}\; \mathrm{MPa} $,泊松比为0.18,密度为2 500 kg/m3,采用等效前、后刚度不变原则,把中柱弹性模量折减为$ 0.857 {\text{×}} 10^{4}\; \mathrm{MPa} $ (陈国兴,2007)。土-地下结构静动力耦合有限元模型宽为120 m,高为42 m,在基岩下20 m处截断,如图5所示。不考虑土与结构之间地震荷载所导致的错动、滑移等现象,仅考虑它们之间完全黏结的情况,即:使用绑定约束(Tie约束)。使用自编的程序VBEA2.0实现地震动从模型底部垂直输入,输入的地震动加速度时程与4.1节一致,如图3所示。

    表  1  大开地铁车站数值模型场地条件参数
    Table  1.  Mechanical parameters of the soil layers for the Dakai metro station
    土质模型各层厚度/m密度/(kg·m−3较软地层模型vS
    /(m·s−1
    真实地层模型vS
    /(m·s−1
    较硬地层模型vS
    /(m·s−1
    泊松比
    人工填土11 900801402900.333
    全新世砂土41 900801402900.320
    全新世砂土31 9001101703200.320
    更新世黏土31 9001301903400.400
    更新世黏土61 9001802403900.300
    更新世砂土52 0002703304800.260
    基岩202 0001 0001 0001 0000.260
    下载: 导出CSV 
    | 显示表格
    图  5  大开地铁车站有限元模型
    Figure  5.  Finite element model of Daikai subway station

    选取车站结构中柱和距车站右侧墙3.5 m处土层的水平峰值加速度和水平最大相对位移,以及车站结构重要位置处的最大主应力来反映结构和地基土的地震响应。上述数据均取地震动力分析时程中绝对值的最大值。结构重点位置监测点设置如图6所示。

    图  6  车站重点位置监测点
    Figure  6.  Monitoring points in important position of subway station

    根据大开车站真实地层模型计算获得的考虑和不考虑静动力耦合作用下的地震动力反应结果对比如图7所示。分析可知,大开车站场地模型在考虑和不考虑静动力耦合作用的情况下,地基土随深度变化的峰值加速度和最大相对位移几乎完全一致,车站中柱的最大相对位移在两种情况下也吻合得很好,而车站结构关键点的最大主应力也有类似的结果,在数值上相差不大。图8给出了真实地层工况下合理进行静动力边界转换后,模型左边界顶、底部节点及模型中部顶点及底边节点的水平方向加速度时程和位移时程曲线,可以看出各时程前段均未出现振荡。软弱地层模型和较硬地层模型计算结果与上述真实地层模型相似,其场地土和中柱的地震动力反应在考虑和不考虑静动力耦合的作用下差别不大,由于篇幅原因不再给出其具体结果。

    图  7  带地下结构的成层场地地震动力反应计算结果
    (a) 场地土PGA;(b) 场地土最大相对位移;(c) 车站中柱最大相对位移;(d) 车站结构重要部位最大主应力
    Figure  7.  Earthquake dynamic response of underground structure in the layered site
    (a) Maximum acceleration of the site;(b) Maximum relative displacement of the site;(c) Maximum relative displacement of the central column;(d) Maximum principal stress at points in the station
    图  8  带地下结构的成层场地左边界及中部节点水平向的加速度时程(左)和位移时程(右)
    (a) 左边界顶部节点;(b) 地表中部节点;(c) 左边界底部节点;(d) 底边界中部节点
    Figure  8.  Horizontal acceleration (left) and displacement (right) time history of nodes at the left boundary and the vertical middle of underground structure in the layered site
    (a) Node at the top of the left boundary;(b) Node at the centre of the surface;(c) Node at the bottom of the left boundary;(d) Node at the centre of the base

    为进一步研究静动力耦合作用的影响对土-地下结构计算模型的影响,利用二阶欧几里得范数(式(6))计算以比较这些模型的边界节点、地基土和结构中柱在考虑静动力耦合作用时得到的加速度时程和位移时程与不考虑静动力耦合作用工况下的误差。三种不同软硬地层情况下的地震动力反应时程误差结果如图9所示。对于三种不同波速的土层模型,静动力耦合作用对土-地下结构计算模型影响不大,对地基土位移、地下结构和边界节点的位移及加速度影响很小,误差几乎均为零;对于软土层而言,不考虑静动力耦合,地基土会产生最大值为1.2%左右的误差,并且22 m以下的基岩层节点的误差仍接近于零。

    图  9  带地下结构的成层场地地震动力反应时程结果误差
    Figure  9.  The error of the earthquake dynamic response of sub-structure in the layered site

    上述结果表明,弹性情况下,静动力耦合作用对土-地下结构计算模型影响很小,仅在结构关键点处有较小的差异,对地基土位移时程以及中柱加速度和位移时程等地震动力反应几乎没有影响。

    综合3.1节和3.2节算例结果,在弹性情况下,如果静动力边界转换正确,那么考虑静动力耦合作用的自由场或土-地下结构计算模型(施加初始应力)与对应的不考虑静动力耦合计算模型动力计算结果基本相同。因此,利用叠加原理给出的以不考虑静动力耦合时的模型的计算结果作为动力反应分析的参考解来验证静动力边界转换正确性,既适用于多土层自由场静动力耦合计算模型又适用于土-地下结构模型的静动力耦合模型。需要注意的是,土-地下结构计算模型的土体和地下结构都需要使用弹性本构关系,且土体与结构之间只能使用绑定约束,这样做的目的是避免模型中产生材料非线性和接触非线性计算,以致影响静动力边界转换合理性的判断。

    在静动力转换时,仅用边界上的反力来替换静力边界条件,不施加单元应力,再进行一次地应力平衡,随后进行动力计算是不合理静动力边界转换的一种。本节用参考解方法来验证这种边界转换方法的正确性和合理性以及分析不合理静动力边界转换前期产生振荡对后续动力分析的影响。首先建立两个100 m×50 m的自由场模型,匀分为5层,材料参数仍使用3.1节中的参数,如图2所示,采用弹性本构模型。其中一个模型不考虑静动力耦合作用(不施加重力引起的初始应力),作为参考解模型,另一个模型施加重力引起的初始应力并进行静动力边界转换来考虑静动力耦合作用(仅施加边界节点反力,不施加单元应力)。

    为了分析模型在不合理静动力边界转换方法下振荡程度以及更好地给出参考解法的使用方法,在本节的自由场模型中设置12个水平加速度监测点,如图2中红色点所示,监测点从模型左边界开始,相邻监测点之间间隔15 m。图10比较了模型左边界及中部的顶部节点的水平加速度时程与对应参考模型的参考值。分析可知,对于模型边界顶部节点和模型顶部中点,不合理的静动力边界转换对初始阶段的水平加速度有一定影响,而后续加速度时程曲线与参考解基本一致。图11给出了12个不同监测点在第一秒内地震荷载下的水平加速度动力响应。分析可知,对于处在同一深度的节点,越靠近左侧边界,不合理静动力转换引起的振荡强度越大,模型中间位置节点的振荡最小;此外,所有边界上的监测点(N1N6)振荡均比其它位置监测点强烈,侧边界顶部节点N1的振荡强度最大。这可能是自由场放大效应和不平衡力共同作用的结果。因此,在使用参考解方法时,建议选取侧边界底部、顶部节点和模型中部顶部和底部的节点的加速度动力响应来比较验证,确认静动力边界转换合理性及黏弹性边界的正确施加和地震动的正确输入。以上结果表明,不施加单元应力的静动力边界转换方法,在弹性条件下,会对动力计算结果加速度的初始阶段产生一定的影响(振荡),而后阶段与参考解基本一致,边界顶部节点的动力响应影响较中间节点更大。

    图  10  模型左边界(a)及中部(b)的顶部节点水平加速度时程
    Figure  10.  Acceleration time history of nodes at the top of the left boundary (a) and the vertical middle model (b)
    图  11  不合理静动力边界转换下模型不同位置节点水平加速度时程曲线
    Figure  11.  Horizontal acceleration time history of different nodes with unreasonable static-dynamic boundary switch method

    黏弹性边界是由模型边界上的并联的弹簧阻尼系统构成,是一种应力型边界。由式(2)可知,在地震输入时,需要输入的等效节点力必须要能够克服弹簧和阻尼器产生的附加应力以及自由场地震动在边界上产生的应力,使得边界与自由场运动一致,才能保证边界上的位移和应力边界条件的连续。由4.1节可知,不合理的静动力边界转换后加速度时程振荡在模型侧面和底面最显著,因此为探究其振荡原因,输出进行合理静动力边界转换和不合理转换两种情况下,多土层自由场模型黏弹性边界节点(N1N6)上弹簧单元在进行静动力边界转换后动力分析步前的弹力绝对值如图12所示。分析可知,不合理静动力边界转换模型边界节点弹簧在动力分析步前已经有很大的弹力,表明静动力边界转换后有力作用在黏弹性边界弹簧单元上,此力称为不平衡力。模型整个体系的地应力平衡是重力及其引起的单元内力和边界反力平衡的结果,因此去掉边界后除了施加边界节点反力外,不施加重力或模型单元内力,必然会使弹簧变形并产生附加应力,意味着动力分析步中地震动(等效节点力)的施加将出现误差。所以在进行土-地下结构静动力耦合动力计算分析时,必须保证静力边界转换成动力边界过程中,整个计算模型处于静力平衡状态,黏弹性边界上不能受到不平衡力的作用,否则会在动力分析步初始时刻引起动力效应,具体表现为模型在未受到地震动力荷载作用前,模型节点动力响应已经出现了较大振荡,最终导致地震动力计算结果不合理。

    图  12  多土层自由场模型黏弹性边界节点弹簧弹力受力情况
    Figure  12.  The magnitude of elastic fore for spring elements located at the viscous-spring boundary

    本文根据弹性叠加原理,提出了一种验证静动力边界转换合理性的方法即参考解法,利用有限元程序ABAQUS,并结合自编等效节点力的计算程序以及黏弹性边界的自动施加程序VBEA2.0,建立了二维成层自由场静动力耦合模型和土-地下结构静动力耦合模型,对比分析了参考解法的适用性;在此基础上,针对黏弹性边界条件下传统静动力边界转换过程出现的振荡反应现象进行了相关研究,探讨了造成振荡的原因和机理。主要结论如下:

    1) 合理考虑其静动力耦合效应需要遵循以下程序:首先进行地应力平衡静力计算;然后进行静动力边界转换;最后进行动力计算。在静动力边界转换过程中,除施加边界节点反力来替代静力边界条件外,还需施加所有单元的应力场以及重力场。

    2) 根据二维成层自由场静动力耦合模型和土-地下结构静动力耦合模型的计算分析结果,以不考虑静动力耦合作用的模型计算结果作为动力反应分析的参考解来判断静动力边界转换的合理性方法(参考解法),能够较好地验证静动力边界转换的合理性,并且能直观地反映出动力计算初期是否有波动以及波动对后阶段动力计算的影响程度。该方法可在土-地下结构的静动力耦合动力分析中为静动力边界转换的合理性检验提供一种有效工具。

    3) 静动力耦合计算模型静力不平衡是导致静动力边界转换出现振荡的物理机制;边界反力、重力、单元应力的施加合理与否是影响静动力边界的不合理转换的主要因素。

  • [1] 何茂富,四川地区区域重力场的基本特征,四川地震,1984, 2: 9——11,

    [2] 赵友年、李春生、赖绍民、咸向明,四川省大地构造若干基本问题(四川省大地构造图简要说明书),31,四川省地质局地质矿产科学研究所出版,1983,

    [3] Mitchell, B. J, and B. M. Hashim, Seismic velocity determinations in the new Madrid seismic zone: a new method using local earthquakes, Bull. Seism. Soc. Am., 67, 413——424, 1977

    [4] Dobrin, M. B,Introduction to Geophysical Prospecting, and ed., 294——303, McGraw——Hill, New York,1976.

    [5] 赵珠,华北地区上地慢速度结构,地球物理学报,超,341——354, 1983.

    [6] 宋仲和,朱介寿,安昌强,张丽娟,北京一萨哈林剖面地馒纵向速度结构,地球物理学报,24,310——318,1981,

    [7] Cart, D. A. and G. A. Bollinger, A regional crustal velocity model for the Southeastern United Staties,Bull. Seism. Soc. Am., 71, 1829——1847, 1981.

    [8] Gupta, I. N., B. W. Barker, J. A. Buretti and Z. A. Der, A study of regional phases from earthquakes and explosions in Western Russia, Bulf. Seism. Soc. Am., 70, 855——872, 1980.

    [9] Bruce, R. J. and D. L. Anderson, Travel——times, apparent velocities and amplitudes of body waves, Bull.Seism. Soc. Am., S8, 339——366, 1968.

    [1] 何茂富,四川地区区域重力场的基本特征,四川地震,1984, 2: 9——11,

    [2] 赵友年、李春生、赖绍民、咸向明,四川省大地构造若干基本问题(四川省大地构造图简要说明书),31,四川省地质局地质矿产科学研究所出版,1983,

    [3] Mitchell, B. J, and B. M. Hashim, Seismic velocity determinations in the new Madrid seismic zone: a new method using local earthquakes, Bull. Seism. Soc. Am., 67, 413——424, 1977

    [4] Dobrin, M. B,Introduction to Geophysical Prospecting, and ed., 294——303, McGraw——Hill, New York,1976.

    [5] 赵珠,华北地区上地慢速度结构,地球物理学报,超,341——354, 1983.

    [6] 宋仲和,朱介寿,安昌强,张丽娟,北京一萨哈林剖面地馒纵向速度结构,地球物理学报,24,310——318,1981,

    [7] Cart, D. A. and G. A. Bollinger, A regional crustal velocity model for the Southeastern United Staties,Bull. Seism. Soc. Am., 71, 1829——1847, 1981.

    [8] Gupta, I. N., B. W. Barker, J. A. Buretti and Z. A. Der, A study of regional phases from earthquakes and explosions in Western Russia, Bulf. Seism. Soc. Am., 70, 855——872, 1980.

    [9] Bruce, R. J. and D. L. Anderson, Travel——times, apparent velocities and amplitudes of body waves, Bull.Seism. Soc. Am., S8, 339——366, 1968.

  • 期刊类型引用(1)

    1. 薛景宏,陆大伟,李东瑞,程安顿,韩通. 地震波斜入射下埋地腐蚀管道地震响应分析. 河北工程大学学报(自然科学版). 2024(04): 18-27+35 . 百度学术

    其他类型引用(1)

计量
  • 文章访问数:  1397
  • HTML全文浏览量:  43
  • PDF下载量:  135
  • 被引次数: 2
出版历程
  • 发布日期:  2011-09-01

目录

/

返回文章
返回