基于GMS软件的海砂矿体三维地质建模和资源量估算

习 龙，倪玉根，何 健，梁 开，陈科衡，夏 真，陈 梅，薛 峭

(1.中国地质调查局广州海洋地质调查局，广东广州 510075；2.南方海洋科学与工程广东省实验室(广州)，广东

广州 511458)

0 引言

砂是世界上开采和交易量最大的矿产资源之一，甚至远远超过化石燃料(Torres et al.，2017)。近年来，由于城市化和人口增长导致日本、美国，尤其是中国、印度和非洲海砂的需求激增，砂的使用量已远超过了其自然形成的速度(Sverdrup et al.，2017；Bendixen et al.，2019)。全球每年大约需消耗320～500亿吨的砂，其主要用于制造混凝土、玻璃、光学仪器、润滑剂和电子工业重要部件等(Bendixen et al.，2019)。砂资源紧缺已经成为关系我国乃至全球经济、社会发展和环境保护的共性问题(Torres et al.，2017；Bendixen et al.，2019)。在河砂开采被限制，产量供不应求的背景下，“向海要砂”，合理合规开发海砂，是解决砂源短缺，保障重大工程建设和社会发展的大势所趋(倪玉根等，2021)。由于需求量大而供应不足，海砂价格暴涨，珠三角地区海砂价格一度突破300元/立方米。因此，在沿海海域开展系统的海砂资源调查，合理规划海砂资源开发利用，显得尤为迫切。海砂资源量估算是海砂资源调查工作中最核心的内容之一，是海砂资源开发的基础，因而精准估算海砂资源量非常重要。

目前，在海砂地质勘查工作中多使用地质块段法来估算海砂资源量，即根据矿床地质特点将矿体划分为若干个地质块段，将它们看作理想的板状体，板状体的厚度为块段内所有见矿工程中矿体的平均厚度，面积为块段平面投影面积，体积为块段体积。该方法的优点是简明，易于计算，缺点是相对粗略，因为海砂矿体的实际形态并非是理想的平板一块，而是高低不平的形态，因此地质块段法计算出的结果可能与实际情况存在一定的误差。在陆域固体矿产、水资源评价行业中，为更加准确地估算矿产资源的资源量/储量，减少矿产开发成本，提高工作效率，规避开发风险，国内外常使用Gocad、Micromine、3DMine、Surpac和GMS(Groundwater Modeling System)等软件或方法来建立矿体的三维地质模型(Wu，2006；罗周全等，2010；余牛奔等，2015；时昆林等，2015；李晓晖等，2015；李杨，2017；朱鹏程，2017；俞嘉嘉等，2021)。三维地质建模的概念最早由Houlding提出，其对三维地质建模的流程和方法做了系统描述，为三维地质建模理论和方法奠定了基础(Houlding，1994)。经过近几十年的发展，国内外专家围绕基础地质、工程地质、矿产地质等学科问题，在三维建模技术方法和建模软件方面取得了很大进度。目前比较常见的三维地质建模方法包括(1)基于钻孔数据的建模方法(明镜，2012；Akiska，2013)；(2)基于剖面的建模方法(吴志春等，2016)，其可分为平行剖面建模法(王勇等，2003；Whiteaker et al.，2012；Miao et al.，2017)和网状交叉剖面建模法(屈红刚等，2008；郭艳军等，2009)；(3)基于多源交互建模方法(薛林福等，2014；Wang et al.，2015)。这些建模方法的地质数据通过地质调查、矿产勘探、地球物理探测等方法手段获取。GMS软件是由美国Brigham Young大学和美国军方在综合MODFLOW、MT3DMS和FEMWATER等已有地下水模块的基础上，共同开发用于地下水模拟的软件系统(贾瑞亮等，2014)。该软件中的Borehole(钻孔)、2D Scatter Point(二维散点)、TINS(三角不规则网)和Solid(实体)模块是建立三维地层模型的主要功能模块。由于GMS软件具有友好的操作界面、优良的三维可视化和处理功能，其在国内外广泛应用于地下水模拟和三维地层构建(陈正华等，2012；刘丽花和张树清，2015；Bayat et al.，2020)。如潘欢迎等(2010)以海口市水文地质钻孔资料为基础，利用GMS软件呈现了海口地下水系统结构的三维结构。本文以南海北部某海域的海砂矿区为例，在收集、解译综合物探和钻探地质资料的基础上，利用GMS软件建立矿区内海砂矿体的三维地质模型和估算海砂矿体体积资源量，并与利用地质块段法估算的该矿体体积资源量进行对比，以便更准确地估算海砂资源量，可更好地支撑海砂资源开发和管理。

1 矿区概况

矿区位于南海北部陆架区，该地区属于南亚热带，海区开阔平坦，气候受海洋调节，冬季受东北季风影响盛行东北风，且可受到强冷空气的影响；夏季受东南季风影响盛行东南风，常有热带气旋影响。

矿区地层揭露均为第四系，现今研究普遍认为该海域的地层年龄下限在40 ka B.P.左右，据Yim(1984)取得的岩心分析和同位素测年资料研究成果，矿区附近的第四纪沉积从下至上可划分为四层，即底层为砂砾层和粉砂质粘土层的冲积相层，14C测年距今33440±1740至36480±830年；次底层为海相层，其岩性主要为淡灰色粉砂质粘土，含少量细砂，其中含有Ammoniaannectens、Quinqueloculina等有孔虫，14C测年为距今28000至30000年；次顶层为冲积层，岩性主要为灰色及棕黄色砂砾层、粉砂及粘土质砂，14C测年为距今16420±660至27660±590年。顶层为全新世海侵层，钻孔约9 m深处的岩心14C测年距今7920±110至6640±100年，地质时代大致相当于中全新世。由此可知,本区在更新世晚期还是陆地，只是在中全新世时，由于大西洋期气候温暖，冰川大量融化，海面上升，本区的海岛才与大陆分离。

本矿区内共涉及钻孔32口，钻孔的平面位置分布图见图1，部分钻孔的岩性信息详见表1。

表1 钻孔岩性信息表

图1 矿区钻孔平面分布图

2 地质数据库构建

海砂矿区三维地质建模数据库的数据源主要包括地质钻孔资料和综合物探资料(多波束测深、浅地层剖面和单道地震测量资料)。多波束测深资料用于生成矿区的海底地形资料，以作为后续三维地质体上边界。浅地层剖面测量用于排除泥层，单道地震测量用于识别砂层，二者结合使用以确定矿体的边界。本文的目的在于呈现海砂矿体的三维地质模型以及估算矿区的海砂资源量，在地质数据库构建过程中把地层简化为砂层和泥层，其中砂层需同时满足如下条件：(1)粒径≥0.063 mm；(2)含砂量≥55%(粒径≥0.063 mm沉积物的质量比)；(3)厚度≥2 m。沉积物粒径小于0.063 mm，且厚度>1 m的为泥层，即夹石。

(1)矿体的边界确定

矿区共有两个海砂矿体，分别为V1和V2矿体，V1矿体为晚更新世晚期冲洪积砂，V2矿体为晚更新世中期冲洪积砂。根据物探(浅地层剖面和单道地震测量)和钻孔资料综合分析得知：V1矿体在矿区内分为西北部V1-1、中部V1-2和东部V1-3三个块段，面积分别为26.5 km2、220.7 km2和73.7 km2，总面积320.9 km2。V1-1块段平均含砂量68.78%，砂层平均厚度4.7 m；V1-2块段平均含砂量64.74%，砂层平均厚度 9.1 m；V1-3块段平均含砂量68.16%，砂层平均厚度7.6 m。V2矿体仅在矿区西南部分布，面积为163.5 km2，平均含砂量68.70%，砂层平均厚度 5.6 m。V1和V2矿体的分布情况见图2。

图2 V1和V2矿体分布图

续表1

(2)地形数据提取

为了得到更加准确的矿体三维地质模型，本文对野外采集的多波束测深资料进行解译，提取建模区的水深地形数据和钻孔孔口高程数据作为控制点，以生成矿区三维地质模型的上界面。

(3)标准层建立

根据矿区的地层特征以及GMS软件三维建模的特点，本研究矿区的标准层共5个岩性层(下文三维模型中的Materials)，共涉及3个泥层和2个砂层(即V1矿体砂层和V2矿体砂层)。

按照Borehole Date(钻孔模块)要求的数据格式对32个钻孔的SoilID(岩性号)、HGUID(层序号)和HorizonID(地层沉积序列号)从下至上进行编号。HorizonID编号是标准层建立的重点，HorizonID代表的是地质体中按照沉积顺序出现的地层单元的上界面。因此，当某一钻孔的最下面地层为建模区域内沉积最早的地层时，其HorizonID需从0开始编号，而SoilID和HGUID则从1开始编号。钻孔ZK1、ZK2……ZK10的岩性、层序和地层沉积序列数据见表2。从表2中可以看到，部分钻孔贯穿整个矿区的岩性层，部分存在岩性层缺失，如ZK3、ZK5和ZK8贯穿整个矿区的5个岩性层，ZK10仅存在顶部的5号岩性层。

表2 钻孔的岩性、层序和地层沉积序列信息

续表2

3 矿区三维地质模型建立

3.1 导入钻孔数据

将按照标准岩性层编制好的钻孔(表3)批量导入软件，给每一层赋予对应的材料(Materials)属性，如SoilID为1、3和5，其对应的Materials为泥层，即非矿层，分别对应于第1层泥、第2层泥和第3层泥(从下至上的顺序)，SoilID为2和4，对应的Materials为砂层，即矿层，分别为下部的砂层(属于V2矿体)和上部的砂层(属于V1矿体)。最终在矿区内生成的钻孔如图3所示。

表3 三维地质模型法和地质块段法资源量估算结果表Table 3 Resource estimation of three-dimensional geological model method and geological block method

图3 矿区内钻孔三维分布图

3.2 构建钻孔联孔剖面

GMS软件提供自动和手动两种方法来创建钻孔联孔剖面。本文先通过钻孔模块下的自动方式初步建立钻孔联孔剖面，然后利用物探和实际钻孔柱状图资料，根据岩性、地层上下关系手动修正钻孔联孔剖面，最后共确定钻孔联孔剖面76个，连接形成的栅状图见图4，建立的钻孔联孔剖面最终将用于建立三维地质模型。

图4 矿区钻孔联孔剖面栅状图

3.3 矿区三维地质体的上边界和平面范围确定

矿区三维地质体的边界包括三维地质体的上、下边界和平面上矿体范围。三维地质体的上边界确定可通过导入海底地形控制点的二维散点数据来实现。在工程中建一个与矿区形状相同的coverage(图层)文件，然后将此图层文件直接转换为TIN模块下的不规则三角网控制面，在此之后将2D散点数据高程Z值通过合适的插值方法插值到已激活的不规则三角网控制面。生成的矿区海底地形不规则三角网控制面见图5，可见矿区整体相对平坦，呈北高南低的趋势，矿区内高程差为18.9 m，坡度为0.945‰。

图5 矿区海底地形不规则三角网控制面

三维地质体的下边界宜采用钻孔的孔底高程来控制。

平面上矿体的范围是通过钻孔揭露的地层和物探解译结果圈定，本文通过钻孔和物探结果圈定的V1和V2矿体平面范围见图2。V1和V2矿体的面文件由ARCGIS软件生成，然后将面文件导入GMS软件转换为可约束后续三维建模的图层文件。特别需要注意的是，此处的图层文件必须得是在conceptual model(概念模型模块)下的图层文件，并且图层的HorizonID(地层沉积序列号)一定要赋值，且限制不同矿体边界的图层与其对应矿体的HorizonID一致，如用于限制V1矿体边界的图层，其HorizonID要与V1在表2中的HorizonID一致，HorizonID赋值为4。

3.4 生成三维地质体

三维地质体利用Borehole(钻孔模块)下的Horizons to Solids(地层沉积序列生成实体模型)来实现。利用矿区内的所有钻孔、钻孔联孔剖面和步骤(3)生成的不规则三角网控制面，选择Natural neighbor(自然邻点插值法)，构建矿区的三维地质体。本文中三维地质体的构建关键需利用不规则三角网控制面来限定三维模型的上边界以及平面范围，利用钻孔的底高程来控制模型的下边界，以确保构建的三维地质体不失真。生成的矿区三维地质体见图6。图7呈现V1矿体的空间展布情况，图8展现V2矿体的空间展布。图9呈现V1矿体的空间展布情况，图10呈现V2矿体的空间展布。

图6 构建的矿区三维地质体

图7 V1海砂矿体三维斜视图

图8 V2海砂矿体三维斜视图

为了更加清晰直观地展示和了解矿区内矿体的展布情况，可使用软件从任意角度、方向切割三维地质体模型生成剖面。如图9南北向切割形成剖面A-A’、 剖面B-B’和剖面C-C’，东西向切割形成剖面D-D’和剖面E-E’，清晰地呈现出海砂矿层、非矿层空间分布以及顶底板分布情况。三维地质体任意方向的剖面呈现有利于海砂勘查工作者和相关管理部门直观了解海砂矿体的地层结构。同时，对于增强海砂勘查研究人员对海砂地层结构的理解，以及深入研究海砂成矿机制具有重要意义。

图9 矿区三维地质体横切面图

4 矿区资源量估算结果

矿区的海砂资源量可直接通过GMS软件生成的三维地质模型读取。经统计，矿区内V1海砂矿体的体积资源量为26.610亿m3，V2海砂矿体的体积资源量为8.665亿m3，总资源量35.275亿m3。传统地质块段法估算的矿区内V1海砂矿体的体积资源量为26.858亿m3，V2海砂矿体的体积资源量为9.186亿m3，总资源量36.044亿m3。两种方法估算的海砂资源量差值为0.769亿m3，差值率为2.13%。

三维地质模型法和地质块段法资源量估算的结果对比发现，三维地质模型法相比传统的地质块段法能够清晰地表达矿体的形态、空间展布规律和分布范围，使得估算结果更加可靠，更符合实际。该方法能更加快捷直观三维可视化展示矿体，便于研究者识别矿体潜在的延伸方向。除此之外，还可以指导海砂矿体的后期开发，因此可以将三维地质模型法推广使用。

矿区三维地质模型的构建很大程度上依赖于钻孔数据，钻孔数据能真实、准确、详细和直观地表达地层单元的原始信息，在三维地质建模过程中发挥不可替代的作用，钻孔数量越多，钻孔数据信息越准确，构建的三维地质模型越符合实际。受水深、水流和钻探船舶等限制，海上钻探施工相较于陆地要难得多，钻孔数据获取成本高，而且钻孔控制的地层信息仅在一定范围内有效(贺怀建等，2002；刘丽花和张树清，2015)，其表达的地层结构具有局限性，因此在实际钻孔稀疏的情况下，引入由实际钻孔和物探资料综合解译得到的虚拟钻孔会使得建立的矿区三维地质模型会更加符合实际(朱良峰等，2004)。

5 结论

(1)利用综合物探(浅地层剖面和单道地震测量)和已有实际钻孔资料，圈定了V1和V2矿体的分布范围，提高了矿区地层信息的精确度，使得建立的矿区三维地层模型会更加符合实际。

(2)基于GMS软件系统构建了数据库，综合利用钻孔、物探资料，剖面图和海底地形控制面等，建立了三维地质模型，有效实现了海砂矿体三维可视化，海砂矿体体积资源量的估算。三维地质模型法与传统地质块段法估算的海砂体积资源量相差0.769亿m3，差值率为2.13%。GMS软件系统构建海砂矿三维地质模型和资源量估算具有一定的推广意义。