第1章 绪论
稳定性(或强度)是土力学的两大主题(变形和强度)之一。根据传统土力学知识，稳定性问题的主要内容包括地基极限承载力、边坡稳定和侧向土压力等。但是，随着人类建设活动的逐步深入，稳定性问题出现许多新的内容，如隧道开挖面的稳定性以及块体材料的承载力等。目前，绝大多数岩土工程都涉及稳定性问题，而且稳定是保证结构安全的前提，因此对其分析理论和方法的研究具有重要的学术意义和工程实用价值。
对一般弹塑性材料而言，受荷作用直至破坏都要经过弹性、约束塑流和自由塑流三个阶段，而要得到比较完整的解，其计算过程都相当复杂，而且不便于实用。在土力学中，许多问题都是以塑性塌陷引起的破坏作为控制条件的。因此，发展直接而有效的方法来计算破坏荷载，其实际意义是不言而喻的。
破坏荷载可以通过基于塑性理论的各种方法来确定，主要有极限平衡法、滑移线法、极限分析法及弹塑性位移型有限元法等，这些方法假设材料为理想(刚)塑性体，尽管理论上并不严格，但也得到了广泛的研究和应用。极限平衡法在岩土工程中应用*为广泛，尤其在边坡稳定性分析方面。这种方法把稳定性问题简化为静力学问题，没有考虑材料的应力-应变关系，也没有考虑材料的流动法则，所以它是一种近似的解法。极限平衡法应用广泛的主要原因在于其处理不同加载条件、复杂的几何边界和土层等问题时相对简单；另外，各种极限平衡法商业化程序的应用，使得该方法更具吸引力。然而，分析中采用大量的假设条件，导致极限平衡法的精度很难量化。
特征线法，又称为滑移线法，是根据边界条件对土体极限平衡微分方程进行求解，从而导出极限平衡区的滑移线场和应力分布，计算得到的基础范围内的边界应力即为极限承载力(Michalowski and Shi，1993；Hansen and Christiansen，1969；Hansen，1961)。如果滑移线场区域的应力分布能够扩展到全部区域，而且满足静力场的要求，则滑移线法解就是下限解；若能建立与滑移线场对应的速度场，而且满足运动场的要求，则滑移线法解就是上限解；若同时满足上述两种要求，则滑移线法解就是极限荷载的精确解。
滑移线法的数值解和图表解在实际问题中的应用都局限于相当简单的情况。主要问题在于扩展的塑性应力场很难同时满足平衡方程、屈服准则和边界条件。大多数情况下，速度特征线往往不能确定，这样相应的应力特征线本身就不能代表完整解。另外，正如极限平衡法一样，滑移线法也没有考虑材料的本构关系，仅部分区域满足平衡方程。更为重要的是，滑移线法缺乏通用性，对于复杂的几何边界或加载情况问题就无能为力。
随着极限分析理论的发展，尤其是Drucker等(1952)证明了上、下限定理之后，各种稳定性问题都可以严格地通过上、下限定理来界定破坏荷载。其中，下限定理阐述了通过静力容许应力场计算得到的荷载总是小于破坏荷载的真值，静力容许应力场在连续体内满足平衡方程、应力边界条件和屈服准则；而上限定理阐述了通过机动容许速度场计算得到的作用荷载和体力不小于破坏荷载真值，机动容许速度场必须满足相容性、速度边界条件和相关联流动法则。
在稳定性分析问题中，极限定理的有效性是显而易见的。因为稳定性问题的真解采用其他方法是难以确定的，而极限分析的上、下限定理则可以界定真值。根据极限分析上、下限定理计算得到的结果可以作为其他非严格方法计算结果的基准尺度(benchmark)：①解法严格；②从应用角度而言，真值能界定在相当小且准确的范围内；③破坏形式不需要事先假定。但受计算条件和方法的限制，传统意义上的极限分析并未充分发挥作用。
随着计算机技术的发展，弹塑性位移型有限元方法在破坏荷载问题的确定中得到了普遍应用。与极限平衡法和滑移线法不同的是，基于弹塑性理论的有限元法可以提供应力和应变分布的计算结果。在有限元分析中，可以模拟与理想塑性材料不同特性(如应变软化和硬化)、复杂加载条件和几何形状、加载历史以及孔隙水压力的生成和消散等情况。每个有限单元节点的应力和应变对应给定的本构关系，并满足平衡方程、相容性和边界条件。尽管有限元法优势明显，但是也只能提供近似解。另外，有限元法往往需要很大的计算成本和很长的时间，而且复杂的本构关系使得该方法在面向实际运用时存在较大的困难。
尽管极限分析法中的上、下限定理可以严格地界定破坏荷载的真值，是岩土工程稳定性问题的重要计算手段，但是采用人工构造应力场和速度场来获得优化解是相当困难的，再应用于实际工程问题几乎是不可能实现的，这就是极富潜力的极限分析法未得到充分利用的主要原因之一。在进行极限分析时，目标就是要寻求优化解——*大的下限解和*小的上限解。随着计算机技术和数值计算能力，尤其是数学规划理论的快速发展，数值极限分析在工程问题中得到应用成为可能。另外，对于受荷作用的结构，进行设计时总是考虑留有一定的安全余量(safety margin)，因而合理而有效地计算破坏荷载在安全性、经济性和时间方面都可以节省大量成本。
基于塑性力学上限定理的极限分析方法理论基础严密，其得出的破坏模式有助于掌握结构体系的力学性能，因而其数值理论研究和应用具有较强的现实意义。
综上所述，本书主要对数值极限分析上限法用于求解破坏荷载时存在的问题及算法进行研究，并就其用于更多的实际工程问题的可行性和有效性进行研究探讨。下面详细回顾一下数值极限分析方法的理论发展及研究现状。
1.1 数值极限分析方法的理论发展及研究现状
从数学角度来看，极限分析的上、下限定理是互为对偶的命题，其示意图如图1.1所示。因为上限定理和下限定理之间存在必然的联系，因而在介绍数值极限分析上限法的同时，不可避免地会涉及数值极限分析下限法的发展。
图1.1 上限定理与下限定理的相互关系
极限分析的计算起源于结构分析领域，并且是由下限法开始的。Kazinczy和Kist通过对连续梁的极限承载力进行研究，提出只要满足平衡方程且梁的任何截面都没有出现达到使梁破坏的弯矩时，结构便不至于破坏，这实际就是下限定理的雏形。1940年，Van从实际工程的角度强调了结构极限设计的实际意义，并且出版了**本关于极限设计的专著。Drucker和Prager(1952)*早将极限分析原理应用于求解土工问题；而Chen(1975)的专著Limit Analysis and Soil Plasticity则完全奠定了极限分析方法成为岩土工程领域有效计算手段的理论基础，并达到极限分析解析解法的一个巅峰。
数值极限分析法可以追溯到20世纪60年代中期。Koopman和Lance(1965)*先利用下限定理和线性规划技术研究了满足线性化屈服准则条件下方形板的承载力问题，发现静力容许应力场只需在网格的离散点处而非整个物体内得到满足。Ceradini和Gavarini(1968a，1968b)、Laudiero(1972)以及Faccioli和Vitiello(1973)在板的稳定性研究中也采用了下限定理与线性规划相结合的方法。
数值极限分析上限法经历了半个多世纪的发展演变，至今仍是理论研究的热点问题之一。下面分别从数值极限分析上限法的单元离散及变量形式、算法研究、数值极限分析的其他方法及应用等方面进行回顾和总结。
1.1.1 数值极限分析上限法的单元离散和变量形式
数值极限分析上限法的主要特点是：通过检验位移协调条件获得优化问题的解，进而构造机动容许速度场。Hodge和Belytschko(1968)*早采用该思想处理二维平面问题和板问题。以Airy应力函数为基本变量，通过插值，将应力函数以单元节点处的变量来表示；单元边界的连接条件通过消除自由度来实现，单元间可以出现应力间断，并采用非线性规划算法进行*小化求解。Hutula (1976)采用相同的思路，推导出数值极限分析上、下限法统一的有限元格式，然后借助线性规划方法得出问题的真正解答。
Casciaro等(1971)选择二次应力形式进行单元离散，发现需要花费较大的成本来检验每个单元的屈服条件。*后建议采用线性离散，这样仅需在单元节点处检验屈服条件；采用分段线性化屈服准则可以有效地借助线性规划进行求解。Christiansen(1981)采用Casciaro等(1971)的线性离散方式来避免较大的优化成本，以实现极限问题的求解。Turgeman和Pastor(1982)在轴对称问题中也采用了这种方法，其中，Mohr-Coulomb条件采用线段线性化近似，锥体由多面体代替。这一思想仍被研究人员采用，如Sloan(1989)、Yu等(1998)、Ukritchon等(2003a，1998)。
Krabbenhoft和Damkilde(2002)采用高级单元的三角形板形式，利用扭转弯矩来构造速度非连续面，其构造形式非常新颖；但是，仅在板的受荷特性中应用，在其研究板桩侧向受荷特性时也使用该方法的上限形式(Krabbenhoft et al.，2005a)。
由于无厚度速度间断面(非连续面)的引入，上限计算结果受非连续面的分布影响较大。因此，部分学者探讨了线性单元的上限有限元网格自适应技术，如Borges等(2001)、Christiansen和Pedersen(2001)、Lyamin和Sloan(2003)、Franco等(2003)、Lyamin等(2005)、Ciria等(2008)。
在数值极限分析上限法的研究过程中，Bottero等(1980)采用与Lysmer(1970)下限法网格相同的形式，并采用笛卡儿坐标系统下的两个速度分量作为节点变量以及大量的未知塑性算子作为基本未知量。为模拟不可压缩条件下的不排水情况(Nagtegaal et al.，1974)，有限元网格划分所使用的三节点常应变三角形单元必须以四边形和其两对角线形式(混合形式)生成。由于上限定理和下限定理的相互关系，这一要求逐渐成为下限有限元网格的主流形式。
由于单元速度分量线性变化，要满足上限定理条件，仅需在节点处确定相容约束。相似地，在三维极限分析上限有限元方法中，可以采用四节点四面体的简单方式对速度场进行离散，以寻求严格的上限解。尽管四节点四面体阶数比较低，但采用与二维情况同类型的分段式线性应力仍有较大优势：一方面，计算精度能满足工程需要；另一方面，同样只需在节点处确定相容约束。这样，程序的通用性更强，对二维问题和三维问题都能有效地解决。Lyamin和Sloan(2002a，2002b)、Yang等(2003)、Salgado等(2004)和Ganneau等(2007)在这方面进行了详细的验证和研究工作。
由上述研究工作可以看出，数值极限分析上限法应采用有限元形式来实现，节点速度矢量作为未知变量，相邻单元通过速度非连续面(velocity discontinuity)分离；速度非连续面增加了额外的自由度，远胜于完全连续的速度场。有限元法与数学规划相结合，将极限分析问题作为线性规划问题来研究，不但使问题求解更加方便，而且进一步推动了塑性极限分析这门学科的发展。
如今，有限元法已成为学术研究和工业应用中使用*广泛的工具。通过文献分析发现，有限元模型有三种基本类型，即位移形式、平衡形式和混合形式。在极限分析的情况下，Hodge和Belytschko(1968)、Nguyen-Dang(1976)、Krabbenhoft和Damkilde(2002)、Lyamin和Sloan(2002a)、Le(2010)对平衡形式进行了研究。位移形式可以在Hodge和Belytschko(1968)、Anderheggen(1977)、Capsoni和Corradi(1997)、Krabbenhoft等(2005a)、Le等(2010c)、Bleyer和Buhan(2013)的相关文献中找到。混合形式允许直接确定应力和位移，并且可以避免体积锁定
展开
