CC BY
57
Вычислительные технологии
Том 2, № 2, 1997
ЗАЖИГАНИЕ ЛЕСНЫХ МАССИВОВ ПОД ДЕЙСТВИЕМ ВЫСОТНОГО ИСТОЧНИКА
____и ____I
ЛУЧИСТОЙ ЭНЕРГИИ*"
А. М. Гришин, В. А. Перминов Томский государственный университет, Россия
С использованием двумерной осесимметричной постановки дается решение задачи о возникновении массовых лесных пожаров в результате техногенных катастроф. Приводятся результаты численных расчетов, из которых следует, что механизм зажигания в данных случаях такой же, как и для столкновительных катастроф (например, Тунгусский взрыв), но количественные характеристики процессов зажигания (время и предельные условия зажигания, форма зоны зажигания на подстилающей поверхности) значительно отличаются. Сравнение предельных размеров зон зажигания для различных мощностей ядерных зарядов в двумерном и квазиодномерном приближении показало хорошее согласование.
Известны публикации, где обсуждаются проблемы техногенной и космической безопасности Земли, обусловленные возможностью ядерной войны [1, 2] и столкновением с небесными телами [3]. Установлено, что, как правило, крупная техногенная или космическая катастрофа сопровождается возникновением массовых лесных пожаров. В связи с оценкой экологических и климатических последствий сильных пожаров представляет интерес прогнозирование воздействия этого процесса на состояние приземного слоя атмосферы. Ввиду того, что натурные исследования в решении таких проблем невозможны, актуальными являются методы математического моделирования. Этот вопрос изучался в рамках квазиодномерной постановки задачи в [4]. Ниже рассматривается задача о начальном этапе воздействия высотного источника лучистой энергии на подстилающую поверхность Земли, покрытую лесной растительностью, в осесимметричной постановке. Целью данного исследования является определение времени и размеров зоны зажигания и изучение механизма протекающих при этом физико-химических процессов.
1. Физическая постановка задачи
Пусть источник лучистой энергии находится на высоте от поверхности Земли (рис. 1). Так как размер источника мал по сравнению с радиусом Земли, то будем считать его точечным источником излучения. На рис. 1 Д — расстояние от центра источника до текущей точки поверхности лесного массива, к — высота лесного массива, ' — эпицентр, г* — радиус зоны
*© А.М.Гришин, В.А.Перминов, 1997.
^ Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований, грант №93-013-17701, и Международного научного фонда, грант №769100.
Рис. 1.
зажигания. На верхнюю границу £ = К лесного массива действует интенсивный лучистый поток, который ослабляется по мере удаления от эпицентра '. Величина плотности лучистого теплового потока при воздушном ядерном взрыве, на основе литературных данных [4-7], аппроксимируется формулой
( ЬрРт Ь I t/tm , t < tm,
^(г,Ь) = 4пО2 1 ехр(-ко(Ь/Ьт - 1)), Ь > 1т,
Рт = 4<-5, Ьт = 0.032^о0-5, (1)
где Рт — максимальная величина светового импульса в момент времени Ьо, кТ/с; Ьт — время максимума тепловыделения источника излучения, с; О — расстояние от центра источника излучения до полога леса, м; Ьр — коэффициент пропускания атмосферы; Ь — угол между направлением вектора плотности потока излучения и верхней границей растительного покрова; ^о — мощность источника, кТ; ко — аппроксимационный коэффициент (ко = 0.75).
Поступление лучистой энергии в растительный покров (¿о < £ < К) вызывает нагрев лесных горючих материалов, испарение влаги и последующее термическое разложение твердого материала с выделением летучих продуктов пиролиза. Затем последние сгорают в атмосфере, взаимодействуя при этом с кислородом воздуха. Из-за наличия силы тяжести нагретые объемы воздуха начинают всплывать вверх, поэтому процесс зажигания лесной растительности оказывается в общем случае связан с гидродинамикой течения. Ввиду того, что на периферии от эпицентра взрыва интенсивность лучистого потока в полог леса невелика, зажигания в нем не происходит. Таким образом, за время действия источника излучения формируется зона первоначального зажигания лесного массива радиусом г*. В идеальном случае она имеет в плане форму круга. Последующее ее развитие определяется взаимодействием мощных восходящих потоков с полем ветра, поскольку эти потоки выносят твердые горящие элементы в приземный слой атмосферы и осуществляют их разброс по окружающей территории, с лесной растительностью, а также с метеорологическими и географическими условиями в заданном районе. Оценки показывают, что величина конвективного теплового потока мала по сравнению с лучистым тепловым потоком, и поэтому
перенос энергии в слое растительности осуществляется в основном благодаря излучению. Это позволяет рассматривать задачу в осесимметричной постановке. Так как процесс зажигания лесного массива сопровождает совокупность различных физических факторов, то его описание целесообразно проводить на различных уровнях сложности, что позволяет оценить роль отдельных факторов, которые игнорируются с целью упрощения математической постановки задачи.
Физическая и математическая модели тепло- и массопереноса при лесных пожарах подробно были рассмотрены в [8], поэтому здесь лишь кратко укажем основные физические допущения и представления об объекте исследования, необходимые для понимания предлагаемой математической модели. Полог леса (совокупность крон деревьев) моделируется однородной, двухтемпературной, реагирующей, недеформируемой пористой средой [8]. Рассматриваются температуры конденсированной 8 и газовой фаз. К первой относятся сухое органическое вещество, влага, конденсированные продукты пиролиза и минеральная часть лесных горючих материалов. В газовой фазе рассматриваются только необходимые для описания реакций горения компоненты а — массовые концентрации (а =1 — кислород, 2 — горючие продукты пиролиза лесных горючих материалов, 3 — остальные компоненты, включая пары воды). Конденсированная фаза, представляющая горючий материал (хвоя и тонкие веточки до 6 мм), не имеет собственной скорости, и ее объемной долей, по сравнению с газовой фазой, в соответствующих уравнениях [8] можно пренебречь. С точки зрения гидродинамики, эта фаза тем не менее оказывает сопротивление перемещению воздушных масс [8]. Определяющим механизмом переноса энергии в данном случае является излучение. В пологе леса поглощает, отражает и переизлучает в основном твердая фаза. Для описания переноса в такой специфической сплошной среде будем использовать диффузионное приближение [8], что является в данном случае оправданным, так как длина свободного пробега излучения в пологе /д << к0 (ко — характерный масштаб высоты), к0 ~ 10 — 15 м, /д ~ 1 м. Цель исследования — численно найти поля температур, скорости и концентраций компонентов в любой момент времени, а также характеристики зажигания, если оно имело место.
2. Математическая постановка задачи
Направим ось г вертикально вверх, ось г — вдоль поверхности Земли (см. рис. 1). Для описания процессов тепломассопереноса в объеме с лесной растительностью используются общие законы сохранения для многофазной среды [8]. Поскольку любое перемещение воздушных потоков в атмосфере носит турбулентный характер, то для их описания используются уравнения Рейнольдса. Сформулированная выше задача в цилиндрической системе координат сводится в этом случае к решению следующих уравнений:
др 1 д . . д . . т + гдг(гру) + д! (рш)
д ( г>) + 1 д (V у2) + д ( ,ч = др +1 д д* г дг дг дг г дг
—рзс^ул/у2 + ш2 —
т,
ут,
(-+дг (-
ру'ш'
(2)
д 1 д д др (рш) + гдг(груш) + ^ (р™2) = —
д*
дг
+ 1 —
дг г дг
(—гру'ш') + д- (— рш'2^
-pscdwVv2 + w2 — pg — wm, (4)
д 1 д д 1 д _ д _
- {pcPT) + -— {rpvcPT) + — (PwcPT) = Гдг (—rPcPv'T') + -Q-z i—PcPw'T') +
+q5R5 + As(Ts — T), (5)
д 1 д д 1 д _ д _
Oí (Pca) + ГдГ (rPvca) + oz (Pwca) = ГдГ (—rPc'*v') + gz (—Pc'*w')—
—R5a — cam, a = 1, 2, (6)
1 д иГдиА + К— k(CUr — 4^) = 0, (7)
r Or \ 3к дг j Oz \ 3к Oz
dTs 4
pmcpi = k(cUr — 4aTs4) — As(Ts — T) + q3R3 — q2^2, (8)
dí
i=i
T? д^2 о M „ Mc д^4 , ,
Pi^TT = —Ri, P2^rr = —R2, Рз^-т = acRi t"T" R3, p^^— = 0, (9)
Oí OÍ OÍ M3 Oí
33 5
ca
Mn '' ca
Pe = PRTY,Mr > Y,ca = 1' =1'
m = (1 — ac)Ri + R2 + Mi R3, R51 = —R3 — 2M2 R5,
( Ei
R52 = V(1 — ac)Ri — R5, Ri = kipi^i exp —
RTs
R2 = k2P2^2T<0-5 exp ( — RTr j ' R3 = k3P<P3Saci exp f — RT
R5 = «*(ff (f) T-2-25 exp(—Ц). (10)
2
Начальные и граничные условия имеют вид:
í = 0: v = 0, w = 0, T = Te, Ts = Te, ca = cae, fi = fie, (11)
Ow OT Oca OUR . .
r = 0: v = 0, — = 0, — = 0, = 0, —- = 0, (12)
дг дг дг Or
Ow n OT П Oca n c OUR , cir n
r = re : 7— = ^ 7— = — = ° ---- + - UR = 0, (13)
Or Or Or 3k Or 2
П П OT ^ Oca ^ OUR n (лл,
z = z0: v = 0, w = 0, -7— = 0, -7— = 0, —— = 0, (14)
Oz Oz Oz
ду дш дТ дса . ,
г = ге: тт- = 0, — = 0, — = 0, = 0, (15)
дг дг дг дг
г = к: ^^^ + = 2дт(г,*).
3к дг 2
Здесь г, г — координаты, г — отсчитывается от напочвенного покрова, а г от оси симметрии; Л1 — Д3, Д5а — массовые скорости пиролиза лесных горючих материалов (ЛГМ), испарения влаги, горения конденсированных и летучих продуктов пиролиза и образования а-компонентов газодисперсной фазы; ср^ рi, — удельные теплоемкости, истинные плотности и объемные доли г-й фазы (1 — сухое органическое вещество, 2 — вода в жидко-капельном состоянии, 3 — конденсированные продукты пиролиза, 4 — минеральная часть, 5 — газовая фаза); , 8 — температура газовой и конденсированной фаз; са — массовые концентрации (а =1 — 02, 2 — СО, 3 — инертные компоненты воздуха); р — давление; и и — плотность энергии излучения; а — постоянная Стефана — Больцмана; к — коэффициент ослабления излучения; Л3 — коэффициент обмена фаз; к — тепловые эффекты, энергии активации и предэкспоненты реакций пиролиза, испарения, горения кокса и летучих продуктов пиролиза; ва — удельная поверхность элемента лесных горючих материалов; а, с, — молекулярные массы соответственно индивидуальных компонентов газовой фазы, углерода и воздушной смеси; в, с^ — удельная поверхность фитомассы и эмпирический коэффициент сопротивления полога леса; — скорость света; у, ш — проекции скорости на оси г и г; ас, V — коксовое число и массовая доля горючих газов в массе летучих продуктов пиролиза; т — массовая скорость образования газодисперсной фазы; д — ускорение свободного падения. Индексы "0" и " е" относятся к функциям в очаге горения и на большом расстоянии от зоны пожара соответственно, штрих — к пульсационной составляющей данной величины.
Система уравнений (1) — (15) описывает процессы переноса в области лесного массива, который включает в себя пространство между подстилающей поверхностью и нижней границей полога леса, полог леса и пространство над пологом. Термодинамические, теп-лофизические и структурные характеристики здесь соответствуют ЛГМ соснового леса [4, 8] и численно равны: Е1/Д = 9400 К, к1 = 3.36 ■ 104 с-1, q1 = 0, £2/Я = 6000 К, к2 = 6 ■ 105 с-1, 52 = 3 ■ 106 Дж/кг, £з/Я = 104 К, кз = 103с-1, ® = 1.2 ■ 107 Дж/кг, £5/Я = 11500 К, к5 = 3 ■ 1013, д5 = 107 Дж/кг, ср1 = 2000, ср2 = 4180, ср3 = 900, ср4 = 1000, ср5 = 1000 Дж/(кг-К), = 1000 м-1, вс^ = 0.1, ас = 0.06, р4^4 = 0.08 кг/м3, V = 0.7, р1 = 500, р2 = 1000, р3 = 200, ре = 1.2 кг/м3, с2е = 0, р3е = 0, Ре = 105 Н/м2, Те = 300 К, с1е = 0.23.
Для замыкания системы уравнений необходимо выбрать модель турбулентности. В первом приближении можно использовать равновесное приближение К — е модели турбулентности [8] (модифицированную теорию Прандтля). Более точным является полный учет диссипации и генерации турбулентной энергии во фронте пожара с использованием точной К — е модели турбулентности [8]. Очевидно, что колебания элементов ЛГМ и химические реакции во фронте лесного пожара значительно влияют на процессы генерации и диссипации турбулентной энергии, однако создание модели турбулентности, учитывающей эти эффекты, представляет собой сложную самостоятельную задачу и выходит за рамки данной работы.
Граничные условия (12) представляют собой условия симметрии, (13) — совокупность "мягких"граничных условий для у, ш, Т, са и условия прозрачности для , (14) — сово-
О 1000 2000 3000
Рис. 2.
купность условий прилипания и адиабатичности, (15) — мягкие граничные условия для v, w, T, ca и условие на верхней границе полога леса, которое получено методом сферических гармоник.
3. Методика решения задачи и обсуждение результатов
Система уравнений (2)-(10) с начальными и граничными условиями (11)-(15) для численного интегрирования была редуцирована к дискретной форме с помощью метода контрольного объема Патанкара — Сполдинга [9]. Сеточные уравнения, возникающие в процессе дискретизации, разрешались с помощью метода SIP [10]. Согласование полей скорости и давления осуществлялось итерационным образом в рамках алгоритма SIMPLE [9]. В результате численного интегрирования были получены поля температур, массовых концентраций компонентов газовой фазы, объемные доли компонентов твердой фазы, а также векторные поля скорости в различные моменты времени. На рис. 2 представлены распределения полей температур газовой (а) (1 — T = 2.6, 2 — T = 3.5, 3 — T = 7.0, 4 — T = 10.0, T = T/Te, Te = 300 K) и твердой (б) (1 — Ts = 2.6, 2 — Ts = 3.5, 3 — Ts = 7.0, Ts = Ts/Te, Te = 300 K) фаз в момент времени, когда зона зажигания имеет максимальные размеры. Представленными на рисунке кривыми описываются процессы, протекающие в пологе леса. Энергия взрыва в данном расчете составляла 1 мТ и считалось, что его высота = 1 км. Согласно [5], в этом случае время прихода ударной волны в эпицентр взрыва превышает
2 с, тогда как зажигание здесь происходит менее чем за одну секунду. Еще значительнее разность между временем зажигания и временем прихода ударной волны реализуется на границе зоны зажигания. Таким образом, можно сделать вывод о том, что в рассматриваемой ситуации зажигание полога леса осуществляется до прихода ударной волны. Этот результат согласуется с данными наблюдений за последствиями американских ядерных взрывов в Хиросиме и Нагасаки [2, 5].
Рис. 3.
Наибольший прогрев наблюдается в центральной части. Над ней происходит всплытие масс нагретых газов за счет действия архимедовой силы. Векторная картина течения представлена на рис. 3. Следует отметить, что в окрестности очага тепломассовыделения возникает собственный ветер, который направлен с периферии к центру области зажигания. Однако скорость ветра на момент воспламенения невелика и для рассматриваемых времен не превышает 0.5 м/с. С процессом зажигания неразрывно связано распределение массовых концентраций горючих летучих продуктов пиролиза (наиболее характерными из них являются оксид углерода (рис. 4, а) (1 — с2 = 1, 2 — с2 = 10-2, 3 — с2 = 10-3, С2 = С2/С1е) и кислород (рис. 4, б) (1 — С1 = 0.7, 2 — С1 = 10-2, 3 — с1 = 10-3, с1 = 01/с1е), который расходуется в процессе горения. В центральной части рассматриваемой области происходит выгорание СО и связанное с этим уменьшение количества кислорода.
Любопытная картина наблюдается на границе области зажигания (см. рис. 4, ), где имеется немонотонность в распределении СО, связаного с процессом зажигания, в результате которого он расходуется. Следует отметить, что процесс зажигания носит газофазный характер, т. е. вначале происходит прогрев твердой и газовой фаз и испарение влаги, затем начинается процесс разложения на конденсированные и летучие продукты пиролиза, последние из которых воспламеняются на верхней границе полога леса. Зона образования летучих горючих продуктов пиролиза по мере прогрева полога леса продвигается в его глубину. На границе области зажигания образуются также продукты, которые не воспламенились в силу недостаточной плотности лучистого потока. За счет архимедовой силы эти продукты всплывают над пологом леса.
На рис. 5 представлены изменения температур газовой (пунктир) и твердой (сплошная линия) фаз на верхней границе полога леса с течением времени для различных расстояний от центра зажигания. Из графиков видно, что реализуются три режима: первый (1) пред-
Рис. 4.
ставляет собой фактически сжигание, второй (2) — так называемый нормальный режим зажигания и третий (3) — невоспламенение, что согласуется с результатами [4].
На рис. 6 представлены распределения массовых концентраций 02 (сплошные линии) и СО (пунктир) для тех же данных, которые использовались при получении кривых на рис. 5. Режимам сжигания и воспламенения (см. пунктирные линии 1, 2 на рис. 6) характерно образование максимумов выделения газообразных горючих продуктов пиролиза, которые затем резко падают за счет выгорания. Соответственно уменьшается и количество кислорода.
Кроме того, были проведены расчеты по определению размеров зон зажигания лесных массивов в рамках квазиодномерного приближения, т. е. считалось, что на начальной стадии зажигания перенос энергии в вертикальном направлении преобладает над горизонтальным. По этой причине предполагалось, что все параметры зависят от времени £ и вертикальной координаты г, а значение радиуса г использовалось как параметр при определении значения qrn (г, £) [4].
На рис. 7 представлено сравнение размеров зон зажигания в зависимости от мощности взрыва, вычисленных в квазиодномерной (кривая 1) и двумерной (кривая 2) постановках. Следует отметить, что радиусы зон зажигания, полученные в рамках квазиодномерного приближения, приблизительно на 5-10% больше найденных в двумерной постановке. По-
Рис. 5.
Рис. 7.
следнее объясняется тем, что в данном приближении не учитывается перенос энергии в горизонтальном направлении. Однако приведенное выше сравнение оправдывает применение квазиодномерного подхода для определения максимальных размеров зон зажигания лесных массивов. Поскольку при зажигании лесных горючих материалов образующиеся продукты пиролиза полностью не выгорают, то в пологе леса и над ним возникает облако газообразных горючих продуктов, смешанных с воздухом. Это облако может взорваться под действием ударной волны, приходящей при взрыве ядерного заряда [11], что имеет место в случае, если время зажигания меньше времени прихода ударной волны. Оценки показывают, что для ядерного взрыва = 2.5 мТ, произведенного на высоте =6.5 км над лесным массивом, до момента прихода ударной волны в зоне зажигания образуется до 16 • 107 кг газообразных горючих продуктов пиролиза. Необходимо отметить, что приведенная оценка не исчерпывает анализа проблемы взрыва газообразных продуктов пиролиза. Возможно, что их инициирование происходит в результате действия не падающей, а отраженной от поверхности Земли ударной волны. Иными словами, механизм детонации топливно-воздушного облака заслуживает отдельного исследования.
Полученные результаты справедливы только до момента прихода ударной волны на подстилающую поверхность. По оценкам, мощная ударная волна должна взорвать газообразные продукты пиролиза, смешанные с воздухом, в результате чего пламенное горение может прекратиться. Однако согласно [2], данный эффект носит временный характер, так как пламенное горение возобновляется из оставшихся очагов тления. В результате этого возникает массовый лесной пожар.
Таким образом, для создания всеобъемлющей теории возникновения массовых лесных пожаров необходимо помимо упомянутых выше физико-химических процессов учитывать еще и взаимодействие падающей ударной волны с очагом зажигания.
Список литературы
[1] Климатические и биологические последствия ядерной войны (Под ред. Е. П. Велихова), Наука, М., 1987.
[2] Питок Б., Акермен Т., КрутцЕн П. и др. Последствие ядерной войны: Физические и атмосферные эффекты. Мир, М., 1988.
[3] Бронштэн В. А. Физика метеорных явлений. Наука, М., 1981.
[4] Гришин А. М., Перминов В. А. Математическое моделирование зажигания полога леса от Тунгусского метеорита. Изв. СО РАН, Сиб. физ.-техн. журн., №6, 1992, 112— 117.
[5] Действие ядерного оружия (Под ред. П. С. Дмитриева). Воениздат, М., 1965.
[6] ГостинцЕв Ю. А., Махвиладзе Г. М., Новожилов В. Б. Начальная стадия развития большого пожара, инициированного излучением. В "Горение конденсированных систем: Материалы IX Всесоюз. симп. по горению и взрыву", Черноголовка, 1989.
[7] ГостинцЕв Ю. А., Махвиладзе Г. М., Новожилов В. Б. Формирование большого пожара, вызванного излучением. Изв. РАН, МЖГ, №1, 1992, 17-25.
[8] Гришин А. М. Математические модели лесных пожаров и новые способы борьбы с ними. Наука, Сиб. отд-ние, Новосибирск, 1992.
[9] Патанкар С. В. Численные методы решения задач теплообмена и динамики жидкости. Энергоатомиздат, М., 1984.
[10] Stone H. L. Iterative solution of implicit approximations of multi-dimensional partial differential equations. SIAM J. Num. Anal., No. 5, 1965, 530-558.
[11] Гришин А. М., Ефимов К. Н., ПЕрминов В. А. Математическое моделирование зажигания лесных массивов в результате столкновительных и техногенных катастроф. В "Вычислительные технологии", ИВТ СО РАН, Новосибирск, 4, №13, 1995, 156-197.
Поступила в редакцию 15 сентября 1995 г.
