一種基于改進粒子系統的池火災實時渲染方法

蘇 謨,郭銳鋒，王麗麗，王鴻亮，馬元婧，趙玉彬

1(中國科學院大學，北京 100049)

2(中國科學院 沈陽計算技術研究所， 沈陽 110168)

1 引 言

池火災是一種常見的火災形式，是可燃液體的自然燃燒.可燃液體(如汽油、柴油等)泄露后流到地面形成液池，遇到火源燃燒而成池火災[1].隨著石油、化工等行業的發展，可燃性液體的生產、使用、存儲等過程越來越多.可燃性液體一旦發生泄漏，很容易形成池火災并發生嚴重的事故[2].本文以環境空氣質量預警平臺項目為背景，研究池火災實時模擬對事故后果模擬分析及事故風險評估尤為重要，同時實時地模擬真實感火焰一直是計算機圖形學的難點之一.

目前，池火的模擬方法主要有三種途徑，分別是基于粒子系統、數學物理模型和紋理貼圖方式.基于數學物理模型模擬是將火焰視為一種特殊的流體，采用NS(Navier-Stokers)方程求解火焰運動，這類方法需要求解一系列復雜的非線性NS方程，但求解計算量大，難以滿足實時性要求[3].郭燕[4]等人，采用外形輪廓與紋理貼圖相結合的方法生成燃燒的效果，這種方法僅構造一個二維的符合一定規律的火焰效果，但很難與周圍的環境發生交互作用及影響火焰動態變化的呈現.粒子系統理論是迄今為止用于描述不規則物體最成熟的理論之一.許多國內外研究學者采用基于粒子系統的火焰模擬方法，但粒子的控制過程過于隨機，單純使用粒子系統無法精確描述燃燒的池火.張亮[5]等人將火焰運動模擬函數融入到粒子系統中實現較好的火焰模擬效果.但未考慮環境因素且影響火焰模擬效果的物理參數有限.劉剛[6]等人采用CUDA平臺加速粒子屬性的更新速度，降低CPU和GPU之間頻繁通信，提高了渲染速度，但粒子運動模型未考慮邊界碰撞檢測，而這部分需要消耗大量的計算資源，本文采用CUDA的并行計算架構對粒子狀態更新加速并對池火邊界碰撞檢測進行優化加速.

綜上，基于數學物理模型的方法難以保證池火實時模擬，而紋理貼圖技術受紋理選擇的影響易造成模擬缺乏真實感.本文提出一種基于改進粒子系統的池火災實時渲染方法，針對粒子系統模擬火焰過于隨機的缺點，綜合文獻[1,2]構建池火數學模型(考慮了環境影響因素)，借鑒文獻[4]思路，將池火數學模型融入粒子系統，形成了改進的粒子系統模型，引入基于自適應二叉樹的剖分算法對燃液碰撞檢測，通過CPU-GPU協同及GPU并行加速實現實時渲染真實感的池火災，本文為池火災的模擬提供了新思路，具有實際應用價值.

2 池火災模型

在傳統的研究池火災模擬分析過程中，形成了很多理論與經驗的數據，運用池火災數學模型，模擬池火燃燒速率、火災發生時的液池半徑、火焰高度、火焰燃燒時間等物理量隨時間、距離的變化過程.

2.1 池火燃燒速率

Babrauskas在進行大規模池火實驗基礎上，提出了池火的燃燒速率估算公式[1]，考慮了液池大小后如公式(1)所示：

mf=mmax[1-exp(-2σβR)]

(1)

其中：mf—可燃液體的質量燃燒速率，kg/(m2·s)；mmax—可燃液體的最大質量燃燒速率(mmax)，kg/(m2·s)；σ—火焰的吸收衰減系數，(m-1)；β—平均光線長度校正系數；R—液池半徑.幾種常見油品的mmax和σβ值如表1所示.

表1 常見燃液的mmax和σβ值Table 1 Value of mmax and σβ of fuel

2.2 液池半徑

考慮到液體的部分流出，假設液體的泄露點為中心呈扁圓柱形的光滑平面上擴散，這時液池半徑R用一下公式計算：

1)瞬時泄漏(泄漏時間不超過30s)時，半徑見公式(2)：

(2)

2)連續泄漏時半徑如公式(3)所示：

(3)

式中：R—液池半徑，m；M—泄漏的液體量，kg；ρ—液體的密度，kg/m3；g—重力加速度，g=9.8m/s2；t—泄漏時間，s；部分可燃液體密度：

表2 部分可燃液體密度Table 2 Some fuel density

對于無邊界阻擋的連續泄露，隨著液池面積擴大燃燒速度加快，當燃燒速度等于泄露速度時，液池半徑達到最大，具體公式(4)如下：

(4)

其中Rmax—液池半徑,(m)；Q液體泄漏流量，(km/s)；mf液池單位面積燃燒速率， (kg/(m2·s)).

2.3 池火高度

計算池火高度，通常使用Thomas關系式，它是以物理原理和實驗量為基礎的[2]，公式(5)如下：

(5)

其中：L—火焰高度，(m)；R—液池半徑，(m)；ρa—空氣密度，ρ0=1.297(kg/m3)；g—引力常數，g=9.8m/s2；mf—燃料的燃燒速率， (kg/m2·s).

2.4 池火持續時間

(6)

公式(6)中，W為剩余燃料的質量，kg.

3 改進的粒子系統模型

池火災在燃燒過程中與傳統火焰燃燒存在一定程度的差異.在使用傳統粒子系統模擬池火災燃燒過程基礎上，本文結合池火災模型，融入池火災燃燒過程的特點，通過對發射器發射面的形狀與面積、發射器發射粒子的高度、發射器發射粒子的速度變化、火焰的顏色、粒子的初始位置重新定義，并考慮了池火與環境障礙物碰撞檢測，建立了改進的粒子系統模型.

3.1 粒子系統的定義

定義1.單個粒子定義為在實數域的一個s維向量，表示為：

Ps={prop1,prop2,…,propi,…,props|s∈N}

(7)

其中prop1,prop2,…,propi,…,props是粒子的s個屬性，通常粒子的屬性包括空間位置，粒子大小、形狀、運動速度和加速度、顏色、亮度、生存周期等.單個粒子是組成粒子系統的基本元素.

定義2.粒子映射是單個粒子到正整數集的關系映射，每個粒子都具有唯一的索引，表示為It到ps的映射：

S(t)={It→Ps|It?Z+,s∈N,t∈R+}

(8)

其中，W(i)=Ps表示索引為i的粒子的狀態和性質.

定義3.粒子系統是粒子映射集的有限集合，表示為：

F(t)={S(t)|t∈{t0,t1,…,tm}}

(9)

公式(9)中，F為粒子在系統時刻在t0,t1,…,tm的狀態集合，S(t0)是初始狀態下粒子系統狀態[7].

3.2 基于池火災的粒子系統模型

根據粒子系統的特點，結合池火災模型，建立池火燃燒速率、液池半徑、火焰高度及火焰持續時間等物理量與粒子系統渲染池火災關鍵參數的關系，形成基于池火災的粒子系統模型.

中共中央黨校(國家行政學院)教授車文輝以《生態文明視域下的鄉村振興》為題，重點闡述了鄉村振興戰略實施過程中的生態文明建設問題，指出了鄉村建設、經濟開發進程中忽視生態、破壞生態的現象以及帶來的嚴重后果，提出生態建設、生態保護是鄉村振興題中應有之義，生態文明好不好決定著鄉村振興的質量。

3.2.1 發射器發射面的形狀與面積

在池火災發生的時候，整個液池都會持續的燃燒，此時火焰是包圍在整個液池的表面，液池表面的面積和形狀作為發射器發射面的形狀和面積約束依據，基于假設液池周圍沒有狹長縫隙等特殊地形，采用圓柱形液池模擬液池形狀，發射面面積S為公式(10)所示：

S=πR2+Cs(0≤R≤Rmax)

(10)

其中R為池火災模型中液池的半徑，CS液池與周邊環境發生碰撞檢測時所影響液池的面積.

3.2.2 發射器發射粒子的高度

發射器發射粒子的高度如公式(11)所示：

H=βLmaxβ∈(0,1]

(11)

其中Lmax為當液池半徑為R時通過mmax求得的火焰最大高度值.通過β值的分布函數計算火焰在燃池中的高度變化.考慮到火焰的收攏趨勢，即大部分火焰在不受外界干擾的情況下，均向中心收攏且中心點的火焰高度最高，向四周方向高度遞減.本文使用正態分布函數來設定β值，使用正態分布計算β，則β值可以使用公式(12)得到：

(12)

其中xp和yp為粒子在初始化點的坐標x和y值,λ跟火焰的燃燒時間有關，y是正態分布函數.由公式(6)可知，池火繼續燃燒時間跟火焰高度的關系可以使用公式(13)進行控制：

λ=1-e-?t

(13)

3.2.3 發射器發射粒子的速度變化

池火粒子的運動與受力分析基于如下假設：

1) 池火粒子可以被視為理想的質點.

2) 風力在整個模擬過程中保持不變.

3) 池火粒子受到的力僅包括重力、風力及空氣阻力的影響(非主要作用力忽略不計).

根據流體力學和空氣動力學，池火粒子在空氣中受到空氣阻力的作用，包括牛頓阻力、摩擦力和粘滯力三個力，計算空氣阻力的經驗公式(14)如下[10]：

(14)

其中，S為粒子的投影面積(取1cm2)，ρ表示空氣密度(取1.29kg/m3)，Cd表示阻力系數(取0.47)，為方便微分方程求解，采用空氣阻力與速度系數的正比關系式為公式(15)所示.

F=kv

(15)

其中k為空氣阻力系數(通過經驗公式計算空氣阻力的比例系數k≈3.03×10-6)，v為粒子運動速度.

根據牛頓第二定律，分析池火粒子的受力情況和加速度可得公式(16)：

ma=-mg-kv

(16)

分量表達形式：

(17)

通過積分計算的粒子各個分量上的速度變化可得：

(18)

為了體現粒子運動的隨機性，在初始速度基礎是添加一定的隨機偏移：

V0=Vinit·rand(0.9,1.1)

(19)

其中Vinit是通過動量守恒定理求得的速度近似下限，粒子速度在上升的過程中粒子遵循動量守恒原則，粒子的初始速度需要達到一個下限的值之后才能運動到最高點，求解空氣阻力做功較復雜，可以使用近似解來求出空氣阻力做功，具體可由動量守恒公式(20)所得.

(20)

綜上，最終加上風的效果，粒子的速度變化如公式(21)所示結果：

V′(t)=V(t)+Vwind

(21)

3.2.4 池火顏色模型

顏色是一個非常重要的物理屬性，通過建立池火顏色模型來提高池火渲染真實感，采用線性插值和顏色漸變兩種方式混合來設定池火粒子的顏色.選擇粒子上升階段總高度的0→α范圍內使用顏色的線性插值，設定α(0,1)值(常用的α值可根據應用渲染需要形成配置策略).在粒子上升總高度的α→1的范圍內使用顏色漸變，顏色模型使用RGBA模型，其中R、G、B分別代表紅、綠、藍分量，A代表透明度.

在粒子上升總高度0→α的范圍內線性插值計算為公式(22)：

(22)

在公式(22)中，O(r,g,b)池火最底部，即初始池火粒子時的顏色值，M(r,g,b)是在火焰高度α處取得的顏色值減去O(r,g,b)得到的顏色的變化范圍，L為粒子根據位置計算得到的能上升的最大高度.h為粒子當前的高度粒子在上升總高度.α→1的范圍內使用顏色漸變，改變A的值從1→0，從而最終消失計算為公式(23)：

(23)

3.2.5 池火粒子初始數量

池火粒子數量決定了池火的密度，池火粒子數量也是影響可視化實時性的重要因素，池火粒子數量過小池火渲染會失真，本文采用公式(24)來確定池火粒子的數量：

(24)

其中meannum是屏幕單位區域內產生單位粒子數目的平均值，方差為varnum，Area為顯示區域， 為粒子大小.

3.2.6 池火粒子初始位置

池火粒子的初始位置決定了粒子表現火焰的密度分布，經典粒子系統中，粒子位置初始化的方法具有普遍性而不具有針對性，針對池火焰的邊界粒子稀疏，中心密集的特點，并且成正態分布，本文采用高斯分布(正態分布)與隨機數相結合的方法確定新粒子的初始位置.設粒子的初始位置為

I(t0)=Position={Xi(t0),Yi(t0),Zi(t0)},1≤I≤numb0

粒子的分布，燃點位于XOZ平面，從該平面向上燃燒(y正方向)，R′為燃燒區域半徑，則火焰粒子的分布可用公式(25)表示：

(25)

3.3 池火與環境障礙物碰撞檢測

池火與環境障礙物的碰撞檢測需要檢測粒子在運動過程中是否與環境障礙物發生碰撞，該問題可以轉化為粒子在時間間隔內所經過的路徑是否與環境障礙物發生碰撞.由于時間間隔比較短，因此粒子的運動軌跡可近似為一條線段.因此將粒子與環境障礙物的碰撞檢測問題轉化為線段與三角面片的相交檢測問題.假設場景有M個三角面片，N個粒子，其時間復雜度為O(M×N).

圖1 池火與環境障礙物碰撞檢測流程圖Fig.1 Collision test between pool fire and environmental obstacle

為了提高碰撞檢測算法的效率，本文根據粒子與環境障礙物的空間相關性，采用ABT空間剖分算法和層次包圍盒技術，ABT空間剖分算法在分割平面(由一個考慮了多種不同參數如分割面的數量和位置等的評分系統決定)的選擇上更科學、合理，從而減少需要檢測的三角面數，節省內存空間，提高算法的效率.算法流程如圖1所示.

4 基于GPU的粒子系統實現模型

在改進的粒子系統模型的基礎上，本文采用基于GPU的粒子系統實現模型，利用CPU-GPU協同及GPU并行加速滿足實時性渲染的要求.

4.1 基于GPU的改進粒子系統

本文圍繞粒子系統構造過程，即粒子產生、粒子屬性更新、粒子銷毀和粒子繪制4個階段[8]，提出基于GPU的改進粒子系統模型，如圖2所示，通過CPU和GPU協同方式實現池火實時渲染，能夠減少CPU和GPU之間的頻繁通信及系統內存和顯存間的頻繁數據傳送，降低CPU負擔[9].在CPU端主要任務包括粒子系統初始化參數設置，建立基于池火災的粒子系統模型，粒子產生及銷毀.鑒于GPU的高并行性和可編程性，將池火粒子屬性更新、池火與環境障礙物的碰撞檢測和池火粒子繪制在GPU端實現，以滿足實時性渲染要求.

4.2 基于GPU的算法實現框架

為了提高池火實時渲染效率，并結合基于GPU的改進粒子系統模型，本文提出了具體的基于GPU的算法實現框架，減少CPU和GPU的通信，實現GPU通用計算和渲染的結合，達到實時渲染池火的目的.如圖3所示，在CPU端主要實現初始化粒子系統參數、分配存儲空間(標量/矢量等)及圍繞交互應用的粒子系統邏輯控制工作.通過將粒子初始數據傳輸至GPU顯存后，由GPU端完成重要的計算和繪制工作.在GPU端主要實現基于池火災的粒子系統模型的計算、迭代更新池火粒子運行狀態、池火邊界碰撞檢測以及池火繪制工作.將計算的頂點數據直接用于GPU的渲染，利用CUDA架構的并行計算能力，將每個粒子映射到GPU處理器的每個線程，將粒子屬性數據存放在全局存儲器，每個線程塊中的各個線程唯一對應于全局存儲器中對應的一個粒子屬性數據，不同線程塊中的不同線程讀取顯存中全局存儲器中對應粒子的屬性數據進行更新計算[10].通過著色器紋理、頂點、光源等數據整合處理后進行粒子的繪制.

圖2 基于GPU的改進粒子系統模型圖Fig.2 GPU-based improved particle system model

5 實驗結果與分析

本實驗硬件配置：CPU Intel Core(TM)i3-4160 3.6GHz，四核處理器，內存為8GB，顯卡型號為NVIDIA GeForce 9600 GT，顯存為1GB.

軟件平臺為Windows 7旗艦版64位操作系統，開發語言為C++，開發工具為Visual Studio 2013，OpenGL 4.3版本，OSG3.2.3版本，NVIDIA驅動版本號為174.16.

圖4給出了采用本文池火災實時渲染方法的效果圖，采用實驗數據(保留一幀的截圖)，當前池火燃燒速率為0.39(kg/m2·s)、燃燒物質為煤油、池火半徑為17.21m、池火高度12.32m.經典粒子系統模擬火焰的過程，粒子的隨機運動一定程度上影響了火焰的真實感，本文通過引入池火災模型作為粒子系統構造的基礎，通過建立基于池火災的粒子系統模型對粒子運行狀態控制，更貼近池火應用場景.本文基于GPU的粒子系統實現模型，確保池火實時渲染效率，當前渲染的頂點數1.8M，三角形面數1.6M，幀率達到53.6fps.

圖3 基于GPU的算法實現框架圖Fig.3 GPU-based algorithm implementation framework

圖4 池火災模擬效果圖Fig.4 Pool fire simulation results

圖5 采用文獻[5]方法模擬出池火災效果圖Fig.5 Pool fire effect via literature[5]

圖6 改進粒子系統的池火效果圖Fig.6 Improved particle system pool fire effect

圖5是采用文獻[5]方法模擬池火災效果圖，由于未考慮環境因素且影響火焰模擬效果的物理參數有限，火焰效果粗糙.圖6是采用改進的粒子系統模型后池火災效果圖，本文構建的池火災模型引入池火燃燒過程的關鍵物理量，并考慮環境因素影響如風力對池火的影響，在燃燒過程中出現火焰聚攏效應，火焰真實感增強，池火的形態會隨風力變化而實時渲染.

圖7是采用文獻[4]算法模擬池火邊界碰撞效果圖，采用層次包圍盒算法來進行碰撞檢測，圖8為采用本文算法模擬池火邊界碰撞效果圖，由于采用了ABT空間剖分算法和層次包圍盒算法相結合方式，可有效進行邊界檢測并作出碰撞反應.

圖7 文獻[4]算法的池火邊界碰撞效果圖Fig.7 Boundary collision effect via literature[4]

圖8 本文算法的池火邊界碰撞效果圖Fig.8 Boundary collision effectvia our method

在文獻[6,8]僅使用CUDA架構加速粒子屬性的更新，本文在GPU端實現迭代更新池火粒子運行狀態和計算池火邊界碰撞檢測.分別用35000、400000個池火粒子來驗證本文算法同使用 CPU 計算在性能上的差異，渲染速度結果如表3所示.

表3 渲染速度對比圖 Table 3 Rendering speed comparison

從表3可以看出，在幀率和每幀繪制的三角面片數量上，本文提出的池火渲染方法比傳統CPU實現池火渲染效果更好.在渲染35000個池火粒子時傳統CPU平臺下可達到平均 30.4幀/秒，基本滿足實時交互能力.但當池火粒子數增加至400000個時，傳統CPU平臺下的 FPS 減少到了平均14.1幀/秒.而采用本文的算法在CUDA平臺下仍然可以達到平均 64.4幀/秒，計算速度是傳統CPU平臺的4倍多.因此，與傳統CPU實現的算法相比較，本文提出的算法CPU的占用率大大降低，降低CPU計算負載，滿足實時性的要求.

6 結束語

本文提出一種基于改進粒子系統的池火災實時渲染方法，將池火災數學模型融入到粒子系統，建立了改進的粒子系統模型，并綜合考慮外力環境因素及池火邊界的碰撞檢測.采用基于GPU的粒子系統實現模型，利用GPU強大的并行計算能力對粒子狀態更新加速，并使用基于自適應二叉樹的剖分算法對池火邊界碰撞檢測進行優化實現加速，有效減少了內存和顯存之間的數據傳輸，從而實現池火災實時渲染過程優化，使性能得到顯著提高，解決了池火災模擬難以實現實時性和真實感的問題.