岳松梅，楊 蘊(yùn)，吳劍鋒，吳吉春

(1.表生地球化學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，南京大學(xué)地球科學(xué)與工程學(xué)院，江蘇 南京 210023;2.河海大學(xué)地球科學(xué)與工程學(xué)院，江蘇 南京210098)

研究大尺度地下水流、污染物濃度和運(yùn)動(dòng)規(guī)律問(wèn)題時(shí)，區(qū)域性模型構(gòu)建中的一個(gè)突出問(wèn)題是勘探資料分布不均或缺乏，以及含水系統(tǒng)本身的非均質(zhì)性，造成了水文地質(zhì)參數(shù)的空間變異性。在通常的地下水?dāng)?shù)值模擬中只能用簡(jiǎn)單的參數(shù)分區(qū)來(lái)描述參數(shù)的非均質(zhì)性，但由于分區(qū)過(guò)大而忽略分區(qū)內(nèi)部參數(shù)的差異性，往往導(dǎo)致模擬結(jié)果的不確定性。尋求一種盡可能利用有限的勘探資料，對(duì)未知區(qū)域內(nèi)的含水層參數(shù)進(jìn)行合理估值的方法，是目前大區(qū)域地下水流模擬中的關(guān)鍵問(wèn)題［1－2］。

在影響地下水系統(tǒng)的所有不確定因素中，滲透系數(shù)是表征含水層特性的一個(gè)關(guān)鍵要素［3－5］，它在空間上的分布變化及其復(fù)雜，其不確定性主要是由自然界含水介質(zhì)的非均質(zhì)各向異性和取樣、測(cè)試中的失真以及試驗(yàn)量測(cè)誤差造成的。這種空間分布的不確定性使我們將其以隨機(jī)變量來(lái)處理，在不可能測(cè)得含水層每一點(diǎn)的介質(zhì)特征和水力特征的前提下，引入含水層空間隨機(jī)場(chǎng)，將實(shí)際含水層視為隨機(jī)場(chǎng)的一個(gè)實(shí)現(xiàn)。但實(shí)際情況中，滲透系數(shù)不僅具有隨機(jī)性，也具有一定的結(jié)構(gòu)性［6－9］。地質(zhì)統(tǒng)計(jì)學(xué)就是研究這種具有“二重性”的區(qū)域化變量的數(shù)學(xué)工具，其以變差函數(shù)為基本工具來(lái)研究分布于空間并呈現(xiàn)一定結(jié)構(gòu)性與隨機(jī)性的自然現(xiàn)象［10－11］，因而用地質(zhì)統(tǒng)計(jì)學(xué)方法來(lái)研究水文水資源水環(huán)境系統(tǒng)參變量的空間變異性，可以定量地揭示水文水資源水環(huán)境系統(tǒng)參變量在空間不同方向上的變化規(guī)律，辨識(shí)各參數(shù)的最大變異方向，尤其是它能給出區(qū)域化變量的最優(yōu)估計(jì)值及估計(jì)方差［10－11］，而這種特性是其他任何一種估計(jì)方法所不曾擁有的。

本文在前人研究的基礎(chǔ)上［12］，將地質(zhì)統(tǒng)計(jì)學(xué)推廣到大尺度的地下水污染運(yùn)移的研究工作中，結(jié)合GSLIB軟件庫(kù)［18］，利用美國(guó)麻省軍事保留區(qū)(Massachusetts Military Reservation，MMR)的場(chǎng)地實(shí)測(cè)數(shù)據(jù)，采用普通克里格法和指示克里格法、順序高斯模擬法和順序指示模擬法等四種地質(zhì)統(tǒng)計(jì)方法，插值估測(cè)和模擬再現(xiàn)含水層滲透系數(shù)隨機(jī)場(chǎng)，進(jìn)而對(duì)比研究四種滲透系數(shù)場(chǎng)對(duì)大尺度污染物運(yùn)移的影響。

1 四種地質(zhì)統(tǒng)計(jì)方法簡(jiǎn)述

1.1 普通克里格和指示克里格方法

普通克里格法屬線性平穩(wěn)地質(zhì)統(tǒng)計(jì)學(xué)范疇，是地質(zhì)統(tǒng)計(jì)學(xué)中最基本的方法，也是最常用的方法之一。設(shè)區(qū)域化變量Z(u)滿足二階平穩(wěn)假設(shè)，其數(shù)學(xué)期望為未知常數(shù)，協(xié)方差函數(shù)和變異函數(shù)存在且平穩(wěn)，以任一未知點(diǎn)為中心點(diǎn)的某一鄰域范圍內(nèi)n個(gè)樣本點(diǎn)，賦予相應(yīng)權(quán)重計(jì)算未知點(diǎn)的屬性值，形成估計(jì)領(lǐng)域內(nèi)n個(gè)信息值的線性組合。通過(guò)求得一個(gè)無(wú)偏的最優(yōu)估計(jì)，并給出每個(gè)估值的誤差方差，該方法即為普通克里格方法(ordinary kriging，OK)［10－11］。

然而，在實(shí)際研究中還會(huì)碰到采樣數(shù)據(jù)中存在特異值得問(wèn)題，所謂特異值是指那些比全部數(shù)值的均值或中位數(shù)高很多的值，其既非分析誤差所致，也非采樣方法等人為誤差引起，指示克里格法(indicator kriging，IK)就是為解決上述問(wèn)題而發(fā)展起來(lái)的一種非參數(shù)地統(tǒng)計(jì)學(xué)法。該方法在給定條件數(shù)據(jù)的前下可以估計(jì)任意位置上的條件累積分布函數(shù)(conditional cumulative distribution function，CCDF)，適用于連續(xù)型和離散型的信息預(yù)測(cè)［10－11］。

相比于普通克里格，指示克里格法是將對(duì)區(qū)域化變量Z(u)的研究轉(zhuǎn)化為對(duì)其指示函數(shù)的研究，它無(wú)需假設(shè)數(shù)值來(lái)自某種特定分布的總體，也無(wú)需對(duì)原始數(shù)據(jù)進(jìn)行變換，可以估計(jì)落于某些區(qū)間的每一部分的值，也可以估計(jì)比某一門坎值高(或低)的那部分值。

1.2 順序高斯模擬和順序指示模擬方法

順序高斯模擬(sequential Gaussian simulation，SGSIM)也稱序貫高斯模擬，是將順序模擬思想和高斯模擬相結(jié)合的一種隨機(jī)模擬方法。它假設(shè)一個(gè)高斯隨機(jī)域，基于變異函數(shù)特征，對(duì)每一個(gè)待模擬的空間位置利用克里格方法求構(gòu)建高斯分布函數(shù)，之后隨機(jī)在其中抽取一個(gè)數(shù)值作為該點(diǎn)的模擬值。

順序指示模擬(sequential indicator simulation，SISIM)將指示克里格方法和順序模擬算法相結(jié)合，實(shí)現(xiàn)了非參數(shù)的連續(xù)型與離散型的分布模擬。該方法不需要對(duì)原始條件數(shù)據(jù)分布的參數(shù)形式做任何假設(shè)，而是在現(xiàn)有資料的基礎(chǔ)上，通過(guò)一系列的門檻值把條件數(shù)據(jù)轉(zhuǎn)化成指示數(shù)據(jù)，相應(yīng)的變差函數(shù)模擬轉(zhuǎn)換為指示變差函數(shù)模型，并采用指示克里格法對(duì)每個(gè)網(wǎng)格點(diǎn)處的條件概率分布進(jìn)行估計(jì)［14－16］。

2 污染物運(yùn)移及其結(jié)果評(píng)估

研究地下水污染問(wèn)題必須用兩類方程來(lái)描述:第一類方程用來(lái)描述液體的流動(dòng);第二類方程用來(lái)描述溶質(zhì)即污染物質(zhì)在地下水中的運(yùn)移。在一般或者較復(fù)雜的水文地質(zhì)條件下，污染物運(yùn)移只有通過(guò)數(shù)值方法來(lái)實(shí)現(xiàn)，而且無(wú)論是采用有限單元或者有限差分方法，都要通過(guò)地下水流模擬地下水的流速場(chǎng)才能由運(yùn)移模型得到污染物運(yùn)移的模擬結(jié)果。在這里，分別采用常用的軟件水流模擬模型MODFLOW和溶質(zhì)運(yùn)移模型MT3DMS來(lái)進(jìn)行模擬。

2.1 滲透系數(shù)場(chǎng)的產(chǎn)生

2.1.1 原始數(shù)據(jù)的基本統(tǒng)計(jì)特征

包括均值、方差、標(biāo)準(zhǔn)差、變異系數(shù)、偏態(tài)系數(shù)和峰態(tài)系數(shù)的計(jì)算與分析;畫出原始數(shù)據(jù)分布直方圖并分析數(shù)據(jù)是否滿足正態(tài)分布特征。

根據(jù)原始數(shù)據(jù)分布直方圖，依據(jù)“偏度、峰度檢驗(yàn)法”對(duì)原始數(shù)據(jù)進(jìn)行正態(tài)分布檢驗(yàn)，進(jìn)而判斷原始數(shù)據(jù)的正態(tài)性。

2.1.2 原始數(shù)據(jù)的穩(wěn)健性處理

地質(zhì)統(tǒng)計(jì)學(xué)必須滿足數(shù)據(jù)為正態(tài)分布這一要求，如不滿足正態(tài)分布，直接計(jì)算的傳統(tǒng)地質(zhì)統(tǒng)計(jì)學(xué)變差函數(shù)必然有偏，進(jìn)而影響估計(jì)的效果，故需對(duì)數(shù)據(jù)進(jìn)行穩(wěn)健性處理。分別對(duì)數(shù)據(jù)先后進(jìn)行去類分析、特異值處理和正態(tài)變換，使數(shù)據(jù)接近正態(tài)分布。

2.1.3 滲透系數(shù)分布場(chǎng)的產(chǎn)生

本文利用GSLIB軟件庫(kù)，對(duì)滲透系數(shù)進(jìn)行原始數(shù)據(jù)基本參數(shù)統(tǒng)計(jì)和數(shù)據(jù)穩(wěn)健性處理，然后采用普通克里格法和指示克里格法、高斯序列模擬法和指示序列模擬法分別對(duì)數(shù)據(jù)進(jìn)行插值和條件模擬，得到相應(yīng)的滲透系數(shù)(lnK)分布場(chǎng)。

2.2 污染羽的空間矩評(píng)估

污染物在空間上的變化特征通常可以通過(guò)刻畫污染羽隨時(shí)間變化的空間矩來(lái)評(píng)價(jià)。污染羽在空間上的零階矩表示污染物的總質(zhì)量;一階矩表示污染羽的質(zhì)心位置;二階矩表示污染羽圍繞質(zhì)心的分布范圍。

污染羽在任一時(shí)刻t關(guān)于原點(diǎn)的零階矩，即污染物總質(zhì)量，可以表示為

污染羽在空間上關(guān)于質(zhì)心的二階矩可以表示為協(xié)方差變量U:

本文以MMR為實(shí)例，建立水流模型和溶質(zhì)運(yùn)移模型，利用空間矩評(píng)估污染羽的變化范圍。由于模型考慮化學(xué)反應(yīng)的兩種基本類型:固液表面反應(yīng)(吸附作用)和一階速率反應(yīng)，污染物的質(zhì)量(零階矩)在運(yùn)移過(guò)程中會(huì)發(fā)生改變，因此我們考慮污染羽的零階矩、一階矩和二階矩隨滲透系數(shù)場(chǎng)的空間變異特征。

3 實(shí)例研究

3.1 算例概況

美國(guó)麻省軍事保護(hù)區(qū)位于科德角(Cape Cod)的法爾茅斯鎮(zhèn)(Town of Falmouth)附近，始建于20世紀(jì)初，面積大約為89 km2。該區(qū)為第四紀(jì)冰川沉積物覆蓋，這些沉積物組成的非承壓含水層為科德角幾個(gè)社區(qū)的主要供水水源。位于MMR的東南部的美國(guó)空軍的給養(yǎng)基地，由于長(zhǎng)期的軍事活動(dòng)，出現(xiàn)了以三氯乙烯(Trichloroethylene，TCE)為主要污染物的化學(xué)物溢漏區(qū)(Chemical－Spill 10，CS－10)，對(duì)整個(gè)區(qū)域的環(huán)境造成嚴(yán)重破壞，直接影響了當(dāng)?shù)鼐用裾５纳a(chǎn)生活。因此，本次研究就以CS－10區(qū)為模擬計(jì)算區(qū)域(圖1)。

圖1 研究區(qū)位置、CS－10模擬區(qū)域及TCE污染治理系統(tǒng)的平面分布［19］

3.2 模型的建立

本文建立了描述研究區(qū)含水層中TCE污染物運(yùn)移的數(shù)學(xué)模型，采用有限差分法求解。其中只有含水層滲透系數(shù)為一隨機(jī)變量，其他模型參數(shù)均為常數(shù)。水流模型中的滲透系數(shù)K是通過(guò)地質(zhì)統(tǒng)計(jì)學(xué)對(duì)實(shí)測(cè)的值進(jìn)行插值和條件模擬得到的，變化范圍較大，在泥沙含水層中可以小于3 m/d，粗砂含水層中能大于91 m/d。在假定研究區(qū)水流呈穩(wěn)定態(tài)，而溶質(zhì)運(yùn)移呈非穩(wěn)定態(tài)的基礎(chǔ)上，整個(gè)研究區(qū)域剖分為159列、161行和21層，總模擬面積為57 km2，總模擬時(shí)間為30 a，1個(gè)應(yīng)力期。水平方向步長(zhǎng)為34m，越接近邊界步長(zhǎng)會(huì)越大，垂直方向上每一層的厚度變化較大，最小可以小于1.5 m，最大可以大于15 m。側(cè)向邊界為定水頭邊界，由區(qū)域水流模型計(jì)算然后插值得到;上部為定流量邊界，補(bǔ)給量變化范圍為41～86 cm/a;下部隔水為零通量邊界。該區(qū)內(nèi)的地下水水流主要是沿水平方向從南向西南流動(dòng)，水平的水力坡度平均值為0.001。

溶質(zhì)運(yùn)移模型是建立在MT3DMS的基礎(chǔ)上的。該模型采用和MODFLOW水流模型相同的網(wǎng)格剖分，側(cè)向和垂向上都是定通量邊界。在邊界上，不考慮彌散項(xiàng)而只考慮對(duì)流項(xiàng)，如果流量為0，則溶質(zhì)運(yùn)移也等于0。其余參數(shù)如，有效孔隙度為0.3，橫向彌散度、縱向彌散度和垂向彌散度分別為11 m、1.1 m 和 0.11 m［19］。TCE 濃度超過(guò) 5μgL－1的污染羽三維分布圖見(jiàn)圖2所示，

圖2 溶質(zhì)運(yùn)移模型運(yùn)行前CS－10 TCE污染羽濃度的三維分布圖

污染物在含水層中運(yùn)移前初始污染物濃度均值為16 μg/L，最高濃度值遠(yuǎn)超過(guò) 4 690 μg/L［13］。在美國(guó)，建立在保護(hù)人類健康和環(huán)境的安全飲用水法規(guī)的條件下，TCE的最大污染標(biāo)準(zhǔn)(MCL)是5 μg/L，這說(shuō)明研究區(qū)域被高濃度的TCE污染。

3.3 滲透系數(shù)分布場(chǎng)的產(chǎn)生

在對(duì)研究區(qū)已知的375個(gè)觀測(cè)數(shù)據(jù)進(jìn)行穩(wěn)健統(tǒng)計(jì)分析，并在此基礎(chǔ)上構(gòu)造該區(qū)域化變量的幾何各向異性套合模型，最后與GSLIB數(shù)據(jù)庫(kù)相結(jié)合，采用四種不同的地質(zhì)統(tǒng)計(jì)學(xué)方法分別對(duì)數(shù)據(jù)進(jìn)行插值和條件模擬，得到相應(yīng)的滲透系數(shù)(lnK)分布場(chǎng)。對(duì)于克里格插值來(lái)說(shuō)，一旦各種參數(shù)確定(變差函數(shù)模型)以后，滲透系數(shù)(lnK)分布場(chǎng)就確定了，對(duì)應(yīng)的污染物運(yùn)移結(jié)果也確定了;而對(duì)于條件模擬來(lái)說(shuō)，則每一次產(chǎn)生的滲透系數(shù)(lnK)分布場(chǎng)都是隨機(jī)的，相應(yīng)的污染物運(yùn)移結(jié)果也是隨機(jī)的。在這里，分別對(duì)克里格插值和條件模擬得出的結(jié)果做出對(duì)比，條件模擬模擬次數(shù)采用300次。

3.4 模擬結(jié)果與分析討論

普通克里格和指示克里格法所得的滲透系數(shù)(lnK)均值在3.7左右，范圍在－3～7之間，標(biāo)準(zhǔn)偏差在1左右，但指示克里格值(0.607 6)比普通格里格值(1.104 9)偏小。四種地質(zhì)統(tǒng)計(jì)方法得到的滲透系數(shù)場(chǎng)(lnK)估計(jì)方差的統(tǒng)計(jì)參數(shù)詳見(jiàn)表1。

表1 四種地質(zhì)統(tǒng)計(jì)方法得到滲透系數(shù)場(chǎng)(lnK)估計(jì)方差的統(tǒng)計(jì)參數(shù)

由表1可以得出，兩種克里格插值估計(jì)方差均值在1左右，指示克里格法略小于普通克里格法，但標(biāo)準(zhǔn)偏差比普通克里格法略大。圖3分別為污染運(yùn)移5 475 d(15 a)和10 950 d(30 a)后，克里格插值下污染羽濃度的三維分布圖。

圖3 克里格插值生成滲透系數(shù)場(chǎng)所對(duì)應(yīng)的15 a和30 a后污染羽分布圖(μgL－1)

從圖3中不難看出，普通克里格法的污染羽擴(kuò)散更加均勻，指示克里格所對(duì)應(yīng)的污染物濃度極值范圍更大。對(duì)于污染羽的二階空間矩(表2)評(píng)估，普通克里格插值對(duì)應(yīng)的污染羽在空間上的展布范圍(二階矩)明顯大于和指示克里格，這主要受滲透系數(shù)的空間變異方差的影響，滲透系數(shù)方差的偏大，污染羽在空間上的二階矩隨之升高，表明污染物圍繞中心的擴(kuò)散范圍不斷加大。這些不同之處是由于指示克里格法是用個(gè)別范圍個(gè)別變異函數(shù)來(lái)刻畫該區(qū)域化變量的關(guān)聯(lián)性的特性，無(wú)需假設(shè)數(shù)值來(lái)自某種特定分布(如正態(tài)分布)的總體，也無(wú)需對(duì)原始數(shù)據(jù)進(jìn)行變換(如對(duì)數(shù)變換)，因此指示克里格法不必去掉重要而實(shí)際存在的高值數(shù)據(jù)即可處理各種不同現(xiàn)象［12］，而普通克里格法表現(xiàn)得更圓滑效應(yīng)，滲透系數(shù)估計(jì)方差大且平滑，進(jìn)而污染羽擴(kuò)散更加均勻。

估計(jì)法是用加權(quán)平均法對(duì)研究區(qū)域的未知量求得線性、無(wú)偏和最小方差的內(nèi)插估計(jì)量及其相應(yīng)的估計(jì)精度，而順序模擬法是通過(guò)系列隨機(jī)模擬對(duì)地質(zhì)變量進(jìn)行局部估計(jì)，克服了估計(jì)法的平滑，能夠更好的再現(xiàn)真實(shí)曲線的波動(dòng)性。對(duì)于條件模擬每一次產(chǎn)生的滲透系數(shù)(lnK)分布場(chǎng)都是隨機(jī)的，相應(yīng)的污染物運(yùn)移結(jié)果也是隨機(jī)的，這里采用300次模擬次數(shù)。

順序高斯模擬和順序指示模擬(模擬次數(shù)300次)所得的滲透系數(shù)(lnK)均值相同，即3.7997，范圍均在－3～7之間，但標(biāo)準(zhǔn)偏差順序指示模擬值(0.5941)比順序高斯模擬值(0.758 2)偏小。這兩種方法得到滲透系數(shù)場(chǎng)(lnK)估計(jì)方差均值在1左右，標(biāo)準(zhǔn)偏差順序高斯模擬值略小，詳見(jiàn)表1。

對(duì)比克里格插值得出的結(jié)果，四種地質(zhì)統(tǒng)計(jì)法所得到的滲透系數(shù)(lnK)均值在3.7左右，所對(duì)應(yīng)污染羽擴(kuò)散范圍基本一致，污染羽的質(zhì)心坐標(biāo)也基本相同，驗(yàn)證了質(zhì)心位置是由滲透系數(shù)的平均值來(lái)決定的理論。圖4分別為污染運(yùn)移5 475 d(15 a)和10 950 d(30 a)后，條件模擬(模擬次數(shù)300次)下污染羽的濃度分布。表3為條件模擬模擬次數(shù)300次所對(duì)應(yīng)的15 a和30 a后污染羽在空間上的分布特征。

綜合圖4和表3，不難看出，順序指示模擬的高濃度范圍明顯高于順序高斯模擬法，但污染羽的二階空間矩略小于后者。

表2 克里格插值方法所對(duì)應(yīng)的15 a和30 a后污染羽在空間上的分布特征

圖4 條件模擬次數(shù)300次模擬生成滲透系數(shù)場(chǎng)所對(duì)應(yīng)的15 a和30 a后污染羽分布圖

兩種模擬方法在模擬次數(shù)增加的情況下，污染羽的擴(kuò)散范圍與質(zhì)心坐標(biāo)沒(méi)有明顯差別，這有力的驗(yàn)證了質(zhì)心位置是由滲透系數(shù)的平均值來(lái)決定的理論;但污染羽的空間二階矩在某一方向上卻有明顯減小的現(xiàn)象，這與二階空間矩主要受滲透系數(shù)的空間變異方差影響的結(jié)論相呼應(yīng)。單獨(dú)對(duì)比順序高斯模擬與順序指示模擬這兩種方法，在模擬次數(shù)增加的情況下，滲透系數(shù)估計(jì)方差與污染物空間矩評(píng)估值這兩方面，后者的變化程度明顯高于前者。這是由于順序指示模擬采用指示克里金來(lái)估計(jì)，指示克里金不同于其他克里金法，它并不依賴于變量的平穩(wěn)性和要求變量服從某種分布。順序指示模擬法可以模擬復(fù)雜各向異性的地質(zhì)現(xiàn)象及連續(xù)性分布的極值，對(duì)于具有不同連續(xù)性分布的變量，可利用不同的變差函數(shù)進(jìn)行表征，建立各向異性的模擬結(jié)果，進(jìn)而更能體現(xiàn)數(shù)據(jù)的差異性。所以，采用順序指示模擬方法模擬的污染羽分布，其高濃度范圍明顯高于順序高斯模擬得出的結(jié)果，并且空間二階矩的變化程度明顯更大。

表3 條件模擬模擬次數(shù)300次所對(duì)應(yīng)的15 a和30 a后污染羽在空間上的分布特征

4 結(jié)語(yǔ)

(1)本文利用野外場(chǎng)地實(shí)測(cè)數(shù)據(jù)，采用四種地質(zhì)統(tǒng)計(jì)方法，插值估測(cè)和模擬再現(xiàn)隨機(jī)滲透系數(shù)場(chǎng)，進(jìn)而對(duì)比研究四種滲透系數(shù)場(chǎng)對(duì)大尺度污染物運(yùn)移的影響。研究結(jié)果表明，污染羽一階空間矩證明了污染羽的質(zhì)心坐標(biāo)基本相同，而四種地質(zhì)統(tǒng)計(jì)法所得到的滲透系數(shù)(lnK)范圍均在－3～7之間，均值在3.7左右，這有力的驗(yàn)證了質(zhì)心位置是由滲透系數(shù)的平均值來(lái)決定的理論;通過(guò)指示克里格和順序指示模擬法得到的滲透系數(shù)(lnK)標(biāo)準(zhǔn)偏差比其余兩種更小，污染羽空間二階矩偏小，這也與二階矩主要受滲透系數(shù)的空間變異方差影響的結(jié)論相呼應(yīng)。

(2)四種地質(zhì)統(tǒng)計(jì)方法在處理樣本數(shù)據(jù)時(shí)各有利有弊，應(yīng)根據(jù)樣本的特點(diǎn)來(lái)作出相應(yīng)對(duì)策。指示克里格與普通克里格得出的結(jié)果最大不同之處在于，前者污染羽的極值特點(diǎn)更加明顯，后者污染羽擴(kuò)散更加均勻。條件模擬追求的是模擬的真實(shí)性，克服了估計(jì)法的平滑效果，能較好地再現(xiàn)真實(shí)曲線的波動(dòng)性;隨著模擬次數(shù)的增大，滲透系數(shù)(lnK)估計(jì)方差與污染羽空間二階矩呈現(xiàn)變小的趨勢(shì)，并且順序指示模擬程度更加明顯，這是由于順序指示模擬在模擬復(fù)雜各向異性的地質(zhì)現(xiàn)象及連續(xù)性分布的極值上更加突出。

(3)需要說(shuō)明的是，本文進(jìn)行的污染物運(yùn)移模擬是建立在四種地質(zhì)統(tǒng)計(jì)方法估計(jì)模擬滲透系數(shù)場(chǎng)的前提下實(shí)現(xiàn)的，不可避免得會(huì)出現(xiàn)估計(jì)誤差，除此之外，在變差函數(shù)理論模型擬合中存在必然的人為誤差和模型誤差。用數(shù)值法計(jì)算中也存在誤差，模型只能計(jì)算位于網(wǎng)格單元中心的結(jié)點(diǎn)處的濃度，而監(jiān)測(cè)孔的實(shí)際位置大多不是位于有限差分網(wǎng)格單元的中心。

［1］Kai－ Wei Juang，Dar－ Yuan Lee，Timothy R.Ellsworth.Using rank－order geostatistics for spatial interpolation of highly skewed data in a heavy － metal contaminated site［J］.Journal of Environmental Quality，2001，30:894－903.

［2］Bastante F G，Ordóez C，Taboada J，et al.Comparison of indicator kriging，conditional indicator simulation and multiple－point statistics used to model slate deposits［J］.Engineering Geology，2008，98:50－59.

［3］Diana dell＇Arciprete，Riccardo Bersezio，F(xiàn)abrizio Felletti，et al.Comparison of three geostatistical methods for hydrofacies simulation:a test on alluvial sediments［J］.Hydrogeology Journal，2012，20:299－311.

［4］Freeze R A.A stochastic conceptual analysis of one－dimensional groundwater flow in non － uniform homogeneous media［J］.Water Resources Research.1975，11(5):725－741.

［5］施小清，吳吉春，袁永生.滲透系數(shù)空間變異性研究［J］.水科學(xué)進(jìn)展.2005，16(2):210－215.

［6］陳彥，吳吉春.含水層滲透系數(shù)空間變異性對(duì)地下水?dāng)?shù)值模擬的影響［J］.水科學(xué)進(jìn)展.2005，16(4):482－48.

［7］劉文婷，朝倫巴根，劉艷偉，王力，宋君.地質(zhì)統(tǒng)計(jì)學(xué)方法在滲透系數(shù)空間變異性研究中的應(yīng)用［J］.水利科技與經(jīng)濟(jì).2010，16(4):364－366.

［8］李少龍，張家發(fā)，楊金忠.滲透系數(shù)空間相關(guān)性對(duì)滲流場(chǎng)統(tǒng)計(jì)特征影響分析［J］.長(zhǎng)江科學(xué)院院報(bào).2009，26(10):17－20.

［9］閆婷婷，吳劍鋒.滲透系數(shù)的空間變異性對(duì)污染物運(yùn)移的影響研究［J］.水科學(xué)進(jìn)展.2006，17(1):29－36.

［10］張仁鐸.空間變異理論及應(yīng)用［M］.北京:科學(xué)出版社.2004.

［11］劉愛(ài)利，王培法，丁園園.地統(tǒng)計(jì)學(xué)概論［M］.北京:科學(xué)出版社.2012.

［12］劉玲玲，吳劍鋒，吳吉春.不同地質(zhì)統(tǒng)計(jì)方法在確定滲透系數(shù)場(chǎng)中的對(duì)比研究［J］.水文地質(zhì)工程地質(zhì).2009，5:66－71.

［13］吳劍鋒，彭偉，錢家忠.基于INPGA的地下水污染治理多目標(biāo)優(yōu)化管理模型［J］.地質(zhì)論評(píng).2011，57(3):437－443.

［14］李少華，張昌民，王振奇.利用順序指示模擬方法預(yù)測(cè)儲(chǔ)集層巖性［J］.新疆石油地質(zhì).2002，23(1):59－62.

［15］宋永忠，于晶，年靜波.隨機(jī)順序指示模擬技術(shù)與應(yīng)用［J］.大慶石油地質(zhì)與開發(fā).2003，1:52－54.

［16］余振，何靜，魏福吉.序貫指示模擬和序貫高斯模擬在某地區(qū)精細(xì)流體預(yù)測(cè)中的聯(lián)合應(yīng)用［J］.天然氣地球科學(xué).2012，23(6):1170－1174.

［17］施小清，吳吉春，吳劍鋒.多個(gè)相關(guān)隨機(jī)參數(shù)的空間變異性對(duì)溶質(zhì)運(yùn)移的影響［J］.水科學(xué)進(jìn)展.2012，23(4):509－515.

［18］Deutsch C V，Journel A G.GSLIB－Geostatistcal Software Library and User’s Guide［M］.2nd ed.New York:Oxford University Press，1998.

［19］Zheng Chunmiao，Wang P P.2002.A field demonstration of the simulation－ optimization approach for remediation system design.Ground Water.40(3):258 ～265.