基于SD模型的長白山安圖縣礦泉水資源承載力研究

谷志琪，卞建民，馬于曦，孫曉慶，李一涵，李佳林

(1.吉林大學(xué)地下水資源與環(huán)境教育部重點(diǎn)實(shí)驗(yàn)室，吉林 長春 130021；2.吉林大學(xué)新能源與環(huán)境學(xué)院，吉林 長春 130021)

水資源是人類生存和社會發(fā)展不可替代的基本物質(zhì)，是生物生存最重要的自然資源之一[1].隨著人們對水質(zhì)要求的逐漸提高，具備純凈、衛(wèi)生、安全等突出優(yōu)勢的天然礦泉水日益受到廣大消費(fèi)者的歡迎[2]，因此，如何科學(xué)合理地開發(fā)利用礦泉水資源十分重要[3].

目前，國內(nèi)外學(xué)者主要從水資源所能承載的最大人口容量[4-5]、支持社會經(jīng)濟(jì)發(fā)展的能力[6-7]、承載的最大閾值[8-9]及綜合支撐社會發(fā)展的各項(xiàng)指標(biāo)[10-11]等四個(gè)方面對水資源承載力進(jìn)行了定義.據(jù)此，本文提出的礦泉水資源承載力的概念指以可持續(xù)發(fā)展為原則，在已形成一定規(guī)模的區(qū)域內(nèi)礦泉水資源量所能支持的區(qū)域內(nèi)社會體系良好發(fā)展的能力.目前，評價(jià)水資源承載力的方法主要有常規(guī)趨勢法[12-13]、主成分分析法[14-15]、模糊綜合評價(jià)法[16-17]、多目標(biāo)決策法[18-19]以及系統(tǒng)動(dòng)力學(xué)法[20-21]等等.其中，系統(tǒng)動(dòng)力學(xué)法擅長處理高階次、非線性、多變量的問題[22]，且能夠模擬不同情景下的水資源承載力，進(jìn)而確定最合適的發(fā)展方案.近年來許多學(xué)者應(yīng)用系統(tǒng)動(dòng)力學(xué)方法對水資源承載力進(jìn)行了研究[23-24]，并耦合其他方法[25-26]對水資源承載力進(jìn)行了評價(jià).

長白山是我國天然礦泉水集中分布區(qū)之一，安圖縣位于長白山腹地，礦泉水資源儲量巨大、類型多樣、品質(zhì)優(yōu)良.近年來，研究區(qū)礦泉水資源的開發(fā)力度不斷加大，有效地推動(dòng)了地區(qū)經(jīng)濟(jì)增長，但在開采過程中存在水源地環(huán)境保護(hù)壓力大、開發(fā)利用率低等問題，在一定程度上阻礙了礦泉水資源的開發(fā)利用.因此，本文利用系統(tǒng)動(dòng)力學(xué)(SD)方法建立了礦泉水資源承載力系統(tǒng)動(dòng)力學(xué)模型，模擬2015—2030年研究區(qū)礦泉水資源承載力在現(xiàn)狀型、發(fā)展型、節(jié)約型、環(huán)境型和協(xié)調(diào)型5種方案下的動(dòng)態(tài)變化，以尋求提升研究區(qū)礦泉水資源承載力的最佳方案，為礦泉水資源的合理開發(fā)、利用和保護(hù)提供理論依據(jù)和科學(xué)支撐.

1 研究區(qū)概況

安圖縣位于吉林省東部(見圖1)，區(qū)內(nèi)廣泛分布軍艦山組玄武巖，發(fā)育了大量的構(gòu)造斷裂，植被覆蓋率達(dá)87%以上；降水量波動(dòng)較大，多年平均值為691.7 mm，多年平均蒸發(fā)量為1 170 mm，多年平均氣溫4.3℃，屬于北中溫帶季風(fēng)半濕潤氣候，總體地勢南高北低、東高西低，中部夾有山間盆地.該區(qū)地表水系發(fā)達(dá)，屬于圖們江及松花江水系，包括布爾哈通河、頭道白河、二道白河等88條河流.截至2015年，研究區(qū)人口20.52萬人，其中城鎮(zhèn)人口12.18萬人、鄉(xiāng)村人口8.34萬人；水資源總量為13.75億m3，占全州水資源量的15%，人均水資源占有量約為6 700 m3.

圖1 研究區(qū)地理位置

為滿足礦泉水資源的大量需求，經(jīng)過對區(qū)內(nèi)礦泉水水源進(jìn)行勘查鑒定，已勘查評價(jià)礦泉水水源26處，礦泉總流量達(dá)1.88萬m3/d，允許開采總量為145 186.83 m3/d.目前已經(jīng)建成5家礦泉水企業(yè)進(jìn)行開發(fā)，利用率僅為3.11%，具有較大的擴(kuò)展空間.

2 數(shù)據(jù)與方法

2.1 數(shù)據(jù)來源

本文使用的社會經(jīng)濟(jì)、人口、環(huán)境數(shù)據(jù)來源于2010—2015年《吉林省統(tǒng)計(jì)年鑒》；水資源統(tǒng)計(jì)數(shù)據(jù)、城鎮(zhèn)及鄉(xiāng)村生活用水定額參照2010—2015年《延邊朝鮮族自治州水資源公報(bào)》.

2.2 研究方法

使用Vensim軟件構(gòu)建礦泉水資源承載力系統(tǒng)動(dòng)力學(xué)模型，結(jié)合層次分析法確定主要變量權(quán)重，模擬預(yù)測2015—2030年研究區(qū)礦泉水資源承載力在設(shè)定的5種方案下的動(dòng)態(tài)變化情況.

2.2.1 系統(tǒng)動(dòng)力學(xué)原理

系統(tǒng)動(dòng)力學(xué)(SD)模型由Forrester教授于1956年提出，主要通過計(jì)算機(jī)仿真技術(shù)，以定性與定量相結(jié)合的方法來解決人口、產(chǎn)業(yè)、資源環(huán)境這類復(fù)雜的系統(tǒng)問題.系統(tǒng)動(dòng)力學(xué)模型本質(zhì)上是一階的時(shí)滯微分方程組[27]，流率變量方程的一般形式為：

LK=LJ+T(RTJK).

(1)

式中：LK，LJ為K，J時(shí)刻的流率變量；T為時(shí)間間隔；RTJK為時(shí)間間隔下的流率變化向量.(1)式可變形為：

(LK-LJ)/T=RTJK.

(2)

安圖縣礦泉水資源承載力系統(tǒng)動(dòng)力學(xué)模型共涉及了38個(gè)方程、9對狀態(tài)變量和速率變量、5個(gè)表函數(shù).由于模型涉及的方程數(shù)量較多，通過線性插值、表函數(shù)等方法對研究區(qū)歷史統(tǒng)計(jì)資料進(jìn)行計(jì)算得出方程中所涉及的相關(guān)參數(shù).

2.2.2 層次分析(AHP)法

(1)構(gòu)建判斷矩陣.確定水資源供需比、總?cè)丝?、?jīng)濟(jì)需水總量、生活需水量、污水排放量、GDP等主要指標(biāo)權(quán)重，對比每兩個(gè)主要指標(biāo)間的重要性，標(biāo)準(zhǔn)如表1所示.

表1 AHP法的標(biāo)度原則

(2)模擬結(jié)果歸一化處理.由于選取的主要指標(biāo)來自不同的系統(tǒng)，它們的量綱以及數(shù)量級存在非常大的差異，因此先要對不同量綱的主要指標(biāo)進(jìn)行歸一化處理，并計(jì)算判斷矩陣的最大特征值及特征向量.

(3)一致性檢驗(yàn).計(jì)算一致性檢驗(yàn)系數(shù)C=λmax/n-1，當(dāng)判斷矩陣階數(shù)超過2階時(shí)，引入均隨機(jī)一致性指標(biāo)R消除檢驗(yàn)系數(shù)C受到的影響，具體見表2.令CR=C/R，當(dāng)CR<0.1時(shí)，認(rèn)為判斷矩陣具有一致性.

表2 隨機(jī)一致性指數(shù)R值

3 結(jié)果與討論

3.1 系統(tǒng)動(dòng)力學(xué)模型構(gòu)建

3.1.1 模型邊界及基準(zhǔn)年的確定

安圖縣礦泉水資源承載力模型以研究區(qū)的行政界限作為模型的空間邊界，時(shí)間邊界為2010—2030年，其中2010—2014年為歷史檢驗(yàn)時(shí)段，2015—2030年為模擬時(shí)段，模擬的時(shí)間步長為1 a.

3.1.2 模型結(jié)構(gòu)

根據(jù)安圖縣水資源的供需狀況、礦泉水資源開發(fā)利用狀況、GDP、污水排放狀況以及人口數(shù)量等因素，將礦泉水資源承載力系統(tǒng)分為4個(gè)子系統(tǒng)，即水資源、經(jīng)濟(jì)、環(huán)境及人口子系統(tǒng)，依據(jù)系統(tǒng)動(dòng)力學(xué)原理，利用軟件Vensim構(gòu)建礦泉水資源承載力系統(tǒng)動(dòng)力學(xué)模型，結(jié)果見圖2.

圖2 礦泉水資源承載力系統(tǒng)動(dòng)力學(xué)模型總流

水資源子系統(tǒng).水資源供需比指供水總量與需水總量的比值，是水資源子系統(tǒng)的核心.其中，需水總量由生活需水量、生態(tài)需水量和生產(chǎn)需水量組成；供水總量由地表水供水量、地下水供水量及污水再生水量組成.

環(huán)境子系統(tǒng)主要研究定量化難度低的廢水及COD污染.廢水污染由生活、工業(yè)及第三產(chǎn)業(yè)污水組成，受到生活、工業(yè)及第三產(chǎn)業(yè)污水排放系數(shù)的影響；COD污染由生活與工業(yè)COD排放量組成，受到生活及工業(yè)COD排放系數(shù)的影響.

經(jīng)濟(jì)子系統(tǒng).根據(jù)傳統(tǒng)的產(chǎn)業(yè)分類法將經(jīng)濟(jì)子系統(tǒng)分為3個(gè)部分：第一產(chǎn)業(yè)(農(nóng)、林、牧、漁業(yè))，第二產(chǎn)業(yè)(工業(yè))和第三產(chǎn)業(yè)(服務(wù)業(yè)).當(dāng)各產(chǎn)業(yè)產(chǎn)值增加時(shí)，不但會使經(jīng)濟(jì)需水量增長，導(dǎo)致總需水量增加，還會使廢水排放量增長，影響生態(tài)環(huán)境.

人口子系統(tǒng).總?cè)丝谑侨丝谧酉到y(tǒng)的核心，根據(jù)城鎮(zhèn)化率將其劃分為城鎮(zhèn)及鄉(xiāng)村人口.城鎮(zhèn)及鄉(xiāng)村人口的數(shù)量隨著城鎮(zhèn)化率的增長而改變，通過城鎮(zhèn)及鄉(xiāng)村生活需水定額引起生活需水量的變化，進(jìn)而影響需水總量.

3.2 模型有效性檢驗(yàn)

在使用模型模擬前應(yīng)檢驗(yàn)?zāi)Ｐ偷慕Y(jié)構(gòu)和功能，以確保模型的有效性和真實(shí)性.本文將研究區(qū)2010—2014年的歷史數(shù)據(jù)與模擬結(jié)果進(jìn)行對比，如果誤差不超過±5%，那么可以保證模擬結(jié)果的精度.

由表3的檢驗(yàn)結(jié)果可見，由于統(tǒng)計(jì)數(shù)據(jù)的局限性和許多不確定因素的影響，2011和2012年需水總量的模擬值與實(shí)際值差別較大.此外，其他變量的模擬結(jié)果與歷史數(shù)據(jù)的相對誤差均小于±5%，說明模型較為真實(shí)有效，可以進(jìn)行下一步模擬預(yù)測.

表3 歷史檢驗(yàn)結(jié)果

3.3 基于層次分析法確定變量權(quán)重

運(yùn)用層次分析法建立的判斷矩陣見表4.

表4 判斷矩陣

經(jīng)計(jì)算，主要變量權(quán)重為0.398 2，0.152 0，0.074 5，0.071 8，0.253 1，0.050 3.一致性檢驗(yàn)中，I=0.090 2≤0.1，R=0.727<0.1，認(rèn)為判斷矩陣的一致性可以接受.

3.4 方案設(shè)計(jì)

根據(jù)研究區(qū)現(xiàn)狀設(shè)計(jì)不同的發(fā)展方案，選取的決策變量有人口自然增長率、城鎮(zhèn)化率，第一產(chǎn)業(yè)、工業(yè)及第三產(chǎn)業(yè)GDP增長率，引水量、提水量及地下水開采量增長率，生活、工業(yè)及第三產(chǎn)業(yè)污水排放系數(shù)，污水處理率及再生利用率等變量.通過調(diào)整不同方案下決策變量的數(shù)值，模擬礦泉水資源承載力隨時(shí)間的發(fā)展趨勢，并對其進(jìn)行分析，提出以下5種設(shè)計(jì)方案.

方案1：現(xiàn)狀型.為便于與其他方案進(jìn)行比較，此方案保持研究區(qū)發(fā)展現(xiàn)狀，各項(xiàng)決策變量數(shù)值不作大的調(diào)整.

方案2：發(fā)展型.由于研究區(qū)的發(fā)展?fàn)顩r相對落后，提出發(fā)展方案以促進(jìn)經(jīng)濟(jì)發(fā)展.此方案綜合考慮研究區(qū)歷史數(shù)據(jù)、“十三五”規(guī)劃及實(shí)際情況，放大經(jīng)濟(jì)增長速度(使第一產(chǎn)業(yè)、工業(yè)及第三產(chǎn)業(yè)GDP較以往增速約20%)，并略微調(diào)整其他參數(shù)(城鎮(zhèn)化率增長25%)反映發(fā)展型方案的反饋效果.

方案3：節(jié)約型.此方案通過適當(dāng)放緩經(jīng)濟(jì)增長速度(使工業(yè)GDP及第三產(chǎn)業(yè)較以往降速約10%)，減少城鎮(zhèn)生活需水量(使城鎮(zhèn)生活需水定額降低約10%)，降低水資源的需求量，并在此基礎(chǔ)上適當(dāng)?shù)販p少開發(fā)利用水資源量(使地下水開采量增長率、提水量增長率降低約40%)，增加與污水利用相關(guān)的決策變量(提高污水再生利用率)，同時(shí)加強(qiáng)引水工程建設(shè)(加入引水量增長率)，增加水資源的供給量.

方案4：環(huán)境型.此方案重點(diǎn)考慮保護(hù)生態(tài)環(huán)境，通過降低與廢水污染相關(guān)的決策變量(使生活、第三產(chǎn)業(yè)、工業(yè)污水排放系數(shù)降低約10%)，提高污水的處理效率(提高污水再生利用率)，減輕生態(tài)環(huán)境壓力.

方案5：協(xié)調(diào)型.綜合分析以上4種方案的長處及弊端，在考慮研究區(qū)實(shí)際情況的基礎(chǔ)上，提出了注重綜合發(fā)展，在經(jīng)濟(jì)以中高速發(fā)展(使工業(yè)GDP及第三產(chǎn)業(yè)較以往增速約10%)的同時(shí)注重節(jié)約資源(降低地下水開采量增長率、提水量增長率約10%，人口自然增長率降低約5%，并加入引水量增長率)及保護(hù)環(huán)境(使生活、第三產(chǎn)業(yè)、工業(yè)污水排放系數(shù)降低約5%，并提高污水再生利用率)的協(xié)調(diào)型方案.

3.5 模擬結(jié)果分析

對提出的5種方案使用模型進(jìn)行模擬.參考相關(guān)文獻(xiàn)，選取總?cè)丝?、水資源供需比、經(jīng)濟(jì)需水量、生活需水量、污水排放量及GDP等從不同方向反映模型結(jié)構(gòu)的變量作為主要變量.根據(jù)主要變量的模擬結(jié)果(見表5)，使用層次分析法確定各變量權(quán)重并計(jì)算綜合承載力(見表6)，得到研究區(qū)2015—2030年礦泉水資源承載力的發(fā)展趨勢(見圖3)，通過對5種方案主要變量的模擬結(jié)果進(jìn)行衡量，得出最佳方案.

圖3 承載力系統(tǒng)動(dòng)力學(xué)模型在5種方案下的模擬結(jié)果

表5 主要變量模擬結(jié)果

表6 綜合承載力計(jì)算結(jié)果

現(xiàn)狀型：至2030年，綜合承載力增長至0.068 7，研究區(qū)總?cè)丝谶_(dá)到21.79萬人，城鎮(zhèn)人口為14.71萬人，GDP為164.47億元；需水總量達(dá)到10 854.80萬m3，水資源供需比逐漸降低至0.95，污水排放量達(dá)到1 898萬m3，COD排放總量達(dá)到3 622.74 t.此方案下綜合承載力增長緩慢，水資源供需基本能夠保持平衡，可以維持社會經(jīng)濟(jì)的可持續(xù)發(fā)展.

發(fā)展型：此方案下該區(qū)經(jīng)濟(jì)快速發(fā)展，到2030年GDP總量達(dá)到181.42億元，工業(yè)GDP增長值達(dá)到46 086.20萬元.城鎮(zhèn)化率的提高、經(jīng)濟(jì)的快速發(fā)展導(dǎo)致總需水量、污水與COD排放量顯著增加，工業(yè)需水量達(dá)到2 143.60萬m3，生活需水量達(dá)到1 297萬m3，污水排放量達(dá)到2 065.08萬m3，COD總排放量達(dá)到4 011.08 t，水資源供需比快速下降至0.93.綜合承載力的上升趨勢不是來自系統(tǒng)本身的調(diào)節(jié)作用，而是在不合理的經(jīng)濟(jì)發(fā)展速度下產(chǎn)生的理論假象，經(jīng)濟(jì)的快速增長使水資源供需矛盾逐漸加大，生態(tài)環(huán)境污染較為嚴(yán)重，影響社會穩(wěn)定性.

節(jié)約型：綜合考慮研究區(qū)當(dāng)?shù)厮Y源量，減緩經(jīng)濟(jì)增長速度，降低提水、地下水開采速率的同時(shí)加強(qiáng)引水工程建設(shè).到2030年，引水量增加至8 538.16萬m3，提水量與地下水開采量分別為286.89，2 175.84萬m3，工業(yè)與第三產(chǎn)業(yè)GDP分別減緩為48.33，80.89億元.此方案通過放緩經(jīng)濟(jì)增長速度，增加引水量、污水再生利用率的方式使總需水量減少，總供水量增多，到2030年水資源供需比為1.25，獲得較大提升，綜合承載力增長至0.071 7，為前4種方案中的最大值.

環(huán)境型：此方案在現(xiàn)狀型方案的基礎(chǔ)上重點(diǎn)考慮生態(tài)環(huán)境保護(hù)，通過降低生活、工業(yè)及第三產(chǎn)業(yè)污水排放系數(shù)、生活及生產(chǎn)COD排放量，減少污染物排放量，同時(shí)增強(qiáng)對污水的處理效率，提高環(huán)境質(zhì)量，增加供水總量.到2030年，污水排放量降低至1 519.70萬m3，COD排放量降低至3 229.46 t，供水總量達(dá)到10 316.7萬m3，水資源供需比為0.963.此方案使污染物排放量明顯下降，但綜合承載力增長緩慢，為0.688，僅高于現(xiàn)狀型方案.

協(xié)調(diào)型：通過分析以上4種方案的優(yōu)缺點(diǎn)，并在考慮研究區(qū)實(shí)際情況的基礎(chǔ)上，提出了同時(shí)注重經(jīng)濟(jì)發(fā)展、節(jié)約資源及保護(hù)環(huán)境的協(xié)調(diào)型方案.到2030年，第一產(chǎn)業(yè)、工業(yè)及第三產(chǎn)業(yè)GDP增長至30.10，53.33，89.33億元，污水排放量降至1 763.42萬m3，COD排放量為3 602.35 t，供水總量達(dá)到1.45億m3，需水總量達(dá)到1.09億m3，水資源供需比提升至1.33，綜合承載力指數(shù)增長至0.073 0.協(xié)調(diào)型方案同時(shí)強(qiáng)調(diào)節(jié)約資源、保護(hù)環(huán)境及合理提升經(jīng)濟(jì)發(fā)展速度，綜合承載力為所有方案中的最大值，具有很高的可行性.

綜上所述，現(xiàn)狀型方案基本能夠維持研究區(qū)水資源的供需平衡與社會經(jīng)濟(jì)的可持續(xù)發(fā)展；發(fā)展型方案GDP增長速度顯著，但水資源供需比降低速度與污染物排放量高于其他方案，只注重經(jīng)濟(jì)發(fā)展而忽略對資源、環(huán)境的影響，為礦泉水資源承載力帶來了巨大的負(fù)擔(dān)；節(jié)約型方案節(jié)水效果明顯，水資源供需比快速提升，但經(jīng)濟(jì)增長速度最為緩慢；環(huán)境型方案顯著降低了生活及生產(chǎn)污染物排放量，但水資源供需比仍然呈現(xiàn)緩慢下降的趨勢，且綜合承載力增長較??；協(xié)調(diào)型方案綜合承載力增長最為迅速，且在保持水資源供需比增長的同時(shí)使經(jīng)濟(jì)以較快的速度增長，降低了污染物排放量，該方案是篩選出的最優(yōu)方案.

4 結(jié)論

本文利用Vensim軟件構(gòu)建了由4個(gè)子系統(tǒng)組成的礦泉水資源承載力系統(tǒng)動(dòng)力學(xué)模型，結(jié)合層次分析法確定主要變量權(quán)重，設(shè)計(jì)了現(xiàn)狀型、發(fā)展型、節(jié)約型、環(huán)境型和協(xié)調(diào)型5種方案并確定了相應(yīng)方案的決策變量，模擬了2015—2030年礦泉水資源承載力在5種方案下的動(dòng)態(tài)趨勢.通過衡量5種方案下主要變量的模擬結(jié)果，確定了最佳的方案，研究結(jié)果表明：

(1)基于系統(tǒng)動(dòng)力學(xué)理論，構(gòu)建了礦泉水資源承載力系統(tǒng)動(dòng)力學(xué)模型，采用2010—2014年的歷史數(shù)據(jù)與主要變量的模擬結(jié)果進(jìn)行檢驗(yàn)，大部分主要變量的誤差均在±5%之內(nèi)，模型真實(shí)有效.

(2)從5種方案的模擬結(jié)果來看，2030年5種方案的綜合承載力排序?yàn)榫C合型(0.073 0)>節(jié)約型(0.071 7)>發(fā)展型(0.069 7)>環(huán)境型(0.068 8)>現(xiàn)狀型(0.068 7).不同方案下綜合承載力都呈增長的趨勢，但單一的方案只能緩解研究區(qū)某一方面的壓力.現(xiàn)狀型方案基本能夠保持研究區(qū)社會經(jīng)濟(jì)的可持續(xù)發(fā)展；發(fā)展型方案研究區(qū)GDP總量能夠快速增加，但會使水資源供需比快速降低，污染物排放量顯著增加；節(jié)約型方案水資源供需比提升顯著，但放緩了經(jīng)濟(jì)增長速度；環(huán)境型方案污染物排放量顯著下降，生態(tài)環(huán)境壓力降低，但綜合承載力增長緩慢；協(xié)調(diào)型方案使經(jīng)濟(jì)以中高速發(fā)展的同時(shí)注重經(jīng)濟(jì)發(fā)展、節(jié)約資源與保護(hù)環(huán)境，對礦泉水資源承載力的提升最為顯著，為研究區(qū)未來發(fā)展的最佳方案.