土石混合体边坡数值模型的自动生成技术

油新华1 王渭明2 李晓3

（1.北京城建集团有限责任公司 北京 100044 2.山东科技大学土木工程学院 山东 泰安 271019 3.中国科学院地质与地球物理研究所 北京 100029）

摘要 首先介绍了土石混合体边坡的概念，并在分析其特点的基础上提出了用于稳定性分析和计算的土石混合体边坡的数值模型，详细研究了土石混合体边坡数值模型的自动生成技术。

关键词 土石混合体 数值模型 自动生成技术

1 土石混合体边坡的概念

土石混合体边坡在我国广泛分布，如香港、广东、福建等地大面积分布的花岗岩残积土边坡，西南地区、长江三峡库区及黄河中上游广泛分布的古滑坡、崩积层边坡，西藏等地分布的冰碛土边坡等［1］。该类边坡在我国各类工程建设，特别是在水利水电工程、交通工程、铁路工程中经常遇到。由于该类滑坡通常规模较大、影响因素众多、突发性强及滑移条件复杂，常给国民经济建设、人民正常生活与生命安全带来严重危害和巨大的财产损失。因此，对这类边坡进行系统分析、研究和总结具有特殊的理论和实用价值。

土石混合体边坡是由土夹碎石或碎块石以及碎石或碎块石夹土等土石混合体组成的边坡。它一般发育在第四系松散堆积层中，主要是由滑坡堆积、残坡堆积、崩坡堆积、冲洪堆积等形成。

对于土体边坡以及岩体边坡来说，人们已经进行了长期、系统、卓有成效的研究工作，并取得了丰富的实践经验。而对于土石混合体边坡来说，由于物质组成和结构的不同而与其他边坡有着显著不同的变形特点，目前人们对于它的研究还只是处于一种定性分析的阶段，或是通过它与某种因素的相关性分析来探讨滑坡的机理［2］，或是通过模型实验来进行稳定性分析［3］，而缺少像研究土坡和岩体边坡那样的理论和技术方法。

边坡稳定性分析一般采取地质分析与计算分析相结合的方法，而计算分析又常采用极限平衡法和数值计算方法。传统的极限平衡法在理论和实用方面均存在着明显的问题［4］。

其一，计算假定不够严密。常规极限平衡法中的条分法求解是一个高次超静定问题，为简化求解而做了很多不合理的近似假定。例如将滑坡看作刚体滑动，假定滑裂面处于极限平衡状态的法向应力与剪应力由条块自重确定，条块间的内力由人为假定合力传递，致使坡脚处剪应力最小，这与坡体实际应力状态与坡脚应力现象不符。

其二，滑坡机理不能解释。由于极限平衡法不考虑边坡实际应力-应变本构关系，不反映边坡失稳的累进破坏及其应力重分布规律，因而无法解释滑坡的破坏机理和形成过程。

数值分析法弥补了极限平衡法的缺陷，较合理地解释了边坡失稳过程，为边坡应力分析和稳定分析提供了较好的数值解。

土石混合体边坡的物质组成是非均质的，岩土体力学性质是非线性的，边界条件是不规则的。过去按照一般的极限平衡法把整个滑体当作刚体，从而给出整个边坡稳定性系数的做法，在理论上有一定的不足。此时数值方法就是人们首选的方法。数值方法必须要有用于数值计算的计算模型，计算模型的好坏决定了计算的成功与否。本文在分析土石混合体边坡特点的基础上，详细研究了土石混合体边坡数值模型的自动生成技术。

2 土石混合体边坡数值模型的自动生成技术

对于土石混合体边坡来说，其物质组成极其复杂且无规律，采用常规的FLAC方法建模，将是十分复杂或者是不可能的。这时如果能采用一套基于ANSYS软件的数值模型自动生成技术，来代替使用FLAC3D进行数值模型的建立和网格的划分，实现FLAC3D建模的直观化、形象化、自动化，将给人们的工作带来很大的方便，并且对于复杂的系统，能大大节省建模所需要的时间。本节通过长江三峡白衣庵滑坡某一剖面的建模过程，详细介绍了这一方法的原理，结果显示这一方法形象、直观、快速、方便，在岩土工程数值计算方面有着很好的应用前景。

2.1 技术路线

在这套模型自动生成技术中，主要使用了三套软件和相应的接口程序：AutoCAD、ANSYS和FLAC3D，前两个主要利用其强大的绘图功能和划分网格的功能，ANSYS 与FLAC3D之间的接口程序（AFP）用Fortran语言编写。整个建模的技术路线如图1所示。

图1 模型自动生成技术流程图

2.2 模型自动生成的过程

2.2.1 工程地质模型的建立

根据所提供的工程地质模型，考虑数值计算所需要的简化，建立简化后的CAD图形，或对已有的CAD图形进行编辑、修改处理，使要赋予不同力学参数的部分或不规则开挖的部分形成一个闭和的图形块，并输出以DXF（Data eXchange File）为扩展名的文件。生成DXF格式的文件，是为了将这种文件转换为IGES（Initial Graphics Exchange Specification）格式。IGES是一种被普遍接受的中间标准格式，可以用来在不同的CAD和CAE系统之间交换几何模型。由于AutoCAD系列软件中只有R12可以将DXF文件转化成IGES文件，所以将生成的DXF 文件调入R12，使其输出IGES 文件。此过程也可以通过ALGOR软件中SUPERDRAW来实现。

在建立模型的CAD图形时，尽量少采用多义曲线等弧线，以免生成的DXF文件太大，可以用多个直线段来代替曲线；但对于直接在ANSYS所提供的CAD系统中操作时，这一步可以省略，并且不必考虑用多个直线段来代替曲线。

2.2.2 网格划分

将IGES格式文件调入ANSYS，由不同的图形块生成不同的面，然后进行网格的划分。在调入时，尽量不要使用缺省的格式，以免线与线之间不闭和而不能生成面（否则，就得进行拓扑修改）。网格的划分可以采用自动划分的形式，让系统自动判断所需划分的网格数；也可以手工进行操作，根据模型中不同网格密度的需要，进行控制线段所分数目或者网格尺寸的限制。网格生成后，需将面元拉伸成体积元，体元生成过程中要注意以下两点：

①不同的面元要赋予不同的参数加以区分；

②不同的面元要分别拉伸成体元，并控制轴向网格数。最后，再进行节点和单元信息的输出，生成节点信息文件NODE.DAT和单元信息文件ELE.DAT。在做这一步处理时，不能先将面拉成体，然后划分网格，这样会由于图形的不规则，不能形成映射网格，而生成四面体网格，但四面体网格是不被FLAC3D认可的。除非对于很规则的体积元才可以进行映射网格的划分，并生成柱体网格，类似于FLAC3D中的BRICK或WEDGE。

2.2.3 FLAD3D模型数据文件的生成

要将ANSYS所生成的节点和单元信息文件能为FLAC3D利用，中间必须有一个接口程序，以实现两者之间数据的交换，使得ANSYS所提供的节点和单元信息格式能与FLAC3D对应起来，从而建立起FLAC3D的数据模型。这个接口程序，也是实现FLAC3D建模直观化、自动化、快捷化的关键所在。

接口程序采用Visual Fortran 5.0编写，目前已经考虑了以下几点：

（1）实现了FLAC3D中常用的BRICK、WEDGE模块与ANSYS的转换，它们之间节点对应关系如图2所示；

（2）可以实现针对体或面的单个或多个INTERFACE（界面元）的生成，并对不同的INTERFACE赋予不同的ID号；

（3）针对不同力学参数或不同开挖步骤的模块，采用不同的GROUP，赋予不同的GROUP名；

（4）如有必要，也可以直接在接口程序中给模型赋予力学参数、边界条件以及初始条件等。

图2 ANSYS 模块与FLAC3D模块之间的对应关系

2.2.4 数值模型的生成

将接口程序生成的数据文件，调入FLAC3D，并加入边界条件和模型参数，即生成数值模型。

2.3 实例

图3是根据三峡白衣庵地区6号小滑坡的工程地质剖面图所建立的工程地质模型，它包含了下部的基岩和上部的土石混合体，在土石混合体中共有52个不同形状的砾石块体。经过一系列的自动处理后，可以生成直接用于 FLAC 计算的数值模型（图4）。模型共5094个单元、11284个节点。该模型由于地质情况复杂，若用FLAC3D所提供的建模方法，要花1～2 d的时间，若用模型自动生成技术，至多用半天时间即可建成，节省了60%以上的时间，而且生成的网格具有高度智能化。

图3 白衣庵6 号小滑坡工程地质模型

图4 白衣庵6 号小滑坡的数值模型

参考文献

［1］孙广忠.中国典型滑坡.北京：科学出版社，1988

［2］谢守益，徐卫亚.降雨诱发滑坡机制研究.武汉水利电力大学学报，1999，32（1）

［3］刘光华.河岸松散堆积层滑坡机制及防治对策研究——以重庆市为例.见：三峡库区地质环境暨第二届中日地层环境力学国际学术讨论会论文集.北京：煤炭工业出版社，1996

［4］郭俊仃，夏季华.超压密土坡稳定分析的非线性有限元法.见：计算机方法在岩石力学及工程中的应用国际学术讨论会论文集.武汉：武汉测绘科技大学出版社，1994

油新华1 何刚1 李晓2

（1.北京城建集团 北京 100025 2.中国科学院地质与地球物理研究所 北京 100029）

摘要 土石混合体边坡是由滑坡堆积、残坡堆积、崩坡堆积、冲洪堆积等形成的一种非均质边坡。由于其物质组成的不均质性，使其与其他边坡具有显著不同的变形性质和破坏机理。针对其物质结构特点，对土石混合体边坡应该采用细观处理技术，将之分别看成均质连续体、等效的均质连续体、非均质不连续体和结构面控制的非连续体4种情况来处理。

关键词 土石混合体边坡 细观处理技术

土石混合体边坡是由土夹碎石或碎块石以及碎石或碎块石夹土等组成的边坡［1，2］。它一般发育在第四系松散堆积层中，主要由滑坡堆积、残坡堆积、崩坡堆积、冲洪堆积等形成。土石混合体边坡在我国有着广泛的分布，如香港、广东、福建等地大面积分布的花岗岩残积土边坡，西南地区、长江三峡库区及黄河中上游广泛分布的古滑坡、崩积层边坡，西藏等地分布的冰碛土边坡等［3］。由于该类边坡通常规模较大、影响因素众多、失稳突发性强、滑移条件复杂，常给国民经济建设、人民正常生活与生命安全带来严重危害和巨大的财产损失。因此，对这类边坡进行系统分析、研究和总结具有特殊的理论和实用价值。

对于土体边坡以及岩体边坡来说，人们已经进行了长期、系统、卓有成效的研究工作，并取得了丰富的实践经验。而对于土石混合体边坡来说，由于物质组成和结构的不同而与其他边坡有着显著不同的变形特点，目前人们对于它的研究还只是处于一种定性分析的阶段，或是通过它与某种因素的相关性分析来探讨滑坡的机理［4，5］，或是通过模型实验来进行稳定性分析［6，7］，而缺少像研究土坡和岩体边坡那样的理论和技术方法。作者认为首先应该认清问题的本质，即土石混合体边坡与众不同的特点，然后以一种新的角度来分析和研究这些问题。

1 细观处理技术的概念

对岩土体的破坏机理研究依据所研究对象的尺寸大小，可分为宏观、细观和微观三类。关于这三类的划分，目前还没有一个统一的标准。按照谢强等人的观点［8］，将野外普遍发育、影响工程岩土体力学特征的断层、节理归为宏观级；将发育在岩土体结构中，用肉眼或显微镜观测，直接影响岩土体力学特征的裂纹归为细观级；将发育在矿物晶体内部，一般对岩土体的宏观力学性质没有直接影响的那些位错归为微观级。从研究实施的角度看，宏观研究多结合野外测绘及现场实验进行。由于观测手段及耗资等原因，所做的研究总是有限的。细观研究多以室内岩样力学实验与显微镜观测的方式进行。以现有的经济技术手段可行性而论，细观水平上岩土体破裂过程的研究是当前主要的研究方法。

传统的岩石力学研究采用的是宏观力学的实验与分析方法，已取得了许多瞩目的成果。随着众多重大岩土工程要求的不断提高和研究工作的逐步深入，人们逐渐认识到，由于岩土体固有的非均质性，仅从唯像学角度以宏观尺寸探讨其力学机理显然是不够的，1998年孙钧提出，应当进一步从细观甚至微观尺度，“三观”相结合，更加全面、深入地来探讨岩石力学问题。

所谓细观尺度，介于微观和宏观之间，它相对于分子、原子的微观尺度足够大，而相对于所涉及对象的宏观尺度又足够小；对于不同的材料及其不同的工程特征尺寸，它是一种相对的概念。对岩石来说，细观尺度可能落到厘米甚至米的范围内。人们认识自然现象的规律，一般总是要经历宏观—细观—微观—宏观的不断反复过程。高一级尺度的研究总是基于低一级尺度的研究之上。这就意味着宏观现象的发生能够在微观研究中得到解释，而微观现象的发生又能够在细观研究中得到解释。

2 细观处理的必要性

一些非均质土，譬如砾质沉积物和崩积物，可以含有包括漂石在内的范围比较广的不同大小的颗粒，从而给岩土工程师提出了一些特殊的问题。对于这些物质的采样和试验的困难已经有过很好的报道［9］，并且得到了一定程度的解决［10］。但是，由于评价这种边坡的稳定性而产生的方法论问题却没有引起足够的重视。常规土的边坡分析方法实际上是建立在级配相对良好的基础之上的，其性状假定为与单元尺寸无关。但是，对于土石混合体边坡来说，由于其物质组成的不均质性，使得有必要从单元尺寸即细观角度上来处理。

土石混合体的力学性质主要受其粒度组成的影响，即受其组成颗粒的大小、形状、均一情况、颗粒的磨圆度特别是孔隙度的影响。如有粘土、粉土或有机质参与时，性质将更加复杂。

图1 土石混合体边坡破坏的三种方式

当土石混合体边坡中含有少量大的颗粒时，一般可能发生穿过颗粒、偏离颗粒和包含颗粒的破碎区变宽3种破坏情况［11］（图1），其破坏机理符合最小抗剪强度原理。哪种现象发生取决于颗粒和周围基质的相对强度、颗粒和破碎区的相对大小、大颗粒的含量以及应力水平。例如穿过颗粒可能在颗粒较弱的情况下发生。如果仅有少数几个大颗粒，则倾向于偏离颗粒破坏，此时破裂区域的几何形状将由于颗粒的存在而在一定程度上略有改变。Hencher ＆ Martin［12］认为当大颗粒含量超过30%时，如香港的崩积物，通过填充物的破裂面一般都要受到大颗粒的干扰。如果存在许多大的颗粒，将会形成锯齿状破裂面（图2）。如果颗粒与边坡尺寸相比太大，或者填充物倾向于沿薄的破裂面破坏，则破碎区变宽的可能性要小于偏离颗粒式的破坏。偏离颗粒破坏在可能发生剪胀的低应力区比限制剪胀的高应力区更可能发生，而在高应力区则可能导致破碎区变宽。

图2 土石混合体边坡的锯齿状破裂面

在上述的讨论中，似乎大颗粒要么增加边坡的稳定性，要么没有大的影响。然而，大颗粒可能降低边坡的稳定性。当沿大颗粒与充填物的交界面或大颗粒之间的剪切强度小于充填物本身的强度时，这种情况就会发生，而且占有相当大的比例。如果大颗粒占据不太理想的位置，则这种作用更加明显。另外还应注意到大颗粒的存在将十分显著地增加坡体的排水性。

从以上的分析可以看出，为了更准确地研究边坡的稳定性，应从细观角度来进行处理。对于重大的工程，应尽可能地进行详细勘察，以确定土石混合体边坡的不均质地质剖面，将其中影响边坡稳定的大的砾石块体找出来，并在稳定性分析中将这些块体考虑进去。

3 细观处理的方法

对于所要研究的具体工程来说，首先要在相对充分的勘察资料基础上建立尽可能反映工程实际情况的岩土工程模型。其中心思想就是首先应该把它看成一种什么样的材料或介质，因此它的内容应包括与边坡有关的基本地质条件，如边坡的尺寸、所含介质的颗粒大小及其均匀性、岩体的风化程度、节理等不连续面的发育程度、岩层的相互组合等。对于不同的工程条件，应该有不同的处理办法。下面是根据土石混合体材料的力学性质研究提出的针对不同边坡工程的处理办法，即根据工程地质条件将边坡材料看成是均质或非均质、连续还是非连续体。

3.1 均质连续体

图3 可看成均质连续体的边坡

最常用的方法是将土石混合体作为一种实质上的均质连续体，它把最主要物质的参数作为区域代表性的参数。这种方法适合于土石混合体中含有少量的砾石块体而且块体的大小和边坡的尺寸相比特别小不足以影响整体力学性质的情况（图3）。这时，将岩土体看成是均质材料是有效的，其工程参数和纯土体相比没有明显的变化，而且不可能发生块体与土体之间结构控制的行为。在这种情况下，没有必要进行野外原位试验，而只需实验室试验即可获得所需的参数。

3.2 等效的均质连续体

当土石混合体中虽然含有较多的砾石，但是砾石分布比较均匀且大小与边坡尺寸相比也比较小时，可以把它看成一种等效的均质连续体（图4）。由前面的分析可知，由于砾石的存在而使得其力学性质有所增强，但是因为试验方法的困难，许多工程师仍将最软弱成分的性质作为区域的代表性质，这种方法一般是保守的。

图4 可等效成均质连续体的边坡

尽管对之进行了整体简化，即将某种物质的性质代表大区域的性质，但是材料对于扰动的敏感性意味着实验室试验得出的数据是不可信的。因此，许多设计者宁愿依赖于原位试验和经验公式。另外值得注意的是，在土石混合体中，砾石与土体之间的接触面在受力变形时变成了不连续面。在许多情况下，不连续面的出现不能被忽略，即使不存在明显的机理，这时不连续面可能占有主要的地位，它对强度的影响是非常大的。

3.3 非均质、非连续体

当土石混合体边坡中含有和边坡相比尺寸较大的砾石块体，即使含量不高但位置较为重要时，则必须要把它看成是非均质、非连续体（图5）。这时，大的砾石块体的存在必然影响着边坡的破坏形态和机理，在进行边坡稳定性分析时，一定要考虑大块体的影响。

图5 非均质、非连续体的土石混合体边坡

在处理非均质岩土体时，由于软弱基质中含有强度较大的核岩（corestone），而使得确定具有实际意义的参数变得十分困难。这时，不可能以岩土体的代表尺寸来测试试样，因此从核岩到较弱的物质的岩土体强度和模量的分布，要么被忽略（假定在没有结构控制时一般是保守的），要么从理论和模拟上来处理［13］。在这方面，许多学者都对之进行了不同程度的研究。Anon以及West等人［11～14］在对较高砾石含量的堆积体边坡分析时，试图计算由不规则破裂面剪胀引起的额外的剪切强度，其几何形状由交错分布的碎屑来控制。Irfan和Tang［15］给出了由于堆积体中强度较高包裹体引起的附加强度的试验性指导原则。只要它们具有相同的特点，相同的方法或许可以适用于不均质岩土体中。Blight［16］以及随后出版的Jaros的工作［17］报道了在预测包含核岩的地基的沉陷时获取实际参数的困难。Jaros将一建筑的低沉陷归结于建筑物下残积土中石英核岩的影响，然而Blight却认为所观察到的沉陷同样可以用稍微不同的参数归因于母岩残积土。DeMello［18］在预测不同压缩性物质区域中建筑物的沉陷时使用相同的方法考虑了同样的问题，他将之归因于不同物质的比例和不同压缩性。

3.4 结构面控制的不连续体

当边坡中含有很多的块石，比如弱风化岩体、碎裂岩体时，其中土体充填物的含量较少，所起的作用不大，起主导作用的是块体与块体之间的节理等不连续面，这时要把它看成由结构面控制的不连续体（图6）。很清楚，在岩土体行为完全被不连续面控制的时候，可以使用正常的岩体力学方法，但是必须注意到风化对沿潜在滑动面剪切强度的影响。

图6 节理面控制行为的土石混合体边坡

由以上的分析可以看出，对于上述四种方法的区分主要是以土石混合体中砾石或块体的含量以及尺寸大小为标准的。但是，其含量多少以及大小如何还没有一个具体的界限，只能根据现场的具体情况和工程经验来确定。作者给出如下建议值：当含石量少于10%时，可以作为第一种情况考虑；当含石量在10%～25%时，属于第二种情况；当含石量为25%～70%时，属第三种情况；当含石量大于70%特别是达到90%时，则可以作为节理岩体来考虑。对于第一和第四种情况，可以使用土力学和岩体力学的方法来处理。而对于第二种情况，可以根据文献中的研究成果来确定土石混合体的强度参数，以进行相应的计算和分析。对于第三种情况，具体方法有两种：

一是经过详细的勘察，将大的砾石块体的位置探明出来，利用土石混合体的实测结构模型进行数值模拟，当然这是一种理想方法，要想达到是非常困难的；二是根据局部的统计资料利用蒙特-卡洛方法模拟砾石的分布情况，在随机结构模型的基础上进行数值模拟，这是一种比较实用、切实可行的方法。

致谢 在本文的写作过程中，得到了许多现场和科研机构工作人员的大力协助，在此向他们，特别是原中国地质环境监测院殷跃平副总工程师、中国科学院地质与地球物理研究所张年学研究员、曲永新研究员表示衷心的感谢！

参考文献

［1］郭俊仃，夏季华.超压密土坡稳定分析的非线性有限元法.见：计算机方法在岩石力学及工程中的应用国际学术讨论会论文集.武汉：武汉测绘科技大学出版社，1994

［2］湖南省水利水电勘测设计院.边坡工程地质.北京：水利电力出版社，1983

［3］孙广忠.中国典型滑坡.北京：科学出版社，1988

［4］王发读.浅层堆积物滑坡特征及其与降雨的关系初探.水文地质工程地质，1995，（1）

［5］谢守益，徐卫亚.降雨诱发滑坡机制研究.武汉水利电力大学学报，1999，32（1）

［6］刘光华.河岸松散堆积层滑坡机制及防治对策研究——以重庆市为例.见：三峡库区地质环境暨第二届中日地层环境力学国际学术讨论会文集.北京：煤炭工业出版社，1996

［7］新杰罗·科玛达.大坝水库运行诱发的边坡不稳定性及其产生的机理（内部资料）

［8］谢强等.岩石细观力学试验与分析.成都：西南交通大学出版社，1997

［9］油新华.土石混合体的随机结构模型及其应用研究.北方交通大学博士学位论文，2001

［10］ Watkins D J.Problems in measuring the properties of glacial tills.M.Sc.Thesis.University of London，Imperial College，1969

［11］ West L J，Hencher S R.Assessing the stability of slopes in heterogeneous soils.Landslides Balkema：Rotterdam，1991

［12］ Hencher S R，Martin R P.The description and classification of weathered rocks in Hong Kong for engineering purposes.Proceedings of the 7thSouth East Asia.Geotechnical Conference.Hong Kong.1982，1

［13］ Lindquist E S，Goodman R E.Strength and deformation properties of a physical model mélange.Proc.1st NARM Symp.Austin Balkema Rotterdam，1994

［14］ Anon.Effect of the coarse fractions on the shear strength of colluvium.Geotechical Engineering Office Report 23.Hong Kong，1993

［15］ Irfan R Y.Potential construction problems for bridge foundations in areas underlain by complex marble formation in Hong Kong.Engineering Geology of Construction Geological Society.London Engineering Geology Special Publication，1995

［16］ Bligth R D.Construction in tropical soils.Proceedings of 2ndInternational Conference on Geomechanics in Tropical soils Singapore.Balkema，1988

［17］ Jaros M B.The settlement of two multi-storey buildings on residual Venterdorp lava.South African Institution of Engineers Symposium on Soil-Structure Interaction.Durban，S.Africa，1978

［18］ Demello V F B.Thoughts on soil engineering applicable to residual soils.Proceedings of 3rdSoutheast Asian Conference on Soil Engineering.1972