基于高斯模型的氣體泄漏源定位反算(上傳)

基于高斯模型的氣體泄漏源定位反算算法研究摘要：隨著經(jīng)濟(jì)發(fā)展和社會(huì)進(jìn)步，?；返氖褂眯枨笤絹碓酱?，隨之而來的危化品泄漏事故不斷增多，給人民群眾的生命財(cái)產(chǎn)帶來了重大威脅。在此類事故的救援處置中，及時(shí)準(zhǔn)確的泄漏源定位是順利實(shí)施救援行動(dòng)的前提條件。本文利用遺傳算法實(shí)現(xiàn)基于高斯氣體擴(kuò)散模型的反算定位，并利用matlab軟件進(jìn)行仿真試驗(yàn)，得出算法在較高的迭代次數(shù)和獲得較多初始數(shù)據(jù)的基礎(chǔ)上能夠?qū)崿F(xiàn)泄漏源反算定位。關(guān)鍵詞：危化品泄漏源定位遺傳算法仿真試驗(yàn)氣態(tài)?；酚捎谝兹家妆?、有毒腐蝕等特性造成的危害極大，在處置此類事故中，泄漏源的準(zhǔn)確定位是保證滅火救援行動(dòng)順利實(shí)施的前提條件。由于泄漏氣體在大氣中的分布是一個(gè)濃度梯度模型，結(jié)合該模型和遺傳算法能夠較好的反算泄漏源的位置，便于消防指戰(zhàn)員快速采取行動(dòng)，準(zhǔn)確及時(shí)控制事故蔓延。一、擴(kuò)散模型為實(shí)現(xiàn)泄漏源的成功定位，必須選用合適的氣體擴(kuò)散模型［1］。現(xiàn)如今的氣體擴(kuò)散理論已經(jīng)比較成熟，根據(jù)泄漏物質(zhì)狀態(tài)一般可分為重氣泄漏和非重氣泄漏［2］，重氣泄漏是指氣體密度大于空氣密度的氣體泄漏，如二氧化碳、氯氣等［2］；非重氣泄漏是指密度與空氣相差不大或者經(jīng)短時(shí)間空氣稀釋后密度與空氣接近的氣體泄漏，如氨氣等［3］。描述氣體擴(kuò)散的模型發(fā)展至今已經(jīng)有很多，但是常見的模型主要有：高斯模型（高斯煙羽模型、高斯煙團(tuán)模型）、BM模型、Sutton模型、FEM3模型，下面對這幾種模型進(jìn)行簡要的描述和分析［4］。（1）高斯模型高斯模型主要應(yīng)用于非重氣擴(kuò)散，它在二十世紀(jì)五、六十年代就開始被廣泛應(yīng)用，它的理論基礎(chǔ)是是以擴(kuò)散物質(zhì)在空間正態(tài)分布，從統(tǒng)計(jì)理論出發(fā)，進(jìn)而得出擴(kuò)散物質(zhì)的濃度分布，可以分為高斯煙羽模型和高斯煙團(tuán)模型，各自有各自的應(yīng)用范圍，其中高斯煙羽模型適用于連續(xù)泄漏或者泄漏事件相對于擴(kuò)散時(shí)間較長的情況，高斯煙團(tuán)模型適用于瞬時(shí)泄漏或泄漏時(shí)間相對于擴(kuò)散時(shí)間較短的情況［5］。（2）BM模型BM模型是一種經(jīng)驗(yàn)?zāi)Ｐ?，由一系列重氣瞬時(shí)擴(kuò)散和連續(xù)泄漏的試驗(yàn)數(shù)據(jù)繪制成的計(jì)算圖表組成，部分學(xué)者對數(shù)據(jù)進(jìn)行無因次處理擬合成解析公式，該模型只能夠適用于重氣泄漏情況下的模擬計(jì)算［5］，由于該模型是從實(shí)驗(yàn)室數(shù)據(jù)推導(dǎo)出來，因而在實(shí)際使用中誤差較大。（3）Sutton模型Sutton模型使用湍流擴(kuò)散統(tǒng)計(jì)理論來處理湍流擴(kuò)散問題，該模型在模擬可燃性氣體的泄漏擴(kuò)散中，會(huì)產(chǎn)生較大誤差［6］。（4）FEM3模型FEM3（3DFiniteElementModel）模型是1979年提出的三位有限源計(jì)算模型，由于其計(jì)算復(fù)雜因而發(fā)展初期緩慢，但是隨著計(jì)算機(jī)的發(fā)展該模型獲得廣泛應(yīng)用，可以處理重氣擴(kuò)散以及處理連續(xù)源泄漏和有限時(shí)間內(nèi)的泄漏，在復(fù)雜地形條件下也可以使用，模型的特點(diǎn)是計(jì)算量大、求解困難、計(jì)算結(jié)果精度較高。他的理論基礎(chǔ)是使用梯度輸運(yùn)理論和混合長理論（K理論）來解決湍流問題。下面對以上模型進(jìn)行歸納比較分析，如下表所示：表1擴(kuò)散模型介紹擴(kuò)散模型特點(diǎn)局限精度適用于連續(xù)泄漏源或是釋放時(shí)間大高斯煙羽模型于等于擴(kuò)散時(shí)間的情形地勢平坦，風(fēng)速大于lm/s較好適用于瞬時(shí)泄漏源或是釋放時(shí)間相高斯煙團(tuán)模型對擴(kuò)散時(shí)間比較短的情形地勢平坦，風(fēng)速大于1m/s較好BM模型簡單易用計(jì)算規(guī)模大較差Sutton模型簡單易用計(jì)算規(guī)模大較差3維模型，計(jì)算精度高，可模擬復(fù)FEM3模型雜地形條件下的氣體擴(kuò)散求解困難，計(jì)算量大較好根據(jù)以上分析結(jié)果，高斯模型（高斯煙羽模型、高斯煙團(tuán)模型）不僅適用范圍廣、計(jì)算規(guī)模小、精度較好、理論完善、易于理解分析［7］，而且提出時(shí)間較早、實(shí)驗(yàn)數(shù)據(jù)多、計(jì)算方便、理論成熟，因此本文采用高斯模型作為研究對象。在實(shí)際情況下，重性氣體在空間中的擴(kuò)散由于受到大氣的稀釋也能夠很好的滿足高斯模型。本文采用的高斯模型為：Pasquill-Gifford模型。對于高斯模型的提出還需要提出以下適用條件：地勢平坦開闊，性質(zhì)均勻，定常態(tài)；氣體擴(kuò)散過程中沒有發(fā)生沉降、化合、分解及地面吸收的發(fā)生，沒有降水等情況；擴(kuò)散空間的風(fēng)速、大氣穩(wěn)定度都均勻穩(wěn)定；云團(tuán)與空氣相對運(yùn)動(dòng)忽略不計(jì)，即云團(tuán)與大氣一起運(yùn)動(dòng)，平均風(fēng)速不小于lm/s；適用尺度不超過20km。泄漏源在地面的投影點(diǎn)為坐標(biāo)原點(diǎn)，x軸為風(fēng)向方向，z軸垂直地面向上，y軸為側(cè)風(fēng)方向，Pasquill-Gifford擴(kuò)散方程如下所示：①瞬時(shí)泄漏高斯煙團(tuán)公式

Q3(2n)32bbbxyzI1(x-utQ3(2n)32bbbxyzI1(x-ut)2exp(_2［b2xrh\、］〕+exp一(、＞2b2J亡［啊_1(詈)z1.1)連續(xù)泄漏煙羽模式C(x,y,z)「1,(z-H)2、［「1.Q+H)2、］exp-r—-)+exp一(-、2b2z2b2z一exp(-1出)2兀ybb2b2yzy1.2)其中：C(x,y,z,為時(shí)刻t時(shí)(x,y,z)處濃度；Q為源強(qiáng)；為風(fēng)速；H為泄漏源有效高度;”分別為下風(fēng)、側(cè)風(fēng)、垂直方向擴(kuò)散系數(shù)；大氣擴(kuò)散系數(shù)是求解擴(kuò)散模型的關(guān)鍵因素，是表征大氣對漂浮在大氣中的污染物質(zhì)擴(kuò)散能力的物理量，與大氣穩(wěn)定度和地面的有效粗糙度有關(guān)。大氣穩(wěn)定度指空氣中某氣團(tuán)由于與周圍空氣存在密度、溫度和流速強(qiáng)度差而產(chǎn)生的浮力使其產(chǎn)生加速度而上升或下降的程度［8；］地表粗糙度反映地表地形、建筑物風(fēng)流影響程度的參數(shù)，地面障礙物越低、密集度越低，地面有效粗糙度長度越小。二、反算算法算法實(shí)現(xiàn)的基本思路為通過擴(kuò)散模型對事故進(jìn)行模擬，將模型的模擬計(jì)算結(jié)果與測量值進(jìn)行比較，建立目標(biāo)函數(shù)，通過遺傳算法對目標(biāo)函數(shù)進(jìn)行最優(yōu)化求解，最終找到泄漏源位置最優(yōu)解，即所求結(jié)果。首先對問題進(jìn)行簡化，計(jì)算中涉及到四個(gè)未知量：泄漏源三維坐標(biāo)、源強(qiáng)，前人對于源強(qiáng)的處理通常采用賦予一定初值的方法，在實(shí)際的情況下，泄漏源強(qiáng)的數(shù)據(jù)獲取十分困難，在泄漏源位置未知的情況下很難確定，當(dāng)作未知量來處理，則會(huì)大大增加計(jì)算量。以連續(xù)泄漏煙羽模式為例，對方程進(jìn)行降維處理。忽略節(jié)點(diǎn)設(shè)置高度的影響以及利用數(shù)學(xué)方法消去源強(qiáng)未知量：CCx,y,0)=exp(CCx,y,0)=exp(兀ybbyzy2b2)-expH2

_-(—)b2z1.3)假設(shè)在地面任意設(shè)置兩個(gè)點(diǎn)：C(x,y,0)111C(x,y,0)111兀ybby1z1exp(1y2—-2b2y1)-exp1H22(hz11.4)兩式相除得：C(x,y,0)=Qexp(-1~y兩式相除得：C(x,y,0)=Qexp(-1~y2-222和bb2b2y2z2y2)?expH2

(b2z21.5)2b?b~y2Z2b?by1z1H2+_?b2y2z2丄工1H?1—?2b22b2y1z1i1.6)1.7)i+1采用此種方法至少需要三個(gè)有效的偵檢節(jié)點(diǎn)才能進(jìn)行泄漏源定位。（1）目標(biāo)函數(shù)的建立1.8)1.8)min工(T-T)2監(jiān)測模擬i=1其中n為監(jiān)測數(shù)據(jù)個(gè)數(shù)。使用不同的目標(biāo)函數(shù)對計(jì)算結(jié)果的精度會(huì)造成影響，比如均方誤差、開平方誤差、開四次方誤差等等。還有一些目標(biāo)函數(shù)如下所示：幾何方差：VG=exp(logT幾何方差：VG=exp(logT監(jiān)測-lOgT模擬）"1.9)歸一化的均方誤差（NRMSE）：1.10)1.111.10)1.11)l（T-T）2NRMSE=測計(jì)~\T?T測計(jì)均方誤差（RMSE）:一l(T-T)2RMSE=j測計(jì)~~■t2測2)遺傳算法遺傳算法是一種解決優(yōu)化問題的全局優(yōu)化算法，最初是由美國密歇根大學(xué)HollandJ.H.教授提出［9］，其思想基礎(chǔ)是自然界生物進(jìn)化理論適者生存的規(guī)律，在專著《自然與人工系統(tǒng)適應(yīng)性行為》中詳細(xì)介紹了遺傳算法，后人對該算法進(jìn)行了改進(jìn)，使之在更多領(lǐng)域得到了推廣［10］?；具^程為：將解集的群體看成一系列帶有染色體特征的種群，染色體帶有遺傳信息，用符號(hào)串標(biāo)識(shí)，適應(yīng)度高的個(gè)體基因更容易被繼承到下一代，基因通過一系列選擇、交叉、

變異操作進(jìn)行變異，最終找到最優(yōu)個(gè)體，對該個(gè)體解碼得到最優(yōu)解，基本步驟包括以下幾個(gè)I種群初始化與編碼第一代種群隨機(jī)獲取，組成問題的初始解，遺傳算法的編碼方式有二進(jìn)制編碼、實(shí)數(shù)編碼、格雷碼編碼、十進(jìn)制編碼等方式，其中二進(jìn)制編碼比較常見，這種編碼方式簡單、方便[11]。其中(x,y)分別用20位二進(jìn)制字符進(jìn)行編碼，h用10位二進(jìn)制字符進(jìn)行編碼。v=cell(1,m);fori=1:mtemp1=round(rand(1,50));temp2=vec2str(temp1);v{i}=temp2;endII個(gè)體評價(jià)適應(yīng)度函數(shù)是度量某個(gè)物種對于生存環(huán)境的適應(yīng)程度，對于生存環(huán)境適應(yīng)程度較高的物種，獲得的繁殖機(jī)會(huì)將增大變異操作進(jìn)行變異，最終找到最優(yōu)個(gè)體，對該個(gè)體解碼得到最優(yōu)解，基本步驟包括以下幾個(gè)I種群初始化與編碼第一代種群隨機(jī)獲取，組成問題的初始解，遺傳算法的編碼方式有二進(jìn)制編碼、實(shí)數(shù)編碼、格雷碼編碼、十進(jìn)制編碼等方式，其中二進(jìn)制編碼比較常見，這種編碼方式簡單、方便[11]。其中(x,y)分別用20位二進(jìn)制字符進(jìn)行編碼，h用10位二進(jìn)制字符進(jìn)行編碼。v=cell(1,m);fori=1:mtemp1=round(rand(1,50));temp2=vec2str(temp1);v{i}=temp2;endII個(gè)體評價(jià)適應(yīng)度函數(shù)是度量某個(gè)物種對于生存環(huán)境的適應(yīng)程度，對于生存環(huán)境適應(yīng)程度較高的物種，獲得的繁殖機(jī)會(huì)將增大［12］，建立適應(yīng)度函數(shù)，計(jì)算個(gè)體適應(yīng)度值。maxrec中存放最優(yōu)值。%計(jì)算種群適應(yīng)度record=adapt(v);maxrec=maxrecord(record,v);III選擇算子選擇是基于適應(yīng)度函數(shù)的一種優(yōu)勝劣汰的過程，是遺傳算法中產(chǎn)生新的子代種群的第一步，種群中個(gè)體的生存能力與適應(yīng)值的大小密切相關(guān)，染色體個(gè)體的選擇是按照個(gè)體適應(yīng)度值計(jì)算出的選擇概率值進(jìn)行隨機(jī)選取。常用的選擇方法有：輪盤賭選擇、隨機(jī)競爭選擇、確定性選擇等［13］。本文采用輪盤賭選擇方法。某染色體被選擇的概率由以下公式?jīng)Q定：Pc=f(X>)Zf(x)i1.13)工fx詁第i個(gè)個(gè)體適應(yīng)度值，functionv=chfather(v,record)岔種群中所有染色體適應(yīng)度之和。globalm;%計(jì)算適應(yīng)度和f=sum(record);%計(jì)算每個(gè)適應(yīng)度pk=record/f;qk=zeros(m,1,'double');fori=1:mqk(i)=sum(pk(1:i));End%模擬轉(zhuǎn)動(dòng)輪盤m次r=rand(m,1);fori=1:mk=1;while(r(i)>qk(k))k=k+1;endr(i)=k;end%保存新序列temp=v;fori=1:mv{i}=temp{r(i)};endendIV交叉算子種群經(jīng)過選擇后，適應(yīng)度大的個(gè)體會(huì)遺傳到下一代，交叉操作在給定的交叉概率或者自適應(yīng)交叉概率的基礎(chǔ)上可以在種群的染色體個(gè)體之間進(jìn)行基因交換［14，］這樣新產(chǎn)生的個(gè)體就能夠替代原個(gè)體，種群的規(guī)模不發(fā)生變化，交叉是遺傳算法中的“繁殖”過程，保證上一代的優(yōu)秀個(gè)體基因被繼承到下一代。在二進(jìn)制編碼中常用的交叉算子有：單點(diǎn)交叉、兩點(diǎn)交叉、均勻較差。本文采用單點(diǎn)交叉。交叉運(yùn)算的概率為pc。functionv=opcrossover(v)globalm;%單點(diǎn)交叉pc=0.25;l=1;while(l==1)r=rand(m,1);mk=find(r<pc);l=numel(mk);endif(mod(numel(mk),2)==1)mk=mk(1:end-1);endr1=randi([1,49],1,numel(mk)/2);fori=1:numel(mk)/2[v{mk(2*(i-1)+1)},v{mk(2*(i-1)+2)}]=onecross(v{mk(2*(i-1)+1)},v{mk(2*(i-1)+2)},r1(i));endendV變異算子種群在經(jīng)過選擇、交叉之后形成新的種群，為了保證種群的多樣性，需要對種群中的個(gè)體基因進(jìn)行變異操作，可以在設(shè)定的變異概率或者自適應(yīng)概率基礎(chǔ)上進(jìn)行，隨機(jī)選擇變異位置，變異操作推動(dòng)了種群的更新，保證算法最終得到的結(jié)果不停留在局部最優(yōu)解上，變異概率為pm。functionv=variation(v)globalm;%變異算子pm=0.01;fori=1:mr=rand(50,1);k=find(r<pm);forj=1:numel(k)v{i}=vari(v{i},k(j));endendendVI算法結(jié)束種群經(jīng)過多次選擇、交叉、變異后，規(guī)模不發(fā)生變化，但是種群的最優(yōu)目標(biāo)函數(shù)值聚集到一個(gè)較小范圍中，此時(shí)可以說遺傳算法收斂，得到了全局最優(yōu)解。通常用三種方式來判斷：一是將每一代種群中最優(yōu)染色體個(gè)體的目標(biāo)函數(shù)值或者適應(yīng)度值與實(shí)際優(yōu)化問題預(yù)先設(shè)定的最優(yōu)化誤差值比較，若小于這個(gè)誤差值則算法收斂到當(dāng)前代的最優(yōu)個(gè)體，否則繼續(xù)搜索；二是通過預(yù)先設(shè)定的種群遺傳最大迭代數(shù)來終止算法；三是判斷目標(biāo)函數(shù)值或者適應(yīng)度值在連續(xù)幾代種群中是否有改進(jìn)來終止算法，本文采用第二種方法。三、仿真試驗(yàn)

對于模型計(jì)算的結(jié)果，可能在x方向上誤差大，而在y方向上誤差小，因此評價(jià)一個(gè)目標(biāo)函數(shù)的好壞，不能使用在這兩個(gè)方向上的誤差，實(shí)際上在實(shí)際的偵檢行動(dòng)過程中，消防員關(guān)心的也不是在這兩個(gè)方向上的偏差，而只是在乎泄漏源是不是真實(shí)的靠近模型計(jì)算出的位置，所以用距離r的偏差來評價(jià)目標(biāo)函數(shù)更合適。1.14)1.14)r=（x-x）2+（y-y）2'實(shí)計(jì)實(shí)測利用matlab軟件來研究算法、目標(biāo)函數(shù)、布局方式對結(jié)果的影響。假設(shè)泄漏源位置為10.5，12.5，10）。三個(gè)監(jiān)測點(diǎn)位置在（40，10,0）、（50,20，0）、（70,0，0），地面有效粗1.15)糙度Z0>0.1m,大氣穩(wěn)定度為B1.15)（5=axb=0.371x0.866，（J=cxd=0.23x0.85yz1）不同的目標(biāo)函數(shù)對定位結(jié)果的影響。使用不同目標(biāo)函數(shù)，設(shè)置三個(gè)監(jiān)測點(diǎn)。結(jié)果如下：表2目標(biāo)函數(shù)對結(jié)果影響目標(biāo)函數(shù)位置x位置y高度h誤差r計(jì)算值誤差計(jì)算值誤差計(jì)算值最小平方和10.43-0.0713.43-0.939.150.851.26186370VG10.450.0513.81-1.319.130.871.57337217NRMSE10.440.0613.3-0.89.140.861.17609523RMSE10.430.0713.62-1.129.150.851.40776418從模擬結(jié)果中可以看出結(jié)果的誤差接近，任一種目標(biāo)函數(shù)優(yōu)勢并不明顯，最小平方和目標(biāo)函數(shù)由于實(shí)現(xiàn)簡單、占用系統(tǒng)計(jì)算資源小，因而本文采用了最小平方和目標(biāo)函數(shù)就能夠滿足要求。（2）算法規(guī)模對結(jié)果的影響以平方和為目標(biāo)函數(shù)，設(shè)置三個(gè)監(jiān)測點(diǎn)，平行放置，在節(jié)點(diǎn)周圍200mX200mX50范圍內(nèi)搜索。結(jié)果如下：優(yōu)化算法迭代次數(shù)n位置計(jì)算值x誤差位置y計(jì)算誤差誤差值高度h計(jì)算誤差誤差值r偏差時(shí)間t100-99.39109.8944.58-32.0821.84-11.84115.08746286.87500-60.8671.3633.34-20.8418.81-8.8174.8610132235.40遺傳算法（種群規(guī)模40）1000-42.5353.0323.97-11.4717.16-7.1654.7266607870.212000-40.0150.5123.59-11.0928.83-18.8355.03468997135.935000-20.4830.9814.5-210.14-0.1431.04480633344.14100009.940.5613.7-1.29.130.871.584455743689.72100009.940.5613.7-1.29.130.871.584455743689.72表3算法對結(jié)果影響結(jié)論：從以上論述看來，對于遺傳算法，有兩個(gè)原因可能造成結(jié)果不準(zhǔn)確甚至偏差較大：一是遺傳基因二進(jìn)制位數(shù)設(shè)置過高，導(dǎo)致交