一种基于Fortran语言的3DEC离散元数值模型建模方法

摘要

本发明公开了一种基于Fortran语言的3DEC离散元数值模型建模方法，包括：首先获取3DEC离散元数值模型区域的三维等高线地形图，然后对地形图进行插值和栅格化处理得到间距规则排列的离散点坐标数据，随后基于Fortran语言进行建模代码的编写和编译，实现栅格化离散点坐标数据生成为3DEC离散元数值模型命令流，最后采用3DEC程序读取模型命令流运行生成3DEC离散元数值模型。本发明计算高效、操作简单、灵活度高、便捷性好，能够对地形复杂且范围较大的研究区域进行高精度、多岩层、可变块体单元尺寸的三维离散元模型构建。

说明书

技术领域

本发明涉及岩土工程数值模拟研究领域，特别涉及一种基于Fortran语言的3DEC离散元数值模型建模方法。

背景技术

3DEC是世界范围内第一款采用离散单元法作为基本理论进行开发并商业化的高级三维分析程序，特别适用于节理岩体等不连续介质的力学特性和破坏特征研究，目前在岩体工程领域得到广泛应用。岩体由岩块和结构面组成，3DEC程序数值模拟中将岩体视为介质离散的散体集合体，将岩体和结构面定义为两个基本单元；在连续介质力学方法基础上，引入接触(如结构面)的概念，把接触作为块体(如完整岩块)之间的边界，利用连续力学方法描述块体的力学行为、非连续力学方法描述接触的力学行为；采用接触算法搜索块体边界接触条件并计算接触受力状态，当接触出现屈服形成剪切滑动或张开时，块体发生运动位移(平动、转动)甚至破坏现象。因此，3DEC程序能够准确高效地模拟不连续介质在承受静载或动载作用时的力学响应，但是，其较弱的模型构建功能一直制约着程序的使用和普及。

3DEC离散元数值模型的块体单元数量决定了模型精度和计算效率。若块体单元数量较多，则数值模型模拟实际工程的精度较高，但是相应的程序计算时间亦较长，甚至超过计算机计算能力而不能计算。反之，若块体单元数量较少，则模型计算效率较高，但是模型模拟实际工程的准确性则相应降低。四面体和三棱柱块体是使用较多、精度较高的3DEC程序建模块体单元，在相同建模精度条件下，四面体块体单元建模所需数量是三棱柱块体单元的数倍，如图11所示。针对研究区域为小范围的模型，可以采用常用的四面体块体单元进行模型构建，块体适应性好，如采用Rhinoceros软件建立三维实体模型并通过Griddle插件进行网格处理生成3DEC模型，或采用专利号CN109242965B所述的通过Midas软件建模并基于C#语言的转换插件转换生成3DEC模型。但是对于研究区域为大范围的模型，若采用四面体块体单元建模往往导致块体数量较多，超过计算机计算能力而不能计算，所以采用三棱柱块体单元建模是更合理的选择。专利号CN107562833B所述的基于MATLAB软件插值处理地表曲面数据，通过三棱柱块体单元建立三维数值模型，只考虑了地表建模，不能实现多岩层模型的构建。论文《复杂地质体的3DEC快速建模》采用Surfer软件和Fish语言实现了地质体的3DEC快速建模，但是Fish语言具有其独特的源代码表达方式，对于一般工程技术人员来说，理解和操作都比较困难，难以推广使用。

Fortran语言是世界上第一个被正式推广使用的高级编程语言，应用历史长久，使用普及度高，尤其适用于数值计算领域，具有执行效率高、机密性好、标准化程度高等优点，因此，有必要发明一种基于Fortran语言的3DEC离散元数值模型建模方法。

发明内容

本发明的目的是在于提供了一种基于Fortran语言的3DEC离散元数值模型建模方法，具有计算高效、操作简单、灵活度高、便捷性好等优点，能够对地形复杂且范围较大的研究区域开展高精度、多岩层、可变块体单元尺寸的三维离散元数值模型构建。

为进一步实现上述目的，本发明采用以下技术方案：一种基于Fortran语言的3DEC离散元数值模型建模方法，包括以下步骤：

步骤a：获取3DEC离散元数值模型区域的三维等高线地形图；

通过考虑研究对象和边界条件的影响，确定3DEC离散元数值模型的研究区域，裁剪三维等高线地形图得到矩形平面形状的研究区域地形图，并对研究区域地形图的四周角点进行相应的高程点数据补充。

步骤b：对地形图进行插值和栅格化处理得到间距规则排列的离散点坐标数据；

通过Surfer或MATLAB软件选择合适的插值函数(如克里金函数、线性函数或三次函数)，并定义合适的间距进行栅格化处理得到间距规则排列的离散点坐标数据，其中栅格间距确定了建模块体单元的尺寸和数量，即决定了模型精度和计算效率；为了实现离散元数值模型计算既精确又高效，可以将数值模型的关键研究区域栅格间距减小，非关键研究区域的栅格间距适当增大。

若模型中包含多个岩层，则将同一位置离散点坐标数据中增加相应的坐标Z值，形成GRID_MODEL.DAT文件，内容如下：

其中，M表示X方向的离散点数，N表示Y方向的离散点数，I表示岩层面数，I层则有I+1个Z值，即：Z

步骤c：基于Fortran语言进行建模代码的编写和编译，实现栅格化离散点坐标数据生成为3DEC离散元数值模型命令流；

通过Fortran语言编写和编译建模代码，创建3DEC建模程序；采用3DEC建模程序将栅格化离散点坐标数据进行循环遍历和语言编写，生成由三棱柱块体单元组成的包含多个岩层的3DEC离散元数值模型命令流，即3DEC_MODEL.DAT文件，内容如下：

其中，POLY表示多面体块体创建命令，REGION I表示岩层编号，PRISM表示棱柱命令，a表示棱柱下底面，b表示棱柱上底面，(x

步骤d：采用3DEC程序读取模型命令流运行生成3DEC离散元数值模型；

采用3DEC程序读取模型命令流，运行生成由三棱柱块体单元组合的初步3DEC离散元数值模型，通过3DEC程序中JOIN on命令可将三棱柱块体单元合并成一个整体块体，还可以通过3DEC程序中JSET、HIDE、SEEK等命令实现断层、结构面等特殊地质构造的构建，生成最终3DEC离散元数值模型。

通过上述技术方案，采用Fortran语言进行建模代码的编写和编译，解决了建模程序执行效率不高、普及应用不强的技术难题；采用不同研究区域选择不同尺寸建模块体单元的方案，解决了既能保证研究对象建模精度又能实现计算高效的建模难题。

与现有技术相比，本发明具有如下有益效果：

(1)本发明的编程语言采用Fortran语言，应用历史长久，使用普及度高，尤其适用于数值计算领域，具有执行效率高、机密性好、标准化程度高等优点。

(2)本发明的3DEC块体单元采用三棱柱块体单元，既能保证复杂地形的建模精度，又能降低建模块体单元数量，实现对地形复杂且范围较大的研究区域开展高精度三维离散元数值模型构建。

(3)本发明的3DEC离散元数值模型建模方法，能够对数值模型的不同研究区域实现不同建模块体单元的尺寸选择，操作简单、灵活度高、便捷性好，能够建立既地形精确又计算高效的3DEC离散元模型。

(4)本发明的3DEC离散元数值模型建模方法，包含多个岩层的建模信息，可以快速准确地进行岩层划分，实现多岩层的3DEC离散元模型快速建模。

附图说明

为了更清楚地说明本发明的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其它附图。

图1为本发明基于Fortran语言的3DEC离散元数值模型建模方法的流程图；

图2为本发明实施例中数值模型研究区域的三维等高线地形图；

图3为本发明图1补充四周角点高程点数据的地形图；

图4为本发明图2插值和栅格化处理后间距规则排列的离散点示意图；

图5为本发明实施例步骤b中GRID_MODEL.DAT文件的部分内容；

图6为本发明实施例步骤c中建模程序打开界面；

图7为本发明实施例步骤c中建模程序运行成功界面；

图8为本发明实施例步骤c中3DEC_MODEL.DAT文件的部分内容；

图9为本发明实施例的初步3DEC离散元数值模型图；

图10为本发明实施例的最终3DEC离散元数值模型图；

图11为3DEC程序中四面体和三棱柱块体单元的建模对比示意图。

具体实施方式

为了更清楚地说明本发明实施例和/或现有技术中的技术方案，下面将对照附图说明本发明的具体实施方式。显而易见地，下面描述中的附图仅仅是本发明的实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图，并获得其他的实施方式。

一种基于Fortran语言的3DEC离散元数值模型建模方法，如图1所示，包括以下步骤：

步骤a：获取3DEC离散元数值模型区域的三维等高线地形图；

通过考虑研究对象和边界条件的影响，确定3DEC离散元数值模型的研究区域，裁剪三维等高线地形图得到矩形平面形状的研究区域地形图，平面尺寸为1100m×250m，如图2所示，并对研究区域地形图的四周角点进行相应的高程点数据(即坐标(X Y Z))补充，如图3中洋红色高程点所示。

步骤b：对地形图进行插值和栅格化处理得到间距规则排列的离散点坐标数据；

通过Surfer或MATLAB软件选择合适的插值函数(如克里金函数、线性函数或三次函数)，并定义合适的间距(此处的合适的间距指在满足模型计算效率下尽可能保证模型精确，即一般情况下尽可能采用较小的栅格间距，但是需要控制离散元模型建模块体单元数量最好不超过10万个)进行栅格化处理得到间距规则排列的离散点坐标数据，其中栅格间距确定了建模块体单元的尺寸和数量，即决定了模型精度和计算效率；为了实现离散元数值模型计算既精确又高效，将数值模型关键研究区域(即研究对象所处的研究区域，平面尺寸为750m×250m)的栅格间距减小为25m，非关键研究区域(即边界条件的影响区域，平面尺寸为350m×250m)的栅格间距适当增大为50m，如图4所示。

若模型中包含多个岩层，则将同一位置离散点坐标数据中增加相应的坐标Z值，形成GRID_MODEL.DAT文件，内容如下：

其中，M表示X方向的离散点数，N表示Y方向的离散点数，I表示岩层面数，I层则有I+1个Z值(即：Z

步骤c：基于Fortran语言进行建模代码的编写和编译，实现栅格化离散点坐标数据生成为3DEC离散元数值模型命令流；

通过Fortran语言编写和编译建模代码，创建3DEC建模程序，如图6所示；采用3DEC建模程序将栅格化离散点坐标数据进行循环遍历和语言编写，生成由三棱柱块体单元组成的包含多个岩层的3DEC离散元数值模型命令流，建模程序运行成功的界面如图7所示，模型命令流形成3DEC_MODEL.DAT文件，内容如下：

其中，POLY表示多面体块体创建命令，REGION I表示岩层编号，PRISM表示棱柱命令，a表示棱柱下底面，b表示棱柱上底面，(x

步骤d：采用3DEC程序读取模型命令流运行生成3DEC离散元数值模型；

采用3DEC程序读取模型命令流，运行生成由三棱柱块体单元组合的初步3DEC离散元数值模型，如图9所示，通过3DEC程序中JOIN on命令可将三棱柱块体单元合并成一个整体块体，还可以通过3DEC程序中JSET、HIDE、SEEK等命令实现断层、结构面等特殊地质构造的构建，生成最终3DEC离散元数值模型，如图10所示。

在本发明的实施方式中研究区域地形图的平面尺寸为1100m×250m，若整个研究区域采用一种栅格间距和两种栅格间距，则对比结果如下：

若整个研究区域全采用25m的栅格间距，则两个岩层的3DEC离散元模型共880个块体单元；若整个研究区域采用两种栅格间距，即研究对象所处的研究区域(平面尺寸为750m×250m)采用25m的栅格间距，边界条件的影响区域(平面尺寸为350m×250m)采用50m的栅格间距，则两个岩层的3DEC离散元模型共740个块体单元。两种方法均能保证研究对象的模型精度为25m，但是采用两种栅格间距所建立3DEC模型的块体单元数量比采用一种栅格间距时少了20％左右，因此本发明中不同研究区域采用不同栅格间距的建模块体单元，能在保证研究对象建模精度条件下降低块体单元数量，提高3DEC模型的计算效率。

本发明的实施方式中，在研究对象所处研究区域(平面尺寸为750m×250m)的341个栅格化离散点插值数据中，采用克里金插值函数和线性插值函数及三次插值函数所得离散点高程数据之间的平均误差分别为0.11％和0.09％，最大误差分别为0.68％和0.66％。三个插值函数所得栅格化离散点插值数据近似一致，均能达到本发明的目的。

本发明的基于Fortran语言的3DEC离散元数值模型建模方法，具有计算高效、操作简单、灵活度高、便捷性好等优点，能够对地形复杂且范围较大的研究区域开展高精度、多岩层、可变块体单元尺寸的三维离散元数值模型构建。

上述实施例仅为本发明的较佳的实施方式，并非对本发明实施方式的限定。在上述实施例的基础上还可以做出其它不同形式的变化或变动，这里无法穷举所有实施方式。凡是属于本发明的技术方案所引申出的显而易见的变化与改动仍属于本发明的保护范围之内。