王威, 陳勃陽, 劉朝峰, 王志濤*

(1.北京工業(yè)大學城市建設學部, 北京 100124; 2. 河北工業(yè)大學土木與交通學院, 天津 300401)

進入20世紀以來,全球范圍內(nèi)海嘯、洪水、干旱、颶風和地震等自然災害的頻率和強度日益加劇,正不斷地威脅著人類社會的生存安全。自然災害數(shù)據(jù)表明,人類社會正遭受著自然災害越來越嚴峻的威脅,損失也越來越巨大[1]。然而,由于災害時空動態(tài)蔓延復雜且后果嚴重、危險性大,很難進行實驗研究,并且受限于實驗方案的安全性,即使進行災害實驗也很難證明數(shù)據(jù)的可靠性。災害模擬仿真通過計算機數(shù)值模擬避免了災害實驗遇到的問題,不僅可以研究災害觸發(fā)→疊加→耦合→累積演化的蔓延機理,還可以對災害進行動態(tài)預測。元胞自動機(cellular automata, CA)是一種時間、空間、狀態(tài)都離散,空間相互作用和時間因果關(guān)系都為局部的網(wǎng)格動力學模型[2],在災害模擬、應急疏散等[3-4]領(lǐng)域有廣泛應用。元胞自動機應用于模擬森林火災蔓延已經(jīng)取得了很多成果[5-6],在三維森林火災模擬方面,張全文等[7]將地理信息系統(tǒng)(geographic information system,GIS)三維可視化和元胞自動機原理相結(jié)合,模擬了“3·29”溫泉鎮(zhèn)火災在三維空間中的動態(tài)蔓延變化。對于城市火災蔓延,日本學者Ohgai等[8]在2007年提出了日本城市社區(qū)火災蔓延的元胞自動機模型,首次將元胞自動機成功應用到城市火災蔓延模擬領(lǐng)域;對于地下礦火蔓延,黃光球等[9]提出了一個基于粒子系統(tǒng)理論的地下礦巷道系統(tǒng)中火災蔓延的元胞自動機模型。在擁擠踩踏仿真方面,張磊等[10]以安全出口前拱形擁擠疏散行人流為研究對象,建立了基于擁擠力的擁擠致傷元胞自動機模型;岳昊等[11]提出了基于行人相互作用力的步行設施內(nèi)行人流擁擠踩踏的元胞自動機模型。在滑坡運動時空動態(tài)模擬領(lǐng)域,意大利學者Avolio等[12]在2000年首次提出了SCIDDICA模型用于模擬滑坡、泥石流等地質(zhì)災害發(fā)展過程。對于洪水蔓延,蔡新等[13]構(gòu)建了基于體積法思想的洪水淹沒模擬元胞自動機模型;劉恒[14]采用二維圣維南方程和基于數(shù)字高程模型(digital elevation model,DEM)的單元柵格大伙房產(chǎn)流模型建立了洪水演進元胞自動機模型。對于疾病傳播,于鑫等[15]提出了一種基于隨機行走的流行病傳播與控制元胞自動機模型;譚欣欣等[16]提出了一種基于隨機行走的考慮個體異質(zhì)性和移動性的SEIR(susceptible-exposed-infected-recovered)傳染病傳播元胞自動機模型。對于恐怖襲擊,曹淑超等[17]建立了考慮恐怖襲擊模式的行人疏散元胞自動機模型,研究了不同襲擊模式和疏散場景下的疏散時間、人員傷亡和位置分布。

綜上所述,目前雖然許多災害都采用了元胞自動機進行動態(tài)仿真模擬,但尚無一種通用的建模方法?；诖?現(xiàn)提出一種面向災害蔓延仿真的元胞自動機通用模型,用于現(xiàn)有模型的改進和新模型的建立。首先建立災害蔓延的理論模型,然后將理論模型的時間和空間離散,之后尋找離散空間單元的狀態(tài)變化規(guī)則,并確定一個時間步長內(nèi)空間單元的遍歷次序,以及不同狀態(tài)的計算次序,最后,按照元胞自動機的結(jié)構(gòu)將每個時間步的計算連接起來,得到按時間步循環(huán)更新狀態(tài)的災害蔓延仿真元胞自動機模型。

1 災害蔓延仿真通用模型

由于災害本身的發(fā)生和發(fā)展往往是一系列相關(guān)事件,而概率型元胞自動機正是按特定規(guī)則演化的隨機動態(tài)系統(tǒng),對因果關(guān)系、函數(shù)關(guān)系并不確定的各種各樣災害事件來說具有強大的描述能力。而對于某些已經(jīng)被人類深入研究過的災害,其與狀態(tài)改變有關(guān)的函數(shù)和參數(shù)已經(jīng)確定,可根據(jù)已有的災害數(shù)學模型再建立一個時空離散的、自底向上計算的數(shù)學模型,并據(jù)此構(gòu)建元胞自動機模型。

1.1 模型構(gòu)建基本思路

步驟1選取具有代表性的變量作為建模初始變量。災害在空間中擴展時,會將正常區(qū)域納入災害區(qū)域,災害在空間中消退時,災害區(qū)域會變成復原區(qū)域,同時,災害區(qū)域內(nèi)部也會隨時間產(chǎn)生惡化或退化等變化。在災害演化過程中,會有許多變量隨之變化,但為簡化計算模型,只需選取具有代表性的變量作為建模初始變量。例如,火災中某一元胞是否著火以及著火后的溫度或燃燒階段,洪水中某一元胞的水深或流速等。

步驟2建立初始變量的計算方法。在取得建模初始變量后,根據(jù)現(xiàn)有理論或由觀測數(shù)據(jù)得到的經(jīng)驗公式等,建立在離散時空中的初始變量的自底向上的計算方法。在初始變量的計算函數(shù)中,會包含很多參數(shù),這些參數(shù)和初始變量一起都將作為模型的災害量。

步驟3建立災害量中變量的計算方法。在取得災害量后,根據(jù)現(xiàn)有理論或由觀測數(shù)據(jù)得到的經(jīng)驗公式等建立在離散時空中的災害量中變量的自底向上的計算方法,其中新加入的參數(shù)作為模型的災害量。

步驟4重復步驟3,直至模型不再加入新的災害變量為止。

步驟5最后,確定所建立的多自變量、多參變量且多步驟的大型函數(shù)中,所有變量在一個離散時間步內(nèi)的計算次序。

1.2 災害元胞自動機屬性集合

災害發(fā)展變化時不同災害變量計算過程相似,在劃分好元胞空間并確定了邊界條件后,每一種災害都可分為4個集合:災害量集S1、時空鄰域集S2、判定規(guī)則集S3、更新規(guī)則集S4,如表1所示。

表1 災害屬性集合Table 1 Disaster attribute set

1.3 模型參數(shù)解析

依據(jù)CA編程方式,將表1中的4類集合、元胞空間和邊界條件重新組合成一個程序整體。只需輸入數(shù)據(jù),就可計算在時間限內(nèi)各個災害變量的時空變化。模型中存在許多隨時間變化的參數(shù),每1種隨時間變化的參數(shù)都可采取如下方式計算。

(1)

(2)

(3)

(4)

上述為在1個時間步中1個時變變量的更新方法。當所有時變變量都更新完成后,進入下1個時間步循環(huán),直至規(guī)定的條件達成。

為進一步詳細闡述本模型的建模步驟,以城市火災蔓延的元胞自動機模型為例進行解析。

2 城市火災蔓延A.Ohgai模型改進

2.1 模型參數(shù)與火災蔓延規(guī)則

表2 改進A.Ohgai的城市火災蔓延模型屬性集合Table 2 Attribute set of the improved urban fire spread model of A.Ohgai

為區(qū)分不同可燃物燃燒后的火災蔓延能力,采用鄰域可變型元胞自動機,當建筑物的燃燒物較多時,這類建筑物的位置處的元胞鄰域范圍更大。根據(jù)Jirou等[20]提出的火向下風向傳播的極限距離公式計算,公式為

DT=1.15(5+0.5v)

(5)

式(5)中:v為風速,m/s;DT為火勢傳播距離極限,m。當風速v=0 m/s時,DT=5.75 m。為保守起見,將燃燒物較少的建筑的鄰域范圍取2.4 m,而燃燒物較多的建筑的鄰域范圍取4.8 m。模型元胞自動機的元胞大小為2.4 m×2.4 m,即上述兩類建筑物的火災時空鄰域的范圍分別為3×3格、5×5格。2021年6月14日,俄羅斯新西伯利亞市一座加油站爆炸,過火面積約1 200 m2,據(jù)此,將加油站位置處的元胞時空鄰域范圍取為36 m×36 m,即15格×15格。模型與時空鄰域判定有關(guān)的參數(shù)集只有按可燃物分類的建筑類型這一個參數(shù),在程序中,與時空鄰域判定有關(guān)的參數(shù)集可表示為

Ci∈{1,2,3,4},Ci∈S1

(6)

式(6)中:1表示燃燒物較多的建筑;2表示燃燒物較少的建筑;3表示加油站;4表示道路空地。

元胞起燃概率的計算方法如下。

(7)

式(7)中:Ht+1為下一時刻的似起燃概率矩陣。

步驟3計算下一時刻的起燃概率。預定義一個和元胞空間同維度的矩陣SP,SP中的值為對應位置元胞的易燃概率D2i,則有

Ht+1=Ht+1SP

(8)

式(8)中:Ht+1為t+1時刻所有元胞的起燃概率。

2.2 實例仿真模擬

選取北京市白塔寺周邊1.2 km×1.2 km范圍為例,將建筑平面染色,獲得建筑物易燃性的染色圖(圖1);根據(jù)實地調(diào)研,依據(jù)燃燒物類型及數(shù)量染色,獲得建筑物的染色圖(圖2);對建筑物燃燒時間的分類,依據(jù)燃燒時間類型將建筑平面染色,獲得建筑物的染色圖(圖3)。

圖1 易燃性Fig.1 Inflammability

計算機程序結(jié)構(gòu)如下。

步驟1采用MATLAB中的imread()命令讀取三張染色圖圖1、圖2、圖3,用imresize()命令將三張圖劃分為m×n的網(wǎng)格,再用rgb2ind()命令將染色圖轉(zhuǎn)化為索引圖。

步驟2根據(jù)索引圖數(shù)據(jù),對易燃概率矩陣SP、可燃物多少矩陣V、燃燒時間類型矩陣FT賦予初值。

步驟3初始化燃燒時間記錄矩陣HT,初始化火災狀態(tài)矩陣H,在初燃的位置處,H初始化值為2。

步驟4遍歷所有元胞,當遇到元胞狀態(tài)H初始化值大于2時,計算此處的時空鄰域。

步驟5根據(jù)計算所得的所有時空鄰域,計算下一時刻所有元胞的起燃概率。

步驟6當未燃元胞的起燃概率大于0.24時,未燃元胞變?yōu)槌跞紶顟B(tài)元胞。當初燃狀態(tài)元胞的已燃燒時間達到元胞燃燒時間類型規(guī)定的轉(zhuǎn)變時間點時,元胞進入發(fā)展階段狀態(tài)。當發(fā)展階段狀態(tài)元胞的已燃燒時間達到元胞燃燒時間類型規(guī)定的轉(zhuǎn)變時間點時,元胞進入猛烈燃燒階段狀態(tài)。當猛烈燃燒階段狀態(tài)元胞的已燃燒時間達到元胞燃燒時間類型規(guī)定的轉(zhuǎn)變時間點時,元胞進入熄滅狀態(tài)。

步驟7更新元胞處于已燃狀態(tài)的時間矩陣HT。

步驟8將元胞染色并顯示圖像。

步驟9回到步驟4,進行循環(huán)。

進行多情景城市火災蔓延仿真分析。情景一:假定居民區(qū)著火,如圖4所示;情景二:假設即將修建一座加油站,加油站著火,如圖5所示。

圖4 居民區(qū)著火蔓延情景模擬Fig.4 Simulation of fire spread in residential area

可以看出,本文模型可以實時計算不考慮救援力量情況下的火災的范圍及狀態(tài),居民區(qū)著火時火災的蔓延速度比較平穩(wěn),加油站著火時能夠瞬間影響較大范圍的可燃物。本文模型將鄰域范圍作為變量并可取不同值,解決了元胞自動機模型應用于災害模擬時無法根據(jù)實時的環(huán)境變化靈活改變?yōu)暮β臃较颉⒕嚯x、速度的問題。但火災模型的結(jié)構(gòu)和計算過程比較簡單,有待進一步改進,同時,模型中的參數(shù)取值也需進一步研究。

3 結(jié)論

城市火災蔓延仿真模型能夠快速展現(xiàn)真實的火蔓延過程,體現(xiàn)了元胞自動機的火災仿真能力。所采用多層的、鄰域可變的非標準元胞自動機通用模型的應用較為簡便。本文模型也可擴展用于離散時空中多種災害的自底向上建模的仿真模擬,即災害變量的改變描述是災害在時間尺度上的變化,災害的時空鄰域描述的是災害在空間尺度中的變化范圍,這方面還需進一步研究。