动载作用下土工格栅加筋土计算模型的研究

动载作用下土工格栅加筋土计算模型的研究

张兴强¹，邓卫东²，陈文金³，闫澍旺¹
(天津大学建工学院岩土工程研究所；2.交通部重庆公路科学研究所.
海洋石油股份有限公司)

摘要：本文在前人工作的基础上提出了在动力荷载作用下土工格栅与土相互作用的有限元分析模型，该模型采用非线性弹簧-阻尼延迟器-质量块系统模拟土工格栅纵肋、横肋与土动力相互作用，采用集中刚度法形成总刚，可用来研究加筋土在动载作用下的变形特性。非线性弹簧的刚度、阻尼延迟器的阻尼比可以通过振动台试验获得。基于该模型编制了加筋土的动力有限元分析程序，有限元分析结果与实测值一致，表明了计算模型和分析程序的正确性。

关键词：加筋土；土工格栅与土动力相互作用；接触面变形模式；接触面模型

收稿日期：2000-05-25
作者简介：张兴强(1970-)，男，安徽人，天津大学博士生，主要研究方向：格栅与土相互作用。

　　在土工格栅加筋土结构的动力分析时，如何模拟土工格栅与土的相互作用是一关键问题。目前一些学者对动载作用下加筋土的计算模型进行了探索：Bray等(1995,1996)进行对比性研究表明在某种条件下接触面一维模型和二维模型有类似的结果。但是一维模型在实际计算中有一些缺点，如：模拟摩擦接触面的不稳定性。在当前数值计算中，解决此问题的方法一般是在地震响应计算中，忽略土工合成材料的作用。Kavazanjian和Matasovic(1995)研究表明：在计算中忽略土工合成材料的作用将过高地估计加速度。Muthucumarasamy等(1992)基于土体与加筋材料之间完全连结的假设，分别采用等效迭代线弹性法和弹性增量法来分析动载作用下加筋土的响应问题。Yegian等(1998)采用等效弹簧-阻尼延迟器-质量块系统来模拟动载作用下土工合成材料与土的接触面特性，并将土工合成材料与土的接触面、土工合成材料等效为一土层来研究加筋土的动力特性。

　　虽然一些研究者从实验和理论上对加筋土在动载作用下的响应问题进行了探索，但是由于加筋材料与土的接触面结构形式、在动力荷载作用下的加筋材料作用机理的复杂性和试验手段等方面的原因，对加筋土计算方法的研究还很缺乏。Muthucumarasamy等(1992)对美国加州大学UCLA加筋挡土墙的计算结果也证实忽略加筋材料与土之间的相互作用将过高地估计加速度，这与Kavazanjian和Matasovic(1995)得出的结论是一致的，因此分析这种问题时应该采用适当的接触面模型来模拟加筋材料与土之间的动力相互作用。本文基于前人工作的基础对土工格栅与土动力相互作用机理进行了研究，提出了一个用于模拟其接触面变形特性的计算模型，并通过算例验证该模型的合理性。

1 接触面动力模型

　　埋在土中的格栅在动载作用下，其锚固能力是通过摩阻力和支承阻力获得的。沿着拉拔荷载方向的纵肋提供摩阻力，垂直方向的横肋提供摩阻力和支承阻力。纵肋的身长率，横肋抗挠刚度会影响摩阻力和支承阻力的发挥。最终的抗拔力(Koerner R.M,1989)为

F(t)=F₁+F₂+F₃

(1)

式中：F₁为纵肋摩阻力；F₂为横肋摩阻力；F₃为横肋支承力。

　　动载作用下土工格栅与土之间将发生相对位移，从整体上来说是不连续的，采用有限元计算须设置接触面单元，将土工格栅与土之间的接触面内设想由无数非线性弹簧、粘滞器、质量块相联系，即构成非线性弹簧-粘滞器单元模型(图1).为了避免无厚度单元在计算中出现两侧单元重叠的现象，采用Desai薄单元作为土工格栅与土之间的接触面单元。

图1 接触面模型示意图

2.1 土工格栅纵肋与土动力相互作用的模拟 采用1个等效非线性弹簧、1个粘滞器和1个质量块来模拟土工格栅纵肋与土之间的接触面动力摩擦特性。接触面的等效质量、刚度和阻尼分别依赖于接触面的土层厚度、土的模量(土体弹性压缩模量和等效剪切模量)和阻尼比，因此可以用弹簧、粘滞器和质量块来模拟接触面。利用振动台试验(Yegian,1998)得到加速度与接触面滑动位移之间的滞回曲线(图2)来计算接触面的剪切刚度和阻尼比。

图2 剪切刚度和阻尼比的确定

Ke=F_m/S_a

(2)

式中：K_e为接触面的等效剪切刚度；F_m=MA_t=W/gA_t，为通过接触面传递的最大剪力；A_t为接触面传递的最大加速度；M和W分别为接触面的质量和重量；S_a为接触面滑动位移。

　　根据Desai等(1985)提出的刚度和概念

K=F_m/A_c/S_a

(3)

式中：A_c为土工格栅纵肋与土的接触面积。则

K=WA_t/A_cgS_a=A_tσ/gS_a

(4)

式中：σ为作用接触面上的正应力；K为接触面的剪切刚度。

G_e=KH_e

(5)

式中：G_e为等效剪切模量；H_e为接触面厚度。

　　根据接触面加速度与滑动位移滞回曲线(图2)，接触面阻尼比可以通过下式来确定

β_h=A₁/4πA₂

(6)

式中：A₁为循环荷载每作用一次所消散的能量；A₂为弹性应变能。

2.2 土工格栅横肋与土相互作用的模拟 有资料表明抗挠刚度较大的横肋在循环荷载作用下成双曲线型(Raju,1998).其荷载-变位关系可用下式表示。

图3 土工格栅横肋与土接触面模型

P/k_s=f(d/L)

(7)

式中：p=F₂+F₃为作用在横肋上总的力；f()d/L)为变位比d/L的显函数；d为横肋中间部位的变化；L为横肋的长度；k_s为横肋的荷载张拉模量。

　　横肋上各点处的变位是不同的，变位值的大小反映横肋抗挠刚度的大小，横肋上各点处的受力特性也不同。为了确切描述可将横肋看成由许多结点组成，考虑到接触面的厚度效应，可以采用四结点单元来模拟每个结点上的力和变形关系。单元由两组沿法向和切向方向的非线性弹簧、粘滞器、质量块组成。弹簧的劲度系数分别为k_n和k_st。k_n与土体的模量有关，k_st可根据横肋模拟锚定板试验确定。资料表明其应力-应变关系曲线接近双曲线(Ragui,1994)，其切线剪切劲度系数为

k_st=(1-R_bq/σ_nN_q)²k₁p_a(σ_n/p_a)ⁿ

(8)

式中：R_b为破坏比，即实际所具有的支承能力与破坏时的应力之比；q为支承阻力；N_q为支承反力系数，可表示：N_q=tan(π/4+φ′/2)exp[(π/2+φ′)tanφ′],其中，φ′为土工格栅孔径内与土的摩擦角。

　　则：土工格栅横肋与土接触面的剪切模量为：

G_he=k_stH_e

(9)

2.3 土工格栅与土接触面的有限元格式 Desai薄单元的本构关系可表示为

(10)

式中：[D_sn]和[D_ns]为考虑正应力与剪应力耦合效应的分量，在这里取为零。[D_ss]为剪应力分量；[D_nn]为正应力分量，可以用显示表示为：

[D_nn]=λ₁[D_nn]_i+λ₂[D_nn]_g+λ₃[D_nn]_st

(11)

式中：脚标i,g和st分别为接触区材料，岩土材料和土工格栅；λ₁,λ₂和λ₃为系数，对于动力问题，取λ₁=0.75，λ₂=0.25，λ₃=0(Zaman,1985).

　　接触面采用弹塑性本构模型，在塑性屈服之前服从弹性应力-应变关系。塑性屈服后，接触面应力增量dσ和应变增量dε之前的关系为：

dσ=[D_ep]dε

(12)

式中：

(13)

其中[D]为弹性矩阵

(14)

式中：

为纵肋与土相互作用的弹性矩阵；

为横肋与土相互作用的弹性矩阵；E为土体的弹性压缩模量，ν为土体的泊松比。A是反映硬化特性的一个变量，与硬化参数H的选取有关。

H=E_T/1-E_T/E

(15)

其中：E_T为弹塑性切线模量。f和g分别为屈服准则和流动准则，采用Mohr-Coulomb屈服准则，流动准则采用相关联流动准则。

f=σ_msinφ+cosθ-/sinφsinθ-ccosφ

(16)

式中：σ_m=I₁/3;

;θ=1/3arcsin(-3

J′3/2

)为Lode角；-π/6<θ<π/6.I₁、J′₂、J′₃分别为第一应力不变量、第二应力不变量和第三应力不变量；c为粘聚力；φ为土工格栅与土接触面的动摩擦角。

[K]=[K]_z+[K]_h=∫S[B]^T[D]_ep[B]dS

(17)

式中：[K]_z为土工格栅纵肋与土接触面的刚度矩阵；[K]_h为土工格栅横肋与土接触面的刚度矩阵；[B]为接触面的几何矩阵。阻尼矩阵采用Relaigh阻尼。

2.4 土工格栅和土体单元的有限元模拟 土工格栅被假设为只能沿轴向变形的一维杆单元，位移函数是线性的，用增量形式表示的单元刚度矩阵为：

(18)

式中：E_ft为单元的切线弹性模量；A为格栅的横切面积；L为格栅单元的长度；s=sinα，c=cosα；α为单元与X轴的夹角。土工格栅的单元质量采用集中质量。

　　土体单元采用等参四结点单元，其本构关系采用经典弹塑性模型，屈服准则为Mohr-Coulomb屈服准则，流动准则采用相关联流动准则。

2.5 土工格栅与土接触面变形模式的模拟 在动载作用下土工格栅与土之间的接触面存在着相对运动模式，一般情况下接触面有4种变形模式：(1)粘结；(2)错动；(3)脱离；(4)嵌入。图4显示了二维理想化的这4种接触面变形模式。

图4 (a)粘结(b)错动(c)脱离(d)嵌入接触面的变形模式

　　假设接触面最初均以粘结形式存在，基于动力平衡方程，已知t时刻的(σ_nn)^t、(τ)^t、(u_n)^t,需要求在t+Δt时刻时的(σ_nn)^t+Δt、(τ)^t+Δt、(u_n)^t+Δt和Δt时刻内的增加量Δσ_nn、Δτ、Δu_n.其中：u_n为接触面平均正向相对位移；σ_nn为接触面的正应力；τ为接触面的剪应力。

(σ_nn)^t+Δt={(σ_nn)^t+Δσ_nn}	(19)
F_s=(sgn){(τ)^t+Δτ}-[c+(σ_nn)^t+Δσ_nn}tanφ	(20)

式中：压应力为正；c为接触面的粘聚力；φ为接触面的摩擦角；

(21)

接触面的变形可以通过表1进行判断。

表1 接触面变形模式


粘结/嵌入	错动	脱离

(σnn)^t+Δt>0并且F_s<0	(σ_nn)^t+Δt>0并且F_s≥0	(σ_nn)^t+Δt<0

　　当接触面处于脱离状态时，正应力和剪应力都不能通过接触面传递。这时有{(σ_nn)^t+Δσ_nn}=R₁且{(τ)^t+Δτ}=R₂,将残余应力R₁、R₂通过系统迭代修正，转换成自平衡结点力，使R₁≌0，R₂≌0.

　　对于薄单元来说，一般不会出现接触面相邻结点互相贯穿现象。但是如果这种情况出现，采取以下的方法：(1)通过将一个任意小的值ε和单元的厚度t′赋给出现贯穿现象相邻结点，使相邻结点分开；(2)贯穿位移的大小通过计算接触面结点的正向位移来判断，被压缩的正向位移为贯穿位移与被分开位移(即εt′)的差；(3)将压缩位移{Δq_s}同下式转化为等效修正结点荷载矢量{ΔQ}

{ΔQ}=∫_V{D}^t_i[B]{Δq_s}dV

(22)

式中：[D]^t_i为t时刻接触面的本构矩阵；[B]为接触面的几何矩阵；V为接触面的单元体积。如果{Δq_s}大于允许误差，继续迭代，直到{Δq_s}小于允许误差。在迭代过程中，只有残余荷载矢量改变，总刚、质量矩阵、阻尼矩阵在每一个时间步中只形成一次，结点速度和加速度在迭代过程不改变。

3 算例分析

　　基于以上接触面动力计算模型，编制了用于计算平面应变问题的弹塑性动力加筋有限元G-D程序。对土体采用弹塑性模型，土体单?捎?结点等单元，计算中需跟踪每一个高斯点的弹塑性应力路径，以便确定相应的弹塑性矩阵；用弹性杆单元来模拟土工格栅；土工格栅与土体之间设接触面单元，接触面单?捎肈esai薄单元，在程序中分别考虑了土工格栅的纵肋与横肋的作用效应，用集中刚度法来形成接触面的总刚。

　　本文采用Richardson等(1977)报道的加州大学加筋挡土墙(UCLA)承受随机激励的试验数据和Muthucumarasamy Yogendrakumar等(1992)对该试验的计算模拟结果作为G\|D程序验证依据。UCLA挡土墙高6.1m，采用4.88m长的钢条加筋，每个钢条相距0.76m，T3R1冲击机位于挡土墙前9.46m，冲击荷载持续时间为0.75s.在选定的位置测量加速度的时间历程。图5中表示在冲击荷载作用下，测量加速度和加筋层所承受动力大小的位置：加速度测量仪用ACC表示，加筋层所受动力大小测量仪用TF表示。ACC1被放置在加筋挡土墙的顶部，距墙面2.44m，ACC2被放置在高3.43m的墙面处，TF1～TF5被放置在钢加筋条上，其放置高度分别为5.72m、4.96m、4.19m、2.67m和0.38m.图6为冲击加速度谱。


图5 UCLA加筋挡土试验墙有限元模型	图6 UCLA加筋挡土试验墙激励加速度

　　Muthucumarasamy Yogendrakumar等(1992)计算时，土体单元的动力剪切模量G_max采用Seed和Idriss(1970)推荐的公式计算：

G_max=21.7K_2maxP_a(σ_m/P_a)^0.5

(27)

其中：σ_m为有效正应力；P_a为大气压力；K_2max为计算常数。

表2 UCLA模型试验土体单元和加筋钢条材料参数


材料类型	K_2max	摩擦角 /度	密度 /(Kn/m³)	泊松比	厚度 /mm	宽度 /mm	弹性模量 /GPa	屈服强度 Mpa	最终强度 MPa

填充土	56	38	18.84	0.35	-	-	-	-	-
粉土	50	36	17.28	0.35	-	-	-	-	-
砂石	70	38	18.84	0.35	-	-	-	-	-
钢条	--	--	-	0.28	3	60	197	255	345

　　图9(a)和(b)分别为UCLA试验墙在ACC1和ACC2处采用G-D程序计算的水平加速度值、Muthucumarasamy等人的计算值与UCLA加筋试验挡土墙实测值的对比分析。结果表明：采用G-D程序计算ACC1处的水平加速度与实测值吻合得较Muthucumarasamy等人的计算值好。图9(a)表明：用G-D程序预测值与实测值在0～0.15s之间有一定的差别，但在0.15～0.75s之间，预测值与实测值的峰值加速度吻合得很好。图9(b)表明：预测值与实测值在0.10～0.35s和0.50～0.75s之间铀俣壤毯芪呛希?～0.10s和0.35～0.50s加速度历程有一定的差别。通过计算与实测值的比较表明：采用G-D程序能较好地预测加筋体的加速度的历程，G-D程序可以很好用于加筋结构的动力计算。


(a)ACC1处的水平加速度	(b)ACC2处的水平加速度
图7 测量与计算的水平加速度比较

4 结论

　　(1)本文在前人工作的基础上，提出了解决平面问题的土工格栅与土相互作用的动力有限元分析模型，该模型采用不同的本构模型分别模拟土工格栅纵肋、横肋与土动力相互作用，具有计算参数少，易于计算机模拟。可用来研究加筋土在动载作用下的变形特性，有助于理解土工格栅与土动力相互作用机理；(2)基于该模型编制了加筋土的动力有限元分析程序G-D，通过与加州大学UCLA加筋挡土墙实测值对比分析表明了程序的正确性，模型的合理性，可真实地模拟加筋材料与土的动力相互作用及加筋土在动载作用下的变形特性；(3)该模型不但可以模拟土工格栅与土的动力相互作用，还可以用于模拟其他加筋材料与土的动力相互作用，如土工织物、钢筋等；基于该模型编制的程序G-D可用于解决道路、铁路、水利等工程中加筋结构的动力计算，具有广泛的适用性和应用面。

参考文献：

[1] 闫澍旺，Ben Barr.土工格栅与土相互作用的有限元分析[J]。岩土工程学报，1997，(6)：56～61.

[2] 全国第四届土工合成材料学术会议论文集[C]。中国土工合成材料工程协会，1996.

[3] Bray, et al.Seismic Stability Procedures for Solid-waste Landfills[J]。J.Geotech.Engrg.,ASCE,1995,121(2):131～151.

[4] Bray, et al.Closure to ‘Seismic Stability Procedures for Solid-waste Landfills’[J]。J.Geotech.Engrg.,ASCE,1996,122(1):952～954.

[5] Kavazanjian, et al.Seismic Analysis of Solid Waste Landfills[C]。Geoenvironment 2000,Geotech.Spec.Publ.No.46,ASCE,N.Y.,1995.

[6] Yegian Harb,Kadakal.Dynamic Response Analysis Procedure for Landfills with Geosynthetic Liners[J]。J.Geotech.Engr.ASCE,1998:1027～1033.

[7] Ragui F Wilson-Fahmy,Robert M Koerner,Leonard J Sansone.Experimental Behavior of Polymeric Geogrids in Pullout[J],J.Geotech.Engineering,ASCE,1994(4):661～677.

[8] Koerner R M,Wayne M H,Carroll R G,Analytic Behavior of Geogrid Anchorage[C]。Proc.Geosynthetics 89 Conf,IFAI:525～536.

[9] Yegian M K,Harb J N,Kadakal U,Dynamic Response Analysis Procedure for Landfills with Geosynthetic Liners[J]。J.Geotech.Engr.ASCE,1998:1027～1033.

[10] Desai C S,Drumm E C,Zaman M M.Cyclic Testing and Modeling of Interfaces[J]。J.Geotech.Engrg.,ASCE,1985(6):793～815.

[11] Muthucumarasamy Yogendrakumar,Richard J.Bathurst,W.D.Liam Finn,Dynamic Response Analysis of Reinforced-Retaining Wall[J]。J.Geotech.Engrg.,ASCE,1992(8):1158～1167.

[12] Md.Musharraf-uz Zaman,Chandrakent S Kesai,Eric C Drumm.Interface Model for Dynamic Soil-Structure Interaction[J]。J.Geotech.Engrg.,ASCE,1994(9):1257～1273.

[13] Seed H B,Idriss I M.Soil Moduli and Damping Factors for Dynamic Response Analysis[C]。Report No.EERC70-10,Earthquake Engrg.Res.Ctr.,Univ of Califormia,Berkeley,Calif.,Dec,1970.

[14] D M Raju,R J Fannin.Load-strain-displacement Response of Geosynthetics in Monotonic and Cyclic Pullout[J]。Can.Geotech.1998:183～193.