陳小月，黃健民，盧 薇

(1.廣東省廣州市地質(zhì)調(diào)查院，廣東 廣州 510440，2.廣州地理研究所，廣東 廣州 510070)

自2007年以來(lái)，廣州金沙洲陸續(xù)發(fā)生巖溶地面塌陷、地面沉降地質(zhì)災(zāi)害[1-3]，據(jù)已有的相關(guān)研究，上述地質(zhì)災(zāi)害的發(fā)生與地下水的異常波動(dòng)密切相關(guān)，并主要發(fā)生在地下水位劇烈波動(dòng)的下降及恢復(fù)期間。為了更好、更深入分析金沙洲地下水變化情況與地面塌陷、地面沉降地質(zhì)災(zāi)害的關(guān)系，本文根據(jù)金沙洲地帶的水文地質(zhì)條件、鉆孔、地下水的水位動(dòng)態(tài)及抽水試驗(yàn)成果，應(yīng)用 FEFLOW軟件建立金沙洲地下水系統(tǒng)的數(shù)值計(jì)算模擬模型。模擬區(qū)范圍以金沙洲界線為界，將整個(gè)金沙洲作為地下水系統(tǒng)的模擬計(jì)算范圍。

FEFLOW是由德國(guó)柏林水資源規(guī)劃與系統(tǒng)研究所(WASY)開發(fā)的一種對(duì)地下水水量及水質(zhì)進(jìn)行模擬的軟件系統(tǒng)，利用它不僅可以計(jì)算出水位、溶液濃度和溫度等標(biāo)量數(shù)據(jù)，而且還可以模擬降水、地表水、地下水的流動(dòng)與轉(zhuǎn)換，計(jì)算出流速、流線及流徑線等向量數(shù)據(jù)，很適合分析模擬區(qū)地質(zhì)災(zāi)害與地下水流場(chǎng)的關(guān)系[4-6]。

1 模擬區(qū)水文地質(zhì)概況

金沙洲四面環(huán)水，東臨沙貝海—白沙河，西側(cè)約1 km為水口涌，沙貝?！咨澈蛹八谟烤鶠榱飨又Я鳎瑑芍Я髟诮鹕持弈喜繀R入珠江。區(qū)域水文主要受珠江干流控制，局部受流溪河影響。在金沙洲模擬區(qū)內(nèi)，其西側(cè)低丘陵的地層為石炭系下統(tǒng)大賽壩組的砂巖、頁(yè)巖、泥巖等，這些巖層隔水性好;東側(cè)石炭系壺天群、石磴子組灰?guī)r及白堊系大塱山組泥質(zhì)粉砂巖與珠江白沙河及沙貝海相連，地下水與珠江水存在水力聯(lián)系;南面石磴子組灰?guī)r及大塱山組泥質(zhì)粉砂巖與佛山市黃岐相接，與區(qū)外地下水存在水量交換。因此，金沙洲具有完整的補(bǔ)、徑、排條件，形成了一個(gè)相對(duì)獨(dú)立的一級(jí)水文地質(zhì)單元。在自然狀態(tài)下，金沙洲地下水順勢(shì)總體自北西向南東流動(dòng)，最終排泄于珠江水系。

2 地下水系統(tǒng)概念模型

為建立模擬區(qū)地下水系統(tǒng)的概念模型，首先要對(duì)實(shí)際的水文地質(zhì)條件加以概化，簡(jiǎn)化或忽略與系統(tǒng)無(wú)關(guān)的某些系統(tǒng)要素和狀態(tài)，以便于數(shù)學(xué)描述，并建立地下水系統(tǒng)模擬模型。

2.1 模擬區(qū)含水層結(jié)構(gòu)的概化

根據(jù)金沙洲區(qū)內(nèi)的地下水系統(tǒng)含水介質(zhì)的物質(zhì)組成及水文地質(zhì)特性，可將其內(nèi)部結(jié)構(gòu)概化為四層，由上至下分別為:

(1)弱含水層(弱透水層):由人工填土、淤泥質(zhì)粘土、淤泥組成。

(2)第一含水層(第四系含水層):主要由砂礫層夾砂，中粗砂、中砂、粗砂組成。

(3)隔(弱透)水層:由殘積礫質(zhì)粘性土組成。

(4)第二含水層(基巖巖溶裂隙含水層):主要由壺天群和石磴子組灰?guī)r組成。

2.2 模擬區(qū)邊界條件概化

金沙洲模擬區(qū)東面緊鄰珠江水系，并與模擬區(qū)地下水有水量較緩，可視為定水頭邊界;西面丘陵地區(qū)為金沙洲地下水含水系統(tǒng)的分水嶺，其山脊線可視為零流量邊界，即為隔水邊界;在模擬區(qū)南側(cè)，第一層含水層和第二含水層與外圍地下水有水量交換，主要為地下徑流側(cè)向補(bǔ)給，可作為徑流邊界;人類工程活動(dòng)抽排地下水導(dǎo)致的降水中心概化為模型的第四類邊界。潛水含水層自由水面為系統(tǒng)的上邊界，通過(guò)該邊界的潛水與系統(tǒng)外發(fā)生垂向水量交換;下部以基巖中風(fēng)化層底板為底部邊界，按零通量邊界處理。

3 地下水系統(tǒng)數(shù)值模擬模型

3.1 地下水系統(tǒng)數(shù)學(xué)模型

從金沙洲模擬區(qū)的空間上看，地下水流整體上以水平運(yùn)動(dòng)為主、垂向運(yùn)動(dòng)為輔，方向性明顯，且地下水系統(tǒng)的輸入、輸出隨時(shí)間和空間變化，為了準(zhǔn)確模擬各含水層的相互影響，將模擬區(qū)內(nèi)的地下水流作為非均質(zhì)、各向異性介質(zhì)中的三維非穩(wěn)定流處理，建立地下水系統(tǒng)數(shù)學(xué)模型如下:

式中:Ω為滲流區(qū)域;h為含水層的水位標(biāo)高(m);Kx、Ky、Kz為分別為 x、y、z方向的滲透系數(shù)(m/d);Kn為邊界面法向方向的滲透系數(shù)(m/d);S為自由面以下含水層儲(chǔ)水系數(shù)(1/m);μ為潛水含水層在潛水面上的重力給水度;ε為含水層的源匯項(xiàng)(1/d);p為人工開采和降水等(1/d);h0為含水層的初始水位分布(m);φ1為第一類邊界上的已知函數(shù);Γ1為滲流區(qū)域的上邊界，即地下水的自由表面;Γ2為滲流區(qū)域的側(cè)向邊界;Γ4為滲流區(qū)域的下邊界，即承壓含水層底部的隔水邊界;n→為邊界面的法線方向;q(x，y，z，t)為定義為二類邊界的單寬流量(m2/d·m)，流入為正，流出為負(fù)，隔水邊界為0。

在所建立地下水系統(tǒng)概念模型中，所有邊界均按照第二類邊界條件來(lái)處理。

3.2 模擬區(qū)網(wǎng)格劃分

將模擬區(qū)范圍內(nèi)平面(即每一個(gè)模擬層)剖分為9 350個(gè)結(jié)點(diǎn)，14 092個(gè)三角單元。模型垂向上共分4個(gè)模擬層，4層模擬層的總結(jié)點(diǎn)數(shù)為9 350×4=37 400個(gè)，總單元數(shù)為14 092×4=56 368個(gè)。金沙洲模擬區(qū)的三維網(wǎng)格剖分圖見圖1。

圖1 金沙洲地下水?dāng)?shù)值模擬區(qū)網(wǎng)格剖分圖

3.3 模擬區(qū)水文地質(zhì)參數(shù)劃分

金沙洲地下水?dāng)?shù)值模擬模型主要的水文地質(zhì)參數(shù)為滲透系數(shù)、給水度、儲(chǔ)水系數(shù)及彈性釋水系數(shù)等。

為保證對(duì)計(jì)算精度不產(chǎn)生影響，并使水文地質(zhì)參數(shù)分區(qū)簡(jiǎn)化，將每一個(gè)滲透系數(shù)分區(qū)同時(shí)視為給水度、儲(chǔ)水系數(shù)分區(qū)。滲透系數(shù)的初始分區(qū)按地形地貌、沉積類型和地層巖性等進(jìn)行劃分，如表1、表2及圖2所示。各區(qū)的參數(shù)初始值根據(jù)抽水試驗(yàn)、沉積類型和地層巖性、地層時(shí)代等特征進(jìn)行估值，垂向滲透系數(shù)根據(jù)類似地區(qū)的計(jì)算經(jīng)驗(yàn)估計(jì)給值。

表1 含水層給水度和釋水系數(shù)分區(qū)表

表2 第四系砂層滲透系數(shù)分區(qū)表 m/d

圖2 金沙洲地下水?dāng)?shù)值模擬計(jì)算參數(shù)分區(qū)圖

基巖裂隙含水層考慮到斷層對(duì)含水層導(dǎo)水性的影響，沿?cái)鄬拥姆较蛏删彌_區(qū)，加大緩沖區(qū)的滲透系數(shù)，見表3。

3.4 源匯項(xiàng)處理

源匯項(xiàng)主要考慮地下水開采量、垂向入滲補(bǔ)給量和蒸發(fā)量?，F(xiàn)狀條件下，模擬區(qū)地下水開采主要為工程活動(dòng)抽排地下水，概化為井單元處理。垂向入滲主要為大氣降水的入滲補(bǔ)給和地表水的入滲補(bǔ)給，分別取多年加權(quán)平均降雨量和河流補(bǔ)給系數(shù)法確定補(bǔ)給量，模擬過(guò)程中不考慮蒸發(fā)量。

表3 基巖巖溶裂隙含水層滲透系數(shù)分區(qū)表 m/d

3.5 模型檢驗(yàn)與校正

為了驗(yàn)證所建立的數(shù)學(xué)模型和模型識(shí)別后確定的水文地質(zhì)參數(shù)的可靠性，采用試估—校正法，取 K36、K37和 K723三個(gè)鉆孔觀測(cè)水位數(shù)據(jù)與模擬水位數(shù)據(jù)進(jìn)行對(duì)比，對(duì)數(shù)值模擬模型進(jìn)行檢驗(yàn)。根據(jù)收集到的地下水觀測(cè)資料，取模型識(shí)別的時(shí)段為2009年7月13日至2009年12月31日，取時(shí)間步長(zhǎng)為1天，共計(jì)1 171個(gè)時(shí)間步長(zhǎng)，以實(shí)際水位和計(jì)算水位的擬合情況作為調(diào)整參數(shù)的依據(jù)。從圖3中可以看出，經(jīng)過(guò)識(shí)別的參數(shù)計(jì)算的水位線與實(shí)測(cè)水位線趨勢(shì)相同，誤差都在容許范圍之內(nèi)。

圖3 實(shí)測(cè)水位與模擬水位檢驗(yàn)圖

4 地下水流場(chǎng)數(shù)值模擬計(jì)算

模擬區(qū)人類工程活動(dòng)抽排地下水造成地下水位波動(dòng)大致分為如下五個(gè)階段(見表 4)，可概化為 1 200 m3/d、2 300 m3/d和3 200 m3/d等三種典型情況，現(xiàn)以抽排水量3 200 m3/d(第四階段)為例進(jìn)行金沙洲地下水流場(chǎng)數(shù)值模擬計(jì)算，模擬結(jié)果見圖4和圖5。

表4 模擬區(qū)抽排地下水階段劃分及其抽排水量

圖4 抽排水量為3 200 m3/d時(shí)第四系砂層含水層模擬區(qū)流場(chǎng)分布圖

圖5 抽排水量為3 200 m3/d時(shí)基巖巖溶裂隙含水層模擬區(qū)流場(chǎng)分布圖

從地下水流場(chǎng)分布圖中可以看出，模擬區(qū)地下水向施工降水中心(豎井)流動(dòng)，形成以人類工程施工降水點(diǎn)為中心的降落漏斗地段。地面塌陷的分布范圍位于橫沙村—沙貝村—鳳岡隔水條帶以西，說(shuō)明隔水條帶阻擋了工程降水并接受來(lái)自于東側(cè)珠江河水的補(bǔ)給，使金沙中學(xué)附近的地下水沿?cái)鄬訌?qiáng)徑流帶由東北側(cè)向工程施工降水點(diǎn)方向流動(dòng)，導(dǎo)致這些地段的地下水流場(chǎng)動(dòng)態(tài)變化加劇。

野外調(diào)查發(fā)現(xiàn)的T9-T14六宗地面塌陷地質(zhì)災(zāi)害均分布于工程施工降水點(diǎn)附近及其所控制地面沉降范圍內(nèi)(BX1)(圖4、圖5);地面塌陷根據(jù)地下水監(jiān)測(cè)資料，在2008年8月15日至2009年5月8日時(shí)間段內(nèi)，金沙洲地下水出現(xiàn)異常波動(dòng)，總體處于快速下降過(guò)程，部分民井逐漸干涸，在此期間，金沙洲共發(fā)生地面塌陷地質(zhì)災(zāi)害14宗，占塌陷總數(shù)的58.3%(圖6)。因此，不論從空間還是時(shí)間角度來(lái)看，模擬區(qū)地面塌陷地質(zhì)災(zāi)害的發(fā)生都和模型分析結(jié)果吻合，說(shuō)明了該模型的建立是合理的。

圖6 巖溶地面塌陷與地下水位動(dòng)態(tài)變化關(guān)系圖

5 數(shù)值模擬的結(jié)論與認(rèn)識(shí)

金沙洲模擬區(qū)內(nèi)斷裂構(gòu)造發(fā)育，斷層及不同規(guī)模的節(jié)理裂隙密集，造成可溶巖發(fā)育多層裂隙和溶洞，洞隙之間連通性好，構(gòu)成區(qū)域地下水的主要徑流場(chǎng)。人類工程施工降水是引起金沙洲地下水流場(chǎng)變化的主要原因:抽排水期間，地下水的流動(dòng)方向由四周向降水中心流動(dòng)，地下水流場(chǎng)分布圖內(nèi)斷層發(fā)育地段的地下水流線加長(zhǎng)，并形成以施工抽水點(diǎn)為中心的降落漏斗，地下水位變化明顯，易引發(fā)地面塌陷和地面沉降地質(zhì)災(zāi)害。因此，基于上述模擬結(jié)果，應(yīng)嚴(yán)格控制金沙洲區(qū)域的工程活動(dòng)施工降水強(qiáng)度，以有效防止地面塌陷及地面沉降地質(zhì)災(zāi)害的發(fā)生。

[1]黃健民，鄭小戰(zhàn)，陳小月，等.廣州金沙洲廣州大學(xué)附屬實(shí)驗(yàn)學(xué)校地面沉降形成機(jī)理研究[J].熱帶地理.2012，32(04):338-343.

[2]郭宇，黃健民，陳建新，等.廣州市白云區(qū)金沙洲地區(qū)地質(zhì)災(zāi)害風(fēng)險(xiǎn)區(qū)劃[J].熱帶地理.2013，33(6):659-665.

[3]郭宇，黃健民，周志遠(yuǎn)，等.廣州市白云區(qū)金沙洲地區(qū)地質(zhì)災(zāi)害現(xiàn)狀及防治對(duì)策[J].中國(guó)地質(zhì)災(zāi)害防治學(xué)報(bào).2013，24(3):100-104.

[4]趙國(guó)紅，寧立波，王現(xiàn)國(guó).新鄭市淺層地下水流數(shù)值模擬及評(píng)價(jià)[J].地下水.2007，29(6):43-46.

[5]郭曉東，田輝，張梅桂，等.我國(guó)地下水?dāng)?shù)值模擬軟件應(yīng)用進(jìn)展[J].地下水.2010，32(4):5-7.

[6]張立志.基于 FEFLOW的大興區(qū)地下水動(dòng)態(tài)模擬研究[D].地質(zhì)大學(xué)(北京).1-66，2009.