基于GMS軟件的海砂礦體三維地質(zhì)建模和資源量估算

習(xí) 龍，倪玉根，何 健，梁 開，陳科衡，夏 真，陳 梅，薛 峭

(1.中國地質(zhì)調(diào)查局廣州海洋地質(zhì)調(diào)查局，廣東廣州 510075；2.南方海洋科學(xué)與工程廣東省實驗室(廣州)，廣東

廣州 511458)

0 引言

砂是世界上開采和交易量最大的礦產(chǎn)資源之一，甚至遠遠超過化石燃料(Torres et al.，2017)。近年來，由于城市化和人口增長導(dǎo)致日本、美國，尤其是中國、印度和非洲海砂的需求激增，砂的使用量已遠超過了其自然形成的速度(Sverdrup et al.，2017；Bendixen et al.，2019)。全球每年大約需消耗320～500億噸的砂，其主要用于制造混凝土、玻璃、光學(xué)儀器、潤滑劑和電子工業(yè)重要部件等(Bendixen et al.，2019)。砂資源緊缺已經(jīng)成為關(guān)系我國乃至全球經(jīng)濟、社會發(fā)展和環(huán)境保護的共性問題(Torres et al.，2017；Bendixen et al.，2019)。在河砂開采被限制，產(chǎn)量供不應(yīng)求的背景下，“向海要砂”，合理合規(guī)開發(fā)海砂，是解決砂源短缺，保障重大工程建設(shè)和社會發(fā)展的大勢所趨(倪玉根等，2021)。由于需求量大而供應(yīng)不足，海砂價格暴漲，珠三角地區(qū)海砂價格一度突破300元/立方米。因此，在沿海海域開展系統(tǒng)的海砂資源調(diào)查，合理規(guī)劃海砂資源開發(fā)利用，顯得尤為迫切。海砂資源量估算是海砂資源調(diào)查工作中最核心的內(nèi)容之一，是海砂資源開發(fā)的基礎(chǔ)，因而精準(zhǔn)估算海砂資源量非常重要。

目前，在海砂地質(zhì)勘查工作中多使用地質(zhì)塊段法來估算海砂資源量，即根據(jù)礦床地質(zhì)特點將礦體劃分為若干個地質(zhì)塊段，將它們看作理想的板狀體，板狀體的厚度為塊段內(nèi)所有見礦工程中礦體的平均厚度，面積為塊段平面投影面積，體積為塊段體積。該方法的優(yōu)點是簡明，易于計算，缺點是相對粗略，因為海砂礦體的實際形態(tài)并非是理想的平板一塊，而是高低不平的形態(tài)，因此地質(zhì)塊段法計算出的結(jié)果可能與實際情況存在一定的誤差。在陸域固體礦產(chǎn)、水資源評價行業(yè)中，為更加準(zhǔn)確地估算礦產(chǎn)資源的資源量/儲量，減少礦產(chǎn)開發(fā)成本，提高工作效率，規(guī)避開發(fā)風(fēng)險，國內(nèi)外常使用Gocad、Micromine、3DMine、Surpac和GMS(Groundwater Modeling System)等軟件或方法來建立礦體的三維地質(zhì)模型(Wu，2006；羅周全等，2010；余牛奔等，2015；時昆林等，2015；李曉暉等，2015；李楊，2017；朱鵬程，2017；俞嘉嘉等，2021)。三維地質(zhì)建模的概念最早由Houlding提出，其對三維地質(zhì)建模的流程和方法做了系統(tǒng)描述，為三維地質(zhì)建模理論和方法奠定了基礎(chǔ)(Houlding，1994)。經(jīng)過近幾十年的發(fā)展，國內(nèi)外專家圍繞基礎(chǔ)地質(zhì)、工程地質(zhì)、礦產(chǎn)地質(zhì)等學(xué)科問題，在三維建模技術(shù)方法和建模軟件方面取得了很大進度。目前比較常見的三維地質(zhì)建模方法包括(1)基于鉆孔數(shù)據(jù)的建模方法(明鏡，2012；Akiska，2013)；(2)基于剖面的建模方法(吳志春等，2016)，其可分為平行剖面建模法(王勇等，2003；Whiteaker et al.，2012；Miao et al.，2017)和網(wǎng)狀交叉剖面建模法(屈紅剛等，2008；郭艷軍等，2009)；(3)基于多源交互建模方法(薛林福等，2014；Wang et al.，2015)。這些建模方法的地質(zhì)數(shù)據(jù)通過地質(zhì)調(diào)查、礦產(chǎn)勘探、地球物理探測等方法手段獲取。GMS軟件是由美國Brigham Young大學(xué)和美國軍方在綜合MODFLOW、MT3DMS和FEMWATER等已有地下水模塊的基礎(chǔ)上，共同開發(fā)用于地下水模擬的軟件系統(tǒng)(賈瑞亮等，2014)。該軟件中的Borehole(鉆孔)、2D Scatter Point(二維散點)、TINS(三角不規(guī)則網(wǎng))和Solid(實體)模塊是建立三維地層模型的主要功能模塊。由于GMS軟件具有友好的操作界面、優(yōu)良的三維可視化和處理功能，其在國內(nèi)外廣泛應(yīng)用于地下水模擬和三維地層構(gòu)建(陳正華等，2012；劉麗花和張樹清，2015；Bayat et al.，2020)。如潘歡迎等(2010)以?？谑兴牡刭|(zhì)鉆孔資料為基礎(chǔ)，利用GMS軟件呈現(xiàn)了?？诘叵滤到y(tǒng)結(jié)構(gòu)的三維結(jié)構(gòu)。本文以南海北部某海域的海砂礦區(qū)為例，在收集、解譯綜合物探和鉆探地質(zhì)資料的基礎(chǔ)上，利用GMS軟件建立礦區(qū)內(nèi)海砂礦體的三維地質(zhì)模型和估算海砂礦體體積資源量，并與利用地質(zhì)塊段法估算的該礦體體積資源量進行對比，以便更準(zhǔn)確地估算海砂資源量，可更好地支撐海砂資源開發(fā)和管理。

1 礦區(qū)概況

礦區(qū)位于南海北部陸架區(qū)，該地區(qū)屬于南亞熱帶，海區(qū)開闊平坦，氣候受海洋調(diào)節(jié)，冬季受東北季風(fēng)影響盛行東北風(fēng)，且可受到強冷空氣的影響；夏季受東南季風(fēng)影響盛行東南風(fēng)，常有熱帶氣旋影響。

礦區(qū)地層揭露均為第四系，現(xiàn)今研究普遍認為該海域的地層年齡下限在40 ka B.P.左右，據(jù)Yim(1984)取得的巖心分析和同位素測年資料研究成果，礦區(qū)附近的第四紀沉積從下至上可劃分為四層，即底層為砂礫層和粉砂質(zhì)粘土層的沖積相層，14C測年距今33440±1740至36480±830年；次底層為海相層，其巖性主要為淡灰色粉砂質(zhì)粘土，含少量細砂，其中含有Ammoniaannectens、Quinqueloculina等有孔蟲，14C測年為距今28000至30000年；次頂層為沖積層，巖性主要為灰色及棕黃色砂礫層、粉砂及粘土質(zhì)砂，14C測年為距今16420±660至27660±590年。頂層為全新世海侵層，鉆孔約9 m深處的巖心14C測年距今7920±110至6640±100年，地質(zhì)時代大致相當(dāng)于中全新世。由此可知,本區(qū)在更新世晚期還是陸地，只是在中全新世時，由于大西洋期氣候溫暖，冰川大量融化，海面上升，本區(qū)的海島才與大陸分離。

本礦區(qū)內(nèi)共涉及鉆孔32口，鉆孔的平面位置分布圖見圖1，部分鉆孔的巖性信息詳見表1。

表1 鉆孔巖性信息表

圖1 礦區(qū)鉆孔平面分布圖

2 地質(zhì)數(shù)據(jù)庫構(gòu)建

海砂礦區(qū)三維地質(zhì)建模數(shù)據(jù)庫的數(shù)據(jù)源主要包括地質(zhì)鉆孔資料和綜合物探資料(多波束測深、淺地層剖面和單道地震測量資料)。多波束測深資料用于生成礦區(qū)的海底地形資料，以作為后續(xù)三維地質(zhì)體上邊界。淺地層剖面測量用于排除泥層，單道地震測量用于識別砂層，二者結(jié)合使用以確定礦體的邊界。本文的目的在于呈現(xiàn)海砂礦體的三維地質(zhì)模型以及估算礦區(qū)的海砂資源量，在地質(zhì)數(shù)據(jù)庫構(gòu)建過程中把地層簡化為砂層和泥層，其中砂層需同時滿足如下條件：(1)粒徑≥0.063 mm；(2)含砂量≥55%(粒徑≥0.063 mm沉積物的質(zhì)量比)；(3)厚度≥2 m。沉積物粒徑小于0.063 mm，且厚度>1 m的為泥層，即夾石。

(1)礦體的邊界確定

礦區(qū)共有兩個海砂礦體，分別為V1和V2礦體，V1礦體為晚更新世晚期沖洪積砂，V2礦體為晚更新世中期沖洪積砂。根據(jù)物探(淺地層剖面和單道地震測量)和鉆孔資料綜合分析得知：V1礦體在礦區(qū)內(nèi)分為西北部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礦體僅在礦區(qū)西南部分布，面積為163.5 km2，平均含砂量68.70%，砂層平均厚度 5.6 m。V1和V2礦體的分布情況見圖2。

圖2 V1和V2礦體分布圖

續(xù)表1

(2)地形數(shù)據(jù)提取

為了得到更加準(zhǔn)確的礦體三維地質(zhì)模型，本文對野外采集的多波束測深資料進行解譯，提取建模區(qū)的水深地形數(shù)據(jù)和鉆孔孔口高程數(shù)據(jù)作為控制點，以生成礦區(qū)三維地質(zhì)模型的上界面。

(3)標(biāo)準(zhǔn)層建立

根據(jù)礦區(qū)的地層特征以及GMS軟件三維建模的特點，本研究礦區(qū)的標(biāo)準(zhǔn)層共5個巖性層(下文三維模型中的Materials)，共涉及3個泥層和2個砂層(即V1礦體砂層和V2礦體砂層)。

按照Borehole Date(鉆孔模塊)要求的數(shù)據(jù)格式對32個鉆孔的SoilID(巖性號)、HGUID(層序號)和HorizonID(地層沉積序列號)從下至上進行編號。HorizonID編號是標(biāo)準(zhǔn)層建立的重點，HorizonID代表的是地質(zhì)體中按照沉積順序出現(xiàn)的地層單元的上界面。因此，當(dāng)某一鉆孔的最下面地層為建模區(qū)域內(nèi)沉積最早的地層時，其HorizonID需從0開始編號，而SoilID和HGUID則從1開始編號。鉆孔ZK1、ZK2……ZK10的巖性、層序和地層沉積序列數(shù)據(jù)見表2。從表2中可以看到，部分鉆孔貫穿整個礦區(qū)的巖性層，部分存在巖性層缺失，如ZK3、ZK5和ZK8貫穿整個礦區(qū)的5個巖性層，ZK10僅存在頂部的5號巖性層。

表2 鉆孔的巖性、層序和地層沉積序列信息

續(xù)表2

3 礦區(qū)三維地質(zhì)模型建立

3.1 導(dǎo)入鉆孔數(shù)據(jù)

將按照標(biāo)準(zhǔn)巖性層編制好的鉆孔(表3)批量導(dǎo)入軟件，給每一層賦予對應(yīng)的材料(Materials)屬性，如SoilID為1、3和5，其對應(yīng)的Materials為泥層，即非礦層，分別對應(yīng)于第1層泥、第2層泥和第3層泥(從下至上的順序)，SoilID為2和4，對應(yīng)的Materials為砂層，即礦層，分別為下部的砂層(屬于V2礦體)和上部的砂層(屬于V1礦體)。最終在礦區(qū)內(nèi)生成的鉆孔如圖3所示。

表3 三維地質(zhì)模型法和地質(zhì)塊段法資源量估算結(jié)果表Table 3 Resource estimation of three-dimensional geological model method and geological block method

圖3 礦區(qū)內(nèi)鉆孔三維分布圖

3.2 構(gòu)建鉆孔聯(lián)孔剖面

GMS軟件提供自動和手動兩種方法來創(chuàng)建鉆孔聯(lián)孔剖面。本文先通過鉆孔模塊下的自動方式初步建立鉆孔聯(lián)孔剖面，然后利用物探和實際鉆孔柱狀圖資料，根據(jù)巖性、地層上下關(guān)系手動修正鉆孔聯(lián)孔剖面，最后共確定鉆孔聯(lián)孔剖面76個，連接形成的柵狀圖見圖4，建立的鉆孔聯(lián)孔剖面最終將用于建立三維地質(zhì)模型。

圖4 礦區(qū)鉆孔聯(lián)孔剖面柵狀圖

3.3 礦區(qū)三維地質(zhì)體的上邊界和平面范圍確定

礦區(qū)三維地質(zhì)體的邊界包括三維地質(zhì)體的上、下邊界和平面上礦體范圍。三維地質(zhì)體的上邊界確定可通過導(dǎo)入海底地形控制點的二維散點數(shù)據(jù)來實現(xiàn)。在工程中建一個與礦區(qū)形狀相同的coverage(圖層)文件，然后將此圖層文件直接轉(zhuǎn)換為TIN模塊下的不規(guī)則三角網(wǎng)控制面，在此之后將2D散點數(shù)據(jù)高程Z值通過合適的插值方法插值到已激活的不規(guī)則三角網(wǎng)控制面。生成的礦區(qū)海底地形不規(guī)則三角網(wǎng)控制面見圖5，可見礦區(qū)整體相對平坦，呈北高南低的趨勢，礦區(qū)內(nèi)高程差為18.9 m，坡度為0.945‰。

圖5 礦區(qū)海底地形不規(guī)則三角網(wǎng)控制面

三維地質(zhì)體的下邊界宜采用鉆孔的孔底高程來控制。

平面上礦體的范圍是通過鉆孔揭露的地層和物探解譯結(jié)果圈定，本文通過鉆孔和物探結(jié)果圈定的V1和V2礦體平面范圍見圖2。V1和V2礦體的面文件由ARCGIS軟件生成，然后將面文件導(dǎo)入GMS軟件轉(zhuǎn)換為可約束后續(xù)三維建模的圖層文件。特別需要注意的是，此處的圖層文件必須得是在conceptual model(概念模型模塊)下的圖層文件，并且圖層的HorizonID(地層沉積序列號)一定要賦值，且限制不同礦體邊界的圖層與其對應(yīng)礦體的HorizonID一致，如用于限制V1礦體邊界的圖層，其HorizonID要與V1在表2中的HorizonID一致，HorizonID賦值為4。

3.4 生成三維地質(zhì)體

三維地質(zhì)體利用Borehole(鉆孔模塊)下的Horizons to Solids(地層沉積序列生成實體模型)來實現(xiàn)。利用礦區(qū)內(nèi)的所有鉆孔、鉆孔聯(lián)孔剖面和步驟(3)生成的不規(guī)則三角網(wǎng)控制面，選擇Natural neighbor(自然鄰點插值法)，構(gòu)建礦區(qū)的三維地質(zhì)體。本文中三維地質(zhì)體的構(gòu)建關(guān)鍵需利用不規(guī)則三角網(wǎng)控制面來限定三維模型的上邊界以及平面范圍，利用鉆孔的底高程來控制模型的下邊界，以確保構(gòu)建的三維地質(zhì)體不失真。生成的礦區(qū)三維地質(zhì)體見圖6。圖7呈現(xiàn)V1礦體的空間展布情況，圖8展現(xiàn)V2礦體的空間展布。圖9呈現(xiàn)V1礦體的空間展布情況，圖10呈現(xiàn)V2礦體的空間展布。

圖6 構(gòu)建的礦區(qū)三維地質(zhì)體

圖7 V1海砂礦體三維斜視圖

圖8 V2海砂礦體三維斜視圖

為了更加清晰直觀地展示和了解礦區(qū)內(nèi)礦體的展布情況，可使用軟件從任意角度、方向切割三維地質(zhì)體模型生成剖面。如圖9南北向切割形成剖面A-A’、 剖面B-B’和剖面C-C’，東西向切割形成剖面D-D’和剖面E-E’，清晰地呈現(xiàn)出海砂礦層、非礦層空間分布以及頂?shù)装宸植记闆r。三維地質(zhì)體任意方向的剖面呈現(xiàn)有利于海砂勘查工作者和相關(guān)管理部門直觀了解海砂礦體的地層結(jié)構(gòu)。同時，對于增強海砂勘查研究人員對海砂地層結(jié)構(gòu)的理解，以及深入研究海砂成礦機制具有重要意義。

圖9 礦區(qū)三維地質(zhì)體橫切面圖

4 礦區(qū)資源量估算結(jié)果

礦區(qū)的海砂資源量可直接通過GMS軟件生成的三維地質(zhì)模型讀取。經(jīng)統(tǒng)計，礦區(qū)內(nèi)V1海砂礦體的體積資源量為26.610億m3，V2海砂礦體的體積資源量為8.665億m3，總資源量35.275億m3。傳統(tǒng)地質(zhì)塊段法估算的礦區(qū)內(nèi)V1海砂礦體的體積資源量為26.858億m3，V2海砂礦體的體積資源量為9.186億m3，總資源量36.044億m3。兩種方法估算的海砂資源量差值為0.769億m3，差值率為2.13%。

三維地質(zhì)模型法和地質(zhì)塊段法資源量估算的結(jié)果對比發(fā)現(xiàn)，三維地質(zhì)模型法相比傳統(tǒng)的地質(zhì)塊段法能夠清晰地表達礦體的形態(tài)、空間展布規(guī)律和分布范圍，使得估算結(jié)果更加可靠，更符合實際。該方法能更加快捷直觀三維可視化展示礦體，便于研究者識別礦體潛在的延伸方向。除此之外，還可以指導(dǎo)海砂礦體的后期開發(fā)，因此可以將三維地質(zhì)模型法推廣使用。

礦區(qū)三維地質(zhì)模型的構(gòu)建很大程度上依賴于鉆孔數(shù)據(jù)，鉆孔數(shù)據(jù)能真實、準(zhǔn)確、詳細和直觀地表達地層單元的原始信息，在三維地質(zhì)建模過程中發(fā)揮不可替代的作用，鉆孔數(shù)量越多，鉆孔數(shù)據(jù)信息越準(zhǔn)確，構(gòu)建的三維地質(zhì)模型越符合實際。受水深、水流和鉆探船舶等限制，海上鉆探施工相較于陸地要難得多，鉆孔數(shù)據(jù)獲取成本高，而且鉆孔控制的地層信息僅在一定范圍內(nèi)有效(賀懷建等，2002；劉麗花和張樹清，2015)，其表達的地層結(jié)構(gòu)具有局限性，因此在實際鉆孔稀疏的情況下，引入由實際鉆孔和物探資料綜合解譯得到的虛擬鉆孔會使得建立的礦區(qū)三維地質(zhì)模型會更加符合實際(朱良峰等，2004)。

5 結(jié)論

(1)利用綜合物探(淺地層剖面和單道地震測量)和已有實際鉆孔資料，圈定了V1和V2礦體的分布范圍，提高了礦區(qū)地層信息的精確度，使得建立的礦區(qū)三維地層模型會更加符合實際。

(2)基于GMS軟件系統(tǒng)構(gòu)建了數(shù)據(jù)庫，綜合利用鉆孔、物探資料，剖面圖和海底地形控制面等，建立了三維地質(zhì)模型，有效實現(xiàn)了海砂礦體三維可視化，海砂礦體體積資源量的估算。三維地質(zhì)模型法與傳統(tǒng)地質(zhì)塊段法估算的海砂體積資源量相差0.769億m3，差值率為2.13%。GMS軟件系統(tǒng)構(gòu)建海砂礦三維地質(zhì)模型和資源量估算具有一定的推廣意義。