李纯, 修占国, 王斐笠, 王秋鸿
东北大学 资源与土木工程学院, 辽宁 沈阳 110819
收稿日期：2016-12-12
基金项目：辽宁省自然科学基金资助项目(2013010311)。
作者简介：李纯(1971-)，男，安徽潜山人，东北大学副教授。

摘要：基于广义塑性理论, 针对Pastor-Zienkiewicz(P-Z)砂土模型静力参数, 划分为3个方面8个参数, 分别分析了各参数变化时对竖向位移的敏感性.研究表明:模型各参数在-40%~40%变化时, 弹性模量E的敏感性最大, 竖向位移变化率可达66.05%, 其次为内摩擦角φ和泊松比μ, 竖向位移变化率在±10%以内.而与塑性模量相关的参数, 竖向位移变化率均低于±1%, 当缺少实测数据时, 可根据土工实验资料取值.研究结论可为P-Z砂土模型在地基变形计算时参数的选取提供依据.
关键词：Pastor-Zienkiewicz本构模型模型参数竖向位移变化率砂土
Sensitivity Analysis of Static Parameters About Deformation in Pastor-Zienkiewicz Sandy Model
LI Chun, XIU Zhan-guo, WANG Fei-li, WANG Qiu-hong
School of Resources & Civil Engineering, Northeastern University, Shenyang 110819, China
Corresponding author: XIU Zhan-guo, E-mail: xiuzhanguo109@126.com
Abstract: Based on the theory of generalized plasticity, for static Pastor-Zienkiewicz(P-Z) model parameters, it can be divided into eight parameters from three aspects. The sensitivity analysis of vertical displacement varying with each parameter was performed. The results showed that elastic modulus E exhibits the largest sensitivity, and the corresponding maximum rate of vertical displacement is 66.05%, when the every parameter is changed between -40%~40%. Followed by the angle of internal friction φ and Poisson's ratio μ, the rate of vertical displacement was within ±10%. But for the parameters about plastic modulus, the rate of vertical displacement was less than ±1%. When lack of measured data, the values of plastic parameters can come from a great amount of geological data. The conclusions can provide important basis for the selection of P-Z sandy model parameters in the deformation calculation of sub grade.
Key Words: Pastor-Zienkiewicz modelmodel parametersvertical displacementrate of changesand
1984年, Zienkiewicz和Mroz提出了广义塑性理论的基本框架^[1].Pastor-Zienkiewicz模型^[2](P-Z模型)是基于广义塑性理论建立的土体本构模型.P-Z模型可以模拟砂土在单调和循环加载下的变形行为^[3].随后, 在不改变模型参数数量的前提下, Pastor等^[4]提出了一个修正的非线性Pastor-Zienkiewicz-Chen模型, 该模型对不排水条件下饱和砂土的动力液化方面有较好的预测能力^[5].李宏恩等^[2]对P-Z模型的静力部分进行修正, 引入模型参数, 并对Toyoura砂进行试验分析.
P-Z模型大多用于模拟地震以及海浪作用下的土体变形^[6], 而对于地基土变形涉及较少.实际工程中, 地基土受竖向荷载作用产生变形, 一般不存在加载后再卸载的情况, 因此, 可将建筑物对地基土的作用视为单调连续加载.本文针对P-Z模型中有关参数对地基变形的贡献进行敏感性分析, 表现为:与弹性状态相关的2个参数, 弹性模量E和泊松比μ; 定义临界状态的2个参数, 相对密实度D_r和内摩擦角φ; 定义塑性模量的4个参数, ξ, H₀, β₀, β₁.通过讨论各参数在变形计算中的敏感性, 为地基计算参数的选取提供依据.
1 静力作用下的P-Z模型P-Z模型定义于p-q-θ空间中, 为了便于数值计算, 将其推广到三维笛卡尔坐标系^[7].该模型在模拟砂土变形时, 采用非相关流动法则(f≠g), 其屈服面方程(f)和塑性势面方程(g)为

(1)

式中:p=I₁; ; I₁, J₂分别为第一应力不变量和第二偏应力不变量; M_g为临界状态线(CSL)的斜率; α为常数; p_e, p_g为初始标准应力.

(2)

(3)

式中:D_r为土的相对密实度; φ为内摩擦角; θ为Lode角; J₃为第三偏应力不变量.
屈服面(f)和塑性势面(g)的加载方向矢量为

(4)

材料的应力增量dσ与应变增量dε的关系为

(5)

式中弹塑性矩阵C^ep在塑性力学中表示为

(6)

对于初次加载, 塑性模量H_L可表示为

(7)

(8)

式中:β₀, β₁均为材料参数; ξ为累计偏应变; η_f为应力比参数; H_v, H_f为应力比对塑性模量的函数; H_s为偏应变硬化函数.
为便于数值计算, 将式(4)从p-q-θ空间转化为笛卡尔坐标系.n与n_gL的计算方法相同, 只需将f替换为g即可.

(9)

(10)

将式(5)展开可得

(11)

假设土体变形通过z方向体现, 土体压缩模量用增量表示为竖向应力增量与应变增量之比:

(12)

式(12)虽是弹性力学公式, 但通过弹塑性矩阵C^ep结合基本假定(dε_x=dε_y=0)进行修正, 使其能够考虑土体的塑性变形.只要求得弹塑性矩阵C^ep, 便得到修正的压缩模量值.
为了说明静力作用下P-Z模型在计算地基变形方面的合理性, 本文针对数值算例进行简化验证, 计算模型如图 1所示.
图 1(Fig. 1)

图 1 数值算例示意图Fig.1 Sketch of numerical example

基底附加应力p₀=189.4 kPa, 运用平面网格划分技术^[8]计算有效应力的传递和地基土的变形时, 算例的参数基准值为:E=28.8 MPa, μ=0.3, H₀=800, φ=36.5°, D_r=42.9%, ξ=0.7, α=0.45, β₀=2, β₁=0.2, Ψ_s=0.21.土层总厚度为4.6 m, 计算出土层变形值, 经修正系数Ψ_s修正, 最终结果为5.730 3 mm.
规范法计算最终变形值为6.277 mm, 与本文计算结果相差0.546 7 mm, 误差为8.71%, 相差不大.误差原因为：对于低压缩性土, 规范法的计算值远大于实际值.同时, 利用P-Z模型计算时, 需考虑大量的模型参数, 导致误差累计.因此, 本文针对模型参数展开敏感性分析.
2 参数敏感性分析基于参数的基准值, 分别从弹性参数、塑性模量参数、临界状态参数3个方面讨论.分别将参数从-40%变化到40%, 计算土层厚0.5 m的变形值, 该变形值可视土层最终变形指标.
2.1 弹性参数分析弹性模量E为材料抵抗弹性变形的指标, 其值越大, 抵抗变形的能力越强.图 2a为E变化时, 对应的竖向位移, 随E的增大, 竖向位移减小, 且E较小时, 竖向位移变化速率相对更快.当E在-40%~40%变化时, 位移变化率在-28.46%~66.05%.可见, 其值的选取对竖向位移变化率高度敏感, 取值时应倍加谨慎.
图 2(Fig. 2)

图 2 弹性参数与位移的关系Fig.2 Relation between elastic parameters and displacement (a)—参数E与位移的关系；(b)—参数μ与位移的关系.

泊松比μ是反映土体侧向与竖向变形比值的参数.由图 2b知, 土体竖向变形随μ的增大而减小, 当μ→0.5时, 变形趋向于0，其变化呈现明显的线性特征.相对于参数E, 其敏感性小得多.当μ在-40%~40%变化时, 其竖向位移变化率在-5.75%~5.15%.
2.2 临界状态参数分析定义临界状态的2个参数为相对密实度D_r和内摩擦角φ.P-Z模型中, M_g为相应临界状态的模型参数.由式(2)知, M_g为内摩擦角φ和Lode角θ的函数, M_f可由D_r和M_g估计.
相对密实度D_r反映了土体松散度.其值越小, 土体越松散, 受载作用下变形越明显.图 3a为D_r变化时竖向位移的值, 当D_r在-40%~40%变化时, 竖向位移变化率在-2.52%~1.94%.
图 3(Fig. 3)

图 3 临界状态参数与位移的关系Fig.3 Relation between CSL parameters and displacement (a)—参数Dr与位移的关系；(b)—参数φ与位移的关系.

图 3b为内摩擦角φ变化时对竖向位移的影响.随φ的增大, 变形呈线性减小, 其变化符合一般规律, 呈现与参数D_r相似的变化趋势, 但其敏感程度却较D_r大.当φ在-40%~40%变化时, 竖向位移变化率在-8.63%~8.61%.
2.3 塑性模量参数分析4个与塑性模量相关的参数, 分别为ξ, H₀, β₀, β₁, 反映了土体在等向压缩时的塑性模量.
图 4a为累计偏应变ξ对竖向位移的影响.式(8)中, ξ 影响偏应变硬化函数H_s, 随ξ的增大, 竖向位移呈非线性减小, 且变化幅度较小.当ξ在-40%~40%变化时, 其竖向位移变化率在-0.12%~0.06%.文献[9]针对ξ=0.6，0.7，0.8时分别进行砂土实验数据的拟合, 结果表明, 当ξ=0.7时拟合结果更为贴近.图 4b为H₀变化时, 竖向位移的变化曲线, H₀增大, 其竖向位移呈线性减小趋势.其敏感性相对ξ大, 但总体较小.当H₀在-40%~40%变化时, 其竖向位移变化率在-0.23%~0.22%.
图 4(Fig. 4)

图 4 塑性参数与位移的关系Fig.4 Relation between plastic parameters and displacement (a)—参数ξ与位移的关系；(b)—参数H0与位移的关系; (c)—参数β0与位移的关系；(d)—参数β1与位移的关系.

图 4c为参数β₀变化时, 竖向位移变化趋势.文献[10]指出, β₀的取值范围在1.5~5.0.由图知, β₀≤1.5时, 其竖向变形不变, 随β₀增加, 竖向位移增大, 且变化率逐渐增大.β₀=5时, 竖向位移增大0.11%.当β₀在-40%~40%变化时, 竖向位移变化率在-0.01%~0.03%.图 4d为参数β₁变化时, 竖向位移变化趋势.文献[10]指出, β₁的取值范围在0.1~0.2.随β₁的增加, 竖向位移减小, 其变化曲线与参数H₀相似, 当β₁在-40%~40%变化时, 竖向位移变化率在-0.06%~0.06%.由式(7)~式(8)知, 参数ξ, H₀, β₀, β₁通过改变H_L实现对C^ep的影响, 其参数敏感性均较低.
为直观比较各参数的敏感性, 将各参数变化率与位移变化率同时比较.图 5a为8个参数变化时, 对位移变化率的影响.由图知, 参数E的敏感性最大, 最大变化率可达66.05%.其次为φ和μ, 变化率均在±10%以内.其余5个参数的敏感性均很小, 变化率在±3%以内.将敏感性较低的7个参数表示在图 5b中, 需要指出, 与塑性模量相关的参数, 竖向位移变化率均低于±1%.当缺少实测数据时, 可根据实验资料取值.
图 5(Fig. 5)

图 5 各参数与位移的关系Fig.5 Relations between parameters and displacement

依据参数敏感性分析结果, 由于定义塑性模量的相关参数对地基变形的敏感性很小, 本文采用地质资料统计和经验取值.利用平面网格划分技术对大型基础下层状砂土地基变形进行数值计算, 将计算值与规范方法计算值对比分析.
3 计算算例为了进一步说明参数敏感性分析结果在计算中的合理性, 选择《建筑地基技术规范》中给出的计算算例为计算对象, 计算模型如图 6所示.
图 6(Fig. 6)

图 6 变形计算必要数据示意图Fig.6 Necessary datas of calculation deformation

基底表面加载方式为p₀=189.4 kPa均布荷载.基底平面尺寸长L=110.2m, 宽B=17.6m, 基底划分20×125个子域.计算参数如表 1所示, 依据参数敏感性分析结果, 塑性模量相关参数依据文献[9-10]取值.
表 1(Table 1)

表 1 地基土的主要物理力学参数Table 1 Main mechanical parameters of foundation soil 土层 E 泊松比μ 内摩擦角 黏聚力 密实度 MPa (°) kPa % 粗砂 25 0.30 35.0 1.3 39.2 砾砂 30 0.25 37.0 2.5 43.8 卵石 40 0.20 39.5 5.0 39.7 粉质黏土 5 0.35 6.7 32.1

表 1 地基土的主要物理力学参数 Table 1 Main mechanical parameters of foundation soil

运用基底平面网格划分技术, 结合基于P-Z本构模型的压缩模量的计算方法, 对该计算算例进行变形计算, 其计算结果如图 7所示.
图 7(Fig. 7)

图 7 地基土的总变形Fig.7 Total vertical deformation of foundation soil (a)—y轴方向截面曲线；(b)—x轴方向截面曲线；(c)—变形空间分布图.

由图 7知, 即使基底表面作用均匀分布荷载, 其最终变形在空间上非均匀分布, 这是由于附加应力在水平和竖向两个方向同时传递和叠加造成的.采用规范法计算, 将计算结果对比, 如表 2所示.
表 2(Table 2)

表 2 计算结果比较Table 2 Comparison of calculating values 土层 土层厚度 本文模型计算值 规范法计算值 误差 m mm mm % 粗砂 4.6 7.323 7.260 0.063 砾砂 4.4 4.670 4.886 0.216 卵石 19.2 10.732 9.310 1.422 粉质黏土 0.8 13.425 14.697 1.272

表 2 计算结果比较 Table 2 Comparison of calculating values

粗砂层的计算差值为0.063 mm，误差为0.86%, 计算结果非常接近.砾砂层的计算差值为0.216 mm，误差为4.42%, 误差较粗砂层大, 但仍然较小.卵石层的计算差值为1.422 mm，误差为15.27%, 误差较大, 原因在于该层深度为19.2 m, 造成误差的持续积累.粉质黏土层的计算差值为1.272 mm，误差为8.67%, 误差原因在于用砂土模型计算黏土的变形.其计算结果的对比曲线及误差如图 8所示.
图 8(Fig. 8)

图 8 各层地基土竖向变形值Fig.8 Vertical deformation of each foundation soil

4 结论1) 依据P-Z模型本构理论, 给出了弹塑性矩阵的计算方法, 并推导了可用于地基土变形计算的压缩模量计算方法.
2) 对模型参数进行敏感性分析, 结果表明:弹性模量E的敏感性最大, 其次为内摩擦角φ和泊松比μ.而与塑性模量相关的参数, 敏感性很低.在计算分析缺少实测数据时, 可依据资料取值, 并通过算例证明其可行性.
3) 依据敏感性分析结果, 简化了模型参数的标定过程.计算结果表明:本文给出的计算方法与规范法相比, 计算结果相近.且本文的计算方法克服了规范法单点计算的缺陷, 计算结果更加直观.
参考文献

[1]	Zienkiewicz O C, Mroz Z. Generalized plasticity formulation and applications to geomechanics[J].Mechanics of Engineering Materials, 1984, 44(3): 655–680.
[2]	Li H E, He Y J, Fan G Y, et al. Recent developments of generalized plasticity model for saturated and unsaturated soils[J]. Water Science and Engineering, 2011, 4(3): 329–344.
[3]	Pastor M, Zienkiewicz O C. A generalized plasticity, hierarchical model for sand under monotonic and cyclic loading[C]//2nd International Symposium on Numerical Models in Geomechanic. Ghent, 1986: 131-137.
[4]	Pastor M, Zienkiewicz O C, Chan A H C. General plasticity and the modelling of soil behavior[J].International Journal for Numerical and Analytical Methods in Geomechanics, 1990, 14: 151–190.DOI:10.1002/(ISSN)1096-9853
[5]	Iraji A, Farzaneh O, Hosseininia E S. A modification to dense sand dynamic simulation capability of Pastor-Zienkiewicz-Chan model[J].Acta Geotechnica, 2014, 9: 343–353.DOI:10.1007/s11440-013-0291-y
[6]	Pastor M, Zienkiewicz O C, Leung K H. Simple model for transient soil loading in earthquake analysis Ⅱ:non-associative models for sands[J].International Journal for Numerical and Analytical Methods in Geomechanics, 1985, 9(5): 477–498.DOI:10.1002/(ISSN)1096-9853
[7]	Akhaveissy A H. Analysis of tunnel and super structures for excavation[J].Scientia Iranica, 2011, 18(1): 1–8.DOI:10.1016/j.scient.2011.03.001
[8]	李纯, 朱浮生, 张学元. 层状砂土地基静力变形形状测试[J].东北大学学报(自然科学版), 2012, 33(7): 1038–1042. ( Li Chun, Zhu Fu-sheng, Zhang Xue-yuan. Deformation properties of layered sandy foundation soil tested with static additional load[J].Journal of Northeastern University(Natural Science), 2012, 33(7): 1038–1042.)
[9]	李宏恩, 李铮, 徐海峰, 等. Pastor-Zienkiewicz状态相关本构模型及其参数确定方法研究[J].岩土力学, 2016, 37(6): 1623–1632. ( Li Hong-en, Li Zheng, Xu Hai-feng, et al. Pastor-Zienkiewicz based constitutive model and determination method of its state parameters[J].Rock and Soil Mechanics, 2016, 37(6): 1623–1632.)
[10]	Farzaneh O, Iraji A. Two-phase model for nonlinear dynamic simulation of reinforced soil walls based on a modified Pastor-Zienkiewicz-Chan model for granular soil[J].Journal of Engineering Mechanics, 2016, 142(2): 1–17.