9 非关联流动法则

前面几节中的变分原理都至关重要地依赖于关联流动的概念。对于某些材料，尤其是金属，此假设与实验相吻合，因此原理可以直接应用。另外，对于摩擦材料，相关流动的假设通常与实验相矛盾。实际上，与Mohr-Coulomb，Drucker-Prager或其他相关破坏准则相关的膨胀通常远大于实验中的实际观察值。因此，通常使用非关联流动法则，依据指定流量G来实现，它与屈服势F有相同公式形式的。

使用非关联流动法则，即使用线性弹性/完全塑性模型，也可以很好地模拟土壤的真实行为。图9.1概述了这种情况。对于松散的沙子和一般固结的粘土，在典型的三轴或双轴排水试验中，剪应力单调增加到某个极限状态。在此过程中，土壤将经历压实，其速率最终将趋于零。使用OPTUM G2中的Mohr-Coulomb模型，可以考虑初始近似线性的应力应变响应和极限强度。在中间范围内，由于理想的，如非硬化的可塑性假设，Mohr-Coulomb模型的响应往往有些僵硬。这些特性与使用关联还是非关联流动法则无关。但变形却（此处以体积应变与剪切应变的关系表示）高度依赖于流动法则。 实际上，对于当前示例，相关流动法则预测的膨胀太大。另外，非关联流规则直接指定了膨胀量，因此可以非常精确地说明变形。

图9.1：Mohr-Coulomb相关/不相关的模型下的土壤变形能力

注意：非关联响应假定没有局部位移发生。

对于稠密的沙子和超固结的粘土情况大致相同。应力应变响应独立于流动法则，相关流动法则过度预测了膨胀量。我们还注意到出现了峰值后应力应变响应减弱。这与扩张趋势减小，最终达到零值相吻合。

根据这些观察结果，很容易得出结论：屈服函数和流动法则可以独立调整-前者与观察到的强度相匹配，后者与观察到的变形行为相匹配。尽管有一些道理，但却有些天真。实际上，引入非关联流动法则会出现许多现象，而这些现象并不能从上述简单的观察中揭示出来。这是在下一节中讨论。

9.1非关联流动结论

流动法则如何响应塑性或弹塑性边值问题的是一个古老的问题，至今仍未解决。由于流动法则对变形的影响非常明显（至少在定性上），因此重点主要放在流动法则如何影响极限载荷上。根据图9.1所示的变形类型，可以得出结论是根本没有影响：极限载荷与流动法则无关。这种观点在某些圈子中仍然得到坚持。它前提是流动法则只会影响“高阶运动约束”。但是，对于大多数土力学问题来说“低”之外，从未明确阐明“高阶运动约束”是如何量化的。这样的陈述很少得到定量结果的支持，例如以有限元分析的形式来解决具有不同摩擦角和膨胀角的问题。如果是这样，那么问题和/或所选的材料参数集通常没有什么实际意义。另一方面，在进行更严格和实际相关的数值分析，结论将不可避免的是，流动法则对极限极限载荷有相当大的影响。例如，Erickson和Drescher（2002），Loukidis和Salgado（2009）和Krabbenhoft等人。 （2012年）报告说，砂土上条形和圆形地基在的承载力降低了45％。承载能力的这种下降以及在典型数值分析中观察到的许多其他可能出乎意料的现象也获得了丰富的理论支持。部分将在下面讨论。

9.1.1结论1：极限载荷的局域性和非唯一性

解决弹塑性边界值问题时通常会做出一个隐式假设，即响应是唯一的，或至少极限载荷是唯一的。如分别针对硬塑料和弹塑性材料的第6节和第8节所示，关联流动法则时，此假设确实成立。

但是，如果非关联流动法则，则响应通常是不唯一的。 换言之，可能存在多个都满足控制方程的解决方案，但每个解决方案都暗示着不同的极限载荷。早已建立了这一事实以及用于分析关于唯一性的本构模型的工具（Rice 1976），但是在地质力学界却发现相对较少的共鸣。

图9.2：双轴测试易于位移的响应

非关联材料的非唯一性与剪切带的出现密切相关，即剪切带的局部变形带往往会突然形成，并通常导致载荷-位移曲线出现明显的软化。该现象可以描述如下。考虑如图9.2所示的双轴测试。通过使用非关联流动法则的本构模型对测试进行建模。最初，响应将是均匀的，即应力和应变状态在整个样本中都是恒定的。然后在某些时候，可能会出现分叉。响应可能沿着变形保持均匀的路径继续，可能会因样品形成剪切带而从根本上发生变化。变形模式的这种改变将导致轴向力的下降，随着继续负荷，轴向力的持续下降。因此，响应是非唯一的：有两个解决方案，均质解决方案和局部解决方案，都满足所有控制方程。

此外，假设样品在第一个分叉点之后仍保持均匀。这并不妨碍在加载过程中在其他某个位置发生分叉。 并且局部解的几何形状（此处由剪切带倾斜度给出）不一定与第一个分叉点相同。 实际上，如果样品在第一个分叉点之后仍保持均质，则通常会有一系列可能的局部极值。

可以证明，本构模型导致局部极值的趋势度量，可以通过声张量的特征值定义：

(9.1)

其中Dep是弹塑性本构模量，P是通常的投影矩阵（1.8），n是剪切带的法向（见图9.2）。 因此，对于给定的应力状态σ和假定的剪切带向量n，如果A的行列式小于或等于零，则可以进入局部极值（localization ）。如果不是，则样品稳定将继续沿着均匀路径。 Bigoni（2000）详细讨论了这一结果和类似结果。Bigoni (2000); Bigoni and Hueckel (1991); Leroy and Ortiz (1989); Rice (1976); Runesson et al. (1992) 等。

9.1.2结论2：减少极限载荷

非关联流的第二个结论是，通常在数值实验中观察到的极限载荷明显低于在关联流假设下获得的相应的极限载荷。通常参照Mohr-Coulomb模型，极限载荷的减小取决于摩擦角和膨胀角的角度的大小及之间的差值。因此，随着摩擦角与膨胀角差值的加大，极限载荷的减小也加大。当差值一定时，减小量通常随摩擦角和膨胀角绝对值的增大而增大。换言之，参数组（ϕ，ψ）=（40◦，10◦）通常比参数组（ϕ，ψ）=（30◦，0◦）导致轴向承载力的降低幅度更大。

这种降低的原因有两个。首先，在大多数情况下，非关联流动法则引起局部极值的趋势将意味着极限状态将由应力和应变的高度局部化状态控制。其次，流动法则对剪切带的运动学施加了一定的限制，在非关联情况下会导致承载能力的降低。第二个效果可以通过图9.3所示的简单示例说明。这里考虑刚性块处于刚性摩擦表面。滑块承受恒定的法向力P，分析可以承受的最大切向力Q。关于问题的确切几何形状，我们想象块和表面之间的任意薄界面，其中相关屈服条件是库仑的：

(9.2)

其中ϕ是摩擦角。

另外，流量由下式给出

(9.3)

假定为非扩张（non-dilative）行为。

图9.3：刚性摩擦面上的刚性块

因此，塑性应变率由下式给出：

(9.4)

此外，应变-位移关系如下

(9.5)

其中是界面上的速度突变，即不连续的上下界面之间的差值。 类似地，是沿界面的速度突变，即界面两端之间的之差。令δ→0，并且需要最大应变率的有限大小，我们发现=0。换言之，无限薄的界面不会纵向拉伸。 

回到流动法则，施加以下约束：

σx = σy    (9.6)

将此约束代入屈服条件，得出：

F = 2|τxy | + 2σy sin ϕ = 0    (9.7)

从中我们得出解决方案：

Q = P sin ϕ      (9.8)

它不是P = Q tan ϕ的标准解，确是假设关联流动G = F的解。

通过tan ϕ ∗ = sin ϕ引入“有效”摩擦角ϕ ∗，可以将解写为

Q = P tan ϕ∗    (9.9)

对于ϕ = 40°，我们有tan ϕ / tan ϕ ∗ = 1.3，这与许多实践准则（例如，欧洲规范7）所规定的摩擦角标准值与设计值之比相似。因此，非关联流的影响也几乎不能忽略，即使对于最少的“运动限制”的问题。

可以将上述过程推广到经历局部变形的任意三维实体（Krabbenhoft等人，2012年，2004年，2010年）。 结果是，有效强度结果由与原始屈服函数相同的函数定义，但涉及到材料参数的减少，其大小取决于摩擦角和扩张角。因此，对于Mohr Coulomb准则，我们有一个有效强度范围，由

(9.10)

图9.4：根据（9.11）的非关联性影响：摩擦角（左）和摩擦系数（右）

有效强度参数由下式给出

   (9.11)

这些关系首先由Hill（1950）针对特殊情况ψ= 0推导，后由Davis（1968）针对一般情况推导。假设摩擦角和膨胀角之间的依赖性不同，则非关联度对有效莫尔-库仑摩擦角的影响如图9.4所示。 同样，对于Drucker-Prager准则，有效强度范围由（Krabbenhoft等人，2012）给出：

(9.12)

其中

  (9.13)

9.1.3结论3：“数值问题”

例如双轴测试的有限元分析中，通常会观察到某种破坏导致局部极值。这可能以材料轻微不均匀性或几何形状或边界条件的轻微破坏的形式出现。此外，经常观察到影响往往取决于有限元网格。如果使用规则的网格，则存在沿着网格边缘发生极限的趋势。考虑到局部破坏带方向的不唯一性，可以预料由此引起的整体力-位移响应。

图9.5：ϕ = 40°和ψ= 10°的无重度土壤上条形基础的荷载-位移响应[根据Krabbenhoft等（2012）]

其次，更重要的是，经常有报道说，涉及非关联本构模型的边值问题的数值解比与关联流动法则的情况要难得多（Carter et al. 2005; Clausen and Krabbenhoft 2008; Loukidis and Salgado 2009; Manzari and Nour 2000）。 对于较高的（但真实）摩擦角和非关联性，这些问题更加明显。同样，材料参数恒定，通常会观察到性能随着模型中有限元数量的增加而下降。

尽管经常被报告为“数值问题”，但这些特征并非完全出乎意料。因为可能存在一系列可能的解决方案，每个解决方案都与不同的局部破坏模式相关联，并且全部都是完全有效的，可以预料，任何数值解都将对物理缺陷以及舍入误差以及定义解方案的过程的确切顺序都非常敏感。最后，结果是负载位移响应趋于振荡。图9.5给出了一个示例，该示例涉及在粘性摩擦和失重的土壤上的条形基础分析（所谓的Nc问题）。与相应的相关问题的响应相比，我们注意到以下几个特征：i）关联问题的刚性更高，除了在响应的最初阶段是纯弹性的； ii）非关联导致极限强度的显着降低（在这种情况下，几乎是40％），iii）对于非关联的问题，超出一定位移水平的负载-位移响应有些振荡。这些振荡是因为边界值问题的不唯一，在物理上对应于超出结构的承载能力首先耗尽的点。这些点在破坏点不同破坏模式之间的切换。如图9.6所示，可以看出，在不同的载荷步之间，失效模式的变化非常明显，超过了明显达到极限载荷的点（在这种情况下，位移为）。

图9.6：不同载荷阶段的塑性剪切应变率（根据Krabbenhoft等2012）

在这种情况下，基于标准牛顿-拉夫森的程序可能会失效。确实，此类程序依靠的关键假设是结果可估算。由于估算未知解的唯一真实可能性是通过最后的收敛解，然而由于应力和变形的模式可能会显著突变，因此不可能收敛。该功能和荷载-位移曲线的波动通常被组合在一起，并标记为“数值问题”。但是，前者肯定是由于求解方法的破坏引起的，但后者是由于本构模型，尽管可能有些违反直觉，但这是可以预期的。

9.2变分公式

出于传统的基于牛顿-拉夫森的解决方案的缺点的考虑，Krabbenhoft等人。 （2012年）开发了一种新的基于优化的解决方案，专门针对非弹性弹塑性。 下面描述该方案的基本原理。

9.2.1摩擦和可塑性

公认的是，除非G = F，也就除非关联流动法则，否则第2节中所述类型的本构模型不允许采用变分公式。在这种情况下，本构方程可以用类似于von Mises的最大塑性耗散原理的变分方程来表达。此外，如果屈服函数是凸的，则可以根据凸函数的数学性质来转化控制方程。这样允许直接分析与存在和唯一性有关的属性。一般的非关联情况不可以这样表述。尽管此缺憾不会给开发类似于相关情况的基于牛顿-拉夫森的传统求解方法带来任何基本障碍，但仍失去了与关联塑性理想数学特性。此外，尽管关联的塑性涉及对称的切线模量，但不规则的流动法则通常会产生一组不对称的离散有限元方程组。最后，尽管可以使用现代数学程序设计的方法非常有效地解决关联塑性问题，但在这种公式化不适用非关联情况。

由于关联可塑性的数值算法的相对效率和强大功能，下面给出了一种数值公式，该数值公式保留了相关联的可计算性的理想属性，但适用于一般非关联模型。 但是，应该注意的是，在局部变形方面，非关联塑性的特征得以保留。这包括在涉及完全塑性非关联模型的边值问题中经常观察到的明显的整体软化。同样，大多数非关联模型所隐含的解的非唯一性也持续存在，并表现出对有限元网格，边界条件等的强烈敏感性。

9.2.2摩擦微观力学

公式背后的基本思想源自与第2节中所述类型的本构模型相关的内部耗散结构。假设该类型的屈服函数

（9.14）

图9.7：微观上的摩擦起因是粗糙物的塑性剪切，限制压力越大意味着粗糙互锁程度越强，从而具有较高的表面剪切强度。

其中p和q分别是偏应力和平均应力的某种量度，M是摩擦系数，k是内聚力。 此外，让塑料势能由

(9.15)

其中N≤M是膨胀系数。

在p-q空间中，塑性应变率由下式给出：

(9.16)

其中和是分别与p和q共轭的体积应变和偏应变。与模型（9.16）相关的内部耗散由下式给出：

(9.17)

这种内部耗散的表达揭示了几个有趣且众所周知的功能。首先，对于关联材料（N = M），耗散与内部内聚力k成正比。对于纯摩擦材料（k = 0）不会发生内部耗散，这与实验显示明显不同。其次，对于N <M，耗散与表观凝聚力成正比，该表观凝聚力包括两项：内部凝聚力k和因非关联度产生的-（M- N）p。后者作为“表观内聚力”支持摩擦是由接触的固体表面上的微观粗糙物的相互机械作用产生的（Bowden and Tabor 1973）。由于粗糙表面上的应力远大于材料的弹性极限，因此，微观上的塑性变形主要决定了宏观观察到的摩擦阻力。图9.7说明了这一点，该图显示了处于不同限制压力下的两个粗糙表面。可以认为塑性剪切是韧性的，纯内聚。对于诸如沙粒之类的脆性材料，该假设是由于在粗糙程度下具有很高的应力水平而成立的，这有效地使原本脆性的材料具有延展性。因此，每个部分的表观剪切强度仅源自压力引起的几何变形，并且库仑摩擦暗示了表观剪切强度与压力之间的线性关系。这种解释促使将屈服函数（9.14）更改为

(9.18)

其中

(9.19)

内聚力明显与压力有关。该材料参数体现了摩擦界面的实际微观力学的所有复杂性，以及它们对所施加载荷的响应演变。因此它的相对简单性是令人惊讶的，但是仍有很多需要额外注意的，例如Bowden和Tabor（1973）所讨论的，但它仍然适用材料非常广泛。

9.2.3时间离散公式

现在假设表观凝聚力是已知的。相关流动法则产生所需的结果，即塑性应变率（9.16）。 在解决边值问题时，表观内聚力当然不是先验的，因为它与解决方案中确定的压力成正比。但是，假设通过一系列虚拟时间步长逐步解决了这些问题，屈服函数的某些部分上原则可以隐式评估，而有些部分明显式评估。假设在时间tn的状态是已知的。然后，在tn + 1处施加的屈服条件可以近似为

(9.20)

其中 

(9.21)

同样，相关流动法则产生所需的时间离散结果，即：

(9.22)

但对表观凝聚力的显式计算意味着对于tn + 1时的新应力状态可能会超出原始屈服函数。同样，近似值可能暗示应力状态下的塑性屈服，否则将被认为是纯弹性的（见图9.8）。 但对于足够小的增量，即tn + 1-tn→0，可认为对表观凝聚力的显式评估所引入的误差趋于消失。这种假设可通过数字实验证实（Krabbenhoft等，2012）。

因此，任何非关联塑性问题都可以近似等效为相关塑性问题。每个时间步都涉及解决一个近似的相关问题。这样，先前开发的用于相关可塑性的方法几乎不需要修改即可适用。此外，后来可以使用一般的或更专业的数学编程方法对变分方程进行求解。 OPTUM G2中采用了此类方法。

图9.8：表观凝聚力的显式评估：原始屈服函数和近似屈服函数

根据原始屈服函数，（a）表示的面积是不允许的，但根据近似屈服函数，则是允许的。 类似地，根据近似屈服函数，由（b）表示的区域是不允许的，但是根据原始屈服函数，该区域在弹性域内。

9.2.4 Mohr-Coulomb

基于上述考虑，平面应变Mohr-Coulomb屈服函数实际计算采用以下形式：

（9.23）

其中

(9.24)

下标0表示当前的已知状态。换言之，仅修改了常规圆锥，而拉伸和压缩截留值保持不变。

9.2.5 Drucker-Prager 

同样，Drucker-Prager屈服函数由下式给出

(9.25)

其中 

(9.26)

下标0表示当前的已知状态。同样，仅修改了常规圆锥，而拉伸和压缩截留值保持不变。

9.3极限分析的状态

经典极限分析主要依赖于关联流动的概念。确实，如第6节所述，上限和下限原理仅在与流动法则相关的假设下有效 对于非关联流动，唯一可用的结果是ψ = ϕ材料的上限也是ψ<ϕ的相应材料的上限。换言之，根据相关流动法则计算出的极限载荷是不安全的。不安全程度取决于不关联的程度（ϕ和ψ之差）和特定问题。如第9.1节所述，可以通过有效材料参数ψ∗和c ∗（对于Mohr-Coulomb）或M ∗和k ∗（对于Drucker-Prager）定义得出“有效强度域”。

长期以来，人们一直建议将这些有效参数用于标准的关联塑性公式中，以计算“非关联极限载荷”（例如Drescher和Detournay 1993）。对于9.1节的简单示例，有效的材料参数肯定会正确的计算结果。对于一般问题（例如通过有限元方法解决的问题），在相关公式中使用有效的材料参数还可以通过正确的摩擦角和膨胀角，通过完整弹塑性计算，得出最终极限载荷大小的合理估算值。不过，应注意引入非关联流动法则所隐含的非唯一性，即有限元分析得出的确切值将高度依赖于几何形状、边界条件、破坏以及为预期局部变形能力来划分的网格优劣。在下面的示例中对此进行了说明。我们考虑一个内聚摩擦材料的棱柱形样品，如图9.9所示。材质模型为Mohr-Coulomb，其弹性参数E和ν以及内聚力c和摩擦角ϕ都大于零。假定是光滑的顶板逐渐向下移动，而垂直侧的法向应力和剪应力保持为零。不管流动法则即ψ的值如何，达到屈服点（A）为止的载荷-位移曲线将如图9.10所示。此时，垂直应力为：

(9.27)

对应于均质变形的连续状态，也与流动法则无关（尽管均质塑性应变场取决于流动法则）。对于关联流动法则，这也是精确（且唯一）的极限荷载。

图9.9：单轴压缩测试设置（左）和可能的局部破坏模式（右）

图9.10：单轴压缩试验的载荷-位移示意图

图9.11：在图9.10中的点A处，相关和非关联流动法则的声张量的行列式（通过弹性模量的行列式标准化）作为剪切带破裂角的函数

为了深入了解可能的局部变形模式，检查了声张量在点A处的应力状态。结果如图9.11所示，其中ψ= ϕ = 40°的相关流动法则和非关联流动法则 ϕ = 40°和ψ= 10°。如图所示，相关流动法则暗示了剪切带破裂角θ= 45°+1/2ϕ的可能性，该剪切带对应于导致破坏载荷著名的滑移线解（9.27）。

另一方面，非关联流动法则允许沿着一系列方向进行局部破坏可能性。可以看出，该范围是45°+1 /2ψ≤θ≤45°+ 1 /2ϕ。 这里所讨论的角度，意味着由剪切带倾斜度在50°≤θ≤65°范围内定义的任何局部解都是可能的。根据第9.1.2节的分析，可以证明相应的破坏荷载由下式给出：

(9.28)

图9.12：根据（9.28）非关联材料的极限载荷为剪切带倾斜角的函数。关联的极限载荷由虚线表示

其中c ∗和ϕ ∗由（9.11）给出。 θ= 45°+ 1/2ϕ ∗时达到最小值，推出

(9.29)

对整个的范围，非关联极限荷载（9.28）小于关联极限荷载（9.27），而当θ= 45°+ 1 /2ψ，两者重合（见图9.12）。这在图9.10中也有说明，其中应变大小δ原则上为零，但在有限元计算中，将取决于剪切带的厚度，该剪切带的厚度受到网格大小的限制。

有趣的是，图9.12还表明，在可能的剪切带倾角范围之外，即θ> 45°+ 1 / 2ϕ时，非关联情况下的极限荷载低于关联情况下的极限荷载。但是，控制方程式不允许在该范围内从均匀变形模式到局部变形模式切换。

有限元计算中，首先通过减小极限载荷，其次通过对破坏的高度敏感性来观察上述分析结果。 图9.13-9.14中显示示例。图9.13显示了ϕ = 40°的单轴压缩试验的各种模拟的载荷-位移曲线。正如预期的那样，相关的计算导致分析极限载荷（9.27）。 接下来，将扩张角设置为ψ= 10°，并且有限元网格的坐标设置为随机扰动样本宽度的10 ^-3倍。生成一系列的载荷-位移曲线，尽管遵循预期的强度降低趋势，但在峰值之后的演变以及最终的残余载荷方面都有些不同。这说明了上面讨论的非唯一性。

关于极限载荷的大小，可以看出它略高于理论上预测的最小值（9.29）。 这是由于网格无法捕获由θ= 62.7°的剪切带倾角定义的最小解。 实际上，图（9.13）所示的局部变形模式对应于，因此，根据（9.28），极限载荷为

图9.13：单轴压缩试验的载荷-位移曲线

图9.14：有限元网格（六节点单元）和u = 1.0时单轴压缩试验的不同局部变形模式（见图9.13）

在这个例子中观察到的趋势完全代表了具有实际意义的典型问（Krabbenhoft等，2012）。 特别是在极限分析的适用性方面，使用原始c和ϕ会导致高估承载力，而使用有效参数c ∗和ϕ ∗则会导致低估承载力。尽管原则上非关联可塑性的本质使得极限载荷毫无意义，但有效参数通常可以为使用正确的c，ϕ和ψ进行的完整弹塑性分析结果提供合理的估计。

9.4不排水条件

尽管在排水条件下的膨胀角可能会对承载力产生某些影响，但在不排水条件下的膨胀角效果显著得多。 实际上，在不排水条件下，必须使Mohr-Coulomb，Drucker Prager和类似的模型，且其扩张角或系数等于零才能获得合理的结果。 其实任何有限的膨胀都会导致无限的承载能力。该结果有时被表示为数值伪像，但实际上完全符合这些模型，因此没有任何怀疑或神秘之处。

要弄明白这一点，请考虑增量中超静孔隙压力与体积应变之间的关系。

图9.15：由于塑性膨胀而产生超静孔压

其中

(9.30)

显然，体积应变∆εv> 0的增加导致超静孔隙水压力∆pe> 0的增加。当总应力恒定，这反过来又导致有效平均应力的减小：

 (9.31)

对于摩擦屈服准则，实际上意味着一旦屈服开始并产生膨胀，材料就会发生明显的“硬化”（见图9.15）。结果是永远不会达到极限状态。因此，必须使用零体积应变的非关联流动法则。所以，对于Mohr-Coulomb模型，必须使用ψ= 0的非关联流动法则。类似地，对于Drucker-Prager模型，必须使用N = 0。