17-A.1               引言

CivilFEM 岩土与基础模块旨在成为一种强大的工具，使土木工程师能够 简便地处理其专业领域中的问题。

为了构建一个紧凑的程序，并便于在其中集成各种 工具，CivilFEM 配备了一个具有以下重要特征的岩土数据库：

• 包含并维护 ANSYS 未涵盖的土体和岩石固有 属性。
• 支持处理不同来源和含义的物理量之间的 相关关系，这是岩土工程中非常常见的做法（例如弹性模量 E 与标准贯入试验 SPT 之间的相关关系）。
• 该模块包含原生库，其中包括土体和岩石的特征属性，以及从参考文献中 获取的大量相关关系。
• 用户可以创建自己的 库。

在此基础上，该模块增加了以下工具：

• 边坡稳定性， 包括地下水位、加筋土等的定义。
• Hoek-Brown模型，用于模拟岩石中基础的行为；Mohr-Coulomb模型或Cam-Clay模型用于土体研究。
• 确定基础的地基基床系数。其确定基于 Boussinesq 方程和有限元法计算模型。
• 生成并计算 支护结构。
• 渗流网。
• 隧道模型。
• 土体初始应力。
• 计算并施加土压力对结构产生的 作用。

 

17-A.2      分层地层的定义

分层地层功能用于处理在模型中未直接采用有限元 进行建模的地层。

可以定义由用户编号标识的不同模型。程序允许每个地层最多定义 20 层（参见命令 ~TERDEF).

其适用的应用包括：

• 基床系数计算。
• 采用 2D 模型计算支护结构。
• 将土压力施加到模型上。
• 桩承台。

对于所有这些应用，都必须指定要采用的地层编号，并检查其是否已定义 相关计算所需的全部属性。

 

 

17-A.3               基床系数模块

利用此功能，用户可以估算基床系数的值，从而近似模拟土体的弹性模型 （土体由其参数 E 和 n) 通过 Winkler 模型（由其地基反力系数定义）。

弹性模型和 Winkler 模型在概念上是不同的，因此不可能进行精确对应。实际上，根据弹性模型得到的地基反力系数并非仅取决于材料的性质，因为它还与基础的形状以及测量基础变形的点位有关。甚至在基础范围以外，该值也不为零。

获取地基反力系数的适当方法，是在基础平面所占区域施加单位荷载并求得沉降。各点的沉降等于地基反力系数的倒数。

为便于用户使用，所执行的命令除可获得点值外，还可获得平均值、最大值和最小值。

17-A.3.1    布辛奈斯克

在所有情况下，我们均以 Boussinesq 给出的由数值为 Q 的集中荷载在地表产生的沉降解为基础，该解在柱坐标中的表达式为：

其中

该规律在原点处存在奇异性，即荷载正下方。

17-A.3.2    任意基础的沉降

考虑承受均布荷载 q 的基础，根据 Boussinesq 对任意荷载给出的表达式，原点处的沉降可通过积分求得

要计算其他任意点的沉降，需要依次进行坐标轴平移，即：

 

17-A.3.3    特殊情况

这些情况未作为单独情况纳入代码或程序中。仅实现了一般情况。

但是，以下特殊情况已通过与 CivilFEM 的结果进行手工对照加以验证，并已用于校核示例。

 

17-A.3.3.1                圆形表面上的均布荷载

考虑半径为 a:

• 板中心下方的沉降

待计算的积分为：

因此，

• 在边缘处：

此时积分为，

其中边缘处沉降与中心处沉降之比为，

• 在任意点处：

其中 K 和 E 分别为第一类和第二类椭圆积分，

17-A.3.3.2                矩形荷载

对于尺寸为 a 的矩形´b，承受均布荷载 q，Steinbrenner（1936）推导出了角点下方沉降的如下表达式：

其中 n=a/b。

因此，要获得任意点处的沉降，应将矩形划分为四个部分，使该点成为它们的公共点。

17-A.3.4    由不同材料组成的地基

当地基由 n 层厚度为 h 的水平土层组成时_i，弹性模量 E_i，以及 泊松比 n_i，计算采用等效值 E_e. E_e 采用 Thenn de Barros（1966）给出的表达式的一般形式计算：

n_e 是以下各项的平均值：n_i

 

 

17-A.3.5    积分

由于该问题针对一般情况建立，因此待求积分如下：

其中 a y b 为基础上待求地基基床系数点的坐标，且 W 为组成基础的所有单元所形成的面积。

该积分表示为各单元积分之和：

 

这些积分采用高斯点进行数值求解。但是，对于包含点 P(a,b) 的单元存在一个问题。在该点处，函数

 

存在奇异性。尽管积分应为有限值，但进行数值计算时会产生问题，因此对于这些单元，采用极坐标求解积分。为此，将单元分解为具有公共顶点的三角形，该公共顶点为点 P(a,b).

图 17-A.3‑1 单元的三角形分解

 

进行极坐标变量变换后，积分为：

 

对于每个三角形，

 

 

17-A.3.6    基床系数与 EFS 之间的关系

EFS（Elastic Foundation Stiffness，弹性地基刚度）是地基承受竖向荷载时的刚度。因此，EFS 直接就是基床系数。

17-A.4               地下连续墙计算 1 ½ D

17-A.4.1    简介

这是采用随开挖标高变化的二维非线性演化有限元模型进行地下连续墙计算。地下连续墙采用梁单元建模，其与地基的相互作用通过梁和接触进行模拟。

地下连续墙的材料可以是钢材或混凝土，也可以将地下连续墙定义为非线性构件，并采用非线性混凝土模块进行计算。

地下连续墙必须竖直，并可考虑未开挖土体部分存在附加荷载；该部分可以具有倾斜面，如图所示。

图 17-A.4‑1 问题示意图

17-A.4.2    地基

可以考虑分层地基。对具有不同土层的地基进行建模时，必须使用命令 ~TERDEF，通过该命令可以定义计算地下连续墙土压力所需的全部材料属性。

17-A.4.3    建模

挡土墙采用二维梁单元进行建模，并在其上施加所有独立于土压力和地下水压力的边界条件及作用，即在挡土墙上设置所有能够阻止其移动或向其施加力的构件（嵌固、锚杆等）。

与土体的相互作用由两对连接到 gap 的弹簧（LINK1 类型）代替，使其仅受压工作。因此，每对弹簧分别用于模拟各方向上的主动土压力和被动土压力。

 

图 17-A.4‑2 挡土墙建模

每个弹簧材料本构关系的输入采用非线性弹性行为模型。

图 17-A.4‑3 材料本构关系

连接各水平线的直线斜率是根据水平地基反力系数推导出的刚度，该系数可作为土体参数输入（参数 HBM 命令 ~TERDEF），也可由 CivilFEM 根据 Boussinesq 理论（第17-A.3.1章）或 Chadeysson 方法（桩承台相关章节中所述）计算得到。

17-A.4.4    计算过程

该模型通过逐步计算求解，其中每个“Load Step”代表开挖的一个步骤（参见命令 ~WALLSTP，用于定义不同的开挖阶段）。每一步的开挖深度必须大于或等于弹簧之间的距离，即挡土墙梁单元的长度。随着开挖推进，通过单元生死功能移除已开挖部分的弹簧，并更新其余弹簧的属性。

每一步对模型进行以下更新：

• 移除 开挖部分中的弹簧
• 改变 开挖面以下弹簧的刚度
• 通过在该步骤开挖节点上施加力，将移除开挖部分弹簧后失衡的静止土压力施加到模型上
• 通过在每个开挖节点上施加挡土墙两侧静止土压力之差，更新位于开挖节点以下弹簧的静止土压力
• 当开挖超过地下水位时，施加静水压力之差

17-A.4.5    挡土墙土压力计算

某点处土体对挡土墙的土压力，计算为该点以上土柱的贡献、黏聚力项以及土体表面荷载项之和。因此可得到以下适用于主动土压力和被动土压力情况的表达式：

图 17-A.4‑4 土压力

其中：

Kh	由土体自重引起的 水平土压力系数。
Khc	由黏聚力引起的 水平土压力。
Khq	由均布荷载引起的 水平土压力系数。

在 CivilFEM 中定义土体时（命令 ~TERDEF），必须向程序指定用于计算这些系数的理论，可选择 Rankine 理论或 Coulomb 理论。如果未作指定，程序将采用 Rankine 理论。

 K_hc 不是由黏聚力引起的水平土压力系数，而是直接表示水平土压力。

 

 是 由黏聚力引起的水平土压力系数

17-A.4.5.1                根据 Rankine 理论的主动土压力

采用 Rankine 理论时，主动土压力系数按下式计算：

 

17-A.4.5.2                根据 Rankine 理论的被动土压力

应用 Rankine 理论，被动土压力系数 按下式计算：

 

17-A.4.5.3                根据 Coulomb 理论的主动土压力

应用 Coulomb 理论，主动土压力系数 按下式计算：

这些系数乘以 cosd，以获得土压力的水平分量。

 

17-A.4.5.4                根据 Coulomb 理论的被动土压力

应用 Rankine 理论，被动土压力系数 按下式计算：

其中：

b	地表面与水平面所成的 角度。
f	材料的 摩擦角。
d	土体-挡墙摩擦角。

 

这些系数乘以 cosd，以获得土压力的水平 分量。

17-A.4.6    支座和锚固

可在挡墙的任意位置定义三种类型的支座。这些支座可在 开挖过程中的任意时刻创建或删除。

内置锚固如下：

 

Tipo	说明
0	固定锚固。
1	铰接锚固。
2	不恢复位移的固定 锚固。

 

固定锚固是一种边界条件，要求定义该支座的节点 的运动（位移和转角）为零。

铰接锚固表示创建一根具有所需特性的梁（单元类型、横截面 等），其一端锚固于挡墙，另一端嵌固于土体中。

不恢复位移的固定锚固与前述固定锚固的工作方式相同， 但其位移不强制设为零，而是取创建该锚固时挡墙的位移。

此外，还可以在两道挡墙之间创建锚固， 并将其模拟为梁。

 

 

17-A.5               边坡稳定性

17-A.5.1    简介

岩土工程模块中此实用工具的用途是为 CivilFEM 用户提供一种 用于分析边坡抗滑安全性的工具。

该分析通过两种概念上不同的 体系进行：

-          经典破坏方法。

-          根据有限元应力-变形模型提供的结果。

17-A.5.2    经典破坏方法

在所有情况下，均假定可能的滑动面满足Mohr-Coulomb关系。

t	剪应力。
sn	法向应力（总应力）。
u	孔隙水压力。
j’	有效应力摩擦角。
c’	有效应力黏聚力。

（下文中，每当提及 c 或 j时，均指有效值 c’ 和 j’).

 

17-A.5.2.1                Fellenius 法

瑞典法或独立条分法。

该方法严格来说仅适用于圆弧，其基础是计算相对于圆心的力矩平衡。该方法不是迭代过程，因此可用于计算初始值。

Fellenius方法忽略条间力的影响，并得到安全系数的以下表达式

在不考虑地震作用时，其中 a 为条分底面的长度。

 

图 17-A.5‑1 条分上的作用

 

17-A.5.2.2                Bishop 法

该方法同样针对圆弧面建立，忽略条分间作用力产生的力矩，从而得到以下方程

 控制迭代过程的参数（容差、最大迭代次数）通过选项参数（命令 ~SLPOPT).

 

17-A.5.2.3                Janbu 法

与前两种方法不同，该方法基于力平衡。

该方法适用于任何类型的破坏面（通常采用多边形线），应用于圆弧时所得数值与 Bishop 法非常接近。

现在需要求解的递推方程为：

其中

p	P/a为每单位条宽的重量。
t	(VI –VD) /a
Q	可能作用于第一条分上的水平作用力。

 

17-A.5.2.4                修正 Janbu 法

为修正忽略条分间竖向相互作用力差异影响这一问题，引入系数 f_o，取决于比值d/L（见附图）以及构成边坡的材料类型；简化Janbu方法得到的安全系数值乘以该系数。

图 17-A.5‑2 参数 L 和 d 的定义

系数 f 的数值_o 可通过以下形式的多项式方程求得

系数 A、B 和 C 取决于滑动面穿过的材料类型，如附表所示。

 

 

对于材料不均匀的边坡，假定 j’>0 或 c’>0，当相应数值在边坡穿过的滑动面长度至少 50% 的范围内不为零时。

 

17-A.5.2.5                边坡上的作用

• 地震作用：

通过施加在条块重心处的加速度引入，其水平和竖向分量为： seismh*g  和 seismv*g 

总重量为： W*(1+ seismv)

滑条中心处的水平力为： W*( seismh)

符号约定：

竖向加速度：向下为正 （Y轴）。

水平加速度：沿滑动方向为正。

• 直接施加在边坡上的 荷载：

这些荷载取自 ANSYS，包括集中力和压力。集中力为施加在该节点上的 所有力的合力，并作为滑移矢量传递给作用于滑动面上的滑条。压力以 类似方式计算并传递。

仅考虑进入滑动体体积内的 荷载。

如果要输入锚杆，应在其两端分别施加两个力，使得当两者均位于圆内时 相互抵消；当圆与锚杆相交时，则考虑其中一个力。

• 边坡被水淹没时作用于其 表面的水荷载：

用户应将水荷载以压力形式输入。该过程不会自动执行，以便用户 拥有更大的操作自由度。

 

17-A.5.2.6                适用性

Bishop法以及一般基于圆形破坏面的各种方法，特别适用于均质黏性材料。 对于非黏性材料或由不同材料组成的边坡，通常采用折线更为合适。

计算假设中的简化认为其夹角应为中等大小的 a。因此，对于以较大角度与地面相交的圆，计算结果会出现偏差。

 

17-A.5.2.7                计算假设

 

图 17-A.5‑3 边坡滑动计算参数

图中参数如下：

W	滑条总重量（P*(1+ seismv)）。
kW	地震的水平分量（P* seismh）。
N	滑条底面的法向力。
Sm	滑条底面产生的剪力。
D	施加在该滑条上的外力。
a	滑条底面与水平面之间的 夹角。
R	圆的半径。
L	滑条底面的长度。
x	圆心到滑条的水平距离 （力矩臂）。
e	滑条中心到圆心的竖向距离 （力矩臂）。

 

力矩平衡：

注意：当出现类似  上部表示向左滑动，下部表示向右滑动。

条块底部的剪力为：

                                        

相对于圆心的力矩平衡：

                                  

其中，MD 是 D 相对于圆心的力矩（顺时针方向）。

代入 S 的表达式_m:

力平衡：

应用水平方向的力平衡：

 

N 的计算：

条块底部竖直方向的力平衡为：

 

其中 X 为各条块的切向力。代入 S 的表达式_m:

采用 Fellenuis 简化方法（忽略条块间的剪力，并将各力投影到与底面垂直的方向）：

 

Bishop 和 Janbu 计算所采用的假设仅忽略条块间的剪力，因此 N 的表达式为：

 

在上述表达式的分母中出现了项

该项可能取接近于零的值，此时计算结果将不再正确。这发生在以下情况：

 

-          a 为负，且  较大

-          a 较大（cos a < 0)  y  较小

 

这两种情况均对应于“缺乏合理性”的滑动面曲线。

通常，在被动区，滑动面的倾角应满足限制条件

 

在主动区也应满足类似条件。

 

为避免此问题，如果在某种情况下达到状态 m，CivilFEM 将舍弃该滑动面_a < 10^-5.

 

 

17-A.5.3    根据应力-变形模型得到的结果进行计算。

 

如果已经采用有限元法进行了计算，则可通过以下方式定义安全系数

其中：

sn	滑动面上的法向应力。
t	滑动面上的切向应力。
a	条块宽度。

 

需要注意的是，这里所采用的安全系数概念不同于弹性方法中的概念，因为后者要求满足两个条件。

• 所有条块的安全系数相同。
• 所有相关材料的黏聚力分量和摩擦分量所对应的安全系数相同。

为获得这些应力，程序使用ANSYS的path函数，沿用户通过破坏圆定义的滑动面几何形状，透明地为用户创建曲线折线。定义该path的点以与划分条块的点相同的方式创建，因此可通过条块数参数进行修改。path的法向量和切向量通过计算每个点处的应力获得（点数 = 条块数 +1），并对与该点相邻条块中的c、u和fi值进行加权。

显然，模型必须已经完成计算并具有应力。此外，从该计算的公式可以推知，在命令中可以输入的加速度~SLPSOL 不会被考虑。如果用户希望考虑这种可能性，则必须将这些加速度输入 ANSYS 模型中。

 

17-A.5.4    孔隙水压力

17-A.5.4.1                运行工况计算

获得某点孔隙水压力的方法是定义等压线网格。位于两条等压线之间某点的压力可通过线性插值得到。

 

17-A.5.4.2                施工结束时的计算

某点的压力计算为该点以上土体的重量乘以材料特征系数“Ru”。如果参数“Rusi”为零（逻辑值 0 或 1），则不考虑该材料中的孔隙水压力。

出现这种情况是因为在低渗透性材料中，施工期间引入的水可能尚未完全排出。

 

17-A.5.4.3                压力网格计算

使用 CivilFEM 求解渗流模型后，命令~WTSLP 生成 文件“jobname.press”，其中包含模型的孔隙水压力，并可在另一台计算机上读取。如果其他用户要使用该文件中的数据，jobname 必须与该文件的名称一致。此外，用户有责任确保两个模型的几何形状相互一致。

文件“jobname.press”必须存在于工作文件夹中。在~SLPSOL 期间创建文件“jobname.slpspd”，用于提供和显示结果。

 

17-A.5.5    破坏线定义

17-A.5.5.1                圆族

圆定义为一组圆心位置和一组切线。

圆心位置由网格确定，网格由三个点以及每条边上的分段数给出：

图 17-A.5‑4 圆心网格

圆的切线由四个点和切线数量定义，这四个点构成一个四边形。

图 17-A.5‑5 切线多边形

将从网格的每个圆心计算至多边形每条切线的圆。通过多边形定义的切线视为无限延伸，即多边形不限制切线的长度。

定义圆的另一种方法是不使用切线，而是指定一个所有圆都必须经过的点（定点圆）。

 

 

17-A.5.5.2                折线族

根据一系列竖直线段（n_sv) 按以下参数定义。

·         x 位置_i

·         y 标高_ij

图 17-A.5‑6 折线生成

折线将以所有可能的方式连接所定义的点生成。

仅考虑凸折线，所有折线 都必须在地形外开始和结束。折线必须经过 所定义的所有 X 坐标。

 

17-A.5.6    加筋土

可以创建不同的加筋土组。每个 组必须具有有效的几何形状，该几何形状由与边坡表面 交点的 X 坐标，以及这些点处加筋条或加筋面的 长度定义。中间加筋条的长度将通过插值确定。

还必须指定加筋条（或加筋面）的数量、每个加筋 面的表观摩擦系数，以及每个加筋体的最大允许力。

图 17-A.5‑7 加筋土

表观摩擦系数必须考虑每层加筋材料的两个摩擦面。

可以指定加筋体是否与边坡表面锚固。在第一种情况下（锚固），位于 地形非滑动部分的加筋体段会试图 滑动，因此会产生阻止运动的摩擦力。如果加筋体未固定，则与运动相反的 摩擦力将由加筋体中锚固程度较低的部分产生（其 摩擦力较小）。

图 17-A.5‑8 锚固或未锚固的加筋体

如果某根加筋条或某个加筋面的摩擦力大于允许值，则该加筋条或加筋面将被视为 断裂，并在计算中不予考虑（不会产生任何摩擦力），仅由该组其余加筋体共同作用。

 

17-A.5.7    错误代码

 

并非所有定义的滑动面都可以 计算。CivilFEM 将提供有效滑动面的安全系数。如果滑动面无法使用，将返回一个代码（负数），以代替安全系数，说明忽略该滑动面的原因。

 

 

17-A.6                Mohr-Coulomb塑性模型

17-A.6.1             简介

CivilFEM中实现了Mohr-Coulomb塑性模型的两种公式，每种公式适用于不同的单元族。在这两种公式中，屈服面均是对标准Mohr-Coulomb面的修改，用于修正其顶点和棱边处的奇异性，因为这些奇异性会导致收敛问题。第一种公式（以下简称：公式I）基于参考文献[1]和[2]，第二种公式（以下简称：公式II）完全基于参考文献[2]。两种公式相似，区别仅在于偏应力平面内逼近标准曲面的方式不同。

17-A.6.2             屈服函数和塑性势

在两种表述中，屈服函数均可统一表示为（[2]）

 

其中

K 是 一个函数，对于每种表述各不相同，hyp 是一个参数，c 是黏聚力，f 是摩擦角。

两种公式的曲面渐近逼近标准Mohr-Coulomb面。参数 hyp 是修改后屈服面的顶点与标准屈服面顶点之间的距离。对于，修改后的屈服面在顶点处是光滑的；对于，其通过标准屈服面的顶点，从而恢复奇异性。

取值为 可获得极佳的近似效果，该值也是 CivilFEM 默认采用的值。

 

屈服面在偏应力平面中的形状由函数K确定。在公式 I（[1]）中

 

其中e是偏应力平面中的偏心率，而 是洛德角，其定义为

 

其中J₂ y J₃  分别为应力张量偏应力分量的第二和第三不变量。当时，屈服面为凸面。对于  偏应力截面为圆形；对于 偏应力截面退化为三角形，顶点处产生奇异性（CivilFEM 仅允许取值）。 取值为 可获得最佳近似；此外，当时，屈服面通过两个主应力相等时对应的标准屈服面顶点

公式 I。不同偏心率下的偏应力截面

在公式 II（[2]）中

 

此处洛德角定义为

θ_T 是角度 ，用于对偏应力平面中的顶点进行圆角过渡（CivilFEM 采用 25º 的固定值，用户无法修改）。  常数A y B为

公式 II。偏应力截面（ e )

在两种公式中，塑性势函数G与塑性函数具有相同的表达式；对于非关联塑性流动，摩擦角f 被胀角y.

17-A.6.3             输入数据

Mohr-Coulomb模型仅适用于土体和岩石（KPLA = 2，在~CFMP). 弹性性质为各向同性。两种公式的材料性质和参数均相同，唯一区别是偏心率参数，该参数仅适用于公式 I。定义该模型的材料性质和参数见下表

材料性质和参数

流动法则的类型由参数IFLOW指定：若为关联流动法则，则IFLOW = 0；若为非关联流动法则，则IFLOW = 1。关联流动中无需指定剪胀角

不允许将黏聚力设为零；对于无黏性材料，必须为其输入一个较小的值。

17-A.6.4              Mohr-Coulomb模型（公式I）与ANSYS模型Drucker-Prager和扩展Drucker-Prager之间的关系

Mohr-Coulomb模型（MC）的公式I具有很强的通用性，并将ANSYS塑性模型Drucker-Prager（DP）和扩展Drucker-Prager（EDP）（线性形式和双曲线形式）作为特殊情况包括在内。下表给出了DP和EDP模型参数与MC模型参数之间的关系，以使这些模型彼此等效。对于非关联流动情况，这些关系类似，但将摩擦f  由剪胀角 y 在塑性流动势中。

DP 与 MC 模型参数之间的关系 (e = 1, hyp = 0)

EDP 与 MC 模型参数之间的关系 (e = 1)

EDP 模型	参数	表达式 F y G	与 MC 参数的等效关系： c, f y hyp
Lineal
双曲型

17-A.6.5             错误消息

求解过程中可能出现以下错误消息：

·      MC 塑性：塑性算法在单元…的积分点处不收敛…

描述：	无法求解单元某个积分点处塑性流动的控制方程。
程序操作：	分析停止。

·      MC 塑性：在单元…的积分点处发现奇异矩阵...

说明：	在求解塑性流动控制方程时，建立了奇异或近似奇异的方程组。这可能是求解控制方程过程中瞬态效应所致。
程序操作：	无。

17-A.6.6             假设和限制

·         仅适用于平面应变、轴对称和 3-D 问题。

·         各公式适用的单元：

公式 I：	PLANE42, SOLID45, SOLID65, PLANE82, SOLID92 y SOLID95
公式 II：	PLANE182, PLANE183, SOLID185, SOLID186 y SOLID187

·         当塑性流动为非关联流动（IFLOW = 1）时，弹塑性模量是不对称的，因此单元刚度矩阵也不对称；所以采用Newton-Raphson方法时，必须启用非对称方程组求解（NROPT, UNSYM).

·         允许大变形和大位移。

·         初始应力仅适用于公式 II（命令 INISTATE y ~TIS）。不考虑初始弹性变形和塑性变形。

17-A.6.7             参考文献

[1]   P. Menétrey 和 K. J.\r Willam，《混凝土三轴破坏准则及其推广》，ACI\r Structural Journal, 92:311–318, 1995年5/6月。

[3]   A.J. Abbo And S. W. Sloan，"A Smooth Hyperbolic Approximation to the Mohr-Coulomb Yield  Criterion"，Computers & Structures, Vol. 54, No. 3, pp. 427-441 (1995)。

17-A.7                Cam-clay塑性模型

17-A.7.1             模型说明

该实现对应于称为 Cam-clay o 修正Cam-clay 由 Roscoe 和 Burland（1968）提出的模型。该模型的基本要素为：

-     塑性 条件

-     相关 流动法则

-     硬化 定律

-     次弹性 本构定律

塑性条件为

 

其中  y  S 为应力张量的偏量分量 σ),  为先期固结压力，M 为临界状态线（CSL）的斜率。

硬化定律为

其中v 为比容（是 一定体积材料中所含固体颗粒所占的体积V), N 是沿各向同性压缩线（ICL）在单位压力下的比容，λ y κ 是ln平面内ICL线和加载-卸载线（URL）的斜率p-q。 然而，为便于求解，采用等效的增量形式更为合适

其中为塑性体积应变。

本构关系为

其中 为应力张量，e y 分别为总应变张量和塑性应变张量。弹性张量C取决于应力状态；假定泊松比为常数，而其余材料常数：体积模量K、剪切模量G和杨氏模量E随应变变化（但在每个时间步内保持不变）：

17-A.7.2             初始状态的指定

在Cam-clay模型中，必须指定初始应力状态  非平凡（INISTATE）以及初始固结状态，后者可通过以下两种方式设置：

(a)  给定超固结比

 

(b)  给定

OCR = 1 表示正常固结状态，即材料过去的压缩程度不大于当前压缩程度。OCR > 1 表示超固结状态，即材料过去的压缩程度大于当前压缩程度。在超固结情况下，初始应力状态位于屈服面内（），材料将发生弹性变形；而在正常固结状态下，初始应力状态位于屈服面（）上，材料将发生塑性变形。

在Cam-clay模型中，位于屈服面外的初始应力状态不具有物理意义，因此在这种情况下，程序将修改 由用户给定为某一数值，使应力状态位于屈服面内或面上（见下表）。如果 ，则分析停止。

初始状态的确定

17-A.7.3             输入数据

该塑性模型仅适用于土体（KPLA = 3，在~CFMP）。弹性性质为各向同性。定义该模型的材料性质和参数如下表所示

材料性质和参数

初始固结压力通过OCR指定（当KP0 = 0时）；或通过P0指定（当KP0 = 1时）。

17-A.7.4             输出值

求解过程中会计算一些感兴趣的变量，这些变量可用于后处理，并作为状态变量保存（OUTRES, SVAR)

17-A.7.5             错误消息

求解过程中可能会出现以下错误消息：

·       CC plasticity：单元 … 的积分点 … 处没有初始平均压应力

·        

描述：	在单元的积分点处存在上述情况，这在此材料模型中是不允许的。
程序操作：	分析停止。
用户操作：	修改单元积分点处的初始应力状态，使其

·       CC plasticity：积分点 … 处的初始比体积小于或接近零， 属于单元 …

描述：	在单元的积分点处。
程序操作：	分析停止。
用户操作：	修改单元的应力状态 和/或材料参数 , N, λ y κ使

·      CC plasticity：塑性算法在单元 … 的积分点 … 处不收敛

描述：	无法在单元的一个积分点处求解塑性流动的控制方程。
程序操作：	分析停止。

17-A.7.6             假设和限制

·         适用于平面应变、轴对称和 3-D问题。

·         只能与PLANE182、PLANE183、 SOLID185、SOLID186和SOLID187单元一起使用。

·         应力为有效应力。

·         支持大变形和大位移。

·         必须定义初始应力状态，这可以通过先前的分析（例如弹性分析）实现，在该分析中保存最后一个 子步骤 (INISTATE, WRITE），或使用命令 ~TIS。应注意分析开始时的面力、体力和边界条件必须与初始应力相容，以避免收敛问题。

·         弹塑性模块是不对称的，因此必须在单元选项中指明单元刚度矩阵的不对称性 (KEYOPT(5) = 1)。此外，采用 Newton-Raphson 方法时，必须启用非对称方程组求解（NROPT, UNSYM).

 

17-A.8      Hoek-Brown破坏准则

17-A.8.1    简介

此工具旨在使CivilFEM用户能够建立满足Hoek-Brown原始（1980年）或修正（1992年）破坏模型的岩石地基模型。

计算将采用与非线性混凝土模块类似的过程，在每一步中使用 Drucker-Prager 材料，其性质随其中现有的荷载水平而变化。

在岩石力学中，通常根据岩体中存在的缺陷或不连续面，将岩体的结构状况分为四类：

图17-A.8‑1 岩体状态分类

同样，围压也分为低围压和高围压两种情况，其界限将在后文（第17-A.8.2.5节）中明确。

根据这些参数，下面将分别介绍原始和修正形式的Hoek-Brown准则；该准则仅适用于低围压以及各情况下所述的岩体结构状况。

17-A.8.2    Hoek-Brown准则

17-A.8.2.1                适用性

在第I组中，未风化岩石可视为完整岩石；当其行为为各向同性时，Hoek-Brown准则适用。

在第II组和第III组中，岩体受少数几组缺陷影响，因此其行为基本上是各向异性的，Hoek-Brown准则不适用于这些情况。

在第IV组和第V组中，岩体受多组缺陷影响，形成了可称为异质性中的均质性的状况，其行为为各向同性。

在这些情况下，Hoek-Brown准则再次适用，尤其是第V组；此时可能更适合采用第17-A.8.2.2节介绍的修正形式。

 

17-A.8.2.2                Hoek-Brown准则的公式

采用 Mohr-Coulomb 准则时，当 Mohr 圆在平面 s-t，最终与固有强度直线相切。

其中c为黏聚力， j 材料的摩擦角。

Hoek 和 Brown 为阐述其准则，选择 采用主应力的形式s₁ y s₃ （莫尔圆与横坐标轴的交点），其 表达式为：

 

其中：

sc	岩石基质的单轴抗压强度。
m, s	取决于岩石特性及其裂隙发育状况的 常数。

 

为了更好地适用于 裂隙非常发育的岩体（IV组和V组），Hoek 和 Brown 于1992年提出了其 修正模型，现介绍如下。

建议在 RMR 指标小于30时采用该模型。

新模型的表达式与 原模型非常相似。

 

当 n = 0.5 时与其重合

新参数 n 是岩体劣化程度的函数，其取值范围为 0.5～0.65。

可采用线性相关关系

   (RMR < 30)

除新指数 n 外，所有相关变量均与初始模型中采用的变量一致；n 的取值 完全取决于岩体的裂隙发育程度。

 

17-A.8.2.3                系数 m 和 s 的确定。

引入新系数 m_o, 该系数仅取决于岩石类型，在岩石完整（无 裂隙）时具有代表性，其定义见附表（Hoek 和 Brown）。

 

 

下表同样由 Hoek（1999）提出，其中 给出了系数 m 的取值_o，括号内的数字 为估计值。

 

 

系数 m 和 s 可通过以下形式的 相关关系获得

系数 s_o，表示完整岩石，始终取为1。

对于系数 A 和 B，当岩石处于天然状态、未受极端作用影响时，Hoek 和 Brown （1988）建议其取值分别为28和9。

 

17-A.8.2.4                单轴抗压强度的估算。

下表（Hoek，1999）为该参数在典型情况下可能取值的范围提供了参考

 

 

17-A.8.2.5                高围压下的岩体。

在高围压下，大多数岩石会由低围压下的脆性行为转变为 延性行为。

两种状态之间的转变点由 Mogi 进行了研究，他于 1966 年发现，当最大主应力与最小主应力之间满足 某一给定比值时，就会发生这种转变

同一作者根据完整岩石（三轴试验， 第 I 组）推导出该值 a = 3 或 4。在对 IV 组和 V 组岩石进行的试验中，过渡并不十分明确，Hoek 和 Brown (1980) 建议取值 a = 2 作为其破坏准则对岩体的适用界限。

 

图 17-A.8‑2 适用范围

17-A.8.3    主应力表达式

17-A.8.3.1                莫尔-库仑

Mohr Coulomb 准则定义为：

在主应力空间中：

其中：

 

 

图 17-A.8‑3 Mohr-Coulomb 破坏准则

 

17-A.8.3.2                Hoek-Brown

Hoek 准则可表示为：

 

Hoek 所采用的符号约定要求考虑 负拉应力

取下式时，该函数的切线为

抗压强度计算如下：

图 17-A.8‑4 Hoek & Brown 破坏曲线的切线

 

瞬时摩擦角为：

黏聚力：

 

 

17-A.8.3.3                材料破坏曲线的分段

曲线被划分为一系列由两个应力 S1 和 S3 定义的线段，位于 坐标轴 s₃。 用于计算摩擦角和黏聚力的应力取对应线段的中点值，即 (S1 + S3) /2

曲线内考虑两个临界点，第一个是与坐标轴 的交点 s₃，该点对应于岩石的抗拉强度。由于难以计算该点，目前不考虑抗拉 强度，并且 s₀ = 0.

这会导致岩石的极限抗拉强度被高估。直线是点 处的切线 s₃ = 0，这是对负应力所采用的值，会得到远大于实际值的拉应力。

 

图 17-A.8‑5 拉伸区线性化 及分段

第二个是点 s_3c ，在该点岩石表现出延性行为，正如 Mogi 的研究所示。

b_n= A_n . s_c

控制曲线分段的参数是摩擦角增量，d_f，用于确定由下式定义的角度 s_o 直到 s_3c 使用公式：

 其为以下变量的函数：s₃, 通过以下公式求解上、下角：

已知角度后，f_o y f_3c将曲线划分为每隔 d_f度。然而，分段曲线保存的数值为应力，因此，为将角度转换为应力， 需对每个角度计算

 

由此可以计算 通过以下公式：

 

对于每一对角度，生成两个应力s_inf y s_sup 曲线逼近所采用的应力为两者的平均值，即位于线段中点的应力。

利用该应力s₃ 代入：

由此可计算黏聚力。首先，需通过以下公式确定中间点处的摩擦角：

最后，该线段的黏聚力为

17-A.8.4    Hoek-Brown 模型的单元

Hoek-Brown 模块中考虑的单元为：

• PLANE 42
• PLANE 82

·         SOLID 45

• SOLID 92
• SOLID 95

每个子步 中用于 更改二维单元材料的应力，是每个节点最大主应力的平均值

对于三维单元，取

17-A.8.5    已进行的试验：Hoek-Brown 无限介质中的圆柱形孔洞

目的是验证处于 30 MPa 原位压应力场作用下的无限弹塑性介质中的应力分布。

材料属性如下：

E = 10¹⁰ Pa

n= 0.25

完整岩石的抗压强度 150·10⁶ Pa

Hoek-Brown 参数为：

RMR = 30

m = 2.94

s = 0.01

孔洞半径为 10 m，模型延伸至该距离的 10 倍。

 

 

图 17-A.8‑6 有限元模型

图 17-A.8‑7 径向应力

图 17-A.8‑8 切向应力

 

17-A.8.6    理论解

切向应力和径向应力分布的理论解由 Hoek 和 Brown 进行了研究（Hoek & Brown Underground Excavations in Rock, London 1982）。

在弹性区内：

其中

P0	各向同性原位应力的大小
re	塑性区半径
	r = r 处的径向应力e

在塑性区内：

其中，a 为孔洞半径，s_c为岩石的单轴抗压强度。数值 s_re y  r_e通过以下公式定义：

其中

下图给出了 CivilFEM 计算结果与理论解之间的比较：

图 17-A.8‑9 理论解与 CivilFEM 的比较

 

17-A.9               渗流

17-A.9.1    简介

此工具旨在使用 CivilFEM 自有单元，方便地求解流体通过 2D 或 3D 多孔介质的渗流问题。

此外，此工具还可简单、便捷地执行以下操作：

• 获取 2D 问题中的 饱和线。
• 将获得的孔隙水压力传递至边坡稳定性问题中。无需使两者所采用的有限 元网格相互一致。
• 计算通过边界的渗流量。

17-A.9.2    行为假设

假定渗流现象遵循以下工作假设：

• 流动过程中不存在 源项和汇项。
• 不可压缩流体。
• 满足 Darcy 定律，渗透系数可能具有各向异性 （沿 x、y、z 方向的渗透系数 k 不同）。
• 现象处于稳态 （相关变量不随时间变化）。

17-A.9.3    控制该现象的偏微分方程

从物理意义上，在区域内定义两个场：一个是标量场 H（水头或势），另一个是 矢量场 v（渗流速度）。假定在 S {0, xyz} 中建立正交基准，根据上一节所述的 假设可推导出以下条件：

结合这些表达式，可得到各向异性介质中稳态渗流的椭圆型方程。

 

17-A.9.4    边界条件

可能出现的条件有多种 类型：

不透水边界

将其视为一条流线，因此 通过该边界的渗流为零：

 （Neumann类型条件）

其中， 为垂直于该表面的法线。

 

恒定势边界（水头 固定）

该部分边界与水 直接接触：

H = H₀  （Dirichlet类型 条件）

其中 H₀ 是总水头。

 

大气压边界或渗出面

在渗出面上，压力为大气压，因此

H = 标高（Dirichlet 型条件）

二维问题中的高程轴为全局坐标系的 Y 轴，而 三维问题中的高程轴为全局 Z 轴。

 

17-A.9.5    计算结果

如前面各节所述，描述该现象的控制方程 以及边界条件均以水力势 H 表示，因此计算首先为每个节点 给出该值。

但是，对于当前分析的问题，这并不是唯一需要 了解的量。在大多数情况下，还需要确定孔隙水压力以及另外两个由水力势 得到的矢量量：

·         孔隙水 压力

                                                                                                                                       

                                                                                                                                                    

·         场梯度

·         流速

 

17-A.9.6    求通过某一表面的渗流量。

 

对于三维问题，通过表面 S 的流量 或渗流量可通过积分求得。

 

其中， 为表面的法向量。

如果问题为平面问题，则单位流量通过曲线积分求得

 

当弧段 L 为平行于坐标轴的 线段时，该积分可化为简单积分。

-          水平线段 [x₀, x₁]

 

-          竖直线段 [y₀, y₁]

 

17-A.9.7    饱和线（二维问题）

所要解决的问题是求得浸润线。为此，需要定义一组区域，并按照如下所示，将岩土数据库中的参数“rusi”值（材料属性）赋予这些区域：

rusi = 0	无渗流（堆石体）时。
rusi = 1	否则。

此外，用户还必须提供一组 ANSYS 线，使其在初始近似中构成浸润线（程序将以此线为起点，通过迭代优化过程求得最终浸润线）。

浸润线有两个端点，分别可满足两种互斥条件之一。

·         固定：点 A 的情况

·         沿滴水面滑动：点 B 的情况。

 

图 17-A.9‑1 带浸润线的渗流

 

17-A.9.7.1                工作原理

几何定义

必须为每个区域指定一种材料，其定义中将包含常数“rusi”的值（命令 ~CFMP).

初始浸润线必须定义为一组线（折线），并将问题的未知量（浸润线的纵坐标）设置在这些线的顶点处。计算过程中，CivilFEM 将迭代这些中间点的位置，直至获得浸润线的最终几何形状。

必须在浸润线的起点和终点处指定应满足的条件（固定点或滑动点）。

边界 条件

所采用的边界条件分别设置在构成渗流区域的各条线上，这些边界条件可分为三类：

• 恒定水力势 （浸没面）
• 水力势等于标高 （自由面）
• 流量为零（防渗面）

17-A.9.8    结果导出

渗流模块将结果写入文件 文件名.PRESS 其中包含模型各节点的结果。这些结果可用于边坡稳定性计算或有效应力表示。

文件格式如下：

BLOQUE 1

版本

节点数

单元数

BLOQUE 2

模型中的每个节点

X 坐标、Y 坐标、Z 坐标、孔隙水压力

BLOQUE 3

每个单元的 I、J、K 节点。

BLOQUE 4

每个单元的相邻单元（3 个）

修改该文件可能导致不可预知的结果，因此不建议对其进行操作。

17-A.9.9    CivilFEM 渗流单元

此工具旨在为 CivilFEM 提供必要的功能，使用户能够通过 CivilFEM 专用的新型有限元单元（PLANE42-SEEP 和 SOLID45-SEEP 渗流单元）处理多孔介质中的渗流问题。这些单元的内部工作原理与 ANSYS 的结构单元 PLANE42 和 SOLID45 完全 相同，因此可以使用该程序中用于前处理和后处理的全部可用资源。

CivilFEM 单元与 ANSYS 结构单元之间的对应关系如下表所示：

 

 

 

17-A.9.9.1                功能

问题 类型

待求解的问题可以是二维或三维问题。

边界条件

渗流问题中使用的边界条件包括：

-          不透水边界：

-          恒定势边界：

      命令：~DLHEAD y ~DAHEAD.

-          大气压力边界或滴水面：

H = cota

命令 ~DLSEEP y ~DASEEP

 

17-A.9.9.2                结果

CivilFEM 提供的结果包括：

-     所有节点的 势 H。

-     所有点的 孔隙压力 u（命令 ~ISOBAR).

该量按下式获得：u = H – coordenada_cota

-     单元节点处的 势梯度（¶H/¶s）：单元节点处的分量 分别沿 x、y、z 方向。

-     单元节点处的 速度矢量：分量分别沿 x、y、z 方向。

-     通过单元各个面的 流量或流率。

 

17-A.9.9.3                模型

模型的建立使用 ANSYS 结构分析功能，并采用 PLANE42 和 SOLID45 单元类型。

 

 

在求解阶段，CivilFEM 会将其替换为相应的 CivilFEM 单元和自由度。

这样，ANSYS Structural 产品的全部功能均可用于模型网格划分。

                                     

单元类型

 

求解

必须使用 CivilFEM 命令 ~LPSOLVE，而不是命令 SOLVE.

 

后处理

使用 CivilFEM 渗流单元时，将生成结果文件 作业名称。SEEP。因此，要访问这些 结果，必须执行命令 FILE, 作业名称,SEEP ，随后执行命令 SET,1 以访问已求解的载荷步。CivilFEM 渗流单元仅支持一个 荷载步，且 始终为第一个载荷步。

使用 CivilFEM 的 HEAD 自由度，通过 ANSYS 常用命令可获得图形和列表输出。例如 PLNSOL，HEAD 可绘制势能等值线。

命令 ~LPRNSOL 可列出节点结果，使用命令 ~PLSEEP y ~ISOBAR 可绘制 渗流结果。

所有单元结果均通过将在积分点获得的相应结果 插值到节点上计算得到。请参见第 13.6 节 “节点和质心数据评估”，该内容位于手册 ANSYS 理论参考

 

17-A.9.9.4                单元说明

 

·         PLANE42-SEEP

该单元是结构单元 PLANE42 的变体，可用于求解渗流问题。

PLANE42-SEEP 可作为具有二维渗流功能的平面单元使用。 该单元有四个节点，每个节点具有一个自由度，即压力水力势。

该单元适用于二维稳态渗流分析。

图 17-A.9‑2 PLANE42-SEEP 单元

输入 数据

该单元的几何形状、节点位置和坐标系如图 17-A.8-2 所示。该单元由四个节点以及材料属性 KXX、KYY、KZZ 定义。 材料的正交各向异性方向与全局坐标方向一致。

单元名称

PLANE42-SEEP

节点

I, J, K, L

自由度

HEAD（水头）

材料 属性

KXX, KYY, KZZ

 

输出数据

与该单元相关的输出解有两种形式：

·         包含在全局节点解中的节点压力势

·         如《单元输出定义》中所示的附加单元输出

可以使用单元 PLANE42-SEEP 的 HEAD 标签获得节点解。

Element Output Definitions 表使用以下 符号：

参数可通过分量名称方法 [ETABLE, ESOL]. 结果文件名 为 Jobname.SEEP，用户必须使用命令 FILE,文件名,SEEP 读取此结果。

PLANE42-SEEP 的等效列指示获取结果时必须使用的 标签。

这些参数可在结果文件中获得。

 

表 17-A.9‑1 PLANE42-SEEP 单元输出定义

 

下表显示了使用命令获取结果的 数值序列 ETABLE y ESOL.

 

表 17-A.9‑2 PLANE42-SEEP。ETABLE 和 ESOL 命令的标签及序列号

 

Name

“Element Output Definitions” 中定义的 输出量

Item

命令的默认 参数标签 ETABLE.

FCn

面编号 n 的 求解参数数值序列。

 

假设和限制

单元的面积不能为零或负值。单元必须位于 X-Y 平面内。

 

产品 限制

此单元只能与 岩土工程模块配合使用。

 

• SOLID45-SEEP

此单元是结构单元 SOLID45 的一种变体，能够求解渗流问题。 

SOLID45-SEEP 可作为具有 三维渗流功能的实体单元使用。该单元有八个 节点，每个节点具有一个自由度，即压力势。

该单元适用于三维 稳态渗流分析。

图 17-A.9‑3 SOLID45-SEEP 单元

 

输入 数据

该单元的几何形状、节点位置和坐标系如图 17-A.8-3 所示。该单元由八个节点以及材料属性 KXX、KYY、KZZ 定义。材料的正交各向异性方向对应于全局坐标方向。

单元 名称

SOLID45-SEEP

节点

I, J, K, L, M, N, O, P

自由度

HEAD

材料 属性

KXX, KYY, KZZ

 

输出数据

与该单元关联的输出解 有两种形式：

• 节点压力势 包含在全局节点解中
• 附加单元输出 如下所示 单元输出定义

节点解可 使用 SOLID45 单元的 UX 标签获得。

参数可通过组件名称方法 [获得ETABLE, ESOL]. 结果文件名为 Jobname.SEEP，用户必须 使用命令 FILE,文件名,SEEP 读取此结果。

SOLID45 的等效列指出了 获取结果时必须使用的标签。

这些参数位于 结果文件中可用.

 

表 17-A.9‑3 SOLID45-SEEP 单元 输出定义

名称（与 SOLID45 的等效项）	定义
EPELX, EPELY, EPELZ	X、Y 和 Z 方向的梯度
SX, SY, SZ	X、Y 和 Z 方向的速度
PRES	穿过某一面的水流量

 

下表显示了 使用以下命令获取结果的数值序列ETABLE y ESOL

表 17-A.9‑4 SOLID45-SEEP ETABLE 和 ESOL 命令的标签及序列号

Name	Item	FC1	FC2	FC3	FC4	FC5	FC6
PRES	SMISC	1	2	3	4	5	6

 

Name

在以下位置定义的输出单元输出定义

Item

命令的默认参数标签ETABLE.

FCn

面编号 n 的求解参数数值序列。

 

假设和限制

不允许使用零体积单元。单元必须按照 SOLID45 中所示的方式编号。当单元编号不正确时，通常会发生这种情况。

所有单元必须具有八个节点。通过将 K 和 L 的节点编号定义为重复，并将 O 和 P 的节点编号重复，可以形成棱柱形单元（参见三角形、棱柱形和四面体单元）。也可以使用四面体形状（对于四面体单元，多余的形状会自动删除）。

 

产品 限制

此单元只能与 岩土工程模块配合使用。

 

17-A.10         土压力

17-A.10.1                简介

CivilFEM 提供了一系列命令，可将与土的主动土压力、被动土压力和静止土压力以及土重相对应的压力和力直接施加到模型中选定的单元上。

17-A.10.2                单元类型

CivilFEM 可将土压力施加到以下单元上：

二维梁单元

            BEAM3, BEAM23, BEAM54

三维梁单元

            BEAM4, BEAM24, BEAM44, BEAM188, BEAM189

板单元

            SHELL43, SHELL63, SHELL93, SHELL181

二维实体单元

            PLANE2, PLANE42, PLANE82, PLANE145, PLANE146

3D实体单元

            SOLID45, SOLID46, SOLID64, SOLID65, SOLID73, SOLID95, SOLID147

面单元

            SURF153, SURF154

17-A.10.3                土体特性

计算土压力时，需要定义土体（命令 ~TERDEF).

土体可以包含不同的土层。

17-A.10.4                土压力 土体

计算假定

土体主动土压力和被动土压力的计算采用库仑理论（如果在土层定义中指定），或采用材料数据库中的土压力系数（默认采用朗肯土压力）。静止土压力采用材料数据库中的土压力系数。假定该土压力是由于土楔滑动产生的，土楔由墙背或结构构件以及通过该构件墙趾的平面限定。计算基于以下假定：

·         土体无黏聚力。

·         土体表面为平面，与水平面形成角度 b，并在按坡面测量的单位面积上施加均布荷载 q。

·         不考虑因存在限制最大土压力土楔的其他墙体而引起的土压力减小，也不考虑可能的筒仓效应。

·         独立构件上的土压力按构件的实际宽度计算；也就是说，对于宽度较小的构件，不通过增大宽度来放大相应的土压力。

计算过程

选定构件上土压力的计算和施加包括以下步骤：

1.    确定角度 a，即每个选定构件的角度。角度 a 是施加土压力的构件表面与土体水平面之间的夹角，测量方向为施加土压力的相反一侧。

对于梁单元，该角度根据构件的轴线计算，即由两端节点确定的方向（对于3D梁单元，还考虑第三节点或方向角 THETA），因此对于变截面梁，角度 a 与实际值会有轻微差异。

对于板单元，角度 a 根据构件前三个节点形成的平面计算。与梁单元一样，对于变厚度板，该角度会与实际值有所差异。

对于2D和3D实体单元，角度 a 由将施加土压力的面上的前三个节点确定，因此始终为实际角度。

为使角度 a 处于土压力系数公式的有效范围内，并使物理问题具有实际意义，必须满足：

a + j £ 180º

a – d > 0º

a + b < 180º

2.    土压力系数的计算。角度计算完成后，a 并根据为土体定义的数值，可从数据库中获得土压力系数，或者按照第 章所述理论，对库仑土压力情况进行计算，17-A.4.5 墙面上的土压力计算

下文中，无论计算时选用的是主动土压力、被动土压力等哪一种土压力系数，均将其记为 K：

3.    法向压力和切向压力的计算。

对于干燥且无附加荷载的土体，位于土面以下深度 z 处某点的压力值由下式给出：

 

该单位深度压力与其作用面形成角度 d（土体-墙体摩擦角）与单元面。

因此，微分单元所承受的总力为：

因此，沿该单元方向表示的微分单元上的压力为：

图 17-A.10‑1 土压力的施加

将该压力投影到单元上，得到：

 

如果土体上存在沿边坡单位长度、数值为 q 的均布荷载，则位于土面以下深度 z 处某点的压力值为：

 

因此，单元上的法向压力和切向压力表达式为：

对于积水土体，地下水位以下的压力采用重度计算：

设地下水位深度为 f，则深度 z（大于地下水位深度）处的压力为：

因此，单元上的法向压力和切向压力为：

4.    法向压力和切向压力在单元上的施加。

梁单元和 Surf153

对于这些单元，将上一节计算得到的法向压力和切向压力直接施加到各个单元上。对于梁单元，CivilFEM 根据梁属性（Beam Properties）中定义的横截面获取单元宽度。对于截面不恒定的梁，取两个截面宽度的算术平均值。

板单元、实体单元和 Surf154

法向压力直接施加在这些单元上，而切向压力则通过沿切向的节点力施加。节点力根据所选单元，通过将切向压力乘以每个节点所对应的面积计算得到。在与角度变化处重合的节点上，a施加的力在大小和方向上，均为对于每个角度值该节点所对应力的平均值a.

17-A.10.5                土体自重引起的侧压力

计算假定

土体自重引起的侧压力按计算构件以上土体的重量以及地表面上的附加荷载进行计算。计算采用以下假定：

·         在构件上施加法向压力和切向压力，其合力为竖向力；此外还考虑水压力，水压力仅具有垂直于构件的法向分量。

·         地表面为平面，与水平面形成角度 b，并沿边坡表面测得每单位面积的附加荷载为 q。

·         不考虑由于与其他构件之间可能产生的拱效应而导致的侧压力减小。

·         孤立构件上的侧压力按构件的实际宽度计算。

计算过程

选定构件上土体自重引起的侧压力的计算和施加包括以下步骤：

1.    确定角度 a，即每个选定构件的角度，方法同土压力计算中的角度确定。

2.    计算法向压力和切向压力。

土体自重在深度为 z 的点处产生的压力，在干燥且无附加荷载的土体情况下，由下列表达式给出：

将该压力投影到构件上可得：

图 17-A.10‑2 侧压力的施加

如果地表面上存在均布荷载，其值为 q，即每单位边坡长度上的荷载，则地表面以下深度为 z 的点处的压力为：

在饱和土体中，地下水位以下的压力采用重度计算：

令 f 为地下水位深度，则深度为 z（大于地下水位）包括：

3.    在单元上施加法向压力和切向压力。压力的施加方式与土压力相同。

 

17-A.11         土体初始应力

17-A.11.1                应力场描述

无变形初始应力通过命令输入 CivilFEM 中的 ~TIS.

在岩土工程中，通常采用一种初始应力场，其中水平分量与竖向分量成比例。

图 17-A.11‑1 微分单元中的应力

参数 k_o称为静止土压力，也称为水平土压力系数或侧压力系数。其取值取决于材料及其固结程度。

17-A.11.2                系数 K 的值₀

17-A.11.2.1             弹性理论

 

根据关系式

在受均布荷载（例如自重）的 Boussinesq 半空间（土体水平）中，可验证关系式：

 

因此

 

或者等价地，

 

因此

 

附表和图表表示了系数 k_o 随泊松比的变化，结果表明弹性理论对该系数给出的值从不大于 1。

 

 

 

17-A.11.2.2             与摩擦角的关系

对于正常固结土，系数 k_o 可通过下式近似求得

 

Jaky 提出的公式。

角度 j 是材料的有效应力摩擦角。

 

17-A.11.2.3             试验值

最早对该系数进行测量的研究由 Karl Terzaghi 完成，他在1920年至1925年间发表了以下数值。

 

 

砂土的建议值有时似乎过小，而黏土的预测值则过高。 目前，由于缺乏更完善的研究，以下数值被认为更为合适：

 

 

上表列出了 j y  u 为获得静止土压力系数所需的数值。

由于自身性质，这两个方程得到的系数 k_o 小于1。然而，在强固结黏土中，通过间接方法测得的系数 k_o 介于2和4之间（Skempton 1961，Bishop 1965），其他类型材料中也存在这种情况。

这些数值无法通过前面各节介绍的方程复现。

17-A.11.3                初始应力

根据给定的地形，初始应力将通过以下表达式计算：

其中 h 为所考虑点的深度，g 材料的重度，k_o 为静止土压力系数。

如果该点上方的土柱由 n 种材料组成，则上述表达式可写为：

17-A.11.4                获取初始应力文件

命令 ~TIS 生成一个 ASCII 文件，其中包含单元中心处的初始应力值。这些值相对于活动（笛卡尔） 参考系计算。

初始应力将通过读取初始应力文件输入模型。该文件必须与单元类型 兼容，这通过参数进行控制：KTIS 命令的 ~TIS

由于命令 INISTATE 已不再支持， 建议在 182、183、 185、186 和 187 号单元类型上定义初始应力。

建议使用获得的初始应力和重力求解第一个荷载状态。这样， 输入的初始应力将进行调整，以便形成与重力平衡的应力状态。因此，输入的 初始应力与平衡状态下的应力之间可能会出现轻微差异。