庫區洲灘潛流帶溫度示蹤流速計算方法

姬雨雨， 施文卿， 陳求穩， 安 磊

(1.河海大學 水利水電學院， 江蘇 南京 210098； 2.南京水利科學研究院， 江蘇 南京 210024)

1 研究背景

潛流帶是河流地表水和地下水動態交換的飽和過渡區域，是河流中溶質交換、生物新陳代謝與繁殖的重要場所，在維持河流生態系統健康和穩定過程中發揮著重要作用[1-3]。為應對全球能源需求和氣候變暖，世界各國推進水電開發[4]，大壩運行導致水位異常波動和泥沙沖淤平衡的改變，在庫區和壩下形成眾多河灣洲灘潛流帶。水庫晝夜發電引起水位波動，促使洲灘內部潛流交換和物質交換頻繁發生，已逐漸成為國內外生態水文研究的熱點區域[5-8]。潛流交換對潛流帶內物化性質(溶解性氣體、溫度、pH、鹽度、營養鹽等)的空間分布和河流水質具有重要影響[6,9-11]。因此，量化潛流帶的潛流交換流速和交換量成為了解其生物地球化學過程的重要基礎和關鍵。

在潛流交換和太陽輻射的雙重作用下，潛流帶內的溫度分布具有時空差異性，該溫度差異可作為天然的示蹤劑，用于表征地下水的流動。溫度示蹤法具有無污染、易于觀測和準確度較高的優點，在計算潛流交換流速和交換量方面的應用廣泛。Alexander等[12]利用連續監測的動態溫度數據，估算水文地質參數并量化潛流交換量，Constant等[13]利用溫度示蹤確定了潛流交換的深度范圍，Halloran等[14]利用溫度示蹤法計算了土體含水飽和度。然而，國內在潛流帶研究中對于溫度示蹤法的應用仍處于起步階段。

因此，本文以瀾滄江漫灣庫區洲灘潛流帶為研究對象，水庫運行期間實時監測洲灘內部水位和水溫的變化過程，利用溫度示蹤法計算潛流交換流速和交換量，豐富溫度示蹤法在庫區河流潛流帶中的應用。

2 研究區域概況

選取云南省瀾滄江漫灣庫區的一處河灣洲灘為研究對象(圖1)，漫灣水庫位于云南省臨滄市境內，建成于1993年，為瀾滄江上第一座以發電為單一目標的百萬千瓦級電站。河灣目標洲灘(24°43′44″N, 100°23′5″E)位于壩前30 km處，為多年泥沙淤積形成(圖2)，洲灘長、寬和淤積深度分別為182、90和60 m。試驗區具有兩個特點：(1)因水庫晝夜發電，庫區水位波動頻繁，洲灘內部潛流交換活躍；(2)試驗區地處云貴高原，太陽輻射強烈且四季氣溫差異較小。

圖1洲灘位置圖2試驗區沉積物層

3 研究方法

3.1 水位和水溫監測

沿庫區至洲灘方向布置W1 ～ W4四口水位和水溫監測井，距離洲灘邊緣分別為0.5、10.5、20.5、35.5 m(圖3)。監測井中安放長度均為4 m的盲管，側面每隔15 cm鉆有直徑為2 cm的小孔，并外覆以尼龍濾網，以實現井管內腔與外界水流通暢交換。監測井底部安裝有水位自動記錄儀(U2000101, OneSetHoBo, USA)，用于監測庫水位晝夜波動，其精度和適用范圍分別為0.5 cm和0 ～ 9 m，沿監測井底部向上每隔60 cm布置一個水溫自動記錄儀(U22001, OneSetHoBo, USA)，用于監測庫水位波動期間地下水溫隨時間發生的波動，利用水準儀確定監測井的安裝高程以保證具有共同的基準面。監測自2016年12月11日開始至2016年12月16日終止，時間間隔為10 min。

3.2 溫度示蹤流速計算方法

地下淺層溫度場受氣象和水力條件影響表現出周期性波動，熱量作為一種天然示蹤劑，可用于計算地下水流速的大小和方向，熱量在潛流帶中的運移數學模型：

(1)

式中:T為溫度，℃；t為時間，s；z為沉積物深度，m；q為水流沿深度方向在孔隙介質中的流速(向下為正)，m/s；ρw為水的密度，kg/m3；ρ為含水介質的密度，kg/m3；cw為水的比熱容，J/(kg·℃)；c為含水介質的比熱容，J/(kg·℃)；λe為含水介質的熱傳導系數，J/(s·cm· ℃)。

Hatch和Keery分別得到了熱傳導方程的解析解，但Keery法只考慮了相位滯后，其解析解僅是流速大小，無法判定方向。Hatch法則綜合考慮了波動振幅的衰減、相位滯后和熱彌散度，其解析解可判定流速的大小和方向，且具有更高的精確度[15-16]，因此本文采用Hatch法計算流速。Hatch法一維熱傳導方程的求解公式為：

(2)

圖3 研究區域監測井布置圖

本滲流模型采用VFLUX計算包計算，該解析解適用于均質多孔介質。通過參考相關資料和室內實驗，確定計算過程中所需的有效熱傳導系數、熱彌散度、水比熱容、土比熱容及洲灘沉積物飽和層孔隙率參數(表1)。

4 結果與討論

4.1 水位和水溫變化

洲灘監測井的水溫波動數據規律表現一致，僅在大小上存在差異。W1監測井靠近庫區，且W1靠近岸邊沒有草甸覆蓋，受水庫水位波動和太陽輻射的雙重作用更加明顯，因此選擇W1監測井水位和水溫數據進行分析。監測周期內研究區域的氣溫和河流水溫分別在12～30℃ 和15～16℃之間，洲灘內部水溫在垂向上表現出“上涼下暖”。在漲水的0～33 h內，受水庫持續蓄水和晝夜發電用水影響，監測井水位呈現持續波動升高趨勢(圖4(a))，因低溫地表水滲入洲灘內部，表層、中層、底層水溫均隨水位上升而下降，分別從開始時刻的17.9、18.2、18.3℃下降至33 h的17.2、17.8、18.1℃(圖4(b))，該時段W1水位高于W2水位，表現為庫區持續補給洲灘，在33 h后，由于水庫持續泄水，W1和W2監測井水位隨之下降(圖4(a))，低溫地表水滲入量逐漸減小，底層水溫在第40 h后開始逐漸上升(圖4(b))。受太陽輻射和底層水溫傳導擴散的雙重影響，表層和中層水溫分別在第40 h和50 h開始上升(圖4(b))。在63 h后，W2水位高于W1水位(圖4(a))，表現為洲灘補給庫區，低溫地表水逐漸潛出，使得表層和中層水溫在63 h后均呈較快上升趨勢，最終分別穩定在18.4℃和18.5℃(圖4(b))。

4.2 水動力學方法

多點原位深層采樣發現洲灘為河流細顆粒泥沙淤積形成，質地比較均一，因此通過變水頭入滲實驗測得平均滲透系數為8.64×102m/d，并按照達西定律計算W1監測井附近的潛流交換流速和交換量(圖5)。在0～63 h，庫區向洲灘補給，流速呈現先增加后減小趨勢，從開始時刻的7.21×10-5m/s 上升為8 h的13.49×10-5m/s，之后呈現遞減趨勢，在63 h處降至0 m/s，監測井間單位寬度上的補給總量為1.85 m3。在63～120 h，洲灘向庫區排泄，流速呈先增加后穩定趨勢，從63 h的0增加到115 h的5.15 m/s，監測井間單位寬度上的排泄總量為0.76 m3。

4.3 溫度示蹤法流速計算

采用VFLUX計算包中的Hatch法計算潛流交換流速和交換量。振幅法計算結果表明，0～66 h為庫區補給洲灘階段，流速呈先增大后減小趨勢，從開始時刻的3.32×10-6m/s逐漸增大至30 h的 12.81×10-6m/s，之后減小至81 h的0，單位寬度上的補給總量為1.92×10-1m3。81～115 h為洲灘補給庫區階段，流速呈逐漸增大的趨勢，從81 h的0逐漸增加至115 h的2.89 m/s，單位寬度上的排泄總量為0.21×10-1m3(圖6(a))。相位法計算結果表明，0 ～ 81 h為庫區補給洲灘階段，流速呈先增大后減小趨勢，從開始時刻的3.32×10-5m/s逐漸增大至30 h的 12.92×10-5m/s，之后減小至66 h的0 ，單位寬度上的補給總量為1.87 m3。66～115 h為洲灘補給庫區階段，流速呈先增大后緩慢減小至穩定的趨勢，從66 h的0逐漸增加至80 h的5.87×10-5m/s，之后緩慢減小并穩定在115 h的3.92×10-5m/s，單位寬度上的排泄總量為0.79 m3(圖6(b))。在庫區水位漲落周期內，振幅法和相位法均呈現出庫區與洲灘間先補給后排泄的過程，但二者潛流交換流速和交換量的數量級均相差10倍，對比水動力學法計算結果，Hatch相位法的結果準確度更高，所以之后的對比均采用Hatch相位法。

4.4 水動力學方法和溫度示蹤法計算結果對比

將水動力學法和Hatch相位法計算結果進行對比(圖7)。庫區補給洲灘階段，在0～8 h，水動力學法計算的流速逐漸增加，在8 h時達到峰值流速13.49×10-5m/s，Hatch相位法計算的流速在0～30 h呈現逐漸增加的趨勢，在30 h時達到峰值流速12.92×10-5m/s；洲灘補給庫區階段，水動力學法計算的流速在80 h時達到峰值流速5.52×10-5m/s，Hatch相位法計算的流速在83 h達到峰值流速5.87×10-5m/s。

可見，利用溫度示蹤法計算因庫區水位波動引起的潛流交換時，流速的峰值存在一定的延遲，這是因為溫度在多孔介質中的擴散和傳導速度滯后于地下水流速。

對比流速曲線，水動力學計算的流速曲線波動明顯，這主要是因為水庫發電引起頻繁性的水位波動，且具有較為密集的數據采集頻率0.1 次/min，而溫度示蹤Hatch相位法計算的潛流交換流速波動性較小，主要是由于溫度示蹤法有數據重取樣的環節，其數據采集頻率為0.5 次/h，數據量相對較少，因此溫度示蹤計算的流速曲線波動性較小。總體來看，Hatch相位法得到的流速范圍為-5.87×10-5～12.92×10-5m/s，單位寬度上潛流交換總量為2.66 m3，水動力學法計算的流速范圍為-5.52×10-5～13.49×10-5m/s，單位寬度上潛流交換總量為2.61 m3，二者較為接近。因此，溫度示蹤法作為一種新興的計算潛流交換流速和交換量的方法，能夠較好地反映水位波動引起潛流交換的流速趨勢，并具有較高的準確性。

表1 相關計算參數

圖4庫區洲灘水位和水溫變化特征圖5水動力學法流速計算結果

圖6 Hatch振幅法和相位法計算得到的流速曲線圖7水動力學方法和溫度示蹤法計算結果對比

4.5 不同深度流速對比

為比較洲灘不同深度處的流速大小，利用Hatch相位法計算了地面下1.1、1.7 m深度處的流速，1.7 m深度處的流速值q2比1.1 m處的流速值q1小(圖8)，主要是因為水位波動引起的潛流交換在潛入和潛出洲灘的過程中存在水壓消散和水頭損失，即深度越深，地下水流速越小。另外，利用溫度示蹤法計算洲灘不同深度處的流速時，對淺層流速變化過程的刻畫要比深層好，這是因為溫度波動信號在沿地層深度方向的傳導和擴散過程中存在信號衰減。

圖8 不同深度地下水流速對比

5 結 論

(1)溫度示蹤法計算結果中，Hatch相位法和振幅法得到的流速曲線趨勢基本一致，僅在水庫落水的第80 h相差較大。對比計算結果，水動力學法計算的潛流交換流速和交換總量分別為-5.52×10-5～13.49×10-5m/s和2.61 m3，Hatch相位法計算的潛流交換流速和交換總量分別為-5.87×10-5～12.92×10-5m/s和2.66 m3，Hatch相位法具有更高的準確性。

(2)在庫區補給洲灘和洲灘補給庫區兩個階段，利用水動力學法計算得到的流速分別在8和79 h達到正負峰值，溫度示蹤法計算的流速分別在30和83 h達到正負峰值。溫度示蹤法計算流速達到峰值的時間存在一定滯后。

(3)溫度示蹤法計算結果表明洲灘上層的流速比下層大，這說明離地越深，傳導的溫度波動信號越衰減，流速越小。