CN113189301B - 一种弹塑性本构模型的构建方法及系统 - Google Patents

Abstract

本发明公开了一种弹塑性本构模型的构建方法及系统，考虑湿度改变引起的膨胀力效应，包括以下步骤：采用D‑P模型定义的广义应力应变向量、屈服准则和流动法则作为本构模型的基本框架；通过无荷试验，求得并改进膨胀土膨胀系数的表达式；对扩展D‑P模型的应力应变向量、屈服准则、流动法则进行修正；推导出应力应变向量、流动法则的增量形式和改进屈服准则；根据试验建立增湿条件下膨胀土强度参数和变形参数随含水率变化数学表达式；用表达式代替应力应变向量、流动法则增量形式和屈服准则中的对应参数，即可得到浅表层裂隙膨胀土材料湿度应力场弹塑性本构模型。

Description

一种弹塑性本构模型的构建方法及系统

技术领域

本发明涉及土木工程领域，特别是一种弹塑性本构模型的构建方法及系统。

背景技术

膨胀土是自然地质过程中形成的一种地质体，在世界范围内分布十分广泛。我国是世界上膨胀土分布最为广泛的国家之一，先后有20多个省区发现有膨胀土，且多为多为海相或湖相沉积膨胀土。由于地壳的长期运动，这些海相或湖相膨胀土往往会隆起在地球表面。

近地表的浅层膨胀土是一种典型的非均匀三相介质，不仅裂隙发育，而且具有吸水膨胀软化和失水收缩开裂等与正常固结粘土不相同的工程性质。近地表的浅层膨胀土受大气影响显著，土中的含水量波动(干湿循环)和胀缩现象反复发生，导致了裂隙的扩展并向深部土层发展，使该部分土体的完整性降低，强度下降，为强风化带，属于大气显著影响区，厚度一般为2.0-3.0m。该部分土呈疏松、碎散状态，具有强渗透性，由于其裂隙发育，且在雨季时吸水膨胀软化，土体强度往往较低，文献[1]《模拟干湿循环和含低围压条件的膨胀土三轴试验》(肖杰，杨和平，林京松，陈冠一，常锦，倪啸.模拟干湿循环和含低围压条件的膨胀土三轴试验[J].中国公路学报,2019,32(1):21-28.)改进以往对于浅表层土采用高围压三轴试验不当之处，采用低围压条件的三轴试验，得到近地表的浅层裂隙膨胀土干湿循环后，其内摩擦角Φ往往比较小，小于22⁰，如若遇水则内摩擦角Φ会更小。根据文献[2]《中国西部公路建设中膨胀土工程地质问题的初步研究》(杨和平，曲永新，郑健龙，刘龙武.中国西部公路建设中膨胀土工程地质问题的初步研究[J].长沙交通学院学报,2003,19(1):19-24.)，由于近地表的浅层膨胀土裂隙性和显著的胀缩特性，使得膨胀土地区修筑的轻型工程在降雨条件下多有破坏，成为全球性的技术难题。

在膨胀土地区修建大型岩土工程、进行能源及原材料开采以及地质现象分析研究等都需对其在降雨条件下的力学特性开展深入研究。浅表层裂隙膨胀土湿度应力场本构模型的构建是对膨胀土力学特性的深入探索，对分析降雨条件下膨胀土浅表层轻型工程破坏机理具有重要意义。

现有的膨胀土的本构模型有几十种，其构建大多从非饱和土的角度出发进行研究，其本构关系包含多方面内容：两相流体的运动和变形、土骨架的变形和强度等。并大致可分为4类：线性模型、弹塑性模型、弹塑性损伤模型(或结构性模型)及多因素耦合模型。从非饱和膨胀土本构模型发展历程的回顾可知：虽然各国学者已经发展了数十个非饱和膨胀岩土本构模型，但得到工程界普遍认可的极少，并且大多数应用起来比较困难，且均存在一定的局限性，这些膨胀性非饱和土的弹塑性模型的研究并没有直接考虑含水率或饱和度对膨胀土力学特的影响，而是用吸力来表示非饱和的状态，工程应用难度较大，并且大多数本构模型对于低内摩擦角的裂隙膨胀土是否合适值得商榷。

根据文献[3]《基于ABAQUS软件的大直径桩承载力-变形分析》(贺嘉，陈国兴.基于ABAQUS软件的大直径桩承载力-变形分析[J].地下空间与工程学报,2007,3(2):306-310.)可知:(1)M-C模型存在有不能反映中主应力σ₂对土体屈服和破坏的影响以及屈服面有棱角,塑性应变增量计算较难收敛的弊端；扩展D-P模型同时考虑了中主应力σ₂及静水压力对屈服与破坏的影响，且屈服面相对于M-C屈服面光滑没有棱角，有利于塑性应变增量方向的确定和数值计算。(2)虽然大多数文献认为，在模拟岩土材料时，M-C模型较D-P模型更加适合，但是，当摩擦角Φ≦22⁰时，扩展D-P模型模拟土体相对于M-C模型能得到较好的计算结果。

文献[4]《膨胀岩体中的湿度应力场理论》(缪协兴，杨成永，陈至达.膨胀岩体中的湿度应力场理论[J].岩土力学，1993，14(4)：49–55.)定义了膨胀应变的公式。

浅表层裂隙膨胀土在雨水作用下往往内摩擦角Φ≦220，然而，现有的扩展D-P模型并不适应浅表层裂隙膨胀土吸水膨胀的特性，尤其无法考虑膨胀效应，目前尚没有合适的可考虑浅表层裂隙膨胀土材料膨胀效应的弹塑性本构模型，严重制约浅表层轻型建筑物的建设，亟需构建。

发明内容

本发明所要解决的技术问题是，针对现有技术不足，提供一种弹塑性本构模型的构建方法及系统，考虑浅表层裂隙膨胀土材料膨胀效应，提高建模精度。

为解决上述技术问题，本发明所采用的技术方案是：一种弹塑性本构模型的构建方法，本构模型弹性部分的增量表达式为：

其中，S₁、S₂分别表示Δτ′、Δσ′的差分型式；Δτ′、Δσ′分别为修正后剪切应力τ′和修正后平均正应力σ′的增量；G为剪切模量；K为体积模量；Δγ^e为弹性状态下的剪应变增量；Δε^e为弹性状态下的正应变增量；μ为泊松比；E为弹性模量；

由于湿度场中膨胀土的屈服状态可以分为体积屈服、剪切屈服、拉伸屈服三种，所以，本构模型塑性应变的增量表达式也可以分为三种：

本构模型体积塑性应变增量表示为：

本构模型剪切塑性应变增量表示为：

本构模型张拉塑性应变增量表示为：

其中，Δe为塑性应变增量；上标pv、ps、pt分别表示体积塑性变形、剪切塑性变形和张拉塑性变形；g^v为体积塑性势函数；g^s为剪切势函数；g^t为张拉势函数；σ为平均正应力；τ为剪切应力；λ^s为塑性剪切流动因子；λ^v为塑性体积流动因子；λ^t为塑性张拉流动因子；N_ψ为与剪胀相关的材料参数。

剪切势函数g^s一般对应于不相关联的流动法则，其表达式为：g^s＝τ+N_ψσ-N_ψP；P为膨胀应力；

默认值为0，当对应于相关联的流动法则时，

ψ为剪胀角，

为扩展DP准则的材料参数。

体积塑性势函数g^v对应于相关联的流动法则，其表达式为：g^v＝σ。

张拉势函数g^t对应于相关联的流动法则，其表达式为：g^t＝σ＝g^v。

膨胀应力P的计算公式为：

其中，

Δw为含水率的变化量；ρ_d为试件的干密度；ρ_w为水的密度；ε_e为无荷膨胀试验中的竖向应变。

所述扩展DP准则包括：

体积屈服准则：f^v＝σ′＝σ-P；当f^v＞0时，说明膨胀土出现体积屈服；

剪切屈服准则：

当f^s＞0时，此时膨胀土开始进行剪切屈服；其中，σ′为考虑膨胀应力的平均正应力，即修正后平均正应力；τ′为考虑膨胀应力的剪切应力，即修正后剪切应力，其值等于τ；

和

为扩展DP准则的材料参数，

其中，

为内摩擦角，c为粘聚力；

张拉屈服准则：f^t＝σ′-σ^t＝σ-P-σ^t；σ^t为抗拉强度；当f^t＞0时，说明膨胀土发生拉伸屈服。

考虑膨胀应力后，剪切应力τ′和平均正应力σ′的计算公式为：

其中，s_ij为偏应力张量；σ_kk为正应力；σ₁₁、σ₂₂、σ₃₃分别为第一、二、三主应力。

剪切应变增量Δγ和体积应变增量Δε的表达式为：

其中，ε₁₁、ε₂₂、ε₃₃分别为第一、二、三主应变；Δε₁₁、Δε₂₂、Δε₃₃分别为ε₁₁、ε₂₂、ε₃₃的增量；ε₁₂、ε₂₃、ε₁₃分别为剪应变。

本发明还提供了一种弹塑性本构模型的构建系统，其包括计算机设备；所述计算机设备被配置或编程为用于执行上述方法的步骤。

与现有技术相比，本发明所具有的有益效果为：(1)本发明进一步完善了膨胀土本构理论体系；(2)本发明合理反映了膨胀土湿度应力场(膨胀效应)对土体变形和次生应力状态的影响，对于认识浅表层膨胀土及其修建在其上轻型建筑物的破坏机理和影响因素，采取针对性的防治措施，可以减少财产损失和人员伤亡，具有重要意义。

附图说明

图1为边坡数值网格模型；

图2为初始孔压场(kPa)；

图3(a)D-P模型孔隙水压力场分布云图；图3(b)本发明模型的孔隙水压力场分布云图；

图4(a)D-P模型饱和度场云图；图4(b)本发明模型的饱和度场云图；

图5(a)D-P模型水平位移云图；图5(b)本发明模型的水平位移云图；

图6(a)D-P模型塑性剪应变云图；图6(b)本发明模型的塑性剪应变云图。

具体实施方式

本发明实施例提出了一种可考虑浅表层裂隙膨胀土材料膨胀效应的弹塑性本构模型的构建方法，所述方法考虑湿度改变引起的膨胀力效应，包括以下步骤：采用D-P模型定义的广义应力应变向量、屈服准则和流动法则作为本构模型的基本框架；通过无荷试验，并根据弹性力学的有关理论和文献[4]，求得并改进膨胀土膨胀系数的表达式，对文献[4]定义的膨胀应变公式进行处理得到膨胀效应的膨胀应力表达式，利用该表达式对扩展D-P模型的应力应变向量、屈服准则、流动法则进行修正；推导出应力应变向量、流动法则的增量形式和改进屈服准则；根据试验建立增湿条件下膨胀土强度参数和变形参数随含水率变化数学表达式；用表达式代替应力应变向量、流动法则增量形式和屈服准则中的对应参数，即可得到浅表层裂隙膨胀土材料湿度应力场弹塑性本构模型；利用具体算例，验证本发明提出的本构模型的合理性。

下面结合一种弹塑性本构模型的构建方法及系统具体实施方式对本发明做进一步说明。

步骤一：首先根据初始含水率和初始干密度制作环刀试件，并进行无荷膨胀试验，假定试件土体为各向同性的弹性体，根据弹性力学的有关理论及文献[4]，采用柱坐标系对无荷膨胀试验的受力变形特点进行求解，即可得到膨胀土试件的膨胀系数β的表达式，及膨胀应力P的表达式，如下式所示：

式中参数均可通过试验获得，其中：μ为泊松比；Δw为含水率的变化量；ρ_d为试件的干密度；ρ_w为水的密度；ε_e为无荷膨胀试验中的竖向应变；E为弹性模量。

步骤二：剪切应力τ和平均正应力σ进行修正。

引入步骤一中的膨胀应力(公式二)，对扩展D-P模型的广义应力向量进行修正，得到修正后剪切应力τ′和修正后平均正应力σ′，如下式所示：

式中，s_ij为偏应力张量；σ_kk为正应力；σ₁₁、σ₂₂、σ₃₃分别为第一、二、三主应力；σ′、τ′分别为修正后的平均正应力和剪切应力。

同τ′和σ′有关的剪切应变增量Δγ和体积应变增量Δε可表示为：

式中，ε₁₁、ε₂₂、ε₃₃分别为第一、二、三主应变；Δε₁₁、Δε₂₂、Δε₃₃分别为ε₁₁、ε₂₂、ε₃₃的增量；ε₁₂、ε₂₃、ε₁₃分别为剪应变。

应变增量可以分解为：

式中，Δγ^e、Δγ^p分别为Δγ的弹性和塑性分量，Δε^e、Δε^p分别为Δε的弹性和塑性分量。

步骤三：为便于数值计算，本发明的扩展D-P弹塑性本构模型的本构方程和屈服准则、流动法则均采用增量形式。其中本构方程弹性部分的增量格式如下：

式中：Δτ′、Δσ′分别为步骤二中修正后剪切应力τ′和修正后平均正应力σ′的增量；G为剪切模量，

K为体积模量，

步骤四：本构模型弹性部分的增量表达式。

对步骤三中的公式六进行转化，其增量表达式为：

式中，S₁、S₂分别代表Δτ、Δσ的差分型式。

步骤五：对屈服准则进行改进。

膨胀土在膨胀力的作用下会出现体积膨胀变形，因此，在扩展DP模型的基础上加入体积屈服准则，因而将扩展DP模型屈服准则改进为步骤一中膨胀应力引起的体积屈服、剪切应力引起的剪切屈服、张拉应力引起的张拉屈服。

①体积屈服准则

对扩展DP模型进行改进后，需增加膨胀应力引起的体积屈服，其函数形式可表示为：

f^v＝σ′＝σ-P (公式八)

当f^v＞0时，说明膨胀土出现体积屈服。

②剪切屈服准则

对扩展DP模型中剪切屈服准进行改进，其函数形式可表示为：

式中：

和

为扩展DP准则的材料参数，

当f^s＞0时，此时膨胀土开始进行剪切屈服。

③张拉屈服准则

张拉屈服函数未发生改变，仍与扩展DP模型中形式相同，可表示为：

f^t＝σ′-σ^t＝σ-P-σ^t (公式十)

式中：σ^t为抗拉强度。

当f^t＞0时，说明膨胀土发生拉伸屈服。

步骤六：对塑性势函数进行改进。

基于塑性位势流动理论，体积塑性势函数g^v与相关联的流动法则匹配，如下所示：

g^v＝σ (公式十一)

剪切势函数g^s与不相关联的流动法则匹配，可表示为：

g^s＝τ+N_ψσ-N_ψP (公式十二)

式中：P为膨胀应力；

默认值为0，当对应于相关联的流动法则时，

其中，ψ为剪胀角，

为扩展DP准则的材料参数。

张拉势函数g^t与相关联的流动法则匹配，其表达式为：

g^t＝σ＝g^v (公式十三)

此时体积塑性势函数对应的流动法则与张拉塑性势函数对应的流动法则相等，可以认为在D-P模型中引入膨胀应力进行修正后的体积和张拉流动法则是等同的。

步骤七：对流动法则进行改进。

根据步骤六的塑性势函数，对扩展DP模型的流动法则进行改进，相应地将扩展DP模型的流动法则分解为体积流动法则、剪切流动法则及张拉流动法则。

①体积屈服时的流动法则

将体积塑性势函数引入土体发生体积屈服时的流动法则，此时其体积塑性应变增量可表示为

②剪切屈服时的流动法则

将剪切塑性势函数引入土体发生剪切屈服时的流动法则，此时其剪切塑性应变增量可表示为

③张拉屈服时的流动法则

将张拉塑性势函数引入土体发生张拉屈服时的流动法则，此时其张拉塑性应变增量可表示为

可以看出，由步骤四和步骤七分别得到本构模型应变增量的弹性部分和塑性部分，即为本构方程，加上上述的屈服函数和流动法则，即为本发明提出的可考虑浅表层裂隙膨胀土材料膨胀效应的弹塑性本构模型。

步骤八：强度和变形参数随含水率变化的关系式。

本发明提出的可考虑浅表层裂隙膨胀土材料膨胀效应的弹塑性本构模型中E、μ、c、

均与土样的含水率密切相关，其拟合的函数可参考文献[5]《增湿条件下膨胀土路堑边坡稳定性数值分析》(刘义高，周玉峰，郑健龙.增湿条件下膨胀土路堑边坡稳定性数值分析[J].岩土工程学报，2007，29(12)：1870–1875.)，如下：

式中，a、b、f、d、k、t、m、n为与膨胀土体本身性质相关的回归参数。

实施例：

为了验证本发明方案的合理性，通过数值模拟对百色地区某非饱和膨胀土均质边坡的渗流场和应力场进行计算。该边坡高度为15m，宽度为25m，坡顶水平面和坡脚水平面均为8m。由于坡度较缓，边坡采用8节点的六面体单元，单元划分为14620个，节点划分为29808个，建立如图1所示为边坡数值网格模型。

由于膨胀土是典型的非饱和土，其膨胀应力在不同含水率状态下表现的应力大小是不同的，需引入非饱和渗流特性对边坡稳定性的影响，从而考虑因含水率变化而发生不均匀膨胀对边坡稳定性的影响。在验证本发明提出的改进本构模型之前，有必要先引入非饱和渗流场。

均质材料采用Drucker-Prager模型的材料参数，粘聚力μ、内摩擦角

弹性模量E、泊松比μ采用

式中：a、b、f、d、k、m、t、n为拟合参数，w为土样的含水率，本发明中根据试验结果得到的百色膨胀土拟合参数为a＝-941.74、b＝206.20、f＝319.03、d＝20.31、k＝12.35、t＝56.73、m＝35.89、n＝12.73，而本发明的本构模型参数增加了膨胀参数，其取值根据式

式中：μ为泊松比；Δw为含水率的变化量；ρ_d为试件的干密度；ρ_w为水的密度；ε_e为无荷膨胀试验中的竖向应变；E为弹性模量。本发明中根据试验结果得到的百色膨胀土参数为ρ_w＝1.0×10³kg/m³、ρ_d＝1.65×10³kg/m³、ε_e＝0.625，其他参数根据湿度场，按式1进行计算。

渗流计算中，土水特征曲线采用V-G模型拟合，其参数a、n、m分别为25、1.852、0.455，饱和渗透系数采用室内变水头试验测得，其值为4.8e-7m/s。数值模型的参数取值可详见表1所示。

表1参数取值

边界条件东西和南北两边加上对水平方向进行约束，底部施加水平和竖直方向的约束，坡脚水平面为地下水位线，孔压为零，上顶面和下底面的孔压分别为-150kPa和60kPa，生成图2所示为初始孔压场。根据百色地区强降雨天气，将降雨强度设置为1×10^-6m/s，持续时间为70h。

降雨是促使膨胀土边坡产生渐进性失稳的环境因素，通过膨胀土边坡水分迁移的时空分布引起边坡发生宏观变形和塑性破坏,边坡的宏观变形和塑性破坏可以通过降雨结束时边坡水平位移和塑性应变的分布情况，来判断水平位移的逐渐发展过程、滑动面破坏的形状与位置，以及塑性应变区逐渐贯通的情况。由于篇幅所限，以及本发明研究重点之故，以下主要对D-P模型和本发明的模型在降雨70h后图3(a)和图3(b)的孔隙水压力场、图4(a)和图4(b)的饱和度场、图5(a)和图5(b)的水平位移场和图6的塑性应变场的分布结果展开分析。

图3(a)和图3(b)分别对引入膨胀应力前后的Drucker-Prager模型的孔隙水压力场计算结果进行对比。当降雨70h结束后，两种模型计算的孔隙水压力随降雨时间的持续而变大，坡面出现较大范围的暂态水压力，且暂态水压力的范围和深度由坡脚至坡顶方向呈逐渐降低的趋势，地下水位线相较于初始状态有了大幅度提高，两者模型计算的降雨入渗规律相近，但降雨入渗深度和暂态饱和区范围略有不同。比较图3(a)和图3(b)可以发现，本发明的模型，即本发明的模型引入了膨胀应力，膨胀应力的存在促使土体密实，起阻水作用，降雨入渗影响范围变小。

图4(a)和图4(b)所示为改进前后D-P模型在降雨70h后的饱和度场云图，两种模型的数值模拟结果的分布规律与孔隙水压力场相一致，当降雨70h停止时，其饱和度分布相较与初始饱和度有了较大的提升，坡面出现大范围的暂态饱和区，使得膨胀土的含水率增加。且比较图4(a)和图4(b)来看，D-P模型模拟的饱和度场的整体饱和度高于本发明的模型，其暂态饱和区的范围和降雨入渗深度也更大于本发明的模型模拟结果。

图5(a)和图5(b)分别对引入膨胀应力前后的Drucker-Prager模型的水平位移计算结果进行对比。D-P模型计算的水平位移趋势和本发明的模型计算的水平位移趋势相同，均具有渐进性破坏的趋势，且由上图可知，D-P模型最大的水平位移为12.35cm，本发明的模型最大的水平位移为17.49cm，说明本发明的模型能体现膨胀应力对边坡位移的影响。

图6(a)和图6(b)为D-P模型改进前后的塑性剪应变的分布云图，从图中可以明显知道：两种模拟结果的最大塑性剪应变出现在坡脚附近，继而沿着滑动圆弧面不断向坡顶扩展，形成塑性滑动面，随着降雨时间的持续，土体持续软化劣化，应力应变发生重分布，塑性破坏区增加，最终出现贯通的塑性滑动面。但D-P模型和本发明的模型的剪应变模拟数值大小结果不相同，D-P模型的最大剪应变为4.7028×10^-5,本发明的模型最大剪应变为8.7028×10^-5，考虑膨胀应力后的塑性应变显著大于不考虑膨胀应力的塑性应变模拟结果，本发明的模型塑性区贯通坡顶，最大剪应变大于D-P模型，塑性区发展的更饱满，说明了膨胀应力在发生塑性应变到破坏过程中,促使应力场重分布，继而伴随着更大的塑性破坏区产生。因此，考虑膨胀应力的模拟结果更加接近实际的破坏面形状，所以也可以看出在该情况下本发明的本构模型要优于D-P本构模型。

从孔隙水压场、饱和度场、水平位移场和塑性应变场的分布结果中，得到了膨胀应力对改进前后D-P模型的模拟结果，可以看出在该情况下本发明的本构模型要优于D-P本构模型，验证了二次开发程序的合理性。

Claims (9)

1.一种弹塑性本构模型的构建方法，其特征在于，本构模型弹性部分的增量表达式为：

其中，S₁、S₂分别表示△τ′、△σ′的差分型式；△τ′、△σ′分别为修正后剪切应力τ′和修正后平均正应力σ′的增量；G为剪切模量；

K为体积模量；△γ^e为弹性状态下的剪应变增量；△ε^e为弹性状态下的正应变增量；μ为泊松比；E为弹性模量；

本构模型体积塑性应变增量表示为：

本构模型剪切塑性应变增量表示为：

本构模型张拉塑性应变增量表示为：

其中，△e为塑性应变增量；上标pv、ps、pt分别表示体积塑性变形、剪切塑性变形和张拉塑性变形；g^v为体积塑性势函数；g^s为剪切势函数；g^t为张拉势函数；σ为平均正应力；τ为剪切应力；λ^s为塑性剪切流动因子；λ^v为塑性体积流动因子；λ^t为塑性张拉流动因子；N_ψ为与剪胀相关的材料参数；ψ为剪胀角，

为扩展DP准则的材料参数，

为内摩擦角。

2.根据权利要求1所述的弹塑性本构模型的构建方法，其特征在于，

剪切势函数g^s表达式为：g^s＝τ+N_ψσ-N_ψP；P为膨胀应力；

默认值为0，当对应于相关联的流动法则时，

3.根据权利要求1所述的弹塑性本构模型的构建方法，其特征在于，

体积塑性势函数g^v表达式为：g^v＝σ。

4.根据权利要求1所述的弹塑性本构模型的构建方法，其特征在于，

张拉势函数g^t表达式为：g^t＝σ＝g^v。

5.根据权利要求2所述的弹塑性本构模型的构建方法，其特征在于，

膨胀应力P的计算公式为：

其中，

△w为含水率的变化量；ρ_d为试件的干密度；ρ_w为水的密度；ε_e为无荷膨胀试验中的竖向应变。

6.根据权利要求2所述的弹塑性本构模型的构建方法，其特征在于，

所述扩展DP准则包括：

体积屈服准则：f^v＝σ′＝σ-P；当f^v>0时，说明膨胀土出现体积屈服；

剪切屈服准则：

当f^s>0时，此时膨胀土开始进行剪切屈服；其中，σ′为考虑膨胀应力的平均正应力，即修正后平均正应力；τ′为考虑膨胀应力的剪切应力，即修正后剪切应力，其值等于τ；

和

为扩展DP准则的材料参数，

其中，c为粘聚力；

张拉屈服准则：f^t＝σ′-σ^t＝σ-P-σ^t；σ^t为抗拉强度；当f^t>0

时，说明膨胀土发生拉伸屈服。

7.根据权利要求1所述的弹塑性本构模型的构建方法，其特征在于，

考虑膨胀应力后，剪切应力τ′和平均正应力σ′的计算公式为：

其中，s_ij为偏应力张量；σ_kk为正应力；σ₁₁、σ₂₂、σ₃₃分别为第一、二、三主应力。

8.根据权利要求1所述的弹塑性本构模型的构建方法，其特征在于，

剪切应变增量△γ和体积应变增量△ε的表达式为：

其中，ε₁₁、ε₂₂、ε₃₃分别为第一、二、三主应变；△ε₁₁、△ε₂₂、△ε₃₃分别为ε₁₁、ε₂₂、ε₃₃的增量；ε₁₂、ε₂₃、ε₁₃分别为剪应变；

△γ^p为△γ的塑性分量，△ε^p为△ε的塑性分量。

9.一种弹塑性本构模型的构建系统，其特征在于，包括计算机设备；

所述计算机设备被配置或编程为用于执行权利要求1～8之一所述方法的步骤。

Images