背景
2019 年 11 月，澳大利亚东部丛林大火肆虐，灾情不断加剧。截至目前，大火已经燃烧了4个多月，土地的过火面积已经超过 26 万平方千米，整个澳洲大陆三分之一被浓白的烟雾覆盖。近 3000 间房屋坍塌成废墟，数十人丧生火海。
图 1. 烈火下的澳大利亚：(a) 澳大利亚山火
图 1. 烈火下的澳大利亚：(b) 一个月内澳大利亚火灾地点累积图
澳大利亚的生态环境已经遭到了严重破坏，生态学家初步估算至少有 4.8 亿只哺乳动物、鸟类和爬行动物在大火中丧生。考拉是受灾最严重的动物之一，它们运动缓慢且只吃桉树叶，而桉树含油量高，极易燃烧。目前为止，已有约 3.7 万只考拉死于大火[1]。专家表示，幸存的考拉可能已不足万只，考拉可能会因此被列为濒危物种。除了考拉，澳大利亚特有的袋鼠生存环境也受到了威胁，即使袋鼠可以快速逃到未着火的区域，但生存环境和食物缺少也是个严重问题。
图 2: 烈火下的动物：(a) 考拉
https://www.zhihu.com/video/1203842124133036032
图 2: 烈火下的动物：(b) 袋鼠
https://www.zhihu.com/video/1203842504489070592
现在大火还在持续，大火烧了4个月，救火行动也持续了4个月，但是毫无效果。据大模头我的概括总结，这场大火之所以难以扑灭，是因为“天时、地利、人不和”：
- 所谓“天时”，澳大利亚 2019 年全年的降雨量是100多年的历史最低值。[2]。2019 年12月17日，是澳大利亚有史以来最热的一天，全国平均气温高达40.9摄氏度[3]。少降雨和高温天气显然有利于大火产生和蔓延。
- 所谓“地利”，澳大利亚约 70% 的国土都属于干旱或半干旱地带，它是世界上最平坦、最干燥的大陆[4]。此外澳大利亚植被覆盖率高，且天然森林中约半数以上都是桉树，桉树富含油脂，遇火即燃。地势和植被助长了大火的蔓延。
- 所谓“人不和”，一方面是澳大利亚政府救火不力，总理莫里森在森林大火危机期间甚至还带全家去了夏威夷度假。另一方面，澳大利亚的火灾很多是人为的蓄意纵火。甚至第一线救火的消防员居然会变成纵火犯[5]！
森林火灾到底是如何蔓延的呢？澳大利亚山火为什么难以控制？本文将构建森林火灾模型模拟澳大利亚山火的蔓延，并研究山火难以控制的原因。
模型
基础模型
在森林火灾蔓延的研究中，最常用的模型之一就是元胞自动机[6]。元胞自动机是复杂系统研究的一个典型方法，特别适合于空间复杂系统的时空动态模拟研究。在元胞自动机模型中，空间被离散成网格，每一个网格被称为元胞。元胞具有有限的离散状态，遵循同样的作用规则， 依据确定的局部规则作同步更新。大量元胞通过简单的相互作用而构成动态系统的演化。
经典的森林火灾元胞自动机模型是由 Drossel 和 Schwabl 在1992年给提的[7]。森林火灾元胞自动机模型定义在正方形网格上，元胞有三种状态：树（未燃烧的树），火（正在燃烧的树）和空（空地）状态。元胞下一时刻状态的更新规则如下：
- 着火： 若某元胞的状态为“树”，且其最邻近的 4 个邻居元胞中有状态为"火", 则该元胞下一时刻的状态由“树”变为"火", 用以模拟火势的蔓延。此外，状态为“树”的元胞也会以一个低概率
![澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第5张图片]( "澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第5张图片") 变为"火"，用以模拟闪电引起的火灾。
- 烧尽：状态为"火"的元胞，其状态在下一时刻将变成“空”，用以模拟着火的树经过一定时间后被烧完。
- 新生：状态为“空”的元胞, 在下一时刻以一个低概率
![澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第8张图片]( "澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第8张图片") 变为"树"，用以模拟新树的长出。
以上三条规则中，着火规则中的闪电以及新生规则都涉及概率，体现了该模型在演化过程中 出现的随机性和不确定性。图 3 是森林火灾元胞自动机模拟效果图。模拟所用程序见附录。
改进模型
基础模型中给出的经典森林火灾元胞自动机有一些缺陷，例如大火向外漫延时并不是各向同性的，此外也无法考虑风的影响。这里我们引入邻居元胞着火概率，对于经典森林火灾元胞自动机，其邻居着火概率如图 4(a) 所示，即中间元胞状态为“火”，则下一时刻上、下、左、右四个邻居着火的概率均为 1。在此基础上，为了使火以各向同性的方式向外蔓延，我们将斜对角的四个元胞也作为邻居，并规定其着火概率为 ![澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第11张图片]( "澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第11张图片") ，具体如图 4(b) 所示。
图 4. 森林火灾元胞自动机邻居着火概率图。左 (a)，右 (b)，右 (c)
在考虑风的作用时，只需要将顺风方向上的邻居着火概率适当加大，而逆风方向上的邻居着火概率适当减小，例如西风可以考虑为图 4(c) 所示的邻居着火概率。为了对比图 4 所示的三种情况，我们设置一个 ![澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第13张图片]( "澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第13张图片") 的元胞空间，初始时刻所有元胞的状态全为“树”，并从最中间元胞开始点火。三种情况的模拟结果如图 5 所示。
图 5. 经典、改进和有风的三种模型测试。左 (a)，右 (b)，右 (c)
需要说明的是在图 5 的测试模拟中，我们删掉了闪电和生长机制。图 5(a) 表明经典森林火灾元胞自动机以矩形的方式向外蔓延，(b) 表明引入斜对角元胞着火概率能使火灾以圆的方式向外蔓延，(c) 表明通过设置特定的邻居着火概率，可以模拟有风的情况。
结果
澳大利亚
接下来我们将改进的森林火灾元胞自动机模型应用到澳大利亚。我们需要根据澳大利亚植被情况构造出元胞自动机的初始状态。图 6 是 2010 年澳大利亚归一化植被指数（Normalised Difference Vegetation Index, NDVI）图。
图 6. 2010 年澳大利亚归一化植被指数图
澳大利亚植被 NDVI 值一般在 0.1 到 0.7 之间，数值越高表明植物冠层密度越大[8]。裸地和积雪的 NDVI 数值接近于零，而水域的数值为负值。各年各月的澳大利亚归一化植被指数数据都可以从澳大利亚气象局官网[9]上下载到。假设 NDVI 值和植被覆盖率的关系如下
![澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第16张图片]( "澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第16张图片")
根据上式，可以很容易将图 6 转化成元胞自动模型所需要的初始状态。具体做法是：以图 6 中每个格点植被覆盖率作为初始为“树”的概率。例如某格点植被覆盖率为 0.7，则该位置元胞的初始状态为“树”的概率为 0.7。按照这种方式，大模头构造出了图 7。
图 7. 元胞自动模型所需要的初始状态
从上图中我们发现，澳大利亚的东南部植被密度较大，北部和西南部也有小部分植被密度较大的区域，这与图 6 澳大利亚归一化植被指数空间分布一致。我们将图 7 作为初始状态，用元胞自动机很容易模拟出大火在澳大利亚蔓延的动态结果。图 8 是不考虑风的情形。
图 8. 澳大利亚大火在无风情况下的蔓延模拟结果
模拟结果显示，在植被覆盖率较低的着火点在蔓延一小片后自动熄灭了。而在植被覆盖率较高的多个着火点不断蔓延，并汇成一片，最终形成了一场巨大的森林火灾。这表明植被的覆盖率决定了能否形成大规模的火灾。这也非常容易理解，如果没什么植被，大火便断了燃料，当然很难形成大规模的火灾。如果考虑风向的影响，不失一般性，我们以西风为例，模拟结果如图 9 所示。结果显示，大火会倾向于往风的方向蔓延。
图 9. 澳大利亚大火在西风情况下的蔓延模拟结果
不同的“树”密度影响
在对澳大利亚进行模拟时，发现只有植被覆盖率较高时才能形成大规模火灾。这里我们再研究一下“树”密度与大规模火灾形成的关系。我们用火灾能否贯穿整片森林来定义大规模火灾。图 ![澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第20张图片]( "澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第20张图片") 是一片 ![澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第21张图片]( "澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第21张图片") 的森林，其中“树”元胞（绿色）的密度为 0.25。图 10 中的“树”元胞从左到右贯穿，故认为火灾能贯穿这片森林。我们对多种尺寸的森林进行了多次测试：按照不同的“树”密度，随机生成森林并检验其“树”元胞是否左右贯穿。测试结果如图 11 所示，从中可以看出大规模火灾发生的概率存在相变点，各尺寸森林发生大规模火灾的概率曲线约在 ![澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第22张图片]( "澳洲变燠洲，考拉成烤拉！澳大利亚山火为什么难以控制 ... 第22张图片") 附近相交于一点。
图 10. 树密度为 16/64=0.25 的森林
图 11. 火灾能贯穿整片森林的概率
这表明当“树”元胞的密度高于 0.6 时就极易形成大规模的火灾，即图 7 中 NDVI 值高于 0.36 的区域极易形成大规模的火灾。这一研究结果与图 1(b) 基本吻合。
结论
- 通过改进经典的森林火灾元胞自动机，能够有效地模拟火灾在有风和无风情况下的蔓延。
- 通过模拟表明，澳大利亚东部容易发生较大规模的火灾，而中西部则相对不太容易。
- 适当的降低植被密度或修建隔离带能有效防止大规模火灾。
附录
%% 微信公众号：数学模型（MATHmodels）
% 联系方式：mathmodels@outlook.com
m = 300; n = 300; % 表示森林的矩阵行列 m x n
Plight = 5e-6; % 闪电概率
Pgrowth = 1e-2; % 生长概率
% 邻居方位 d 和点燃概率 p
d = {[1,0], [0,1], [-1,0], [0,-1]};
p = [ 1, 1, 1, 1];
% % 改进元胞自动机
% d = {[1,0], [0,1], [-1,0], [0,-1], [1,1], [-1,1], [-1,-1], [1,-1]};
% p = [ones(1,4), ones(1,4)*(sqrt(1/2)-1/2)];
% % 考虑风的情况
% d = {[1,0], [0,1], [-1,0], [0,-1], [1,1], [-1,1], [-1,-1], [1,-1], [0,-2]};
% p = [ 0.80, 0.30, 0.80, 1.00, 0.12, 0.12, 0.30, 0.30, 0.8];
% 空=0, 火=1, 树=2
E = 0; F = 1; T = 2;
S = T*(rand(m,n)<0.5);
imh = image(cat(3, S==F, S==T, zeros(m,n)));
axis image;
for t = 1:3000
% 计算邻居中能传播着火的个数
sum = zeros(size(S));
for j = 1:length(d)
sum = sum + p(j) * (circshift(S,d{j})==F);
end
isE = (S==E); isF = (S==F); isT = (S==T); % 找出三种不同的状态
ignite = rand(m,n)<sum | (rand(m,n)<Plight); % 着火条件
% 规则 1: 着火
Rule1 = T*(isT & ~ignite) + F*(isT & ignite);
% 规则 2: 烧尽
Rule2 = F*isF - F*isF;
% 规则 3: 新生
Rule3 = T*(isE & rand(m,n)<Pgrowth);
S = Rule1 + Rule2 + Rule3;
set(imh, &#39;cdata&#39;, cat(3, isF, isT, zeros(m,n)) )
drawnow
end参考资料
[1] JUST LIKE US Heartbreaking Australia fire pics show devastated koala mourning its dead friend by burying its head in its hands: https://www.thesun.co.uk/news/10764264/australia-fire-koala-mourning-dead-friend/
[2] Annual climate statement 2019: http://www.bom.gov.au/climate/current/annual/aus/#tabs=RAInfall
[3] Australia logs hottest day ever — 40.9 C — for now: https://www.japantimes.co.jp/news/2019/12/18/asia-pacific/australia-hottest-day-ever/
[4] Outback Australia - the rangelands: https://www.environment.gov.au/land/rangelands
[5] 澳大利亚大火可能永远无法熄灭！消防员报复社会：不灭火，还放火: http://3g.163.com/news/article/F3441CLC05351EQY.html
[6] 周吕文，元胞自动机简介及其在数学建模中的应用（大学生数学建模竞赛指南 4.3 节）: https://github.com/mathor/book
[7] Wikipedia, Forest-fire model: https://en.wikipedia.org/wiki/Forest-fire_model
[8] Normalised Difference Vegetation Index: http://www.bom.gov.au/climate/austmaps/about-ndvi-maps.shtml
[9] Monthly NDVI Average for Australia: http://www.bom.gov.au/jsp/awap/ndvi/index.jsp
上一篇：21-22财年，澳洲境内州担保特点总结
下一篇：澳洲投资移民——已经做过888A的生意，还能再用来申请 ... |
窥视卡
雷达卡
提升卡
置顶卡
沉默卡
喧嚣卡
变色卡
千斤顶
照妖镜
