【正文】 for k=1:n yc0(k+1)=ycl(k+1)ycl(k)。 ycl(1)=x0(1)。 y(i)=x0(i+1)。n=length(x0)。hold on。,linspace(0,18449),39。textdecoration:none} 27 附錄 7 由于 8 種重金屬的空間分布 圖的程序都差不多，下面以 As 的空間分布圖為例 . x=data(:,1)。fontfamily:Monospaced。D:\課 件 \畢業(yè)論文 \A\濃度表 .xls39。fontsize:13pt。 /SHEET=name 39。 P綜 =sqrt((P ^2+ maxiP ^2)/2) 公園路地區(qū)： P = mean(a(5,:))。 P綜 =sqrt((P ^2+ maxiP ^2)/2) 山區(qū)： P = mean(a(3,:))。) 25 附錄 4 矩陣 a 表示下表中各種重金屬在不同區(qū)域中的單項(xiàng)污染指數(shù) 元素 As Cd Cr Cu Hg Ni Pb Zn 生活區(qū) 工業(yè)區(qū) 山區(qū) 交通區(qū) 公園綠地區(qū) 生活區(qū)： P = mean(a(1,:))。 y5=data5 (:,2)。 z4=data4(:,3)。 scatter3(x3,y3,z3,39。ro39。),hold on x2= data2(:,1)。 y1=data1 (:,2)。v439。 scatter3(x5,y5,z5,39。mp39。) ,hold on x4= data4(:,1)。 y3=data3 (:,2)。 z2=data2 (:,3)。 contourf(X,Y,Z) ,hold on scatter3(x1,y1,z1,39。 [X,Y,Z]=griddata(x,y,z,linspace(0,27823)39。)。 z=data (:,3)。三次多項(xiàng)式插值法分別對(duì) 8種重金屬污染物的濃度進(jìn)行數(shù)據(jù)插值，編程實(shí)現(xiàn) (附錄 7).畫出圖像后，觀察圖像的極值 點(diǎn)，根據(jù)附件 2的數(shù)據(jù)，利用 MATLAB的搜索函數(shù) find，先從重金屬的最高濃度開始查找，然后找到該點(diǎn)對(duì)應(yīng)的 x， y坐標(biāo)，最后把找到的坐標(biāo)跟圖像對(duì)比，符合的就是圖像的極值點(diǎn)，找到第一點(diǎn)后，再用次高濃度查找，直到找到所有的極值點(diǎn)為止 . 2. 模型的結(jié)果 圖 5 As的濃度圖 從圖 5中，可以看出 As的濃度圖中分別有 2個(gè)極值點(diǎn)，可以大概確定極值點(diǎn)的位置為： 14 (18134， 10046)、 (12696， 3024). 圖 6 Cd的濃度圖 從圖 6中，可以看出 Cd的濃度圖中分別有 1個(gè)極值點(diǎn)，可以大概確定極值點(diǎn)的位置為：(21439， 11383). 圖 7 Cr的濃度圖 從圖 7中，可以看出 Cr的濃度圖中分別有 2個(gè)極值點(diǎn)，可以大概確定極值點(diǎn)的位置為：(3299， 6018)、 (4592， 4603). 15 圖 8 Cu的濃度圖 從圖 8中，可以看出 Cu的濃度圖中分別有 2個(gè)極值點(diǎn)，可以大概確定極值點(diǎn)的位置為：(2383， 3692)、 (3299， 6018). 圖 9 Hg的濃度圖 從圖 9中，可以看出 Hg的濃度圖中分別有 4個(gè)極值點(diǎn)，可以大概確定極值點(diǎn)的 位置為：(2383， 3692)、 (2708， 2295)、 (13694， 2357)、 (15248， 29106). 16 圖 10 Ni的濃度圖 從圖 10中，可以看出 Ni的濃度圖中分別有 1個(gè)極值點(diǎn)，可以大概確定極值點(diǎn)的位置為： (3299， 6018). 圖 11 Pb的濃度圖 從圖 11中，可以看出 Pb的濃度圖中分別有 2個(gè)極值點(diǎn)，可以大概確定極值點(diǎn)的位置為： (1647， 2728)、 (4777， 4897). 17 圖 12 Zn的濃度圖 從圖 11中，可以看出 Zn的濃度圖中分別有 5個(gè)極值點(diǎn)，可以大概確定極值點(diǎn)的位置為： (9328， 4311)、 (13797， 9621)、 (12696， 3024)、 (4948， 7293)、 (2383， 3692). 根據(jù)各種重金屬濃度圖中的極值點(diǎn)，可以確定各種重金屬污染源的位置，如表 9. 表 9 各種重金屬污染源的位置 元素 污染源位置 /m As (18134， 10046)、 (12696， 3024) Cd (21439， 11383) Cr (3299， 6018)、 (4592， 4603) Cu (2383， 3692)、 (3299， 6018) Hg (2383， 3692)、 (2708， 2295)、 (13694， 2357)、 (15248， 29106) Ni (3299， 6018) Pb (9328， 4311)、 (13797， 9621) Zn (9328， 4311)、 (13797， 9621)、 (12696， 3024)、 (4948， 7293)、(2383， 3692) 由上表可以看出， As和 Zn有共同的污染源，來自交通區(qū)的 (12696， 3024)； Cr、 Cu、Ni有共同的污染源，來自交通區(qū)的 (3299， 6018)； Cu、 Hg、 Zn有共同的污染源，來自工 18 業(yè)區(qū)的 (2383， 3692)； Pb、 Zn有共同的污染源，來自生活區(qū)的 (9328， 4311)和交通區(qū)的(13797， 9621). [5] 根據(jù)重金屬在土壤中的擴(kuò)散特點(diǎn)，該問題的解決還可以選用高斯擴(kuò)散模型進(jìn)行建模 . 污染源選為坐標(biāo)原點(diǎn)，污染處重金屬元素的濃度記為 ( , , ,0)C x y z ，則時(shí)刻 s 處的無窮空間中的任一點(diǎn)重金屬濃度為 ( , , , )i i i iC x y z s ，單位距離內(nèi)通過單位法向面積的流量與濃度梯度成正比，則 iq gradC??? ? ? ?=1,2,3,i i? ??? 是擴(kuò)散系數(shù)， grad表示梯度，負(fù)號(hào)表示由濃度高向濃度低的地方擴(kuò)散 . 考察空間域 ? ，其體積為 V ，包圍 ? 的曲面為 1S ， 1S 為一規(guī)則的球面， 1S 外法線向量為 = , ,1xynzz?????? . 則在 ? ?,s s s?? 內(nèi)通過 ? 的流量 ? ?11 . 1sss sQ q n d d s???? ? ?? ? 內(nèi)重金屬元素的增量為 ? ? ? ? ? ?2 , , , , , , 2VQ C x y z s s C x y z s d V? ? ???????? 又由于污染物釋放的重金屬的總量為 ? ?00 3sssQ p d V d t?? ?? ? ??? 由質(zhì)量守恒定律，有 ? ?0 1 2 4Q Q Q?? 19 根據(jù)曲面面積的 Gauss積分公式有 ? ?1 . 5sVq n d d iv q d V? ??? ??? 其中 div 表示散度記號(hào)，由 ??1 式至 ??5 式并利用積分中值定理整理得 ? ? ? ?2 2 22 2 2 , 0 , , , 6i i y zC C C Cd iv g r a d C s x y zs x y z? ? ? ???? ? ? ?? ? ? ? ? ? ? ? ???? ? ? ??? 這是無界區(qū)域的偏微分方程，根據(jù)假設(shè) (5)，初始條件為作用在坐標(biāo)原點(diǎn)的污染 源函數(shù)，可記作 ? ?0( , , , 0 ) ( , , ) 7C x y z Q x y z?? 0Q 表示釋放的重金屬元素總量， ? ?,x y z? 是單位強(qiáng)度的擴(kuò)散系數(shù) .方程 ??6 滿足方程 ??7 的解為 ? ? ? ? ? ?2 2 20 3 / 2, , , , , ,44iiQ x y zC x y z s e i x y zksks ??? ???? 其中 k 為 t 時(shí)間內(nèi)傳播的速度 . 問題 4 的分析與求解 (1)模型一 8種重金屬元素空間分布圖比較清晰、直觀，分析也比較詳細(xì)；