基于GMS的某典型石化園區地下水污染運移研究

劉景蘭，李立偉，朱明奕，葛菲媛，秦 磊

(1.天津市地質研究和海洋地質中心，天津 300170； 2.天津市生態環境科學研究院，天津 300191)

自1984年以來，我國石化園區經歷了從有到無，再發展壯大的過程。1990年全國僅有省級以上園區13個，2014年則發展至近700家，且其中半數以上分布于我國沿海地區。隨著我國石油化工行業的迅速發展，石化園區對周邊水土環境構成的污染風險也成為了管理部門和學界關注的焦點。石化園區中的生產活動，往往會造成較為嚴重的水土污染問題。近些年，隨著我國石化園區環境管理制度的不斷完善，園區生產所產生的廢水、廢氣和廢渣逐步實現了統一管理，導致污水處理廠和固廢暫存場地等區域成為了石化園區的潛在污染源。園區污水處理設施在非正常狀況下可能發生滲漏事故，在設施檢修的間隔期，污染物滲入地下水對周邊環境所構成的風險尤其值得關注，針對非正常狀況下污染物在地下水中遷移過程的預測對園區環境管理和園區周邊敏感目標的保護具有重要的現實意義。

隨著計算機模擬技術的發展，數值模擬方法被廣泛用于地下水污染物遷移的模擬預測。目前國內外最為常用的數值模擬軟件包括GMS(Groundwater Modeling System)、Visual MODFLOW和FEFLOW等，其中GMS和Visual MODFLOW是集成MODFLOW、MODPATH和MT3DMS等模擬程序開發的商業軟件，這些模擬程序都是基于有限差分方法的，而FEFLOW軟件則是基于有限元方法開發的，這些模擬軟件為地下水污染物遷移轉化的預測研究提供了有效工具。翟勇等[1]利用GMS軟件對我國東北某化工場地地下水中揮發酚類物質在未來10年中的遷移和污染物濃度衰減過程進行了模擬預測。何華燕等[2]利用FEFLOW軟件對氯代烴類污染物在地下水中的稀釋擴散過程進行了模擬，發現鹵代烴類污染物遷移能力較弱，擴散范圍較小，對周邊環境影響有限。魏傳云等[3]利用Visual MODFLOW軟件的MT3DMS模塊，對某化工園區污染物在不同情景下的運移規律進行了模擬預測，并重點針對連續源強和瞬態源強條件下污染物遷移的差異性進行了對比分析。姜鳳成等[4]采用數值模擬方法對某化工場地石油類污染物和重金屬Zn在地下水中的時空演變規律進行了對比分析，發現石油類污染物泄漏后在含水層中可發生自然衰減，污染物濃度和范圍隨時間逐漸減小，而重金屬Zn會隨地下水不斷向下游遷移。王平等[5]針對某有機污染場地地下水中污染物對周邊擬建水庫的潛在影響進行了數值模擬分析，發現水庫建設會改變周邊地下水流場，從而使得地下水污染羽范圍有所擴大，但并不會對水庫水質造成影響。張鳳娟等[6]對淄博市張店化工園區建設對地下水環境的影響進行了模擬預測，發現在20年中石油類污染物在地下水中遷移了350 m的距離，影響面積達到3萬m2，對地下水環境影響嚴重，有必要采取有效的防滲和監測措施加以防范。鄧一榮等[7]利用GMS軟件的MT3D模塊對粵港澳大灣區某典型化工污染場地苯系污染物在地下水中的遷移過程進行了模擬分析，發現溶解態苯系污染物沿地下水流向往下游遷移，通過識別檢驗發現運移過程主要受對流彌散作用的控制。

本文以我國北方沿海某石化園區為研究區，以園區內污水處理廠作為潛在污染源，以石油烴為特征污染物，利用GMS軟件的MODFLOW和MT3DMS模塊建立石油烴三維遷移數值模型，對非正常狀況下污染物滲入地下水后在未來30年中的遷移過程進行模擬預測，研究可為我國沿海地區石化園區環境管理提供科學依據。

1 材料與方法

1.1 研究區概況

研究區范圍如圖1所示，該石化產業園區建立于2003年7月，產業結構以精細化工和醫藥制造為主，占地面積約746.84 hm2。研究區地貌類型屬海積沖積低平原，地勢平坦，地面坡度小于1/10 000。該石化園區被河流圍繞，北側為城排河，南側依次為荒地排河和獨流減河，西側城排明渠自北向南溝通城排河和獨流減河，東側是長青河。

圖1 研究區范圍和水系分布示意Fig.1 Schematic of range and water system distribution of the study area

1.2 研究區水文地質條件

研究區場地位于我國北方海積低平原的咸水及鹽鹵水分布區。根據區域水文地質調查成果[8-10]，研究區400 m以淺的平原松散地層孔隙水可劃分為4個含水巖組：第Ⅰ含水組屬于淺層地下水系統；第Ⅱ—Ⅳ含水組屬深層地下水系統。第Ⅰ含水組為咸水及鹽鹵水，底板埋深約100 m，其中頂部潛水含水層隔水底板埋深在18 m左右；第Ⅱ含水組在調查評價區內為咸水，咸水底界埋深在200 m以深；深部的Ⅲ、Ⅳ含水組均為承壓淡水。第Ⅰ含水組屬濱海平原沖海積層咸水及鹽鹵水區，該部分地下水無開發利用情況，不作為居民生活飲用水使用[11-12]。因此，研究將地下水環境保護目標設置為第Ⅰ含水巖組的頂部潛水含水層，以防止污染物進入該潛水含水層向下游遷移擴散對園區外部地下水和地表水體造成污染，并對下游敏感目標構成環境風險。

第Ⅰ含水巖組頂部含水層地下水水位埋深在0.24～1.29 m，平均水位埋深為0.70 m。包氣帶主要為人工填土層素填土、粉質黏土；含水層為第四系松散巖類，巖性以粉質黏土、淤泥質黏土為主。潛水含水層厚度為15.01～18.71 m，淤泥質黏土單層厚度2.90～9.60 m，粉質黏土單層厚度0.90～12.70 m。含水層滲透系數為0.05～0.20 m/d，單位涌水量為0.015～0.270 m3/(h·m)，富水性弱。地下水化學類型為Cl-Na型水，pH值為7.45～7.66，礦化度為15～54 g/L。

1.3 概念模型建立

1.3.1 含水層結構

模擬區內主要分布第四系松散巖類孔隙水，含水巖組由上至下分別為人工堆積層、全新統上組陸相沖積層和全新統中組淺海相沉積層，巖性由上至下依次為素填土、粉質黏土、粉土和淤泥質黏土。該潛水含水巖組富水性較差，其下伏粉質黏土厚度較大，分布連續，且水平滲透系數和垂直滲透系數在10-7～10-8cm/s，構成了地下水潛水良好的隔水底板。因此，本次將第四系松散巖類孔隙水概化為統一的潛水含水層，同時考慮垂向上的巖性變化，構成非均質各向異性的多層結構地下水潛水三維模型。

1.3.2 邊界條件

根據模擬區水文地質條件及地下水流場特征，模擬區內地下水由區域中心分別向北西、南西和南東方向的河流進行徑流排泄，幾乎無地下水的側向流入補給，因此將模型側向邊界均概化為定水頭邊界。模型的上邊界為潛水面，整個含水層系統通過該邊界接收大氣降水補給和蒸發排泄，從而與外界進行垂向上的水量交換。

1.3.3 水文地質參數

由于區內潛水富水性較差，且水質普遍咸化，園區內企業幾乎不對區內地下水潛水進行開發利用。因此，模擬過程中無需考慮人工開采地下水所導致的地下水潛水位的動態變化。基于上述原因，此次所構建地下水流數值模型為穩定流模型。

地下水穩定流模型所需率定的水文地質參數主要為滲透系數K，根據所收集的區內抽水試驗計算結果，結合場地含水層結構特征，確定潛水含水層水平滲透系數Kh為0.05～0.20 m/d，將垂直滲透系數Kv設置比水平滲透系數小一個數量級。根據經驗值，有效孔隙度設置為0.3。由于水動力彌散尺度效應的存在，難以通過野外或室內彌散試驗獲得真實的彌散度。因此，模擬參照經驗值，按照保守的評價原則，取縱向彌散度值為15 m，橫向彌散度值為1.5 m，垂向彌散度為0.15 m。

1.3.4 源匯項

模擬區內地下水主要接收降水補給，側向流入補給幾乎可以忽略。根據包氣帶巖性和地下水埋深特征，結合經驗值，確定大氣降水入滲系數為0.11。由于地下水潛水水位埋深較淺，一般為0.5～2.0 m，蒸騰作用較為強烈，區內地下水主要通過蒸發排泄，少部分地下水通過側向流出進行排泄。此外，根據地下水潛水監測流場形貌判斷，研究區內存在局部小規模地下水降落漏斗，推測該處可能存在少量地下水人工開采。

1.3.5 源強概化

根據該化工產業園區污水處理廠設計進水水質指標限值以及各污染因子進水濃度和標準限值所計算的標準指數，發現石油類標準指數最高，對地下水所構成的潛在威脅最大。因此，將研究中污水處理設施滲漏情境下對地下水可能造成影響的預測因子確定為石油類。

根據項目工程設計，污水處理各構筑物的池體均為鋼筋混凝土。池體的泄漏量參考《給水排水構筑物工程施工及驗收規范》(GB 50141—2008)中關于滿水試驗驗收的要求，鋼筋混凝土池體滿水試驗驗收標準為2.0 L/(m2·d)，假設項目在非正常狀況下池底由于地面不均勻沉降或地下水對池體的腐蝕等多種因素影響下，出現防滲層破裂情況，破裂程度引起的地下水滲漏量按照驗收標準的10倍計算[13-17]，即20 L/(m2·d)。

沉淀池池體長9.7 m、寬8.5 m、高5.0 m，本次采用池體的表面積作為泄漏面積，即264.45 m2，結合單位面積滲漏量20 L/(m2·d)，可計算得到非正常狀況下沉淀池池體滲漏量為5.289 m3/d。由工程分析可知石油類污染物的進水濃度為20 mg/L，假定污水處理設施的地下水防控或檢漏周期為120 d，即發生非正常狀況后120 d發現并進行修復切斷滲漏源，因此進入含水層中石油類污染物滲漏量為12 213.6 g。

1.4 數學模型和軟件選取

1.4.1 地下水流模型

針對本次所建立的非均質各向異性三維穩定流地下水流模型，選取下列數學模型和定解條件進行求解：

(1)

式中，h為地下水潛水位；Ω為地下水滲流區域；Γ1和Γ2分別為模型第一類邊界和第二類邊界；Kxx、Kyy和Kzz分別為x、y和z主方向的滲透系數；λ為通過單位水平面積潛水面注入含水層的水量；ω為源匯項，包括降水入滲補給和井抽水量等；h0(x，y，z)為第一類邊界的地下水水頭函數；q(x，y，z)為第二類邊界單位面積流量函數。

1.4.2 溶質運移模型

地下水水質模型是定量評價和研究地下水水質問題的數學模型。建立地下水水質模型時，將溶質運移過程中的對流、彌散和吸附等反應項考慮在內，用于描述地下水中溶質的運移狀態。模型的偏微分方程式及定解條件如下：

(2)

式中，方程左側前2項為彌散項，之后2項為對流項，最后一項為由化學反應或吸附解吸所產生的溶質增量；Dxx和Dyy分別為x和y兩個主方向的彌散系數；μx和μy為x和y方向上的實際水流速度；c和c0分別為溶質濃度和初始濃度；Ω為溶質運移的區域；Γ為三類邊界；φ為邊界溶質通量；v為滲流速度；gradc為濃度梯度。

1.4.3 軟件選取

研究所選用的GMS軟件綜合了MODFLOW、MODPATH和MT3DMS等已有的地下水模型基礎上研發的。GMS軟件相較于Visual MODFLOW等地下水模擬軟件最大的特點是提出了概念模型建模的方法，這為GMS與GIS軟件間提供了良好的數據接口，使得建模過程更為靈活簡便。除此之外，GMS因其強大的前處理和后處理功能以及優良的三維可視化功能，目前已經成為了國際上最受歡迎的地下水模擬軟件之一。

1.5 網格剖分和應力期劃分

1.5.1 網格剖分

空間網格剖分是建立地下水水流和溶質運移數值模型的基礎。考慮到需要利用所建立的模型對重點區域污染物遷移過程進行模擬預測，因而在相應區域進行了網格的加密剖分。東西向剖分網格64個，南北向剖分網格74個，垂向上剖分網格5層，共剖分網格23 680個，其中有效計算單元15 960個。具體空間網格剖分如圖2所示。

圖2 數值模型空間網格剖分示意Fig.2 Schematic of the spatial grid profiles of the numerical model

1.5.2 應力期劃分

此次模擬將地下水非穩定流模型和溶質運移模型模擬期劃分為2個應力期：①第1個應力期設置為120 d，對應污染物入滲階段；②第2個應力期設置為11 000 d(30年)，對應污染物擴散階段，最長時間為污水處理廠最長服務年限(30年)。

2 結果與討論

2.1 地下水流模型識別檢驗

2.1.1 地下水位擬合

通過模型參數識別，得出模型計算地下水流場和實際監測地下水流場擬合情況如圖3所示。

圖3 計算水位與實際監測水位擬合Fig.3 Fitting of calculated water level and actual monitoring water level

從圖3中可以看到，模型計算地下水流場和實際監測地下水流場擬合情況較好，二者之間的水位誤差整體上小于±0.1 m，說明所建立的地下水流模型能夠較為準確地反映評價區地下水潛水的真實流動狀態，所建立的模型可用于評價區地下水中污染物遷移過程的模擬預測。

2.1.2 地下水均衡分析

通過模型參數識別，得出模型驗證期地下水水量均衡結果見表1。從表1中可以看出，模擬區地下水潛水主要接受大氣降水入滲補給，并主要通過潛水蒸騰方式進行排泄，并有少量地下水通過側向流出進行排泄。

表1 模擬區地下水均衡分析Tab.1 Groundwater balance analysis in simulated area

2.2 污水處理設施滲漏模擬預測

2.2.1 地下水流場動態模擬

滲漏事故情景下，污水滲入地下會構成新的地下水潛水補給來源，除此之外，模擬區地下水系統的其他源匯項基本保持不變。因此，非正常狀況下污染物在地下水中的遷移預測，有必要在已有地下水流模型的源匯項條件和含水層特征的基礎上，增加入滲補給地下水潛水的污水水量，重新進行地下水潛水非穩定流模型計算。污水入滲時間設置為120 d，計算得到的污水處理廠周邊地下水流場變化情況如圖4所示。

圖4 非正常狀況下污水處理廠周邊地下水流場變化情況Fig.4 Changes of groundwater flow field around sewage treatment plant under abnormal conditions

從圖4中可以看到，隨著滲漏時間延長，在污水處理設施周邊開始形成地下水水丘，且規模不斷增大，在120 d時達到最大，且在模擬1 000 d時也未能完全恢復。地下水水丘的形成使得泄漏點周邊地下水水力梯度增大，這會加快污染物向周邊的擴散速率。

2.2.2 地下水溶質運移模擬

污水處理設施泄漏對地下水污染預測結果如圖5所示。

圖5 污水處理設施泄漏對地下水污染預測結果Fig.5 Prediction results of groundwater pollution caused by leakage of sewage treatment facilities

從圖5中可以看出，在模擬期內石油類對潛水含水層將會造成污染，并存在超標現象。石油類在潛水含水層中不同時間的超標范圍及最大運移距離見表2。

表2 非正常狀況下污水處理設施泄漏對含水層影響范圍Tab.2 List of influence range of aquifer caused by leakage of sewage treatment facilities under abnormal conditions

從表2中可以看到，1 000 d污染物最遠遷移41 m左右，超標范圍為4 879 m2左右；30年(11 000 d)污染物最遠遷移56 m左右，超標范圍為6 482 m2左右。根據模擬結果，在污染物發生泄漏30年后，地下水中污染物運移距離有限，污染羽前緣距離下游地表水體還有120 m左右的距離，這主要與研究區地下水水力梯度較小、同時含水層滲透系數較低有關。模擬結果表明，污水處理設施泄漏對周邊地表水體影響較小，但由于周邊土壤和地下水遭受污染，對泄漏點周邊敏感目標可能構成環境風險，有必要在污水處理設施下游進行地下水水質定期監測。

3 結論

(1)利用GMS構建了研究區三維非穩定地下水流和地下水溶質運移數值模型，通過對污水處理設施非正常狀況下廢水的泄漏進行模擬預測可知：池體泄漏條件下，泄漏點周邊會形成地下水水丘，在120 d時規模達到最大，加速了污染物向周邊的遷移速率；研究區屬低滲透性地層，預測期內污染物運移距離較小，尚未對下游地表水體產生影響。

(2)為防控地下水污染，避免對泄漏點周邊地下水環境造成影響，結合預測結果，建議在污水處理設施下游30～50 m處布設跟蹤監測井，定期開展地下水監測工作，以便及時發現污染并截斷污染源。

(3)此次模擬遵從保守性原則，未考慮污染物的吸附和降解等生物地球化學過程，后期可根據管理需要，通過補充相關反應動力學和熱力學參數，進一步提高模型模擬預測精度。