城市表層土壤重金屬污染分析數(shù)學(xué)建模國家二等獎

2011高教社杯全國大學(xué)生數(shù)學(xué)建模競賽承諾書我們仔細閱讀了中國大學(xué)生數(shù)學(xué)建模競賽的競賽規(guī)則.我們完全明白，在競賽開始后參賽隊員不能以任何方式（包括電話、電子郵件、網(wǎng)上咨詢等）與隊外的任何人（包括指導(dǎo)教師）研究、討論與賽題有關(guān)的問題。我們知道，抄襲別人的成果是違反競賽規(guī)則的,如果引用別人的成果或其他公開的資料（包括網(wǎng)上查到的資料），必須按照規(guī)定的參考文獻的表述方式在正文引用處和參考文獻中明確列出。我們鄭重承諾，嚴格遵守競賽規(guī)則，以保證競賽的公正、公平性。如有違反競賽規(guī)則的行為，我們將受到嚴肅處理。我們參賽選擇的題號是（從A/B/C/D中選擇一項填寫）： 我們的參賽報名號為（如果賽區(qū)設(shè)置報名號的話）：所屬學(xué)校（請?zhí)顚懲暾娜簠①愱爢T(打印并簽名)：1.2.3.指導(dǎo)教師或指導(dǎo)教師組負責(zé)人(打印并簽名)：日期：年月日賽區(qū)評閱編號（由賽區(qū)組委會評閱前進行編號）：2011高教社杯全國大學(xué)生數(shù)學(xué)建模競賽編號專用頁賽區(qū)評閱編號（由賽區(qū)組委會評閱前進行編號）：賽區(qū)評閱記錄（可供賽區(qū)評閱時使用）：評閱人評分備注全國統(tǒng)一編號（由賽區(qū)組委會送交全國前編號）：全國評閱編號（由全國組委會評閱前進行編號）：城市表層土壤重金屬污染分析摘要本文以某城區(qū)為例，對其土壤地質(zhì)環(huán)境重金屬含量進行取樣分析，找出取樣點被污染的主要原因并尋找到取樣點的污染源。為此，對于第（1）和第（2）小題我們對取樣點的數(shù)據(jù)統(tǒng)計分析，用MATLAB軟件畫出8張元素的分布圖，各代表同一重金屬元素在不同功能區(qū)不同含量范圍內(nèi)的分布點，統(tǒng)計出點的數(shù)目，對均值和方差進行對比得出該城區(qū)內(nèi)不同區(qū)域重金屬的污染程度和主要污染原因。針對第（3）小題，我們建立了2個模型來研究，干塵擴散模型——模型一，高斯水流擴散模型——模型二。模型一是為了求解重金屬元素附著在干塵粒上的擴散方式和范圍，運用煙霧擴散模型，結(jié)論是根據(jù)取樣點的某金屬含量程度來判斷該點距污染源的距離，則污染源在以該距離為半徑的圓的圓周上，在根據(jù)取樣點的金屬含量高低分布情況，找到污染源。模型二是研究重金屬在土壤中的擴散情況，在相間物質(zhì)交換為平衡的條件下,可用阻滯系數(shù)來表示其影響，可得到一個微分方程帶入?yún)?shù)，可以找到某金屬的污染源坐標(biāo)。關(guān)鍵詞：煙霧擴散模型圖例統(tǒng)計污染物運移物質(zhì)交換MATLAB1問題的提出隨著城市經(jīng)濟的快速發(fā)展和城市人口的不斷增加，人類活動對城市環(huán)境質(zhì)量的影響日顯突出。對城市土壤地質(zhì)環(huán)境異常的查證，以及如何應(yīng)用查證獲得的海量數(shù)據(jù)資料開展城市環(huán)境質(zhì)量評價，研究人類活動影響下城市地質(zhì)環(huán)境的演變模式，日益成為人們關(guān)注的焦點。按照功能劃分，城區(qū)一般可分為生活區(qū)、工業(yè)區(qū)、山區(qū)、主干道路區(qū)及公園綠地區(qū)等，分別記為1類區(qū)、2類區(qū)、……、5類區(qū)，不同的區(qū)域環(huán)境受人類活動影響的程度不同。本文根據(jù)某個城市的5個區(qū)域的采樣數(shù)據(jù)分析出8種主要重金屬元素在該城區(qū)的空間分布，并分析該城區(qū)內(nèi)不同區(qū)域重金屬的污染程度。通過數(shù)據(jù)分析，說明重金屬污染的主要原因。分析重金屬污染物的傳播特征，由此建立模型，確定污染源的位置。附件1列出了采樣點的位置、海拔高度及其所屬功能區(qū)等信息，附件2列出了8種主要重金屬元素在采樣點處的濃度，附件3列出了8種主要重金屬元素的背景值?，F(xiàn)要求通過數(shù)學(xué)建模來完成以下任務(wù)：(1)給出8種主要重金屬元素在該城區(qū)的空間分布，并分析該城區(qū)內(nèi)不同區(qū)域重金屬的污染程度。(2)通過數(shù)據(jù)分析，說明重金屬污染的主要原因。(3)分析重金屬污染物的傳播特征，由此建立模型，確定污染源的位置。(4)分析你所建立模型的優(yōu)缺點，為更好地研究城市地質(zhì)環(huán)境的演變模式，還應(yīng)收集什么信息？有了這些信息，如何建立模型解決問題？2問題的分析對于問題和問題，主要采用圖形對照法，根據(jù)附件中的數(shù)據(jù)使用MATLAB圖形處理軟件圖形化，分析出同一重金屬元素在不同功能區(qū)及不同金屬元素在同一功能區(qū)的含量分布情況，找出污染的主要原因；而對于問題，為了更好的研究重金屬傳播特性和確定污染源位置，我們建立了兩種模型：模擬重金屬在空氣中的傳播，建立干塵擴散模型——模型一；考慮到金屬在土壤里的傳播，建立土壤擴散模型——模型二。在問題中，對模型的優(yōu)缺點進行了評價，提出了其它在模型中沒有考慮到的一些會影響重金屬污染程度的因素，并在模型一和模型二基礎(chǔ)上進一步優(yōu)化。問題一：金屬元素符號說明：砷鎘鉻銅汞鎳鉛鋅重金屬含量分段說明：、、、、、、、含量安全范圍分別為：、、、、、、、，我們分段時從均值開始到安全范圍上限，以后的范圍斷是人為分段（由于相關(guān)數(shù)據(jù)難以查到），具體分布可以參照表1經(jīng)過參考附錄中的坐標(biāo)點及海拔我們用軟件實現(xiàn)了該城區(qū)的三維模擬再根據(jù)已知數(shù)據(jù)用畫出俯瞰的功能區(qū)分布圖如圖2示，經(jīng)過圖形處理，得到該市大概的地區(qū)分布圖形如示：由砷的相關(guān)數(shù)據(jù)畫出各個區(qū)砷含量分布：由和知砷在部分工業(yè)區(qū)含量在以上，屬于重度污染，還有部分工業(yè)區(qū)砷含量在，交通區(qū)砷含量普遍在內(nèi)，屬于輕度污染；由鎘的相關(guān)數(shù)據(jù)得出該城市鎘的含量分布用比，可以得出該城市在山區(qū)鎘含量為，屬于輕度污染；在公園和生活區(qū)污染相對要嚴重，含鎘量在；在交通區(qū)含鎘量在；而工業(yè)區(qū)含鎘量則在之間，污染相當(dāng)嚴重。由鉻元素的數(shù)據(jù)可以得到鉻元素的分布：通過和對比，該城市普遍鉻元素含量在，在生活區(qū)含量在，而工業(yè)區(qū)和交通區(qū)污染則更嚴重，含量在之間；由附件中銅元素的數(shù)據(jù)可以畫出銅元素在該城市的分布：通過和的對比，發(fā)現(xiàn)在工業(yè)區(qū)和交通區(qū)銅元素含量在之間，其它功能區(qū)含量普遍低于；根據(jù)附件中汞元素的分布數(shù)據(jù)得出：與對比后可以看出工業(yè)區(qū)﹑交通區(qū)以及地處工業(yè)區(qū)附近的公園綠地區(qū)含量偏高，在之間，生活區(qū)和山區(qū)汞含量則低于；從附件中鎳元素的數(shù)據(jù)得到了鎳元素在該城區(qū)的分布：通過和比較，工業(yè)區(qū)和交通區(qū)以及生活區(qū)鎳含量偏高，在之間，而其它功能區(qū)則普遍低于；同樣的，由附件中的數(shù)據(jù)可以得到鉛和鋅元素在城市的分布和：通過和圖3對比鉛以及鋅在交通區(qū)﹑業(yè)區(qū)以及生活區(qū)含量高，鉛含量在80μg/g～480μg/g之間，鋅含量則是500μg/g～3800μg/g的范圍內(nèi)；上述圖對比為宏觀上的判斷，針對此，我們列出了同一重金屬元素在不同含量段在整個城區(qū)的分布，并統(tǒng)計出每個含量段的取樣點數(shù)，如表1As(μg/g)Cd(ng/g)Cr(μg/g)Cu(μg/g)含量范圍樣點數(shù)含量范圍樣點數(shù)含量范圍樣點數(shù)含量范圍樣點數(shù)3.6-5.488130-1906431-4914813.2-20.4585.4-1513190-50015949-1009920.4-15021115-30.135500-100039100-2005150-500101000-20007200-5003500-10000500-100021000-26002Hg(ng/g)Ni(μg/g)Pb(μg/g)Zn(μg/g)含量范圍樣點數(shù)含量范圍樣點數(shù)含量范圍樣點數(shù)含量范圍樣點數(shù)35-516112.3-19.916831-438669-977151-1509719.9-507043-8010997-16068150-5004150-100180-15053160-50094500-10007100-1501150-2509500-1500131000-16003250-48031500-38005表1通過表1的分析，、、、、在該城區(qū)屬于重度污染；、達到污染；分析，、、、、在該城區(qū)屬于重度污染；、達到污染；屬于輕度污染綜上，得出重金屬分布情況：八種重金屬元素在工業(yè)區(qū)、交通區(qū)含量較高，污染較嚴重；在生活區(qū)和公園綠地區(qū)含量相對要少，屬于輕度污染；而山區(qū)各種金屬元素含量除了接近工業(yè)區(qū)和交通區(qū)的區(qū)域，重金屬含量都在安全范圍內(nèi)。問題2：主要思路：分別對五個功能區(qū)重金屬元素的分布情況進行分析，找出該區(qū)某重金屬含量超過安全范圍的主要原因。通過均值和偏差的對比得討論如下：功能區(qū)均值砷（As）鎘（Cd）鉻（Cr）銅（Cu）汞（Hg）鎳（Ni）鉛（Pb）鋅（Zn）生活區(qū)7.48502.97118.16174.8953.6623.37116.22501.39工業(yè)區(qū)4.21323.8350.7950.4778.9814.6960.96264.86山區(qū)5.11224.9539.5925.0890.5214.9444.81104.51交通區(qū)5.83286.7243.9641.61312.8117.7157.92162.75公園綠地區(qū)5.39236.0538.9118.4147.4614.8741.0592.31標(biāo)準(zhǔn)偏差砷（As）鎘（Cd）鉻（Cr）銅（Cu）汞（Hg）鎳（Ni）鉛（Pb）鋅（Zn）生活區(qū)3.93281.54166.97415.523081.3519.4996.14547.46工業(yè)區(qū)1.48162.9745.2452.79146.886.1225.11607.79山區(qū)3.71249.7422.1321.75261.599.7933.31165.03交通區(qū)2.46173.1614.7829.891719.685.3531.21159.96公園綠地區(qū)2.39189.9016.338.0448.645.8812.4267.64對照圖和表所列數(shù)據(jù)，得以下的結(jié)論：1、生活區(qū)鉛和鋅的含量相對其他重金屬含量較高，考慮到人們的生活環(huán)境，如人們所用到的家具和房間的涂料，以及家用汽車的增加導(dǎo)致尾氣的排放、生活垃圾等都含有大量的鉛；含鋅化合物主要用作油漆、油墨、塑料、廣告顏料、化妝品、蓄電池的原料或填充劑。所以可以看出該市的人們在生活區(qū)使用的涂料、化妝品、電池和相關(guān)商品較多并且處理不的不完善，導(dǎo)致含鋅量較高。2、工業(yè)區(qū)所有重金屬的含量都很高。長期以來我國由于科技和經(jīng)濟欠發(fā)達，居民的環(huán)境保護意識薄弱，工業(yè)生產(chǎn)過程中產(chǎn)生的廢水、廢氣、廢渣未經(jīng)處理直接排放。其中含有大量的重金屬，經(jīng)過自然的沉降雨水的淋溶等途徑進入土壤導(dǎo)致Hg、Pb、Cu、Cr、As、Ni、Zn、Cd含量逐年增加導(dǎo)致工業(yè)區(qū)所有重金屬含量都很高。3、山區(qū)各重金屬含量大部分都在正常范圍內(nèi)，部分地區(qū)屬于輕度污染。主要原因是受到工業(yè)區(qū)和生活區(qū)的影響。工業(yè)區(qū)通過空氣傳播和水流傳播以及生活區(qū)的生活垃圾污染導(dǎo)致部分山區(qū)的輕度污染。4、交通區(qū)鉛和鎘、銅、鋅的含量較高。交通區(qū)汽車燃油、潤滑油和鍍金部分的燃燒或磨損釋放出大量重金屬元素。一般以道路為中心成帶狀分布，導(dǎo)致路邊土壤中這幾種重金屬元素明顯高于非交通線路兩側(cè)的土壤。5、公園綠地區(qū)汞、鋅、鉛的含量相對較高。首先考慮到該區(qū)所處的地理位置,地處生活區(qū)、工業(yè)區(qū)、交通區(qū)附近。所以公園較嚴重的重金屬污染可能來源于交通、工業(yè)污染和居民生活廢棄物等污染。問題3：一模型的假設(shè)模型一干塵擴散模型1、假設(shè)無風(fēng)情況下，重金屬看作在地面某一點向四周等強度的在無窮空間擴散，不計地形影響。2、重金屬物質(zhì)附著在干沉上的擴散服從擴散定律，即單位時間通過單位法平面面積的流量與它的濃度梯度成正比。3、重金屬物質(zhì)在穿過降雨時，其強度由于雨水的吸收而減少，減少比率為常數(shù)。4、風(fēng)的影響只考慮重金屬污染源的某一風(fēng)向、某一風(fēng)速的影響。風(fēng)的影響范圍，隨著時間而減小。模型二土壤水流擴散模型1、重金屬在土壤多孔介質(zhì)中有三種存在形式，假設(shè)相間物質(zhì)交換為平衡2、延滯系數(shù)唯一。二、模型的建立1.模型概述重金屬擴散的范圍和強度主要受到污染源的自然擴散、風(fēng)以及降雨的影響。自然擴散符合煙霧擴散模型，在給定污染源總含量等數(shù)據(jù)時，可以求出靜態(tài)擴散范圍。風(fēng)力具有加大擴散面積的作用，而降雨具有吸收稀釋重金屬元素含量的作用。2.參數(shù)示意污染源濃度污染源重金屬擴散的速度表示某點的濃度污染源重金屬擴散的總量擴散的時間降雨對重金屬的吸收率有風(fēng)時候的速度有風(fēng)時候的方向風(fēng)力影響隨時間減小量泄漏點距有效地面的高度(m)分別為用濃度標(biāo)準(zhǔn)差表示的軸上的擴散參數(shù)3.模型的求解求解主要思想：在無降雨情況下污染源重金屬總量等于擴散后的總量，所以有，當(dāng)有降雨的情況下，擴散后金屬含量有變化，吸收后剩余量為.3.1模型一建立過程建立空間坐標(biāo)系：將重金屬開始擴散時設(shè)為，以污染源為坐標(biāo)原點，建立以正東為x正方向，正北為y正方向，建立三維直角坐標(biāo)系，如下：3.1.1推導(dǎo)重金屬隨干沉擴散范圍半徑（無風(fēng)無降雨的擴散模型）時刻t無窮空間中任意一點的重金屬含量記為P。根據(jù)假設(shè)2，單位時間通過單位法相面的流量為：·····················s是擴散系數(shù)，grad表示梯度，負號表示由濃度高向濃度低的地方擴散?？疾炜臻g域Ω，Ω的體積為V，包圍Ω的曲面為S，S的外法線向量為n,則在[t,t+Δt]內(nèi)通過Ω的流量為：··················而Ω內(nèi)重金屬物質(zhì)含量的增量為：·····················由擴散前后金屬總量不變得：································根據(jù)曲面面積分的奧氏公式：·······················其中是散度記號。由可得：無界區(qū)域的拋物型偏微分方程······················根據(jù)假設(shè)1，初始條件為作用在坐標(biāo)原點的點源函數(shù)，可以記作：···························得：·······這個結(jié)果表明。對于任意時刻t，重金屬物質(zhì)濃度C的等值面是球面x2+y2+z2=R2,并且隨著球面半徑R的增加C的只是連續(xù)減少的；當(dāng)R—>∞或t—>∞時，C—>0通過將x2+y2+z2=R2代入方程（8）可以解出R：··············R為在給定擴散系數(shù)s，影響時間t，重金屬含量為C，以及污染源總含量W的情況下的受自然擴散影響范圍的半徑，即在不考慮風(fēng)力、降雨情況下的影響范圍。通過matlab，可以實現(xiàn)該模型的圖像生成，如圖一。在此假定：s=0.00001(擴散系數(shù))，W=3000（μg/g）（污染源總量，假定），C=103.20（重金屬在檢測點物質(zhì)濃度，這里取鉛），t=1000000(秒)（擴散時間）r=0.5(假設(shè)雨水對于重金屬吸收率為50%)ws=0.04(公里每小時)（假設(shè)作用在擴散源的風(fēng)速為40米/小時）wa=0.785（假設(shè)風(fēng)向與x0y面平行，且與x軸成正45度角）wd=0.3（假設(shè)風(fēng)力影響范圍從R開始，每天減弱300米）我們調(diào)用函數(shù)ZJ_Shu(k,Q,C,t,r,wv,we,wj)（函數(shù)本身可參考附錄），加入?yún)?shù)為ZJ_Shu(0.0001,3000,103.2,1000000,0.5,0.04,0.785,0.3)，符號意思看參數(shù)說明表，結(jié)果為：經(jīng)過去不同重金屬的均值，我們得到元素取值（平均值）半徑RAs3.624.9564×10^3Cd13015.3267×10^3Cr3118.1306×10^3Cu13.220.4128×10^3Hg3517.8477×10^3Ni12.320.6261×10^3Pb3118.1306×10^3Zn6916.4310×10^33.1.2考慮風(fēng)和降雨影響的模型建立降雨因素：假設(shè)降雨對于放射性物質(zhì)的吸收率為r,則方程（9）變?yōu)椋猴L(fēng)力因素：只考慮產(chǎn)生在輻射原點某一角度，某一強度的風(fēng)力影響。風(fēng)力的影響范圍從R開始每天減少，并且對于影響范圍內(nèi)每一點的作用相同。3.2模型結(jié)果分析模型的分析結(jié)果說明在給定的條件下，擴散的最遠距離為24.95公里。由于現(xiàn)實中風(fēng)向的不確定性，出于安全考慮，在污染源周邊24.95公里范圍的地域內(nèi)都有污染，而根據(jù)距離，可以判斷出污染源在以距離為半徑的圓周上，尋找到工業(yè)區(qū)或者是生活區(qū)。模型二：重金屬元素在土壤中的傳播模型分析：重金屬元素在土壤多孔介質(zhì)中有三種可能的存在形式：:溶于水中、揮發(fā)為氣體及吸附于固體顆粒。揮發(fā)性有機污染物在水、氣、固體顆粒三相間的物質(zhì)交換與分配是決定其運移的重要因素,在相間物質(zhì)交換為平衡的條件下,可用阻滯系數(shù)來表示其影響。污染物與土壤固體顆粒間的非平衡吸附解吸是相間交換中影響污染物運移的最重要的機制,由于天然土壤具有固有的不均質(zhì)性,必須用多個反應(yīng)系數(shù)才能準(zhǔn)確描述污染物與固體顆粒間的非平衡吸附解吸。在土壤環(huán)境中,一系列的機制控制著污染物的運移:一、地下水流決定了污染物的運動方向和速率;二、擴散使污染物產(chǎn)生縱向及橫向的轉(zhuǎn)移;三、污染物與土壤顆粒中有機質(zhì)及礦物質(zhì)之間的吸附解吸、污染物在土壤包氣帶中的水氣界面處的物質(zhì)交換使污染物的運移受到阻滯作用;四、由于具有揮發(fā)性,污染物還隨氣體遷移、擴散;五、土壤中的生物與化學(xué)作用使污染物降解,生成無害物質(zhì)或其他有害物質(zhì)。要預(yù)測污染物的運移和其歸宿,必須對土壤—污染物—空氣—水這一復(fù)雜的相互作用的系統(tǒng)及污染物在其中的諸多遷移機制有充分理解,并把該系統(tǒng)模型化從而得到污染源的位置。模型的建立符號說明氣體中的污染物濃度亨利系數(shù)水中的污染物濃度顆粒表面的污染物濃度污染物在水與固體顆粒間的分配系數(shù)分配系數(shù)土壤含水率土壤含氣率土壤干容重孔隙水流速度水平距離擴散系數(shù)比例系數(shù)距離傳播速度時間某已知點的橫坐標(biāo)某已知點的縱坐標(biāo)污染源的橫坐標(biāo)污染源的縱坐標(biāo)1氣—水—固三相間的物質(zhì)交換在有水氣共存的土壤中,即土壤非飽和帶中,污染物隨水運移時,水中的污染物會與氣體及土壤顆粒發(fā)生物質(zhì)交換。水氣兩相間的物質(zhì)交換平衡可用亨利定律來表示式中,為氣體中的污染物濃度,為水中的污染物濃度,H為亨利系數(shù)。不同物質(zhì)的H值可在文獻中查到。對于水中污染物與固體顆粒間的物質(zhì)交換,經(jīng)常用到以下這一簡單的線性方程:=式中,為固體顆粒表面的污染物濃度,為污染物在水與固體顆粒間的分配系數(shù)。有機物質(zhì)在水與固體顆粒間的分配與土壤顆粒中的有機質(zhì)含量、土壤顆粒的大小與表面積等因素直接相關(guān),而在天然土壤中,土壤的這些性質(zhì)往往不是均一的,用線性的方程來描述水—固間的物質(zhì)分配與實測數(shù)據(jù)往往不很吻合,而非線性方程則更好地表示了其間的物質(zhì)交換。方程是一個經(jīng)常被使用的非線性吸附等溫線方程,其表達式為:式中,為分配系數(shù),表示水—固吸附非線性指數(shù)。方程、都屬經(jīng)驗公式,其中的參數(shù)值可通過專門的實驗測出。污染物由氣相吸附在固體顆粒表面的量一般很小,可以忽略不計。假設(shè)氣相為靜止,且不考慮污染物在氣相中的擴散,污染物在非飽和帶中隨水運移可由下列一維對流—擴散方程來表示:式中,，分別為土壤含水率與含氣率,為土壤干容重,為孔隙水流速度,為水平距離,為擴散系數(shù)。方程左邊第二、第三項分別表示污染物在固相與氣相中的分配。在水—氣、水—固之間的物質(zhì)分配為平衡的條件下,可把方程、(或)代入,且方程兩邊同除以得:為延滯系數(shù),延滯系數(shù)Rt代表污染物運移時由于發(fā)生了水—氣及水—固體顆粒之間的物質(zhì)交換而產(chǎn)生的延滯作用。在不同的土壤含水率條件及不同的污染物運移承載相(氣或水)情況下,延滯系數(shù)的表達式是不同的,如考慮污染物隨氣相在非飽和帶中運移,則可由污染物通過土柱的時間與不發(fā)生水—固或氣—固吸附時非揮發(fā)性試劑通過土柱的時間的對比實驗來求得。近來有實驗表明,實驗測得的污染物在非飽和土柱中隨氣相運移時的延滯系數(shù)大于由求得的預(yù)測值。這部分額外的延滯效應(yīng)則應(yīng)歸于污染物在水—氣界面上的積累。污染物在水—氣界面上的積累不僅使污染物運移受阻,而且其影響還相當(dāng)大,如等在實驗中發(fā)現(xiàn)(三氯乙烯)的滯后是由于在水—氣界面處的積累造成的。因此,揮發(fā)性有機污染物在水—氣—固三相的分配可能包括以下幾個方面:污染物在水—氣兩相間的分配、污染物在水—固體顆粒兩相間的分配及污染物在水—氣界面處的累積。同時考慮這三個方面后的污染物隨氣體運移的延滯系數(shù)為:式中,為污染物在水氣界面處累積的吸附系數(shù),為有效水氣界面面積。由所求出的延滯系數(shù)可以看出和即濃度有關(guān)，我們可以假設(shè)延滯系數(shù)即擴展系數(shù)，顯然延滯系數(shù)和傳播速度成正比而從以污染源為起點延滯系數(shù)將逐漸降低，設(shè)污染源的位置為任意一地點（采樣的點）位置為。單位時間傳播的距離的積分就是污染源到這點的距離從而由得根據(jù)和已知數(shù)據(jù)便可以求出污染源的位置。用MATLAB求解并取代表性坐標(biāo)如元素坐標(biāo)As(μg/g)（5013，7202）Cd(ng/g)（2724，2576）（4034，6005）Cr(μg/g)（4632，4315）（3504，6013）Cu(μg/g)（2367，3443）（3845，6012）Hg(ng/g)（4812，1845）（4865，7011）Ni(μg/g)（4055，6089）Pb(μg/g)（2002，2951）（4990，5310）Zn(μg/g)（14024，9529）（3970，5978）（5372，7537）經(jīng)過數(shù)據(jù)分析，主要污染源在該城區(qū)的西南方向，集中了工業(yè)區(qū)，交通區(qū)，生活區(qū)，海拔在4——8之間。問題44模型的評價對于干塵擴散模型——模型一優(yōu)點:考慮了無風(fēng)無雨時候的理想情況和有風(fēng)有雨的情況，把現(xiàn)實生活中的大部分情況都包括在內(nèi)，并且得出了不同情況時候的傳播關(guān)系。缺點：考慮的因素較少，忽略了人為，濕度，自然環(huán)境的吸收能力等因素。得到的結(jié)果的準(zhǔn)確性自然較低。對于土壤擴散模型——模型二優(yōu)點：通過文獻，得到了直觀且簡單的計算方法和關(guān)系式。缺點：也只簡單的求出延滯系數(shù)，通過所求的式子得到各個量之間的關(guān)系，得到濃度和位置的一個關(guān)系式。對不同的重金屬的傳播途徑的分析和利用不夠透徹。為更好地研究城市地質(zhì)環(huán)境的演變模式，還應(yīng)收集的信息有:本地區(qū)不同地方的降水量、地基的復(fù)雜程度、人們的活動范圍以及將來社會的發(fā)展趨向（如哪里新建住房等）等。有了以上這些信息。我們認為應(yīng)綜合的根據(jù)以下各方面的原因分析。（1）降雨——降雨的多少必定影響這個地區(qū)的土壤和污染等；（2）地基的復(fù)雜程度——地基的不同也對這個地區(qū)的承受能力和變化有影響；（3）人們的活動范圍——人們的活動對這個地方的演變起到重要的作用，人多的地方演變的必定更快。（4）將來社會的發(fā)展趨向——社會的發(fā)展趨向?qū)Νh(huán)境的影響更是至關(guān)重要，直接影響人們活動的地域，從而影響環(huán)境。5模型的推廣1、本模型可以擴展到其他污染物的污染，只需檢測出相應(yīng)相應(yīng)的污染物在各個區(qū)域的含量分布，就可以找出污染源；2、本模型對擴散型污染物都有一定的適用性，可以參照此解決一些擴散性污染的問題；3、本模型具有較強的應(yīng)用價值。環(huán)境的污染容易破壞生態(tài)鏈，從而影響人類的生活，重金屬的污染，甚至可以直接危害人類的健康。該模型可以通過控制一些因素找出被污染城市的污染源，從而及時制定應(yīng)急措施，減小損失。參考文獻楊桂元，數(shù)學(xué)模型應(yīng)用實例，合肥：工業(yè)大學(xué)出版社，2007年6月第一版汪國強，數(shù)學(xué)建模優(yōu)秀案例選編，華南：華南理工大學(xué)出版社，2001年8月第一版周博謝東來，MATLAB科學(xué)計算，三河市：機械工業(yè)出版社，2010年5月第一版附錄：圖形實現(xiàn)主代碼A=[......];%此矩陣數(shù)據(jù)太多，在此省略bianhao=A(:,1);%編號As=[A(:,2)];Cd=[A(:,3)];Cr=[A(:,4)];Cu=[A(:,5)];Hg=[A(:,6)];Ni=[A(:,7)];Pb=[A(:,8)];Zn=[A(:,9)];x=[A(:,10)];y=[A(:,11)];z=[A(:,12)];gnqu=[A(:,13)];X_hang=length(gnqu);%scatter(x,y,5,z)%散點圖[X,Y,Z]=griddata(x,y,z,linspace(0,28654)',linspace(0,18449),'v4');%插值%pcolor(X,Y,Z);%shadinginterp%偽彩色圖%所有監(jiān)測點網(wǎng)絡(luò)圖fig1=figure(1);set(fig1,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdonplot(x,y,'-*')title('檢測點網(wǎng)絡(luò)圖');xlabel('x');ylabel('y');fig2=figure(2);set(fig2,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdon%Ask=0;fori=1:X_hangif3.6<=As(i)&As(i)<=5.4plot(x(i),y(i)+250,'^');holdonelseif5.4<As(i)&A(i)<=15plot(x(i)+250,y(i),'s','markerfacecolor','g');holdonelseifAs(i)>15plot(x(i),y(i)-250,'d','markerfacecolor','r');text(x(i)+200,y(i)+200,sprintf('(%0.1f)',As(i)),'color','r','fontsize',10,'fontweight','bold');holdonendendendendtitle('As的分布','backgroundcolor','y','color','b','fontsize',12,'fontweight','bold');xlabel('x');ylabel('y');text(1000,17900,'

-重度污染As(μg/g)>15','backgroundcolor','k','color','w');text(1000,16900,'△-輕度污染3.6<=As(μg/g)<=5.4','backgroundcolor','k','color','w');text(1000,15900,'□-污染5.4<As(μg/g)<=15','backgroundcolor','k','color','w');%Cdfig3=figure(3);set(fig3,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdonfori=1:X_hangifCd(i)<130&Cd(i)<190plot(x(i),y(i),'v','markerfacecolor','k');%text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Cd(i)),'color','k','fontsize',10,'fontweight','bold');holdonelseif190<=Cd(i)&Cd(i)<=500plot(x(i),y(i),'p','markerfacecolor','g','markersize',10);%text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Cd(i)),'color','k','fontsize',10,'fontweight','bold');holdonelseif500<Cd(i)&Cd(i)<=1000plot(x(i),y(i),'o','markerfacecolor','r');%text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Cd(i)),'color','r','fontsize',10,'fontweight','bold');holdonelseif1000<Cd(i)&Cd(i)<2000plot(x(i),y(i),'d','markerfacecolor','r','markersize',12);text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Cd(i)),'color','r','fontsize',10,'fontweight','bold');holdonendendendendendtitle('Cd的分布','color','b','backgroundcolor','y','fontsize',12,'fontweight','bold');xlabel('x');ylabel('y');text(100,17900,'▲-130<Cd(ng/g)<190','color','w','backgroundcolor','k');text(100,16900,'★-190<=Cd(ng/g)<=500','color','w','backgroundcolor','k');text(100,15900,'●-500<Cd(ng/g)<=1000','color','w','backgroundcolor','k');text(100,14900,'◆-1000<Cd(ng/g)<2000','color','w','backgroundcolor','k');%Crfig4=figure(4);set(fig4,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdonfori=1:X_hangif31<=Cr(i)&Cr(i)<=49plot(x(i),y(i),'v','markerfacecolor','k');holdonelseif49<Cr(i)&Cr(i)<=100plot(x(i),y(i),'p','markerfacecolor','g');holdonelseif100<Cr(i)&Cr(i)<=200plot(x(i),y(i),'o','markerfacecolor','b');holdonelseif200<Cr(i)&Cr(i)<=500plot(x(i),y(i),'d','markerfacecolor','r','markersize',10);holdonelseif500<Cr(i)&Cr(i)<=1000plot(x(i),y(i),'s','markerfacecolor','r','markersize',12);text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Cr(i)),'color','r','fontsize',10,'fontweight','bold');holdonendendendendendendtext(100,17900,'▲-32<Cr(μg/g)<49','color','w','backgroundcolor','k');text(100,16900,'★-49<=Cr(μg/g)<=100','color','w','backgroundcolor','k');text(100,15900,'●-100<Cr(μg/g)<=200','color','w','backgroundcolor','k');text(100,14900,'◆-200<Cr(μg/g)<500','color','w','backgroundcolor','k');text(100,13900,'■-500<Cr(μg/g)<1000','color','w','backgroundcolor','k');title('Cr的分布','color','b','backgroundcolor','y','fontsize',12,'fontweight','bold');xlabel('x');ylabel('y');%Cufig5=figure(5);set(fig5,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdonfori=1:X_hangif13.2<=Cu(i)&Cu(i)<=20.4plot(x(i),y(i),'v','markerfacecolor','k');holdonelseif20.4<Cu(i)&Cu(i)<=150plot(x(i),y(i),'p','markerfacecolor','g');holdonelseif150<Cu(i)&Cu(i)<=500plot(x(i),y(i),'o','markerfacecolor','b');holdonelseif500<Cu(i)&Cu(i)<=1000plot(x(i),y(i),'d','markerfacecolor','r','markersize',10);holdonelseif1000<Cu(i)&Cu(i)<=2600plot(x(i),y(i),'s','markerfacecolor','r','markersize',12);text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Cu(i)),'color','r','fontsize',10,'fontweight','bold');holdonendendendendendendtext(100,17900,'▲-13.2<Cu(μg/g)<20.4','color','w','backgroundcolor','k');text(100,16900,'★-20.4<=Cu(μg/g)<=150','color','w','backgroundcolor','k');text(100,15900,'●-150<Cu(μg/g)<=500','color','w','backgroundcolor','k');text(100,14900,'◆-500<Cu(μg/g)<1000','color','w','backgroundcolor','k');text(100,13900,'■-1000<Cu(μg/g)<2600','color','w','backgroundcolor','k');title('Cu的分布','color','b','backgroundcolor','y','fontsize',12,'fontweight','bold');xlabel('x');ylabel('y');%Hgfig6=figure(6);set(fig6,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdonfori=1:X_hangif35<=Hg(i)&Hg(i)<=51plot(x(i),y(i),'v','markerfacecolor','k');holdonelseif51<Hg(i)&Hg(i)<=150plot(x(i),y(i),'p','markerfacecolor','g');holdonelseif150<Hg(i)&Hg(i)<=500plot(x(i),y(i),'o','markerfacecolor','b');holdonelseif500<Hg(i)&Hg(i)<=1000plot(x(i),y(i),'d','markerfacecolor','r','markersize',10);holdonelseif1000<Hg(i)&Hg(i)<=1600plot(x(i),y(i),'s','markerfacecolor','r','markersize',12);text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Hg(i)),'color','r','fontsize',10,'fontweight','bold');holdonendendendendendendtext(100,17900,'▲-35<Hg(ng/g)<51','color','w','backgroundcolor','k');text(100,16900,'★-51<=Hg(ng/g)<=150','color','w','backgroundcolor','k');text(100,15900,'●-150<Hg(ng/g)<=500','color','w','backgroundcolor','k');text(100,14900,'◆-500<Hg(ng/g)<1000','color','w','backgroundcolor','k');text(100,13900,'■-1000<Hg(ng/g)<=1600','color','w','backgroundcolor','k');title('Hg的分布','color','b','backgroundcolor','y','fontsize',12,'fontweight','bold');xlabel('x');ylabel('y');%Nifig7=figure(7);set(fig7,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdonfori=1:X_hangif12.3<=Ni(i)&Ni(i)<=19.9plot(x(i),y(i),'v','markerfacecolor','k');holdonelseif19.9<Ni(i)&Ni(i)<=50plot(x(i),y(i),'p','markerfacecolor','g');holdonelseif50<Ni(i)&Ni(i)<=100plot(x(i),y(i),'o','markerfacecolor','b');holdonelseif100<Ni(i)&Ni(i)<=150plot(x(i),y(i),'d','markerfacecolor','r','markersize',10);text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Ni(i)),'color','r','fontsize',10,'fontweight','bold');holdonendendendendendtext(100,17900,'▲-12.3<Ni(μg/g)<19.9','color','w','backgroundcolor','k');text(100,16900,'★-19.9<=Ni(μg/g)<=50','color','w','backgroundcolor','k');text(100,15900,'●-50<Ni(μg/g)<=100','color','w','backgroundcolor','k');text(100,14900,'◆-100<Ni(μg/g)<150','color','w','backgroundcolor','k');title('Ni的分布','color','b','backgroundcolor','y','fontsize',12,'fontweight','bold');xlabel('x');ylabel('y');%Pbfig8=figure(8);set(fig8,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdonfori=1:X_hangif31<=Pb(i)&Pb(i)<=43plot(x(i),y(i),'v','markerfacecolor','k');holdonelseif43<Pb(i)&Pb(i)<=80plot(x(i),y(i),'p','markerfacecolor','g');holdonelseif80<Pb(i)&Pb(i)<=150plot(x(i),y(i),'o','markerfacecolor','b');holdonelseif150<Pb(i)&Pb(i)<=250plot(x(i),y(i),'d','markerfacecolor','r','markersize',10);holdonelseif250<Pb(i)&Pb(i)<=480plot(x(i),y(i),'s','markerfacecolor','r','markersize',12);text(x(i)+200,y(i)+200,sprintf('(%0.1f)',Pb(i)),'color','r','fontsize',10,'fontweight','bold');holdonendendendendendendtext(100,17900,'▲-31<Pb(μg/g)<43','color','w','backgroundcolor','k');text(100,16900,'★-43<=Pb(μg/g)<=80','color','w','backgroundcolor','k');text(100,15900,'●-80<Pb(μg/g)<=150','color','w','backgroundcolor','k');text(100,14900,'◆-150<Pb(μg/g)<250','color','w','backgroundcolor','k');text(100,13900,'■-250<Pb(μg/g)<=480','color','w','backgroundcolor','k');title('Pb的分布','color','b','backgroundcolor','y','fontsize',12,'fontweight','bold');xlabel('x');ylabel('y');%Znfig9=figure(9);set(fig9,'menubar','none','position',[2551001000600]);contourf(X,Y,Z);%等高線圖holdonfori=1:X_hangif69<=Zn(i)&Zn(i)<=97plot(x(i),y(i),'v','markerfacecolor','k');holdonelseif97<Zn(i)&Zn(i)<=160plot(x(i),y(i),'p','markerfacecolor','g');holdonelseif160<Zn(i)&Zn(i)<=500plot(x(i),y(i),'o','markerfacecolor','b');holdonelseif500<Zn(i)&Zn(i)<=1500plot(x(i),y(i),'d','markerfacecolor','r','markersize',10);holdonelseif1