«ГАЗОВЫХ И ДВУХФАЗНЫХ ВЫБРОСОВ В ОТКРЫТОЙ АТМОСФЕРЕ ...»

Глава 2. Подъем термика и перенос примесей 2.3. Начальное состояние запыленного термика Как было показано в Главе 1 при анализе приближенной модели термика в стратифицированной по плотности атмосфере, тепловая энергия термика представляет собой интеграл уравнения притока тепла и сохраняется неизменной на начальной и автомодельной стадиях движения, до тех пор, пока не начнет сказываться температурная стратификация среды (см. (1.9)). Покажем, при каких условиях тепловая энергия может считаться интегралом движения в сжимаемой среде. Для этого рассмотрим случай чисто газового термика, переход к которому в уравнениях (2.10)–(2.17) осуществляется при Yp = 0 (в этом случае p = 0). Введем приращения давления и температуры по отношению к их значениям в невозмущенной температуре p = P Pa, = T Ta, преобразуем уравнение притока тепла (2.12) при помощи уравнений состояния и неразрывности (2.10) к дивергентному виду, после чего заменим производные от Pa и Ta по вертикальной координате z через соотношения (2.18). В результате получим уравнение притока тепла в стратифицированной сжимаемой атмосфере совпадающее с уравнением (1.6) Главы 1 при p 0 и a. Интегрирование (2.22) по полупространству z 0 с учетом граничных условий на оси симметрии, подстилающей поверхности и на бесконечности (2.15)–(2.17) дает баланс тепла Q = dV в термике:

где dV = 2 r d d — элемент объема. Первый член в правой части (2.23), по поr z рядку величины равный O(Ma2 ), описывает расход энергии на генерацию акустиГлава 2. Подъем термика и перенос примесей ческих возмущений, второй и третий — охлаждение (при J > 0) или нагрев (при J < 0) объема газа при его вертикальном перемещении в атмосфере. Из физических соображений ясно, что вклад первого и третьего членов мал для рассматриваемых здесь существенно дозвуковых течений. Второй член дает существенный вклад только на стадии зависания, когда стратификация начинает играть определяющую роль.

Поэтому на начальном и автомодельном участках подъема тепловая энергия сохраняется постоянной: Q(t) = const = Q0.

В работах [16, 17] на основе анализа экспериментальных данных по динамике турбулентных термиков показано, что на автомодельном участке движения термика справедлив закон подъема (см. (1.1)) где zt — координата верхней кромки облака, z — координата виртуального источника, B0 — полный начальный запас плавучести облака, который в случае несжимаемой среды является интегралом движения наряду с полным запасом тепла Q0. Поскольку для сжимаемой атмосферы сохраняющейся величиной является лишь интеграл тепла, удобно представить эту экспериментальную зависимость в безразмерном виде где I0 с точностью до множителя 2 представляет собой безразмерную тепловую энергию термика Q0, являющуюся интегралом движения.

Поскольку величина тепловой энергии термика определяется мощностью взрыва и, как правило, хорошо известна, целесообразно использовать ее в качестве определяющего параметра задачи, а начальный радиус термика R0, входящий в гауссово распределение (2.14), определять численно исходя из соотношения I 0 = d d. В случае чисто газового облака начальный перегрев 0 задавался, исходя из максимальной температуры в центре термика 6000 K, характерной для ядерных взрывов [18, 108].

Глава 2. Подъем термика и перенос примесей При интерпретации результатов расчетов применительно к термикам ядерных взрывов полезна связь безразмерной начальной энергии I0 c мощностью взрыва W, выраженной в мегатоннах. Для используемого во всех расчетах набора характерных масштабов (2.19), (2.20) эта связь имеет вид I0 = 0,676W. При этом полная выделяющаяся энергия Q = 4,18 · 1015 W Дж, а при вычислении переводного коэффициента учитывалось, что в тепловую энергию переходит лишь часть общей энергии Q0 = 0,35Q [16, 18]. Представляет интерес сравнить получаемые таким образом размеры термика с имеющимися в литературе данными по термикам ядерных взрывов. Так, для облака, соответствующего взрыву с W = 1 на высоте H0 = 1,5, характерная полуширина гауссового распределения температуры (2.14) оказалась равной R0 = 0,59, а радиус облака, соответствующий избыточной температуре на уровне 1% от максимальной, составил 0,95, (в размерных единицах 950 м), что хорошо согласуется с приведенной в литературе величиной 1000 м [18, 108].

При наличии в теле термика распределенных по объему дисперсных частиц выделяющееся при образовании термика тепло Q0 распределяется между газовой и дисперсной фазами:

что в безразмерном виде выражается как Перераспределение энергии может вызвать изменение как температуры, так и размера облака. Для нахождения начальных распределений параметров в термике требуется детальное рассмотрение процесса формирования термика при точечном взрыве с учетом вовлечения дисперсной фазы, что само по себе представляет весьма сложную задачу. Поэтому в настоящей работе предполагается, что радиус запыленного термика R0 равен радиусу чисто газового облака, которое образовалось бы при выделении суммарного тепла I0 в отсутствие дисперсных частиц, а температура в термике изменяется по сравнению с чисто газовым случаем за счет поглощения части Глава 2. Подъем термика и перенос примесей тепла частицами. После нахождения радиуса R0, по заданной безразмерной полной массе примеси Mp = Mp /0 L3 из распределения плотности примеси (2.14) можно определить максимальную плотность 0, а затем, разрешая численно интегральное соотношение (2.26), определить максимальную начальную избыточную температуру 0.

Начальное состояние термика с дисперсными частицами удобно характеризовать двумя параметрами, описывающими тепловой баланс в термике и действующие силы. Способность примеси запасать часть выделившегося тепла будем описывать параметром = Cp Mp T0 /Q0 = p Mp /2I0, равным отношению тепла, поглощаемого примесью при повышении ее температуры на величину T0, к полной тепловой энергии термика Q0. Этот параметр позволяет единым образом описать распределение тепла между фазами. На рис. 2.1 кривыми 1, 2 показаны доли энергии g = Ig /I0 и p = Ip /I0, запасенной в газовой и дисперсной фазах соответственно.

Данные зависимости получены в широком диапазоне начальных параметров термика 0 Mp 10, 0 p 10, 0,1 I0 2,7 при H0 = 1,56, 0 = 21. Здесь же кривой 3 показано отношение максимальной начальной температуры 0 в запыленном термике к соответствующей величине 0 в чисто газовом облаке, демонстрирующее понижение температуры за счет передачи части энергии дисперсной фазе.

Рассмотрим теперь силы, действующие на нагруженный частицами газ в термике, выделив отдельно суммарный вес частиц F и суммарную подъемную силу F+, действующую на газ. В безразмерных переменных F = Mp, тогда как Если начальная высота термика не слишком велика, можно положить в знаменателе Ta 1 (для наибольшей использовавшейся в расчетах высоты H0 = 6 это выполняется с погрешностью порядка 15%), что дает F+ 2Ig = 2I0g. Весовое влияние примеси будем характеризовать параметром, равным отношению веса частиц Глава 2. Подъем термика и перенос примесей к подъемной силе Связь между величинами и имеет вид p = /g (). На рис. 2.1 показана зависимость функции = p от параметра (кривая 4). Рост загрузки облака примесью (увеличение ) приводит к уменьшению подъемной силы F+ вследствие уменьшения запасаемого газом тепла и к одновременному росту суммарного веса частиц. При достаточно большой загрузке, силы F+ и F сравниваются ( = 1), при этом критическое значение параметра = можно найти из решения уравнения p = /g ( ). При p = 1 это значение равно 0,34. Заметим, что если частицы не влияли бы на подъемную силу F+ за счет отбора тепла у газа, критическое значение параметра было бы равным = 1 = 0,71. Следовательно, тепловое воздействие примеси снижает критическую загрузку примерно вдвое.

2.4. Автомодельный подъем термика, весовое и тепловое влияние Рассмотрим сначала подъем чисто газового термика, возникающего при мощном взрыве вблизи земной поверхности. Это позволит выделить особенности, связанные со сжимаемостью среды и рассмотреть отдельно эволюцию облака в отсутствие влияния дисперсной фазы. Для того, чтобы характеризовать динамику подъема облака на стадии, соответствующей корневому закону роста координаты его верхней кромки, определим автомодельную координату верхней кромки t :

где производную следует брать на прямолинейном участке зависимости zt (t1/2 ) (см.

рис. 1.1). Экспериментальная величина автомодельной координаты верхней кромки, равная числовому множителю в законе подъема (2.24), составляет в среднем texp 4,35 [16, 17]. Как показывает анализ размерностей (см. Главу 1), в качестве основного определяющего параметра задачи на автомодельной стадии движения Глава 2. Подъем термика и перенос примесей Рис. 2.1. Начальное распределение тепла между газовой и дисперсной фазами в термике, содержащем дисперсные частицы: 1 — доля тепла в газовой фазе g, 2 — доля тепла в дисперсной фазе p, 3 — отношение начальных температур запыленного и газового термиков 0 / 0, 4 — отношение веса примеси к выталкивающей силе p Глава 2. Подъем термика и перенос примесей выступает число Грасгофа, характеризующее соотношение между силами плавучести и вязкости:

Для того, чтобы найти зависимость автомодельной координаты верхней кромки термика от этого параметра, а также от параметров задачи, определяющих начальное состояние термика, была проведена серия расчетов в диапазоне I0 = 0,34 2,71, 0 = 5 21, H0 = 1,56 7,42, Re = 20 60, J = J1, что соответствует варьированию числа Грасгофа в пределах 200 Gr 2400. Для приведенного выше набора определяющих параметров (2.19), (2.20) это соответствуют термикам, образующимся при взрывах мощностью W = 0,5 4,0 (Q = (2,1 16,7) · 1015 Дж) на высотах 1560–7420 м и поднимающимся в тропосфере с типичным значением параметра стратификации.

Результаты расчетов, проведенных при 0 = 21, I0 = 0,68 представлены на рис. 2.2 в виде зависимостей автомодельной координаты верхней кромки от числа Грасгофа t (Gr), полученных при начальных высотах облака H0 = 1,56, 3,32, 5, и 7,42 (точки 1–4). Чтобы убедиться в том, что параметры 0 и I0 влияют на динамику подъема не индивидуально, а лишь через число Грасгофа, были проведены соответствующие расчеты. Так, при I0 = 0,68 и H0 = 1,56 варьировалась начальная температура термика 0 : во всем диапазоне 5 0 21 отличия в значениях автомодельной координаты при этом не превышали 2%. Варьирование мощности взрыва I0 в диапазоне I0 = 0,342,71 также дало слабые отклонения: на рис. 2.2 результаты таких расчетов при H0 = 1,56 нанесены точками 5–7, которые практически ложатся на кривую, соответствующую базовому значения I0 = 0,68 (точки 1). Таким образом, можно выделить два основных параметра, определяющих величину автомодельной координаты верхней кромки t — число Грасгофа и начальную высоту термика H0.

Функция двух параметров t (Gr, H0 ), соответствующая кривым рис. 2.2, может быть представлена в виде произведения двух функций, каждая из которых зависит Глава 2. Подъем термика и перенос примесей Рис. 2.2. Зависимости автомодельной координаты верхней кромки термика t от числа Грасгофа Gr при различных начальных высотах (точки 1–4) и тепловых энергиях (точки 5–7) облака Глава 2. Подъем термика и перенос примесей лишь от одного параметра:

t = (Gr)(H0 ).

Вид функций и может быть восстановлен по результатам численных расчетов, представленных на рис. 2.2. Каждая из этих функций определяется с точностью до постоянного множителя. Будем считать, что функция, описывающая зависимость автомодельной координаты верхней кромки от высоты облака, удовлетворяет условию (0) = 1. Тогда (Gr) представляет собой единую кривую, описывающую подъем термика при единичном начальном давлении и температуре, а функция (H0 ) дает поправку на реальные условия на высоте взрыва.

Обе функции и, позволяющие в широком диапазоне параметров описать подъем мощных турбулентных термиков в неоднородной сжимаемой атмосфере, представлены на рис. 2.3. Функция (Gr) представлена на рис. 2.3а сплошной кривой, здесь же показаны расчетные точки, соответствующие данным рис. 2.2 (разброс точек не превышает 2%). Функция построена на рис. 2.3б в зависимости от приведенной начальной высоты термика H0 = Ma2 H0, равной отношению высоты H0 к характерному вертикальному масштабу изменения плотности невозмущенной атмосферы 1/Ma2 (кривая 1). Здесь же представлены зависимости от безразмерного давления PH и плотности атмосферы H на высоте взрыва (кривые 2 и 3 соответственно). По результатам численных расчетов найдены следующие аппроксимационные формулы, описывающие функции и :

Найденные зависимости позволяют по заданной тепловой энергии термика и его начальной высоте определить значения турбулентных коэффициентов переноса, обеспечивающих совпадающий с наблюдаемым в экспериментах закон подъема термика на автомодельной стадии, т. е., равенство рассчитанной и наблюдаемой в экспериментах координат верхней кромки термика. Для этого по начальной высоте центра термика следует, воспользовавшись рис. 2.3б или аппроксимацией (2.31), найти знаexp exp чение функции. Тогда из условия t = t (где t 4,35 — см. 2.24) получим, Глава 2. Подъем термика и перенос примесей что требуемое значение функции составит 4,35/, после чего из рис. 2.3а или формулы (2.31) определяется число Грасгофа Gr. Далее, зная тепловую энергию термика I0, можно определить число Рейнольдса Re = (Gr/I0 )1/2 (см. (2.29)).

Перейдем теперь к анализу подъема запыленного термика. Динамика подъема облака, нагруженного дисперсной примесью, может отличаться от чисто газового термика за счет активного воздействия частиц. Сравним скорость подъема запыленного и газового термиков на автомодельном участке движения. Для этого, аналогично (2.28), (2.29), введем автомодельную координату верхней кромки термика t и турбулентный аналог числа Грасгофа Gr, определив их через безразмерное количество тепла, запасенного газовой фазой Ig :

При уменьшении количества примеси в термике величины t и Gr стремятся к соответствующим значениям t и Gr в чисто газовом облаке с тем же суммарным запасом тепла I0, поэтому, сравнивая автомодельную координату верхней кромки газового и запыленного термиков можно судить о степени влияния дисперсной фазы в процессе подъема. Результаты расчетов движения термика, переносящего активную примесь, представлены на рис. 2.4 в виде зависимости t Gr. Термик первоначально располагался на высоте H0 = 1,56 (при этом безразмерная начальная высота составила H = Ma2 H0 = 0,2); отношение теплоемкостей фаз полагалось равным p = 1, поскольку теплоемкость многих твердых веществ (пыль, сажа, песок) близка к теплоемкости воздуха при нормальных условиях.

Сплошная кривая на рис. 2.4 отвечает подъему чисто газового термика и описывается соотношениями (2.30). Точки 1–6 получены при суммарной энергии облака I0 = 0,68 для = 0,1; 0,2; 0,3; 0,5; 0,7 и 0,95 соответственно. При малой загрузке, когда величина не превышает значения 0,4, автомодельные координаты верхней кромки чисто газового и нагруженного термиков в диапазоне 250 Gr с хорошей точностью совпадают: точки 1–3 ложатся на сплошную кривую. Дальнейшее увеличение количества примеси ( > ) приводит к замедлению подъема Глава 2. Подъем термика и перенос примесей Рис. 2.3. Функции и, задающие единую зависимость автомодельной координаты верхней кромки облака в сжимаемой атмосфере. Точки 1–4 соответствуют данным рис. 2. Глава 2. Подъем термика и перенос примесей Рис. 2.4. Зависимость автомодельной координаты верхней кромки t термика, нагруженного дисперсной примесью, от числа Грасгофа Gr.

Глава 2. Подъем термика и перенос примесей облака, вызванному весовым, инерционным и тепловым воздействием примеси (точки 4–6). Аналогичные результаты, полученные при загрузках = 0,2 и 0,5 и при других значениях начального запаса тепла I0 представлены на рис. 2.4 точками 7, 8 (I0 = 0,34) и 9, 10 (I0 = 2,7). Во избежание неправильной интерпретации этих результатов подчеркнем, что равенство автомодельных координат t запыленного и чисто t газового термика при малых загрузках означает, что термик с частицами поднимается на автомодельном участке с той же скоростью, что и чисто газовый термик, имеющий тепловую энергию Ig, меньшую, чем полная выделившаяся энергия I0. Такой термик поднимается медленнее, чем чисто газовый термик, обладающий полной тепловой энергией I0, что является проявлением активной роли примеси.

Перейдем теперь к анализу механизма воздействия примеси на движение газа в процессе подъема запыленного термика. Для этого сопоставим результаты расчетов термика с дисперсной примесью, проведенных при одинаковых начальных условиях, но с изменением параметров, характеризующие степень взаимодействия газовой и дисперсной фаз. В качестве базового варианта выбрано облако с параметрами I 0 = 0,68, p = 1, H0 = 1,56, = 0,5 (при этом = 0,2, g = 0,57, Ig = 0,39 — см.

точки 4 на рис. 2.4). В основном варианте загрузка термика немного превышает, отклонения от единой зависимости (2.30) только начинают проявляться. При такой загрузке для малых Gr термик всплывает медленнее чисто газового, тогда как при больших Gr — быстрее (ср. на рис. 2.4 сплошную кривую со штриховой, проходящей через точки 4).

На рис. 2.5 представлены результаты расчетов движения чисто газового и запыленного термика при последовательном учете весового и теплового взаимодействия фаз, проведенных при Gr = 1650. Использованы безразмерные координаты = t/tJ, zt = (t H0 )/Ig tJ, где tJ = J 1/2 — характерное время подъема термика, определяемое параметром стратификации в нижних слоях атмосферы (тропосфере). Эти координаты не зависят от выбора линейного масштаба L. Динамика подъема чисто газового термика с энергией Ig, содержащего пассивную в тепловом и Глава 2. Подъем термика и перенос примесей весовом отношении примесь, показана кривой 1. Тангенс угла наклона прямолинейного участка равен величине автомодельной координаты верхней кромки t, которая при данном числе Грасгофа составляет 4,45. В этих расчетах для исключения взаимодействия фаз во всех уравнениях плотность примеси p полагалась нулевой, а для отслеживания динамики подъема дополнительно вводилась пассивно переносимая субстанция, для концентрации которой решалось уравнение вида (2.4).

Чтобы учесть только весовое воздействие примеси, рассматривалось облако с начальными параметрами, соответствующими основному варианту, однако теплоемкость частиц полагалась нулевой (что соответствует p = 0 в уравнении (2.3). В этом случае дисперсная фаза пассивна в смысле обмена теплом, но активна с точки зрения ее весового воздействия. Динамика подъема такого термика представлена на рис. 2. кривой 2, автомодельная координата верхней кромки составляет t = 3,40. Кривая получена расчетом с параметрами основного варианта, но с малым начальным содержанием частиц ( = 3,3 · 102) и большой теплоемкостью частиц. При этом весовое воздействие примеси несущественно, но начальная эффективная теплоемкость среды остается неизменной за счет специального выбора величины p, обеспечивающей то же значение параметра p, что и в основном варианте. Таким образом, примесь активна с точки зрения ее теплового влияния, но пассивна в весовом отношении. Для кривой 3 автомодельная координата верхней кромки составляет t = 5,10.

Наконец, результаты расчета основного варианта, в котором частицы оказывают как весовое, так и тепловое воздействие, описывается кривой 4 с t = 4,45 (верхняя точка 4 на рис. 2.4).

Из сопоставления кривых 1–4 на рис. 2.5 можно оценить относительный вклад различных механизмов межфазного взаимодействия, определяющих закономерности подъема нагруженного термика. На начальном этапе движения примесь сосредоточена в центре облака, ее концентрация велика, что препятствует развитию восходящего конвективного движения и ведет к затягиванию стадии разгона запыленного термика (кривые 2, 4) по сравнению с чисто газовым облаком (кривая 1). На автомоГлава 2. Подъем термика и перенос примесей Рис. 2.5. Сравнение динамики подъема чисто газового термика (кривая 1) и запыленного термика при учете только весового (2) или теплового (3) воздействия примеси, а также обоих механизмов взаимодействия (кривая 4) Глава 2. Подъем термика и перенос примесей дельной стадии примесь, не обменивающаяся теплом с газом (кривая 2) значительно замедляет подъем, тогда как при наличии аккумулирующих тепло, но не обладающих массой частиц (кривая 3) развивается более быстрое восходящее течение. Указанные факторы влияют на скорость подъема противоположным образом: с одной стороны, примесь утяжеляет облако и замедляет его движение, с другой же стороны, запасая часть общего тепла и сообщая его газу в процессе подъема, примесь препятствует остыванию термика и способствует тем самым его более быстрому всплытию.

Суммарная скорость движения запыленного термика, таким образом, определяется противодействием этих факторов. При малой загрузке ( < ) они взаимно компенсируются, что и объясняет совпадение на рис. 2.4 точек 1–3 со сплошной кривой.

С увеличением загрузки ( > ) тепловое воздействие примеси становится недостаточным для компенсации ее весового влияния, облако замедляется (точки 1– 3 на рис. 2.4). В переходной же области (, точки 4 на рис. 2.4) преобладание весового или теплового механизма воздействия зависит от интенсивности конвекции, характеризуемой числом Грасгофа Gr. Если для Gr = 1650 скорость подъема запыленного термика на автомодельном участке оказывается выше, чем чисто газового (ср. наклоны кривых 1 и 4 на рис. 2.5, то при Gr = 460 наблюдается обратная картина: в перечисленных выше случаях поочередного «включения» механизмов межфазного взаимодействия автомодельная координата верхней кромки равна соответственно t = 3,55, 2,20, 4,05 и 3,20. Это означает, что с уменьшением вязкости среды (с ростом Gr) эффекты, связанные с передачей тепла от примеси несущему газу на автомодельной стадии движения облака проявляются более ярко.

Из приведенных выше результатов следует, что для не слишком больших загрузок динамика подъема запыленного термика в неоднородной сжимаемой атмосфере описывается зависимостью (2.30), полученной для чисто газовых облаков. Поэтому, в отсутствие экспериментальных данных о воздействии примеси на динамику движения, представляется возможным при помощи соотношений (2.30) определять турбулентные коэффициенты переноса и для термика с частицами. По известной полной Глава 2. Подъем термика и перенос примесей тепловой энергии облака I0, массе находящихся в нем частиц Mp и их относительной теплоемкости p можно определить энергию, запасенную газом Ig. Затем, зная начальную высоту облака H0 и воспользовавшись функциями и (см. (2.30)) можно найти величину числа Грасгофа, при котором автомодельная координата верхней кромки принимает наблюдающееся в экспериментах значение t 4,35. При таком подходе обеспечивается непрерывность перехода от запыленных к чисто газовым термикам при уменьшении суммарного количества примеси в плавучем облаке.

2.5. Структура всплывающего термика при различной загрузке примесью Перейдем теперь к анализу внутренней структуры термика на различных стадиях его эволюции. Представляет интерес рассмотреть сначала подъем чисто газового термика, после чего, варьируя загрузку облака, выявить особенности подъема запыленного термика. В качестве базового варианта примем следующий набор параметров, определяющих начальное состояние термика: I0 = 0,68, H0 = 1,57, 0 = 21. Безразмерные параметры, задающие начальное состояние атмосферы и свойства газа, были фиксированы согласно (2.19), высота тропопаузы принята равной HT = 10.

При размерных масштабах (2.20) такой термик соответствует облаку взрыва с полной энергией Q = 4,18 · 1015 Дж на высоте H0 = 1560 м, поднимающемуся в стандартной двухслойной атмосфере с тропопаузой на уровне 10 км.

Выбор турбулентных коэффициентов переноса осуществлялся в соответствии с описанной в разделе 2.4 процедурой: числа Прандтля и Шмидта полагались равными единице, величина числа Грасгофа составила Gr = 1350, что для указанного запаса тепла I0 дает число Рейнольдса Re = 44,5. Все стадии эволюции рассчитывались сквозным образом. Для отслеживания эволюции термика введена пассивно переносимая примесь, для концентрации Y которой решается конвективно-диффузионное уравнение вида (2.13) с теми же начальными и граничными условиями (2.14)–(2.17), что приведены выше для активной примеси. Максимальная начальная «плотность»

Глава 2. Подъем термика и перенос примесей примеси 0 полагалась равной единице.

Динамика подъема чисто газового термика представлена на рис. 2.6 в виде зависимости от t1/2 координаты верхней кромки zt (кривая 1), горизонтальной пунктирной линией обозначено положение тропопаузы. Штриховой кривой нанесена зависимость от t1/2 вертикальной координаты центра тяжести примеси При подъеме термика последовательно наблюдаются все основные стадии движения.

На начальной стадии (t 4,5) термик «разгоняется» и выходит на автомодельную стадию (4,5 t 15), которая характеризуется линейным участком зависимости zt (t1/2). Затем наступает стадия зависания (рассчитано первое колебание термика), облако «прошивает» тропопаузу, частично проникая в стратосферу.

Перейдем к рассмотрению структуры термика на основных стадиях эволюции.

На начальном этапе нагретый газ, всплывая, вовлекает в движение окружающий воздух, что приводит к формированию вихревого тороидального течения. Сферической вначале облако приобретает характерную грибовидную форму. Области максимальной температуры и концентрации примеси находятся на оси симметрии. Максимальная скорость подъема газа составляет Umax = 1,6 1,8. Структура термика в момент времени t = 1,6 представлена на рис. 2.7а, где показаны изолинии массовой концентрации примеси Y (слева) и избыточной температуры (справа). Для каждой переменной 10 линий уровня построены с равными интервалами между минимальным и максимальным значениями соответствующей функции: минимальная и максимальная концентрации на изолиниях Y составляют (Y )min = 1,83 · 102, (Y )max = 0,183, минимальная и максимальная избыточная температура равна min = 0,18, max = 1,8. В размерных единицах максимальный перегрев газа в термике при этом составляет 550 К. В правой части рисунка стрелками показано поле скорости, масштаб которой указан отрезком. Такая же схема представления структуры облака использована и ниже.

На автомодельной стадии подъема облако увеличивается в размерах, при этом Глава 2. Подъем термика и перенос примесей Рис. 2.6. Динамика подъема термика при различной загрузке дисперсной примесью: 1 — чисто газовый термик, 2 — умеренно запыленное облако, = 0,56, 3 — сильно запыленное облако, = 0,96. Сплошные кривые — координата верхней кромки zt, штриховые — координата центра тяжести примеси zc Глава 2. Подъем термика и перенос примесей Рис. 2.7. Эволюция чисто газового термика (I0 = 0,68, H0 = 1,57, 0 = 21): а — начальная стадия (t = 1,6), б — автомодельный подъем (t = 6,1), в — выход в стратосферу, зависание (t = 13,8), г — колебания облака (t = 35,4) Глава 2. Подъем термика и перенос примесей форма облака практически остается неизменной, что хорошо согласуется с результатами расчета автомодельного термика в атмосфере с переменной плотностью, полученными в Главе 1. Тороидальный вихрь, сопровождающий подъем термика, также увеличивается в размерах. Области максимальной температуры и концентрации примеси сходят к оси, что наблюдалось и в расчетах других авторов, причем положение максимума температуры не совпадает с центром завихренности (см., например, [76]).

Хотя газ в термике остывает, его переохлаждения не происходит и избыточная температура остается положительной. Структура течения в момент t = 6,1 представлена на рис. 2.7б, (Y )min = 3,33 · 103, (Y )max = 3,33 · 102, минимальная и максимальная избыточная температура равна min = 2,1 · 102, max = 0,22. В размерных единицах максимальный перегрев газа в облаке к этому моменту составляет 70 К.

Автомодельная стадия продолжается в течение временного интервала 4,5 t 15. При этом до момента t 10 тепловая энергия термика I0 сохраняется с точностью до 5%, а к концу автомодельной стадии I0 уменьшается на 30% по сравнению со своим начальным значением, что вызвано влиянием устойчивой стратификации атмосферы. Интенсивное подъемное течение газа, формирующееся у оси симметрии в течение первых двух стадий, на стадии зависания вызывает переохлаждение газа по отношению к окружающей атмосфере. Особенно интенсивно переохлаждение происходит при выходе термика в стратосферу, которая характеризуется более высоким параметром стратификации J2 > J1. Структура термика в момент t = 13,8 приведена на рис. 2.7в, (Y )min = 8,4 · 104, (Y )max = 8,4 · 103, min = 1,67 · 102, max = 2,62 · 102 (в размерных единицах отклонения температуры газа от окружающей заключены в пределах ±8 К). При дальнейшем подъеме происходит инверсия температурного поля, когда вся область у оси симметрии имеет температуру ниже окружающей, а на периферии находится более теплый газ.

Колебания зависающего облака сопровождаются периодической перестройкой вихревой структуры течения. Точке максимального подъема, достигаемой к моменту времени t = 28, соответствует система из двух противоположно закрученных вихрей Глава 2. Подъем термика и перенос примесей приблизительно равной интенсивности, нижний из которых закручен в том же направлении, что и вихрь, сопровождавший подъем термика до зависания, а верхний — в противоположном направлении. Далее переохлажденное облако начинает опускаться, моменту прохождения положения равновесия соответствует основной вихрь с опускным движением газа у оси. Структура термика на этой стадии представлена на рис. 2.7г, t = 35,4, (Y )min = 2,5·104, (Y )max = 2,5·103, min = 1,87·102, max = 1,65 · 102. Нижней точке вновь отвечает система двух противоположно закрученных вихрей. Облако перегрето по отношению к атмосфере, оно начинает снова подниматься и т. д. В целом стадия зависания характеризуется термогравитационными колебаниями термика с периодически меняющейся вихревой структурой.

Представляет интерес сравнить полученные результаты с данными расчетов других авторов, в которых использовались более сложные модели турбулентности. В работе [81] использовались уравнения, записанные в рейнольдсовых напряжениях, позволяющие учесть анизотропию турбулентности, вызванную влиянием закрутки в вихревом кольце. В работах [83,84] исследовались термики воздушных и приповерхностных взрывов, при этом использовалась алгебраическая модель и дифференциальная k модель турбулентности. Сравнение результатов с представленными выше расчетами позволяет заключить, что структура термика воспроизводится качественно одинаково при использовании всех моделей турбулентности: термик быстро принимает характерную «грибовидную» форму, поднимаясь затем в автомодельном режиме. Количественные отличия состоят в том, что термик, рассчитанный на основе дифференциальных моделей турбулентности, более «выгнут», а максимум температуры дальше отходит от оси симметрии, чем в случае постоянных коэффициентов переноса. Что касается динамики подъема, то приведенные в [83, 84] координаты верхней кромки от времени для разных моделей турбулентности практически совпадают между собой на начальном и автомодельном участках подъема. На стадии зависания дифференциальные модели дают несколько большую высоту подъема облака.

Следует иметь в виду, однако, что результаты расчетов по k модели турбулентГлава 2. Подъем термика и перенос примесей ности весьма чувствительны к выбору начальных распределений кинетической энергии турбулентности и скорости ее диссипации. В частности, представленные в [84] результаты показывают, что высота зависания приповерхностного термика весьма существенно меняется при изменении начального уровня турбулентности. Поэтому вопрос о наиболее адекватной модели термика представляется открытым.

Перейдем теперь к анализу структуры запыленного термика, несущего активную дисперсную примесь. При умеренных загрузках ( < ) процесс подъема запыленного облака аналогичен рассмотренному выше случаю чисто газового термика, несущего пассивную примесь. При большой начальной загрузке, как уже было видно из кривых на рис. 2.5, 2.4, подъем термика по сравнению с соответствующим газовым облаком замедляется. На начальной стадии влияние примеси может быть особенно сильным вследствие ее высокой концентрации. Поскольку температура и концентрация распределены неравномерно, даже при равенстве веса частиц и суммарной выталкивающей силы локальной компенсации сил не происходит, в запыленном термике имеется тяжелое центральное ядро, обладающее меньшей (при большой загрузке — отрицательной) плавучестью, и более легкая периферийная область. При загрузках, близких или превышающих критическую величину =, облако на начальной стадии может сильно деформироваться или даже делиться, когда более тяжелое ядро опускается на поверхность, а периферийная часть поднимается вверх.

Подобные процессы наблюдались при расчетах неизотермических дисперсных облаков, изучавшихся численно в работах [98–101]. Отметим, что в случае достаточно больших загрузок начальная стадия движения характеризуется малыми скоростями газа. В этом случае односкоростная и однотемпературная модель среды, принятая в настоящей работе, может оказаться применимой лишь для очень мелких частиц, скорость витания которых должна быть меньше характерной скорости конвекции.

Поэтому при > эволюция термика с частицами не изучалась, а при полученные результаты следует рассматривать как качественные.

Рассмотрим более подробно структуру запыленного термика в тех случаях, когда Глава 2. Подъем термика и перенос примесей влияние примеси только начинает сказываться ( = 0,56) и когда оно становится существенным ( = 0,96). Расчеты проводились при I0 = 0,68, H0 = 1,56, p = 1, т. е., для термиков с суммарным запасом тепла и высотой, совпадающими с базовым вариантом, использованным выше при исследовании чисто газового термика.

При загрузке примесью = 0,56 параметр, описывающий распределение тепла между фазами, равен = 0,225, доли тепла, запасенного газом и дисперсной фазой, составляют g = 0,54, p = 0,46, максимальная температура в центре облака равно p = 5,9. Динамика движения термика показана на рис. 2.6 кривыми 2: сплошная кривая показывает координату верхней кромки облака zt, штриховая — координату центра тяжести примеси, вычисляемую аналогично (2.33). На рис. 2.8 представлены изолинии избыточной температуры и линии равной плотности примеси p = Yp, построенные с равными интервалами между минимальным и максимальным значениями соответствующей функции, а также поле скорости.

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

На рис. 2.8а структура термика показана в момент t = 1,6 (максимальные и минимальные значения функций на изолиниях max = 1,9, min = 0,19, p,max = 0,65, p,min = 6,5 · 102). В отличие от случая пассивной примеси или малой загрузки термика, в центральной части конвективное движение газа практически отсутствует, газ обтекает центр облака.

По мере того, как концентрация примеси в центре термика уменьшается, вихревое течение охватывает всю приосевую область. Пространственные распределения избыточной температуры и плотности примеси некоторое время сохраняют еще отличия, порожденные на начальном этапе. Так, при t = 6,1 (рис. 2.8б, max = 0,85, min = 8,5 · 102, p,max = 0,17, p,min = 1,7 · 102) облако вытянуто в вертикальном направлении. В дальнейшем такие различия постепенно исчезают, термик, соГлава 2. Подъем термика и перенос примесей провождаемый тороидальным вихрем, поднимается вверх — см. рис. 2.8в (t = 13,8, max = 0,12, min = 7,7 · 103, p,max = 3,7 · 102, p,min = 3,7 · 103). На стадии зависания концентрация частиц столь мала, что не оказывает влияния на термогравитационные колебания термика около уровня стабилизации и его температурный распад. На рис. 2.8г приведена структура течения в момент t = 26 (max = 8,8 · 103, min = 1,4 · 102, p,max = 4,8 · 103, p,min = 4,8 · 104). Термик находится вблизи точки зависания, у оси сформировалась переохлажденная область, над основным вихрем возник вторичный, закрученный в противоположную сторону. Дальнейшая эволюция термика аналогична рассмотренному выше чисто газовому термику.

Указанные особенности проявляются еще сильнее при увеличении загрузки термика. Структура облака, рассчитанная при I0 = 0,68, H0 = 1,56, p = 1, = 0,96, = 0,32 представлена на рис. 2.9а, зависимости от времени координаты верхней кромки термика zt и центра тяжести примеси zc нанесены на рис. 2.6 сплошной и штриховой линиями 3 соответственно. В этом случае загрузка близка к критической = 1, вес частиц в термика велик настолько, что его центральное ядро опускается на поверхность — см. рис. 2.9а (t = 3,0, max = 2,0, min = 0,2, p,max = 0,83, p,min = 8,3 · 102. При t = 6,9 облако вытянуто по вертикали, максимальная концентрация достигается вблизи поверхности (рис. 2.9б, max = 0,62, min = 6,3 · 102, p,max = 0,43, p,min = 4,3 · 102). Лишь после осаждения некоторой доли примеси вихревое тороидальное течение начинает увлекать частицы вверх, термик принимает грибовидную форму, а примесь из приповерхностной области выносится на большие высоты — см. рис. 2.9в (t = 13,8, max = 0,18, min = 1,8 · 102, p,max = 0,17, p,min = 1,7 · 102) и рис. 2.9г (t = 35,4, max = 2,3 · 102, min = 7,5 · 105, p,max = 2,5 · 102, p,min = 2,5 · 103).

Таким образом, наличие большого количества активной примеси приводит к замедлению подъема запыленного термика и к изменению его внутренней структуры.

Наиболее сильно эти отличия проявляются в том случае, когда суммарный вес примеси превышает 40% от действующей на газ выталкивающей силы.

Глава 2. Подъем термика и перенос примесей Рис. 2.8. Эволюция запыленного термика при умеренной загрузке = 0,56: а — t = 1,6, б — t = 6,1, в — t = 13,8, г — t = 35, Глава 2. Подъем термика и перенос примесей Рис. 2.9. Эволюция запыленного термика при большой загрузке = 0,96: а — t = 3,0, б — t = 6,1, в — t = 13,8, г — t = 35, Глава 2. Подъем термика и перенос примесей Изложенные выше результаты получены для термиков, которые образуются при воздушных взрывах и в начальный момент не касаются подстилающей поверхности.

При наземных взрывах источник тепловыделения находится на поверхности, а термик имеет форму, близкую к полусфере [108]. Как отмечается в [92], именно такие термики могут эффективно захватывать и выносить пыль на большие высоты, поскольку в этом случае в объем облака попадает большое количество примеси.

Внутренняя структура приповерхностного термика при отрыве от подстилающей поверхности рассматривались при I0 = 0,68 для различной загрузки облака примесью. Начальные распределения величин при z 0 задавались в виде (2.6) c H0 = 0, так что в момент t = 0 термик имел форму полусферы, опирающейся на подстилающую поверхность.

Нагретый газ, поднимаясь вверх вблизи оси симметрии, формирует в атмосфере вихревое течение, при этом вдоль подстилающей поверхности возникает направленный к центру горизонтальный поток холодного воздуха, который «отсекает» область максимальной температуры от земной поверхности. Структура приповерхностного чисто газового термика, переносящего пассивную примесь, представлена на рис. 2.10а в момент t = 2,6 (max = 1,6, min = 0,16, p,max = 0,19, p,min = 1,9 · 102, при задании начальных условий концентрация пассивной примеси Y в центре облака полагалась равной единице). В дальнейшем облако принимает характерную грибовидную форму — см. рис. 2.10б (t = 8,5, max = 0,19, min = 1,9 · 102, p,max = 3,4 · 102, p,min = 3,4 · 103), после достаточного удаления от поверхности дальнейшая эволюция термика практически не отличается от термика воздушного взрыва. Экспериментально процесс отделения термиков от поверхности изучался, например, в [34], а расчеты подобного класса течений, проведенные в работах [74,84], дали качественно близкие результаты. По такому же механизму происходит отрыв от поверхности и слабо запыленных ( < ) термиков.

Глава 2. Подъем термика и перенос примесей Рис. 2.10. Начальная стадия подъема приповерхностного термика: а, б — газовый термик, в, г — запыленный термик ( = 0,7) Глава 2. Подъем термика и перенос примесей При большой загрузке, как и в случае термиков воздушных взрывов, начальная стадия подъема запыленного термика подвержена влиянию тяжелого ядра примеси.

На рис. 2.10в, г представлена структура термика при I0 = 0,68, = 0,7, p = 1.

В начале процесса тяжелое ядро, обладающее отрицательной плавучестью, начинает «растекаться» по поверхности z = 0, тогда как на периферии, где концентрация частиц меньше, возникает вихревое течение газа — рис. 2.10в (t = 2,6, max = 2,0, min = 0,2, p,max = 1,13, p,min = 0,11). По мере уменьшения концентрации примеси растекание по поверхности прекращается и восходящее течение охватывает всю приосевую область — рис. 2.10г (t = 8,5, max = 0,7, min = 7,2·102, p,max = 0,54, p,min = 5,4 · 102). Как и в случае термика воздушного взрыва, при большой загрузке процесс подъема запыленного термика характеризуется большей (по сравнению с чисто газовым облаком) продолжительностью начальной стадии и меньшей скоростью подъема на автомодельном участке движения.

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

Для описания вертикального распределения примеси в атмосфере будем использовать приведенную вертикальную координату где Bg = g B0 — суммарный запас плавучести газа в дисперсном облаке. Координата z имеет смысл отношения расстояния между данной точкой и высотой взрыва к характерному расстоянию Ig ts, проходимому термиком за время ts = J представляющему собой характерную высоту зависания термика. Положение троГлава 2. Подъем термика и перенос примесей попаузы по отношению к точке образования термика характеризуется приведенной высотой тропопаузы zT = (HT H0 )/Ig t1/2. Введенные координаты не зависят от линейного масштаба L и потому их использование предпочтительнее исходных безразмерных координат.

В качестве характеристик загрязнения верхних слоев атмосферы будем рассматz ривать безразмерную линейную плотность примеси (), определяемую интегрированием создаваемого термиком распределения концентрации примеси по радиальной координате, а также долю примеси, выносимую за тропопаузу:

Функция () нормирована так, что интеграл от нее по высоте z равен единице. Использование в определении z величины Ig позволяет единым образом рассматривать как запыленный, так и чисто газовый термик, переносящий пассивную примесь. Аналогичные параметры применялись для описания выноса примеси в стратосферу, например, в [16, 50, 102].

Рассмотрим сначала случай пассивно переносимой примеси (как и ранее, при моделировании такой добавки в определяющих уравнениях полагалось p = 0, а концентрация пассивной субстанции Y находилась решением уравнения конвективной диффузии (2.4)). Расчеты проводились в диапазоне 1,3 zT 5,8, что достигалось варьированием величин 1,5 H0 6,0, 9 HT 16, 0,13 I0 2,71 (для указанного выше конкретного примера (2.19), (2.20) это соответствует облакам, возникающим при взрывах мощностью W = 0,2 4 (Q = 8,4 · 108 1,7 · 1010 МДж) при высоте тропопаузы HT = 9–16 км). Полученная в расчетах зависимость относительного количества вещества, выносимого в стратосферу p, от параметра zT представлена на рис. 2.11а. Точки 1–5 отвечают значениям I0 = 0,13, 0,34, 0,68, 1,35, 2,71. При малых значениях тепловой энергии расчетные точки достаточно хорошо описываются единой кривой, нанесенной сплошной линией. С увеличением параметра I0 при фиксированном zT наблюдаются отклонения от полученной зависимости — Глава 2. Подъем термика и перенос примесей доля выносимой примеси несколько уменьшается (точки 4, 5).

Представляет интерес сравнить данные по выносу пассивной примеси, полученные на основе модели сжимаемой среды, с имеющимися в литературе. Штриховой линией на рис. 2.11а нанесены результаты, полученные в работе [50] путем численного интегрирования уравнений несжимаемой среды в приближении Буссинеска.

Штрихпунктирная линия соответствует аналитическому решению задачи о подъеме термика, найденному в работах [16, 41] (в указанных работах для величины выброса получены верхняя и нижняя оценки, которые оказались весьма близкими, кривая на рис. 2.11а отвечает максимальной оценке p). Количественное отличие результатов, полученных в данной работе (сплошная кривая) и в [16, 41, 50], объясняется тем, что используемая здесь модель среды учитывает уменьшение плотности атмосферы с высотой. При выходе термика в разреженные слои атмосферы снижается действующая на газ выталкивающая сила, а следовательно, уменьшается высота зависания облака и доля выносимой примеси.

Обратимся теперь к вопросу о влиянии начальной загрузки термика на вынос в стратосферу активной примеси. Данные расчетов, проводившихся при I 0 = 0,68, H0 = 1,56, 9 HT 16 представлены на рис. 2.11б. Отношение теплоемкостей фаз p, как и ранее, полагалось равным единице. Степень загрузки облака варьировалась в широких пределах = 0,07, 0,2, 0,5, 0,56, 0,7, 0,95 (точки 1–6). Для сравнения на рис. 2.11б сплошной кривой нанесена зависимость p (T ), полученная выше для пассивной примеси (см. рис. 2.11а). При < = 0,4 термик с частицами поднимается подобно чисто газовому, имеющему запас тепла I g < I0. Поэтому вынос активной примеси описывается той же зависимостью p (T ) — точки 1, 2 попаz дают на кривую рис. 2.11а. При > начинает сказываться воздействие примеси, облако поднимается медленнее, что приводит к уменьшению выброса в стратосферу (точки 3–6).

Более детальную информацию о распределении примеси по высоте, создаваемом в результате подъема термика, можно получить, сравнивая между собой проГлава 2. Подъем термика и перенос примесей Рис. 2.11. Доля примеси, выносимой термиком в стратосферу: а — пассивная примесь, штриховая кривая — численное решение [50], штрих-пунктирная кривая — аналитическое решение [16], б — активная примесь Глава 2. Подъем термика и перенос примесей фили () для разных степеней загрузки, при этом целесообразно сопоставлять результаты, отвечающие на рис. 2.11б одному значению абсциссы zT. На рис. 2. построены функции (), полученные при zT = 3, = 0, 0,2, 0,5, 0,7, 0,85 (кривые 1–5), пунктирной линией обозначена тропопауза. Результаты для пассивной примеси (кривая 1) и активной в случае малой загрузки (кривые 2, 3) различаются незначительно. При большой загрузке (кривые 4, 5) максимум линейной плотности примеси смещается вниз, в сторону меньших z, при этом в стратосферу попадает меньшее количество примеси.

Чтобы исследовать применимость полученных выше результатов к приповерхностным термикам, была проведена серия расчетов с I0 = 0,68 (что соответствует Q = 4,18 · 1015 Дж) при различных загрузках облака. Расчеты показали, что зависимости p (T ), получающиеся при малой загрузке для приповерхностных и возz душных термиков, близки между собой. Так, при высотах тропопаузы HT = 9, и 12 (при этом zT = 3,3, 3,7 и 4,4) вынос пассивной примеси приповерхностным термиком равен p = 0,62, 0,46 и 0,22 соответственно. На рис. 2.11а эти результаты нанесены точками 6, которые попадают на сплошную кривую, полученную для термиков воздушных взрывов. Аналогичное сопоставление выноса активной примеси (p = 1) приповерхностными и воздушными термиками проводилось для = 0, и 0,7. В первом случае приповерхностный термик выносит в стратосферу такое же количество примеси, что и соответствующий ему при том же zT воздушный (см. точки 7 на рис. 2.11б). При большой загрузке вынос примеси за тропопаузу уменьшается для обоих термиков: точки 8 ( = 0,7) для приповерхностного и точки 5 для соответствующего ему воздушного термика лежат ниже сплошной кривой.

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

Глава 2. Подъем термика и перенос примесей Рис. 2.12. Распределения линейной плотности примеси, создаваемые в атмосфере при подъеме термика: кривые 1–5 соответствуют = 0, 0,2, 0,5, 0,7, 0,85. Пунктирная линия zT = 3 обозначает положение тропопаузы Глава 2. Подъем термика и перенос примесей Приведенные на рис. 2.11 и 2.12 результаты получены при неизменном значении параметра p = 1. Чтобы установить степень «чувствительности» данных по выносу аэрозоля к выбору конкретного значения p, на рис. 2.13а проводится сравнение результатов для фиксированных I0 = 0,68, H0 = 1,56, HT = 10, полученных при p = 0,25, 1,0, 2,0 (точки 1–3) и различных начальных загрузках. По оси абсцисс отложена величина 0 = Mp /2I0 = g, которая, как следует из соотношения F+ = 2Ig, имеет смысл отношения веса примеси к максимально возможной при данном количестве выделившегося тепла подъемной силе, т. е., силе, действующей на чисто газовый термик с энергией I0. Важный вывод, следующий из рис. 2.13а, состоит в том, что при небольших загрузках вынос примеси практически не зависит от величины p. По-видимому, происходит взаимная компенсация двух механизмов воздействия частиц на газ: при одном и том же количестве примеси термик, соответствующий меньшему значению p, обладает большей подъемной силой, но эффективная теплоемкость среды при этом меньше, что ведет к более быстрому его остыванию.

Напротив, при большем p подъемная сила меньше, но газ остывает медленнее. При дальнейшем увеличении загрузки весовой механизм начинает преобладать над тепловым, поэтому кривые расходятся.

Для сравнения штриховой линией на рис. 2.13а нанесена зависимость доли пассивной примеси (расчет чисто газового термика). В этом случае под 0 следует понимать такую величину загрузки, что при p = 1 «газовая» часть Ig тепловой энергии запыленного термика с загрузкой 0 совпадает с энергией рассматриваемого газового термика. Видно, что все полученные результаты при малых 0 очень близки, поэтому расчет выноса активных частиц в стратосферу для таких загрузок можно производить на основании модели пассивной примеси, подправив соответствующим образом тепловую энергию термика.

Воспользовавшись результатами рис. 2.13а, можно определить зависимость полной массы частиц, попадающих в стратосферу s = 0 p от степени загрузки облака примесью 0. Соответствующие данные приведены на рис. 2.13б; точки 1– Глава 2. Подъем термика и перенос примесей и штриховая кривая имеют тот же смысл, что и на рис. 2.13а. Видно, что увеличение загрузки термика примесью ведет к возрастанию массы вынесенной примеси лишь до определенного предела. Дальнейшее увеличение загрузки 0 вызывает замедление подъема термика и, как следствие, уменьшение величины s. Максимально возможный вынос зависит от параметра p : чем меньшая часть тепла идет на нагрев частиц при образовании термика, тем большее количество дисперсной примеси может попасть в стратосферу.

В заключение проведем сравнение полученных выше результатов с экспериментальными данными. В [92] обсуждаются вопросы выноса пыли в стратосферу термиком, образующимся в результате приземного ядерного взрыва. При этом указывается, что наиболее вероятный диапазон загрузки термика на единицу мощности взрыва составляет (0,2–0,5)·109 кг/Мт. Приняв p = 1, получим, что величина при этом изменяется в пределах (3–7,5)·102, тогда = 0,05–0,18 (см. рис. 2.1). Таким образом, загрузка термика примерно на порядок меньше критической, что согласуется с приведенными в [92] оценками. Термик, образующийся при взрыве с тепловыделением Q = 4,18 · 109 МДж и касающийся в начальный момент подстилающей поверхности, при высоте тропопаузы 12 км, согласно нашим расчетам, выносит в стратосферу 20–25% содержащейся в нем пыли. Это примерно вдвое меньше оценок [92], сделанных в предположении о равномерном распределении пыли по объему зависшего облака.

2.7. Метод расчета течений сжимаемого газа Численное решение уравнений вязкого теплопроводного сжимаемого газа (2.1)–(2.4) осуществлялось при помощи неявной трехслойной экстраполяционной схемы покоординатного расщепления [109,110], модифицированной с целью улучшения консервативных свойств. Эта схема хорошо зарекомендовала себя при решении различных задач теории горения, механики дисперсных сред, расчете конвективных течений как в замкнутых областях, так и в открытой атмосфере [98–102, 111–113]. В ее основе Глава 2. Подъем термика и перенос примесей Рис. 2.13. Влияние отношения теплоемкостей газовой и дисперсной фаз p на долю p и полную массу s = 0 p выносимой в стратосферу примеси. Абсцисса 0 = g равна отношению веса примеси к максимально возможной при данной мощности взрыва подъемной силе, сплошные кривые соответствуют чисто газовому термику Глава 2. Подъем термика и перенос примесей лежит неявная аппроксимация пространственных производных с последующей линеаризацией и факторизацией разностного оператора. Близкая по реализации схема для расчета сжимаемых течений была предложена в [114].

При решении задачи о подъеме термика повышенные требования предъявляются к консервативным свойствам схемы, которая должна обеспечивать выполнение баланса тепловой энергии облака. Исходная разностная схема [109, 110] является неконсервативной, что может приводить к возникновению «источников тепла», вычислительного происхождения. Для улучшения консервативных свойств разностной схемы была осуществлена модификация, близкая по сути к подходу, используемому в консервативных при установлении схемах [115]. В этих схемах правая часть разностных уравнений, вычисляемая по известным с предыдущего слоя значениям функций, выражается через разностный аналог дивергенции вектора потоков, а левая (неявная) часть остается в недивергентной форме. Полученные таким образом схемы консервативны при установлении, т. к. при выходе на стационар выполняются разностные аналоги законов сохранения. Поскольку недивергентный оператор действует не на сами функции, а на их приращения между соседними временными слоями, можно ожидать хорошие консервативные свойства схемы и при расчете нестационарных течений, по крайней мере при достаточно малых шагах по времени.

Продемонстрируем основные особенности схемы [109, 110] и произведенные для улучшения ее консервативных свойств модификации на примере расчета чисто газового термика. Уравнения (2.1)–(2.4) могут быть записаны в операторном виде как где f — вектор искомых функций, — дифференциальный оператор, равный сумме операторов r, z, действующих по соответствующим координатным направлениям, и rz, содержащего смешанные производные, F — вектор правой части, содержащий компоненты массовых сил и объемные источники тепла. Численная аппроксимация производной по времени в уравнении (2.34) осуществляется по схеме первого порядка точности, а линеаризация конвективных членов производится относительно Глава 2. Подъем термика и перенос примесей значений искомых функций, полученных экстраполяцией относительно текущего и предыдущего временного слоя. В качестве новой переменной используется разность приращений функции на соседних временных слоях где — экстраполированное значение искомой функции. После аппроксимации вреf менной производной по схеме первого порядка «назад по времени» и линеаризации оператора относительно получим [I + t] = tR, где I — единичный оператор, = ( — линеаризованный оператор. Пренебреf) гая смешанными производными, факторизуем оператор в левой части: I + t I + t [r + z ] [I + tr ] [I + tz ]. В результате получим систему одномерных уравнений, которая при аппроксимации пространственных производных центральными разностями имеет блочную трехдиагональную матрицу, для обращения которой используются продольно-поперечные векторные прогонки:

(здесь — вспомогательная функция, под понимаются соответствующие разностные операторы). При решении уравнения (2.36) от системы уравнений отщепляется уравнение для v, так что векторная прогонка используется для совместного решения разностных уравнений для u, P и T вдоль сеточных линий в радиальном направлении. Соответственно, при решении (2.37) отщепляется уравнение для u, а методом векторной прогонки решаются уравнения для v, P и T вдоль вертикальных сеточных линий. Таким образом, размер блоков в каждой точке сетки при векторных прогонках составляет 3 3. Отметим, что при консервативной записи решаемых совместно уравнений для компонент скорости, температуры и давления в двумерном случае возникает необходимость обращения матриц размера 4 4, а в трехмерном Глава 2. Подъем термика и перенос примесей размер матриц возрастает до 5 5, что приводит к значительному увеличению вычислительных затрат.

Модификация численной схемы [109, 110], позволяющая значительно улучшить ее консервативные свойства, состоит, во-первых, в переходе к сетке второго рода и, во-вторых, в применении дивергентной записи оператора в правой части R уравнений (2.35), (2.36). На сетке второго рода значения искомых переменных определены в центрах ячеек, а для аппроксимации граничных условий вводятся фиктивные ячейки, расположенные симметрично к приграничным ячейкам относительно соответствующих границ области. Для перехода к дивергентным аппроксимациям в правой (явной) части разностных уравнений (2.35), (2.36) заметим, что R представляет собой невязку исходных уравнений, в которые вместо точного значения искомых функций на верхнем слое f n+1 используется экстраполированное по двум предыдущим слоям значение. Это позволяет преобразовать правую часть к дивергентному виду при помощи уравнения неразрывности. Суть вносимых изменений продемонстрируем на примере уравнения притока тепла (2.12), которое для чисто газового термика в безразмерных переменных имеет вид Соответствующая компонента вектора правой части в исходной схеме [109,110] имеет вид Переходя к дивергентной записи конвективных членов в уравнении притока тепла, преобразуем это выражение к виду Из-за наличия последнего члена в квадратных скобках выражение для RT не является полностью дивергентным, однако для дозвуковых течений следует ожидать, что этот член мал. Аналогичным образом могут быть преобразованы и оставшиеся Глава 2. Подъем термика и перенос примесей компоненты вектора R, причем в правых частях уравнения импульса Ru и Rv к дивергентному виду возможно привести все члены, а в выражении для RP, как и выше, к дивергентной форме записи приводятся конвективные и диффузионные члены.

Для иллюстрации эффективности модификации приведем данные о сохранении тепловой энергии (2.22) газового термика. При недивергентной записи правой части (схема [109, 110]) накопление ошибок аппроксимации приводит к сильному (в несколько раз за время подъема термика) росту интеграла тепла Q. Некоторого улучшения выполнения балансов можно добиться путем искусственного зануления возмущений температуры вдали от термика. Типичное поведение интеграла тепла Q(t), отнесенного к начальному значению Q(0), представлено на рис. 2.14. Кривая 1 отвечает расчету со срезанием возмущений температуры. К моменту наивысшего подъема термика (t 25) дисбаланс энергии достигает порядка 50%. Переход к модифицированной схеме позволяет значительно улучшить выполнение баланса тепла, соответствующая зависимость интеграла энергии от времени показана на рис. 2. кривой 2: дисбаланс не превышает 5%, при этом процедура срезания возмущений не применялась.

Основная часть расчетов подъема газовых и запыленных термиков выполнена на ортогональных сетках, содержащих 40 50 ячеек в радиальном и вертикальном направлениях соответственно. Сетки сгущались у оси симметрии и у подстилающей поверхности. Поскольку на начальной стадии термик локализован в относительно малой области у земной поверхности, а по мере подъема его размеры значительно возрастают, в ходе расчета периодически проводилось расширение расчетной области с соответствующим перестроением сетки и квадратичной переинтерполяцией решения со старой сетки на новую. Расширение области производилось когда верхняя кромка термика достигала 2/3 высоты расчетной области. Типичный расчет начинался в области размером 5 6 (в безразмерных единицах, где масштаб длины приблизительно равен начальному диаметру термика), тогда как на стадии зависания область имела размер 20 30. Расчеты велись с переменным шагом по времени. ЕсГлава 2. Подъем термика и перенос примесей Рис. 2.14. Выполнение баланса тепла при использовании исходной схемы [109,110] (кривая 1) и модифицированной схемы с дивергентными аппроксимациями правой части (кривая 2) Глава 2. Подъем термика и перенос примесей ли на начальной стадии число Куранта, рассчитанное по характерной скорости звука Ku = t/(h·Ma), составляло 0,25–0,5 (здесь h — минимальный пространственный шаг сетки), то к концу расчета оно возрастало до Ku = 3–4.

Для оценки сходимости решения по сетке был проведен ряд тестовых расчетов на сетках с удвоенным числом узлов 80100. Отличия в динамике подъема (зависимости координаты верхней кромки от времени) составили порядка 5%, что подтверждает достаточность использования указанных сеток. Программа была реализована на языке программирования FORTRAN-77.

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

Модель учитывает весовое и тепловое влияние примеси на течение несущего газа.

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

Проведены сквозные расчеты всех стадий эволюции воздушных и приповерхностных термиков с различной загрузкой, определены временные границы каждой стадии и показана внутренняя структура термика и генерируемого им вихревого течения. Показано, что примесь оказывает двоякое воздействие на динамику движения облака: с одной стороны, частицы утяжеляют газ и замедляют подъем, с другой же стороны, тепло, запасенное примесью и отдаваемое затем газу, замедляет остывание Глава 2. Подъем термика и перенос примесей облака, способствуя развитию более интенсивного конвективного течения.

При малых начальных загрузках (если общий вес частиц не превышает 40% от начальной подъемной силы, действующей на термик) динамика и структура запыленного термика аналогичны случаю чисто газового облака с запасом тепла, равным доле полной энергии, идущей на нагрев газовой фазы. С ростом количества находящейся в термике примеси сверх указанной границы начинает сказываться влияние дисперсной фазы на структуру термика, наиболее сильно проявляющееся на начальной стадии процесса. Наличие примеси приводит к ослаблению, а при значительной загрузке — и к подавлению восходящего конвективного течения газа в ядре термика. При больших загрузках происходит сильная деформация термика и затягивается начальная стадия подъема.

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

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

Глава Образование и горение газовых облаков при аварийных выбросах в атмосферу. Физическая теория 3.1. Аварийные выбросы, их особенности и методы изучения Одним из основных источников опасности в современной химической промышленности является аварийный выброс пожаровзрывоопасных веществ в атмосферу [116–118]. Типичными причинами такого выброса служат разгерметизация или разрушение производственных резервуаров, трубопроводов, средств транспортировки и хранения. В случае воспламенения образующаяся топливовоздушная смесь может сгорать как во взрывном режиме (в этом случае главную опасность представляют генерируемые ударные волны), так и в режиме медленного горения, когда наиболее существенным поражающим фактором является тепловое излучение пламени. Реализация того или иного режима сгорания облака зависит от реакционной способности топлива, внешних условий (например, степени загроможденности пространства) и от условий зажигания (задержки воспламенения, энерговклада при инициировании реакции и т. п.).

В настоящее время разработан и применяется ряд инженерных методик оценки опасностей химических производств (см., например, [119–122]). Эти методики Глава 3. Теория выбросов основаны на схематизации физических процессов, протекающих при том или ином сценарии развития аварии, и описании каждого процесса простыми математическими моделями с привлечением эмпирической информации для выбора возникающих констант. Применительно к описанию последствий разгерметизации резервуаров со сжатым газом обычно рассматриваются два предельных случая — непрерывное истечение из небольшого отверстия в оболочке резервуара, либо мгновенное высвобождение всей массы газа при полном разрушении сосуда. Результатом непрерывного истечения, как правило, является стационарная струя газа, тогда как мгновенный выброс приводит к формированию компактного облака. Соответственно, зажигание выброшенного газа от расположенного вблизи места аварии источника приводит в первом случае к возникновению горящего факела, а во втором — к формированию огненного шара.

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

Обзор современных достижений в исследованиях этих объектов можно найти в литературе [116–122].

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

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

Только после того, как газ разбавляется воздухом до низких концентраций, примесь может считаться нейтрально плавучей (пассивной).

Струйные течения, возникающие при непрерывном выбросе газов в атмосферу изучены весьма подробно, что обусловлено широкой распространенностью этого класса потоков в природе, технологических установках, двигателях и т. д. Имеются многочисленные измерения скоростных и концентрационных полей в инерционных и плавучих струях, получены корреляционные соотношения для скорости вовлечения, изучены свойства турбулентности на различных участках струи (см., например, [9, 12, 123–129]). Большой объем эмпирических данных накоплен и по реагирующим струям, включая горящие факелы как в неподвижной атмосфере, так и в спутном и поперечном сносящем потоках [130–141]. В результате были разработаны модели, позволяющие описать геометрию струи, оценить количество газа, концентрация которого находится между нижним и верхним концентрационными пределами воспламенения, построить траекторию факела в сносящем поле ветра, найти длину турбулентного пламени, оценить тепловые потоки, создаваемые горящим факелом.

Если расход газа в источнике меняется со временем, в атмосфере развивается нестационарное течение. Одним из возможных видов таких истечений является струя, развивающаяся (стартующая) при резком включении источника, расход газа через который затем поддерживается постоянным. Как показывают экспериментальные исследования (см., например, [8,142–144]) и численные расчеты [145], такая струя состоит из головной части, где течение носит вихревой характер, и следующей за ней квазистационарной части, распределения скорости и концентрации в которой практически совпадают с соответствующими распределениями в установившейся наГлава 3. Теория выбросов порной струе. Установление стационарной струи, таким образом, происходит сразу за движущимся передним фронтом, координата которого в осесимметричном случае растет со временем по корневому закону t1/2 [6, 8]. Аналогичный характер имеет и развитие положительно плавучей струи.

Весьма важен в практических приложениях случай, когда действие источника ограничено по времени. Именно такой характер носит истечение из резервуаров конечного объема, при котором время истечения ограничено выбросом всей массы заключенного в резервуаре газа. Как правило, расход газа через выходное отверстие источника при таких истечениях зависит от времени, что обусловлено понижением внутреннего давления как за счет уменьшения массы газа в резервуаре, так и вследствие резкого падения температуры, вызванной расширением газа [119, 146–148]. В то время как процессы истечения газа из разгерметизированного резервуара (расчет расхода газа как функции времени) изучены достаточно подробно [149], лишь в отдельных работах предпринимались попытки описания распределений концентрации газа в атмосфере, создаваемых источниками конечного времени действия [150–152].

Если масса газа выбрасывается в атмосферу достаточно быстро (в пределе — мгновенно), дальнейшая эволюция облака во многом определяется его интегральными параметрами, такими как суммарный начальный импульс и начальный запас плавучести. В случае чисто инерционных направленных выбросов возникающее конусообразное облако распространяется, смешиваясь с окружающим воздухом, расширяясь и замедляясь. В осесимметричном случае координата его переднего фронта растет со временем как t1/4 [6, 153]. Наоборот, если выделение газа происходит так, что начальный импульс незначителен (например, как в экспериментах с мыльными пузырями), эволюция облака определяется только плавучестью газа. Такое облако представляет собой концентрационно-тепловой термик, рассмотренный в Главе 1.

Как правило, возникающее при мгновенном выбросе облако в течение достаточно длительного времени переобогащено топливом, поэтому при немедленном зажиГлава 3. Теория выбросов гании оно может гореть в диффузионном режиме, поднимаясь под действием сил плавучести и образуя огненный шар [116, 117]. Основным отличием огненного шара от тепловых термиков, рассмотренных в Главах 1–2, является наличие в облаке химических реакций, сопровождающихся непрерывным выделением тепла и поддерживающих высокую температуру на протяжении всего времени жизни огненного шара.

Термик же, образующийся в результате мгновенного выделения тепла, быстро остывает, смешиваясь с холодным внешним воздухом. Горение огненного шара продолжается до тех пор, пока не будет израсходовано все топливо, при этом огненный шар образует грибовидное облако, «ножка» которого — это сильное восходящее конвективное течение. В типичном для аварий с участием углеводородных топлив огненном шаре масса топлива порядка 100 тонн сгорает в атмосфере за время порядка 10–20 с, выделяя энергию 5 · 1012 Дж, причем около 20–30% этой энергии выделяется в виде теплового излучения [141]. Таким образом, для огненных шаров характерно наличие мощного, но кратковременного импульса излучения.

Источником огненного шара, как правило, служит выброс газа или двухфазной горючей смеси в атмосферу, хотя при определенных условиях огненные шары могут образовываться и при зажигании паров жидкого криогенного топлива, пролитого на поверхность [154–156]. Важным случаем выброса, ведущего к возникновению парожидкостного облака, является так называемый взрыв расширяющихся паров вскипающей жидкости, или BLEVE (Boiling Liquid Expanding Vapor Explosion). Это явление характерно для аварий, в которых происходит разрушение оболочки резервуара со сжиженным газом. Резкий сброс давления приводит к объемному вскипанию жидкости, в результате чего часть ее переходит в пар, увлекая оставшуюся жидкость в виде мелкодисперсного аэрозоля. BLEVE считается наиболее опасным исходом аварии со сжиженными углеводородами, поэтому его изучение в последнее время вызывает повышенный интерес [157–163].

Экспериментальное изучение огненных шаров началось в середине 60-х годов при исследовании последствий взрывов ракетных топлив [164–167], содерГлава 3. Теория выбросов жащих в качестве компонент керосин и кислород, а также водород и кислород.

Важность исследования огненных шаров с точки зрения обеспечения промышленной безопасности стимулировала дальнейшие усилия по изучению этого явления [21, 155, 168–184]. Диапазон экспериментально изученных огненных шаров простирается от лабораторных опытов с массой топлива порядка нескольких грамм или килограмм [21, 155, 168–172, 179, 180], до крупномасштабных полигонных испытаний, в которых исследовались огненные шары с массой топлива, исчисляемой единицами и десятками тонн [158,183]. Как правило, в качестве топлива использовались углеводороды, широко применяемые в качестве топлив или сырья и, следовательно, транспортируемые и хранимые в больших количествах — метан, пропан, бутан, пентан [155, 158, 169–172, 179]. Горение облаков водорода в условиях открытой атмосферы исследовано в [21, 167, 168, 173, 174]. Горение и детонация двухфазных облаков керосина и дизельного топлива исследовались в крупномасштабных экспериментах [182–184], а распространение пламени по двухфазному облаку изучалось в работе [185].

В подавляющем большинстве имеющихся на настоящий момент экспериментальных работ основное внимание уделялось изучению интегральных характеристик огненных шаров — максимального диаметра горящего облака DF B, времени жизни огненного шара tF B, высоты его подъема за время горения HF B. Перечисленные величины, как правило, описываются корреляционными соотношениями, имеющими вид степенных функций массы топлива M0 в огненном шаре:

В таблице 3.1 представлены данные по характеристикам огненных шаров, приведенные к единицам CИ, использовались как оригинальные экспериментальные работы, так и обзоры [120, 121, 178, 186–188].

Сравнительный анализ моделей огненных шаров, представленных различными авторами, показывает, что, несмотря на большие различия в массах топлива, типах горючего и условиях проведения экспериментов, практически все имеющиеся корГлава 3. Теория выбросов Таблица 3.1. Экспериментальные данные по характеристикам огненных шаров Глава 3. Теория выбросов реляции дают зависимость максимального диаметра огненного шара D F B от массы топлива M0 в виде степенной функции (3.1) с показателем b1, близким к 1/3. Причина этого заключается в том, что при горении углеводородных топлив объем необходимого для полного сгорания воздуха значительно превосходит объем самого горючего [169, 170]. Следовательно, с точки зрения газового состава, горящий огненный шар в основном состоит из нейтрального газа и продуктов горения, нагретых за счет выделяющегося тепла. Углеводородные газы имеют близкие теплоты сгорания и температуры пламени и, поэтому, весьма близкие степени расширения при горении. Таким образом, объем огненного шара определяются лишь полным количеством выделяющегося тепла, пропорционального полной массе топлива. Диаметр же огненного шара, пропорциональный кубическому корню из объема нагретой области, практически не зависит от конкретных параметров топлива. С этой точки зрения теоретическая зависимость диаметра огненного шара от массы топлива, приведенная в [177], согласно которой диаметр обратно пропорционален кубическому корню из молекулярной массы горючего, представляется не вполне верной.

Времена жизни огненного шара tF B, полученные разными авторами, согласуются значительно хуже: как видно из Таблицы 3.1, показатель степени b2 заключен в диапазоне от 1/6 до 1/3, большой разброс имеет и коэффициент a2. Существуют практические рекомендации, согласно которым для огненных шаров малого масштаба (с массой топлива менее 10 кг) время жизни шара лучше описывается степенной функцией с показателем 1/6, а крупномасштабные шары соответствуют степени 1/3 [116]. Такие «правила» являются слишком грубыми, поскольку они не отражают реальных причин расхождений в экспериментальных зависимостях. Эти расхождения отчасти порождаются тем, что эксперименты ставятся в различных условиях и для регистрации времени жизни огненного шара применяются различные методики.

Однако более существенными являются различия в начальных условиях — образован ли огненный шар при зажигании неподвижного объема газа, либо горючее выбрасывалось в атмосферу в виде направленной струи. При направленном выбросе Глава 3. Теория выбросов топлива в атмосферу процессы смешения интенсифицируются, что ведет к увеличению массовой скорости горения и сокращает время существования огненного шара.

Исследованию влияния начальной скорости на время жизни огненного шара посвящена экспериментальная работа [180], где методом высокоскоростной киносъемки измерялись времена сгорания выбросов метана и пропана с массой топлива от 1,5 до 13 г при скоростях истечения до 88 м/с. Исследовались направленные горизонтальные и вертикальные выбросы топлива из отверстий различного диаметра и с различными выходными насадками. В результате были получены следующие зависимости времени сгорания огненного шара от скорости выброса. В случае горизонтального направленного выброса топлива тогда как при вертикальных выбросах где характерные масштабы времени и скорости введены как t = (M0 /a )1/6 g 1/2, U = (M0 /a )1/6 g 1/2.

Обработка экспериментальных данных и их представление в виде (3.1) тесно связана с получением соответствующих теоретических зависимостей методами подобия и размерности [190, 191]. Приведем основные результаты такого анализа в условиях, когда течение в огненном шаре носит развитый турбулентный характер, так что скорость турбулентного горения определяются интенсивностью смешения топлива с окислителем из окружающего воздуха, а кинетические константы реакции горения не входят в число определяющих параметров. Как показывают эксперименты, вязкие эффекты становятся несущественными при массах топлива, превышающих 0,5 г [179]. В работе [180] выделено 17 определяющих параметров, в число которых входят характеристики окружающей атмосферы (плотность a, удельная теплоемкость Cp, температура Ta, коэффициент диффузии, кинематическая вязкость и Глава 3. Теория выбросов коэффициент теплопроводности), характеристики топлива (теплота сгорания Hc, массовый стехиометрический коэффициент в реакции горения, температура пламени, плотности исходного топлива и продуктов горения), параметры выброса (общая масса топлива M0 и скорость истечения U0), а также ускорение силы тяжести g и характеристики возникающего при горении огненного шара - его диаметр DF B, высота подъема HF B и время горения tF B. Поскольку все величины выражаются через четыре независимых размерности (масса, длина, время, температура), то, согласно Пи-теореме, связь между 17 размерными величинами может быть выражена в виде соотношения между 17 4 = 13 безразмерными комплексами. При рассмотрении огненных шаров углеводородных топлив многие из этих безразмерных комплексов (отношение температуры пламени к внешней температуре, степень расширения газа, безразмерная теплота сгорания) практически являются константами и, поэтому, могут быть исключены из дальнейшего рассмотрения. Кроме того, для развитых турбулентных течений несущественными оказываются и параметры, основанные на ламинарных коэффициентах переноса — числа Грасгофа, Прандтля и Шмидта.

Определим в качестве масштаба длины величину L = (M0 /a )1/3. В задаче имеются две характерные скорости, связанные с начальным импульсом газа U 0 и с действием сил плавучести U = (L g)1/2 = (M0 /a )1/6 g 1/2. Соответствующие характерные времена вводятся как t = L /U0 и t = L /U = (M0 /a )1/6 g 1/2.

Тогда в общем случае зависимость диаметра огненного шара от времени имеет вид В случае чисто плавучих выбросов (U0 U ) время развития конвекции t много меньше времени t, поэтому существенной оказывается только зависимость от t/t. Соответственно, время жизни огненного шара определяется моментом достижения максимума функцией (3.4), что дает tF B /t = const или tF B t. Если же скорость выброса велика, так что преобладающее влияние имеют силы инерции (U0 U), в качестве определяющего времени выступает t. В этом случае время горения tF B t. Таким образом, для чисто плавучих выбросов время жизни огненного Глава 3. Теория выбросов шара пропорционально M0, тогда как для инерционных выбросов это время пропорционально M0 и обратно пропорционально скорости U0. Заметим, что в качестве U0 может выступать не только скорость направленного выброса, но и, например, характерная скорость расширения газа при взрыве резервуара высокого давления.

При данном анализе размерностей в число определяющих параметров не были включены величины, характеризующие тепловое излучение огненного шара. Связано это с тем, что процессы радиационного теплопереноса не влияют напрямую на скорость турбулентного горения, определяемую только процессами смешения в турбулентных вихрях. Характерное время процессов остывания огненного шара за счет высвечивания лучистой энергии, как показывают оценки [192], более чем на порядок величины превышает характерное время смешения газа в огненном шаре с окружающим воздухом. Радиационные теплопотери ведут к снижению температуры огненного шара, однако размеры горящего облака и время его существования при этом изменяются незначительно [180]. С этой точки зрения не выдерживает критики анализ размерностей, проведенный в [116, 193, 194] и основанный на введении лишь одного характерного времени, определенного как отношение полной тепловой энергии горения к суммарной мощности теплопотерь облака излучением. В то же время, в число определяющих параметров не была включена величина ускорения силы тяжести, являющаяся фундаментальной при рассмотрении плавучих облаков, к которым относятся горячие огненные шары. Хотя авторы и получили зависимость времени горения tF B M0, наблюдающуюся в экспериментах, аргументация такой зависимости не может считаться удовлетворительной.

Анализ размерностей дает две различные зависимости времени горения огненного шара от массы топлива в зависимости от характеристик выброса, что и объясняет большой разброс экспериментальных зависимостей, приведенных в Таблице 3.1.

Экспериментальные зависимости (3.2), (3.3) дают экстраполяцию между предельными случаями чисто инерционных и чисто плавучих выбросов на тот диапазон параметров, где влияния обоих факторов сравнимы между собой.

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

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

3.2. Классификация выбросов конечной продолжительности Рассмотрим горючий газ, истекающий в атмосферу через круглое отверстие диаметром D со скоростью U0 в течение времени tr. Смешение выбрасываемого газа с окружающим воздухом может быть охарактеризовано некоторым характерным временем tmix, зависящим параметров источника, свойств вещества, геометрии выброса Глава 3. Теория выбросов и т. д. Распределения концентраций, возникающие в атмосфере в результате выброса, зависят от соотношения времен истечения и смешения, возможные типы истечений схематически изображены на рис. 3.1а-в. Если выброс происходит практически мгновенно (за время tr tmix ), в атмосфере возникает переобогащенное топливом облако, сгорающее при зажигании в диффузионном режиме в виде огненного шара (рис. 3.1а). Напротив, в если продолжительность выброса значительно превосходит характерное время турбулентного смешения газа с окружающим воздухом t mix, в атмосфере сформируется квазистационарная струя, зажигание которой ведет к образованию горящего факела (рис. 3.1в). Между этими двумя предельными случаями может наблюдаться множество промежуточных, соответствующих различным соотношениям характерных времен tr и tmix (рис. 3.1б).

Чтобы получить количественный критерий, позволяющий классифицировать выбросы конечной продолжительности, рассмотрим два типичных вида выброса при разгерметизации газовых емкостей — облако конечной массы, выпущенное с направленным начальным импульсом, и развивающуюся нестационарную струю, возникающую при резком включении источника массы и импульса. Считается, что силы плавучести пренебрежимо малы (чисто инерционное течение), не учитывается влияние ветра и не рассматриваются препятствия, ограничивающие выброс. В каждом случае в качестве характерного времени смешения принимается время, за которое объемная концентрация выброшенного газа падает до верхнего концентрационного предела горения (ВКП), что дает временной масштаб процессов перемешивания, делающих газовую смесь горючей. Отметим, что используется не нижний (НКП), а верхний предел горения, поскольку целью данной главы является описание переобогащенного топливом облака, которое при зажигании горит в диффузионном режиме. Расчет линий уровня концентрации, соответствующих НКП, обычно производят для оценки количества газа, находящегося в пределах горения, который может в случае зажигания сгорать во взрывном режиме или детонировать. Движение газа считается осесимметричным, кроме того, принимается, что время установления авГлава 3. Теория выбросов Рис. 3.1. Выбросы конечной продолжительности: а) приводящие к образованию облака, б) промежуточные, в) струйные Глава 3. Теория выбросов томодельных распределений скорости и концентрации намного меньше характерного времени разбавления газа до ВКП. Это обусловлено тем, что при струйном истечении автомодельные распределения устанавливаются на расстояниях порядка нескольких диаметров выходного отверстия [8, 12, 145], тогда как разбавление газа до ВКП, как будет показано ниже, происходит на значительно больших расстояниях. Следовательно, основная часть времени разбавления приходится на автомодельную стадию.

Эволюция облака мгновенного выброса может быть описана интегральной моделью [8, 153, 176, 195]. Облако аппроксимируется конусом, размеры которого возрастают по мере движения из-за смешения с окружающим воздухом. Динамика движения облака и изменение его объема находятся интегрированием законов сохранения массы, энергии и количества движения. Коэффициент вовлечения C, равный среднему тангенсу угла расширения облака, находится, согласно [8], в диапазоне от 0,16 до 0,58, причем большой разброс данных обусловлен внутренне присущей нерегулярностью турбулентного облака. Ниже используется наиболее типичное значение коэффициента вовлечения, принятое в литературе C 0,25 [8, 153]. Увеличение объема облака за счет смешения с воздухом вполне удовлетворительно описывается формулой где V — текущий, V0 — начальный объемы облака, — плотность газа, индексы a и r здесь и ниже относятся к окружающему воздуху и газу на уровне выходного отверстия соответственно. Поскольку общая масса горючего газа в облаке есть величина постоянная, его средняя объемная концентрация падает с увеличением объема как C = b V0 ma /a V mg (где ma и mg — молекулярные массы воздуха и газа), отсюда время, необходимое для разбавления газа до ВКП оценивается как где CU F L — объемная концентрация, соответствующая ВКП.

Процесс смешения газа с воздухом в нестационарной (развивающейся) Глава 3. Теория выбросов струе при резком включении источника может быть описан следующим образом.

Как показывают эксперименты, развивающаяся струя состоит из головной части и следующей за ней конической части, которая подобна установившейся напорной струе [8,142]. Аппроксимируем распределение концентрации в развивающейся струе модельным, заменив головную часть резким фронтом, находящимся на расстоянии HJ от виртуального источника, причем перед фронтом концентрацию будем считать нулевой, а за фронтом будем использовать поле концентрации в стационарной осесимметричной струе с гиперболическим законом затухания осевой концентрации и гауссовым распределением концентрации в радиальном направлении при z HJ :

где Cax — осевая объемная концентрация, BJ — эмпирическая константа, значение которой по данным различных экспериментов находится в диапазоне от 4.0 до 5.9 [9, 12, 125, 196–198], ниже используется приведенное в [12] значение BJ = 4,75.

Тангенс угла раскрытия струи J, согласно данным [9, 198], составляет 0,10 – 0,16, в дальнейшем принимается рекомендованное в [9] значение J = 0,132.

Зависимость координаты фронта от времени можно найти из закона сохранения массы, приравняв массу газа, выброшенную источником к моменту t, массе газа, находящейся в атмосфере и определяемой интегрированием концентрации (3.6), домноженной на плотность газа при нормальных условиях g,a.

Интегрирование, производимое в пределах от 0 до по радиальной координате и от 0 до HJ по осевой координате, дает зависимость где значение множителя равно 2,45 при J = 0,132 и BJ = 4,75.

Интересно, что тот же функциональный вид зависимости HJ (t) может быть получен с использованием интегральной модели [8, 176, 195], в которой модели струя аппроксимируется конусом с высотой HJ и радиусом основания RJ = J HJ. Приравнивая текущее количество движения конуса 3/2a V U (где V = RJ HJ /3 — Глава 3. Теория выбросов объем конуса, U = dHJ /dt — скорость фронта) полному импульсу (D 2 /4)r U0 t, созданному источником к моменту t, и интегрируя полученное дифференциальное уравнение, придем к соотношению (3.7), в котором числовой множитель должен быть заменен на (J )1/2. Значения тангенса угла раскрытия внешней границы струи, приводимые в литературе, составляют J = 0,166 – 0,25 [8, 142, 153], что дает значения константы от 2,0 до 2,4, довольно близкие к величине, полученной выше.

Характерное время смешения в развивающейся струе введем как время, в течение которого передний фронт струи достигает точки, в которой осевая концентрация равна CU F L. Из (3.6), (3.7) следует, что это время равно Выразим теперь время истечения tr как отношение полной массы выброшенного газа M0 = 0 V0 к полному массовому расходу газа через выходное отверстие G:

где r, U0 — плотность и скорость газа на уровне выходного отверстия, Cd = 0,6–0,9 — коэффициент истечения, учитывающий особенности выходного отверстия [119–121]. Ниже используется значение Cd = 0,85, рекомендованное для газовых истечений в [199, 200].