0.   引言

在对含饱和砂土层场地地震反应研究中,孔隙水压力上升引起的砂土液化是一种典型的地震灾害现象。饱和砂土动力有效应力本构关系是否能真实的表征饱和砂土的物理力学特性,对场地地震液化反应模拟的准确性至关重要,因此可液化土层的动力响应模拟成为国内外学者的研究热点。合适的饱和砂土动力有效应力本构不仅仅需要模拟液化过程中土体的软化现象,还需要对液化过程中的土体应力应变关系的变化细节进行描述。

选择一个合理的土体动力非线性本构模型是液化土层地震动模拟计算分析的前提,以双曲模型为主的土体弹塑性动力本构模型广泛应用于土层的地震反应分析中。基于Masing法则^[1]的双曲模型是国内外常用的基于骨架曲线的土体动力非线性本构,Rasenblueth^[2]和Newmark^[3]考虑到可能出现应力值超过土体极限应力值的不合理现象,扩展了Masing法则。赵丁凤等^[4]在Davidenkov本构模型基础上,提出代替“扩展Masing法则”的修正加载准则,解决了该准则在程序实现中需要记录大量数据的问题。为使应力应变关系符合阻尼比试验结果,王志良等^[5]提出了“阻尼比退化系数”对Masing关系进行了修正。李小军在王志良模型的基础上提出“基于骨架曲线修正”^[6]处理的“动态骨架曲线”^[7]的概念,给出了能考虑土体试验特性且形式简单的一维土体动力本构关系。等效线性化是土层一维地震反应分析中常采用的方法,但该方法是平均滞回曲线意义上的线性方法,无法体现土体的塑性变形,且难以考虑孔压上升对饱和砂土层剪切强度的影响。笔者在李小军模型的基础上提出了以对数函数为骨架曲线^[8],可拟合阻尼比、刚度比试验曲线的土层时域非线性地震反应分析方法,并开发了相应的一维地震反应计算程序。该方法得到的本构模型滞回圈较宽,能更好的体现土体的非线性和塑性特性。国内外专家学者根据循环三轴或剪切试验提出了不同的孔压增长模型,如以seed^[9]为代表的根据孔压比和振次比关系提出的孔压增长公式和Ishibash等^[10]根据孔压比与动应力及振次的关系确定的孔压增长模型,但是这些模型不能用于不均等固结土体的孔压增长计算,且不适用于不等幅荷载情况下的孔压计算。孙锐等^[11]、陈龙伟等^[12]在seed模型基础上提出了适用于水平场地的不均等固结饱和砂土孔压发展增量计算模型。付海清等^[13]以加速度、埋深、砂土密实度等现场参数为基本指标,构建直接适用于水平场地饱和砂土孔压增量模型。上述孔压模型只能模拟在应力循环作用下的孔压增长,无法反映应力卸载条件下的实际孔压值的减小。实际孔压模型由不可逆孔压项和可逆孔压项组成,文中所用的孔压计算模型是在陈龙伟的不可逆孔压增量的基础上,加入笔者提出的可逆孔压模型,两者结合得到的可体现孔压波动的新型孔压模型。该模型能较好的模拟应力循环作用下饱和砂土的实际孔压变化情况。

此外,孙锐^[14]通过高精度动三轴的试验结果发现砂土循环最大剪切模量$G_{\max ,N}$、极限剪应力${\tau }_{\text{ult},N}$和孔压比的关系与砂土类型和相对密度无关,且可以描述成统一的线性关系。基于此,本文提出了以对数动骨架曲线本构土层时域非线性地震反应分析方法为基础,结合能计入可逆超静孔隙水压的孔压计算模型,通过饱和砂土的有效应力对其剪切模量和抗剪强度进行修正,得到可模拟含饱和砂土场地液化地震反应的土层时域非线性有效应力本构。

在自行编制的土层地震反应分析程序Soilresp1D的基础上,加入基于对数动骨架曲线并考虑可逆超静孔隙水压力的有效应力本构,实现了可液化场地动力反应分析。通过对含饱和松砂层场地数值分析,并与基于不同孔压增长模型的有效应力本构的时域分析结果和实际地震动观测记录对比表明,本文提出的基于对数动骨架曲线并考虑可逆超静孔隙水压的孔压计算模型的有效应力本构的模拟结果较为合理,与实际地震动记录相符,且反映了液化对场地地表加速度峰值、反应谱及饱和砂土层抗剪强度的影响规律。

1.   有效应力分析方法

1.1   考虑可逆超静孔隙水压的孔压计算模型

含饱和砂土场地地震反应模拟结果是否能反映饱和砂土液化对场地地震反应的影响取决于合理的孔压计算模型,这一模型应能真实模拟饱和砂土层在每一应力循环后的孔压值。

常用的孔压模型只能反映孔压在应力循环加载作用下的增长,无法体现应力卸载时孔压的减小,如图1（a）中的红线所示,而饱和砂土不排水循环三轴试验实测孔压值如图1（a）中的蓝线所示。据此可将饱和砂土不排水循环三轴试验实测孔压值分为可逆孔压$u_{\text{re}}$（图1（a）黑线）和不可逆孔压$u_{\text{ir}}$（图1（a）红线）,如此可知,不可逆孔压$u_{\text{ir}}$随着循环次数的增加单调增长,本文将依据陈龙伟^[12]提出的适用于水平成层场地的不均等固结饱和砂土的孔压发展增量计算模型表示,这一模型是通过振动台和振动离心机土层液化试验得到的,能真实反映实际水平成层场地在水平剪切地震动作用下的土体单元水平往返剪切力,孔压发展增量计算模型：

图  1  饱和砂土不排水循环三轴试验过程中孔压变化

Figure  1.  Variation of pore pressure of saturated sand by undrained cyclic triaxial tests

下载: 全尺寸图片 幻灯片

式中$u_{\text{ir,0}}^{*}$为初始不可逆孔压比；$\Delta u_{\text{ir},N}^{*}=\Delta u_{\text{ir},N}/$ ${\overline{\sigma }}_0$,$\Delta u_{\text{ir},N}$为第N次应力循环中的不可逆孔压的孔压增量；${\overline{\sigma }}_0$为初始有效围应力；$\Delta u_{\text{ir},N-1}^{*}=\Delta u_{\text{ir},N-1}/{\overline{\sigma }}_0$,$u_{\text{ir},N-1}$为第N-1个应力循环中的不可逆孔压；$c_{1,0},$ $c_{1,a},c_{1,b},A_{4,0}$为试验参数,与土的密实程度有关；$K_{\text{c}}$为固结比；$N_{\text{ep}}$为等效应力作用次数,${\tau }_i$为第i次循环剪应力幅值（$1\le i\le N$）；${\tau }_N$为第N次循环剪应力幅值,$\alpha$为材料参数。式中相关参数取值详见文献[12],在此不再赘述。

由图1（a）可知,可逆孔压$u_{\text{re}}$是非正的,对于每一循环加载过程,逐渐增大至零,随卸载逐渐减小,反向加载时,可逆孔压亦逐渐增大。根据如图1（b）所示的饱和砂土不排水循环三轴试验数据分析可得,每一应力循环的可逆孔压与该循环的应力呈线性关系,且每一循环的可逆孔压初始值为$-\Delta u_{\text{ir}}$。可逆孔压$u_{\text{re},N}$定义为

在对饱和砂土层进行时域非线性地震反应有效应力分析时,根据上述两式即可得到循环应力作用下饱和砂土的实际孔压值（图1（a）中的u线）。

1.2   基于孔压的剪切模量和抗剪强度修正方法

饱和砂土层每一循环的最大剪切模量和抗剪强度受这一循环孔压变化的影响显著,最大剪切模量和抗剪强度随着孔压比的上升不断减小。为考虑在循环荷载作用下孔压的上升对土体软化过程的影响,则每经过一次应力循环需对砂土的最大剪切模量和抗剪强度进行修正。孙锐等^[14]通过拟合高精度动三轴的试验数据发现饱和砂土最大剪切模量$G_{\max ,N}$、极限剪应力${\tau }_{\text{ult},N}$和孔压比的关系几乎不受砂土类型和相对密实度的影响,且能简单描述成统一的线性关系,其具体函数表达式为

式中,$G_{\max ,N}$,${\tau }_{\text{ult},N}$分别为N-1次应力循环后砂土的最大剪切模量和极限剪应力,$U_N=u_N/{\overline{\sigma }}_0$。

1.3   基于对数动骨架曲线含饱和砂土层场地时域非线性本构有效应力分析方法

在文献[8]提出的基于对数动骨架曲线时域非线性土体动本构基础上,考虑可逆超静孔隙水压的孔压计算模型及剪切模量和抗剪强度随孔压变化的修正关系,提出了一种新的有效应力本构：

式中,相关参数的物理含义及其取值详见文献[8]及前文相关定义,在此不再赘述。

基于文献[8]中的式（18）,（19）,（21）结合式（4）所表述的有效应力本构,可形成一种基于对数动骨架曲线时域非线性土体动本构并考虑可逆超静孔隙水压力的含饱和砂土场地时域非线性地震反应分析方法,这一方法可以计算不均等固结条件下的水平成层场地在不等幅的地震荷载作用下的孔压变化情况,且可模拟饱和砂土液化过程中的砂土剪切模量和抗剪强度软化特征。将这一方法引入笔者研发的土层场地时域非线性动力反应分析程序Soilresp1D,实现了含饱和砂土场地时域非线性地震反应分析。

图2为依据本文有效应力分析方法得到的饱和松砂土体试验模型等幅循环荷载作用下的应力应变曲线及孔压变化,数值模拟中当孔压比达到1时,砂土层发生液化,且液化后土层的剪切模量取0.0125 $G_{\max }$^[15]。由图可知本文有效应力本构所反映的应力应变关系合理,能实现在动荷载作用下孔压上升引起的饱和砂土剪切模量和抗剪强度软化特征。

图  2  基于本文本构模拟试验模型中饱和砂土的动力反应

Figure  2.  Stress-strain curves and simulated results of variation of pore pressure based on constitutive model for effective stress

下载: 全尺寸图片 幻灯片

2.   可液化土层场地非线性地震反应对比分析

利用本文方法,对含饱和松砂层场地和实际观测台阵场地进行地震液化反应分析,并通过基于丰万玲孔压增长模型的有效应力本构的计算结果和实际地震动观测记录对比分析,验证了基于对数动骨架曲线并考虑可逆超静孔隙水压力的孔压模型的有效应力本构的合理性。同时,讨论了饱和砂土层液化对场地地表加速度峰值、反应谱及饱和砂土体剪切模量的影响。

2.1   基于不同孔压模型的有效应力本构计算结果对比

数值分析中输入图3（a）所示时长82 s、峰值加速度2 m/s²、离散时间步距为0.02 s的地震动时程。按照时域非线性分析稳定性的要求并结合地震动截止频率（15 Hz）确定数值模拟时间步距为0.0025 s,为实现数值分析将输入地震动时程通过插值成离散散步距为0.0025 s。计算场地模型包含黏土层、饱和砂土层和下伏基岩,总厚32 m,各岩土层模型计算参数见表1。根据石兆吉给出的液化判别式得到的该计算模型中砂土的液化临界剪切波速为267.1 m/s^[16],各砂土层的实测剪切波速均小于该值,所以该模拟场地为液化场地。按照土层时域动力反应分析方法精度要求确定了沿场地深度的离散步距$\Delta z$为1 m,如此将计算模型划分成32层,其中层7～层9为饱和砂土层。在对模型进行动力反应计算时,每一次应力循环后（应力两次过零线）需对$G_{\max ,N}$和${\tau }_{\text{ult},N}$修正,当饱和砂土层发生液化后土层参数保持不变。数值模拟结果如图3（b）～（d）、图4所示,图3（b）～（e）为本文方法分析结果和基于丰万玲孔压增长模型的有效应力本构分析结果,图4为本文方法分析得到的含饱和砂土场地模型在地震动作用下液化和未液化砂土层的应力应变曲线。

图  3  输入地震动时程和两种有效应力本构的场地地震反应

Figure  3.  Time histories of input ground motion and seismic responses of two constitutive models for effective stress

下载: 全尺寸图片 幻灯片

表  1  场地计算模型参数

Table  1.  Model parameters for site salulation

土类	深度/m	剪切波速/(m·s-1)	密度/(t·m-3)	$b/a$	$a_1$ /(10-3)	$b_1$
黏土	6	120.0~142.2	1.95	1160	0.82	5.9
细砂	9	142.2~153.3	1.49	1934	2.24	5.2
黏土	30	153.3~231.0	1.95	1160	0.82	5.9
基岩	32	511.0	2.65

下载: 导出CSV 

| 显示表格

图  4  本文方法模拟液化与未液化饱和砂土层应力-应变关系曲线

Figure  4.  Simulated results of stress-strain relationship curves between liquefied and unliquefied saturated sand

下载: 全尺寸图片 幻灯片

数值模拟结果显示,基于本文分析方法得到的地面加速度峰值为1.62 m/s²,基于丰万玲孔压增长模型有效应力本构计算得到的地面加速度峰值为1.5 m/s²。根据孙锐对比液化与非液化场地加速度反应谱,可知液化场地减少地震动高频分量但同时对低频分量显著放大^[17]。由图3（d）两种本构地面加速度反应谱,本文本构计算场地在周期为1 s内的高频地震反应明显减小,在周期大于1 s的地震反应有明显放大,丰万玲本构计算结果则无该现象。根据反应谱对比图判别本文本构计算该模型场地发生液化,丰万玲本构计算该场地未发生液化,该结果与本文孔压计算结果一致。由图3（b）,（c）可以看出,丰万玲有效应力本构分析所得地面加速度时程与入射加速度时程无明显差异,本文分析结果在饱和砂土液化后高频地震动得到了明显的衰减。图3（e）所示结果显示,丰万玲孔压增长模型分析孔压增长缓慢,最大孔压比仅为0.1,未发生液化；而本文分析结果表明,层7于16.55 s孔压比达1.0,即开始液化,总体表现为孔压比随作用时间逐渐增加,当加速度接近峰值时,上层饱和砂层孔压迅速上升,并发生液化；对于下层饱和砂土,孔压比随作用时间逐渐增加,达0.54后,趋于平稳,不难发现其平稳段对应于上层饱和砂土液化阶段；同时可以看出本文方法模拟结果能较真实的反映饱和砂土地震作用下的实际孔压变化特征,即孔压随地震荷载作用有增有减；相比于本文模拟结果,丰万玲有效应力本构虽能描述孔压增长过程,但不能合理的模拟饱和砂土孔压在地震荷载作用下孔压变化特征,未能反映含饱和砂土液化对场地地震反应的影响。

此外,由图4可以看出,层7饱和砂土层在液化前应变较小,土体剪切模量较大,应力应变曲线较陡,液化后应变增大,剪切模量减小,应力应变曲线迅速变缓。层9饱和砂土层孔压比未达1.0,未发生液化,其应变较小,剪切模量的折减较为缓慢。

综上所述,本文方法能再现饱和砂土层在强地震作用下孔压变化特征,能反映液化饱和砂土层的减震效应和高频地震动的衰减特征。

2.2   Fire Station场地模型计算结果与实测结果对比

本文以加州Corralitos地区的Fire Station场地影响观测台阵的场地资料为基础,建立了含饱和砂土层场地地震液化反应分析模型,以地表下102 m基岩处实测加速度时程为作用地震动,加速度时程如图5（a）所示,其峰值加速度为0.75g。按照一维土层场地时域非线性计算精度和稳定性要求确定了各土层离散步距和计算时间步距,采用本文方法模拟得到了场地地表加速度反应,并与实际地表观测记录对比,验证了本文方法数值模拟的可靠性。分析模型厚度为102 m,由松细砂、中密砂及黏土组成,图5（b）给出了模型剪切波速随深度变化趋势,模型相关参数如表2所示。同样根据石兆吉给出的液化判别式得到该实际场地模型上层厚度为26 m细砂的液化临界剪切波速为371.7 m/s,上层细砂的实测剪切波速均小于该值,该场地为液化场地。

图  5  实测基岩地震动加速度时程和场地剪切波速随深度的变化

Figure  5.  Time histories of measured bedrock acceleration and variation of shear wave velocity with depth

下载: 全尺寸图片 幻灯片

表  2  场地计算模型参数

Table  2.  Model parameters for site calculation

土类	层厚/m	深度/m	密度/(t·m-3)	$b/a$	$a_1$ /(10-3)	$b_1$
松细砂	26	26	1.92	1711	7.18	3.717
黏土	12	38	1.68	1527	3.73	4.331
松细砂	22	60	1.92	1051	8.36	4.245
中密砂	40	100	1.68	709	10.53	4.081
基岩	2		1.92

下载: 导出CSV 

| 显示表格

由表2可以看出,分析模型覆盖土层主要为砂和黏土,浅层砂土为饱和松细砂。图6（a）,（b）所示实测和计算所得的地表加速度时程均显示10s左右地震动达到峰值,土层发生液化,该结果与图6（d）所示计算孔压比一致,验证了本文计算方法的可靠性与有效性。图6（d）所示不同深度饱和砂土层孔压比表明,浅层饱和松细砂层孔压比随时间快速增加,且埋深越浅发生液化越快,而下层饱和松细砂层孔压增长较为缓慢,未发生液化。由此可知本文有效应力本构模型能真实反映饱和砂土层液化特征,亦表明地表下20～30 m深度内的饱和砂土能产生液化,这一结果的可信性得到了2011年新西兰Mw6.3地震中发现20～30 m深度内饱和砂土产生液化现象的证实^[18]。图6（c）所示为地表与基岩加速度反应谱,根据孙锐提出的液化场地地表反应谱特征,可知实测地表与计算地表反应谱均显示该场地发生液化。在约0.3 s周期内,计算地表加速度反应谱相对基岩加速度反应谱明显减小,与实测地表加速度反应谱相近,变化趋势相似。

图  6  数值模拟结果与观测结果

Figure  6.  Simulated and observed results

下载: 全尺寸图片 幻灯片

由此不难看出,饱和砂土层液化对短周期地震动具有明显的减震效应,本文方法模拟结果可信,能模拟饱和砂土层孔压变化规律和液化对地表加速度峰值及反应谱的影响特征。

3.   结论

从对数骨架曲线出发并考虑可逆超静孔隙水压的孔压模型,结合每一加载循环饱和砂土体的剪切模量、抗剪强度修正方法,提出了一种基于对数动骨架曲线时域非线性土体动本构并考虑可逆超静孔隙水压力的含饱和砂土场地时域非线性地震反应分析方法。通过典型算例验证了本文方法的可行性和模拟结果的合理性,得到以下4点结论。

（1）基于饱和砂土不排水循环三轴试验结果,提出了考虑可逆超静孔隙水压的孔压计算模型,该模型不仅能模拟应力循环加载作用下孔压的增长,还能体现应力卸载时孔压的减小。

（2）数值模拟结果对比分析表明,本文方法能再现饱和砂土层在强地震作用下孔压变化特征,能反映液化饱和砂土层的减震效应和高频地震动的衰减特征。

（3）场地影响观测台阵Fire Station分析模型数值结果显示,本文有效应力本构模型能真实反映饱和砂土层液化特征,分析方法是可行的,且结果可信。

（4）饱和砂土液化对地震动峰值加速度与加速度反应谱等地表地震动参数,以及饱和砂土体剪切模量和抗剪强度具有显著影响。