一种基于元胞自动机的森林火灾蔓延模拟方法

1.本发明涉及森林火灾蔓延模拟技术领域，具体涉及一种基于元胞自动机的森林火灾蔓延模拟方法。


背景技术：

2.森林火灾，是指失去人为控制，在林地内自由蔓延和扩展，对森林、森林生态系统和人类带来一定危害和损失的林火行为。森林火灾是一种突发性强、破坏性大、处置救助较为困难的自然灾害，大空间森林火灾常带来巨大的生命财产损失和自然生态破坏。深入研究森林火灾发生、发展与演化规律，提升对林火行为的理解与掌握，对森林火灾预防、监测、控制有着重要作用，能够为森林火灾综合防灾减灾提供更加科学有效的辅助决策支撑。
3.到目前为止，国内外学者对森林火灾蔓延过程进行了大量的科学研究并提出了许多科学有效的森林火灾模拟模型：火灾蔓延经验模型主要是依据火灾历史统计数据和火灾实验模拟数据分析火灾蔓延规律，但计算复杂，模型的迁移适用性不高；火灾蔓延概率模型主要使用数学概率方法对火灾蔓延过程进行统计分析和描述，对火灾蔓延这一随机模拟过程进行数学概率化建模，从概率论的角度分析了火灾蔓延的内在规律，但不适用于大场景的模拟。
4.元胞自动机是一种时空离散的局部动力学模型，是复杂系统研究的一个典型方法，特别适合用于空间复杂系统的时空动态模拟研究，能够比较容易地整合大气环境、地形地貌等影响火灾蔓延过程的动态变量，并支撑对火灾蔓延过程的可视化模拟，能够比较直观地演示火灾蔓延发生的过程。但目前的元胞一般只包含了未燃和燃烧两种燃烧情况，不能完全表达元胞的状态，导致林火蔓延模型的准确性和实用性不高。


技术实现要素：

5.根据现有技术的不足，本发明的目的是提供一种基于元胞自动机的森林火灾蔓延模拟方法，将元胞状态划分为未燃、初燃，全燃并能够点燃周围元胞、渐灭和全灭五种燃烧状态，提高林火蔓延模型的准确性和实用性。
6.为了解决上述技术问题，本发明采用的技术方案为：
7.一种基于元胞自动机的森林火灾蔓延模拟方法，包括：
8.步骤1、获取模拟区域的数字高程数据、植被覆盖数据和自然气象数据；
9.步骤2、将模拟区域划分为规则的格网空间，每个格网单元作为一个元胞，每个元胞包含一个状态空间矢量，状态空间矢量包括元胞的燃烧状态、数字高程数据、植被覆盖数据和自然气象数据；
10.步骤3、计算邻域元胞对中心元胞的蔓延速度；
11.步骤4、以实际树木从初燃到熄灭的过程设计状态转变规则，将元胞的燃烧状态划分为未燃、初燃，全燃并能够点燃周围元胞、渐灭和全灭，基于元胞自动机根据林火蔓延速度、元胞的大小和即时元胞的燃烧状态计算下一时刻元胞的元胞的燃烧状态，实现元胞的
燃烧状态更新；
12.步骤5、输出不同时间段不同元胞的燃烧状态；
13.步骤6、渲染不同时序的元胞的燃烧状态，实现森林火灾蔓延的可视化展示。
14.进一步地，在步骤1中，植被覆盖数据包括植被覆盖类型，自然气象数据包括风力等级、风速、风向、温度和湿度。
15.进一步地，在步骤3中，根据初始林火蔓延速度、风调整系数、地形坡度调整系数和可燃物指数获取邻域元胞对中心元胞的蔓延速度。
16.进一步地，初始林火蔓延速度与风力等级、温度和湿度相关，初始林火蔓延速度r0的计算公式为：
17.r0＝0.03t+0.05w+0.01h-0.3
18.其中，t为温度，单位为℃，h为日最小湿度，单位为rh％，w为蒲福风级；
19.蒲福风级w的计算公式为：
[0020][0021]
其中，v为风速，单位为m/s，int表示取整数；
[0022]
通过风调整系数kw获取风速和风向对林火蔓延速度的增益作用，风调整系数kw的计算方式为：
[0023]kw
＝e
0.1783v
′
[0024]
其中，v
′
为邻域元胞的风速；
[0025]
通过地形坡度调整系数k
φ
获取地形对林火蔓延速度的增益作用：
[0026]
通过可燃物指数ks获取可燃植被类型对林火蔓延速度的影响；
[0027]
火从中心元胞向邻域元胞的蔓延速度的计算公式为：
[0028][0029]
进一步地，平铺松针的可燃物指数ks为0.8，红松和华山松的可燃物指数ks为1.0，枯枝落叶的可燃物指数ks为1.2，茅草和杂草的可燃物指数ks为1.6，莎草和矮桦的可燃物指数ks为1.8，牧场和草原的可燃物指数ks为2.0。
[0030]
进一步地，根据火从中心元胞向邻域元胞的蔓延方向分为西北方向、北向、东北方向、东向、东南方向、南向、西南方向和西向，将中心元胞北向、东向、南向和西向的领域元胞作为中心元胞的邻胞，将西北方向、东北方向、东南方向和西南方向的领域元胞作为中心元胞的次邻胞，其中为中心元胞向周围8个方向的领域元胞蔓延的速度分量方向，表示四个象限中的任意风向，则中心元胞西北方向、北向、东北方向、东向、东南方向、南向、西南方向和西向的8个领域元胞的速度分量分别表示为：
[0031][0032]
进一步地，邻域元胞(k，l)相对于中心元胞(i，j)的可以表示为：
[0033]
或
[0034]
其中，h
k，1
和h
i，j
表示邻域元胞(k,l)和中心元胞(i,j)中心位置的高度值，假设在一个元胞内高度值都是相同的，a表示元胞的边长大小，代表元胞的对角线长度，当元胞(k,l)是中心元胞(i,j)邻胞时使用计算公式当邻域元胞(k,l)是中心元胞(i,j)次邻胞时使用计算公式当时，g值为0，表示上坡对蔓延速度的增强作用；当时，g值为1，表示下坡对蔓延速度的抑制作用,当时，g值为0，表示上坡对蔓延速度的增强作用；当时，g值为1，表示下坡对蔓延速度的抑制作用。
[0035]
进一步地，将风方向投影到中心元胞8个领域方向上，可得西北方向、北向、东北方向、东向、东南方向、南向、西南方向和西向火从中心元胞向8个领域元胞的蔓延速度分别表示为：
[0036][0037][0038][0039][0040][0041]
[0042][0043][0044]
其中，θ为风向与正北方向的夹角。
[0045]
进一步地，所述步骤4中，元胞自动机表示为{zn，s，n，f}，其中zn表示n维元胞空间，s表示状态，n表示邻域元胞的数量，f表示为状态转换规则，林火蔓延地理元胞自动机。
[0046]
进一步地，令s＝0表示未燃，s＝1表示初燃，s＝2表示全燃并能够点燃周围元胞，s＝3表示渐灭，s＝4表示全灭；
[0047]
遍历火场中燃烧状态可能改变的元胞，满足：
[0048][0049]
其中，&&表示且，表示存在，||表示或，当元胞状态为s＝1且其邻域内存在s＝2的元胞，或者元胞的状态为s＝2或s＝3时，元胞状态可能发生改变；
[0050]
若此时s＝0，可燃且邻域元胞存在状态为s＝1的初燃元胞，则计算得到下一时刻该中心元胞的元胞值，应用公式：
[0051][0052][0053]
其中，对和取整数，l表示元胞的边长，δt为时间步长，表示元胞燃烧状态更新的时间间隔；t为当前时刻，t+δt为下一时刻；(i,j)代表中心元胞地理位置；表示即时中心元胞(i,j)的状态，表示下一时刻中心元胞(i,j)状态；表示当前时刻火从中心元胞向西北方向、东北方向、东南方向、和西南方向的4个领域元胞的蔓延速度；表示当前时刻火从中心元胞向北向、东向、南向和西向的4个领域元胞的蔓延速度；r
max
为中心元胞周围8个领域元胞的最大蔓延速度；m为步长系数；
[0054]
若此时s＝1恰好初燃，则下一时刻s＝2，全燃；
[0055]
若此时s＝2且邻域元胞s≥2或不可燃，则s＝3，逐渐熄灭；
[0056]
若此时s＝3，则下一时刻s＝4，完全熄灭。
[0057]
与现有技术相比，本发明具有以下优点和有益效果：
[0058]
本发明所述的一种基于元胞自动机的森林火灾蔓延模拟方法，充分结合元胞自动机的特点，将元胞的状态划分为未燃、初燃，全燃并能够点燃周围元胞、渐灭和全灭未燃，充分考虑了元胞燃烧的各种状态，划分5个状态更精确，具有计算量小、计算速度快、模型精度
高、可视化效果好、广泛适用于大小场景尺度下的森林火灾模拟的特点，保持了算法效率和模拟精度的平衡，为森林火灾蔓延模拟领域提供了一种有效的方法。
附图说明
[0059]
图1为本发明实施过程的流程图。
[0060]
图2为本发明将模拟区域划分为规则的格网空间示意图。
[0061]
图3为本发明元胞八邻域示意图。
[0062]
图4(a)为本发明以(i,j)为火点，其全燃并具有点燃周围元胞的能力。
[0063]
图4(b)为本发明t+δt时刻(i,j)以一定速度点燃周围8邻域的示意图。
[0064]
图4(c)为本发明t+2δt时刻8邻域全燃的示意图。
[0065]
图4(d)为本发明t+3δt时刻(i,j)渐灭，8邻域以一定速率点燃各自8邻域的示意图。
[0066]
图5(a)为本发明以(i,j)为火点可视化展示森林火灾蔓延的模拟示意图。
[0067]
图5(b)为本发明t+δt时刻可视化展示森林火灾蔓延的模拟示意图。
[0068]
图5(c)为本发明t+2δt时刻可视化展示森林火灾蔓延的模拟示意图。
[0069]
图5(d)为本发明t+3δt时刻可视化展示森林火灾蔓延的模拟示意图。
具体实施方式
[0070]
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述。
[0071]
一种基于元胞自动机的森林火灾蔓延模拟方法，包括：
[0072]
步骤1、获取模拟区域的数字高程数据、植被覆盖数据和自然气象数据；
[0073]
步骤2、将模拟区域划分为规则的格网空间，每个格网单元作为一个元胞，每个元胞包含一个状态空间矢量，状态空间矢量包括元胞的燃烧状态、数字高程数据、植被覆盖数据和自然气象数据；
[0074]
步骤3、计算邻域元胞对中心元胞的蔓延速度；
[0075]
步骤4、以实际树木从初燃到熄灭的过程设计状态转变规则，将元胞的燃烧状态划分为未燃、初燃，全燃并能够点燃周围元胞、渐灭和全灭，基于元胞自动机根据林火蔓延速度、元胞的大小和即时元胞的燃烧状态计算下一时刻元胞的元胞的燃烧状态，实现元胞的燃烧状态更新；
[0076]
步骤5、输出不同时间段不同元胞的燃烧状态；
[0077]
步骤6、渲染不同时序的元胞的燃烧状态，实现森林火灾蔓延的可视化展示。
[0078]
在步骤1中，植被覆盖数据包括植被覆盖类型，自然气象数据包括风力等级、风速、风向、温度和湿度。
[0079]
在步骤3中，影响森林火灾的影响参数包括初始林火蔓延速度、风调整系数、地形坡度调整系数以及可燃物指数，得到邻域元胞对中心元胞的蔓延速度。
[0080]
具体地，初始林火蔓延速度与风力等级、温度和湿度相关，初始林火蔓延速度r0的计算公式为：
[0081]
r0＝0.03t+0.05w+0.01h-0.3
[0082]
其中，t为温度，单位为℃，h为日最小湿度，单位为rh％，w为蒲福风级；
[0083]
蒲福风级w的计算公式为：
[0084][0085]
其中，v为风速，单位为m/s，int表示取整数；
[0086]
通过风调整系数kw获取风速和风向对林火蔓延速度的增益作用，风调整系数kw的计算方式为：
[0087]kw
＝e
0.1783v
′
[0088]
其中，v
′
为邻域元胞的风速；
[0089]
通过地形坡度调整系数k
φ
获取地形对林火蔓延速度的增益作用：
[0090]
通过可燃物指数ks获取可燃植被类型对林火蔓延速度的影响；
[0091]
火从中心元胞向邻域元胞的蔓延速度r的计算公式为：
[0092]
r＝r0·kw
·kφ
·ks
[0093]
其中，不同可燃植被类型的可燃物指数ks不同，常见的可燃物指数ks如下表所示：
[0094][0095]
本发明中，根据火从中心元胞向邻域元胞的蔓延方向分为西北方向、北向、东北方向、东向、东南方向、南向、西南方向和西向，将中心元胞北向、东向、南向和西向的领域元胞作为中心元胞的邻胞，将西北方向、东北方向、东南方向和西南方向的领域元胞作为中心元胞的次邻胞，其中为中心元胞向周围8个方向的领域元胞蔓延的速度分量方向，表示四个象限中的任意风向，则中心元胞西北方向、北向、东北方向、东向、东南方向、南向、西南方向和西向的8个领域元胞的速度分量分别表示为：
[0096][0097]
由于坡度变化会导致地形差异的变化。因此在元胞空间中，的表示方法需要改进，对于8个领域元胞中的任何一个元胞(k，l)都有各自的相对于中心元胞的值，也都有一个各自的坡度值因此，邻域元胞(k，l)相对于中心元胞(i，j)的可以表示为：
[0098]
或
[0099]
其中，h
k，1
和h
i，j
表示邻域元胞(k,l)和中心元胞(i,j)中心位置的高度值，假设在一个元胞内高度值都是相同的，a表示元胞的边长大小，代表元胞的对角线长度，当元胞(k,l)是中心元胞(i,j)邻胞时使用计算公式当邻域元胞(k,l)是中心元胞(i,j)次邻胞时使用计算公式当时，g值为0，表示上坡对蔓延速度的增强作用；当时，g值为1，表示下坡对蔓延速度的抑制作用,当时，g值为0，表示上坡对蔓延速度的增强作用；当时，g值为1，表示下坡对蔓延速度的抑制作用。
[0100]
将风方向投影到元胞8个领域方向上，得到如表1所述的邻胞方向上的夹角余弦。
[0101]
表1 邻胞方向上的夹角余弦
[0102][0103]
次邻胞方向的夹角余弦
[0104][0105][0106]
根据表1，将风方向投影到中心元胞8个领域方向上，可得西北方向、北向、东北方向、东向、东南方向、南向、西南方向和西向火从中心元胞向8个领域元胞的蔓延速度分别表示为：
[0107]
[0108][0109][0110][0111][0112][0113][0114][0115]
其中，θ为风向与正北方向的夹角。
[0116]
林火蔓延地理ca算法可描述为：由不同地理位置燃烧状态，不同火情元胞组成的地理元胞空间，按照林火蔓延的规则，随着离散的时间推进，每棵树木燃烧状态不断变化，且每棵树木的状态变化只与其所处状态以及摩尔邻域树木所处状态有关的林火蔓延系统，其中转换规则驱动整个系统运行的核心。由于在算法中，状态变化是只与它自己所处的状态和摩尔邻域的树木所处状态相关，采用8邻域相较于采用4邻域所得到的状态更准确。
[0117]
8个邻域元胞中的任何一个元胞(k,l)都有各自的相对于中心燃烧元胞的(k
φ
)
k，l
值，也都有一个各自的坡度值因此，邻域元胞(k,l)相对于中心燃烧元胞(i，j)的k
φ
需要改变传统的计算方式，可以提高模拟精度。
[0118]
元胞自动机是定义在一个由离散、有限状态的元胞组成的元胞空间，按照一定的局部规则在离散的时间上演化的动力学系统，能够抽象为{zn，s，n，f}，即{n维元胞空间，状态，邻域数，状态转换规则}。林火蔓延地理元胞自动机算法可描述为：由不同地理位置燃烧状态(s)，不同火情元胞组成的地理元胞空间(zn)，按照林火蔓延的规则(f)，随着离散的时间推进，每棵树木燃烧状态不断变化，且每棵树木的状态变化只与其所处状态以及摩尔邻域树木(n)所处状态有关的林火蔓延系统，其中转换规则s
t+δt
＝f(s
t
，n)是驱动整个系统运行的核心，n＝8为摩尔邻域。
[0119]
在步骤4中，元胞的状态划分为未燃、初燃，全燃并能够点燃周围元胞、渐灭和全灭，s＝0表示未燃，s＝1表示初燃，s＝2表示全燃并能够点燃周围元胞，s＝3表示渐灭，s＝4表示全灭；
[0120]
遍历火场中燃烧状态可能改变的元胞，满足：
[0121][0122]
其中，&&表示且，表示存在，||表示或，当元胞状态为s＝1且其邻域内存在s＝2的元胞，或者元胞的状态为s＝2或s＝3时，元胞状态可能发生改变。
[0123]
若此时s＝0，可燃且邻域元胞存在状态为s＝1的初燃元胞，则计算得到下一时刻
该中心元胞的元胞值，应用公式：
[0124][0125][0126]
其中，对和取整数，l表示元胞的边长，δt为时间步长，表示元胞燃烧状态更新的时间间隔；t为当前时刻，t+δt为下一时刻；(i,j)代表中心元胞地理位置；表示即时中心元胞(i,j)的状态，表示下一时刻中心元胞(i,j)状态；表示当前时刻火从中心元胞向西北方向、东北方向、东南方向、和西南方向的4个领域元胞的蔓延速度；表示当前时刻火从中心元胞向北向、东向、南向和西向的4个领域元胞的蔓延速度；r
max
为中心元胞周围8个领域元胞的最大蔓延速度；m为步长系数。
[0127]
若此时s＝1恰好初燃，则下一时刻s＝2，全燃；
[0128]
若此时s＝2且邻域元胞s≥2或不可燃，则s＝3，逐渐熄灭；
[0129]
若此时s＝3，则下一时刻s＝4，完全熄灭。
[0130]
m越大，算法效率越高，模拟精度越低，在这里取0.125。
[0131]
现有的基于元胞自动机的森林火灾蔓延模拟方法通常将元胞的燃烧状态划分为未燃、燃烧和熄灭三个状态，但是未燃到燃烧的状态转换并不符合实际情况，未考虑是否具有点燃周围元胞的能力，可能在其中的某一时期并不能点燃周围元胞，若直接进入燃烧状态，则不能精确展现森林火灾蔓延情况，误差较大。而本发明将元胞的状态划分为未燃、初燃，全燃并能够点燃周围元胞、渐灭和全灭未燃。充分考虑了元胞燃烧的各种状态，划分5个状态更精确。
[0132]
本发明的一个实施例中，图4(a)-图4(d)展示了元胞自动机算法模拟林火蔓延的过程。图4(a)表示(i,j)为火点，其全燃并具有点燃周围元胞的能力；图4(b)为t+δt时刻(i,j)以一定速度点燃周围8邻域的示意图；图4(c)为t+2δt时刻8邻域全燃的示意图；图4(d)是t+3δt时刻(i,j)渐灭，8邻域以一定速率点燃各自8邻域的示意图。
[0133]
图5(a)-图5(d)展示了元胞自动机模拟林火蔓延的各时刻状态。在模拟中，将气温设置为25℃、相对湿度为30％，风速为5m/s，风向为正北方向。采用本发明提供的一种基于元胞自动机的森林火灾蔓延模拟方法，能够清楚的展示森林火灾蔓延的情况。
[0134]
综上所述，本发明基于元胞自动机原理，针对现有森林火灾蔓延模型存在的不足，充分考虑了可燃物、风、温度、湿度和坡度等因素对林火蔓延速度的影响，结合以遥感与地理信息系统核心的空间信息技术，利用计算机实现林火实时动态蔓延模拟：首先输入模拟区域的数字高成数据、植被覆盖数据、自然气象数据，然后对模拟区域进行规则格网划分，并包含各单元的燃烧状态、数字高程、植被覆盖类型、湿度、温度、风速、风向、风力等级状态
值，再利用上述状态值和经验模型，分别计算初始林火蔓延速度、风调整系数、地形坡度调整系数，根据植被覆盖类型查表得出可燃物指数，并利用四个所得值和森林火灾蔓延改正模型计算得到邻域元胞对中心元胞的蔓延速度，接着根据元胞状态转变规则，遍历火场中燃烧状态可能改变的元胞，对元胞状态进行更新，其中针对s＝0、可燃且邻域存在s＝1的初燃元胞的元胞，需要根据林火蔓延速度、元胞的大小、目前的状态来计算得到下一时刻的状态，最后输出不同时间段不同点燃烧等级，并进行渲染，实现森林火灾蔓延的可视化展示。
[0135]
本发明充分结合元胞自动机的特点，充分利用风向、风速、温度、湿度、坡度等数据，具有计算量小、计算速度快、模型精度高、可视化效果好、广泛适用于大小场景尺度下的森林火灾模拟的特点，保持了算法效率和模拟精度的平衡，为森林火灾蔓延模拟领域提供了一种有效的方法。
[0136]
综上所述，本发明可靠实用，对于森林火灾蔓延的模拟具有较好的准确性和实用性，与实际进行对比验证，具有很好的效果，实用性和可行性较好。
[0137]
显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。