0.   引言

近年来，随着基础设施建设需求的日益增长，城市地下空间被进一步开发与利用，而既有建筑设施的存在对地下工程施工提出了较高的变形控制要求。尤其是对于上海地区的深基坑施工，市区内轨道交通设施纵横交错、老旧保护建筑星罗棋布，这些重要设施对于地基土体变形极为敏感，这就要求深基坑开挖过程需要控制周边土体变形处在厘米甚至毫米级的微变形水平。严苛的变形控制措施也必然导致周围土体的应变处在较小水平（0.01%~0.1%）内^[1]。

大量室内土工试验表明，当土体应变水平处在非常小应变时（不大于0.001%），其剪切模量远高于土体处在大应变（大于0.1%）状态下的值^[2]，而当土体的应变处在以上两个阈值之间时，剪切模量随着应变的增大呈现非线性减小的趋势。深基坑开挖过程中，围护结构外的土体大多处在此范围内，目前现行上海市《基坑工程技术标准》^[3]建议宜采用能够考虑土体小应变刚度的本构模型进行分析。近些年来，土体小应变刚度特性及相关的本构模型逐渐得到关注^[4-5]，其中Benz^[6]在硬化土HS模型的基础上提出了小应变硬化（HSS）模型，该模型已被Plaxis和Z-Soil等商业有限元软件采用，并且在基坑开挖问题分析中表现出了较好的变形预测能力^[7-8]。

目前上海等软土地区对于土体HSS参数取值开展了较为系统的研究^[9-12]，但HSS模型在工程应用中较为复杂。此外，作为一种半经验半理论的模型，不能从机理上反映土体的临界状态特性，且难以反映土体处在超固结状态的应力应变关系^[13]。如图 1所示，深基坑开挖卸荷过程中，不同位置的土体将经历不同的应力路径，挡墙外侧A点的土单元应力变化方向朝向相应的屈服面外部，即${n_A} \cdot {\text{d}}{\sigma _A} > 0$，随着基坑开挖的进行，应力点向外扩张并始终处在正常固结状态。而坑底B，C两点的应力变化方向朝向屈服面内部，应力点进入历史最大屈服面的内部，即进入超固结状态。大量室内试验与现场测试表明，与正常固结土相比，超固结土呈现出高刚度的弹塑性力学响应^[14]。然而传统本构模型如莫尔-库仑（MC）模型、修正剑桥（MCC）模型等以及HSS模型由于采用单屈服面理论进行计算，认为屈服面内的应力点移动只产生纯弹性变形响应^[15]。这导致此类模型均难以准确反映基坑周围不同位置处土体经历的复杂应力路径和超固结特性。

图  1  基坑开挖周边不同位置土体应力路径示意图

Figure  1.  Stress paths of soils at different locations around excavation

下载: 全尺寸图片 幻灯片

为考虑复杂加卸载路径下土体产生的超固结特性，已有一些基于多屈服面理论框架提出相应本构模型的研究^[16-19]。其中，边界面模型^[16-17]是近年来提出的一种基于临界状态理论的先进土体本构模型，通过对塑性模量进行插值计算，使得当前应力点在边界面内部移动时仍能产生弹塑性响应，克服了传统模型的缺点。为了进一步使土体本构模型能够考虑小应变刚度特征，现有研究采用了不同的方式，例如提出了“滑块”模型和“气泡”模型，此类方法的介绍详见文献[20]。Yao等^[20]通过提出与应变状态相关的土体弹性回弹参数，以评估小应变范围内土体刚度非线性变化，并提出小应变统一硬化（SSUH）模型；张硕等^[21]基于上下负荷面模型并对土体刚度参数计算方法进行了修正，形成了考虑小应变下土体刚度衰减的上海软土统一本构模型。

上述的边界面模型能够反映超固结土的力学响应，但并未考虑土体的小应变刚度特征，需要对其进行必要的修正。此外，在对实际问题开展有限元数值模拟时，需要通过合理的应力更新算法将此类本构模型开发进商业软件的材料库。而目前常使用的显式或完全隐式格式的应力更新算法在计算效率、资源占用以及程序实现的复杂性方面各有不足之处^[22-23]。因此有必要基于更合理的应力更新算法将此类先进本构模型嵌入数值模拟软件中以实现其工程应用。

基于岩土边界面本构理论，同时考虑小应变状态下的土体刚度参数演化特征，本文提出了一种能够考虑土体应变刚度特征的改进边界面模型。基于ABAQUS有限元二次开发平台，使用半隐式应力更新算法实现了模型的有限元应用，并通过试验和理论预测结果对比验证了模型与应力更新算法的合理性。最后在有限元软件中使用本文模型进行基坑开挖数值模拟，通过与现场实测值和其他本构模型计算结果对比，论证了本文模型在分析基坑开挖问题时的合理性。

1.   考虑小应变刚度的边界面模型

边界面模型的特点是通过引入映射法则和塑性模量插值函数，合理反映了应力点在边界面内部移动过程中塑性模量的演化，与MC、MCC等模型相比能够较好地反映土体的超固结特征以及经历复杂应力路径时的变形响应^[16-17]。且边界面模型基于土体临界状态理论提出，与HSS模型相比能够更加严谨地反映土体的力学特征。本节针对Dafalias提出的边界面模型展开改进，引入了小应变状态下土体剪切模量随应变非线性衰减规律，对边界面模型的弹塑性本构关系进行了修正，使模型能够更加合理地反映土体小应变刚度特性，并且对原本的边界面形状函数进行了简化，使其更易在实际工程中应用。本文涉及的表达式均在三维应力空间中给出，模型涉及的变量采用张量形式，应力均采用有效应力进行计算。根据弹塑性力学理论，应力增量导致的总应变增量由以下两部分组成：

式中：$\dot \varepsilon _{ij}^{\text{e}}$，$\dot \varepsilon _{ij}^{\text{p}}$分别为弹性和塑性应变增量。弹性部分采用非线性弹性关系进行计算，其增量型应力应变关系为

式中：${\boldsymbol{C}_{ijkl}}$为四阶弹性刚度张量，

式中：K，G分别为体积模量和剪切模量；${\boldsymbol{\delta} _{ij}}$为Kronecker张量。

1.1   考虑小应变刚度的剪切模量修正

大量室内试验表明，土体的刚度与应变紧密相关，当土体处在小应变范围内，即切应变$\gamma < 0.1\%$时，土体的刚度与其当前应变大小呈现出明显的非线性关系^[24]，为表征小应变范围内土体模量的变化特征，Santos等^[25]给出了小应变状态下土体当前剪切模量随剪应变的变化关系：

式中：G₀为初始剪切模量；${\gamma _{0.7}}$为模型参数，代表剪切模量降至其初始值70%时对应的切应变；$a$为系数，$a = 3/7$。

因此本文为考虑小应变状态下土体刚度随应变的非线性变化关系，采用上式对剪切模量的计算展开修正，即当土体应变处在小应变状态时采用修正的剪切模量，而当土体应变超出小应变范围后（$\gamma \geqslant 0.1\%$），继续采用原先本构模型中的剪切模量进行计算，即

式中：${\gamma _{\text{c}}}$为小应变与常规应变状态的分界点，本文按照Atkinson等^[2]的建议取0.1%；E为弹性模量，按照非线性弹性关系根据泊松比$\nu$、回弹线斜率$\kappa$、初始孔隙比${e_0}$和当前应力状态进行计算：

1.2   边界面方程

Dafalias等^[26]于1986年提出的边界面方程在p-q空间中由两段椭圆弧和一段双曲线组合而成，分段的边界面函数导致模型参数高达14个，且部分参数无明确的物理含义，使得模型标定过于复杂^[27]。因此本文对边界面形状进行合理简化以减少模型参数。简化后模型在三维应力空间中的示意图如图 2所示，图 2中${\bar \sigma _{ij}}$和${\sigma _{ij}}$分别表示边界面上的像应力点和加载面上的当前应力点；${n_{ij}}$为曲面的单位外法向量。边界面$\bar F$和内部的加载面f在应力空间中为同源相似的椭球体，

图  2  边界面模型示意图

Figure  2.  Diagram of bounding surface model

下载: 全尺寸图片 幻灯片

其中

式中：带有上标的变量表示与边界面上像应力点相关的量，而不带有上标的变量表示与加载面上当前应力点相关的量；R为边界面形状参数；M为临界状态线的斜率；${\bar P_C}$和${P_A}$分别为边界面和加载面在p-q空间中与p轴的交点。

1.3   硬化法则

本文模型采用独立的硬化参数${\bar P_C}$以反映边界面的硬化规律，由黏土的等向压缩规律可知，${\bar P_C}$的演化与土体塑性体变的累积相关，即

式中：$\lambda$，$\kappa$分别为e-lnp坐标下土体压缩线和回弹线的斜率；$\dot \varepsilon _{\text{v}}^{\text{p}}$为土体的塑性体应变增量，

式中：$\left\langle {\;} \right\rangle$代表Macaulay括号，$\varLambda$为塑性增量乘子，基于相关联流动法则，${n_{ij}}$为边界面或加载面上某点的单位外法向量，

式中：${\bar K_{\text{p}}}$，${K_{\text{p}}}$分别为边界面塑性模量和实际塑性模量。

1.4   映射法则及塑性模量

本文模型采用广泛使用的径向映射法则^[17]，将映射中心定为坐标原点，如图 2所示，边界面与加载面上的应力关系可为

式中：b为映射系数，且$b \geqslant 1$。随着当前应力状态与像应力状态之间相对关系的变化，b值也随之变化。

将式（15）代入式（13）的第二部分可得

同时，考虑边界面的一致性条件可得映射系数b的增量计算方程为

式中：边界面塑性模量${\bar K_{\text{p}}}$的计算为

而当前加载塑性模量根据边界面与加载面之间的相对距离进行插值计算^[27]：

式中：h，m为模型参数；p_a为标准大气压，取101.325 kPa；$\eta$为应力比；$\delta$为当前应力点与像应力点之间的距离，如图 2所示，${\delta _0}$为其初始值，在计算中可取${\bar P_C}$的初始值。可以看到，当土体处在正常固结初始状态时则$\delta = 0$，即真实的塑性模量与边界面塑性模量相等；而当土体处在超固结状态时，式（19）对真实塑性模量的插值计算能够较好地反映土体的高刚度弹塑性力学响应，这避免了常规单面本构模型的不足。

1.5   弹塑性本构关系

根据弹塑性力学理论，在应力增量的作用下，塑性应变增量大小为

由式（1），（2），同时考虑式（5）提出的考虑小应变的剪切模量计算形式，可得增量型本构关系为

1.6   模型参数的确定

本模型共计10个模型参数，包含临界状态线斜率$M$，初始孔隙比${e_0}$，泊松比$\nu$，土体压缩线斜率$\lambda$，回弹线斜率$\kappa$共5个经典的剑桥模型参数；边界面形参数$R$、形状硬化参数$h$和插值函数控制参数$m$共3个边界面控制参数；以及初始剪切模量${G_0}$和参考剪应变${\gamma _{0.7}}$共2个新添加的小应变参数。与拥有13个参数的HSS模型相比，本文模型的所需参数更少。

以上参数中的$\lambda$，$\kappa$，$M$，${e_0}$可通过一维压缩试验和常规的三轴试验确定，对于实际工程应用也可采用下式根据现场测试土体性质参数进行计算得到^[28]

式中：${C_{\text{c}}}$，${C_{\text{s}}}$分别是土体的压缩指数与回弹指数；$\varphi '$为土体的有效内摩擦角；参数$R$，$h$，$m$可通过室内三轴试验确定，Dafalias等^[26]建议$R$取值2.0或2.72，建议$m$取值为0.02。

而对于小应变参数，${G_0}$与${\gamma _{0.7}}$均可通过共振柱试验获取，为方便工程应用，可采用式（5），利用参考切应变直接进行对${G_0}$进行计算^[21]：

而参考切应变${\gamma _{0.7}}$根据课题组前期针对上海软土的研究经验，其取值范围在$1.6 \times {10^{ - 4}}\sim 4.2 \times {10^{ - 4}}$ ^[9-10]。

2.   本构模型有限元模拟的二次开发

针对上述考虑小应变刚度修正的软土边界面本构模型展开二次开发，基于一种半隐式应力积分算法将上述模型嵌入到ABAQUS的本构模型材料库中。

2.1   ABAQUS二次开发概述

ABAQUS计算过程采用Newton-Raphson迭代解法进行求解，每次求解均会按照图 3所示的逻辑调用UMAT子程序。主程序传入一个适当的应变增量$\Delta \varepsilon$，并由UMAT程序迭代计算对应的应力增量和Jacobian矩阵直至达到收敛。而上述计算过程的核心是用户通过编程完成本构关系计算与各类变量更新，所采用的应力更新算法是否高效合理对此有很大的影响^[29]。

图  3  ABAQUS二次开发UMAT子程序调用原理

Figure  3.  Calling mechanism of UMAT subroutine in ABAQUS

下载: 全尺寸图片 幻灯片

2.2   半隐式应力更新数值算法

应力更新算法的作用是在UMAT调用过程中，由传入的应变增量计算得到应力增量$\Delta \sigma$和其他待更新的内变量。常用的全隐式应力积分算法在求解应力增量$\Delta {\sigma _{n{\text{ + 1}}}}$时，均采用第n+1步的流动方向与塑性模量值，即将它们视作未知量参与迭代求解，这样导致需要求解大量的非线性方程组，占用大量运算资源^[22]。而本文采用半隐式应力更新算法对改进边界面模型进行二次开发，将上述两个变量按照显式格式进行更新，以实现高效的数值计算。下文出现的变量下标n与n+1分别代表上一个分析步末和当前分析步的变量，变量上标(k)代表塑性修正的第k次迭代过程。半隐式应力更新算法的示意图如图 4所示，计算过程如下：

图  4  半隐式应力更新算法示意图

Figure  4.  Diagram of semi-implicit stress updating algorithm

下载: 全尺寸图片 幻灯片

（1）初始化变量

（2）由式（14），（19）计算第n步的流动方向与塑性模量，并计算硬化参数的相应模量：

（3）假定主程序传入的应变增量全部产生弹性变形，按照下式进行弹性试探：

（4）使用式（7）校验边界面屈服方程$\bar F_{n{\text{ + 1}}}^{(k)}$，若$\bar F_{n{\text{ + 1}}}^{(k)} < {\text{ITOL}}$，则迭代计算收敛，跳出循环，ITOL为迭代允许容差，本文取$1.0 \times {10^{ - 7}}$；如果弹性试探后不满足上式则进行塑性修正过程。

（5）基于Newton-Raphson方法对应力更新中的本构关系进行求解，按照下式计算每个迭代步的塑性乘子增量$\Delta\varLambda _{n{\text{ + 1}}}^{(k)}$：

（6）由式（20）计算应力增量，并获得其他内变量的增量，同时完成变量更新：

应力应变及内变量更新后令$k = k + 1$，返回第（4）步进入下一次迭代，直至应力点被回拉至屈服面上，结束计算。

应力更新过程中需为ABAQUS主程序提供Jacobian矩阵以保证半隐式算法的二次收敛性，本文采用一致切线刚度矩阵作为Jacobian矩阵：

3.   本构模型与应力更新算法的验证

为验证改进后本构模型的合理性，将改进后的本构模型预测结果与试验结果进行对比。张硕等^[21]使用安装有局部位移传感器的改造三轴仪针对上海地区软土开展了三轴试验，重点关注了在小应变范围内的土体刚度演化规律。采用此试验结果进行对比验证，土体参数如表 1所示，其中剑桥模型参数采用张硕等针对上海地区黏土的取值，其余参数按照上述模型参数取值方法确定。

表  1  土体参数取值

Table  1.  Soil parameters

土体参数	13 m (4层土)	39 m (6层土)
$M$	1.10	1.15
$\nu$	0.35	0.35
$\lambda$	0.155	0.07
$\kappa$	0.02	0.013
${e_0}$	0.86	0.66
$R$	2.72	2.72
$h$	10	10
$m$	0.02	0.02
${\gamma _{0.7}}$	1.8×10-4	2.8×10-4

下载: 导出CSV 

| 显示表格

图 5为本文本构模型针对上海地区软土在小应变范围内的应力-应变关系的预测曲线及其与三轴试验的对比结果，图中切应变$\gamma = \sqrt {2{e_{ij}}{e_{ij}}/3}$，${e_{ij}}$为偏应变张量。同时也给出了没有经过刚度参数修正的原始边界面模型的计算结果。如图 5所示，采用改进后边界面模型得到两个不同深度土层的应力-应变曲线，曲线斜率随着切应变的增大呈现出迅速减小的趋势，与三轴试验结果吻合，而未经修正的原始模型难以反映出土体的这一特征。

图  5  小应变下模型预测应力-应变关系与试验对比结果

Figure  5.  Comparison of stress-strain curves between predicted and test results under small-strain state

下载: 全尺寸图片 幻灯片

随后，为验证本构模型UMAT程序的正确性，在ABAQUS有限元软件中建立了如图 6（a）所示的模型模拟三轴排水剪切试验，采用单网格的立方体模型反映土单元应力状态，对其中3个不相互平行的表面施加外法线方向的位移约束，其余的3个表面施加围压约束边界条件。剪切过程为单元体试样施加逐渐增大的轴向应变同时控制围压保持不变。将UMAT的计算结果与梁发云等^[10]针对上海软土开展的三轴试验对比，如图 6（a）所示。随后将使用UMAT的有限元计算结果与本文边界面模型的理论解进行对比，并对不同超固结比的土体进行了计算，初始围压均为60 kPa。如图 6（b），6（c）所示，采用本文半隐式应力更新算法的有限元数值计算结果与本构模型理论解十分吻合。

图  6  三轴排水剪切试验预测结果对比

Figure  6.  Comparison between predicted and test results of drained triaxial shear tests

下载: 全尺寸图片 幻灯片

4.   基坑开挖数值模拟的应用

4.1   工程概况

选取上海轨道交通13号线某地铁车站的基坑开挖工程作为计算案例。此地铁车站基坑呈长条形，总长约188 m，宽约20~26 m，开挖深度在17.50~26.22 m。整个基坑开挖分东、西两个标准段，东、西两个竖井及中部的换乘竖井5个部分。其中，位于东部的标准段1开挖宽度20 m，长约60 m，开挖平均深度为17.5 m。将针对此标准段1展开数值模拟研究，此段基坑开挖断面及土层分布如图 7所示，基坑工程所处场地属典型的上海地区土层，详细场地资料见Tan等^[30]的研究。

图  7  基坑开挖施工剖面及土层分布示意图

Figure  7.  Diagram of profile and soil layers of excavation

下载: 全尺寸图片 幻灯片

基坑施工过程采取分层支护开挖，整个开挖过程共分6层，每层开挖标高如图 7所示，基坑支护方式选取地下连续墙结合内支撑的形式，除第一道支撑采用现浇混凝土支撑外，其余4道支撑均采用钢管支撑。

4.2   模型建立

建立了东侧标准段基坑开挖的三维数值模型，考虑模型边界效应对结果的影响，建立有限元计算模型尺寸为28 m×100 m×60 m，由土体，地连墙围护结构和支撑3部分组成。土体和地下连续墙采用8节点实体单元模拟，支撑采用2节点梁单元模拟。土层按照地勘报告进行划分，沿深度方向共分6层土，土体本构模型使用上述改进的边界面模型，通过ABAQUS自定义材料输入模型参数，通过Input文件为每个单元赋予本文边界面模型中状态变量的初值。每层土的参数均按照上述的模型参数确定方法和现场地勘报告实测数据综合确定，如表 2所示。

表  2  各土层模型参数取值表

Table  2.  Model parameters of soil layers

土层\参数	$M$	$\nu$	$\lambda$	$\kappa$	${e_0}$	$R$	$h$	$m$	${\gamma _{0.7}}$/10-4
①填土，②黏土	1.20	0.3	0.044	0.009	0.807	2.72	50	0.02	3.2
②3砂质粉土	1.31	0.3	0.031	0.014	0.728	2.72	50	0.02	3.2
④淤泥质黏土	0.78	0.35	0.081	0.020	1.315	2.72	50	0.02	3.2
⑤1黏土	1.20	0.3	0.041	0.014	0.882	2.72	50	0.02	3.2
⑤2粉质黏土	1.29	0.3	0.017	0.003	0.982	2.72	50	0.02	3.9
⑦粉细砂，⑧粉质黏土	1.38	0.3	0.010	0.002	0.965	2.72	50	0.02	3.9

下载: 导出CSV 

| 显示表格

地下连续墙和内支撑均采用线弹性本构模型进行模拟，地下连续墙和第一道现浇支撑视为混凝土材料，弹性模量定为20 GPa，泊松比0.17；钢支撑的弹性模量定为200 GPa，泊松比0.3。支撑布置形式参照实际施工过程中采用的布置方式^[30]。

为方便网格划分和提高计算效率，使用ABAQUS前处理软件中的布尔运算功能将挡墙与支撑建立成为一个整体部件，分层支护工况使用Model Change功能实现，整个有限元模型共划分21812个单元。

4.3   基坑开挖模拟及结果分析

数值模拟分析步计算过程与实际基坑开挖施工步骤相同，如表 3所示，数值模拟计算过程包括Initial初始分析步、激活地应力场并平衡土体应力、激活地下连续墙及相应接触、依次激活支撑结构以及开挖土体。计算过程中，支护结构的激活以及土体的逐层开挖使用生死单元法实现。此外，已失效单元对应的接触对也需要依次取消激活，否则会造成模型难以收敛。

表  3  基坑开挖分析步设置

Table  3.  Construction steps of excavation

分析步	分析步工况描述
Initial	ABAQUS内置初始分析步
1	激活地应力场，地应力平衡计算
2	施工围护结构并激活相应接触对
3	激活第一道砼支撑，完成第一层土体开挖
4	完成第二层土体开挖，激活第二道钢支撑
5	完成第三层土体开挖，激活第三道钢支撑
6	完成第四层土体开挖，激活第四道钢支撑
7	完成第五层土体开挖，激活第五道钢支撑
8	完成最后一层土体开挖，结束计算

下载: 导出CSV 

| 显示表格

图 8为基坑开挖至不同深度时地下连续墙水平位移沿深度的分布曲线，分别绘制了现场实测值和采用不同土体本构模型的数值模拟结果，其中HSS计算结果由吴瑞拓等^[31]使用Plaxis有限元软件模拟得到。如图 8（a）所示，使用本文模型计算基坑开挖后的挡墙最大水平位移出现在当前开挖面处附近，并且地下连续墙根部的变形很小，与现场实际情况相同。本文模型计算得到的最大挡墙位移为15.3 mm，与HSS模型和MCC模型的结果相比更接近于现场实测值。

图  8  开挖至不同水平时挡墙水平位移沿深度分布曲线

Figure  8.  Lateral displacements of retaining wall versus depth after different excavation steps

下载: 全尺寸图片 幻灯片

图 8（b），（c）分别为开挖第五层完成和开挖至坑底完成时的挡墙水平位移曲线。可以看出，随着开挖的进行，挡墙最大位移出现的位置向下移动，并且总是出现在开挖面的附近，这符合徐中华等^[32]提出的深基坑开挖挡墙变形规律，并且与HSS模型的计算结果相比，使用本文模型能够更加准确地模拟开挖面之上挡墙的水平变形。

在基坑开挖数值模拟分析时，所选用的土体本构模型对计算结果有着巨大的影响，从图 8中可以明显看出，由于MCC模型不能反映土体在小应变状态和超固结状态下表现出的高刚度特征^[33]，导致严重高估基坑挡墙的水平变形。大量研究也表明，使用传统土体模型计算基坑开挖时，由于不能考虑土体小应变刚度特征，围护结构变形计算结果与实际情况严重不符，并且挡墙的变形模式也与实际监测差距较大^[34-35]。而使用本文模型得到的深基坑围护结构变形模式与实际情况吻合，这是由于本文使用的改进边界面模型能够合理考虑土体的硬化特征，并且也对小应变状态下的土体剪切模量进行了修正，因此对于挡墙外侧土体变形的计算更加适用^[36]。此外，由于边界面模型更适合模拟经历复杂应力路径和处在超固结状态时的土体变形特性，与其他本构模型得到的结果相比，本文模型计算得到的挡墙最大水平位移不论是出现位置还是数值均更接近于实际工程监测结果。

5.   结论

考虑小应变状态下土体剪切模量的非线性衰减，对边界面模型的刚度矩阵进行改进，形成了考虑土体应变刚度特征的软土边界面本构模型，并基于ABAQUS有限元软件的用户开发接口，实现了模型的二次开发。通过有限元方法模拟实际地铁车站基坑开挖过程，对基坑围护结构水平位移进行了计算，得到3点主要结论。

（1）在小应变范围内，本文模型的预测结果与室内三轴试验结果吻合，土体的应力-应变曲线的斜率随着当前切应变增大而迅速减小，而未经改进的原始模型难以反映此特点。

（2）在有限元软件中对三轴剪切试验进行了建模计算，其数值模拟结果与理论预测结果吻合的较好，说明本文针对改进模型提出的半隐式应力更新算法能够准确计算土体的变形响应特征。

（3）基坑开挖过程中，挡墙的最大位移出现在开挖面附近，并且数值模拟得到的挡墙变形规律以及最大水平位移与工程实测值吻合。与使用其他本构模型进行的有限元模拟结果相比，本文模型的计算结果更接近工程实测。