某制藥廠廢水排放對(duì)地下水污染的數(shù)值模擬預(yù)測(cè)研究

楊 鑫,寧立波

(1.中國地質(zhì)大學(xué)(武漢)地球科學(xué)學(xué)院，湖北 武漢 430074;2.中國地質(zhì)大學(xué)(武漢)環(huán)境學(xué)院，湖北 武漢 430074)

制藥工業(yè)關(guān)系著國計(jì)民生，是經(jīng)濟(jì)發(fā)展的支柱行業(yè)之一。但是，制藥行業(yè)屬于精細(xì)化工，制藥廢水的組成十分復(fù)雜，具有有機(jī)污染物種類多、濃度高，COD和BOD5值高，NH3-N濃度和固體懸浮物SS濃度高、色度深、毒性大等特征[1-3]。據(jù)統(tǒng)計(jì)，我國制藥行業(yè)廢水的排放量較大，約占全國工業(yè)企業(yè)污水排放量的2%[2]。制藥廢水在排放過程中污水管道發(fā)生滲漏可能造成的地下水污染問題一直備受關(guān)注，也是目前地下水污染中較難解決的問題之一，其中的關(guān)鍵問題在于地下水中污染物的運(yùn)移特征及其范圍較難確定。數(shù)值模擬方法是解決該類問題的途徑之一，前人已做了大量的研究工作[4-6]，主要運(yùn)用的數(shù)值模擬軟件包括Modflow、GMS、Tough2等[6-9]，但針對(duì)制藥廢水泄漏而產(chǎn)生的地下水污染范圍進(jìn)行數(shù)值模擬的研究較少。新鄭市屬于鄭州市規(guī)劃的航空港區(qū)，其功能龐大，對(duì)環(huán)境質(zhì)量的要求較高，但尚未進(jìn)行地下水污染方面的環(huán)境影響評(píng)價(jià)工作?；诖?，本文以擬進(jìn)駐該區(qū)的某制藥廠為研究對(duì)象，運(yùn)用GMS軟件對(duì)該制藥廠廢水排放過程中污水管道發(fā)生滲漏可能造成的地下水污染范圍進(jìn)行了數(shù)值模擬研究，以為該地區(qū)地下水資源的合理開發(fā)與保護(hù)提供科學(xué)依據(jù)。

1 研究區(qū)概況

該制藥廠選址于河南省新鄭市規(guī)劃新區(qū)(見圖1)，該地區(qū)屬暖溫帶大陸性季風(fēng)氣候，多年平均蒸發(fā)量為1 859.7 mm，多年平均降雨量為683.64 mm。研究區(qū)地處黃河和雙洎河沖洪積平原區(qū)，地形平坦，地表徑流遲緩，地下水埋深較淺，且包氣帶巖性大部分為粉土，局部為粉砂，結(jié)構(gòu)松散，極有利于大氣降水的入滲補(bǔ)給。該地區(qū)地層以第四系為主，主要含水巖組為松散巖類孔隙含水巖組，其補(bǔ)給來源主要為大氣降水的入滲補(bǔ)給，地下水的水化學(xué)類型較簡單，以HCO3-Ca、HCO3-Ca·Mg型水為主。

圖1 研究區(qū)位置平面圖Fig.1 Location of the study area

2 研究區(qū)數(shù)值模型的建立

2.1 研究區(qū)水文地質(zhì)條件與模型概化

水文地質(zhì)概念模型是將含水層實(shí)際的邊界性質(zhì)、內(nèi)部結(jié)構(gòu)、滲透性質(zhì)、水力特征和補(bǔ)給排泄等條件進(jìn)行概化，便于進(jìn)行數(shù)學(xué)與物理模擬[4-5]。水文地質(zhì)概念模型的核心要素是邊界條件、內(nèi)部結(jié)構(gòu)和地下水流態(tài)，通過對(duì)調(diào)查區(qū)的巖性構(gòu)造、水動(dòng)力場(chǎng)、水化學(xué)場(chǎng)的分析，可以確定水文地質(zhì)概念模型的要素[6]。

2.1.1 研究區(qū)水文地質(zhì)條件

研究區(qū)地下水屬松散巖類孔隙水，含水層由第四系上更新統(tǒng)和全新統(tǒng)的細(xì)砂、中砂、砂礫石層構(gòu)成。研究區(qū)雙洎河河谷區(qū)的含水層巖性為砂、砂礫石，厚度為15～20 m，地下水水位埋深淺。雙洎河河谷區(qū)的外圍地帶和雙洎河北部支流的河谷地帶，含水層巖性主要為粉土，局部夾少量薄層粉細(xì)砂和卵礫石，厚度為20～30 m，地下水水位埋深為25 m左右；研究區(qū)內(nèi)的其他地區(qū)地下水水位埋深大，含水層巖性為黃土，夾少量的砂層，富水性差。研究區(qū)的淺層含水巖組主要由黃河多次遷徒、改道、泛濫沖積而成，在含水層的顆粒粗細(xì)、埋藏分布、滲透性能等方面具有明顯的分帶性；含水層大體由西南向東北方向延伸，呈單層出現(xiàn)，局部地段夾有粉土或粉質(zhì)黏土薄層透鏡體，厚度為20～40 m，個(gè)別地段厚度小于20 m；含水層頂板埋藏深度為5～50 m不等,巖性以粉土、粉砂為主，從上游至下游顆粒稍有變細(xì)，整體厚度約為60 m，屬潛水含水層，下伏粉質(zhì)黏土、黏土隔水層。

研究區(qū)淺層地下水主要補(bǔ)給來源是大氣降水的入滲補(bǔ)給，此外雙洎河的側(cè)向滲漏也直接補(bǔ)給該區(qū)淺層地下水；淺層地下水的徑流方向大致與地形變化趨勢(shì)一致，總體方向由西南、西北山前向東、東南方向徑流，淺層地下水的排泄方式主要有人工開采、越流排泄和徑流排泄等，地下水埋深基本都大于4 m，蒸發(fā)排泄量很小，可以忽略。

2.1.2 研究區(qū)水文地質(zhì)概念模型的建立

根據(jù)研究區(qū)的水文地質(zhì)條件，構(gòu)建了研究區(qū)的水文地質(zhì)概念模型，見圖2。

圖2 研究區(qū)水文地質(zhì)概念模型圖Fig.2 Conceptual model map of hydrogeology in the study area

(1) 模擬范圍：研究區(qū)模擬范圍為40 km2，由于該地區(qū)淺層地下水含水層的埋深約60 m，結(jié)合鉆孔資料，模擬垂直深度為60 m。

(2) 邊界條件：①水平邊界，即模擬區(qū)西與南左側(cè)邊界為定流量補(bǔ)給邊界，東、北與南右側(cè)邊界為定流量排泄邊界；②垂直邊界，模擬區(qū)垂向地下水補(bǔ)給包括大氣降水入滲補(bǔ)給、河流側(cè)向入滲補(bǔ)給、灌溉入滲補(bǔ)給；地下水排泄主要為側(cè)向徑流排泄與人工開采。

(3) 含水層結(jié)構(gòu)特征：含水層大體由西南向東北方向或由西向東延伸，本次模型含水層為層1(粉細(xì)砂)，隔水層為層2(粉質(zhì)黏土、黏土)。

研究區(qū)水文地質(zhì)概念模型的水力特征參數(shù)見圖3，污染源為某藥業(yè)污水處理廠管道破裂滲漏形成，為一線形連續(xù)性污染源。

圖3 研究區(qū)水文地質(zhì)概念模型的水力特征參數(shù)圖 (負(fù)值表示排泄)Fig.3 Hydraulic characteristic parameters of the study area hydrogeological conceptual model (negative values represent excretion)

2.2 數(shù)學(xué)模型的構(gòu)建

有機(jī)污染物在地下水中的運(yùn)移非常復(fù)雜，影響因素除對(duì)流、彌散作用以外，還存在物理、化學(xué)、微生物等作用，這些作用常常會(huì)使地下水中污染物總量減少，運(yùn)移擴(kuò)散速度減慢。本次模擬中假設(shè)污染物質(zhì)在地下水運(yùn)移過程中不與含水層介質(zhì)發(fā)生反應(yīng)，可以被認(rèn)為是保守型污染物質(zhì)。所以，本文選擇地下水流模型和溶質(zhì)運(yùn)移模型兩個(gè)模型，聯(lián)合求解水流方程和溶質(zhì)運(yùn)移方程，即可得到污染物質(zhì)的空間分布。

2.2.1 地下水流模型

三維、非均質(zhì)、各向異性的層流、非穩(wěn)定潛層地下水流數(shù)學(xué)模型可表示為

式中：Kx、Ky、Kz分別為含水層在x、y、z方向的滲透系數(shù)(m/d)；Kxx為Kx在x方向上的投影；Kyy為Ky在y方向上的投影；Kzz為Kz在z方向上的投影;x，y，z方向?yàn)楦飨虍愋越橘|(zhì)滲透系數(shù)的主方向;h為地下水水位(m);h0為地下水初始水位(m)；Kn為邊界上垂直于等水位線方向的滲透系數(shù)(m/d);Ω為研究區(qū)范圍；Γ2為邊界；w為降雨入滲補(bǔ)給量(m3/a)。

2.2.2 溶質(zhì)運(yùn)移模型

不考慮污染物在含水層中的吸附、交換、揮發(fā)、生物化學(xué)反應(yīng)，地下水中溶質(zhì)運(yùn)移的數(shù)學(xué)模型可表示為

式中：α為含水層的彌散度(m)；vi，vj分別為i和j方向上的速度分量(m/d)；|v|為速度模；C為模擬污染物的濃度(mg/L)；ne為有效孔隙度；i,j為1，2，3，x1=x，x2=y，x3=z；t為模擬時(shí)間(d)；C′為源匯的污染物濃度(mg/L)；W為源匯單位面積上的通量；Dij為彌散系數(shù)(m2/d)。

2.3 模型識(shí)別與參數(shù)確定

2.3.1 模擬流場(chǎng)及初始條件

本次模擬以2016年5月的地下水流場(chǎng)作為初始流場(chǎng)(見圖4)，以2017年5月統(tǒng)測(cè)的地下水流場(chǎng)作為模擬流場(chǎng)。

圖4 模型初始地下水流場(chǎng)Fig.4 Initial groundwater flow field of the model

2.3.2 模擬區(qū)網(wǎng)格剖分

模擬區(qū)網(wǎng)格剖分單元格為50 m×50 m，其網(wǎng)格剖分見圖5。

圖5 模擬區(qū)網(wǎng)格剖分圖Fig.5 Mesh subdivision map of simulation area

2.3.3 模型識(shí)別與參數(shù)確定

模型的識(shí)別與驗(yàn)證是整個(gè)模擬中極為重要的工作，通常需要在反復(fù)地調(diào)整參數(shù)和某些源匯項(xiàng)的基礎(chǔ)上才能達(dá)到較為理想的擬合結(jié)果[10]。本次模型識(shí)別與驗(yàn)證過程采用試估-校正法，該法屬于反求參數(shù)的間接方法之一[11]。本次模擬時(shí)間為2018年1月到2022年12月，每個(gè)時(shí)間段內(nèi)包括若干時(shí)間步長，時(shí)間步長由模型自動(dòng)控制，嚴(yán)格控制每次迭代的誤差。模型最終識(shí)別的研究區(qū)含水層水文地質(zhì)參數(shù)分區(qū)見圖6，模型識(shí)別的降雨入滲系數(shù)、含水層滲透系數(shù)和給水度見表1，灌溉入滲系數(shù)為0.10，該區(qū)地下水水位埋深普遍大于5 m，潛水蒸發(fā)排泄小可忽略不計(jì)。

圖6 研究區(qū)含水層水文地質(zhì)參數(shù)分區(qū)圖Fig.6 Zoning maps of hydrogeological parameters of aquifers in the study area

編號(hào)降雨入滲系數(shù)水平方向滲透系數(shù)/(m·d-1)給水度層1層2層1層210.150.80.008 60.0470.03420.202.10.001 00.0470.03430.151.00.001 00.0470.03440.153.50.001 00.0340.03450.101.00.001 00.0340.03460.121.00.001 00.0400.03470.101.50.008 60.0470.034

研究區(qū)不同巖性介質(zhì)彌散度的初值依據(jù)前人研究成果及經(jīng)驗(yàn)值確定，最終確定研究區(qū)不同巖性介質(zhì)的平均縱、橫向彌散度，見表2。

表2 研究區(qū)不同巖性介質(zhì)彌散度經(jīng)驗(yàn)參數(shù)表

3 模擬結(jié)果與分析

污染源位置為遂成藥業(yè)污水處理廠，假定污水管道出現(xiàn)長10 m、寬2 cm的裂縫，管道天然基礎(chǔ)層的滲透系數(shù)取值為0.047 m/d，則滲漏量約為10×0.02×0.047×1 000=9.4 kg/d，管道在無防滲設(shè)置的情況下發(fā)生滲漏，即為連續(xù)的線性污染源。根據(jù)污水中主要污染物監(jiān)測(cè)指標(biāo)，考慮污染物與巖土之間的反應(yīng)和自身降解性，選取CODCr為地下水污染預(yù)測(cè)因子，根據(jù)遂成藥業(yè)污水處理廠進(jìn)水水質(zhì)狀況，進(jìn)水中CODCr濃度為419 mg/L。

本文運(yùn)用GMS軟件對(duì)擬建的遂成藥業(yè)污水處理廠廢水排放過程中污水管道發(fā)生滲漏可能造成的地下水污染范圍進(jìn)行了數(shù)值模擬預(yù)測(cè)研究，模型運(yùn)行時(shí)間起始于2018年1月1日，模擬預(yù)測(cè)了100 d、730 d 和1 825 d(2022年12月31日)三個(gè)時(shí)間段內(nèi)研究區(qū)地下水中污染物CODCr的最大運(yùn)移距離和污染暈的運(yùn)移分布情況，其預(yù)測(cè)結(jié)果見表3和圖7。

污染組分的擴(kuò)散具有彌散特征，初期是以污染源為中心向周圍彌散，之后隨地下水的流向污染暈的主軸方向與其趨向一致(向南運(yùn)移)，所以，污染組分CODCr的檢出范圍以整個(gè)擴(kuò)散范圍進(jìn)行統(tǒng)計(jì)，最大運(yùn)移距離以污染源(見圖7中顏色最深的中心點(diǎn))為起點(diǎn)向東南方向的運(yùn)移距離為準(zhǔn)。預(yù)測(cè)結(jié)果表明：該制藥廠廢水排放過程中管道滲漏發(fā)生100 d后，地下水中污染物CODCr的檢出范圍為148 520 m2，超標(biāo)范圍為0 m2，地下水中污染物CODCr的最大運(yùn)移距離為238 m；滲漏發(fā)生730 d后，地下水中污染物CODCr的檢出范圍為165 500 m2,超標(biāo)范圍為0 m2，地下水中污染物CODCr的最大運(yùn)移距離為250 m；滲漏發(fā)生1 825 d后，地下水中污染物CODCr的檢出范圍224 000 m2，超標(biāo)范圍為0 m2，地下水中污染物CODCr的最大運(yùn)移距離為362 m(見表3)。受該地區(qū)西側(cè)與西南側(cè)兩個(gè)地下水降落漏斗的影響，該污染物CODcr滲漏后傾向于西南側(cè)運(yùn)移且運(yùn)移速度緩慢，至規(guī)劃年2022年末距東側(cè)地表水體距離較遠(yuǎn)，對(duì)周邊水體的影響較小。

表3 研究區(qū)地下水中污染物CODCr的預(yù)測(cè)結(jié)果

圖7 研究區(qū)地下水中污染物CODCr的污染暈運(yùn)移分布 預(yù)測(cè)圖Fig.7 Prediction map of CODCr pollution halo transport in groundwater in the study area

由圖7可見，研究區(qū)地下水中污染物CODCr的污染暈隨地下水的流向而擴(kuò)散，其擴(kuò)散方向與地下水的流向一致；地下水中污染物CODCr的運(yùn)移距離隨時(shí)間的增加而增加，但整體上看在預(yù)測(cè)期內(nèi)地下水中污染物CODCr的運(yùn)移距離仍處在該制藥廠的建設(shè)區(qū)范圍內(nèi)，不會(huì)對(duì)區(qū)域地下水造成明顯的污染。

4 結(jié)論與建議

(1) 制藥廠建設(shè)所排放的廢水含有復(fù)雜且大量的有機(jī)污染物，可能會(huì)對(duì)廠區(qū)范圍內(nèi)的地下水造成污染，所以必須在制藥廠建設(shè)前對(duì)廢水排放可能造成的地下水污染范圍進(jìn)行預(yù)測(cè)研究。

(2) 本文運(yùn)用GMS軟件對(duì)擬建的遂成藥業(yè)污水處理廠廢水排放過程中污水管道發(fā)生滲漏可能造成的地下水污染范圍進(jìn)行了數(shù)值模擬研究，模型運(yùn)行時(shí)間起始于2018年1月1日，模擬預(yù)測(cè)了100 d、730 d和1 825 d(2022年12月31日)三個(gè)時(shí)間段內(nèi)研究區(qū)地下水中污染物CODCr的最大運(yùn)移距離。模擬結(jié)果表明：滲漏發(fā)生100 d后，地下水中污染物CODCr的最大運(yùn)移距離為238 m；滲漏發(fā)生730 d后，地下水中污染物CODCr的最大運(yùn)移距離為250 m；滲漏發(fā)生1 825 d后，地下水中污染物CODCr的最大運(yùn)移距離為362 m。

(3) 研究區(qū)地下水中污染物CODCr的污染暈隨地下水的流向而擴(kuò)散，其擴(kuò)散方向與地下水的流向一致；地下水中污染物CODCr的運(yùn)移距離隨時(shí)間的增加而增加，但整體上看，在預(yù)測(cè)期內(nèi)不會(huì)對(duì)區(qū)域地下水造成明顯的污染。但是，鑒于該制藥廠污染物的極大危害性，建議在廠區(qū)范圍內(nèi)及其周邊設(shè)置監(jiān)測(cè)井，定期檢測(cè)地下水中污染組分的變化，一旦發(fā)現(xiàn)地下水受到污染，應(yīng)立即采取相應(yīng)的措施進(jìn)行治理。