CN109918841A - 基于有限单元与极限平衡的土体稳定性评价方法 - Google Patents

Abstract

本发明属于岩土工程设计与施工中的土体稳定性评价技术领域，具体涉及一种基于有限单元与极限平衡的土体稳定性评价方法专利申请事宜。所述土体稳定性包括边坡稳定性、地基承载力和支挡结构土压力三类问题，具体包括：1）地形测绘、地质勘察与有限单元模型构建，2）总体约束方程构建，3）自变量，4）目标函数，5）评价与结论等步骤。本申请以有限单元与极限平衡理论为基础，能够适应土体性质的不连续、不均匀、在破坏时呈现的大变形等特性，且不需要引入最危险滑裂面假定、并能够同时适用于土体稳定性评价中的边坡稳定性、地基承载力和支挡结构土压力三类问题，其评价方法在理论上更为严格，结论与工程实践也更为符合。

Description

基于有限单元与极限平衡的土体稳定性评价方法

技术领域

本发明属于岩土工程设计与施工中的土体稳定性评价技术领域，具体涉及一种基于有限单元与极限平衡的土体稳定性评价方法专利申请事宜。

背景技术

土体稳定性评价包括三类问题：边坡稳定性、地基承载力和支挡结构土压力。对于这三类问题，Columb和Rankine的方法被广泛应用到边坡稳定性和地基承载力评价中，属于一种极限分析的方法。在土体稳定性评价中，尽管还有有限元等其它方法，但极限平衡法因其理论简单、概念清晰、求解方便等优点而被广泛使用。

经典极限平衡法的计算步骤可归纳为两步：首先是假定一个滑裂面，计算该面发生失稳的安全系数；然后是寻找最小的安全系数，也就是寻找临界滑裂面。为了对土体的稳定性进行评价，首先需要假定土体因某种扰动因素从目前的稳定状态进入极限平衡状态，此时，土体内出现一个假想的滑裂面，在该滑裂面上，每一点的法向应力和切向应力都满足Mohr-Columb强度准则。

Fellenius于1927 年提出了边坡稳定分析的圆弧滑动分析法，也就是瑞典圆弧法。该方法中，滑裂面假定为圆弧，土条底面的法向应力被看作是土条重力在土条底面法线方向的投影，因此对圆心取矩时不出现，使计算工作量大大简化。但该法忽略了条间力的作用，不能满足所有静力平衡条件。当滑裂面不是圆弧时，该法不再适用。

Bishop法在瑞典圆弧法的基础上进行了发展，该法仍然假定滑裂面为圆弧，但考虑了条间力的影响。通过土条的竖向力平衡及力矩平衡推求出安全系数。Bishop法考虑了条间力，比瑞典圆弧法前进了一大步，对于均质土坡，目前普遍认为该法可以得到比较合理的解。但是，滑裂面为圆弧的假定限制了该法的应用范围。

Janbu法假定滑裂面为任意形状，且假定条间力的作用点位置（即推力线位置）为己知。对每一土条，可列出两个静力平衡方程和一个力矩平衡方程，属于严格方法。但其推力线的假定必须符合条间力的合理性要求。尽管该法的滑裂面为任意形状，但必须事先假定，因此该法仅适用于滑裂面位置、形态均基本确定的边坡稳定分析。

不平衡推力传递法是按照折线滑动面将滑动土体分块，假定块间力的合力与上一土块底面平行。该法对滑体仅进行分块，而没有分条，因而计算工作量大大减小。同样，该法也仅适用于滑裂面位置、形态均基本确定的边坡稳定性分析。

上述所列仅是极限平衡方法范畴中用于边坡稳定性评价的几种典型方法。总体上，由上述方法可知：在经典极限平衡法中，无论是假定滑裂面为圆弧形、还是任意形状，也无论是严格方法、还是非严格方法，都必须首先假定滑裂面的形态（如 Fellenius法和Bishop法），或者同时假定滑裂面的形态和位置（如Janbu法和不平衡推力传递法）。但这种假定在实际应用中并不完全准确，使得土体极限平衡分析法在实际应用中受到一定限制。

随着计算机及工程实践的发展，更多的工程尝试采用有限元法对土体稳定性进行分析，与传统的极限平衡法相比，有限元法的优点有滑裂面的形状或位置不需要事先假定，并能够提供应力与变形的全部信息，理论体系相对严密，但有限元法需要引入土体中最不易弄清的本构关系，使得有限元法在实际应用时受到了很大限制。

因此，设计一种不需要引入最危险滑裂面假定、不需要引入土体中最不易弄清的本构关系的土体稳定性评价方法，对岩土工程中边坡稳定性、地基承载力和支挡结构土压力工程的设计与施工都有着重要意义。

发明内容

本申请目的在于提供一种基于有限单元与极限平衡的土体稳定性评价方法，该方法不需要引入最危险滑裂面假定、不需要引入土体中最不易弄清的本构关系、并能够同时适用于土体稳定性评价中的边坡稳定性、地基承载力和支挡结构土压力三类问题，其评价方法在理论上更为严格，结论与工程实践也更为符合。

本申请所采取的技术方案详述如下。

一种基于有限单元与极限平衡的土体稳定性评价方法，所述土体稳定性包括边坡稳定性、地基承载力和支挡结构土压力三类问题，所述方法具体包括以下步骤：

（一）地形测绘、地质勘察与有限单元模型构建

针对待评价土体，进行地形测绘、地质勘察和土体物理力学参数测定，并建立土体稳定性评价的有限单元模型；

需要简要解释强调的是：

地形测绘应覆盖土体稳定性评价的范围；

地质勘察应包括：查明评价范围内土体的地下水位线、具有不同物理力学性质的各类土体的空间分布等，并为评价范围内各类土体物理力学参数的测定提供必要的试样试件；

需要测定的土体物理力学参数包括：土体内摩擦角、粘滞系数、天然容重、干容重、饱和容重、含水率等；

根据各类土体的空间分布、地下水位线以及工程实践对计算精度的需求，将评价范围内的土体划分为有限个单元，形成土体稳定性评价的有限单元模型；实际划分中：

对平面二维问题，可采用三角形和/或多边形单元进行划分；

对空间三维问题，可采用四面体和/或多面体单元进行划分；

（二）总体约束方程构建

需要说明的是，为了便于直接采用弹性力学中的有关结论，本申请中节点应力的大小与方向均采用弹性力学中的有关规定，即正应力以拉应力为正，切应力以使微分体产生逆时针转动趋势为正；

所述总体约束方程由屈服准则约束方程、节点平衡约束方程和应力已知边界条件约束方程共同构成；具体而言：

（2.1）构建屈服准则约束方程

根据测定出的土体物理力学参数及相应的土体强度理论，建立评价范围内各类土体以应力形式表达的屈服准则；

以平面问题、Mohr-Columb屈服准则为例，设有限单元模型的节点编号依次为，第个节点所对应的节点应力分量分别为、和，该点在某一方向上的正应力和切应力分别为和，则有：

， 公式(1)；

， 公式(2)；

进一步，设测定出的土体内摩擦角为，并记，根据Mohr-Columb屈服准则，其节点应力在任一方向上的正应力和切应力均应满足Mohr-Columb准则，可得到：

， 公式(3)；

而公式(3)等价于

， 公式(4)；

若记为：

， 公式(5)；

， 公式(6)；

则可得到公式(3)等价于；

进一步将公式(1)与公式(2)代入公式(5)的左边表达式，则可得到：

， 公式(7)；

公式(7)中：

；

根据公式(7)，可得到公式(5)等价于：

， 公式(8)；

根据Mohr-Columb准则，对任意的 QUOTE ，公式(8)均应成立，进一步可得到：

， 公式(9)；

由于土体不能承受拉应力，因此可限定及，此时可得到公式(9)等价于

， 公式(10)；

据上可得到，在限定及的条件下，公式(5)等价于公式(10)；

类似地，将式(1)与式(2)代入公式(6)的左边表达式，并限定及，可得到公式(6)等价于

， 公式(11)；

综合公式(10)与公式(11)，对平面问题，Mohr-Columb准则可综合地写为

， 公式 (12)

公式(12)中，；

也即，公式(12)即为本申请的土体稳定性评价方法的屈服准则约束方程；

需要说明的是，在公式(12)中，对于每个节点，Mohr-Columb准则有4个不等式约束方程；若节点总数为，则共有个屈服准则不等式约束方程；

对三维空间问题、或其它以应力形式表达的屈服准则问题，均可采用上述方法，类似地建立屈服准则约束方程；

（2.2）节点静力平衡约束方程构建

对于平面问题，单个节点的静力平衡方程只包括 和 两个方向的平衡方程；对于三维空间问题，单个节点的静力平衡方程则包括、和 三个方向的平衡方程；

以平面问题为例，可采用如下方法建立节点静力平衡约束方程：

如图1所示，设围绕节点的所有单元按逆时针方向编号依次为，围绕节点的所有节点按逆时针方向编号依次为；设节点 的坐标为 ，节点的坐标为，节点的坐标为，节点的应力为、和；设相邻节点之间的正应力及切应力均按直线变化；此时可得到如下结论：

节点的正应力对节点在方向力的贡献为：

， 公式(13)；

节点的切应力对节点在方向力的贡献为：

， 公式(14)；

节点的正应力对节点在方向力的贡献为：

， 公式(15)；

节点的切应力对节点在方向力的贡献为：

， 公式(16)；

对平面问题，取单元厚度为单位厚度，设土体在和方向上的单位体力分别为与，单元的面积分别为，可得到：

单元的体力对节点在方向的总合力为：

，公式(17)；

单元的体力对节点在方向的总合力为：

，公式(18)；

据上，可得到关于节点的静力平衡约束方程为：

， 公式(19)

公式(19)中，；

其中：

为节点的正应力对节点在X方向力的贡献；

为节点的切应力对节点在X方向力的贡献；

为节点的正应力对节点在Y方向力的贡献；

为节点的切应力对节点在Y方向力的贡献；

为单元的体力对节点在X方向的总合力；

为单元的体力对节点在Y方向的总合力；

也即，公式(19)即为本申请的稳定性评价方法的静力平衡约束方程；

需要说明的是，在公式 (19)中，对每个节点有2个静力平衡等式约束方程；若单元模型的节点总数为，则共有个静力平衡等式约束方程；

对三维空间问题，可采用上述方法，类似地建立静力平衡等式约束方程；

（2.3）应力已知边界条件约束方程构建

以平面问题为例，对边坡稳定性和支挡结构土压力问题，在应力已知边界上，令其位于应力已知边界上节点的应力分量等于其已知量；

对地基承载力问题，则需要把所需的地基承载力设定为应力已知边界，并与其它应力已知边界共同构成地基承载力问题的应力已知边界；令其位于应力已知边界上的节点的应力分量等于其已知量；

在应力已知边界上，设已知的节点正应力和切应力分别为和 ，节点处的边界方向为，则有：

， 公式(20)；

公式(20)中，与应力已知边界条件有关；

即，公式(20)即为本申请的稳定性评价方法的应力已知边界条件约束方程；

需要说明的是，公式(20)中，应力已知边界条件约束方程的个数与应力已知边界条件及单元划分有关；

对三维空间问题，可采用上述方法，类似地建立应力已知边界条件约束方程；

（三）自变量

以所有节点的应力分量为自变量；

以平面问题为例，自变量为：

， 公式(21)；

公式(21)中，；

对三维空间问题，可采用上述方法，类似地建立自变量；

（四）目标函数

平面问题中，对边坡稳定性问题，目标函数是判断满足步骤（二）的总体约束方程的自变量是否存在；若存在，判定该土体为稳定的；若不存在，判定该土体为不稳定的；

平面问题中，对地基承载力问题，目标函数是判断满足步骤（二）的总体约束方程的自变量是否存在；若存在，判定该地基承载力满足要求；若不存在，判定该地基承载力不满足要求；

平面问题中，对支挡结构的土压力问题，设待求边界上以节点应力为变量表示的土压力的大小为：

， 公式(22)；

公式(22)中， QUOTE 与待求边界上的应力有关；

对极限最小土压力，目标函数是求解满足步骤（二）的总体约束方程条件下的P的最小值；

对极限最大土压力，目标函数是求解满足步骤（二）的总体约束方程条件下的P的最大值；

对三维空间问题，可采用上述方法，类似地建立本申请的稳定性评价方法的目标函数；

（五）评价与结论

由步骤（二）的总体约束方程、步骤（三）的自变量、步骤（四）的目标函数共同构成稳定性评价方法的数学模型；

借助计算机程序情况下，采用最优化求解的方法，对上述数学模型进行计算求解；

对边坡稳定问题，是判定数学模型是否有解；若数学模型的解存在，判定该土体为稳定的；若解不存在，判定该土体为不稳定的；

对地基承载力问题，是判定数学模型是否有解；若数学模型的解存在，判定该地基承载力满足要求；若解不存在，判定该地基承载力不满足要求；

对支挡结构土压力问题，是求出土压力P的最小或最大值；P的最小值即为极限最小土压力，P的最大值即为极限最大土压力。

总体而言，本申请所提供的土体稳定性评价方法，以有限单元与极限平衡理论为基础，避开了土体工程中最不易弄清的本构关系，能够适应土体性质的不连续、不均匀、在破坏时呈现的大变形等特性，且不需要引入最危险滑裂面假定、并能够同时适用于土体稳定性评价中的边坡稳定性、地基承载力和支挡结构土压力三类问题，其评价方法在理论上更为严格，结论与工程实践也更为符合。在工程实践上具有更好的通用性，具有较好的实用价值。

附图说明

图1是围绕节点i的单元及节点局部编号示意图；

图2是支挡结构的土压力问题示意图；

图3是支挡结构的土压力问题的有限单元模型；

图4是围绕内部节点的节点局部编号示意图；

图5是围绕拐角节点的节点局部编号示意图；

图6是围绕上、下边界中间节点的节点局部编号示意图；

图7是围绕左、右边界中间节点的节点局部编号示意图。

具体实施方式

下面结合实施例对本申请做进一步的解释说明，以使本领域技术人员能够更加清晰地了解本申请的设计方案。

实施例

本实施例仅以土体稳定性评价中的支挡结构土压力问题为例，具体例为竖直支挡结构单位宽度内的水平向最小土压力问题，就本申请的土体稳定性评价方法介绍如下。

（一）地形测绘、地质勘察与有限单元模型构建

如图2所示，地形测绘：设定土体上表面水平，右侧范围无限大；地质情况：设定土体性质均匀单一，无地下水分布；土体物理力学参数：设定容重𝛾为18.5kN/m³，内摩擦角𝜑为37°；土体稳定性评价范围：设定竖向深度5m、水平长度为10m。

由于所需要解决的属平面二维问题，为便于计算，采用水平矩形单元对评价范围内的土体进行单元划分，将稳定性评价范围内的土体沿竖直方向划分为5个单元、6个节点，沿水平方向划分为10个单元、11个节点，共计50个单元、66个节点，并形成有限单元模型（如图3所示）。

（二）总体约束方程构建

需要说明的是，为了便于直接采用弹性力学中的有关结论，节点应力的大小与方向均采用弹性力学中的有关规定，即正应力以拉应力为正，切应力以使微分体产生逆时针转动趋势为正。

所述总体约束方程由屈服准则约束方程、节点静力平衡约束方程和应力已知边界条件约束方程共同构成。

设定所划分单元模型的节点编号依次为，其第个节点所对应的节点应力分量分别为、和。

（2.1）屈服准则约束方程构建

采用Mohr-Columb屈服准则，设第个节点的平面应力为为、和，则屈服准则约束方程为

， 公式 (23)

公(23)中，。

（2.2）节点静力平衡约束方程构建

如图4所示，设围绕内部节点的按逆时针方向所有节点的局部编号依次为，设相邻节点之间的正应力及切应力均按直线变化；可得到：

X方向的静力平衡方程为：

， 公式(24)；

Y方向的静力平衡方程为：

， 公式(25)。

对四个拐角节点（即，如图5所示），

X方向的静力平衡方程为：

， 公式(26)；

Y方向的静力平衡方程为：

， 公式(27)。

对上、下边界中间节点（即，如图6所示），

X方向的静力平衡方程为：

， 公式(28)；

Y方向的静力平衡方程如下：

， 公式(29)。

对左、右边界中间节点（即，如图7所示），

X方向的静力平衡方程为：

， 公式(30)；

Y方向的静力平衡方程为：

， 公式(31)。

综合公式（24）~公式（31），可得到节点静力平衡约束方程如下：

， 公式(32)；

公式(32)中，。

需要说明的是，公式（32）中，有66个X方向节点静力平衡方程和66个Y方向节点静力平衡方程，共计有132个节点静力平衡等式约束方程。

（2.3）应力已知边界条件约束方程构建

在本实施例中，土体上表面的Y向应力及上表面上的切应力应为零，因此可得到

， 公式(33)；

公式(33)中，。

也即，公式(33)即为本实施例的应力已知边界条件约束方程。

需要说明的是，公式（33）中共有22个应力已知边界条件约束方程。

综上，总体约束方程包括屈服准则约束方程公式（23）、节点静力平衡约束方程公式（32）、和应力已知边界条件约束方程公式（33）。

（三）自变量

以所有节点的应力分量为自变量；具体有：

， 公式(34)；

公式(34)中，，共计198个自变量。

（四）目标函数

针对支挡结构土压力问题，待求边界上以节点应力为变量表示的单位宽度的土压力P可表示为

， 公式 (35)；

公式(35)中，P采用土力学中的习惯，土压力以压力为正；另由前面的设定可知，式中dy为1.0m。

P即为本实施例的目标函数。

（五）评价与结论

由步骤（二）的总体约束方程、步骤（三）的节点应力自变量和步骤（四）的目标函数共同构成稳定性评价的数学模型；

借助计算机程序情况下，采用最优化求解的方法，可得目标函数P的最小值为44.4966kN。

作为对比，在上述设定条件下，根据经典的库仑主动土压力理论，计算可得到本实施例对应的库仑主动土压力为57.4850kN。从这一结果对比可以看出，经典库仑主动土压力所得结果约为本申请评价方法所得的极限最小土压力的1.29倍。

对此结果差异，分析可知：

基于相关假定的经典库仑主动土压力对应的土体应力状态为：①土体应力满足静力平衡方程；②土体应力满足Mohr-Columb屈服准则；③土体某点的应力分量是土体埋深的函数，而与该点离开支挡结构的距离无关。根据本申请所提供的评价方法，土体内某点应力分量不仅与土体埋深有关，而且与该点离开支挡结构的距离有关。因此，从理论上来说，本申请提供的土体稳定性评价方法在理论上更为严格，与实际也更为符合。

实践上，本申请提供的土体稳定性评价方法可对某些实际问题给出很好的解释，如堆料仓内的堆料对仓壁产生的侧压力问题，当仓内堆料足够高时，仓壁所受到的侧压力的大小并不随着堆料高度呈直线变化。这一现象不能从经典库仑主动土压力理论得到合理解释，而本申请所提供的评价方法则可较好地解释这一现象。

综上，本申请所提供的土体稳定性评价方法，在理论基础和工程实践上均具有优越性。

Claims (6)

1.一种基于有限单元与极限平衡的土体稳定性评价方法，其特征在于，所述土体稳定性包括边坡稳定性、地基承载力和支挡结构土压力三类问题，所述方法具体包括以下步骤：

（一）地形测绘、地质勘察与有限单元模型构建

针对待评价土体，进行地形测绘、地质勘察和土体物理力学参数测定，根据各类土体的空间分布、地下水位线以及工程实践对计算精度的要求，将评价范围内的土体划分为有限个单元，建立土体稳定性评价的有限单元模型；实际划分中：

对平面二维问题，采用三角形和/或多边形单元进行划分；

对空间三维问题，采用四面体和/或多面体单元进行划分；

（二）总体约束方程构建

所述总体约束方程由屈服准则约束方程、静力平衡约束方程和应力已知边界条件约束方程共同构成；

（三）自变量

以所有节点的应力分量为自变量；

（四）目标函数

对边坡稳定性问题，目标函数是判断满足步骤（二）的总体约束方程的自变量是否存在；若存在，判定该土体为稳定的；若不存在，判定该土体为不稳定的；

对地基承载力问题，目标函数是判断满足步骤（二）的总体约束方程的自变量是否存在；若存在，判定该地基承载力满足要求；若不存在，判定该地基承载力不满足要求；

对支挡结构的土压力问题，设待求边界上以节点应力为变量表示的土压力的大小为P；对极限最小土压力，目标函数是求解满足步骤（二）的总体约束方程条件下的P的最小值；对极限最大土压力，目标函数是求解满足步骤（二）的总体约束方程条件下的P的最大值；

（五）评价与结论

由步骤（二）的总体约束方程、步骤（三）的自变量、步骤（四）的目标函数共同构成稳定性评价方法的数学模型；

采用最优化求解的方法，对上述数学模型进行计算求解；

对边坡稳定问题，是判定数学模型是否有解；若数学模型的解存在，判定该土体为稳定的；若解不存在，判定该土体为不稳定的；

对地基承载力问题，是判定数学模型是否有解；若数学模型的解存在，判定该地基承载力满足要求；若解不存在，判定该地基承载力不满足要求；

对支挡结构土压力问题，是求出土压力P的最小或最大值；P的最小值即为极限最小土压力，P的最大值即为极限最大土压力。

2.如权利要求1所述基于有限单元与极限平衡的土体稳定性评价方法，其特征在于，步骤（一）中，地形测绘应覆盖土体稳定性评价的范围；

地质勘察应包括：查明评价范围内土体的地下水位线、具有不同物理力学性质的各类土体的空间分布等，并为评价范围内各类土体物理力学参数的测定提供必要的试样试件；

需要测定的土体物理力学参数包括：土体内摩擦角、粘滞系数、天然容重、干容重、饱和容重、含水率。

3.如权利要求1所述基于有限单元与极限平衡的土体稳定性评价方法，其特征在于，步骤（二）中，对平面问题

采用Mohr-Columb屈服准则，具体屈服准则约束方程构建为：

，

公式中，设单元模型的节点编号依次为，第个节点所对应的节点应力分量分别为、和；设测定出的土体内摩擦角为，；

需要说明的是，对于每个节点，Mohr-Columb准则有4个不等式约束方程；若节点总数为，则共有个屈服准则不等式约束方程。

4.如权利要求1所述基于有限单元与极限平衡的土体稳定性评价方法，其特征在于，步骤（二）中，节点静力平衡约束方程构建过程中：对于平面问题，单个节点的静力平衡方程只包括 和 两个方向的平衡方程；对于三维空间问题，单个节点的静力平衡方程则包括、和 三个方向的平衡方程；

针对平面问题，设围绕节点的所有单元按逆时针方向编号依次为 ，围绕节点的所有节点按逆时针方向编号依次为；设节点 的坐标为 ，节点的坐标为，节点的坐标为，节点的应力为、和；设相邻节点之间的正应力及切应力均按直线变化；最终可得到对节点的静力平衡约束方程为：

，

；

其中：

为节点的正应力对节点在X方向力的贡献；

为节点的切应力对节点在X方向力的贡献；

为节点的正应力对节点在Y方向力的贡献；

为节点的切应力对节点在Y方向力的贡献；

为单元的体力对节点在X方向的总合力；

为单元的体力对节点在Y方向的总合力；

需要说明的是，针对平面问题，对每个节点有2个静力平衡等式约束方程；若单元模型的节点总数为，则共有个静力平衡等式约束方程。

5.如权利要求1所述基于有限单元与极限平衡的土体稳定性评价方法，其特征在于，步骤（二）中，应力已知边界条件约束方程中，

对边坡稳定性和支挡结构土压力问题，设定为：在应力已知边界上，令其位于应力已知边界上节点的应力分量等于其已知量；

对地基承载力问题，则需要把所需的地基承载力设定为应力已知边界，并与其它应力已知边界共同构成地基承载力问题的应力已知边界，令其位于应力已知边界上的节点的应力分量等于其已知量；

在应力已知边界上，设已知的节点正应力和切应力分别为和 ，节点处的边界方向为，则有：

；

此式即为应力已知边界条件约束方程；式中，与应力已知边界条件及单元划分有关。

6.如权利要求1所述基于有限单元与极限平衡的土体稳定性评价方法，其特征在于，该方法具体用于解决竖直支挡结构单位宽度内的水平向最小土压力问题或者堆料仓内的堆料对仓壁产生的侧压力问题。