朱建杰，周晅毅，顧 明

（同濟(jì)大學(xué)土木工程防災(zāi)國(guó)家重點(diǎn)實(shí)驗(yàn)室，上海 200092）

當(dāng)公共場(chǎng)所發(fā)生毒氣襲擊或者化工廠發(fā)生?；沸孤r(shí)，需要快速對(duì)泄漏源進(jìn)行定位和強(qiáng)度評(píng)估，從而給疏散民眾和制定應(yīng)急措施提供指導(dǎo)。由于事故發(fā)生時(shí)，流場(chǎng)和濃度場(chǎng)往往具有隨時(shí)間變化的特性，因此在時(shí)變流場(chǎng)下進(jìn)行污染源參數(shù)反演具有重要的應(yīng)用價(jià)值。

污染源參數(shù)反演的主要方法有直接求解法［1-2］、優(yōu)化方法［3-4］和概率方法［5-6］。直接求解法通過(guò)構(gòu)造控制方程的反問(wèn)題，直接求解源參數(shù)。Zhang等［1］利用Tikhonov 正則化法［7］反演了污染源強(qiáng)度，得到了污染物釋放速率的變化曲線。Wei 等［2］將Tikhonov正則化法和貝葉斯方法相結(jié)合，得到了多個(gè)污染源的污染物釋放速率和污染源位置的概率分布。優(yōu)化方法通過(guò)構(gòu)造目標(biāo)函數(shù)，求解使目標(biāo)函數(shù)最小的源參數(shù)。江思珉等［3］利用模擬退火算法反演了地下水污染源的源強(qiáng)度。曾令杰和高軍［4］通過(guò)遺傳算法對(duì)空調(diào)系統(tǒng)引發(fā)的污染進(jìn)行了快速溯源。貝葉斯概率方法考慮了測(cè)量濃度和模擬濃度之間的誤差，將誤差量化為似然概率，通過(guò)貝葉斯公式計(jì)算相應(yīng)源參數(shù)的概率分布。Keats 等［5］利用伴隨方程求解傳感器的伴隨濃度場(chǎng)，使模擬濃度的計(jì)算量大大降低，得到了非時(shí)變濃度場(chǎng)下污染源位置和強(qiáng)度的概率分布。Guo 等［6］考慮了污染物濃度場(chǎng)從零濃度狀態(tài)到穩(wěn)定狀態(tài)的變化過(guò)程，利用MCMC（Markov chain Monte Carlo）抽樣算法直接對(duì)符合后驗(yàn)概率的樣本參數(shù)進(jìn)行抽樣，獲得了污染源的位置和強(qiáng)度信息，然而該研究的污染源所在的流場(chǎng)仍然是不隨時(shí)間變化的。

本文在Keats［5］、Guo［6］的研究基礎(chǔ)上，將基于貝葉斯推斷的污染源參數(shù)反演方法擴(kuò)展至?xí)r變流場(chǎng)的場(chǎng)景?？紤]了流場(chǎng)變化對(duì)濃度場(chǎng)的影響，求解伴隨方程得到了不同污染源參數(shù)下傳感器位置處的模擬濃度，再通過(guò)貝葉斯公式計(jì)算相應(yīng)污染源參數(shù)的后驗(yàn)概率。比較了不同位置污染源的參數(shù)反演效果，還討論了利用污染物擴(kuò)散不同階段的測(cè)量濃度數(shù)據(jù)時(shí)，不同測(cè)量數(shù)據(jù)對(duì)反演效果的影響。

1 貝葉斯推斷

1.1 貝葉斯公式

貝葉斯公式可表示為

式中：P（m|d）為后驗(yàn)概率，是已知測(cè)量濃度為d條件下，污染源參數(shù)的概率分布；P（m）為先驗(yàn)概率，是未知測(cè)量信息情況下的污染源參數(shù)概率分布；P（d|m）為似然概率，表示在污染源參數(shù)為m的條件下，模擬濃度c與測(cè)量濃度d之間的相似程度。測(cè)量濃度d為傳感器測(cè)量得到的濃度，模擬濃度c通過(guò)污染物擴(kuò)散數(shù)值模擬獲得。

對(duì)于本文的污染源參數(shù)反演問(wèn)題，污染源參數(shù)為m=｛xs，ys，zs，qs｝，其中xs，ys，zs為污染源的空間坐標(biāo)，qs為污染物的釋放速率；測(cè)量濃度d=｛d1，1，d1，2，…，di，j，…，dN，M｝，其中di，j表示第i個(gè)傳感器第j個(gè)測(cè)量濃度值，M為每個(gè)傳感器記錄數(shù)據(jù)的總數(shù)，N為傳感器的數(shù)量。

1.2 先驗(yàn)概率

先驗(yàn)概率是在未得到任何測(cè)量信息情況下污染源參數(shù)的概率。這里認(rèn)為先驗(yàn)概率服從均勻分布，即：

式中：C為常數(shù)。

1.3 似然概率

似然概率表示污染源參數(shù)為m的情況下，模擬濃度c與測(cè)量濃度d之間的相似程度，通常用它們之間誤差的概率來(lái)表達(dá)。許多學(xué)者［5-6，8-13］認(rèn)為模擬濃度c與測(cè)量濃度d之間的誤差符合高斯分布，也有一部分學(xué)者［14-16］考慮到不同傳感器之間的濃度值存在量級(jí)差，認(rèn)為利用對(duì)數(shù)高斯分布描述誤差可以降低量級(jí)差對(duì)反演準(zhǔn)確性的影響?？紤]到高斯分布在污染源反演研究領(lǐng)域應(yīng)用更加廣泛，并且具有良好的反演效果，故本文采用高斯分布描述模擬濃度與測(cè)量濃度之間的誤差，似然概率為

式中：ci，j（m）表示污染源參數(shù)為m的情況下第i個(gè)傳感器第j個(gè)模擬濃度值；σi，j為模擬濃度ci，j（m）和測(cè)量濃度di，j之間誤差的標(biāo)準(zhǔn)差。假設(shè)誤差與測(cè)量濃度處于同一量級(jí)［11-13］，可以將誤差標(biāo)準(zhǔn)差σi，j設(shè)為與測(cè)量濃度di，j相同的值，即：

由于反演需要利用不同時(shí)刻多個(gè)傳感器的濃度信息，假設(shè)不同測(cè)量值之間相互獨(dú)立，可以將所有測(cè)量濃度對(duì)應(yīng)的似然概率相乘。于是，式（1）中的似然概率P（d|m）表示為

2 模擬濃度計(jì)算方法

利用貝葉斯推斷解決污染源參數(shù)反演問(wèn)題，往往要計(jì)算特定源參數(shù)條件下的模擬濃度ci，j，傳統(tǒng)的模擬濃度計(jì)算方法需要解對(duì)流擴(kuò)散方程。為了獲得完整的污染源參數(shù)后驗(yàn)概率分布，須求解所有可能污染源參數(shù)條件下的對(duì)流擴(kuò)散方程。這種方法的計(jì)算量十分龐大，反演的效率非常低。Keats 等［5］提出利用伴隨方程計(jì)算模擬濃度，使得模擬濃度計(jì)算問(wèn)題變得高效且結(jié)果精確。下面將對(duì)時(shí)變流場(chǎng)中的伴隨方程進(jìn)行推導(dǎo)。

2.1 對(duì)流擴(kuò)散方程

對(duì)流擴(kuò)散方程為

式中：Cm為源參數(shù)為m條件下的模擬濃度場(chǎng)；U為污染源所在的流場(chǎng)；K為擴(kuò)散系數(shù)；S（m）為源項(xiàng)。當(dāng)污染物釋放速率恒定且污染源為點(diǎn)源時(shí)，S（m）=qsh（x-xs，y-ys，z-zs），其中h（·）為delta函數(shù)，表示僅在坐標(biāo)為xs，ys，zs的位置存在釋放速率為qs的污染源。式（6b）和式（6c）分別為對(duì)流擴(kuò)散方程的邊界條件和初始條件。模擬濃度場(chǎng)Cm與傳感器模擬濃度ci，j（m）的關(guān)系為

式中：xi=｛xi，yi，zi｝表示第i個(gè)傳感器的空間坐標(biāo)；tj表示第j個(gè)測(cè)量數(shù)據(jù)的測(cè)量時(shí)刻。

由于污染物擴(kuò)散對(duì)流場(chǎng)的影響較小，故對(duì)流擴(kuò)散方程（6a）的流場(chǎng)U不會(huì)受到源項(xiàng)變化的影響。在這種情況下，濃度場(chǎng)Cm與污染物釋放速率qs呈線性關(guān)系；而濃度場(chǎng)Cm與污染源位置xs，ys，zs呈非線性關(guān)系。當(dāng)污染源位置發(fā)生變化時(shí)，需要求解新的對(duì)流擴(kuò)散方程，再通過(guò)式（7）計(jì)算傳感器位置處的模擬濃度。對(duì)流擴(kuò)散方程的求解是非常耗時(shí)的，由于在一次反演中污染源往往存在非常多的潛在位置，若對(duì)每一個(gè)潛在位置都求解對(duì)流擴(kuò)散方程，反演計(jì)算的效率將變得非常低。

2.2 伴隨濃度

為了快速計(jì)算任意參數(shù)條件下傳感器位置處的模擬濃度，引入了伴隨濃度C*這一變量。將伴隨濃度C*與源項(xiàng)S（m）相乘，將其乘積在時(shí)間和空間上積分得到：

利用分部積分法和高斯公式，可以將式（8）轉(zhuǎn)化為

伴隨方程（10a）可以理解成以傳感器為源項(xiàng)，在tj時(shí)刻存在脈沖的污染物釋放，同時(shí)將時(shí)間推進(jìn)方向和流場(chǎng)速度方向進(jìn)行反向處理后得到的對(duì)流擴(kuò)散方程。當(dāng)對(duì)流擴(kuò)散方程（6a）滿足邊界條件（6b）和初始條件（6c）并且伴隨方程（10a）滿足邊界條件（10b）和初始條件（10c）時(shí)，式（9）中ndSdt=0，可得到伴隨濃度場(chǎng)C*與傳感器模擬濃度ci，j（m）的關(guān)系式：

以上為利用伴隨方程計(jì)算模擬濃度的推導(dǎo)過(guò)程。由伴隨方程（10a）可知，伴隨濃度場(chǎng)C*i，j跟源參數(shù)m無(wú)關(guān)，它僅跟第i個(gè)傳感器的坐標(biāo)和測(cè)量時(shí)刻tj相關(guān)。當(dāng)傳感器位置xi和測(cè)量時(shí)刻tj不變時(shí)，計(jì)算不同源參數(shù)條件下的模擬濃度只需將相應(yīng)源項(xiàng)S（m）代入式（11）即可，式（11）只涉及積分運(yùn)算，而利用對(duì)流擴(kuò)散方程時(shí)需要求解偏微分方程，前者的計(jì)算量遠(yuǎn)遠(yuǎn)小于后者，這樣就可以極大地加快反演效率。設(shè)每個(gè)傳感器測(cè)量的次數(shù)為M，傳感器的數(shù)量為N，則所需求解伴隨方程的總數(shù)為M×N。

3 算例分析

3.1 污染物擴(kuò)散數(shù)值模擬

3.1.1 模擬工況

如圖1 所示，建筑模型尺寸為0.1 m×0.1 m×0.2 m（L·B·H），其周?chē)嬖贏、B、C、D、E 5 個(gè)不同位置的污染源，均位于地面，用空心方塊表示，對(duì)5個(gè)污染源分別做污染物擴(kuò)散數(shù)值模擬。來(lái)流剖面遵循指數(shù)律，建筑頂部高度處風(fēng)速為uH=4.2 m·s-1，指數(shù)α=0.25。10個(gè)傳感器位于建筑背風(fēng)側(cè)，用實(shí)心黑點(diǎn)表示，編號(hào)為1～10，位于距地面0.062 5H高度處。

圖1 污染源與傳感器布置方式Fig.1 Layout of sources and sensors

3.1.2 計(jì)算域與模擬設(shè)置

根據(jù)日本建筑學(xué)會(huì)的計(jì)算流體動(dòng)力學(xué)（computational fluid dynamics，CFD）指南［17］，計(jì)算域（圖2）尺寸設(shè)為3.4 m×2.0 m×0.8 m（17H·10H·4H）；計(jì)算域入口與建筑模型的距離為4.5H；計(jì)算域的側(cè)面和出口與建筑模型的距離分別為4.75H和12H。

圖2 計(jì)算域Fig.2 Computational domain

流場(chǎng)模擬以及污染物擴(kuò)散模擬參數(shù)設(shè)置如表1所示，湍流模型采用剪切應(yīng)力輸運(yùn)（shear stress transfer，SST）k-ω；流場(chǎng)的計(jì)算時(shí)間為50 s（1 050H/uH）。當(dāng)流場(chǎng)已經(jīng)達(dá)到穩(wěn)定（t=630H/uH）時(shí)，污染物開(kāi)始釋放。由于該模擬工況與日本東京工藝大學(xué)于2006年做的污染物擴(kuò)散風(fēng)洞試驗(yàn)工況類(lèi)似，故來(lái)流平均速度u和湍動(dòng)能k根據(jù)試驗(yàn)數(shù)據(jù)確定，釋放速率為qs=5.83×10-6m3·s-1。

表1 邊界條件與數(shù)值模擬設(shè)置Tab.1 Boundary conditions and solver settings for simulation case

3.1.3 網(wǎng)格無(wú)關(guān)性檢驗(yàn)

為了確保模擬結(jié)果不受網(wǎng)格劃分方案影響，需要進(jìn)行網(wǎng)格無(wú)關(guān)性檢驗(yàn)（表2）。計(jì)算域采用結(jié)構(gòu)化網(wǎng)格劃分，建筑模型周?chē)约皞鞲衅鞲浇鼌^(qū)域進(jìn)行了網(wǎng)格加密，粗、中、細(xì)三種網(wǎng)格的最小單元尺寸分別為H/20、H/30和H/40，網(wǎng)格的伸縮比設(shè)置為1.05，總網(wǎng)格數(shù)分別為571 400、1 254 123和2 105 408。

表2 網(wǎng)格無(wú)關(guān)性檢驗(yàn)參數(shù)設(shè)置Tab.2 Parameter settings of grid independent analysis test

圖3 為污染源C 處釋放的污染氣體在建筑模型下風(fēng)向0.5H和H位置且y=0 處的模擬濃度剖面，可見(jiàn)僅在建筑下風(fēng)向0.5H且高度大于1.2H的位置處，粗網(wǎng)格條件下模擬得到的氣體濃度小于中網(wǎng)格和細(xì)網(wǎng)格條件下的模擬結(jié)果，其他位置處三種網(wǎng)格的模擬濃度曲線重合度較高，故本文選取中網(wǎng)格對(duì)污染物擴(kuò)散進(jìn)行數(shù)值模擬。

圖3 不同網(wǎng)格劃分方案的平均濃度剖面Fig.3 Profile of mean concentration for different meshing schemes

3.1.4 自保持性驗(yàn)證

為了保證計(jì)算域入口至建筑模型處的來(lái)流剖面不會(huì)發(fā)生顯著變化，需要進(jìn)行自保持性驗(yàn)證。圖4a為空風(fēng)場(chǎng)條件下，數(shù)值模擬和風(fēng)洞試驗(yàn)中的x向平均風(fēng)速剖面從計(jì)算域入口至建筑位置處的變化情況。在數(shù)值模擬結(jié)果中，x向平均風(fēng)速?gòu)挠?jì)算域入口至建筑位置處的變化較小，與風(fēng)洞試驗(yàn)的測(cè)量結(jié)果也較為吻合。圖4b為空風(fēng)場(chǎng)條件下，數(shù)值模擬和風(fēng)洞試驗(yàn)中的湍動(dòng)能剖面從計(jì)算域入口至建筑位置處的變化情況。數(shù)值模擬和風(fēng)洞試驗(yàn)得到的湍動(dòng)能均出現(xiàn)了一定的衰減，相比風(fēng)洞試驗(yàn)結(jié)果，數(shù)值模擬的湍動(dòng)能衰減程度更大。Kim和Baik［18］研究了數(shù)值模擬中湍動(dòng)能對(duì)污染物濃度的影響，根據(jù)研究結(jié)果可知本文數(shù)值模擬中的湍動(dòng)能衰減對(duì)污染物濃度場(chǎng)的影響不大。

圖4 空風(fēng)場(chǎng)條件下的x向平均風(fēng)速與湍動(dòng)能剖面（y/H=0）Fig.4 Profile of mean velocity in x direction and turbulent kinetic energy in empty wind field (y/H=0)

3.1.5 流場(chǎng)

圖5 為3 號(hào)傳感器位置處的x向的風(fēng)速（u）時(shí)程，流場(chǎng)穩(wěn)定后的風(fēng)速存在有規(guī)律的周期性波動(dòng)。

圖5 3號(hào)傳感器處x向瞬時(shí)風(fēng)速時(shí)程Fig.5 Temporal velocity in x direction for Sensor 3

3.1.6 對(duì)流擴(kuò)散方程計(jì)算測(cè)量濃度

當(dāng)流場(chǎng)已達(dá)穩(wěn)定狀態(tài)（t=630H/uH）時(shí)，釋放污染物氣體。圖6 給出了分別由污染源A～E 釋放的污染氣體，擴(kuò)散至3 號(hào)傳感器位置處隨時(shí)間變化的濃度，圖中的濃度值經(jīng)過(guò)了量綱一處理，其中c0=qs/（uH·H2）=34.7×10-6。從污染物釋放初期（t=630H/uH）開(kāi)始，濃度由零逐漸發(fā)展至穩(wěn)定，受到周期脈動(dòng)的流場(chǎng)影響，污染物的濃度在發(fā)展過(guò)程中和穩(wěn)定后均表現(xiàn)出有規(guī)律的周期波動(dòng)特性。這里通過(guò)對(duì)流擴(kuò)散方程（6a）對(duì)濃度場(chǎng)進(jìn)行求解，將式（7）的計(jì)算結(jié)果作為傳感器的測(cè)量濃度d。

圖6 3號(hào)傳感器濃度時(shí)程Fig.6 Concentration history for Sensor 3

為了在后續(xù)研究中探究利用污染物擴(kuò)散不同階段的測(cè)量數(shù)據(jù)反演源參數(shù)時(shí)不同測(cè)量數(shù)據(jù)對(duì)反演結(jié)果的影響，以t=660H/uH為分界將污染物擴(kuò)散過(guò)程劃分成兩個(gè)階段。在t=660H/uH之前為發(fā)展階段，在t=660H/uH之后為穩(wěn)定階段。由于3號(hào)傳感器的濃度時(shí)程最具有代表性，其他傳感器的濃度變化與之類(lèi)似，故不再畫(huà)出其他傳感器的濃度時(shí)程圖。

3.1.7 伴隨方程計(jì)算模擬濃度

2.2 節(jié)介紹了在計(jì)算任意源參數(shù)條件下某一傳感器位置處某一時(shí)刻的模擬濃度時(shí)，基于伴隨方程的計(jì)算方法在節(jié)省計(jì)算量方面具有優(yōu)勢(shì)，故本文通過(guò)求解伴隨方程（10a），再利用式（11）計(jì)算模擬濃度。圖7 對(duì)比了污染源分別為A～E 時(shí)，3 號(hào)傳感器位置處的測(cè)量濃度和模擬濃度的計(jì)算結(jié)果。當(dāng)污染源為A時(shí)，模擬濃度與測(cè)量濃度相比偏??；當(dāng)污染源為B時(shí)，模擬濃度與測(cè)量濃度相比偏大；對(duì)于其他位置的污染源，兩種濃度的計(jì)算結(jié)果基本保持一致。雖然污染源分別為A 和B 時(shí)，測(cè)量濃度和模擬濃度在3 號(hào)傳感器位置的結(jié)算結(jié)果存在差異，但是兩者均表現(xiàn)出了一致的周期性波動(dòng)趨勢(shì)，數(shù)值上的誤差可能是由網(wǎng)格劃分、微分方程離散和四舍五入等因素導(dǎo)致的，這表明基于伴隨方程的模擬濃度方法可以在忽略數(shù)值誤差的情況下達(dá)到與基于對(duì)流擴(kuò)散方程時(shí)相同的效果。在反演的過(guò)程中，這里的誤差可以作為數(shù)值模擬與實(shí)際測(cè)量的差異。值得注意的是，在基于對(duì)流擴(kuò)散方程計(jì)算5 個(gè)不同位置污染源產(chǎn)生的污染物擴(kuò)散至3號(hào)傳感器位置處在某一時(shí)刻的模擬濃度時(shí)，須求解5 種源參數(shù)條件下的對(duì)流擴(kuò)散方程；而基于伴隨方程的方法只需針對(duì)3 號(hào)傳感器在同一測(cè)量時(shí)刻求解一次，即可達(dá)到相同效果，大大節(jié)省了計(jì)算量。

圖7 3號(hào)傳感器位置處測(cè)量濃度與模擬濃度比較Fig.7 Comparison of measurements and simulated concentrations at location of Sensor 3

3.2 污染源參數(shù)的后驗(yàn)概率

假設(shè)污染源位于地面，即污染源的豎向坐標(biāo)zs=0，因此僅對(duì)污染源的xs、ys坐標(biāo)和污染源的污染物釋放速率qs進(jìn)行反演。

3.2.1 污染源參數(shù)的傳感器與污染源相對(duì)位置的影響

圖8、圖9 給出了當(dāng)污染源分別位于A～E 時(shí)通過(guò)貝葉斯推斷得到的污染源位置和強(qiáng)度的后驗(yàn)概率。在污染源位置的后驗(yàn)概率云圖中，星號(hào)處為污染源的真實(shí)位置，空心圓為傳感器；在污染源強(qiáng)度的后驗(yàn)概率柱狀圖中，虛線表示污染源的真實(shí)強(qiáng)度。污染源A位置的后驗(yàn)概率分布在x方向呈較寬的分布，后驗(yàn)概率均值為-1.29H，與真實(shí)值之間的誤差為0.54H，后驗(yàn)概率極值點(diǎn)為-0.75H，與真實(shí)值一致，但其后驗(yàn)概率標(biāo)準(zhǔn)差高達(dá)0.40H，表明污染源位置x坐標(biāo)的反演結(jié)果具有較大的不確定性；污染源A的強(qiáng)度后驗(yàn)概率均值和極值點(diǎn)分別為0.72qs和0.62qs，反演結(jié)果與真實(shí)值相比偏小。雖然污染源B位置的反演誤差和不確定性與污染源A相比有所減少，但是在其強(qiáng)度反演結(jié)果中，后驗(yàn)概率分布較為平坦，后驗(yàn)概率均值為1.87qs，與真實(shí)值之間的誤差高達(dá)87%，后驗(yàn)概率標(biāo)準(zhǔn)差為0.80qs，表明污染源B強(qiáng)度的預(yù)測(cè)可靠性較低。

圖8 污染源A～E的位置概率云圖Fig.8 Probability contour of location for pollution source A-E

圖9 污染源A～E的強(qiáng)度概率密度Fig.9 Probability density of strength for pollution source A-E

污染源A 和污染源B 的反演誤差較大，是因?yàn)橐訟 和B 為污染源做污染物擴(kuò)散數(shù)值模擬時(shí)，傳感器位置處的測(cè)量濃度和模擬濃度之間存在較大的誤差，所以當(dāng)數(shù)值模擬不能準(zhǔn)確反映真實(shí)的測(cè)量結(jié)果時(shí)，會(huì)導(dǎo)致反演出現(xiàn)較大的偏差。同時(shí)，在反演結(jié)果中可以發(fā)現(xiàn)，污染源A 和污染源B 的源參數(shù)后驗(yàn)概率標(biāo)準(zhǔn)差也大于其他位置的污染源。這是因?yàn)檫@兩個(gè)污染源與傳感器之間的距離較遠(yuǎn)。由似然概率式（3）可知，當(dāng)污染源參數(shù)xs，、ys和qs能正確描述真實(shí)污染源時(shí)，似然概率達(dá)到最大值，即概率極值。在計(jì)算距傳感器較遠(yuǎn)的污染源所產(chǎn)生的污染物在傳感器位置的模擬濃度c（xs，ys，qs）時(shí)，模擬濃度值對(duì)于源參數(shù)的變化不敏感，即和較小，后驗(yàn)概率在概率極值點(diǎn)附近的變化不明顯，導(dǎo)致后驗(yàn)概率標(biāo)準(zhǔn)差較大，因此污染源A 和污染源B 的反演結(jié)果具有較大的不確定性。為了減少反演結(jié)果的不確定性，可以根據(jù)歷史信息（曾經(jīng)發(fā)生過(guò)泄露的位置）和經(jīng)驗(yàn)設(shè)計(jì)傳感器布置方式。

在計(jì)算污染源C～E 產(chǎn)生的污染物擴(kuò)散至傳感器位置的濃度值時(shí)，測(cè)量濃度和模擬濃度的計(jì)算結(jié)果吻合較好，因此源參數(shù)反演的誤差較小。同時(shí)，傳感器的濃度值對(duì)于污染源C、污染源D 和污染源E的參數(shù)變化靈敏度較高，導(dǎo)致了源參數(shù)的后驗(yàn)概率標(biāo)準(zhǔn)差較低，反演的不確定性較小。

3.2.2 利用不同擴(kuò)散階段測(cè)量數(shù)據(jù)對(duì)反演結(jié)果的影響

如圖6 所示，以t=660H/uH為分界線將污染物擴(kuò)散過(guò)程分為兩個(gè)階段。第一個(gè)階段為發(fā)展階段，包含了污染物濃度從零濃度達(dá)到穩(wěn)定狀態(tài)的發(fā)展過(guò)程；第二個(gè)階段為穩(wěn)定階段，完全由周期性波動(dòng)的測(cè)量數(shù)據(jù)構(gòu)成。下面分別利用污染源在C處時(shí)傳感器位置處發(fā)展階段、穩(wěn)定階段和所有時(shí)刻的測(cè)量數(shù)據(jù)對(duì)污染源參數(shù)進(jìn)行反演，反演結(jié)果分別如圖10、圖11所示。相比利用穩(wěn)定階段的測(cè)量數(shù)據(jù)的污染源位置反演結(jié)果，利用發(fā)展階段的測(cè)量數(shù)據(jù)可以降低位置誤差和減小后驗(yàn)概率標(biāo)準(zhǔn)差，但穩(wěn)定階段數(shù)據(jù)提升的效果并不顯著。圖表中的結(jié)果也反映出利用所有時(shí)刻的測(cè)量數(shù)據(jù)時(shí)，其反演結(jié)果的精確程度提高非常有限。

圖10 利用不同階段測(cè)量數(shù)據(jù)的污染源A～E的位置概率云圖Fig.10 Probability contour of location for pollution source A-E by using measurement from different stages

圖11 利用穩(wěn)定階段測(cè)量數(shù)據(jù)的污染源A～E的強(qiáng)度概率密度Fig.11 Probability density of strength for pollution source A-E by using measurement from different stages

4 結(jié)論

本文基于貝葉斯推斷對(duì)時(shí)變流場(chǎng)下的污染源參數(shù)進(jìn)行了反演，比較了不同位置污染源的參數(shù)反演結(jié)果，利用不同擴(kuò)散階段的測(cè)量數(shù)據(jù)探討了其對(duì)反演結(jié)果的影響，得到了以下結(jié)論：

（1）模擬濃度與測(cè)量濃度之間的誤差是影響反演精確性的主要因素，當(dāng)數(shù)值模擬能準(zhǔn)確地反映污染物擴(kuò)散時(shí)，污染源反演的準(zhǔn)確性能夠得到極大的提高。

（2）當(dāng)污染源與傳感器距離較近時(shí)，污染源參數(shù)的后驗(yàn)概率分布比較集中，反演結(jié)果的不確定性較小，真實(shí)的污染源參數(shù)可能出現(xiàn)在較窄的區(qū)間，污染源所需搜索的范圍較??；當(dāng)污染源與傳感器距離較遠(yuǎn)時(shí)，污染源參數(shù)的后驗(yàn)概率分布較寬，反演結(jié)果的不確定性較大，真實(shí)的污染源參數(shù)可能出現(xiàn)在較寬的區(qū)間，需要較大范圍搜索污染源。

（3）不同階段的污染物擴(kuò)散測(cè)量數(shù)據(jù)對(duì)反演結(jié)果不會(huì)造成顯著的影響。

本研究的不足之處在于污染源參數(shù)反演是在已知流場(chǎng)的前提下進(jìn)行的，流場(chǎng)的未知性會(huì)影響模擬濃度的計(jì)算，從而可能造成反演不準(zhǔn)確。如何減少流場(chǎng)未知性對(duì)反演的影響，需要進(jìn)一步深入研究。