Численное моделирование динамики пожара с учетом рельефа местности и внешнего поля скоростей Текст научной статьи по специальности «Физика»

Л. Ю. КАТАЕВА, д-р физ.-мат. наук, профессор Нижегородского государственного технического университета им. Р. Е. Алексеева, г. Нижний Новгород, Россия

Д. А. МАСЛЕННИКОВ, аспирант Нижегородского государственного технического университета им. Р. Е. Алексеева, г. Нижний Новгород, Россия И. Е. БЕЛОЦЕРКОВСКАЯ, ассистент Нижегородского государственного технического университета им. Р. Е. Алексеева, г. Нижний Новгород, Россия

УДК 519.63,536.46

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ДИНАМИКИ ПОЖАРА С УЧЕТОМ РЕЛЬЕФА МЕСТНОСТИ И ВНЕШНЕГО ПОЛЯ СКОРОСТЕЙ

Излагаются результаты численного моделирования динамики распространения лесного пожара на местности заданной геометрии, равномерно покрытой лесной растительностью, под воздействием внешнего поля скоростей. Предложен подход к численному моделированию данной задачи на основе модификации метода крупных частиц на разнесенной сетке и коррекции скоростей химических реакций с учетом расхода вещества за время расчетного шага. Проведен анализ полученных результатов и их сопоставление с известными экспериментальными и теоретическими данными.

Ключевые слова: математическое моделирование; лесной пожар; подстилающая поверхность; численное моделирование; рельеф.

Введение

Известно, что внешнее поле скоростей оказывает существенное влияние на динамику развития лесного пожара, а рельеф в свою очередь является важным фактором, формирующим данное поле. Еще Кюри и Фонс [1,2] указывали, что наличие склона сказывается на увеличении интенсивности переноса тепла от пламени к горючим материалам перед фронтом пожара и что эффект, оказываемый наличием склона, являясь относительно низким при отсутствии ветра, в сочетании с последним может быть существенным. Спустя более 10 летБарроуз [3] также пришел к выводу, что склон способствует увеличению среднего размера пожара. Эксперименты Мерфи [4] продемонстрировали, что "головное" и "хвостовое" распространение огня наблюдалось на горизонтальных поверхностях и плоскостях с уклоном 14 и 27° при скоростях ветра от 0 до 4,4 м/с с шагом 0,4 м/с. Для более высоких скоростей ветра скорость движения передней кромки пожара значительно больше, чем задней. Макартур в 1968 г. [5] также получил результаты, показывающие, что склон может значительно влиять на скорость распространения огня, особенно сразу после возгорания. При этом по сравнению с плоским рельефом скорость распространения огня может возрастать в 2 раза при уклоне 10° и в 4 раза — при уклоне 20°.

В 1972 г. Роттермелем [6] была создана полуэмпирическая модель, прогнозирующая скорость распро-

странения лесного пожара в зависимости от типа растительности, погоды и геометрии рельефа. Конев в 1984 г. [7] получил зависимость скорости низового пожара. Почти 10 лет спустя Вейз [8] представил результаты серии 65 экспериментов с уклонами от -30 до 30° и при скоростях ветра до 1,1 м/с. Виегас в 2005 г. [9] исследовал связь скорости распространения огня с уклоном местности в контексте модели распространения огня на восходящих холмистых ущельях.

Среди работ, посвященных численному моделированию динамики лесных пожаров, можно выделить работы [10-12], однако в работе [10] рельеф местности никак не учитывался, а в [11, 12] рассматривалось распространение пожара по бесконечной наклонной плоскости.

Постановка задачи и особенности численного моделирования динамики распространения лесного пожара

Как уже было отмечено выше, характер распространения пожара по местности с переменным рельефом существенно отличается от его динамики на ровной или наклонной поверхности. Таким образом, при моделировании более сложного рельефа требуется учитывать эффекты обтекания, не сводящиеся к зависимости скорости ветра от угла наклона подстилающей поверхности.

© Катаева Л. Ю., Масленников Д. А., Белоцерковская И. Е., 2012

В качестве основной системы уравнений использовалась полная математическая модель пожара [11] в предположении независимости процесса от координаты у. На рис. 1 представлена геометрия задачи и указаны направления движения фронта пожара и поле скоростей. Приняты следующие обозначения: h — высота (глубина) препятствия; R — полуширина большего основания препятствия; г — полуширина меньшего основания препятствия в форме равнобедренной трапеции в сечении плоскости xOz. Для более наглядного отображения направлений потоков газовой фазы длина стрелок на рис. 1 непропорциональна скорости газовой фазы. Задача рассматривается в декартовой системе координат; z = 0 соответствует уровню моря. Высота растительности считается постоянной относительно рельефа. Начальный вид очага задается бесконечным вдоль оси у параллелепипеда заданных размеров с поддерживаемой температурой в течение заданного времени.

Для численного моделирования использовался метод частиц в ячейках Харлоу [13] на разнесенной сетке. При этом применялся метод расщепления по координатам и физическим процессам. На первом этапе определялось стационарное поле скоростей с учетом рельефа местности, аэродинамического сопротивления и вязкости слоя растительности заданной высоты. На втором этапе полученное поле скоростей использовалось в качестве начальных условий для моделирования динамики распространения пожара.

Математическая модель получена на основе общей математической модели лесных пожаров [11], включающей учет физико-химических процессов, протекающих при сушке и горении конденсированных и летучих продуктов пиролиза. Используется модель излучения [14], основанная на добавлении фиктивной диффузии в дифференциальную модель [11]. Математическая модель задачи представляется в виде законов сохранения в дифференциальной форме и дополнительных соотношений и включает в себя:

• уравнение неразрывности газовой фазы:

др5 , 3(Р5^) , 5(Р5 W)

дt

дл

дг

_ е;

(1)

• уравнение сохранения количества движения для проекций скорости и:

д(Р5 и) , д(р5и2) , д(р5UW)

дл

+ * ди | + ^ I- (2)

дt дР дл

дл У дл -Р5 Беи^й

дг

д ( ди ^*

дг

дг

W2

уравнение сохранения количества движения для проекций скорости W:

д( р 5 W) д( р д(р 5 W2)

дt дР

дг дл * дл

дл дW

дг

д ( дW

+ | + |- (3)

дг I ^ * дг

-Р5 Беё

^л/й

2 + W2;

уравнение сохранения энергии:

5

Е (ргФгСрг )Т

1_ г _ 1

д*

дл 1 ' дТ

дл

д (р5 ер5) д (р5 ер5)

дл

д дг

дТ

На

45 Я5

дг Ср5 те

(4)

+ дзЯз - 42Я2 + к, (еик - 4стТ4); уравнение распространения излучения:

д ик д*

д и я

дл У Зк Е дл

д и я

ду у ЗкЕ ду

- к, (е ик - 4стТ4);

(5)

уравнение сохранения концентраций компонентов газовой фазы:

д(р 5Са) , д(р 5Саи) д(р 5Са W) _

д*

дл [р 5 * дл

дл

дг

з- 1 Р5 * дг

дг . а _ 1, 2;

уравнение состояния газовой фазы:

Р _РаЯТ. _1 _ Е•

М ' М аЕ1 Ма ;

3 5

Е са _ 1; Еф. _ 1;

а _ 1 г _ 1

+ Я5а , (6)

(7)

уравнения сохранения объемных долей компонентов твердой фазы:

дФ 2

р дФ1 _ Я

1;

р2

ЗГ ^

р3 % _ аеЯ1 -М-Яз;

(8)

д*

М1

массовые скорости реакций соответственно пиролиза сухого органического вещества РГМ, испарения влаги, горения конденсированных продуктов пиролиза, горения летучих продуктов пиролиза:

Я1 _ к1р 1Ф1 ехр |- Ет I;

Я2 _ к2р 2ф2Т~0,5 ехр |- Ят

(9)

Я3 _ к3Бстр5Ф3^1 ехр | - —

г, м А

Земля

г, М А

Очаг / -

■ Лес (

Земля

20

40

60

80

2, М А

Очаг

->->->->->->->->->->->->->->->->->->->->->->->->-5

10

о

20

40

60

80

Рис. 1. Геометрия расчетной области (слева) и стационарное поле скоростей в отсутствии очага горения (справа): а — холм; б — ровная местность; в — овраг заданных размеров

Я5 =

к5М2,25 ехр |- х^25Х2, х > 0,05;

к5М2Т~2,25 ехр | - | Х1Х2, х < 0,05;

Е,

С

X; =

3 С

X—М

Д /Г 1

(10)

к=1 Мк

• массовые скорости образования компонентов газовой фазы:

Я 51 = -Я 3 -

Я 5 М1

2 М2

Я 52 = (1 с г - Я 5,

(11)

где ср5, р5, ф5 — удельная теплоемкость, истинная плотность и объемная доля газовой фазы; и, Ж — горизонтальная и вертикальная составляющие скорости газовой фазы; ? — время; Q — массовая скорость образования газовой фазы; Q = (1- ас)Я1 + Я2 + + (Мс /М1)Я3; ас, уг — коксовое число РГМ и массовая доля горючего газа в общей массе летучих продуктов пиролиза соответственно; Я1, Я2, Я3, Я5 — массовые скорости реакций соответственно пиролиза сухого органического вещества, испарения влаги, горения конденсированных продуктов пиролиза и горения летучих продуктов пиролиза; Ма, Мс, М—

молярные массы соответственно индивидуальных компонентов, углерода и смеси в целом; Р—давление; Цр Хр —коэффициенты соответственно динамической вязкости, турбулентной теплопроводности и турбулентной диффузии; Б — удельная площадь поверхности фитомассы полога; са — эмпирический коэффициент аэродинамического сопротивления; р;, ф; — плотности и объемные доли компонентов твердой фазы (; =1 — сухое органическое вещество,; = 2 — влага,; = 3 — конденсированные продукты пиролиза,; = 4 — зола); Т — температура среды; д2, д3, д5 —тепловые эффекты реакций соответственно испарения, горения кокса и окисления летучих продуктов пиролиза; кБ — спектральный коэффициент поглощения; с — скорость света; иЯ — плотность потока излучения; а — постоянная Стефана-Больцмана; V — величина, обратно пропорциональная скорости искусственной диффузии; къ — коэффициент ослабления; Са — массовые концентрации компонентов газовой фазы (а = 1 — кислород, а = 2 — горючие компоненты продукта пиролиза, а = 3 — инертные компоненты газовой фазы, нереагирующие компоненты продукта пиролиза и водяного пара); Е и к — энергия активации и пред-экспонента реакций соответственно пиролиза сухого органического вещества РГМ, испарения влаги,

горения конденсированных продуктов пиролиза, горения летучих продуктов пиролиза; Е1 = 9400Я К; Е2 = 6000Я К; Е3 = 10000Я К; Е5 = 11500Я К; к1 = = 36300 с-1; к2 = 600000 с-1К0,5; к3 = 1000 кг-с-1/м2;

мольс 1/м3; р1 =360 кг/м3; 5а

к5 = 31013 К2'25 удельная площадь поверхности элемента горючих материалов; Ск, Мк — концентрация и масса к-й компоненты газовой фазы; к — индекс суммирования в формуле (10).

Для получения стационарного поля скоростей с учетом заданной геометрии на свободных границах задавалась постоянная скорость ветра, аддитивно скомбинированная с условием проницаемой границы, основанной на скорости распространения бегущей волны. При этом температура всей системы была равна температуре окружающей среды, поэтому химико-физические процессы не учитывались.

Граничные условия на твердых границах: скорость равна нулю (условие прилипания); твердая граница не пропускает тепловую энергию и вещество (нулевая динамическая вязкость и теплопроводность на границе). Верхняя и боковые границы являются проницаемыми. Описанные выше граничные условия представляются следующими соотношениями:

г hmяY:

Т _ Т

и _ ие;

С а _ С а

W _

Р5 - Р5

а _ 1, 2;

(12)

х = 0:

Т _ Те; и _

ие -

Р5 - Р5е

W _ 0; Са _ Сае; а _ 1,2;

_'М (13)

X лп

Т _ Те; и _

и

Р5 - Р5е

Ь ' е ' Р

" тах 1 5

W _ 0; Са _ Сае; а _ 1,2;

_М (14)

О:

В = 0; ц* = 0; = 0; и = 0; W = 0, (15)

где О — твердая граница; е — индекс, соответствующий невозмущенному состоянию среды. Для решения полученной системы уравнений используется метод частиц в ячейках [13]. Особенностью реализации является использование разнесенной сетки [15], что обеспечивает более точную аппроксимацию граничных условий как на проницаемой, так и на твердой границе. Схема разнесенного шаблона представлена на рис. 2, где и и W — горизонтальная и вертикальная проекции скоростей соответственно; Ф — остальные переменные; — значения расчетных величин за границами области, используемые для расчета входящих потоков. Индекс '^г" означает, что данная величина рассчиты-

Ф,

Ж -

Ф,

7«*

Ж -V

Ф,

Ж -V

6 Фгш и^ Ф и

и Ф Цг 6

1 2 3 4 5 6

Рис. 2. Структура разнесенной сетки

вается на основании граничных условий (12)-(14) для проницаемых границ и (15) — для земли.

Результаты численного моделирования динамики распространения низового и верхового пожаров и их анализ

На рис. 3 представлено распределение температур для различной геометрии на моменты развития пожара 6, 10, 14 и 18 с. Скорость ветра растет с высотой, вследствие чего наблюдается наклон фронта пожара по направлению ветра. При наличии холма наклон фронта менее выражен, что связано с увеличением горизонтальной компоненты скорости над поверхностью плато вследствие его обтекания. В случае геометрии в форме оврага вихрь перенаправляет большую часть потока нагретого воздуха вверх, и фронт пожара продвигается несколько выше и дальше. Рельеф также существенно влияет на ширину фронта пожара: при распространении пожара на подъеме и ровной возвышенности имеет место расширение фронта, тогда как на спуске — сужение.

На рис. 4 представлено распределение концентраций кислорода для различной геометрии на моменты развития пожара 6, 10, 14 и 18 с. Наличие неоднородного рельефа приводит к увеличению размеров области, содержащей концентрацию кислорода менее 1 %. На рис. 4, в данная изолиния выпукла в отрицательном направлении, что связано с формированием вихря.

На рис. 5 представлено распределение концентраций летучих продуктов пиролиза для различной геометрии на моменты развития пожара 6, 10, 14 и 18 с. Препятствие в форме оврага способствует вдуву кислорода в области активного сгорания продуктов пиролиза за счет вихря, что приводит к их активному выгоранию и, как следствие, к понижению концентрации. Данный факт находит подтвержде-

тах

О 20 40 60 80 х, м 0 20 40 60 80 х, м 0 20 40 60 80 х, м

г

О 20 40 60 80 х, м 0 20 40 60 80 х, м 0 20 40 60 80 х, м

Рис. 3. Распределение температур для холма (слева), горизонтальной плоскости (в центре) и оврага заданных размеров (справа): а — 6 с; б — 10 с; в — 14 с; г — 18 с;--500 К;----- 1000 К;--1500 К

ние на рис. 3, в, из которого видно, что более высокие температуры характерны для геометрии в форме оврага. Для горизонтальной поверхности наблюдается наибольшая концентрация летучих продуктов пиролиза. Анализируя взаимосвязь данных, представленных на рис. 3 и 5, можно заметить, что тем-

пература порядка 2000 К достигается в зонах с минимальными концентрациями летучих продуктов пиролиза, что, с одной стороны, связано с необходимостью высокой температуры для существенной интенсивности реакции, а с другой — с выделением большого количества тепла в ходе реакции.

О 20 40 60 80 х, м 0 20 40 60 80 х, м 0 20 40 60 80 х, м

20 40 60 80 х, м 0 20 40 60 80 х, м 0 20 40 60 80 х, м

20 40 60 80 х, м 0 20 40 60 80 х, м 0 20 40 60 80 х, м

Рис. 4. Распределение концентраций кислорода для холма (слева), горизонтальной плоскости (в центре) и оврага заданных размеров (справа): а — 6 с; б — 10 с; в — 14 с; г — 18 с;--20 %;-----10 %;--1 %

На рис. 6 представлено поле скоростей газовой фазы для различной геометрии рельефа на моменты развития пожара 6,10,14 и 18 с. При отсутствии не-однородностей рельефа пожар вносит незначительные возмущения в поле скоростей, тогда как в случае

оврага и холма можно видеть большую интенсивность движения газовой фазы, а также более выраженные турбулентные эффекты.

На рис. 7 представлено поле скоростей собственного ветра (т. е. движения газовой фазы вслед-

О 20 40 60 80 х, м 0 20 40 60 80 1,м 0 20 40 60 80 х, м

20 40 60 80 х, м О

20 40 60 80 х, м 0 20 40 60 80 х, м

80 х, м О

80 х, м О

80 х, м

Рис. 5. Распределение концентраций летучих продуктов пиролиза для холма (слева), горизонтальной плоскости (в центре) и оврага заданных размеров (справа): а — 6 с; б — 10 с; в — 14 с; г — 18 с;--1 %;-----10%;--30 %

ствие образования газообразных продуктов горения) для различной геометрии подстилающей поверхности на моменты 6, 10, 14 и 18 с.

Под собственным ветром в данной работе понимается разность между полями скоростей, сформировавшимися на текущий и начальный моменты

времени, т. е. искажение поля скоростей за счет избыточного давления, создаваемого пожаром.

Влияние конвекции, излучения и диффузии в сочетании с внешним полем скоростей приводит к образованию "подсоса" свежего воздуха в зону горения.

20 40 60 80 х, м 0 20 40 60 80 х,н 0 20 40 60 80 х, м

б

О 20 40 60 80 х, м 0 20 40 60 80 1,м 0 20 40 60 80 х, м

в

О 20 40 60 80 х, м 0 20 40 60 80 х, м 0 20 40 60 80 х, м

80 х, м

80 х, м О

80 х, м

Рис. 6. Поле скоростей газовой фазы для холма (слева), горизонтальной плоскости (в центре) и оврага заданных размеров (справа): а — 6 с; б — 10 с; в — 14 с; г — 18 с

Работа выполнена при финансовой поддержке РФФИ 11-08-97074-р_Поволжъе_а "Математическое моделирование, мониторинг и прогноз динамики лесного пожара (на примере Нижегородской об-

***

ласти)", 13-03-91164-ГФЕН_а "Эксперименталъное исследование кинетики и механизма термического разложения лесных горючих материалов и процессов распространения пламени по их слою".

80 х, м О

80 х, м О

80 х, м

б

О 20 40 60 80 х, м 0 20 40 60 80 1,м 0 20 40 60 80 х, м

в

О 20 40 60 80 х, м 0 20 40 60 80 х, м 0 20 40 60 80 х, м

80 х, м О

80 х, м О

80 х, м

Рис. 7. Поле скоростей собственного ветра для холма (слева), горизонтальной плоскости (в центре) и оврага заданных размеров (справа): а — 6 с; б — 10 с; в — 14 с; г — 18 с

СПИСОК ЛИТЕРАТУРЫ

1. Curry J.R., Fons W. L. Rate of spread of surface fires in the ponderosa pine type of California // J. Agr. Res. — 1938. —No. 57. — P. 239-267.

2. Curry J. R., Fons W. L. Forest fire behavior studies // Mech. Engng. — N. Y., 1940. — P. 219-225.

3. Barrows J. S. Fire behavior in northern Rocky Mountain forests. Pub. No. 29. 1951. Ogden, UT: U. S. Department of Agriculture, Forest Service, Intermountain Forest and Range Experiment Station.

4. Murphy P. J. Rates of fire spread in an artificial fuel, MSc. : Thesis. — Bozeman, MT: Montana State University, 1963.

5. McArthur A. G. Fire Behaviour in Eucalypt Forests. Leaflet No. 107. 19th Commonwealth Forestry Conference, India, 1968.

6. Rothermel R. C. A Mathematical Model for Predicting Fire Spread in Wildland Fuels. Intermountain Forest & Range Experiment Station, Forest Service, US Dept. of Agriculture, 1972.

7. Конев Э. В.,ДоррерГ. А. Анализ процесса распространения лесных пожаров и палов // Теплофизика лесных пожаров. — Новосибирск : ИТФ СО АН СССР, 1984. — С. 99.

8. Weise D. R. Modelling wind and slope-induced wildland fire behavior. Ph. D. Dissertation. Berkeley : University of California at Berkeley, 1993. — 130 p.

9. Viegas D. X. A mathematical model for forest fires blowup // Combust. Sci. and Tech. — 2005. — No. 177. — Р. 27-51.

10. ПерминовВ. А. Математическое моделирование возникновения верховых и массовых пожаров : дис. ... д-ра физ.-мат. наук. —Томск, 2010.

11. Катаева Л. Ю. Анализ динамических процессов аварийных ситуаций природного и техногенного характера : дис. ... д-ра физ.-мат. наук. — Нижний Новгород, 2010.

12. Катаева Л. Ю., Белоцерковская И. Е., Масленников Д. А., Куркин А. А. Сравнение аналитического и численного решения математической модели низового пожара с учетом влияния угла наклона подстилающей поверхности // Пожаровзрывобезопасность.—2010. — Т. 19, № 11. — С. 25-31.

13. Масленников Д. А., Катаева Л. Ю., Белоцерковская И. Е. Об особенностях моделирования излучения при пожарах // Тр. НГТУ. — 2012. — № 2(95). — С. 66-75.

14. Бабкин А. В., КолпаковВ. И., ОхитинВ. Н., СеливановВ. В. Численные методы в задачах физики быстропротекающих процессов : учебник для втузов: в 3 т. Т. 3. Прикладная механика сплошных сред / Науч. ред. В. В. Селиванов. — 2-е изд., испр. — М.: Изд-во МГТУ им. Н. Э. Баумана, 2006. — 520 с.

15. Чернов В. В. Моделирование трехмерной конвекции в мантии Земли с применением неявного метода расщепления по физическим процессам / Институт геологии и минералогии СО РАН // Вычислительные технологии. — 2006. — Т. 11, № 4. — С. 73-86.

Материал поступил в редакцию 28 сентября 2012 г.

Электронные адреса авторов: kataeval2010@mail.ru; dmitrymaslennikov@rambler.ru.

Из пожарно-технического энциклопедического словаря

МОДЕЛИРОВАНИЕ ПОЖАРА — исследование процессов возникновения и развития пожара и борьбы с ним путем построения и изучения модельного очага пожара. Моделирование пожара — одна из категорий теории процессов возникновения, развития и борьбы с пожаром, на которой базируется любой метод научных исследований, как теоретический, при котором применяются различные знаковые, абстрактные модели (математический аналог физико-химических процессов возникновения и развития пожара), так и экспериментальный, использующий простейшие модели, лабораторные, крупномасштабные и натурные.

МОДЕЛЬНЫЙ ОЧАГ ПОЖАРА — искусственный очаг пожара установленной формы и размеров, предназначенный для оценки характеристик различных огнетушащих веществ и средств их подачи. В России в соответствии с международными нормами приняты модельные очаги пожаров классов А и В. Модельный очаг пожара класса А представляет собой кубический штабель из брусков древесины. В качестве горючего материала используют древесину хвойных пород определенного размера и влажности. В зависимости от размера брусков и их количества в штабеле модельный очаг пожара может иметь обозначения: минимальный — 1 А, максимальный — 20 А. Модельный очаг пожара класса В представляет собой круглый противень из листовой стали. В качестве горючего материала применяют бензин марки А-76 летний по ГОСТ 2084.