，，Jun Shimada，

(1.長江科學院 水土保持研究所，武漢 430010； 2.中山大學 地理科學與規(guī)劃學院，廣州 510000； 3.熊本大學 自然科學研究院，日本 熊本 8608555； 4.長江勘測規(guī)劃設計研究院 市政與交通院，武漢 430010)

1 研究背景

地下水資源作為重要的戰(zhàn)略儲備，在供水安全、健康及生態(tài)保護方面意義重大；地下水資源相關研究成果對科學制定水資源管理措施，嚴格實行水資源管理制度具有重要支撐[1]。近年來，隨著全球氣候變化、國家城鎮(zhèn)化進程加快和人口的急劇增長，地下水資源危機日益嚴重[2-3]，嚴重威脅國民用水安全、破壞生態(tài)平衡、引發(fā)一系列鏈生災害，制約我國社會經濟的可持續(xù)發(fā)展。地下水超采和海水入侵造成地下水采補失衡，可供水資源量減少[4]；城市、工農業(yè)污水和廢水任意排放，引起地下水硝酸鹽和重金屬離子超標，威脅城市用水安全[5-6]；被污染的地下水物質通量排泄入海，引發(fā)近海岸的生態(tài)、環(huán)境和漁業(yè)問題[7]。上述由地下水資源危機引發(fā)的一系列問題對經濟社會發(fā)展和生態(tài)環(huán)境保護帶來極不利的影響。

對地下水資源進行準確評價和科學管理是解決我國社會經濟發(fā)展過程中面臨的地下水資源和環(huán)境問題的有效途徑。準確評價地下水資源，是對地下水資源的數量、質量、時空分布特征和可開采量進行科學、全面的分析和估計[8]。地下水資源評價的常規(guī)手段(原位觀測、模型模擬和同位素示蹤等方法)具有一定局限性，迫切需要一種成本低、易操作、無污染、能連續(xù)操作的天然示蹤方法。受地下水、地質條件和外部因素影響，地質體溫度具有一定時空差異，這些差異能夠表征地下水活動[9]。地下溫度場不僅與滲流場產生耦合作用，還受溫度場邊界條件影響，因此，地下水溫度能夠有效記錄地表溫度變化，反演地表溫度時間序列[10]。

本文綜述了地溫示蹤技術在地下水科學中的應用，從地質體能量變化特征及其影響因素出發(fā)，分析“地表溫度-地下溫度場-滲流場”系統(tǒng)中的科學問題。在已有綜述研究[9,11]的基礎上，將溫度示蹤從“工程地下水滲漏探測” 和“地表水-地下水交換研究”領域中的應用擴展到“地表暖化-地下水流動對地溫影響剖析”、“含水層參數識別”和“濱海含水層水循環(huán)過程解析”，闡述地溫示蹤方法在水文學研究領域的應用研究。

2 地溫示蹤原理與方法

2.1 地溫分布的影響因素

地溫隨深度增加而增加是一個普遍規(guī)律，但其幅度有很大差異，這與地質構造和區(qū)域的自然條件、水文地質條件及各種要素有關[12]。區(qū)域性地溫分布差異一般由地質構造決定，一定深度的地溫垂向分布特征主要受地表溫度變化和地下水活動的影響。地表溫度變化和地下水流動過程對天然地熱梯度產生干擾，這種干擾通常強烈而迅速，并顯示為清晰的溫度變化信號，使溫度隨深度的變化曲線發(fā)生偏離[13]，溫度信號的偏離程度可用于研究地表溫度變化和地下水活動。

地表溫度變化通過地殼巖石的熱傳遞向地底深處傳播，溫度變化幅度隨深度的增加呈指數遞減[14]。地殼巖石作為一個低頻濾波器，使長時地表溫度變化信號比短時信號傳遞得更深[14]。地表短時日變化和季節(jié)變化信號只能傳遞到地下十幾米深處，地表溫度的長期變化趨勢會影響上百米深處的地溫分布，地表溫度變化對地溫的擾動量可以用來精確反演地表溫度變化。Pollack等[15]利用全球358個鉆孔地溫數據重構地表氣溫變化，結果表明過去5個世紀全球氣溫增加了約1 ℃，僅20世紀全球氣溫就增加了約0.5 ℃。Huang等[16]的研究基于全球區(qū)域616個鉆孔溫深曲線，構筑了地表過去5個世紀氣溫變化，結果表明北半球氣溫增加幅度比南半球高。利用地溫擾動重構地表溫度變化方法簡單易行[17-18]，但其計算精確性受地溫偏移量提取和地殼巖層熱擴散系數率定方法的影響[19-20]，相關問題的解決是當前研究的重點。

圖1 地下水流動和地表暖化影響下溫度-深度曲線(據Anderson[13]修改)Fig.1 Temperature-depth profiles influenced by groundwater flowing and surface warming (revised from Anderson[13])

溫度在地殼巖層中的熱傳遞作用包括熱傳導和熱對流，在地下水活躍的含水層，由水體流動引起的熱對流作用占主導地位。如圖1所示，地表溫度恒定、無地下水流動條件下，均質地層地溫隨深度呈線性遞增[17]；地表溫度恒定、地下水補給和排泄條件下，均質地層的地溫-深度曲線向上分別呈現“凹”和“凸”不同形態(tài)，入滲和排泄速度值越大，曲線的彎曲程度越大[13]；受地表溫度增加和地下水流動同時影響時，溫度-深度曲線在淺層部分向右偏移[10,21-22]。受地表溫度變化和地下水流動的影響，地層溫度受到擾動，溫度-深度曲線的形態(tài)是表征這種擾動的信號，能夠有效反演地下水流動和地表溫度變化。建立在地下滲流場和溫度場的耦合作用上，不同邊界條件下的水-熱流運移控制方程描述了地層溫度的分布特征和演化過程，其解析解或數值解能揭示地表溫度變化和地下水流動對地層能量特征的影響。

2.2 地溫示蹤方法

利用溫度示蹤法探明地表溫度變化和地下水流動的研究始于20世紀50年代左右，其理論基礎可以歸結為Stallman[23]提出的各向同性、均質、飽和、多孔介質三維非穩(wěn)定水-熱流運移的控制方程?？刂品匠叹_描述了地下水流場和溫度場的耦合過程，是利用地溫示蹤地表溫度變化和地下水流動的理論基礎。

水-熱流運移控制方程描述了含水層中水體運移過程的能量特征，當地下水流速=0時，反映地表溫度變化在地殼巖層中的熱傳導過程；地下水流速≠0時，則描述地殼巖層中的熱傳導和熱對流傳熱過程。因此，該方程在不同邊界條件下的解析解或數值解能有效揭示地表溫度變化和地下水流動對地溫分布情況的影響。

解析解使用具體函數的形式，描述地下溫度與時間、地下水流速、介質物理特性等自變量之間的關系。解析解對地溫分布和變化情況的理解有益，但其推求過程受限于當前數學物理發(fā)展水平，僅適合簡單形式問題。利用解析解法對水-熱運移控制方程進行求解能夠有效探究和剖析地表氣溫變化-地下水流場-地下熱場之間的交互作用和機理。

國外眾多學者研究了不同邊界條件下一維和二維水-熱穩(wěn)定和非穩(wěn)定流問題(表1)，采用數學物理方法導出其解析解，并用于重構地表溫度演變過程、估算含水層及不同界面間水流滲流速度、反演率定地殼巖層熱力學參數等。數值解法通常將偏微分方程簡化為線性方程，并以數值方法為媒介，求出其近似解。數值解雖然是控制方程的近似結果，但對大多數情況而言，近似解足夠準確，能夠滿足應用要求。同時，數值解法適用范圍廣，不同邊界和初始條件的水-熱流運移控制方程均能通過編程求出數值解。

20世紀80年代以后，數值模擬被廣泛應用到水-熱流運移的研究中。對一維非穩(wěn)定流水-熱運移問題來說，數值解通常用于驗證推導的解析解的正確性[27,33]。數值模擬方法更為廣泛應用于二維和三維復雜邊界非穩(wěn)定流水-熱運移問題：Ferguson和Woodbury[34]對加拿大Winnipeg地區(qū)地下溫度進行監(jiān)測，并使用二維數值模型進行模擬，結果表明由于城市建筑釋放熱的原因，建筑區(qū)下的地溫相對其他區(qū)域高5 ℃；Ferguson等[35]從理論上探索了地下水流動對地溫的影響，指出在利用地溫重構地表氣溫變化的研究中，必須剝離區(qū)域地下水流動對地溫的影響；Bense和Beltrami[36]基于二維水-熱數值模型，對森林砍伐地區(qū)，地下水區(qū)域流動對地溫的影響進行模擬和分析，結果表明地下水側流對地溫擾動較大，在實際研究工作中應予以考慮；Schardtj等[37]使用三維數值模擬方法研究火山硫化物礦區(qū)的水-熱運移過程；Yang[38]使用三維數值方法模擬了斷陷沉積巖盆地的水熱運移過程；水-熱流運移控制方程的數值解法還用于估算復雜水文地質條件下含水層的地下水流速[39]和分析受抽水影響含水層的熱力學擾動[40]。

表1 水-熱運移控制方程解析解及其應用實例總結Table 1 Summary on analytical solutions of water-heat transport governing equations and their applications

隨著計算機技術的發(fā)展，出現了眾多針對水-熱流運移控制方程數值計算的軟件[13]，其中大部分軟件都擁有可視化界面，有利于軟件的推廣和應用。根據軟件側重點不同進行分類，TOUGH2主要應用于對地熱系統(tǒng)的模擬和解析，SHEMAT和BASIN2側重于對復雜地質過程涉及的水-熱運移過程的剖析，VS2DH，SUTRA和FEFLOW用于淺層含水層地下水動力學和熱力學過程模擬。

3 地下水科學中的應用

3.1 地面暖化-地下水流動對地溫影響的研究

地下水流動過程對天然地熱梯度產生干擾，這種干擾通常強烈而迅速，并顯示為清晰的溫度變化信號，使溫度隨深度的變化曲線發(fā)生異常[13]。地溫梯度不僅記載了區(qū)域地下水流動信息，還能有效揭示地表氣溫突變、氣候變暖等區(qū)域氣候變化特征[14,16]。根據Taniguchi等[21]和Anderson[13]的研究，地面暖化情形下，“凹”和“凸”不同形態(tài)的溫深曲線在淺層部分向右偏移，地表溫度增加的速率和幅度決定了偏移程度的大小(圖1)?；谏鲜隼碚摚琓aniguchi團隊在UNESCO GRAPHIC (Groundwater Resources Assessment under the Pressures of Humanity and Climate Change)項目的支持下，分別對亞洲大型城市東京[21]、大阪[41]、曼谷[42]城市化過程中，地下水流動對地溫的影響進行了評價。Taniguchi等[31]和Abe等[43]還分別對日本熊本和奈良平原地溫進行監(jiān)測和研究，利用地溫探明了城市化過程中區(qū)域地下水資源的時空分布特征。

地面暖化-地下水流動對地溫影響研究是在地下水流活躍的含水層，探究由地表溫度變化驅動的地殼巖體溫度場與地下水滲流場之間的耦合作用。在氣候變化和人類活動的雙重影響下，地表溫度持續(xù)增加，地溫擾動因素由傳統(tǒng)的“一元”轉變?yōu)椤岸?。相關研究由地溫示蹤地下水流速[24-25]或地表升溫速率[14,16]轉向剖析地表溫度變化-地下水流動對地溫梯度的影響[10]。區(qū)域地下水資源量通常由地下水的垂向“補給”和“排泄”量決定[8]，相關學者[10,32-33]側重使用解析方法剖析地表溫度波動在一維非穩(wěn)定流中的傳播過程。Bense和Beltrami[36]基于二維水-熱數值模型，對森林砍伐地區(qū)，地下水區(qū)域流動體系對地溫的影響進行模擬和分析，結果表明地下水測流對地溫擾動較大，在實際研究工作中應予以考慮；筆者在博士論文[44]中，使用二維水-熱數值模型對日本熊本地區(qū)城鎮(zhèn)化過程影響、地表溫度線性增加情景下，區(qū)域地下水流動對地溫的擾動情況進行模擬，并依據計算結果繪制了1930—2010年間日本熊本地區(qū)地溫分布情況。成果(圖2)表明：①地面暖化對地溫的影響多集中在近地面淺水層；②由于地下水入滲作用的影響，地溫的暖化效應在地下水入滲區(qū)傳導較深；③受地下水擾動影響，排泄區(qū)地溫梯度遠大于入滲區(qū)。

圖2 1930—2010年日本熊本地區(qū)地表暖化情形下二維地下水系統(tǒng)地溫分布(來自Dong[44])Fig.2 Subsurface temperature distribution in Kumamoto plain, Japan from 1930 to 2010 under the condition ofsurface warming (from Dong[44])

3.2 含水層水文地質參數計算研究

含水層水文地質參數是反映含水層或透水層水文地質特性的指標，是進行地下水資源評價的重要基礎資料，也是水文地質計算時不可或缺的數據。常規(guī)的含水層參數計算方法從地下水動力學角度出發(fā)，基于現場抽水試驗資料，采用傳統(tǒng)優(yōu)化算法、高斯-牛頓法和遺傳算法，對水文地質參數進行反演[45-46]。但是，單一的水頭分析方法在含水層參數率定過程中存在不確定性，不同降水入滲系數和含水介質導水系數組合可以模擬得到相同的水頭分布[47-48]，因此無法準確獲得含水層參數信息；同時，描述水頭分布的數學模型在某些參數估算過程中并不收斂，導致無法率定含水層參數[46]。鑒于此，相關研究[45-48]同時通過比較不同水文地質參數下溫度模型模擬值和實測值，實現對含水層水文地質參數的率定。

利用地溫示蹤技術計算含水層水文地質參數具有一定的局限性。水文地質參數主要包括含水介質的滲透系數和導水系數、承壓含水層的儲水系數、潛水含水層的重力給水度、降水入滲系數等參數，溫度示蹤法通常用于率定含水層的滲透系數[47-48]。含水層水力學和熱力學參數敏感性分析研究[49]表明：相較于熱力學參數而言，含水層滲透系數和孔隙率是含水層熱運移過程的主控因素。在理想化的條件下，含水層中溫度的時空分布特征主要受地下水流速的影響，溫度信號可以反演地下水流速，含水層的滲透系數則可以通過地下水流速和水力梯度推算得到。但是，在各向異性含水層中，熱傳遞過程較復雜，無法準確量化地下水流動產生的溫度變化[49-51]，地溫示蹤技術僅能用于率定各向同性含水層的水文地質參數[51-52]。

3.3 濱海含水層地下水循環(huán)過程解析

濱海地區(qū)是人類經濟和社會活動最為活躍的地區(qū)，人類活動極易破壞濱海地區(qū)含水層系統(tǒng)的平衡狀態(tài)，造成地面沉降、海水入侵以及海水污染等環(huán)境地質問題[53]。受當前監(jiān)測技術的限制，濱海含水層水文地質結構和參數獲取較為困難，對濱海地區(qū)地下水循環(huán)過程研究造成一定影響。地溫示蹤技術具有低成本、易實現、能連續(xù)監(jiān)測的優(yōu)點，在研究濱海地區(qū)地下水循環(huán)領域具有一定應用前景。當前利用溫度解析濱海含水層地下水循環(huán)過程主要集中在海水入侵和海底地下水入海通量計算。

海水入侵地下水是由于人為超采地下水造成水動力平衡的破壞，海水入侵會引起地下水水質變咸、土壤鹽漬化和水田面積減少等一系列環(huán)境問題。相關研究[54-55]通常通過對地下水水頭和鹽度的監(jiān)測和分析來研究海水-淡水界面特征，受海岸帶復雜水文地質條件限制，需要相對經濟和可靠的新方法，實現對海水入侵地下水的研究。Taniguchi[56]首次嘗試使用溫度監(jiān)測手段代替常規(guī)分析方法對海水入侵地下水狀況進行解析：對日本礪波平原觀測井的溫-深曲線進行測量，根據垂向一維水-熱運移原理，將地溫梯度產生變化的深度確定為咸淡水界面位置。上述方法與常規(guī)水力學方法分析的結果比較符合，具有較強實用性。

海底地下水入海通量(Submarine Groundwater Discharge, SGD)對濱海生態(tài)系統(tǒng)的營養(yǎng)平衡起著十分重要的作用，對其時空分布規(guī)律的研究是當前水文水資源領域的熱點[53]。目前國際上采用原位觀測、水文模型和同位素示蹤的方法對SGD進行計算[7]。溫度作為一種低成本、穩(wěn)定和無污染的示蹤劑，可應用于SGD相關研究中。Taniguchi等[26]使用地溫示蹤方法對西澳海岸帶的SGD進行估算，分析結果與滲流儀的觀測數據具有一定差異性，可能是溫度只能對來自內陸的地下水補給量進行示蹤，忽略了海水進入陸地后經過循環(huán)又排入海水的水量。

4 研究展望

(1)傳統(tǒng)的地溫示蹤技術主要用于探測工程地下水滲漏、估算地下水滲流速度和地表水-地下水交換研究中，隨著溫度監(jiān)測儀器的精細化和數學物理方法的發(fā)展，地溫示蹤技術被應用于多邊界、復雜水文地質環(huán)境問題的剖析中，如剖析地下水流動-地面暖化對地溫的影響，估算水文地質參數，解析濱海含水層地下水循環(huán)過程。

(2)當前，國際水文地質學的研究熱點已從單一的含水層水動力學問題向過渡帶含水層水流、氣流、鹽分和污染物動態(tài)及交互作用方面轉移[53]。地下水溫度作為能量的表征指標，與含水層水體的滲流過程聯系較為緊密，相關研究大多建立在各向同性、均質、飽和流的基礎上，利用溫度示蹤技術解析地下水科學問題目前僅適用于各向同性含水層，關于地下水氣相和固相的交互過程中溫度變化的原理有待進一步深化。

(3)當前，我國在地溫示蹤技術應用側重于“工程地下水滲漏探測” 和“地表水-地下水交換研究”領域，利用溫度示蹤復雜水文地質過程方面的研究遠落后于美國、日本等發(fā)達國家，相關理論及在國內研究區(qū)域的應用亟需加強。