利用應力場預測熱液區域
——以TAG區為例

陳欽柱, 陶春輝, 廖時理, 李懷明, 鄧顯明

(1. 溫州科技職業學院 園林與水利工程學院，浙江 溫州325000；2. 國家海洋局第二海洋研究所 海底科學重點實驗室，浙江 杭州310012)

?

利用應力場預測熱液區域
——以TAG區為例

陳欽柱1, 陶春輝2*, 廖時理2, 李懷明2, 鄧顯明2

(1. 溫州科技職業學院 園林與水利工程學院，浙江 溫州325000；2. 國家海洋局第二海洋研究所 海底科學重點實驗室，浙江 杭州310012)

本文根據TAG區的鉆探資料及巖心測試結果，建立了雙層地質模型,在此基礎上利用ANSYS應力軟件并結合TAG熱液區的地形數據對該區進行應力模擬。結果表明:熱液噴口區域與最大水平應力低值區有較好的對應關系。其中仍處于活動狀態的TAG丘狀體區呈現明顯的局部最大水平應力低值;而已經停止活動并且不具有典型噴口地形的MIR丘狀體區域則處于最大水平應力的非封閉低值區域。據此，本文在TAG丘狀體區域圈定了5個噴口可能區域，鉆探結果揭示區內存在較好的礦化和蝕變現象，表明應力場預測法可能是一種有效的成礦預測方法。

海底地形；TAG熱液區；應力場；成礦預測

1 引言

地應力是地質力學與巖體力學研究的基本內容之一。現今對地應力的研究主要集中在巖層中各部位的應力狀態的模擬及推測。通過分析應力分布，可以為工程建筑區域穩定性提供依據，同時還能由應力場的分布推測可能產生破裂的區域[1]。Sch?pa等[2]通過分析巖體自重引起的應力場推測武爾長諾島上熱噴口的位置。

地應力形成的原因十分復雜，至今仍是未完全解決的問題。地應力的形成主要與地球的各種動力運動有關，溫度、水壓梯度也可引起相應的應力場。其中，構造應力場和巖體自重應力場是現今地應力場的主要組成部分[3]。陸地上對應力場的研究大多數是基于實測應力數據，通過有限元數學模型回歸分析或者邊界載荷調整等方法尋找最優的邊界條件，使模型在應力實測點的應力狀態與實測應力值最接近，從而模擬得到整個研究區的應力狀態[4—6]。但在有些情況下，也會在只考慮巖體自重條件下計算整個研究區域的應力場[2]。

現有關于海底應力場的研究主要集中在洋中脊形成過程中的應力場變化。Buck 等[7]通過考慮巖層密度和局部構造對應力的影響來研究洋中脊地形形成的原因及過程。然而，隨著陸地上礦產資源的不斷被消耗利用，海底賦存的資源日漸受到重視，其中海底多金屬硫化物礦產以其巨大的潛在價值尤其受到科學家的關注[8]。在此背景下，科學家希望通過研究海底熱液區應力場的特點為圈定海底多金屬硫化物富集區提供一些依據。但是由于海上作業成本及難度的原因，至今還難以如陸地上一樣獲得實測點的應力數據。因此，海底熱液區的應力場研究一直是一個難點。最近，在僅有多波束海底地形數據的情況下， Petukhov等[9]嘗試通過對不同熱液區的海底應力場進行模擬，最后總結發現海底切應力場與熱液噴口存在一定的關系。本文在此基礎上，融合了Sch?pa等[2]利用應力場尋找陸地熱噴口的理論，提出利用最大水平應力場尋找熱液噴口可能存在區域的方法。

2 方法介紹

為了盡可能接近真實地形的應力場的模擬，本次研究主要采取有限元的方法求解目標區的應力場。有限元應力求解法具有處理復雜多變模型的能力，適合處理不規則的地形模型[4]。現今市場上存在多種可選的有限元應力分析軟件，比如ANSYS和FLAC 3D，ANSYS在靜力分析中具有支持大應變、大變形且可選分析單元種類豐富的特點。本文主要采用ANSYS軟件分析巖體受力變形情況。具體模擬分析過程分以下3個步驟：

首先，利用FLAC 3D和surfer軟件建立符合真實海底地形的地質實體模型。針對大西洋的巖性特征，模型采用雙層地質體結構，并根據大西洋TAG區鉆探結果確定相應巖層的楊氏模量、泊松比和密度等物性參數[10]。

然后，根據實際需要確定相應的模型邊界條件。據前人研究表明，區域構造穩定地區，淺層地應力主要受到起伏地形影響[11]。考慮模型邊界條件的時候，力求不對目標區產生太大的邊界效應。因此，在只考慮巖體自重引起的應力條件下，為了使下底面邊界不對海底目標區域的應力結果產生太大影響，模型中相對加大下底面邊界的截止深度。對于其他幾個面的邊界條件的處理方式具體如下：固定模型前后和左右邊界面的x，y方向的位移；底面邊界則施加x,y,z三方向邊界約束；上表面則不做任何約束，使其自由活動[12]；施加的重力加速度g，使用以下計算公式獲得[13]:

g=9.780 49(1+0.005 288sin2φ

-0.000 006sin22φ-0.000 308 6Z) m/s2，

其中，φ為緯度值，Z為海拔高度值，獲得TAG區重力加速度為9.802 m/s2。

最后，通過ANSYS應力分析軟件計算海底應力場并分析結果。與Petukhov利用模擬的切應力場結果研究熱液噴口的位置相比，本文認為最大水平主應力場更能體現一個區域的應力場的特點，而且更易于結合陸地上尋找熱噴口的相關理論。因此，在結果討論部分主要研究最大水平主應力場與噴口位置的關系。

3 地質背景

Trans-Atlantic Geo-Traverse(TAG)熱液區是在大西洋中脊最早發現的大型熱液活動區，該區位于慢速擴張的中大西洋洋中脊26.08°N附近，靠近洋中脊裂谷東壁的谷底中[14](圖1)。ODP158航段的鉆探結果表明，TAG熱液區的基底巖石主要由玄武巖枕狀熔巖組成，這些巖石通常遭受了不同程度的熱液蝕變作用[15]。據估計，該區的熱液活動歷史有20 000年之久，并具有多期次熱液活動的特征[15]。

圖1 TAG區域地形圖，其中兩個方框的區域為應力模擬模型的范圍Fig.1 Topographic map of TAG hydrothermal field, the boxes show modeling regions

TAG熱液區的面積大約為25 km2，水深在2 300～4 000 m之間，在此熱液區上發現3處熱液帶：兩個非活動的熱液丘狀體殘留區MIR(26°08.7′N，44°48.4′W)和AIVIN區(26°09.54′～26°10.62′N，44°48.89′～44°48.50′W)以及現今仍然處于活動期的TAG丘狀體 (26°08.21′N, 44°49.57′W)[16]。此外，根據大西洋洋中脊深鉆巖石測試結果，推斷在地層表面至150 m范圍內為固結玄武巖以及礫狀沉積物的混合物，150 m以下則主要是固結玄武巖[10]。

4 應力模擬及結果

4.1 模型的建立

文中采用的多波束地形數據主要是通過geomap軟件獲取的公開數據，其中包括TAG熱液丘狀體的局部高精度地形圖以及Tivey等[17]在TAG區獲取的高精度多波束數據。借助于sufer的強大網格劃分能力和flac5.0對數據整合能力，實現了多波束數據和ANSYS實體模型的銜接。

因為ALVIN殘留熱液區的分布面積較分散，不利于體現熱液集中區域的應力場的特點，并且綜合考慮計算機的計算能力和模型網格劃分精度，本文建立了兩個實體模型，其中一個為3 500 m×2 340 m的TAG熱液區小比例尺模型(圖2)，另一個則為230 m×265 m的局部大比例尺的TAG丘狀體模型(圖3)。兩個模型的地形都是根據網格化的多波束數據確定，并把經緯度坐標轉換成相應的距離坐標。根據大西洋洋中脊的巖層特性，兩個模型均為雙層結構，同時按照Deep-Sea Drilling Project(DSDP)第37航段的巖心測量結果，賦予兩層地層相應地層物性參數[10](詳見表1)。

表1 巖石物性參數(據Hyndman和Drury [10]修改)

圖2 TAG熱液區小比例尺模型Fig.2 The small scale model of TAG hydrothermal field

圖3 TAG丘狀體大比例尺模型Fig.3 The large scale model of TAG doom model

4.2 計算結果

在計算的過程中，模型統一采用solid185單元類型，使用ANSYS對模型進行自由網格劃分，劃分的結果見圖2和圖3。其中TAG區模型共劃分了2 526 983個單元和453 029個節點，TAG丘狀體模型共劃分了4 237 477單元和848 739個節點。在靜態應力條件下，采用Pre-Condicton (CG)求解器對節點結果進行計算。通過對比節點z方向的位移結果和模型區域的地形(圖4，圖5)，可以看出在地勢較高的部位節點z方向位移值較大，在地勢低洼的部位節點z方向位移值較小，這一結果符合重力變形基本理論。且在模型的底面，z方向位移為0，結果和最初施加的邊界條件約束相吻合。該模型z方向位移結果的數量級和前人重力模擬位移結果的相對數量級相似，可知本次模擬結果是可信的[18]。

圖4 小比例尺模型z方向位移結果圖Fig.4 The z direction displacement of small scale model

圖5 TAG丘狀體模型z方向位移結果Fig.5 The z direction displacement of TAG doom model

5 分析討論

海底多金屬硫化物是海底熱液活動的產物。海底熱液循環系統的形成與海底斷裂和裂隙等密不可分。海底裂隙為熱液循環系統提供了熱液流通和礦物質交換沉淀的場所，為熱液系統循環的產生提供了條件。因此，可以通過尋找裂隙的可能區域，來間接尋找海底多金屬硫化物。

早在1951年 Anderson就發現巖石中的3個主應力提供了水壓破裂方向的相關信息，即裂隙將沿著最大主應力的方向發育，同時裂隙垂直于最小主應力方向[19]。 Sch?pa等在研究陸地火山熱噴口與巖體自重引起的應力的關系時發現高溫熱噴口主要出現在水平最大主應力低值區[2]。因此研究水平和垂直應力的大小和方向有助于分析流體的運動趨勢，同時也能為推測海底熱液噴口位置提供一定依據。

據此，提取兩個模型海底面的最大水平主應力結果信息，分析最大水平主應力與噴口斷裂的關系。在圖6小比例尺模型結果圖中，最小應力值為0.1 MPa，最大應力值為1.48 MPa，符合巖體重力應力模擬的數量級。圖中可以看出主應力的高值主要出現在地形變化較大的地方，應力高值相對集中在所選區域的中部和中東部位，其中中部高值呈明顯的南北向分布。此外，還能看到應力高值區域附近的相對應力低值主要分布在山脊一側，這一現象和Sch?pa等的發現一致。另外，TAG丘狀體正好處在應力低值區，且在一個小范圍內被應力高值包圍，這也和Sch?pa等文章中地熱溢出點多處于最大水平應力低值區相符合。由于MIR區是非活動的殘留熱液區，其地形經過后期的改造已經不具有典型的凸起的似火山地形特征，而是處于一塊地形變化較大的山坡上。相應的在最大水平應力值圖中其應力值處于應力高值與應力低值的交界處，且靠近應力低值區。總的來說，小比例尺的應力模擬結果，再一次確認了噴口易出現在最大水平應力低值的觀點，也就是說，裂隙趨向于阻力小的方向發育。這一經驗性結論，可以為未來在較大范圍內圈定海底熱液硫化物可能區域提供一定依據。

圖6 小比例尺模型最大水平應力結果等值線圖Fig.6 The maximum horizontal stress result of the small scale model

圖7 TAG丘狀體最大水平主應力模擬結果圖Fig.7 The maximum horizontal stress result of TAG doom model

對于大比例尺的TAG丘狀體區域的模擬結果(圖7)，我們采用最大水平應力低值來預測丘狀體的熱液易噴發區域。在應力結果圖中，根據鉆探資料標出煙囪及鉆孔的位置，此外，根據Humphris[20]研究結果表明TAG-1、TAG-2、TAG-5 3個區域富含大量的硬石膏。從結合了3D地形效果的應力結果圖中可以看出，TAG丘狀體出現了4個最大水平應力低值區，這些相應的低值區可以認定為是熱液噴發易發生區域。首先，在TAG丘狀體頂端出現小范圍的低值區，此處的低值與熱液容易在地形最高點噴發的事實相一致，同時據已有資料顯示該區域內有煙囪體分布，表明現今該區域還有熱液噴出現象[20]。此外，按地形從高到低往下的第二層與第三層地勢較平緩區域也分別出現多處低值區，且出現的低值區呈現環帶狀分布并與山脊位置一致；據前人鉆探研究表明在圈定的應力低值區域內發現大量的硬石膏等礦物，表明相應區域曾經可能經歷較強的熱液活動。最底層由于地勢變化小，出現最大主應力為低值，這一結果與理論相符。

6 結論

本文利用海底地形數據模擬TAG熱液區自重應力場，從應力角度分析熱液噴口與地應力關系，較簡單地形預測方法能更綜合反映地形及巖性信息。小比例尺模擬結果表明，現今仍處于活動狀態的火山狀TAG丘狀體對應于最大水平主應力局部低值區，且在丘狀體周圍區域出現小范圍的應力高值區。而對于已經停止熱液活動且并不具有典型的火山地形的MIR區，在應力結果圖中依然對應最大水平應力低值區，和TAG丘狀體區域相比只是沒有相應的應力高值包圍。此外，根據局部的TAG丘狀體區域應力模擬結果，推測了5處熱液活動易發區域，并且與已有的鉆孔資料進行對比，發現最大水平應力低值區與熱液產物富集區存在一定關系。因此，該發現可以為今后進一步縮小熱液噴口尋找范圍提供參考依據。

致謝：本研究在模型模擬及分析討論過程中得到全俄地質與世界海洋礦產資源科學研究院的S.I. Petukhov研究員的寶貴建議，在此深表謝意。

[1] Martin C D, Kaiser P K, Christiansson R. Stress, instability and design of underground excavations[J]. International Journal of Rock Mechanics and Mining Sciences, 2003, 40(7): 1027-1047.

[2] Sch?pa A, Pantaleo M, Walter T R. Scale-dependent location of hydrothermal vents: Stress field models and infrared field observations on the Fossa Cone, Volcano Island, Italy[J]. Journal of Volcanology and Geothermal Research, 2011, 203(3): 133-145.

[3] 侯明勛, 葛修潤, 王水林. 水力壓裂法地應力測量中的幾個問題[J]. 巖土力學, 2003, 24(5): 840-844.

Hou Mingxun, Ge Xiurun, Wang Shuilin. Discussion on application of hydraulic fracturing method to geostressmeasurement[J].Rock and Soil Mechanics, 2003, 24(5): 840-844.

[4] 何江達, 謝紅強, 王啟智, 等. 官地水電站壩址區初始地應力場反演分析[J]. 巖土工程學報, 2009, 31(2): 166-171.

He Jiangda, Xe Hongqiang, Wang Qizhi, et al. Inversion analysis of initial geostress in dam site of Guandi Hydropower Project[J]. Chinese Journal of Geotechnical Engineering , 2009，31(2)：166-171.

[5] 于崇, 李海波, 李國文, 等. 大連地下石油儲備庫地應力場反演分析[J]. 巖土力學, 2010, 31(12): 3984-3990.

Yu Chong, Li Haibo, Li Guowen, et al. Inversion analysis of initial stress field of Dalian underground oil storage cavern[J]. Rock and Soil Mechanics, 2010, 31(12): 3984-3990.

[6] 金偉良, 喻軍華, 鄒道勤. 模擬高邊坡初始地應力場的邊界位移法[J]. 土木工程學報, 2003, 36(10): 72-75.

Jin Weiliang, Yu Junhua, Zou Daoqin. A new method of simulating initial ground stress field in high slope[J]. China Civil Engineering Journal, 2003, 36(10): 72-75.

[7] Buck W R, Lavier L L, Poliakov A N B. Modes of faulting at mid-ocean ridges[J]. Nature, 2005, 434(7034): 719-723.

[8] 吳仲瑋, 孫曉明, 王琰. 中印度洋脊 Edmond 熱液區多金屬硫化物中超顯微金銀礦物的發現及其成礦意義[J]. 礦物學報, 2013(2): 670-671.

Wu Zhongwei, Song Xiaoming, Wang Yan, et al. The discovery of native gold in massive sulfides from the Edmond hydrothermal field, Central Indian Ridge and its significance[J]. Acta Petrologica Sinica, 2013(2): 670-671.

[9] Petukhov S I, Anokhin V M, Mel'nikov M E, et al. Geodynamic features of the northwestern part of the Magellan Seamounts, Pacific Ocean[J]. Journal of Geography and Geology, 2015, 7(1): 35.

[10] Hyndman R D, Drury M J. The physical properties of oceanic basement rocks from deep drilling on the Mid-Atlantic Ridge[J]. Journal of Geophysical Research, 1976, 81(23): 4042-4052.

[11] 白世偉, 李光煜. 二灘水電站壩區巖體應力場研究[J]. 巖石力學與工程學報, 1982, 1(1): 45-56.

Bai Shiwei, Li Guangyu. Research on stress field around dam area of Ertan Hydropower Station[J]. Chinese Journal of Rock Mechanics and Engineering, 1982, 1(1): 45-56.

[12] 柴賀軍, 劉浩吾, 王明華. 大型電站壩區應力場三維彈塑性有限元模擬與擬合[J]. 巖石力學與工程學報, 2002, 21(9): 1314-1318.

Cai Hejun, Liu Haowu, Wang Minghua. Stress field simulation and fitting by 3D elastoplasticity FEM for large hydropower project[J]. Chinese Journal of Rock Mechanics and Engineering, 2002, 21(9): 1314-1318.

[13] Jacobs J A, Russell R D, Wilson J T. Physics and Geology[M]. New York: McGraw-Hill, 1959.

[14] Evans R L. A seafloor gravity profile across the TAG hydrothermal mound[J]. Geophysical Research Letters, 1996, 23(23): 3447-3450.

[15] Humphris S E, Cann J R. Constraints on the energy and chemical balances of the modern TAG and ancient Cyprus seafloor sulfide deposits[J]. Journal of Geophysical Research: Solid Earth (1978—2012), 2000, 105(B12): 28477-28488.

[16] Rona P A, Hannington M D, Raman C V, et al. Active and relict sea-floor hydrothermal mineralization at the TAG hydrothermal field, Mid-Atlantic Ridge[J]. Economic Geology (plus the Bulletin of the Society of Economic Geologists), 1993, 88(8):1989-1989.

[17] Tivey M K, Humphris S E, Thompson G, et al. Deducing patterns of fluid flow and mixing within the TAG active hydrothermal mound using mineralogical and geochemical data[J]. Journal of Geophysical Research: Solid Earth, 1995, 100(B7): 12527-12555.

[18] Martel S J, Muller J R. A two-dimensional boundary element method for calculating elastic gravitational stresses in slopes[J]. Pure and Applied Geophysics, 2000, 157(6/8): 989-1007.

[19] Anderson E M. The dynamics of faulting and dyke formation with applications to Britain[M]. Edinburgh:Oliver & Boyd, 1951.

[20] Humphris S E. Rare earth element composition of anhydrite: implications for deposition and mobility within the active TAG hydrothermal mound[C]//Proceedings-Ocean Drilling Program Scientific Results. National Science Foundation, 1998: 143-162.

Analyzing the gravitational stress field to forecast hydrothermal field- A case study of TAG hydrothemal field

Chen Qinzhu1, Tao Chunhui2, Liao Shili2，Li Huaiming2, Deng Xianming2

(1.CollegeofLandscapeArchitectureandWaterConservancyEngineering,WenzhouVocationalCollegeofScienceandTechnology,Wenzhou325000,China；2.KeyLabofSubmarineGeosciences,SecondInstituteofOceanography,StateOceanicAdministration，Hangzhou310012,China)

According to drilling data, two layer’s model is established by corresponding physical properties. Based on topographic data, here we use ANSYS software to simulate the gravitational stress field of TAG hydrothermal field. The result shows that, the hydrothermal vent locations are corresponding to the low stress area of the maximum horizontal stress. The stress of the active TAG doom area are characterized by local stress low surrounding by relative stress high, and the inactive MIR doom area without typical hydrothermal vent topography features located on unclosed area of the maximum horizontal stress low. Based on these conclusions, we predicted five places where the hydrothermal liquid probably flew out. According to drilling result, significant anhydrite was recovered in these areas, indicating that gravitational stress could be used for hydrothermal activity prediction.

topography; TAG hydrothermal area; stress field; hydrothermal activity prediction

10.3969/j.issn.0253-4193.2017.01.005

2016-04-24；

2016-06-28。

國家重點基礎研究發展計劃(2012CB417305)；國際海域資源調查與開發“十二五”重大項目(DY125-11-R-01, DY125-11-R-05)；國際海底管理局(ISA)捐贈基金。

陳欽柱(1990—),男，溫州市蒼南縣人，主要從事海底硫化區域預測研究。E-mail:chenqinzhu@aliyun.com

*通信作者：陶春輝(1968—)，男，研究員，主要從事海底熱液硫化物研究。E-mail:taochunhuimail@163.com

P738.6

A

0253-4193(2017)01-0046-06

陳欽柱, 陶春輝, 廖時理, 等. 利用應力場預測熱液區域——以TAG區為例[J]. 海洋學報, 2017, 39(1): 46-51,

Chen Qinzhu, Tao Chunhui, Liao Shili, et al. Analyzing the gravitational stress field to forecast hydrothermal field- A case study of TAG hydrothermal field[J]. Haiyang Xuebao, 2017, 39(1): 46-51, doi: 10.3969/j.issn.0253-4193.2017.01.005