砂土材料状态相关临界状态各向异性模型

李学丰，袁琪，王 兴

（宁夏大学 固体力学研究所，银川 750021）

砂土的力学响应与其材料状态紧密相关［1－3］。传统本构理论往往忽略了材料状态对应力－应变关系的影响，将不同孔隙比同一种砂视为不同的材料，使得同一种砂必须使用多组材料参数才能描述其力学特性，相对繁杂，存在一定的局限性。材料状态相关临界状态理论则能有效克服这些局限性［4－7］。

砂土材料状态不仅要考虑其相对密度、有效围压，还要考虑其细观各向异性和应力加载方式。这几个因素的相互影响造成了砂土强度－变形的复杂性。在这个复杂变化过程中实质变化的是材料细观状态，因此，结合其细观状态研究力学特性逐步引起关注［8－14］。Mooney等［8］用摄影技术研究剪切带的演化规律时发现了临界状态线在e-p平面不唯一。Yoshimine［9］、Nakata等［10］的大量试验逐步也证实了临界状态线不唯一，其不唯一的根本原因是砂土初始各向异性的存在。Dafalias［11］最早应用该特性建立粘性土的本构关系。对于砂土，Wan等［12］、Li等［13］和Dafalias等［14］许多学者都根据砂土临界状态的不唯一性建立了本构关系。考虑临界状态不唯一建立本构的关键是对土体细微观结构与材料状态关系的合理描述。最近，Li等［15］、Zhao等［16］、Abouzar等［17］对砂土细观组构对临界状态不唯一关系进行了研究，为合理描述砂土各向异性力学特性奠定了基础，也为材料状态相关临界状态理论提供了新思路。

为合理描述各向异性砂土的强度－变形特性，在黄茂松等［6－7］模型及砂土各向异性分析［18］的基础上，在材料状态相关临界状态理论的框架下，基于砂土临界状态的不唯一，考虑其细观组构对材料临界状态、峰值状态、相变状态、剪胀性及其硬化规律的影响，建立各向异性砂土的本构模型，期望逐步实现砂土细观结构对宏观力学特性影响的定量描述。

1 砂土各向异性对材料状态的影响

临界状态理论基于粘性土试验提出，很好描述了粘性土的临界特性。对于砂土，许多试验表明其临界应力比与孔隙比关系不唯一，根本原因是其各向异性的存在。

1.1 砂土各向异性临界状态描述

试验表明砂土的临界状态线在e-lnp平面不是直线，这与粘性土有很大的差别。Li等［4］通过对砂土试验的归一化得到表达式为

式中：eΓ、λc和ζ为决定临界状态线的材料参数；p为当前平均应力；pa为大气压力，如图1。式（1）是围压的函数，不同围压下，各向同性砂有唯一的临界状态线，不同的相变线。临界状态是不同围压下相变状态的临界值［2］。

图1 Toyoura砂相变状态线与临界状态线 （文献［2］）Fig.1 Phase transformation and Critical-state lines of Toyoura

当各向异性存在时，砂土临界状态不唯一。考虑各向异性对砂土临界状态线影响的时，式（1）中eΓ、λc和ζ都要受各向异性的影响。然而，同时考虑3个参数的影响，将会极大增加本构描述的复杂性。Poulo等［19］研究认为砂土颗粒排列主要对其临界状态线斜率影响最大，基于这个研究成果Li等［13］提出的临界状态线只考虑eΓ影响，Yang等［20］通过细观试验又做了相应的改进。Abouzar等［17］针对组构对π平面临界状态特性的研究很好解决了这个问题。基于以上分析和笔者的各向异性破坏特性［19］方面的工作，将文献［19］新定义的状态变量引入式（1）建立了临界状态线

式中：t是状态参数，式（2）在ec－（p／p）ζ平面内变为相互平行的直线，t可以由常规三轴压缩、伸长试验由式（2）直线截距差得到，本文使用Yang等［21］结合Toyoura砂细观定量检测的三轴试验确定；A为是笔者等在文献［19］中定义一个新的各向异性状态变量

图2 状态变量A与中主应力系数的关系Fig.2 Anisotropic state variable VS Intermediate principal stress coefficient.

图3 各向异性临界状态线Fig.3 Anisotropic Critical-state lines

式（2）是A的函数，A受应力状态、各向异性程度和应力与土体组构角度3个因素的影响，因此，随3个因素的变化会导致式（2）的变化，直接导致临界状态线的不唯一。如图3所示，其它材料参数保持不变时，材料为各向同性时，A恒等于零，A对式（2）没有影响；当a＝0.2时（A的变化如图2），临界状态线随中主应力系数b的变化在一定区间变化，这个区间为Verdugo等［2］通过试验得到Toyoura砂的准临界状态和临界状态区间。

1.2各向异性临界状态线的特点

现有试验成果表明各向异性对砂土临界状态的影响主要表现在e-（p／pa）ζ平面、π平面上和摩尔圆平面上的影响规律。

1）根据式（2）可以直观地看出，在e-（p／pa）ζ平面上，临界状态线是各向异性状态变量的函数，其斜率受砂土各向异性程度、应力状态以及组构与应力的几何关系3个主要因素的影响，临界状态线是相互平行的直线。

2）根据材料状态相关临界状态理论，考虑各向异性后，临界状态应力比的一般表达式为

式中：Mcs为常规压缩三轴试验确定的临界状态应力比；ψ（A）是考虑各向异性的状态参量。根据Been等［1］定义状态参数，增加了各向异性对状态参数的影响，其表达式为

式中：e为当前孔隙比为ec（A）为式（2）定义的各向异性临界状态线。土体达到临界状态时，必须满足孔隙比和应力比分别同时达到临界状态孔隙比和临界应力比两个基本条件。则式（5）可以表示为

式（6）和文献［19］破坏准则形式一致，因此，在π平面上表现的各向异性特性和各向异性破坏规律相似。如图4所示，在π平面上，临界状态线的形状和位置都受各向异性的影响。各向异性时，其中心自然偏离了静水压力轴，而且各向异性越大，偏离静水压力轴中心越远，同时其形状改变也越大。在常规三轴压缩点（图中与σ1轴的交点）各向同性和各向异性值相同，定义A采用了基准点的思想很好描述了这个特性，这与文献［21］的试验结论一致。图4中峰值线和相变状态线也具有以上特点。结合图1，临界状态是相变状态的极限表现，合理描述临界状态的各向异性对砂土本构的峰值状态和相变状态有重要作用。

图4 π平面上的状态线Fig.4 State line on the octahedral plane.

3）在同一孔隙比和应力状态的条件，式（2）定义可以较好描述随砂土沉积角度变化土体临界状态值的试验规律。如图5所示，当砂土的沉积面方向和小主应力轴为90°时，临界状态强度最大，0°时到最小值，随沉积面角度变化临界强度单调递减。这个规律和Oda等［22］等的真三轴试验结果一致。

综上所述，临界状态线能较好描述各种空间、应力状态和各向异性程度砂土临界状态变化规律。

图5 摩尔圆平面上的各向异性临界状态线Fig.5 Anisotropic Critical-state lines on the Mohr plane.

2 模型框架

2.1 屈服面

屈服面采用式（7）所示形式。

式中：M为硬化内变量；g（θσ）为罗德角的插值函数，采用 William等［23］建议的椭圆形插值函数。

2.2 剪胀性

剪胀性不仅与孔隙比和围压有关，同时，也受砂土初始各向异性的影响。参考Gajo等［24］提出的剪胀方程并考虑各向异性的影响，将方程修正为

式中：d0为模型参数；Md为状态转换应力比。试验结果表明：Md不是常量，而是与材料状态相关。Md参考Li等［13］建议的形式，本文Md三维空间的表达式为

式中：kd模型参数，同样从式（3）可以看出Md是各向异性状态变量的函数。如图1和图5所示，Md通过状态参数将材料的当前状态和临界状态紧密联系起来，式（8）描述的剪胀性也是A的函数。

当砂土为各向同性时，临界状态性唯一，其剪胀和剪缩关系也唯一。如图6（a）所示，当前孔隙比小于临界应力比砂土特性表现为剪胀，反之为剪缩，其关系唯一；当砂土为各向异性时，由于临界状态线不唯一，其剪胀和剪缩特性不唯一。如图6（b）所示，考虑各向异性后，临界状态线受各向异性大小、应力状态和应力与组构几何关系3个关键因素的影响，临界状态线不唯一，为相互平行的直线。对于图6（a）中状态A′和状态B′，在各向同性条件下分别表现为减缩和剪胀的点，然而图6（b）则分别表现为剪胀和减缩。可以看出，当考虑各向异性对材料状态的影响后，剪胀性关系的判断将更加全面、细致。当ψ（A）＞0时表示砂土处于松散状态，受力会发生剪缩现象，ψ（A）＜0则表示砂土处于密实状态，受力会发生剪胀现象，而且｜ψ（A）｜大小可以反映剪胀性的大小。ψ（A）通过各向异性状态变量，也将砂土的各向异性程度、应力状态和砂土组构与应力之间的关系对砂土剪胀性的影响进行了量化。因此，引入各向异性状态变量，既丰富了砂土材料状态描述，又丰富了其剪胀性描述。

图6 状态参量与剪胀的关系Fig.6 Relation between the state parameter and dilatancy

2.3 硬化规律

关于内变量M的硬化规律，采用Li等［25］建议的增量双曲线形式

式中：G为弹性剪切模量；hs为模型参数，hs＝（h1－h2e），h1和h2为模型常数；Mp为峰值应力比；hs和Mp都是各向异性状态变量的函数。模型参数会受各向异性影响的表达式为

式中：kh为模型参数。

峰值应力比受各向异性影响的表达式为

式中：kp是模型参数。从式（12）可以看出峰值应力比也是应力状态、各向异性参量以及应力和组构张量的几何关系3个因素的函数。3个因素中各向异性程度对另外两个因素有实质性影响，各向同性时，A≡0，式（12）只是应力状态的函数，自然退化为各向同性的形式；各向同性时，各向异性影响规律与本文各向异性对砂土临界状态影响相似。对于同一孔隙比、相同各向异性的砂，在同一围压条件下，峰值应力随砂土组构方向与主应力方向的几何关系相关，随几何关系的变化，峰值应力也不断变化，这种变化可以描述考虑主应力轴旋转的硬化内变量规律。可见，考虑各向异性后，砂土的峰值应力比描述更加丰富，同时硬化规律的描述也更全面。

综合式（6）、式（9）和式（12）可以看出：砂土的临界状态、峰值状态和相变状态都是A的函数。通过A将砂土的细观结构和宏观现象联系起来，砂土细观参量、细观组构与应力状态几何关系的变化直接描述了宏观力学特性的变化。如图4所示，材料为各向同性时，3种状态在π平面上为轴对称的各向同性规律；当考虑砂土的各向异性时，各向异性明显影响了3个状态线在π平面上的形状和位置。

3 模型参数确定与验证

3.1 模型参数的确定

各向异性模型参数主要涉及两个参数，即各向异性参量a和模型参数t。这两个参数都反映了砂土的细观组构特性对宏观力学特性影响，需要用宏细观结合的试验确定。Yang等［20］的Toyoura砂细观定量分析表明：落砂法观测的幅值参量平均值为a＝0.214、湿捣法为a＝0.091。参照该试验成果，本文取a＝0.12。参数t为临界状态的试验参数，它的确定需要结合细观实验，用常规三轴压缩、伸长试验临界状态线确定，如式（2）。根据 Yang等［20］试验，两种制样方法得到的值取平均值后得到t＝0.26。其余模型参数与黄茂松等［6］模型参数一致，模型验证参数是根据Verdugo等［2］Toyoura砂试验确定，另外13个模型参数和文献［7］的模型参数表中参数一致。

3.2 模型参数验证

对于a和t两个参数的验证见文献［6］。本文用这组参数对砂土不同沉积面角度三轴试验强度进行验证。图7为Lam等［26］真三轴试验结果，试验采用Toyoura砂，撒砂法制样，e＝0.71，围压98 k Pa，用相同围压对3种试样进行试验，试样分别为H／W＝0.25、H／W＝1.00和H／W＝2.00。该试验没有砂土细观观测结果，参考Yang等［20］细观试验，取a＝0.12。如图7所示，3种试样尺寸在相同围压条件下，强度随沉积面角度的变化规律都相似，即随沉积面角度变化强度单调减小，试验得到90°和0°时强度相差4°左右。图7中的模拟可以看出各向异性参数可以较好模拟随砂样沉积面角度不同强度的变化规律，模拟的最大值和最小值的变化略小于试验值。通过试验模拟再次验证了各向异性参数的有效性。

图7 各向异性砂土的内摩擦角试验模拟Fig.7 Simulation of the peak friction angle for anisotropic sand

3.3 各向异性对强度的影响

图8 三轴压缩不排水试验与模型模拟对比Fig.8 Comparison between undrained triaxial compression test results and model responses

图9 三轴伸长不排水试验与模型模拟对比Fig.9 Comparison between undrained triaxial extension test results and model responses

图10 三轴压缩不排水试验与模型模拟对.Fig.10 Comparison between undrained triaxial compression test results and model responses

图11 三轴伸长不排水试验与模型模拟对比Fig.11 Comparison between undrained triaxial extension test results and model responses

图8至图11是Yang等［17］试验结果和本文模型模拟的对比。通过模拟可以看出，模型可以较好描述试验初期砂土表现的剪缩性、接近状态转化线则表现剪胀性特点。图10和图11是三轴拉伸试验和模型模拟结果。模型对落砂法试样表现的剪缩软化特性、湿捣法试样表现的流动特性、最终出现颈缩特性都能较好描述。图8至图11的模拟情况表明，模型可以用一套模型参数，结合细观定量得到的关系异性参数即可较好模拟各向异性砂土的试验应力－应变－强度特性，初步验证了宏细观结合方法描述砂土特性的合理性。

4 结语

考虑砂土各向异性对其临界状态的影响，建立了砂土的各向异性本构模型，并做了模型试验验证。

将新的各向异性状态变量引入砂土的临界状态方程后，可以综合考虑材料的孔隙比、围压、应力状态和各向异性等状态相关量对材料状态的影响。扩展了砂土材料状态相关的描述范围。

引入状态变量后，模型能够较好描述π平面上各向异性砂的临界状态、相变状态和峰值状态的变化规律。各向异性越大，形状变化也越大，状态面偏离的静水压力轴也越远；各向同性时，3个状态面自然退化到各向同性状态。模型可以自然描述砂土细观状态量的引起宏观剪胀特性的变化，量化了砂土物理状态量变化对剪胀性及硬化规律的影响。

［1］Been K，Jefferies M G.A state parameter for sands［J］.Geotechnique，1985，35（2）：99-112.

［2］Verdugo R，Ishihara K.The steady state of sandy soil［J］.Soils and Foundations，1996，36（2）：81-91.

［3］Cai Z Y，Li X S.Deformation characteristics and critical state of sand ［J］. Chinese Journal of Geotechnical Engineering，2004，26（5）：697-701.

［4］Li X S，Dafalias Y F，Wang Z L.State-dependent dilatancy in critical-state constitutive modeling of sand［J］. Canadian Geotechnical Journal，1999，36：599-611.

［5］Yao Y P，Sun D A，Luo T.A critical state model for sands dependent on stress and density ［J］.International Journal for Numerical and Analytical Methods in Geomechanics，2004，28：323-337.

［6］黄茂松，李学丰，贾苍琴.基于材料状态相关理论的砂土双屈服面模型［J］.岩土工程学报，2010，31（11）：1764-1771.

Huang M S，Li X F，Jia C Q.A double yield surface constitutive model for sand based on state-dependent critical state theory ［J］. Chinese Journal of Geotechnical Engineering，2010，31（11）：1764-1771.（in Chinese）.

［7］黄茂松，李学丰，钱建固.各向异性砂土的应变局部化分析［J］.岩土工程学报，2012，34（10）：1885-1892.

Huang M S，Li X F，Qian J G.On Strain localization of anisotropic sands ［J］. Chinese Journal of Geotechnical Engineering，2012，34（10）：1885-1892.（in Chinese）.

［8］Mooney M A，Finno R J，Viggianii G.A unique critical state for sand［J］.Journal of Geotechnical and Geoenvironmental Engineering，1998，124 （11）：1128-1138.

［9］Yoshimine M，Ishihara K，Vargas W.Effects of principal stress direction and intermediate principal stress on undrained shear behavior of sand ［J］.Soils and Foundations 1998，38（3）：179-188.

［10］Nakata Y， Hyodo M， Murata H，et al.Flow deformation of sands subjected to principal stress rotation［J］.Soils and Foundations，1998，38 （2）：115-128.

［11］ Dafalias Y F. Bounding surface plasticity I：Mathematical foundation and hypoplasticity ［J］.Journal of Engineering Mechanics，ASCE，1986，112（12）：966-978.

［12］Wan R G，Guo P J.A pressure and density dependent dilatancy model for granular materials［J］.Soils and Foundations，1999，39（6）：1-12.

［13］Li X S，Dafalias Y F.Constitutive modelling of inherently anisotropic sand behavior ［J］.Journal of Geotechnical and Geoemvironmental Engineering，ASCE，2002，128（10）：868-880.

［14］Dafalias Y F，Papadimitriou A G，Li X S.Sand plasticity model accounting for inherent fabric anisotropy ［J］.Journal of Engineering Mechanics，2004，130（11）：1319-1333.

［15］Li X S，Dafalias Y F.Anisotropic critical state theory：the role of fabric ［J］.Journal of Engineering Mechanics，ASCE，2012，138（3）：263-275.

［16］Zhao J，Guo N.Unique critical state characteristics in granular media considering fabric anisotropy ［J］.Geotechnique，2013，63（8）：695-704.

［17］Sadrekarimi A，Olson S M.Residual state of sands［J］.Journal of Geotechnical and Geoemvironmental Engineering，ASCE，2014，128（10）：868-880.

［18］李学丰，黄茂松，钱建固.宏细观结合的砂土各向异性破坏准则［J］.岩石力学与工程学报，2010，29（9）：1885-1892.

Li X F，Huang M S，Qian J G.Failure criterion of anisotropic sand with the method of macro-micro incorporation［J］.Chinese Journal of Rock Mechanics and Engineering，2010，29 （9）：1885-1892. （in Chinese）

［19］Poulos S J，Castro G，France J W.Liquefaction evaluation procedure ［J］.Journal of Geotechnical Engineering，ASCE，1985，111（6）：772-792.

［20］Yang Z X，Li X S，Yang J.Quantifying and modelling fabric anisotropy of granular soils［J］.Geotechnique，2008，58（4）：237-248.

［21］Ishihara K.Liquefaction and flow failure during earthquakes［J］.Geotechnique，1993，43（3）：351-415.

［22］Oda M，Koishikawa I，Higuchi T.Experimental study on anisotropic shear strength of sand by plane strain test［J］.Soils and Foundations，1978，18（1）：25-38.

［23］William K J，Warnke E P.Constitutive model for the triaxial behavior of concrete ［C］／／International Association for Bridge and Structure Engineering Proceedings，Bergamo，Italy，1975，19：117-131.

［24］Gajo A，Muir W D.A kinematic hardening constitutive model for sands：the multiaxial formulation ［J］.International Journal for Numerical and Analytical Methods in Geomechanics，1999，23：925-965.

［25］Li X S，Daflias Y F.Dilatancy for cohesionless soils［J］.Geotechnique，2000，50（4）：449-460.

［26］Lam W K，Tatsuoka M.Effect of initial anisotropic fabric andσ2on strength and deformation characteristics of sand ［J］.Soils and Foundations，1988，28（1）：89-106.

（编辑王秀玲）