DOI: 10.11817/j.issn.1672-7207.2017.07.019

基于广义位势理论的非共轴本构模型验证

温勇^{1, 2}，杨光华^{2, 3, 4}，汤连生¹，钟志辉⁵，姚捷⁶，张玉成^{2, 3}

(1. 中山大学 地球科学与工程学院，广东 广州，510275；

2. 广东省水利水电科学研究院，广东 广州，510610；

3. 广东省岩土工程技术研究中心，广东 广州，510610；

4. 武汉大学 土木建筑工程学院，湖北 武汉，430072；

5. 上海大学 上海市应用数学和力学研究所，上海，200072；

6. 中铁第四勘察设计院集团有限公司，湖北 武汉，430063)

摘要：基于广义位势理论提出的考虑拟弹性塑性变形的弹塑性模型，把塑性应变增量分解为满足弹性分解准则的拟弹性部分和符合传统塑性理论假设的纯塑性部分，这样得到的塑性应变增量方向不仅与应力状态有关，而且与应力增量相关，为解决土的非共轴性问题提供一种有效的方法。通过单剪试验结果及含主应力轴旋转的土体平面问题数值模拟对模型的合理性进行检验。研究结果表明：本文模型计算结果与试验结果吻合较好，且能够描述单剪试验过程中的非共轴现象；此外，与共轴模型数值模拟结果相比，本文模型能够考虑主应力轴旋转产生的土体塑性变形，计算结果更符合实际。

关键词：非共轴；本构模型；广义位势理论；主应力轴旋转；数值模拟

中图分类号：TU443        文献标志码：A         文章编号：1672-7207(2017)07-1817-07

Verification of non-coaxial constitutive model based on generalized potential theory

WEN Yong^{1, 2}, YANG Guanghua^{2, 3, 4}, TANG Liansheng¹, ZHONG Zhihui⁵, YAO Jie⁶, ZHANG Yucheng^{2, 3}

(1. School of Earth Science and Engineering, Sun Yat-Sen University, Guangzhou 510275, China;

2. Guangdong Research Institute of Water Resources and Hydropower, Guangzhou 510610, China;

3. Geotechnical Engineering Technology Center of Guangdong Province, Guangzhou 510610, China;

4. School of Civil and Architectural Engineering, Wuhan University, Wuhan 430072, China;

5. Shanghai Institute of Applied Mathematics and Mechanics, Shanghai University, Shanghai 200072, China;

6. China Railway Siyuan Survey and Design Group Co. Ltd., Wuhan 430063, China)

Abstract: An elastic-plastic model considering quasi-elastic-plastic deformation was proposed based on the generalized potential theory, in which the plastic strain increment was decomposed into quasi-elastic part and pure-plastic part, and the quasi-elastic part obeyed elastic rule and the pure-plastic part obeyed the traditional plasticity theory. The direction of plastic strain increment obtained by this proposed model not only depends on total stress but also the stress increment, which gives an effective method to solve the problem of the non-coaxiality of soil. The results of a simple shear test and numerical simulation of a plane problem of soil involving the rotation of principal stress axes were used to verify the rationality of the proposed model. The results show that the calculation results of the proposed model agree well with the test ones and the non-coaxial phenomenon during the simple shear test can be described reasonably. What is more, the plastic deformation of soil caused by the rotation of principal stress axes can be considered by the proposed model, which makes the calculation results of the proposed model conform to reality better than the coaxial model.

Key words: non-coaxiality; constitutive model; generalized potential theory; rotation of principal stress axes; numerical simulation

当土体主应力方向发生旋转时，塑性主应变增量方向与主应力方向之间存在着非共轴现象，这已被很多室内试验所证实^[1-6]，例如单剪试验、空心扭转试验等。同时，在实际岩土工程问题中(例如偏荷载、交通及波浪荷载作用等)，由于主应力轴发生旋转所引起的非共轴现象也是普遍存在^{[5, 7-8]}。然而，由于传统的本构模型大多建立在双轴和三轴试验的基础上，这些试验受加载系统的限制，在加载过程中荷载的增量始终限制在应变的主方向，使得本构模型的屈服面和塑性势面的数学表达式往往定义为当前的主应力或应力不变量的函数，再结合传统塑性理论的正交流动法则假设，导致得到的塑性应变增量的主方向和主应力的方向始终保持一致(即共轴)。因此，基于传统塑性理论的本构模型无法客观地描述土体的非共轴现象。针对这个问题，国内外学者提出了基于多面概念、边界面概念、亚塑性概念以及中间面概念等非共轴模型来预测主应力和塑性主应变增量的变化趋势。但这些模型的建立没有基于人们熟悉的弹塑性理论，因而在实际应用中(如有限元计算)受到一定的限制^[9-10]。广义位势理论^[11]直接从数学原理出发建立土的本构模型，其特点是数学原理明确，何处引入何种假设清晰，且新理论可表达更为一般的土的本构关系，包含传统塑性本构理论作为其特例，从而为研究土的本构模型提供了更为广阔的理论基础。本文作者基于广义位势理论模型对土的非共轴特性进行分析，提出基于广义位势理论的非共轴本构模型，并通过土体单剪试验结果及含主应力轴旋转的土体平面问题数值模拟对模型的合理性进行检验。

1  传统塑性位势理论的数学实质及其共轴假定

Mises在1928年类比弹性应变增量可用弹性势函数对应力求偏导的表达式，提出了塑性势的概念，其数学表达式为

               (1)

式中：为塑性应变增量；为应力分量；g为塑性势函数；为待定系数，由主空间的本构方程确定。

式(1)即为传统的塑性位势理论，其实质是规定了塑性应变增量各分量之间的比例关系。为了更清楚地了解传统塑性位势理论的数学实质，从一般情况进行探讨^[11]。设主空间上得到的塑性应变增量与应力增量的本构关系为：

           (2)

式中：为塑性主应变增量；为主应力增量。

将主空间上的塑性主应变增量看作1个矢量，若刚好为1个有势场的矢量，则从数学的矢量场理论可知，存在1个势函数g，使得

             (3)

若进一步假设和的3个主方向相同，则由导数表示的坐标变换关系可得

            (4)

将式(3)代入式(4)，有

             (5)

式(5)与式(1)相同，皆为传统的塑性位势理论公式。由以上分析可知：传统塑性位势理论的数学实质是假定了主空间上的本构关系矢量为一有势场矢量，且和的3个主方向一致。因此，传统的塑性位势理论隐含了塑性主应变增量方向和主应力方向共轴的假定。

2  基于广义位势理论的非共轴本构模型

2.1  广义位势理论及其简化弹塑性模型

杨光华等^[11]从数学角度出发建立了更一般的岩土材料本构理论——广义位势理论，它把本构模型的研究分为2个主要部分：1) 用线性无关的势函数梯度矢量来表述主空间的应力应变关系；2) 通过导数表示的坐标变换方法，把主空间上的应力应变关系转换为六维空间的应力应变关系，形成以势函数表述的本构关系。传统的位势理论都可归结为由主空间上的本构关系通过数学坐标变换而得到，不同的只是坐标转换过程中采用的数学假设。

根据广义位势理论可建立简化的弹塑性模型。以应力空间为例，忽略Lode角等影响，塑性应变增量与应力增量的关系一般可表示为

             (6)

式中：A，B，C和D为塑性系数；p和q分别为平均主应力和广义剪应力；和分别为塑性体应变和塑性剪应变。

类似传统弹塑性理论，把总的应变增量分解为弹性部分和塑性部分。其中，弹性部分应变根据弹性理论求得，即

             (7)

式中：为弹性柔度矩阵，可由广义胡克定律确定。

塑性部分应变则根据广义位势理论可得

             (8)

式中：为塑性柔度矩阵，其表达式为^[11]

      (9)

因此，由式(7)~(9)可得到用于有限元等数值分析的一般应力空间本构方程，即

         (10)

由此可见：基于广义位势理论的简化弹塑性模型只需确定4个系数，即可获得计算所需的弹塑性柔度矩阵；广义位势理论模型对A，B，C和D关系不进行规定，且A，B，C和D可通过常规三轴试验、等向压缩试验及等p试验确定^[11]。

进一步分析可知^[11]，传统塑性位势理论实质上规定了A，B，C和D这4个系数满足一定的数学条件，即

              (11)

且当传统塑性位势理论模型为非关联模型时，有

             (12)

当传统塑性位势理论模型为关联模型时，有

             (13)

此时，A，B，C和D可以通过邓肯-张模型中的切线模量E_t、切线泊松比μ_t来表示^[11]，也可以由剑桥模型中的各向等压固结参数λ、回弹参数κ及破坏常数M来表示^[12-13]。

因此，传统塑性位势理论模型可以看作为广义位势理论模型在一定假设条件下的特例。

2.2  基于广义位势理论模型的土的非共轴性分析

传统的塑性位势理论隐含了塑性主应变增量方向和主应力方向共轴的假定。塑性主应变增量和主应力的共轴性与传统的试验观察是一致的，因为在双轴、三轴或真三轴试验的加载过程中应力的主方向始终保持不变，而应变的增量也限制在应力的主方向。然而，含主应力轴旋转的试验研究表明^{[2, 7, 14]}：土的塑性应变增量流动方向不仅仅取决于当前的应力状态，而且取决于当前的应力增量，即塑性应变增量主方向与主应力方向并不始终保持一致；且即使在应力不变量保持不变的条件下，单纯的主应力方向旋转也会引起塑性变形。

针对土的塑性应变增量方向存在非唯一性问题，杨光华等^[15-16]基于广义位势理论提出了考虑拟弹性变形的土体弹塑性模型，其主要思想是把传统不可恢复的塑性应变增量分解为具有弹性应变特性的拟弹性部分和纯塑性部分，如图1所示。其中：拟弹性部分遵循弹性法则，并与应力增量有相同的方向，采用弹性模型表示；纯塑性部分遵循传统塑性理论的假设，方向具有唯一性，可以采用符合塑性理论的假设来建模。

根据以上建模思想，若进一步从p-q平面来分析，式(6)可写成为

         (14)

式中：A_pe和D_pe为拟弹性参数；A_pp，B，C和D_pp为纯塑性系数。

图1  塑性应变增量分解

Fig. 1  Decomposition of plastic strain increment

拟弹性塑性应变部分遵从弹性法则，假设符合广义虎克定律，则有

          (15)

纯塑性应变部分遵从传统塑性理论，即

           (16)

式中：塑性系数A_pp，B，C和D_pp满足传统塑性理论假设的数学条件，即式(11)。

为了方便计算，文献[15]建议：把纯塑性部分直接假设符合关联模型，并将关联模型计算结果与试验结果之间的差值部分归为拟弹性塑性部分，以此确定参数A_pe和D_pe，而此时纯塑性系数A_pp，B，C和D_pp可通过邓肯-张模型参数或剑桥模型参数来确定^[11-12]。经这样处理后的模型不但计算结果与试验结果更符合，而且保持了关联模型弹塑性矩阵对称的优点，利于数值计算。此外，由此得到的塑性应变增量方向不仅与应力状态有关，而且与应力增量相关，这为解决土的非共轴性问题提供了一种有效的方法。下面以土体单剪试验及含主应力轴旋转的土体平面问题数值模拟对模型进行检验。

3  模型的试验及数值验证

3.1  试验验证

单剪试验和空心圆柱扭剪试验是典型的主应力轴旋转试验，但由于空心圆柱扭剪试验中的荷载及位移边界条件比较复杂，而单剪试验可以看成平面问题，其荷载及位移边界条件相对简单。因此，本文选取单剪试验来检验模型的合理性。

在单剪试验加载过程中，主应力和塑性主应变增量方向都会发生偏转，如图2所示^[17]，两者的偏转角计算公式由弹性力学可知，即

           (17)

图2  单剪试验中主应力和塑性主应变增量方向的变化^[17]

Fig. 2  Directions of principal stress and principal plastic strain increment in simple shear test^[17]

根据前面基于广义位势理论模型的土的非共轴性分析，将塑性应变增量分解为具有弹性应变特性的拟弹性部分和纯塑性部分，且和分别采用弹性模型和符合传统塑性位势理论假设的模型来计算。对于平面问题，纯塑性部分的计算可写成为

          (18)

式中：为塑性柔度矩阵，即式(9)，其中的塑性系数A，B，C和D通过邓肯-张模型参数或剑桥模型参数确定。

拟弹性部分的计算可写成为

 (19)

式中：E_pe和m_pe分别为拟弹性弹性模量和泊松比，通过将关联模型计算结果和试验结果之间的差值部分归为拟弹性塑性部分，根据式(15)确定。

采用文献[14]给出的单剪试验结果进行验证。试样分别由直径1.6 mm和3.0 mm、长50 mm的铝棒混合堆积而成，通过一个二维应力装置进行单剪试验，并测定相应的应力，和以及应变，和。该试验材料的初始孔隙比为0.23，起始应力为，，，试验过程中保持不变，剪应力步长。为了进行对比分析，选取类剑桥模型^[12]作为共轴模型，模型参数如表1所示。

表1  模型参数^[14]

Table 1  Parameters of model^[14]

根据前面介绍的方法及式(18)和式(19)，计算出拟弹性部分和纯塑性部分，两者相加得到总的塑性应变增量，将其代入式(17)即可得到塑性主应变增量方向角的变化情况。与共轴模型计算结果进行对比，结果如图3所示。

图3  主应力及塑性主应变增量方向角计算结果与试验结果比较

Fig. 3  Comparison between calculation results and test results of directions of principal stress and plastic principal strain increment

由图3可知：共轴模型计算得到的塑性主应变增量方向角和主应力方向角始终保持相等，即共轴；而本文模型计算得到的塑性主应变增量方向角和主应力方向角在开始时偏差较大，但随着剪应变增加，两者之间的偏差逐渐减小，最终趋于相等，这与试验结果较吻合，说明本文模型能够较好地描述土体单剪试验过程中的非共轴现象。

3.2  含主应力轴旋转的土体平面问题数值模拟

为了进一步检验本文模型的合理性，基于FLAC3D对本文模型进行二次开发，并将其应用于分析一类含主应力轴旋转的土体平面问题，即偏荷载作用问题，数值计算模型及边界条件如图4所示。该模型长10 m，高5 m，土体采用文献[18]中的砂土材料。为了便于说明，考虑以下2种荷载工况：1) 顶部施加均布荷载q₀=500 kPa；2) 顶部左起0~3 m段再施加偏荷载q₁=1 500 kPa，以使土体产生主应力轴旋转。

图4  含主应力轴旋转的土体平面问题计算模型

Fig. 4  Calculation model of plane problem of soil involving rotation of principal stress axes

同样，为了与本文模型进行对比，选取(考虑剪胀性的)类剑桥模型^[13]作为共轴模型，模型参数^[18]如表2所示。本文模型则根据前面介绍的方法求得。利用这2个模型分别计算不同荷载工况下土体的主应力方向和竖向位移。

表2  模型参数^[18]

Table 2  Parameters of model^[18]

图5和图6所示分别为均布荷载作用下共轴模型与本文模型的计算结果。由图5和图6可知：在均布荷载作用下，2个模型计算得到的主应力方向都未发生旋转，最大竖向位移分别为13.9 mm和13.5 mm，基本一致。因此，在均布荷载作用下，2个模型计算得到的土体主应力方向和竖向位移都基本一致，这主要是均布荷载作用下土体为均匀沉降，此时土体不存在主应力轴旋转，故共轴模型与非共轴模型的计算结果基本一致。

图5  均布荷载作用下土体的主应力方向

Fig. 5  Direction of principal stress of soil under uniform load

图6  均布荷载作用下土体的竖向位移

Fig. 6  Vertical displacements of soil under uniform load

图7和图8所示分别为偏荷载作用下共轴模型与本文模型的计算结果。由图7和图8可知：在偏荷载作用下，2个模型计算得到的主应力方向均发生了旋转，最大竖向位移分别为42.5 mm和47.1 mm，本文模型计算结果比共轴模型的计算结果大一些。这主要是偏荷载作用下土体主应力轴发生了旋转，而主应力轴旋转会产生塑性变形^{[2, 7, 14]}，共轴模型由于隐含了共轴假定，不能反映主应力轴发生旋转时产生的塑性变形，而本文模型能够反映土体的非共轴，因此，导致其计算结果比共轴模型的计算结果大一些，这更符合实际，也与文献[19]得到的结论定性一致，从而进一步验证了本文模型的合理性。

图7  偏荷载作用下土体的主应力方向

Fig. 7  Direction of principal stress of soil under eccentric load

图8  偏荷载作用下土体的竖向位移

Fig. 8  Vertical displacements of soil under eccentric load

4  结论

1) 土体的非共轴现象是存在的，而传统的塑性位势理论隐含了塑性主应变增量方向和主应力方向共轴的假定，因而无法客观地描述土体的非共轴现象。

2) 基于广义位势理论提出的考虑拟弹性变形的弹塑性模型，把总的塑性应变分解为满足弹性分解准则的拟弹性部分和符合传统塑性理论假设的纯塑性部分，这样得到的塑性应变增量方向不仅与应力状态有关，而且与应力增量相关，从而为解决土的非共轴性问题提供了一种可能有效的方法。本文模型计算结果与试验结果吻合较好，且能够较好地描述土体单剪试验过程中的非共轴现象。

3) 与共轴模型相比，本文模型能够较好地反映主应力轴发生旋转时产生的塑性变形，计算结果更符合实际，从而进一步验证了其合理性。

参考文献：

[1] ROSCOE K H, BASSETT R H, COLE E R. Principal axes observed during simple shear of a sand[C]// Proceedings of the Geotech Conference. Oslo, Norway, 1967: 231-237.

[2] PRADEL D, ISHIHARA K, GUTIERREZ M. Yielding and flow of sand under principal axes rotation[J]. Soils and Foundations, 1990, 30(1): 87-99.

[3] JOER H A, LANIER J, FAHEY M. Deformation of granular material due to rotation of principal axes[J]. Geotechnique, 1998, 48(5): 605-619.

[4] 沈扬. 考虑主应力方向变化的原状软粘土试验研究[D]. 杭州: 浙江大学建筑工程学院, 2007: 61-128.

SHEN Yang. Experimental study on effect of variation of principal stress orientation on undisturbed soft clay[D]. Hangzhou: Zhejiang University. School of Civil Engineering and Architecture, 2007: 61-128.

[5] 童朝霞. 应力主轴循环旋转条件下砂土的变形规律与本构模型研究[D]. 北京: 清华大学水利水电工程系, 2008: 37-80.

TONG Zhaoxia. Research on deformation behavior and constitutive model for sands under cyclic rotation of principal stress axes[D]. Beijing: Tsinghua University. Department of Hydraulic and Hydropower Engineering, 2008: 37-80.

[6] 严佳佳. 主应力连续旋转下软粘土非共轴变形特性试验和模型研究[D]. 杭州: 浙江大学建筑工程学院, 2014: 65-111.

YAN Jiajia. Experimental research on non-coaxial deformation behavior and modeling of soft clay under principal stress rotation[D]. Hangzhou: Zhejiang University. School of Civil Engineering and Architecture, 2014: 65-111.

[7] ISHIHARA K, TOWHATA I. Sand response to cyclic rotation of principal stress directions as induced by wave loads[J]. Soils and Foundations, 1983, 23(4): 11-26.

[8] 罗强. 非共轴本构模型的数值应用及离心机试验研究[D]. 大连: 大连理工大学土木工程学院, 2012: 74-115.

LUO Qiang. Numerical application and centrifugal modeling study of a non-coaxial constitutive model[D]. Dalian: Dalian University of Technology. School of Civil Engineering, 2012: 74-115.

[9] 黄茂松, 孙海忠, 钱建固. 粗粒土的非共轴性及其离散元数值模拟[J]. 水利学报, 2010, 41(2): 172-181.

HUANG Maosong, SUN Haizhong, QIAN Jiangu. Non-coaxial behavior of coarse granular aggregates simulated by DEM[J]. Journal of Hydraulic Engineering, 2010, 41(2): 172-181.

[10] YANG Y, YU H S. Numerical simulations of simple shear with non-coaxial soil models[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2006, 30(1): 1-19.

[11] 杨光华, 李广信, 介玉新. 土的本构模型的广义位势理论及其应用[M]. 北京: 中国水利水电出版社, 2007: 30-71.

YANG Guanghua, LI Guangxin, JIE Yuxin. The general potential theory for soil constitutive model and its application[M]. Beijing: China Water Power Press, 2007: 30-71.

[12] 杨光华, 温勇, 钟志辉. 基于广义位势理论的类剑桥模型[J]. 岩土力学, 2013, 34(6): 1521-1528.

YANG Guanghua, WEN Yong, ZHONG Zhihui. Similar Cam-clay model based on the generalized potential theory[J]. Rock and Soil Mechanics, 2013, 34(6): 1521-1528.

[13] WEN Yong, YANG Guanghua, ZHONG Zhihui, et al. A similar cam-clay model for sand based on the generalized potential theory[C]// Proceedings of IACGE 2013-Challenges and Recent Advances in Geotechnical and Seismic Research and Practices, ASCE GSP. Chengdu, 2013: 162-169.

[14] MATSUOKA H, SAKAKIBARA K. A constitutive model for sands and clays evaluating principal stress rotation[J]. Soils and Foundations, 1987, 27(14): 73-78.

[15] 杨光华, 姚捷, 温勇. 考虑拟弹性塑性变形的土体弹塑性本构模型[J]. 岩土工程学报, 2013, 35(8): 1496-1503.

YANG Guanghua, YAO Jie, WEN Yong. Elastic-plastic model for soils considering quasi-elastic-plastic deformation[J]. Chinese Journal of Geotechnical Engineering, 2013, 35(8): 1496-1503.

[16] 温勇, 杨光华, 傅旭东, 等. 基于广义位势理论的岩土塑性应变增量方向非唯一性问题的研究[J]. 中南大学学报(自然科学版), 2015, 46(4): 1428-1435.

WEN Yong, YANG Guanghua, FU Xudong, et al. Study on the problem of non-uniqueness of the direction of plastic strain increment of geomaterials based on the generalized potential theory[J]. Journal of Central South University (Science and Technology), 2015, 46(4): 1428-1435.

[17] 黄茂松, 孙海忠, 钱建固. 粗粒土的非共轴性及其离散元数值模拟[J]. 水力学报, 2010, 41(2): 173-181.

HUANG Maosong, SUN Haizhong, QIAN Jiangu. Non-coaxial behavior of coarse granular aggregates simulated by DEM[J]. Journal of Hydraulic Engineering, 2010, 41(2): 173-181.

[18] YAO Yangping, SUN Dean, MATSUOKA H. A unified constitutive model for both clay and sand with hardening parameter independent on stress path[J]. Computers and Geotechnics, 2008, 35(2): 210-222.

[19] 刘元雪, 郑颖人. 含主应力轴旋转的土体平面应变问题弹塑性数值模拟[J]. 计算力学学报, 2001, 18(2): 239-241.

LIU Yuanxue, ZHENG Yingren. Elastoplastic numerical simulation of plane strain problems of soils involving rotation of principal stress axes[J]. Chinese Journal of Computational Mechanics, 2001, 18(2): 239-241.

(编辑  杨幼平)

收稿日期：2016-08-07；修回日期：2016-10-14

基金项目(Foundation item)：国家自然科学基金资助项目(51378131, 41572277)；广东省自然科学基金资助项目(2015A030313118)；广东省水利科技创新项目(2009-25)；广东省岩土工程技术研究中心开放基金资助项目(201404002) (Projects(51378131, 41572277) supported by the National Natural Science Foundation of China; Project(2015A030313118) supported by the Natural Science Foundation of Guangdong Province；Project(2009-25) supported by the Technological Innovation of Guangdong Water Conservancy; Project(201404002) supported by the Geotechnical Engineering Technology Center of Guangdong Province)

通信作者：汤连生，教授，博士生导师，从事岩土力学、地质工程研究；E-mail: eestls@mail. sysu.edu.cn