地质力学学报  2008, Vol. 14 Issue (2): 149-157
引用本文
汪吉林, 丁陈建, 吴圣林. 基于FLAC3D的复杂地貌三维地质建模[J]. 地质力学学报, 2008, 14(2): 149-157.
WANG Ji-lin, DING Chen-jian, WU Sheng-lin. 3D GEOLOGICAL MODELING OF COMPLEX LANDFORMS BASED ON FLAC3D[J]. Journal of Geomechanics, 2008, 14(2): 149-157.
基于FLAC3D的复杂地貌三维地质建模
汪吉林1 , 丁陈建2 , 吴圣林1     
1. 中国矿业大学资源与地球科学学院, 江苏 徐州 221006;
2. 徐州中国矿大岩土工程新技术发展有限公司, 江苏 徐州 221006
摘要:针对复杂地貌和地质构造条件下数值模拟的建模问题, 基于FLAC3D平台提出了从局部到整体的建模思想, 采用Brick单元和Tetrahedron单元通过编程生成复杂地貌条件下的三维地质体模型。利用趋势面分析法拟合出层状连续地质体, 将其定义为不同的"群", 建立了多层状三维地质体结构模型。采用interface界面模拟出了断层面, 并与趋势面相结合, 构建了断层上下盘中的地层。认为通过数学分析结合计算机编程, 可以实现直接在FLAC3D中进行复杂地貌和地质构造条件下的地质体建模。
关键词数值模拟    建模    FLAC3D    趋势面分析    复杂地貌    
文章编号:1006-6616(2008)02-0149-09
3D GEOLOGICAL MODELING OF COMPLEX LANDFORMS BASED ON FLAC3D
WANG Ji-lin1 , DING Chen-jian2 , WU Sheng-lin1     
1. College of Resources and Earth Science, China University of Mining and Technology, Xuzhou 221008, Shandong, China;
2. CUMT Geotechnical Engineering & New Technology Development Co., Ltd., Xuzhou 221008, Shandong, China
Abstract: With respect to the problem of numerical modeling for simulation in conditions of complex lan dforms and geological structure, the idea of modeling from the part to whole was put forward based on the FLAC3D platform.Using the primitive mesh shapes of brick and tetrahedron, a 3D model of the geologic body under complex landform conditions was constructed by means of computer programming.A lay ered continuous geological body was fitted up using the trend surface analysis and defined as different "groups".Thus a multi-layered 3D structure model of the geologic body was constructed.The fault plane in the geologic body was simulated by using interface.On that basis, combined with the trend surface analysis, the strata of the hanging wall and footwall of a fault were simulated.The authors consider that a 3D model of the geologic body with complex landforms and geological structure can be constructed in FLAC3D directly by using mathematical analysis and computer programming.
Key words: numerical simulation    modeling    FLAC3D    trend surface analysis    complex landform    

自然界中地质体几何形态是及其复杂的, 描述地质体形态的传统方法是采用平面图、切面图、等高线图、剖面图等多种图件相结合。随着计算机技术的进步, 三维地质体可视化的技术得到了很大的发展[1~6], 与此同时, 数值模拟方法在采矿工程和地质工程等领域也得到了广泛应用[7~8]。由于数据结构的差异, 地质体三维可视化的网格与用于数值模拟计算的网格之间存在本质差异, 而在大多数数值模拟分析中, 模型建立、空间单元划分等前处理工作繁琐、甚至难于进行, 妨碍了数值模拟技术的应用, 若简化地质模型则又会影响模拟效果, 因此有研究者采用数据转换的方法在“可视”与“可算”之间建立桥梁[9~10]

数值模拟的理论方法包括有限单元法、离散单元法、边界单元法、快速拉格朗日分析等及其相互耦合[7], 其中快速拉格朗日分析(FLAC)可以较准确地模拟材料的屈服、塑性流动、软化乃至大变形, 是岩土工程领域的通用软件[10], 但FLAC3D不能直接进行图形处理, 其建模难度更大, 有研究者利用其他有限元程序(ANSYS)中进行复杂地质体建模, 通过接口程序转换为FLAC3D软件平台中的模型[11], 这需要熟悉ANSYS的文件结构, 且难以甄别和控制程序执行中产生的差错。

无论从可视化三维地质体数据库中转换, 还是在不同的数值模拟软件之间转换, 都会增加三维地质建模的工作量, 同时对使用者提出更多的技术要求和操作技巧。针对FLAC3D建模难度大的问题, 特别是存在复杂地貌和地质构造的条件下, 笔者提出通过数据点采集、趋势面拟合和计算机编程, 直接在FLAC3D中建立复杂地形地质条件的三维地质体模型, 较好地解决了FLAC3D中地质建模难的问题。

1 地形地质条件概述

本文以山西307国道复线阳泉市天子庙隧道地段作为模拟对象, 说明基于FLAC3D软件平台建立三维地质体模型的方法和过程。天子庙隧道为一条双线并行隧道, 左线长888.10m, 右线长819.05m, 地形较复杂, 地表高差较大, 在研究区范围内为900~1054m (图 1)。研究区场地下部开采过4层煤(3#、8#、12#煤和15#煤), 累计煤层厚度最大达11m。研究区地表露头为上石盒子组(P2s)下段~下石盒子组(P1x), 主要为红褐色砂泥岩及灰白色砂岩, 黄色或黄绿色泥岩及砂岩等。160~180m勘探深度范围内的地层为山西组(P1s) ~太原组(C3t)的煤系地层。

图 1 天子庙隧道地形图 Figure 1 Topographic map of the Tianzi temple tunnel

研究区范围内地层总体上较平缓, 以NW走向为主, 倾向SW, 局部有一定的变化, 倾角5~10°。节理裂隙发育, 倾向以SE和SW向为主, 一般倾角较大。研究区内的岩石组合在力学强度上具有层状分布的特征, 故根据研究区的岩石组合特征, 将其按工程地质特性划分为砂岩岩组、泥岩岩组、灰岩岩组、砂泥岩互层岩组和破碎岩组(采空区堆积岩体)等5个岩组类型。

2 复杂地貌的三维建模 2.1 地形数据采集与处理

图 1所示为采空区场地的研究范围, 其中NW—SE方向包含隧道的图框为模拟的地质体范围。为了建立模拟对象的三维立体几何模型, 首先需要采集地表的高程数据, 通过手工统计写入文本文件以便于编程调用, 模拟范围长×宽=800×300m, 以20m间距划分为16行、41列, 共656个地表高程数据点, 其立体示意图如图 2 (a)所示。为了使模拟的地形更符合实际, 将20m间距的数据进行样条插值、加密处理, 可以得到更加细致的、10m间距三维地形模型(图 2 (b)), 数据加密可以通过编程实现。

图 2 天子庙隧道地形立体示意图 Figure 2 3D sketch map of the Tianzi temple tunnel
2.2 FLAC3D空间单元生成

数值模拟建模首先要将模拟对象的几何体划分为若干微小单元体。FLAC3D采用的标准单元体主要有六面体块体(Brick)、楔形体(Wedge)、角锥体(Pyramid)、四面体(Tetrahedron)、圆柱体(Cylinder)等多种单元体, 以适应不同模拟对象的需要[12]。对于不规则的地质体而言, 综合考虑单元划分的有效性和可操作性, 采用六面体和四面体比较适宜(图 3), 前者划分简便, 文件容量小, 但模型不够精细; 后者划分繁琐, 文件容量大, 模型精细程度高。

图 3 FLAC3D中的六面体块体单元和四面体单元 Figure 3 Brick and tetrahedron in FLAC3D

FLAC3D中一般地网格划分方法是先定义单元体类型和x, y, z各方向的分割数, 再定义各个节点的坐标, 即可生成一个网格化的整体模型。这种方法对于形状规则的几何体是可行的, 比如机械零部件、工程建筑等, 但对于形态复杂、没有规律的地质体, 即使将模型几何形态简化也非常难以进行划分, 导致了FLAC3D在岩土工程、地质工程领域的推广应用受到限制。

从整体到局部的建模方法有局限性, 如果改为从局部到整体进行建模, 则可以有效地解决复杂地质体建模的问题。具体思路是:(1)根据建模对象的平面网格划分(mn列), 得到x, y坐标集合; (2)建模对象的地表高程数据采集点与平面网格节点是对应的, 根据最大高差确定垂向的网格划分数k; (3)将平面网格每个节点处的地形垂向高差进行k等分, 因此得到一个三维的数组(m×n×k); (4)逐行逐列、从底向上, 生成一个个微小的六面体块体单元, 这些小六面体单元堆砌而成一个整体的模型; (5)从图 3可以看出一个六面体单元实际上可以由6个四面体单元组成, 如果需要对不规则地质体进行更精细地划分, 可以采用四面体单元堆砌生成地质体模型。

由于微小单元体很多, 必须通过编程生成FLAC3D可执行的命令文件, 以上建模思想可以通过任何一种编程语言实现。根据以上的建模思想, 采用20m间距的平面网格, 垂向高差10等分, 对研究区隧道底板水平以上的地质体(900~1054m)进行建模, 可得到由六面体单元和四面体单元生成的三维地质体模型(图 4)。对于本文中的模型, 平面网格若采用10m间距、六面体单元, 则生成的模型文件较大, 导致生成速度缓慢; 采用10m间距、四面体单元时, 因生成速度过于缓慢而不得不中止。为了取得较好的模拟效果并保证可操作性, 建模中的单元类型、网格尺寸应当根据模型大小和计算精度共同确定。

图 4 900m高程以上三维地质体模型 Figure 4 3D geological model at an elevation of >900m
3 层状连续地质体的拟合 3.1 地层底板几何形态的拟合

建立三维的地质体模型, 必须先准备好模拟对象的空间形态数据。在天子庙隧道的研究区场地, 模拟的地质体为多层状结构, 且地层具轻微褶皱的形态。通过收集采矿资料、现场钻探、山体工程地质测绘等工作, 得到了地下煤层的底板等高线图, 并根据岩层的工程地质特性自15#煤层底板向上共划分了19个岩组。为了在建模中使程序能够自动识别不同的岩组, 对高程资料充分的3#、8#、12#、15#煤底板和地表第14~19岩组的底板进行趋势面拟合, 其他高程资料不够充分的岩组底板, 则以相邻或相近的已知岩组的底板近似递推。

趋势面分析是应用数理统计学中的多元回归分析原理研究地质变量的空间分布和变化规律的一种数学方法。在地质工作中应用较多的是多项式趋势面分析, 即用多项式函数来拟合地质变量的空间变化趋势。地质上常用的是二维趋势面分析, 当多项式的次数增高时, 多项式能逼近任意连续函数[13]。一般说来, 多项式的次数越高, 趋势值与实际值的偏差越小, 但次数越高, 计算越复杂, 同时造成趋势面的过多曲折, 反而降低了效果。研究区地层总体上比较平缓, 略具褶皱形态, 通过试算, 选择2次趋势面拟合即可得到满意的效果(图 5)。2次趋势面的拟合方程为

图 5 12煤底板等高线原图(左)与拟合图形(右) Figure 5 Contour map of No.12 coal seam floor (the left is the original map, and the right is the fitting map)

式中Zi为拟合的高程, xi为横坐标, yi为纵坐标, (b0, b1, …, b5)为待定系数。各煤层和部分岩组拟合方程的待定系数及拟合效果检验如表 1所示。

表 1 各煤层和部分岩组底板趋势面拟合成果 Table 1 Fitting parameters of trend surface analy sis for various coal seams and some floors
3.2 利用内嵌Fish语言生成地层

模拟对象除了地上的山体部分, 还有地下的煤层、岩层, 将地下720~900 m水平的地质体进行网格划分, 并与上部山体的模型拼接, 就可以建立其研究对象的整个立体的模拟, 由于煤层、岩组的厚度不一, 为了生成精确的地层, 采用了10m间距的平面网格, 垂向上在不同高程段的分割密度随地层的厚薄而疏密不同, 研究对象的地质体共划分为237600个Brick单元(图 6)。

图 6 研究对象整体的三维初始模型 Figure 6 Original 3D model of the whole study object

以上生成的模型尚未进行地层(岩组)的划分, 可以通过Fish语言简单编程而实现。Fish语言是FLAC3D中内嵌的编程语言, Fish语言小巧灵活、语法结构简单, 用户通过它自定义变量和函数, 以及进行逻辑计算, 可以大大地扩展FLAC3D的功能。程序设计中, 对整个模型的单元体进行搜索, 获取各单元体的x, y, z坐标, 并代入前文拟合的趋势面方程中计算, 根据计算结果的z坐标值, 判断该单元处于哪个层位。将不同的煤、岩组定义为不同的“群”, 就可得到近似于实际的三维地质体模型结构(图 7)。

图 7 研究对象的三维地质体模型结构 Figure 7 Structure of the 3D geological model of the study object
4 断层面的模拟

建立三维地质体模型时, 将不可避免遇到断层问题。在天子庙隧道研究区建模中, 由于地下采煤的缘故, 地表多年来塌陷严重, 山顶发育大量裂缝, 而研究区断层稀少且落差小, 所以, 建模中将研究对象视为层状散体结构, 未考虑断层影响, 应该近似符合实际情况。

为了完善三维地质体建模问题的探讨, 本文针对研究区的地质条件, 假定存在一规模较大的断层, 讨论存在断层的条件下如何建立三维地质体模型, 因主要探讨断层问题, 故将地层简化, 仅考虑煤层(并适当放大尺寸)被断开的情形。

FLAC3D中的interface界面可用来模拟体与体之间节理、断层或摩擦界面的影响。界面可以附着于一个子网格上或多个子网格之间, 也可以放置在空间内的任何位置。当子网格接触界面时, 产生剪切力和法向力。FLAC3D允许沿界面是子网格之间相互作用的平面, 这些界面可以滑动、分离。

模拟中选择断层面上三个已知点, 计算出断层面的平面方程, 定义紧邻空间几何平面位置的网格分界面为interface, 以模拟断层面。与前文描述的趋势面相结合, 再通过Fish语言编程, 就可以实现在研究对象的模型中构建模拟的断层面和上下断盘中的地层(图 8)。

图 8 研究对象三维地质体中断层的构建 Figure 8 Fault in the 3D geological model of the study object

由于模型的单元体是Brick单元, 所以模拟的断层面为锯齿状, 若模型的单元体为Tetrahedron单元, 则可以较好的消除锯齿的影响, 但本文的模拟对象规模太大, 不适宜采用Tetrahedron单元, 若模拟对象规模较小, 则可以采用Tetrahedron单元建模并构建断层面。

5 结论

(1) 基于FLAC3D平台可以采用如下方法实现复杂地貌条件的地质体建模:将模拟对象划分为mn列的平面网格; 根据最大高差将垂向划分为k层网格; 将平面网格每个节点处的地形垂向高差进行k等分, 得到一个三维的数组(m×n×k); 逐行逐列、从底向上, 生成一个个微小的块体单元, 堆砌而成一个整体的模型。

(2) 一般情况下模型单元体可采用Brick单元, 如果需要对不规则地质体进行更精细地划分, 可以采用Tetrahedron单元。为了取得较好的模拟效果并保证可操作性, 建模中的单元类型、网格尺寸应当根据模型大小和计算精度共同确定。

(3) 应用趋势面分析法拟合层状连续地质体, 通过FLAC3D中内嵌的Fish语言编程, 搜索整个模型单元体的坐标, 并代入趋势面方程, 根据计算出的z坐标值, 判断该单元体所处层位, 将不同的地层定义为不同的“群”, 得到近似于实际的三维地质体模型结构。

(4) 根据断层面的平面方程, 定义紧邻空间几何平面的网格分界面为interface界面, 可以模拟断层面。与趋势面相结合, 可以实现在模型中构建模拟的断层面和上下断盘中的地层。Tetrahedron单元可以消除Brick单元的锯齿效应。

参考文献
[1]
Michaud LH. Development of a 3D geologic model and database for dragline mine planning on a microcomputer[J]. International Journal of Surface Mining & Reclamation, 1992, 6(4): 179-185.
[2]
曹代勇, 李青元, 朱小弟, 等. 地质构造三维可视化模型探讨[J]. 地质与勘探, 2001, 37(4): 60-62.
[3]
王纯祥, 白世伟, 贺怀建. 三维地层可视化中地质建模研究[J]. 岩石力学与工程学报, 2003, 22(10): 1722-1726. DOI:10.3321/j.issn:1000-6915.2003.10.026
[4]
何满潮, 李学元, 刘斌, 等. 侵入型岩体三维可视化构模技术研究[J]. 煤田地质与勘探, 2004, 32(4): 29-32. DOI:10.3969/j.issn.1001-1986.2004.04.010
[5]
武强, 关文革, 贾丽萍, 等. 面向矿区复杂地质体的四面体生成算法[J]. 中国矿业大学学报, 2005, 34(5): 617-625. DOI:10.3321/j.issn:1000-1964.2005.05.016
[6]
ZHOU Ye, LIN Ge, GONG Fa-xiong, et al. Application of 3D visual techniques in Daliuta coal mine[J]. Journal of China University of Mining and Technology, 2006, 16(1): 77-81.
[7]
唐辉明, 晏鄂川, 胡新丽. 工程地质数值模拟的理论与方法[M]. 武汉: 中国地质大学出版社, 2001.
[8]
汪吉林, 姜波, 王超勇. 矿井构造应力场的模拟研究[J]. 地质力学学报, 2007, 13(3): 239-246. DOI:10.3969/j.issn.1006-6616.2007.03.006
[9]
侯恩科, 吴立新, 李建民. 三维地学模拟与数值模拟的耦合方法研究[J]. 煤炭学报, 2002, 27(4): 388-392. DOI:10.3321/j.issn:0253-9993.2002.04.012
[10]
王明华, 白云. 层状岩体三维可视化构模与数值模拟的集成研究[J]. 岩土力学, 2005, 26(7): 1123-1126. DOI:10.3969/j.issn.1000-7598.2005.07.024
[11]
廖秋林, 曾钱帮, 刘彤, 等. 基于ANSYS平台复杂地质体FLAC3D模型的自动生成[J]. 岩石力学与工程学报, 2005, 24(6): 1010-1013. DOI:10.3321/j.issn:1000-6915.2005.06.019
[12]
Itasca Consulting Group, Inc. FLAC user' s manual [M]. Minnesot a: State University of Minnesota, CSA, 1997.
[13]
韩金炎. 数学地质[M]. 北京: 煤炭工业出版社, 1987.