鐘振楠,康鳳新,宋明忠,郎旭娟,柳祿湧,傅朋遠,李志杰
1) 山東省地質礦產勘查開發局第六地質大隊,山東威海,264209;2) 山東省地質礦產勘查開發局,濟南,250013;3) 河北地質大學水資源與環境學院,石家莊,050031;4) 山東省深部金礦探測大數據應用開發工程實驗室,山東威海,264209;5) 河北省水資源可持續利用與產業結構優化協同創新中心,石家莊,050031;6) 河北省水資源可持續利用與開發重點實驗室,石家莊,050031;7) 河北省水文學及水資源重點學科,石家莊,050031
內容提要: 魯東地熱區在地質構造單元上位于沂沭斷裂帶昌邑—大店斷裂以東,地熱資源豐富。本文采集了魯東地熱區招遠地熱田內一眼2000 m深地熱井(DRZK01)中的40塊巖芯樣品進行巖石熱導率、巖石生熱率測試及分析,結合測溫資料及收集資料對該區地熱通量構成及分層熱流進行了分析研究;采集區內典型地熱流體樣品進行水化學分析并采用合適的地熱溫標估算了該區熱儲溫度;綜合研究成果建立了該區地熱成因概念模型。研究結果顯示,該區巖石熱導率數值較高,測量值在2.8~5.7 W/(m·K)之間,普遍高于上地殼平均熱導率,地溫梯度較高,為31.8℃/km;利用熱導率和地溫梯度計算的地熱通量102 mW/m2中熱傳導分量為(73.2±6.18) mW/m2,對流分量為(28.76±6.18) mW/m2,其中熱傳導分量中地殼熱流為22.5 mW/m2,地幔熱流為(50.74±6.18) mW/m2,二者比值為1∶1.98~1∶2.52,為“熱幔冷殼”型熱結構。石英溫標計算熱儲平均溫度為128.6 ℃,熱循環深度約3634 m。研究結果豐富了該區地熱系統理論并為該區地熱資源開發利用提供一定的理論支撐。
在全球氣候變化的大背景下,地熱資源由于具備清潔、方便、蘊藏量大且可再生等優點而日益受到世界各國重視(王貴玲等,2000;龐忠和等,2012,2014;袁利娟等,2020; 李泓泉等,2020; 原若溪等,2021)。招遠地熱田位于魯東地熱區,魯東地熱區是山地陸塊地殼演化最為復雜的地區,尤其是中生代以來,地殼活動尤為頻繁劇烈,地熱資源豐富(徐希強等,2015)。為了科學高效利用該區地熱資源,查明該區“源、通、蓋、儲”因素,開展地熱資源賦存成因機理研究是不可缺少的關鍵環節(金秉福等,2000;邱楠生,2001;張濤,2011;史猛等,2019;趙輝等,2019)。大地熱流、地溫梯度、巖石熱導率、巖石生熱率作為基本的地熱學熱參數,為地球動力學演化以及油氣生成研究提供地熱背景和基礎參數(邱楠生,1998;葉正仁等,2001;藺文靜等,2016;孫旭東,2020);李學倫等(1997)利用地下水中SiO2含量計算山東半島硅熱流值,認為膠北隆起在大地構造上屬于較高熱異常區。Jiang Guangzheng等(2016)依托一口位于膠東半島萊州灣附近的金礦勘探井的相關資料,計算得到一高質量大地熱流數據。趙輝等(2019)對膠東地熱田地熱流體的補徑排特征及地熱田的地熱資源進行了評價研究。前人研究為本文研究提供了豐富的資料支撐但是同時仍然存在很多問題,如該區巖石圈熱結構方面的研究基礎比較薄弱,地熱熱量構成及地熱成因等方面尚存在很多需要深入研究的問題。因此本文在前人研究的基礎上,以招遠地熱田內一眼2000 m深地熱井為研究對象,采集巖芯樣品、典型地熱流體樣品、收集相關數據,分析巖層熱導率、生熱率及地溫分布特征,計算該區地熱通量;估算熱儲溫度及熱循環深度,并建立了招遠地熱田地熱成因概念模型,豐富了該區地熱系統的相關理論,為該區地熱資源可持續、高效開發利用提供理論支撐。
魯東地熱區以牟平—即墨斷裂為界分為東西兩大構造單元,招遠地熱田位于魯東地熱區西部華北板塊膠遼隆起區的膠北隆起北部。魯東地熱區主要出露太古代、元古代、中生代、新生代的地層(圖1)(崔煜烽等,2018)。該區中生代印支期由于地殼發生強烈的脆性張裂變形作用,形成了一系列近SN向張性斷裂以及NE向、NW向和近EW向扭性斷裂,并導致下部地殼部分巖石熔融,從而產生大規模的花崗質巖漿,沿斷裂入侵形成花崗巖體,有的沿NE向斷裂入侵形成一些淺成脈巖,由于上述地質演變,區域上具有良好的地熱形成條件(田禹,2015)。
本次研究中,巖芯樣品采集工作在招遠地熱田一眼深2000 m,地理坐標為120°24′41″,37°21′38″的地熱鉆孔(DRZK01)中進行(地熱井具體位置見圖1),50 m取芯一次,共采集了40塊巖芯樣品,所取巖芯樣品巖性為中生代玲瓏序列大莊子單元的二長花崗巖(巖性柱狀圖見圖2)。巖芯熱導率數據是由中國科學院地質與地球物理研究所巖石熱物性與地熱測量實驗室檢測,使用德國生產的TCS(Thermal Conductivity Scanning)熱導率自動掃描儀,其測量范圍為0.2~25 W/(m·K),測量精度為±3%。樣品規格依照實驗室要求進行樣品預處理,樣品長度大于3 cm,至少保留了一個平整面,平整面上沒有裂隙、空洞,并且未遭受蝕變。巖芯放射性元素(鈾、釷、鉀含量)數據是由澳實分析檢測(廣州)有限公司檢測,其中鈾、釷元素含量采用ICP-MS(電感耦合等離子體質譜儀)分析,分析儀器是由美國生產的Agilent(型號7700X),測試范圍0.05×10-6~1000×10-6,精度控制相對偏差<10%,準確度控制相對誤差<10%。鉀含量采用X射線熒光光譜儀分析,分析儀器是荷蘭生產的PANalytical(型號PW2424),測試范圍0.01~20%,精度控制相對偏差<5%,準確度控制相對誤差<2%。放射性含量樣品制樣方式相同,樣品烘干后,破碎過篩,保留直徑1~2 mm樣品300 g做正樣,其余樣品為副樣,副樣留存備用。將正樣用無污染缽振動碾磨至200目;取30 g正樣加入偏硼酸鋰—四硼酸鋰熔劑,混合均勻,于熔爐中在1025 ℃熔化。熔液冷卻后,用包含了硝酸、鹽酸和氫氟酸的混酸消解定容,再用等離子體質譜儀定量分析。另取30 g正樣,均分為兩份:一份試樣中加入含硝酸鋰的硼酸鋰—硝酸鋰熔融助熔劑,充分混合后,高溫熔融。熔融物倒入鉑金模子形成扁平玻璃片后,再用X射線熒光光譜儀分析。同時另一份試樣放入馬弗爐中,于1000 ℃灼燒。冷卻后稱重。樣品加熱前后的重量差即是燒失量。本方法所測的鉀元素含量為燒失量的結果和XRF測得的元素氧化物結果相加。檢測數據真實可靠。
水質全分析樣品用500 mL塑料瓶采集并由中國地質科學院水文地質環境地質研究所國土資源部地下水礦泉水及環境監測中心進行檢測,水質全分析是按照中華人民共和國國家標準引用天然礦泉水檢驗方法(GB/T8538-2008)檢測。
本次研究中井溫測量工作使用Server6000便攜式數控測井系統完成,測量范圍0~230℃,終孔后靜井18 d 測溫,可保證井溫恢復至穩態或似穩態。并且所測溫度及對應深度都經過Matrix Logging System處理后獲取,溫度測量精確度為±0.1 ℃。野外測溫時,探頭下降速率約為6 m/min,探頭每隔0.1 m記錄一個溫度數據,抽析間隔2 m。
鉆孔測溫曲線及地溫梯度曲線見圖2,由圖2可以看出,測溫曲線整體為“上凸”型,說明該鉆孔溫度受到區域上升流的影響,地下水作為載體將熱流攜帶至地表淺部(王一波等,2019),因此直接利用鉆孔中巖石熱導率數據和地溫梯度數據通過公式

圖2 招遠地熱田鉆孔柱狀圖及測溫曲線圖Fig. 2 Drilling column diagram and temperature measurement curve of geothermal well in the Zhaoyuan geothermal field
q=-K(dT/dH)

的計算值不能稱為大地熱流,可稱為地熱通量。由圖2看出,該區地溫增長分為5個層段:① 0~500 m層段地溫增長迅速,平均地溫梯度為159.3 ℃/km,地溫梯度過大,這可能是由地表溫度影響強烈原因造成的,因此進行溫度場分析時,為了減小誤差地層表層的溫度數據一般會被舍棄(胡圣標等,1994;Huang Shaopeng et al.,2000;Pollack et al.,2003);② 500~946 m層段溫度隨深度平穩增加,表明熱量傳遞主要以熱傳導方式進行,地溫梯度在該層段較為穩定,平均值為54.7 ℃/km;③ 1300~1800 m層段由于受到上升流影響,地溫梯度基本不變;④ 946~1171 m層段和1300~1809 m層段地溫梯度出現負增長,可能是由于局部破碎導致冷水混入造成的;⑤ 1171~1300 m層段溫度比較穩定,當溫度曲線出現分段現象時,地溫梯度可采用上下兩段的溫度曲線來進行校正(徐明等,2011)。綜合上文對于地溫梯度的分段分析,該區地溫梯度由500~1300 m層段的測溫數據計算獲取為31.8 ℃/km。
巖石熱導率是表征巖石傳熱特性的物理量,它是研究巖石圈熱結構和地球深部熱狀態的重要參數之一(倪守斌等,1999;欒錫武等,2003;章邦桐等,2010)。本次研究所取巖芯熱導率測量值在2.8~5.7 W/(m·K)之間,平均值為4.0 W/(m·K)(圖3)。

圖3 招遠地熱田鉆孔熱導率值統計直方圖Fig. 3 Statistical histogram of thermal conductivity in the Zhaoyuan geothermal field
巖石本身特性(巖石礦物成分、顆粒粒度、孔隙度、滲透率等)、溫度、壓力等因素會影響巖石熱導率數值,對于花崗巖等致密巖石,影響其熱導率的因素主要是溫度(Brigaud et al.,1989)。由于本次研究中取到的巖石樣品均為致密的二長花崗巖,巖石本身特性較一致,因此僅對巖石熱導率測試數據進行了溫度校正,校正過程是根據Sass等通過實驗得出的經驗公式進行的(Sass et al.,1992):
(1)
(2)
式中,λ(0),λ(25)分別為0℃和25℃時巖石的熱導率,T為巖石原位形成溫度,λ(T)為其原位溫度對應下的巖石熱導率,即為校正后巖石熱導率。按照上述方法對研究區測量的40個巖石熱導率數據進行了校正,校正后的巖石熱導率數據見圖4。由圖4b可以看出,本次研究獲得的巖芯熱導率值可以分為三部分(圖4b中分別由三個橢圓圈出),每一部分熱導率值基本都呈現出隨著深度增加而增大趨勢;而且相同巖性的巖芯熱導率數值存在一定變化范圍,這可能是由相同巖性之間物質成分存在差異造成的,也可能是受物質成分分布樣式等因素影響(崔景偉等,2019)。校正后熱導率數值分布在2.5~5.3 W/(m·K)之間,均值為3.5 W/(m·K),所取巖石樣品熱導率測量值都高于上地殼平均熱導率2.5 W/(m·K),反映出該區具有良好的地熱地質背景。

圖4 招遠地熱田校正后鉆孔巖芯熱導率統計直方圖(a)及熱導率與深度關系圖(b)Fig. 4 Statistical histogram of thermal conductivity after correctionand relationship between thermal conductivity and depth of the Zhaoyuan geothermal field
巖石的放射性生熱率(A)是單位體積巖石中所含放射性元素在單位時間由衰變所釋放的能量,單位為μW/m3,它表征著巖石自身生熱能力的高低。巖石放射性生熱率是描述地球內部熱狀態的一個非常重要的熱物性參數,可以為研究地球深部熱分布狀態、熱流的構成、探討地球動力學過程提供重要的數據支撐(胡圣標等,1994,2001;何麗娟等,2006;Tang Boning et al.,2019)。
巖石中所含的放射性元素很多,但是對生熱產生貢獻的元素主要是U、Th、K 三種元素,巖石放射性生熱率(A)就是通過測量計算巖石中U、Th、K 這3種元素的含量獲取(李文慶,2015;郎旭娟,2016)。本文巖石放射性生熱率數據是按照Rybach于1976年提出的計算公式獲取(Rybach,1976;邱楠生,2002):
A=0.01ρ(9.52CU+2.56CTh+3.48CK)
(3)
其中,A為巖石放射性生熱率(μW/m3),ρ為巖石密度(g/cm3),CU、CTh、CK分別為巖石中鈾(μg/g)、釷(μg/g)、鉀的含量(%)。
本次研究中巖石放射性生熱率計算結果見圖5。巖石生熱率數值統計直方圖(圖5a)表明,該區巖石生熱率值分布在0.4~2.2μW/m3之間,計算數據多數集中在1.2~1.6μW/m3之間,平均值為1.3μW/m3。圖5b顯示,生熱率數值和深度之間存在弱線性關系,生熱率計算結果總體上隨著深度增加而減小。

圖5 招遠地熱田鉆孔生熱率統計直方圖(a)及鉆孔生熱率隨深度變化圖(b)Fig. 5 Statistical histogram of heat production and relationship between heat production and depth of the Zhaoyuan geothermal field
地層深部的放射性生熱率的分布是解決該區地殼上地幔熱流配分及殼內熱結構、熱狀態等問題的關鍵(趙平,1995)。由于實際情況限制,本次研究中只取到了該區0~2000m層段的巖石樣品,對于2000m以下地層深部巖石放射性生熱率數據利用Rybach通過實驗得出的VP—A關系公式獲取(Rybach et al.,1984;趙平,1995;劉峰等,2020):
50 MPa壓力條件下:
(4)
100 MPa壓力條件下:
(5)
200 MPa壓力條件下:
(6)
其中,A為巖石放射性生熱率,VP為地震波速。
利用上述公式,代入收集的研究區地震波速數據(圖6),得出區內5000~36000 m(36000 m深度約為該區莫霍面深度(方寶明,2006))層段的巖石生熱率數據,詳細數據見表1。

圖6 華北克拉通超長觀測距探測剖面不同塊體二維速度結構圖(據王帥軍等修改,2014) Fig. 6 Two-dimensional velocity structure diagram of different blocks in the North China Craton ultra-long observation range detection profile(after Wang Shuaijun et al, 2014&)
由表1可以看出,利用VP—A關系公式計算的巖石生熱率值隨著深度的增加而遞減,對比0~2 km層段實測生熱率曲線(表1),可以看出公式5和公式6計算的結果偏大,利用公式4計算的生熱率結果最接近實際測量值。因此將公式4(50 MPa壓力條件下)的計算結果作為該區深部巖層生熱率結果。

表1 招遠地熱田 5.0~36.0 km(莫霍面)生熱率計算結果及0~2 km實測值
地熱流體中,某些決定熱水中溶解組分比例的水巖反應依賴于溫度,因此可用這些組分的比值來估算熱儲溫度,這種熱儲溫度估算方法稱為地熱溫標法,常用的地熱溫標有二氧化硅地熱溫標計、Na—K地熱溫標計、Na—K—Ca地熱溫標計等(Fournier,1977;Arnorsson,1985;劉昭等,2014)。根據研究區地熱流體的水化學特征及地熱地質情況,結合各地熱溫標適用條件,選取了合適的地熱溫標對該區6口地熱井進行熱儲溫度估算,估算結果見表2。

表2 招遠地熱田深部熱儲溫度估算結果Table 2 Estimated results of reservoir temperature in the Zhaoyuan thermal field
由計算結果可以看出,在多種地熱溫標中,玉髓地熱溫標計算結果偏低,因為鉆孔DRZK01井深1000 m以上時,實測溫度達120 ℃左右,故玉髓溫標計算結果不是熱儲溫度的真實反映;Na—K及Na—K—Ca地熱溫標計算出的結果偏高,石英溫標和Na—K—Ca溫標計算結果相差較大,且Na—K—Ca溫標計算結果偏高,遠高于石英溫標計算結果,而石英溫標計算結果高于K—Mg溫標計算結果,故最終結果選用石英溫標計算的熱儲溫度(102.88℃~148.84℃)為該區熱儲溫度的真實反映,均值為128.63℃。
地熱田內地熱流體補給來源主要為大氣降水(Sanliyuksel et al.,2011;Dotsika,2012),大氣降水入滲后經過深循環上升至地表形成地熱水,循環深度越大,地熱水溫度越高。地下熱水循環深度可由下列公式推算出來:
(7)
其中,H為地下熱水循環深度,m;t1為地下熱水深部最高溫度,℃,取上文計算的熱儲溫度均值128.63 ℃,t2為該區恒溫帶溫度,取當地常年平均氣溫14 ℃;h為恒溫度厚度,m,取值為30 m;I為地溫梯度,取上文分析的0.0318 ℃/m。將以上數值代入公式(7)中得出該區地下熱水循環深度約為3634.7m。
7.1.1地熱田地熱通量中對流分量
該區地熱通量利用公式
來計算,q為地表熱流(mW/m2),k為巖石熱導率[W/(m·K)],dT/dH為地溫梯度(℃/km),負號表示熱流傳導方向。式中地溫梯度取值為上文利用鉆井(DRZK01)500~1300 m 層段測溫曲線分析的31.8 ℃/km,巖石熱導率數據采用鉆井(DRZK01)500~1300 m層段巖芯熱導率檢測平均值3.2 W/(m·K),最終利用公式q=-k(dT/dH)計算該區地熱通量約為 102 mW/m2,具有較高的地熱背景。由于該區內未開展過系統的大地熱流測量工作,該區大地熱流值可取利用距離該地熱田約40 km的萊州三山島黃金科鉆的數據計算的熱流值 73.24±6.18 mW/m2(Jiang Guangzheng et al.,2016),由此可推算出本次研究中計算的地熱通量102 mW/m2中對流分量為(28.76±8.76) mW/m2。
7.1.2地熱田地熱通量中傳導分量
7.1.1 節中采用的大地熱流值73.24±6.18 mW/m2即為該區地熱通量中的傳導分量,即在地表觀測到的大地熱流,主要由地殼中放射性元素衰變產生熱能(地殼熱流qc)和地幔深部熱能(地幔熱流qm)兩部分組成(Jaupart et al.,2007;Czechowski et al.,2012;姜光政等,2016;劉紹文等,2017)。地殼中放射性元素衰變產生的熱能
qc=Σqa
式中:qa=AD, 為地殼中各結構層由放射性元素衰變產生的熱流,A為生熱率,D為相應層段厚度。地幔熱流由公式
qm=q-qc
獲取(Birch et al.,1968;Roy et al.,1968;徐青等,1992;Seipold et al.,1998 ;Ranalli,2005;Duchkov et al.,2009)。進而可獲得一個地區地殼熱流與地幔熱流配分比,分析其巖石圈熱結構類型(Huang Shaopeng et al.,1992)。
本次研究中,0~2 km層段巖石生熱率數值采取本文實際測算的巖石生熱率均值1.3 μW/m3;2 km以下每個層段的生熱率數值用表1中利用公式4(50 MPa壓力條件下)計算的結果,取值為各層段上下層面生熱率的均值。利用公式分別計算每一層段的地殼熱流值、地幔熱流值和熱流值,具體計算結果見表3。

表3 招遠地殼熱結構計算結果表Table 3 The results of thermal structure of the lithosphere in Zhaoyuan
由計算結果可以看出,招遠地區放射性元素衰變產生的熱能和為22.5 mW/m2,即地殼熱流為22.5 mW/m2,區內大地熱流采用73.24±6.18 mW/m2,就可估算出區內地幔熱流為(50.74±6.18) mW/m2,進一步估算出區內殼、幔熱流比為1∶1.98~1∶2.52,巖石圈熱結構屬于“熱幔冷殼”型,與汪集旸等分析的中國東部地區熱結構都為為“熱幔冷殼”型研究結果一致(汪集旸等,1986;唐顯春等,2020)。
由上述分析可知,該區計算的地熱通量102 mW/m2由三部分熱量構成,對流分量為28.76±8.76 mW/m2; 傳導分量中地殼放射性元素衰變產生熱量為 22.5 mW/m2; 深部地幔熱量為 50.74±6.18 mW/m2,三者比例約為1∶0.8∶1.76,地熱田內熱量約有二分之一是來自深部地幔。
綜合上文研究,研究區地熱成因概念模型如圖7所示:在較高的區域熱背景下(73.24 mW/m2),地熱田周圍山區的降水入滲后通過斷層或構造有利地段進入熱田地下水流動系統,然后沿著熱田東西兩側2 km處出露的深大斷裂—招平斷裂和玲瓏斷裂向深部運移,在運移過程中,不斷接受來自地幔深部的熱量(50.74±6.18 mW/m2)、巖石放射性元素衰變產生熱量(22.5 mW/m2)和對流熱量(28.76±8.76 mW/m2),溫度不斷升高(被加熱至約120~148 ℃),在約3.6~3.8 km深的部位兩條斷裂相交,從而使得地熱水沿著兩組深大斷裂交匯處形成的良好通道向上運移近地表處形成地熱異常區(地熱田)。兩條深大斷裂——招平斷裂和玲瓏斷裂對該區地下熱水的形成起到了重要作用:① 兩條斷裂是大氣降水(冷水)徑流到深部從而進行加溫的主要通道;② 兩條斷裂在約3.8km深處交匯,從而形成深部熱流上涌的有利通道。兩處斷裂多期活動,尤其是玲瓏斷裂在第四系以來仍有活動的痕跡,致使處于兩條斷裂中間位置的地熱田內構造發育,形成以NNE向為主的構造群,對其深入調查發現,NNE向F1斷裂(圖7)為主要構造,破碎帶寬15~20 m,發育深度大于2 km,斷裂北部有礫石充填,南段張性較好。從構造行跡及附近節理互切關系分析,F1斷裂為玲瓏斷裂的伴生斷裂,至少經歷了先壓后張兩次地殼運動,是一條長期活動至第四系且仍有活動的斷裂構造。依據水文地質資料顯示,地熱田內F1斷裂處富水性最強,大于10 m3/h,遠離F1斷裂,富水性逐漸減弱,進一步推斷F1斷裂為地熱流體上涌的主要通道。

圖7 招遠地熱田地熱成因概念模型Fig. 7 The conceptual model of geothermal genesis in Zhaoyuan thermal field
(1)通過測試分析顯示,該區巖石熱導率數值呈現出隨著深度加深而增大的規律,數值分布在 2.5~5.3 W/(m·K)之間,平均值為4.0 W/(m·K),高于上地殼平均熱導率2.5 W/(m·K);依據鉆孔測溫資料計算,該區地溫梯度值為31.8 ℃/km,結合熱導率測試數據,估算該鉆孔地熱通量為102 mW/m2,反映出該區良好的地熱地質背景。
(2)在多種地溫溫標中,石英溫標最適合用來計算該區熱儲溫度,估算熱儲溫度約為128.6 ℃;地下熱水循環深度約3634 m。
(3)利用地震波推算出的該區深部巖石生熱率分布規律與實測情況基本一致,并以此為依據推算出該區各層段的熱流分布規律。通過分析,該區地熱通量中既有來自傳導熱量(73.24 mW/m2)的貢獻也有對流熱量(28.76±8.76 mW/m2)的貢獻,其中傳導熱量中地殼熱流與地幔熱流配分比約為1∶1.98~1∶2.52,巖石圈熱結構屬于“熱幔冷殼”型。
(4)該區地熱系統屬“斷控型”熱儲,地下熱水補給通道及上涌通道主要是招平斷裂和玲瓏斷裂兩條深大斷裂,并且玲瓏斷裂至今仍在活動;物探資料顯示兩條深大斷裂在地下約3.8 km處交匯,這也正是本文計算的該區地下熱水的熱循環深度范圍,驗證了本文研究的正確度及可信度。