數(shù)學(xué)建模論文-城市重金屬污染及地質(zhì)環(huán)境演變分析

1、<p>　　2011高教社杯全國大學(xué)生數(shù)學(xué)建模競賽</p><p><b>　　承 諾 書</b></p><p>　　我們仔細閱讀了中國大學(xué)生數(shù)學(xué)建模競賽的競賽規(guī)則.</p><p>　　我們完全明白，在競賽開始后參賽隊員不能以任何方式（包括電話、電子郵件、網(wǎng)上咨詢等）與隊外的任何人（包括指導(dǎo)教師）研究、討論與賽題有關(guān)的問題。&

2、lt;/p><p>　　我們知道，抄襲別人的成果是違反競賽規(guī)則的, 如果引用別人的成果或其他公開的資料（包括網(wǎng)上查到的資料），必須按照規(guī)定的參考文獻的表述方式在正文引用處和參考文獻中明確列出。</p><p>　　我們鄭重承諾，嚴格遵守競賽規(guī)則，以保證競賽的公正、公平性。如有違反競賽規(guī)則的行為，我們將受到嚴肅處理。</p><p>　　我們參賽選擇的題號是（從A/B/C

3、/D中選擇一項填寫）： </p><p>　　我們的參賽報名號為（如果賽區(qū)設(shè)置報名號的話）： </p><p>　　所屬學(xué)校（請?zhí)顚懲暾娜褐袊痛髮W(xué)（北京） </p><p>　　參賽隊員

4、(打印并簽名) ：1.崔天賜 </p><p>　　2. 李彥奇 </p><p>　　3. 劉逸飛 </

5、p><p>　　指導(dǎo)教師或指導(dǎo)教師組負責(zé)人 (打印并簽名)： </p><p>　　日期： 年 月 日</p><p>　　賽區(qū)評閱編號（由賽區(qū)組委會評閱前進行編號）：</p><p>　　2010高教社杯全國大學(xué)生數(shù)學(xué)建模競賽</p><p>&

6、lt;b>　　編 號 專 用 頁</b></p><p>　　賽區(qū)評閱編號（由賽區(qū)組委會評閱前進行編號）：</p><p>　　賽區(qū)評閱記錄（可供賽區(qū)評閱時使用）：</p><p>　　全國統(tǒng)一編號（由賽區(qū)組委會送交全國前編號）：</p><p>　　全國評閱編號（由全國組委會評閱前進行編號）：</p><

7、;p>　　城市重金屬污染及地質(zhì)環(huán)境演變分析</p><p><b>　　摘要</b></p><p>　　在問題1中，我們首先根據(jù)離散的300余個樣本的坐標及其相應(yīng)的重金屬濃度，用matlab分別擬合繪制8種重金屬元素濃度隨坐標位置變化的濃度分布漸變圖，用不同顏色代表不同的污染程度，定性直觀地看出8種重金屬元素在該城區(qū)的空間分布。同時，我們引入污染指數(shù)，利用地累

8、積指數(shù)法對8種元素在5塊區(qū)域的污染程度進行定量計算，并且對污染程度進行分級，不僅定量分析了不同區(qū)域的重金屬污染情況，而且對第一部分的擬合圖做了檢驗，兩部分吻合很好。</p><p>　　在問題2中，我們建立主成分分析模型，運用主成分分析法，用較少的新變量（主成分）代替原來較多的舊變量（8種重金屬元素），既減少了變量個數(shù)，又保證了新變量能夠覆蓋大部分信息。這樣便可以比較容易地通過主成分挖掘出各污染元素與區(qū)域功能之間

9、的關(guān)系，進而確定重金屬污染物的來源。 在問題3中，我們分析了重金屬污染物在地層中的擴散屬于“分子擴散”而無“對流擴散”，符合菲克擴散定律的傳播特征。我們根據(jù)空間中任一體積元內(nèi)金屬污染物的物料守恒，結(jié)合反映金屬元素傳播特征的菲克擴散定律，建立金屬元素在地層中的擴散模型，得到濃度關(guān)于坐標的拉普拉斯方程。我們采用有限差分法進行數(shù)值求解，在問題一濃度分布圖的基礎(chǔ)上，利用直觀濃度最大點附近的若干已知樣本數(shù)據(jù)作為邊界條件，借助Excel進行

10、迭代運算，使所得到的數(shù)值解與真實情況最大誤差為0.001。再用matlab篩選出濃度最大點的坐標即為要求的污染源位置。由此分別求得8種重金屬的污染源。</p><p>　　在問題4中，我們首先分析了研究地質(zhì)環(huán)境轉(zhuǎn)換的具體內(nèi)容，也就是既要研究不同區(qū)域地質(zhì)環(huán)境狀況的“橫向?qū)Ρ取保忠芯坎煌瑫r間地質(zhì)環(huán)境狀況的“縱向?qū)Ρ取?，對前三問建立的模型進行了優(yōu)缺點分析。并且提出要從物質(zhì)組成、地質(zhì)結(jié)構(gòu)、動力作用三個方面補充收集14

11、類信息，建立層次分析模型，得到了一套地質(zhì)環(huán)境狀況評價指標，也就是各因素所占權(quán)重（權(quán)重的變化可以看出演化的過程）及綜合評價指標。再將不同時期的評價體系進行對比，即可得出城市地質(zhì)環(huán)境在的演變模式。</p><p>　　關(guān)鍵詞：地污染指數(shù)；主成分分析；菲克定律；層次分析</p><p><b>　　問題重述</b></p><p>　　為了對某城市

12、城區(qū)土壤地質(zhì)環(huán)境進行調(diào)查，研究人類活動影響下城市地質(zhì)環(huán)境的演變模式，現(xiàn)將城區(qū)按功能劃分為生活區(qū)、工業(yè)區(qū)、山區(qū)、主干道路區(qū)及公園綠地區(qū)等，分別記為1類區(qū)、2類區(qū)、……、5類區(qū)，進行研究，不同的區(qū)域環(huán)境受人類活動影響的程度不同。</p><p>　　將所考察的城區(qū)劃分為間距1公里左右的網(wǎng)格子區(qū)域，按照每平方公里1個采樣點對表層土進行取樣、編號，并用GPS記錄采樣點的位置。應(yīng)用專門儀器測試分析，獲得了每個樣本所含的多種

13、化學(xué)元素的濃度數(shù)據(jù)。另一方面，按照2公里的間距在那些遠離人群及工業(yè)活動的自然區(qū)取樣，將其作為該城區(qū)表層土壤中元素的背景值。</p><p>　　這樣便得到了采樣點的位置、海拔高度及其所屬功能區(qū)；8種主要重金屬元素在采樣點處的濃度；8種主要重金屬元素的背景值。</p><p>　　依據(jù)以上信息，通過數(shù)學(xué)建模解決以下問題：</p><p>　　（1）給出8種主要重金屬元

14、素在該城區(qū)的空間分布，并分析該城區(qū)內(nèi)不同區(qū)域重金屬的污染程度。</p><p>　　（2）通過數(shù)據(jù)分析，說明重金屬污染的主要原因。</p><p>　?。?）分析重金屬污染物的傳播特征，由此建立模型，確定污染源的位置。</p><p>　?。?）分析所建立模型的優(yōu)缺點，為更好地研究城市地質(zhì)環(huán)境的演變模式，還應(yīng)收集什么信息？有了這些信息，如何建立模型解決問題？<

15、/p><p><b>　　問題分析</b></p><p>　　對于問題1，分析題目要求可以看出它有兩部分構(gòu)成：重金屬元素在各區(qū)域的空間分布分析，以及不同區(qū)域重金屬污染程度評價。要分析重金屬元素在各區(qū)域的空間分布，結(jié)合題目所給出的離散的樣品點，我們想到用matlab擬合出每種重金屬元素的分布圖，可以直觀地看出重金屬的分布特征。對于第二部分，要想對不同區(qū)域重金屬的污染程度進

16、行評價，只要引入一個污染評價指標便可以定量評價每個區(qū)域的污染程度。</p><p>　　對于問題2，要分析重金屬污染的原因，也就是各種金屬的來源，需要找出區(qū)域的功能與不同污染元素之間的關(guān)系。但由于污染元素一共有8種，變量因子過多，且原始數(shù)據(jù)紛雜，不利于找出內(nèi)在聯(lián)系，于是想到用較少的新變量代替原來較多的舊變量，而且保證新變量盡可能多的保留原來變量所反映的信息。這樣就可以通過分析少量新變量之間的關(guān)系，挖掘出內(nèi)在聯(lián)系。

17、主成分分析法便可以實現(xiàn)這一目的。</p><p>　　對于問題3，要得到污染源，必須要得到地層中重金屬元素的濃度分布，以求濃度最大點。首先我們要明確重金屬元素在地層中的擴散規(guī)律，然后根據(jù)這一規(guī)律建立擴散模型。通過查閱大量文獻，我們總結(jié)出重金屬元素在地層中的擴散只有“分子擴散”，而無“對流擴散”，滿足菲克定律。這樣我們根據(jù)特定地層區(qū)域內(nèi)重金屬“物料守恒”建立了擴散模型，求得了地層中的濃度分布，進而得到濃度最大點也就

18、是污染源。</p><p>　　對于問題4，要研究地質(zhì)環(huán)境演變時，既要研究不同區(qū)域地質(zhì)環(huán)境狀況的“橫向?qū)Ρ取?，又要研究不同時間地質(zhì)環(huán)境狀況的“縱向?qū)Ρ取?，并且要明確地質(zhì)環(huán)境這一概念所包含的因素（重金屬污染只是其中一種）。明確了以上兩點也就明確了我們研究地質(zhì)環(huán)境演變規(guī)律所需要的信息。而在“橫向?qū)Ρ取边^程中，我們必須要對地質(zhì)環(huán)境狀況有一套評價指標，也就是各因素所占權(quán)重（權(quán)重的變化可以看出演化的過程）及綜合評價指標?！?/p>

19、縱向?qū)Ρ取敝恍枰獙Ρ炔煌瑫r間同一區(qū)域的地質(zhì)環(huán)境狀況即可。</p><p><b>　　模型假設(shè)</b></p><p>　?、僦亟饘傥廴疚餄舛仍诳臻g的分布連續(xù)；</p><p>　?、谥亟饘傥廴疚镌诘貙又械臄U散系數(shù)為常數(shù)；</p><p>　?、墼摮菂^(qū)的重金屬污染物分布在短時間內(nèi)不隨時間變化，是穩(wěn)定狀態(tài)。</p&g

20、t;<p><b>　　模型建立與求解</b></p><p>　　4.1重金屬元素分布與不同區(qū)域污染程度（問題1）</p><p>　　4.1.1各重金屬元素在不同區(qū)域的分布圖</p><p>　　要分析8種重金屬元素在該城區(qū)的空間分布，我們想到了根據(jù)給出的樣本測量數(shù)據(jù)，把8種元素分開討論，剔除樣本中的特異點后，對于每一種元素都

21、用matlab軟件的“griddata”語句繪制出該種元素在5個區(qū)域內(nèi)的連續(xù)分布（在不同區(qū)域內(nèi)的濃度大?。┑木W(wǎng)格圖，進而直觀地看出不同區(qū)域每一種重金屬污染的程度。如下圖，是通過matlab繪制的各種重金屬元素的污染分布圖（繪圖源代碼見 附錄—Matlab程序—第1問）。</p><p>　　由于圖片所占頁面空間太大，這里只貼出As、Cd元素的分布圖作為代表，剩下6種元素的分布圖見附件。</p>&l

22、t;p>　　對于8種重金屬元素的分布圖來說，紅色表示濃度最大，紫色表示濃度最小。顏色由暖色至冷色表示濃度遞減。</p><p>　　圖1：As元素在各區(qū)域的分布</p><p>　　從圖1中可以直觀地看出As元素主要分布在工業(yè)區(qū)、生活區(qū)與綠地周圍，而在交通區(qū)分布較少。在山區(qū)大部分地區(qū)分布較少，少部分地方有分布。</p><p>　　圖2：Cd元素在各區(qū)域的分布

23、</p><p>　　對于圖2，可以直觀地看出Cd元素主要分布在工業(yè)區(qū)與交通區(qū)周圍，次之分布在生活區(qū)與綠地，在山區(qū)分布較少。</p><p>　　另外，由其它6種重金屬元素的分布圖（見附件—圖像—第1問）可以看出Cr元素主要分布在生活區(qū)與交通區(qū)，在工業(yè)區(qū)、綠地與山區(qū)分布較少。</p><p>　　對于Cu元素，發(fā)現(xiàn)工業(yè)區(qū)與交通區(qū)分布較多，其它地方分布較少。</

24、p><p>　　對于Hg元素，發(fā)現(xiàn)只有少部分生活區(qū)與交通區(qū)分布較多，其它部分基本上沒有受到Hg元素的污染。</p><p>　　對于Ni元素，發(fā)現(xiàn)在山區(qū)、生活區(qū)與工業(yè)區(qū)分布較多，在交通區(qū)與綠地分布較少。</p><p>　　對于Pb元素，發(fā)現(xiàn)在交通區(qū)分布最多，在部分生活區(qū)與部分工業(yè)區(qū)也有所分布，而在其它區(qū)域都分布較少。</p><p>　　對于Z

25、n元素，發(fā)現(xiàn)在生活區(qū)與工業(yè)區(qū)分布較多，其它區(qū)域分布較少。</p><p>　　4.1.2 地累積指數(shù)法評價不同區(qū)域重金屬污染程度</p><p>　　地累積指數(shù)法是一種常用的研究沉積物中重金屬污染的定量指標，其計算公式為： （1）</p><p>　　式

26、中：是實測元素在地層中的含量；為元素在地層中的背景值；是考慮了各地巖石差異可能會引起的變動而取得變動系數(shù)（一般情況下取值為1.5）?！?】</p><p>　　得到的地積累指數(shù)可分為不同的等級，分級情況與污染程度的對應(yīng)關(guān)系見下表：</p><p>　　表1：重金屬污染級別、地積累指數(shù)（）和分級比較</p><p>　　對于本題，是題目中所給的各種元素的背景值。是通過

27、題目所給的數(shù)據(jù)統(tǒng)計得到的各區(qū)域中各種重金屬元素的平均含量，如下表：</p><p>　　表2：各區(qū)域中各種重金屬元素的平均含量</p><p>　　把代入到式（1）后，得到5塊區(qū)域內(nèi)各種重金屬元素的地積累指數(shù)值（），如下表：</p><p>　　表3：各區(qū)域內(nèi)各種重金屬元素的地積累指數(shù)值</p><p>　　通過上表各區(qū)域的污染級別可以看出：

28、</p><p>　　對于生活區(qū)，Cd、Cr、Cu、Hg、Pb、Zn這幾種重金屬元素屬于偏中度污染；As、Ni兩種重金屬元素屬于輕度污染。</p><p>　　對于工業(yè)區(qū)，Hg是重污染；Cu是偏重污染；Zn是中度污染；As、Cd、Pb屬于偏中度污染；Cr、Ni屬于輕度污染。</p><p>　　對于山區(qū)，所有的重金屬元素都是輕度污染。</p><

29、p>　　對于交通區(qū)，Hg是偏重污染；Cu是中度污染；Cd、Pb、Zn均屬于偏中度污染；As、Ni屬于輕度污染。</p><p>　　對于綠地區(qū)，Cd、Cu、Hg、Zn屬于偏中度污染；As、Cr、Ni、Pb均屬于輕度污染。</p><p>　　另外，為了對比各區(qū)域的總污染程度，我們對各個區(qū)域內(nèi)的所有重金屬元素的地積累指數(shù)求和，得到各區(qū)域的總地積累指數(shù)（），結(jié)果如下表：</p>

30、;<p>　　表4：各區(qū)域的總地積累指數(shù)（）</p><p>　　通過上表比較可以發(fā)現(xiàn)不同區(qū)域的污染程度為：工業(yè)區(qū)>交通區(qū)>生活區(qū)>綠地區(qū)>山區(qū)。</p><p>　　通過對比4.1.1與4.1.2的結(jié)論，發(fā)現(xiàn)重金屬元素的分布情況與各區(qū)域內(nèi)的重金屬污染程度相互呼應(yīng)，吻合的很好。該城區(qū)的重金屬污染特點，我們可以總結(jié)為：該城區(qū)的主要污染元素為Hg、Cu、Z

31、n；且污染主要集中在工業(yè)區(qū)與交通區(qū)。</p><p>　　4.2主成分分析模型（問題2）</p><p>　　要研究重金屬污染的主要原因，就要發(fā)現(xiàn)不同區(qū)域內(nèi)主要污染物與該區(qū)域功能的關(guān)聯(lián)性，進而確定各區(qū)域污染原因之間的關(guān)聯(lián)性，以確定整個城區(qū)污染的原因。</p><p>　　但由于導(dǎo)致每個區(qū)域污染的重金屬元素眾多，沒法直接通過重金屬元素之間的關(guān)系確定各區(qū)域污染原因之間的

32、關(guān)聯(lián)性。于是想到用較少的新變量代替原來較多的舊變量，而且保證新變量盡可能多的保留原來變量所反映的信息。主成分分析法便可以實現(xiàn)，主成分分析法是把原來多個變量劃為少數(shù)幾個綜合指標的一種統(tǒng)計分析方法。 </p><p>　　4.2.1主成分分析模型的建立</p><p>　　假定有個地理樣本，每個樣本共有個變量描述，這樣就構(gòu)成了一個階的地理數(shù)據(jù)矩陣：</p><p>&l

33、t;b>　?。?）</b></p><p>　　要想從這么多變量的數(shù)據(jù)中抓住地理事物的內(nèi)在規(guī)律，自然要在維空間中加以考察，這是比較麻煩的。為了克服這一困難，就需要用較少的幾個綜合指標來代替原來較多的變量指標，這些綜合指標（即新變量）選取最簡單的形式就是取原來變量指標的線性組合，適當(dāng)調(diào)整組合系數(shù)，使新的變量指標之間相互獨立且代表性最好。</p><p>　　若原來的變量指標

34、為，它們的綜合指標——新變量指標為，則</p><p><b>　?。?）</b></p><p>　　在（3）式中，系數(shù)由下列原則來決定：</p><p><b>　　①（）相互無關(guān)；</b></p><p>　?、谑堑囊磺芯€性組合中方差最大者；是與不相關(guān)的的所有線性組合中方差最大者；……；是與都

35、不相關(guān)的的所有線性組合中方差最大者。</p><p>　　這樣決定的新變量指標分別稱為原變量指標的第一，第二，…，第主成分。其中，在總方差中占的比例最大，的方差依次遞減。在實際問題的分析中，常挑選前幾個最大的主成分，這樣既減少了變量的數(shù)目，又抓住了主要矛盾，簡化了變量之間的關(guān)系。</p><p>　　從以上分析可以看出，找主成分就是確定原來變量()在諸主成分（）上的載荷（），從數(shù)學(xué)上容易知

36、道，它們分別是的相關(guān)矩陣個較大的特征值所對應(yīng)的特征向量。</p><p><b>　　模型計算步驟：</b></p><p>　?、儆嬎阆嚓P(guān)系數(shù)矩陣：</p><p><b>　　（4）</b></p><p>　　為原變量與的相關(guān)系數(shù)，，計算公式為：</p><p><

37、;b>　?。?）</b></p><p>　?、谟嬎闾卣髦蹬c特征向量。</p><p>　　先解特征方程，常用雅克比法求出特征值，并使其按大小順序排列：。</p><p>　　然后，分別求出對于特征值的特征向量，要求，即，其中表示向量的第個分量。</p><p>　?、塾嬎阒鞒煞重暙I率及累計貢獻率。</p>&

38、lt;p>　　貢獻率： （6）</p><p>　　累計貢獻率： （7）</p><p>　　一般累計貢獻率達到80%以上的特征值所對應(yīng)的為第1、第2、…、第個主成分。</p&

39、gt;<p><b>　　④計算主成分載荷</b></p><p><b>　?。?）</b></p><p><b>　　4.2.2模型求解</b></p><p><b>　?、儆嬎阆嚓P(guān)系數(shù)矩陣</b></p><p>　　對所給的數(shù)據(jù)進

40、行分區(qū)域標準化處理，利用式（5）得到5個區(qū)域的相關(guān)系數(shù)矩陣如下：（具體計算源代碼見 附錄—matlab—相關(guān)系數(shù)）</p><p><b>　　生活區(qū)：</b></p><p><b>　　工業(yè)區(qū)：</b></p><p><b>　　山區(qū)：</b></p><p><b

41、>　　交通區(qū)：</b></p><p><b>　　綠地區(qū)：</b></p><p>　?、谟嬎闾卣髦?、特征向量，求貢獻率與累計貢獻率</p><p>　　根據(jù)式（6）（7）可以求得各個區(qū)域的主成分的貢獻率與累計貢獻率，見下表：（具體計算源代碼見 附錄—matlab程序—第2問）</p><p>　　表

42、5：生活區(qū)特征值及貢獻率 表6：工業(yè)區(qū)特征值及貢獻率</p><p>　　表7：山區(qū)特征值及貢獻率 表8：交通區(qū)特征值及貢獻率</p><p>　　表9：綠地區(qū)特征值及貢獻率</p><p>　　一般情況下認為前幾個較大的主成分累計貢獻率達到65%便可以代替原來的多個指

43、標，能基本覆蓋其中的信息。</p><p>　　所以，根據(jù)以上5個區(qū)域的貢獻率表格可以看出，對于生活區(qū)只需求出前3個主成分即可；對于工業(yè)區(qū)只需求出前2個主成分即可；對于山區(qū)，只需求出前3個主成分即可；對于交通區(qū)，只需求出前3個主成分即可；對于綠地區(qū)，只需求出前3個主成分即可。</p><p><b>　?、矍笾鞒煞州d荷</b></p><p>

44、　　對于5個區(qū)域，利用式（8），分別求出原變量在各個選定的較大主成分上的載荷。如下表：（具體算法見 附錄—matlab—第2問；詳細數(shù)據(jù)見 附錄—附表—第2問）</p><p>　　表10：各區(qū)域主成分載荷</p><p>　　（說明：由于計算載荷時我們只精確到小數(shù)點后一位，所以，有的載荷值經(jīng)過“四舍五入”，為0.0）</p><p><b>　?、芙Y(jié)果分

45、析：</b></p><p>　　對于生活區(qū)，第一主成分的特點表現(xiàn)在因子變量在Cr、Cu、Hg的濃度上有較高的正載荷，反映了Cr、Cu、Hg的富集程度；第二主成分的特點表現(xiàn)在因子變量As的濃度上有很大的正載荷，反映了As的富集程度；第三主成分的特點表現(xiàn)在因子變量Hg上有很大的正載荷，反映了Hg的富集程度。</p><p>　　對于工業(yè)區(qū)，第一主成分的特點表現(xiàn)在因子變量Cu、Pb

46、上有較高的正載荷，反映了Cu、Pb的富集程度；第二主成分的特點表現(xiàn)在因子變量As、Ni上有較高的正載荷，反應(yīng)了As、Ni的富集程度。</p><p>　　對于山區(qū)，第一主成分的特點表現(xiàn)在因子變量Cu、Ni、Zn上有較高的正載荷，反應(yīng)了Cu、Ni、Zn的富集程度；第二主成分的特點表現(xiàn)在因子變量Cd上具有很高的正載荷，反應(yīng)了Cd的富集程度；第三主成分的特點表現(xiàn)在因子變量As、Cd、Cu、Hg上具有很高的正載荷，反映了

47、As、Cd、Cu、Hg。</p><p>　　對于交通區(qū)，第一主成分的特點表現(xiàn)在因子變量Pb上具有較高的正載荷，反映了Pb的富集程度；第二主成分的特點表現(xiàn)在因子變量As上具有很高的正載荷，反映了As的富集程度；第三主成分的特點表現(xiàn)在因子變量Cd、Hg上具有較高的正載荷，反映了Cd、Hg的富集程度。</p><p>　　對于綠地區(qū)，第一主成分的特點表現(xiàn)在因子變量Cd、Hg上具有較高的正載荷，

48、反映了Cd、Hg的富集程度；第二主成分的特點表現(xiàn)在因子變量As、Cu上具有較高的正載荷，反映了As、Cd的富集程度；第三主成分的特點表現(xiàn)在因子變量As、Hg上具有較高的正載荷，反映了As、Hg的富集程度。</p><p>　　總結(jié)上面的分析可以發(fā)現(xiàn)：對于As，在5個區(qū)域的第一主因素上都沒有反應(yīng)，在生活區(qū)、工業(yè)區(qū)、交通區(qū)以及綠地區(qū)的第二主因素上有所反應(yīng)，說明As元素的污染主要來源于人類生活、工業(yè)排放、汽車尾氣排放以

49、及大氣沉降等；對于Cd，只在綠地區(qū)的第一主因素上，可能來源與汽車尾氣或大氣沉降；對于Cr，只在生活區(qū)的第一主元素上有所反應(yīng)，表明起污染主要來源于人類生活；對于Cu，在生活區(qū)、工業(yè)區(qū)以及山區(qū)的第一主元素上均有所反應(yīng)，說明其污染來源于人類生活、工業(yè)排放和土壤本身；對于Hg，在生活區(qū)與綠地區(qū)的第一主因素上均有所反應(yīng)，說明其主要來源于人類生活和大氣沉降；對于Ni，只在山區(qū)的第一主因素上有所反應(yīng)，說明其來源于土壤本身；對于Pb，發(fā)現(xiàn)其在交通區(qū)與工

50、業(yè)區(qū)的第一主因素上有所反應(yīng)，說明其來源于汽車尾氣排放與工業(yè)排放；對于Zn，發(fā)現(xiàn)其只在山區(qū)的第一主因素上有所反應(yīng)，說明其污染主要來源與土壤本身。</p><p>　　4.3重金屬在地層中的擴散模型（問題3）</p><p>　　4.3.1重金屬元素在地層中的擴散特點</p><p>　　由于我們研究表層土壤中重金屬元素的分布，所以其擴散規(guī)律不涉及到隨流體流動而產(chǎn)生的“

51、對流擴散”，也就是說其擴散情況只有“分子擴散”。分子擴散完全是因為組分分布不均勻，即在空間中存在濃度梯度，導(dǎo)致這些組分依靠分子熱運動從高濃度帶擴散到低濃度帶，最后趨于平衡。根據(jù)描述物質(zhì)“分子擴散”現(xiàn)象的規(guī)律——菲克擴散定律，在單位時間內(nèi)通過垂直于擴散方向的單位截面積的擴散物質(zhì)流量（稱為擴散通量）與濃度沿該截面外法向的方向?qū)?shù)成正比，也就是說，濃度沿截面法向單位距離的變化率越大，擴散通量越大。</p><p>　　

52、4.3.2擴散模型的建立</p><p>　　在地層中，若各點金屬元素濃度不同，就會有重金屬元素由高濃度帶向低濃度帶擴散。確定污染源的位置，就是要找到重金屬元素的濃度最高點，故問題歸結(jié)為研究地層中重金屬元素的濃度分布，即根據(jù)重金屬元素在地層中的擴散特點建立濃度關(guān)于位置的函數(shù)關(guān)系。

53、

54、 </p><p>　　假設(shè)地層中任意一點在時刻t的重金屬污染物濃度為。假設(shè)在時間區(qū)間內(nèi)，流過點處面積為的曲面重金屬流量為，則由菲克定律的得：</p><p><b>　?。?）</b&g

55、t;</p><p>　　其中比例系數(shù)稱為重金屬的擴散系數(shù)（大于0），在一般情況下都把它作為常數(shù)處理。是曲面的外法線單位向量。負號表示流量方向與濃度梯度方向相反。</p><p>　　在地層中任取一封閉曲面，記它所包圍的區(qū)域為。于是，從時刻t1到時刻t2，通過曲面流出該區(qū)域的重金屬污染物質(zhì)量為：</p><p><b>　?。?0）</b>&l

56、t;/p><p>　　若內(nèi)有污染源不斷產(chǎn)生污染物，設(shè)污染源密度為（在時刻t于任一點處單位時間、單位體積所產(chǎn)生的污染物質(zhì)量），則在時刻t1到t2內(nèi)產(chǎn)生的污染物質(zhì)量為：</p><p><b>　?。?1）</b></p><p>　　由于物質(zhì)擴散，在時刻t1到t2內(nèi)，中的污染物質(zhì)量的改變量為：</p><p><b>

57、;　?。?2）</b></p><p><b>　　由物料守恒定律有：</b></p><p><b>　　，即</b></p><p><b>　?。?3）</b></p><p>　　式（13）就是我們建立的重金屬污染物在地層中的擴散模型。根據(jù)這個方程便可以求得

58、。</p><p>　　由于我們研究的是該城區(qū)地層在這一確定時刻的污染情況，所以與時間無關(guān)，即的值為0。另外，在研究該問題的時候，我們在短時間內(nèi)不考慮新的污染物的注入，所以。而又為常數(shù)，這樣我們的模型可以簡化為即：</p><p><b>　?。?4）</b></p><p>　　由于方程（14）是偏微分方程,只有在問題具有高度對稱的情況下,才

59、能求出解析解,但我們研究的重金屬元素在地層中的擴散受到地貌、濃度分布、以及其他復(fù)雜因素的影響，顯然不能符合高度對稱這一條件。對于這些情況,只能尋求它的數(shù)值解?！?】</p><p>　　4.3.3模型求解（拉普拉斯方程的數(shù)值求解）</p><p>　　有限差分法是求解拉普拉斯數(shù)值解的一種常用算法【3】，其思想是用差分代替導(dǎo)數(shù)，用網(wǎng)格將求解區(qū)域覆蓋，對于平面拉普拉斯方程，第i行，第j列小格的

60、值由拉普拉斯五點差分格式給出：</p><p><b>　　（15）</b></p><p>　　我們根據(jù)這個算法來求上述擴散模型的拉普拉斯方程的數(shù)值解，以確定As元素的一個污染源為例，具體操作步驟如下：</p><p>　?、僭谟蓡栴}1得到的As元素濃度分布圖上，大致確定濃度最高點的位置，取以（12380，10043），（12380,3532

61、），（15107，10043），（15107，3532）為頂點的矩形框?qū)⑵浼案浇?個已知濃度點圈在內(nèi)；</p><p>　?、趯⒕匦慰騽澐譃?25×136個邊長為20m的小正方形，得到126×137個交點，將8個已知濃度點的位置分別取為離其最近的交點處；</p><p>　?、塾肊xcel中的單元格代表交點，單元格中的數(shù)值代表As元素濃度，建立126×137

62、的矩陣，在相應(yīng)單元格中輸入8個已知的濃度值，在其他未知濃度值的單元格中輸入0，將Excel表初始化；</p><p>　?、苋鐖D3所示，利用Excel作為計算工具利用公式（15）對各濃度值進行32767次迭代運算，使其最大誤差為0.001。</p><p>　　根據(jù)以上步驟，我們求得了矩形區(qū)域內(nèi)126×137個點處的As元素濃度值，即求得了該范圍內(nèi)滿足8個已知點濃度條件的拉普拉斯

63、方程的數(shù)值解。利用matlab編程找到其中的濃度最大點，即為該區(qū)域As元素污染源。根據(jù)以上步驟，我們在Excel中輸入相關(guān)數(shù)據(jù)以及公式，進行迭代計算，得到了8種元素的濃度最高的坐標，也就是污染源的位置。下圖是As元素一個污染源的計算截圖。</p><p>　　圖3：Excel執(zhí)行有限差分法截圖</p><p>　　這樣，我們便得到了8種重金屬元素的污染源位置坐標，如下表：</p>

64、;<p>　　表11：8種重元素的污染源位置坐標</p><p>　?。ㄕf明：每種元素的含量均有多個峰值，也就是有多個污染源。序號1,2…表示污染源的編號）</p><p>　　4.4城市地質(zhì)環(huán)境演變模式分析（問題4）</p><p>　　4.4.1所建模型局限性</p><p>　?、僭谇?問建立的所有模型中，城區(qū)地層評價指標

65、都針對于城區(qū)以及不同區(qū)域的重金屬污染物而建立的，這些評價指標都只能評價重金屬污染狀況，而重金屬污染狀況只是地質(zhì)環(huán)境的很小的一個方面，所以我們前面所建立的模型雖然能比較準確的評價地層重金屬污染的情況，但卻不足以全面地反應(yīng)地層的地質(zhì)環(huán)境狀況，模型反應(yīng)的信息太過于局限；</p><p>　?、谝芯吭摮菂^(qū)地質(zhì)環(huán)境的演化模式，必須要對該城區(qū)的地質(zhì)狀況進行縱向?qū)Ρ龋簿褪遣煌瑫r間的地質(zhì)環(huán)境狀況對比。而我們在前3問建立的模型

66、，都是反映該城區(qū)地層在這一特定時間的地質(zhì)環(huán)境狀況，不包含時間這一參數(shù)，也就不能對地質(zhì)環(huán)境狀況進行縱向?qū)Ρ取７从巢怀龅刭|(zhì)環(huán)境狀況隨時間的變化情況。</p><p>　　4.4.2研究地質(zhì)環(huán)境狀況演變所需要的信息</p><p>　　地質(zhì)環(huán)境指由巖石圈、水圈和大氣圈組成的環(huán)境系統(tǒng)。城市發(fā)展對環(huán)境的改造與擾動主要通過三個方面, 即: 大規(guī)模的地面和地下工程建設(shè)；城市及其鄰區(qū)地下資源的開發(fā)；城市垃

67、圾污染?；?.4.1的局限性分析，我們可以總結(jié)出研究地質(zhì)狀況演變模式所需要的信息：</p><p>　　地基巖土的類別與結(jié)構(gòu)；巖土體承載力；特殊類巖土體引起的地面變形；地表水位置；地下水水位及侵蝕性；土壤厚度；水土流失程度；平均侵蝕模數(shù)；其他地質(zhì)資源；地質(zhì)構(gòu)造；地面坡度；地震及強度；洪水和河流侵蝕；新構(gòu)造活動和風(fēng)化作用【4】。</p><p>　　以上信息均需要有不同年份的對比。<

68、/p><p>　　4.4.3基于“層次分析法”模型的城市地質(zhì)環(huán)境演變模式分析</p><p>　　我們在研究城市地質(zhì)環(huán)境演變模式時，主要有兩部分工作：①根據(jù)“層次分析法”，結(jié)合反映地質(zhì)環(huán)境的主要信息確定一個評價某一區(qū)域地質(zhì)環(huán)境狀況的評價指標體系，也就是地質(zhì)環(huán)境狀況各因素所占權(quán)重（權(quán)重的變化可以看出演化的過程）及綜合評價指標，對該城區(qū)的不同區(qū)域進行橫向?qū)Ρ龋虎诟鶕?jù)不同年份得到的該城區(qū)的評價指標體

69、系來進行縱向?qū)Ρ龋治銎涞刭|(zhì)演化過程。</p><p>　　①層次分析法建立評價體系</p><p>　　層次分析法是一種系統(tǒng)分析綜合方法，用于求解層次結(jié)構(gòu)或網(wǎng)絡(luò)結(jié)構(gòu)等復(fù)雜系統(tǒng)的評估問題。它的基本原理是對評價系統(tǒng)的有關(guān)方案的各種要素分解成若干層次，并以同一層次的各種要求按照上一層要求為準則，進行兩兩判斷比較和計算，求出各要素的權(quán)重，進而確定綜合評價指標。</p><p&

70、gt;　　具體建模過程如下：【5】</p><p>　　1、假設(shè)要比較某一層n個元素對上一層一個因素O的影響。取兩個元素，用表示對O影響之比，全部比較結(jié)果組成成對比較矩陣：</p><p><b>　?。?6）</b></p><p>　　2、采用“方根法”解特征方程，求得特征向量。其中即為這一層n個評價元素的權(quán)重。</p>&l

71、t;p>　　3、按照上述評價因子及其權(quán)值，確定綜合評價指標：</p><p>　　②根據(jù)不同年份得到的該城區(qū)的評價指標體系來進行縱向?qū)Ρ?，分析其地質(zhì)演化過程。</p><p><b>　　模型評價</b></p><p><b>　　5.1模型優(yōu)點</b></p><p>　　問題1，我們分兩

72、部分解決，第一部分是畫出各元素的分布圖，從感性上認識其分布特征，直觀；第二部分運用地積累指數(shù)法定量計算各區(qū)域的污染情況，既把各區(qū)域的污染狀況進行了量化又對第一部分進行了驗證。</p><p>　　問題2，主層次分析模型中運用主要因素對復(fù)雜的變量因子進行歸類，挖掘出了各變量之間并不明顯的內(nèi)在聯(lián)系，進而確定出了各重金屬污染的可能來源。</p><p>　　問題3，根據(jù)重金屬元素擴散規(guī)律建立的擴

73、散模型，能真實地反映地層中重金屬元素的擴散，并且應(yīng)用有限差分法巧妙地對沒有解析解的模型進行了數(shù)值求解。</p><p>　　問題4，在研究地質(zhì)環(huán)境演變時，既研究了不同區(qū)域地質(zhì)環(huán)境狀況的“橫向?qū)Ρ取保盅芯苛瞬煌瑫r間地質(zhì)環(huán)境狀況的“縱向?qū)Ρ取?，并且對地質(zhì)環(huán)境狀況的各要素做了層析分析，確定了其權(quán)重以及綜合評價指標。模型考慮全面，且邏輯清晰。</p><p><b>　　5.2模型缺點

74、</b></p><p>　　問題2模型雖然能夠清晰地反映變量之間的內(nèi)在聯(lián)系，很好的解釋污染物的來源，但是計算量較大，數(shù)據(jù)處理繁瑣。</p><p>　　問題3中的擴散模型，為了計算簡便，把擴散系數(shù)規(guī)定為常數(shù)，而實際它是隨不同位置地層微觀上的變化而變化的，會對結(jié)果造成一定誤差。</p><p><b>　　參考文獻</b></

75、p><p>　　祝云龍 姜加虎 孫占東 黃群 王紅娟 周云凱，洞庭湖沉積物中重金屬污染特征與評價，湖泊科學(xué)，2008,20(4): 477-485</p><p>　　曾紹標 熊洪允 毛云英，應(yīng)用數(shù)學(xué)基礎(chǔ)，天津：天津大學(xué)出版社，1993.3</p><p>　　賈新民 嚴文，有限差分法求解拉普拉斯方程，昌吉學(xué)院學(xué)報，2009年第5期：105-109</p>

76、<p>　　倪師軍 魏倫武 張成江等，城市地質(zhì)環(huán)境風(fēng)險性分區(qū)評價體系，地質(zhì)通報，2006.11：1279-1286</p><p>　　蔣啟源 謝金星 葉俊，數(shù)學(xué)模型，北京：高等教育出版社，2003.8</p><p><b>　　附錄</b></p><p>　　下面六幅圖從左至右分別為Cr，Cu，Hg，Ni，Pb，Zn六種元素的