基于K-均值算法的煙氣層高度判定方法

李 昂，侯 岳，任 凱，王康勃

(海軍工程大學動力工程學院，湖北 武漢 430033)

煙氣層高度值是表征室內火災煙氣特性的一個重要參數，是人員疏散時間、建筑防火安全性評估、防排煙設計過程中的關鍵計算依據[1]，在火災煙氣試驗中，基本都要對煙氣層高度予以判定和計算。對煙氣層高度的判定，一般是將溫度、CO體積分數、CO2體積分數[2]中的一種作為分層表征參數，許多相關學者通過研究形成了一些煙氣層高度判定方法，如：Cooper等[3]提出的N-百分比法；Quintiere等[4]提出的基于上層區域平均和質量等效的方法；Janssens等[5]提出的基于最大梯度的判定方法；He等[6]提出的積分比值法；毛少華等[7]針對煙氣層溫度較低、溫度數據梯度較小的問題，提出了采用方差方法分析每層內部溫度均一性和層間溫度的最大差異性，在此基礎上構建一個溫度閾值與方差之間的關系式，最小溫度閾值對應的高度即是煙氣層高度。

目前應用較多的煙氣層高度判定方法是積分比值法。但應用積分比值法確定煙氣層高度存在溫度數據處理較為繁瑣的問題，由于受算法的限制，該法僅能分析對象某一垂直方向上的溫度數據，而實際計算過程中，假如要想得知整個空間內的煙氣層高度變化情況，基本需要采用兩種方式：一是沿煙氣蔓延的方向一一處理每個垂直方向上的溫度數據；二是先將一定區間內的溫度數據進行預處理，將若干個熱電偶數據合并為一個數據，但煙氣層高度在縱向上的變化受到通風、空間結構、火源功率等因素的影響，很難準確判斷合理的區間劃分位置和數量，而且一般來說由于煙氣具有卷吸效應[8]，越向下游，煙氣層的高度會變得越低，空間內整個煙氣層高度較為一致的情況很難出現，因此無論采取何種溫度數據處理方式，當煙氣溫度數據集規模較大時，采用積分比值法確定煙氣層高度均存在工作量較大且過程復雜的問題。

數據挖掘理論中的K-均值算法[9]是利用一定的算法步驟自動抽取數據特征，并將不同數據進行分類，被廣泛應用于數據分類問題。因此，本文將K-均值算法用于判定煙氣層高度，并通過火災煙氣試驗，對比分析了積分比值法、視覺法和K-均值算法3種煙氣層高度判定方法所得到的結果，以驗證K-均值算法的有效性。

1 K-均值算法用于煙氣層高度判定計算模型

K-均值算法已被廣泛應用于用戶特征分類、故障診斷、模式分析等人工智能領域[10]，是較為成熟的數據挖掘計算方法，具有模型數學思想較為簡單、無需人工標記數據特征、自動進行類別劃分的優點，因此非常適用于數據分類問題。其基本思想是通過初始化類別中心點，對樣本數據與中心的相似性進行測度，對類間進行測度，從而自動將數據進行分類。

1.1 K-均值算法的基本步驟

K-均值算法是一種空間數據劃分或分組的重要方法[11]，該方法主要是將研究對象的空間距離指標按照相似性準則劃分到若干個子集中，使得相同子集中各元素間差別最小，而不同子集中各元素間差別最大。最常用的是歐幾里得距離，其表達式為

(1)

K-均值算法的計算過程如下：

(1) 由用戶確定所要聚類的準確數目K，并隨機選擇K個對象(樣本)，每個對象稱為一個種子，代表一個簇(類)的均值或中心。

(2) 對剩余的每個對象，根據其與各簇中心的距離將它賦給最近的簇。

(3) 重新計算每個簇內對象的平均值，形成新的聚類中心。

(4) 這個過程重復進行，直到下式準則函數收斂為止：

(2)

式中：E為所有研究對象的平方誤差總和；p為空間的點，即數據對象；mi為簇Ci的平均值。

按照這個準則生成的結果簇趨向于獨立和緊湊，見圖1。

圖1 空間聚類過程Fig.1 Spatial clustering process

1. 2 簇數和聚類中心的確定

根據K-均值算法的基本計算思路，進行煙氣層高度的判定主要涉及到兩個核心問題：一是設定最佳的類別數量(K值)，數據挖掘中稱之為總的簇數，確定了簇數就確定了煙氣層的縱向劃分區間；二是設定合理的初始聚類中心。這是兩個影響到后續數據處理和計算量的主要問題。

在經典算法中初始聚類中心是隨機選擇的，這樣的方式有可能會使聚類結果陷入局部最優，一般來說對于這一問題的處理有多種方法，如CCIA算法、模擬退火算法、遺傳算法等[13]，但這些算法的計算過程過于復雜，存在計算量較大的問題。

如何選擇較為合理的初始聚類中心從而避免聚類結果陷入局部最優，其基本思路是先設定一個大致合理的初始聚類中心，之后根據算法再確定準確的初始聚類中心。火災煙氣由于卷吸效應的存在，煙氣層與空氣層之間存在一個“過渡層”，過渡層溫度大致恰好位于高低溫中間狀態，確定了過渡層位置后，根據文獻[14]，由于各層內溫度大致呈現由下向上增加的趨勢，因此可以合理假定其聚類中心為幾何中心,這樣設定的聚類中心就大致吻合了實際的溫度分布規律，從而避免了設定的初始聚類中心不準確導致聚類結果陷入局部最優的問題。

圖2 煙氣層聚類中心選定Fig.2 Selection of clustering center of smoker layer

如圖2所示，假設已知p1,p2,…,pn點的數據，高度為z1,z2,…,zn，溫度為T1,T2,…,Tn，縱向坐標為x1,x2,…,xn，確定煙氣層和空氣層初始聚類中心的基本步驟如下：

(1) 劃分區間。在縱向上劃分區間，根據每一步所設定的聚類數K值，設定劃分的區間數Rn=K/2，在縱向上將熱電偶所占據的空間均等劃分為Rn個。

(2) 數據點分類。根據所劃分的區間空間縱向坐標，將熱電偶數據按區間劃分，設區間的分界面縱向坐標為xRi,Ri+1，區間內任何一點pi的縱向坐標xi≤xRi,Ri+1。

(3) 估算過渡層位置。設某一區間內數據點共有n個，因此過渡層豎向坐標和縱向坐標為

(3)

(4) 確定煙氣層和空氣層的初始聚類中心。設煙氣層內的數據點為m個，空氣層內的數據點為l個，因此兩者的初始聚類中心由下式計算可得：

(4)

式中：zsj、xsj、Tsj分別為煙氣層豎向坐標、縱向坐標和溫度值；zaj、xaj、Taj分別為空氣層豎向坐標、縱向坐標和溫度值。

1.3 K-均值算法確定煙氣層高度的計算步驟

按照前述方法可以設定初始的聚類中心和聚類數，某一煙氣層的溫度數據應當是在一定區域內聚集，意為針對一個類別而言，除了應使溫度數據間的差異最小外，還應使一個類別的數據總體空間幾何距離的差異最小，因此樣本的向量空間應當包括縱向坐標值、豎向坐標值、溫度值3個向量。

按照K-均值算法的步驟，得到各個簇Ci的分類結果后，需要得到煙氣層的高度值，即求解出煙氣層與空氣層分界面的坐標，具體的計算步驟如下：

(1) 確定樣本集合。設樣本集合內包含m個樣本，即D={p1,p2,…,pm}，其中pi(i=1,2,…,m)代表每一個熱電偶采集的數據，每個樣本有3個向量，即pi=(xi;zi;Ti),其中，xi為第i個熱電偶的縱向坐標;zi為第i個熱電偶的豎向高度坐標;Ti為第i個熱電偶的溫度數據。

(3) 設定初始聚類中心。初始聚類中心共有K個,即(μ1,μ2,…,μK)，初始聚類中心的選定按照第1.2節中所述方法。

(4) 計算樣本與聚類中心的距離。按照下式計算每一個樣本pi(i=1,2,…,m)與每一個聚類中心μp(1≤p≤k)的距離：

(5)

(5) 進行樣本歸集。根據歐式距離最小值，將相應的樣本歸入相應的簇，即Ci=Ci∪{pi}。

(6)

(7) 評價聚類結果最優性。對于第(2)步中給出的每一個K值，在聚類算法中，一般采用Silhouette指標評價哪個聚類結果是最優的[14-15]，該指標反映了聚類結構的類內緊密性和類間分離度，按照下式計算：

(7)

式中：a(i)為樣本i與類內所有其他樣本的平均距離;b(i)為樣本i到其他每個類中樣本平均距離的最小值。

得到了最佳聚類數kopt，即得到了所劃分的煙氣層縱向區間數，同時根據最佳聚類數，可以得到最優聚類中心的集合：

(8)

式中：Cμp為對應kopt的最優聚類中心集合。

(9)

將每一個簇的計算結果進行比較，如果xi=xj、Ti>Tj,則簇Ci為

{(p1,y1),(p2,y2),…,(pk,yk)},yi∈{1}

(10)

反之，則簇Ci為

{(p1,y1),(p2,y2),…,(pk,yk)},yi∈{-1}

(11)

式中：yi為樣本標記，1代表煙氣層，-1代表空氣層。

(9) 確定煙氣層高度值。數據樣本歸類到相應的簇內，從縱向上就劃分了區間，對煙氣層高度的求解采用一個煙氣層區間的平均高度來表征高度，這樣可以突出表征煙氣層的區間特征，此時只需要找到煙氣層聚類中心與樣本形成一定夾角范圍內的邊緣樣本[16]，求取其對應高度的平均值即可，其計算公式如下：

(12)

式中:Hs為某一區間煙氣層平均高度(m)；pci為簇Ci的某一個樣本;zq為pci的豎向坐標值(m);z1,z2,…,zr為符合條件的各個樣本的豎向坐標值(m)；zj為得到的某一位置的豎向坐標最小值(m)；xμ′和zμ′分別為簇Ci聚類中心的縱向和豎向坐標值(m)；θ1,θ2為簇Ci邊緣分布的特征角(°)。

2 K-均值算法用于確定煙氣層高度的火災煙氣試驗

火災煙氣試驗在海軍工程大學艦船火災模擬試驗艙進行，試驗艙共有3個部分，分別為模擬通道、狹小艙室和典型艙室，其中本次試驗選取在典型艙室進行，其尺寸為5 m×5 m×2.5 m(長×寬×高)，整體為鋼制框架結構，試驗艙室外觀如圖3所示。

圖3 試驗艙室外觀Fig.3 Appearance of the test chamber

2. 1 試驗測量儀器布置和工況設置

2.1.1 試驗測量儀器的布置

試驗火源采用柴油作為燃燒材料，容器為水浴鋼制油盤，形狀為方形，尺寸為0.6 m×0.6 m×0.2 m(長×寬×深)，油盤內部采用鋼制格柵進行分割，用于改變燃燒面積；試驗布置了熱電偶樹，從1號門至2號門每間隔0.5 m布置一個熱電偶樹，其布置了8個熱電偶樹，每個熱電偶樹從底部至頂部每間隔0.2 m布置一個熱電偶，共布置了108個熱電偶；此外，試驗艙室布置激光片源提供縱剖面光場，利于觀察煙氣層高度在整個試驗過程中的變化情況，同時在光場中布置激光片源和標尺，通過查看標尺的視覺方法確定煙氣層高度。試驗測量儀器的布置，見圖4。

圖4 試驗測量儀器布置示意圖Fig.4 Schematic diagram of test measuring instruments

2.1.2 試驗工況的設置

對于試驗工況的設定，根據美軍所完成的艦艇火災煙氣特性試驗表明[17]，火源面積位于艙室空間尺寸(除設備占用面積的有效尺寸)的3/50～7/50區間內時，試驗結果基本能夠表征典型火災時煙氣的運動特征，火源面積占艙室空間尺寸的比例過小，則煙氣產生速率不足，火源面積占艙室空間尺寸的比例過大時，其試驗結果與典型火源面積比例區間內的試驗結果并無本質區別。因此，為了保證試驗的典型性，將火源面積占艙室空間尺寸的比例設定為3/50和7/50兩種比例，即火源面積設定為0.16 m2和0.36 m2。另外，油盤底部安裝了4個液壓傳感器，通過監測壓力變化換算為液面高度變化，從而計算油盤燃料的燃燒速率，大致確定火源不同的燃燒階段，并選取不同燃燒階段的典型時間點，進行后續數據分析。試驗過程中燃料的燃燒速率按照下式計算：

為了保證試驗工況的典型性，設置了兩種工況，兩種工況的油盤面積(即火源面積)和通風模式都不相同，油盤置于試驗艙中心位置，試驗工況見表1。

表1 試驗工況

2. 2 試驗結果與討論

按照前述的K-均值算法判定煙氣層高度的計算步驟，采用MATLAB軟件編制數據處理程序，用于處理溫度數據，并針對兩種試驗工況，分別采用K-均值算法、積分比值法、視覺標尺法計算煙氣層高度，通過對比分析來驗證K-均值算法的有效性。

2.2.1 C1工況試驗結果及討論

圖5和圖6顯示了C1工況時燃料的燃燒速率和距離火源中心0.5 m處不同典型高度處煙氣層溫度的分布情況。

圖5 C1工況時燃料的燃燒速率Fig.5 Fuel burning rate under working condition C1

圖6 C1工況時距離火源中心0.5 m處不同典型 高度處煙氣溫度分布圖Fig.6 Smoke layer temperature distribution at different typical heights 0.5 m away from the fire source center under working condition C1

由圖5和圖6可以看出：在0～120 s區間內燃料的燃燒速率上升較快，在120～200 s區間內燃料的燃燒較為穩定，200 s后燃料的燃燒進入衰減階段，故選取不同燃燒階段的時間點，即選定50 s、170 s、260 s三個時間點作為溫度數據處理的時間點。

利用K-均值算法分別對50 s、170 s、260 s時間點的溫度數據進行了迭代，圖7和圖8分別顯示C1工況時通過K-均值算法、積分比值法計算所得到的煙氣層高度值。

圖7 C1工況時基于K-均值算法的煙氣層高度Fig.7 Smoke layer heights obtained by K-means algorithm under working condition C1

由圖7可以看出：K-均值算法將50 s時煙氣層高度劃分為了6個區間，距離火源最遠的位置，煙氣層高度為2.38 m，火源中心煙氣層高度最低，其值為0.33 m；將170 s時煙氣層高度劃分為了6個區間，其最高點為2.2 m，最低點為0.4 m；將260 s時煙氣層高度劃分為了3個區間，其最高點為1.5 m，最低點為1.3 m。

圖8顯示了積分比值法得到的煙氣層高度，積分比值法是按照不同的豎向高度來計算其煙氣層高度的。

圖8 C1工況時基于積分比值法的煙氣層高度Fig.8 Smoke layer heights obtained by internal ratio method under working condition C1

由圖7和圖8可知：在50 s、170 s時間點，由K-均值算法計算得到的煙氣層高度值在距離火源較遠的位置要高于積分比值法的計算結果，而在火源中心位置兩者計算得到的煙氣層高度較為接近，原因在于K-均值算法計算得到的煙氣層高度是一個區間內的高度，而積分比值法本質上是最小二乘法，其計算得到的煙氣層高度是單個縱坐標點所對應的煙氣層高度，由于算法的差異性，在遠離火源的位置，由于卷吸效應較為明顯，煙氣層溫度數據梯度較大，因此兩者的計算結果有著明顯的差異；通過比較兩種算法260 s時煙氣層高度的計算結果，可以明顯看出，無論是距離火源較遠的位置還是火源中心位置附近，兩者計算結果的差異并不明顯，主要原因在于此時煙氣層溫度較為均一，因此兩者計算結果的差異不大。

2.2.2 C2工況試驗結果及討論

圖9和圖10顯示了C2工況時燃料的燃燒速率和距離火源中心0.5 m處不同典型高度處煙氣層溫度的分布情況。

圖9 C2工況時燃料的燃燒速率Fig.9 Fuel burning rate under working condition C2

圖10 C2工況時距離火源中心0.5 m不同典型 高度處煙氣層溫度分布圖Fig.10 Smoke layer temperature distribution at different typical heights 0.5 m away from the fire source center under working condition C2

由圖9和圖10可知：在0～140 s區間內燃料的燃燒速率上升較快，在140～200 s區間內燃料的燃燒大體保持穩定，在240～340 s區間內燃料的燃燒進入衰減階段，故選取典型燃燒階段的時間點，即選定100 s、220 s、330 s三個時間點作為溫度數據處理的時間點。

利用K-均值算法和積分比值法，對C2工況時煙氣層高度進行了計算，其計算結果見圖11和圖12。

圖11 C2工況時基于K-均值算法的煙氣層高度Fig.11 Smoke layer heights obtained by K-means algorithm under working condition C2

由圖11可以看出：K-均值法將100 s時煙氣層高度劃分為了7個區間，其最高點為2.05 m，最低點為0.31 m；將220 s時煙氣層高度劃分為了6個區間，其最高點為1.6 m，最低點為0.31 m；將330 s時煙氣層高度劃分為了5個區間，其最高點為1.13 m，最低點為0.32 m。

圖12 C2工況時基于積分比值法的煙氣層高度Fig.12 Smoke layer heights obtained by integral ratio method under working condition C2

通過比較兩種算法的計算結果(見圖11和圖12)可知：100 s時積分比值法計算得到的煙氣層高度在側壁附近接近K-均值算法的計算結果，而在火源與側壁的中間位置附近，K-均值算法的計算結果顯示出明顯的區間特征，且出現“振蕩”現象；由于單側通風存在，最低位置并不位于火源中心，而是在中心靠后的空間區域，因此220 s時同樣出現了“振蕩”現象，兩者的計算結果在煙氣層高度值最低點較為接近；均值法計算結果顯示出，330 s時煙氣層高度仍然有較為明顯的區間特征。

利用視覺標尺法觀察煙氣層高度，可以進一步確定煙氣層高度，表2顯示了不同工況、不同燃燒時間時，通過3種方法計算得到的煙氣層高度值經數學平均處理后的結果。

表2 不同煙氣層高度判定方法計算結果的比較

3 結 論

通過火災煙氣試驗，將K-均值算法與積分比值法進行對比，驗證了K-均值法判定煙氣層高度的計算步驟，得到如下結論：

(1) K-均值算法與積分比值法相比較，兩種方法對煙氣層高度值的計算結果在各個位置點大體上一致，但是具體數值上存在差異，尤其是煙氣層溫度數據出現較大梯度時，兩者的差異較為明顯，但總體趨勢和數值范圍兩者基本保持一致，故以K-均值算法為基礎構建的煙氣層高度判定數學模型是有效的。

(2) 通過計算結果表明：K-均值算法可以自動將煙氣層進行縱向區間的識別和劃分，從而得到一個區間的煙氣層高度，而積分比值法只能得到單個位置的煙氣層高度，并且K-均值算法的計算結果能夠更加突出表征煙氣層高度的縱向區間分布特征。