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

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

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

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

1 災(zāi)害蔓延仿真通用模型

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

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

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

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

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

步驟4重復(fù)步驟3,直至模型不再加入新的災(zāi)害變量為止。

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

1.2 災(zāi)害元胞自動機屬性集合

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

表1 災(zāi)害屬性集合Table 1 Disaster attribute set

1.3 模型參數(shù)解析

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

(1)

(2)

(3)

(4)

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

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

2 城市火災(zāi)蔓延A.Ohgai模型改進

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

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

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

DT=1.15(5+0.5v)

(5)

式(5)中:v為風(fēng)速,m/s;DT為火勢傳播距離極限,m。當風(fēng)速v=0 m/s時,DT=5.75 m。為保守起見,將燃燒物較少的建筑的鄰域范圍取2.4 m,而燃燒物較多的建筑的鄰域范圍取4.8 m。模型元胞自動機的元胞大小為2.4 m×2.4 m,即上述兩類建筑物的火災(zāi)時空鄰域的范圍分別為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計算下一時刻的起燃概率。預(yù)定義一個和元胞空間同維度的矩陣SP,SP中的值為對應(yīng)位置元胞的易燃概率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,初始化火災(zāi)狀態(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)。

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

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

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

3 結(jié)論

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