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

Российская Академия Наук

Институт Проблем Механики

На правах рукописи

ЯКУШ СЕРГЕЙ ЕВГЕНЬЕВИЧ

УДК 536.46

ГИДРОДИНАМИКА И ГОРЕНИЕ ГАЗОВЫХ И ДВУХФАЗНЫХ

ВЫБРОСОВ В ОТКРЫТОЙ АТМОСФЕРЕ

Специальность 01.02.05 — Механика жидкостей, газа и плазмы

ДИССЕРТАЦИЯ

на соискание ученой степени доктора физико-математических наук

Москва — 2000 Оглавление Введение 4 1. Автомодельный турбулентный термик в атмосфере с переменной плотностью 1.1. Основные закономерности эволюции плавучих облаков и методы их исследования.............................. 1.2. Математическая модель термика в сжимаемой атмосфере....... 1.3. Численное моделирование термика в экспоненциальной изотермической атмосфере.............................. 1.4. Приближенное аналитическое решение................. 1.5. Выводы................................. 2. Подъем термика и перенос дисперсных примесей в атмосфере 2.1. Математическое моделирование термиков на основе уравнений сжимаемого газа. Прикладные задачи.................... 2.2. Уравнения движения термика с дисперсной примесью......... 2.3. Начальное состояние запыленного термика............... 2.4. Автомодельный подъем термика, весовое и тепловое влияние примеси 2.5. Структура всплывающего термика при различной загрузке примесью. 2.6. Вынос дисперсной примеси в стратосферу............... 2.7. Метод расчета течений сжимаемого газа................ 2.8. Выводы................................. 3. Образование и горение газовых облаков при аварийных выбросах в атмосферу. Физическая теория 3.1. Аварийные выбросы, их особенности и методы изучения....... 3.2. Классификация выбросов конечной продолжительности........ 3.3. Истечение из резервуаров низкого и высокого давления........ 3.4. Безразмерный вид критерия....................... 3.5. Сравнение с экспериментом....................... 3.6. Влияние задержки воспламенения на режим сгорания топлива.... 3.7. Коэффициент участия топлива при горении огненного шара...... 3.8. Примеры использования модели.................... 3.9. Выводы................................. 4. Огненный шар при горении выбросов газового углеводородного топлива 4.1. Математические модели огненных шаров................ 4.2. Постановка задачи............................ 4.2.1. Основные уравнения....................... 4.2.2. Модель образования и выгорания сажи............. 4.2.3. Модель переноса излучения................... 4.3. Определяющие параметры....................... 4.4. Горение огненного шара: расчет без учета излучения.......... 4.5. Горение излучающих углеводородных шаров.............. 4.6. Структура радиационного поля в огненном шаре............ 4.7. Тепловые потоки и оценка воздействия огненного шара........ 4.8. Метод расчета существенно дозвуковых течений............ 4.8.1. Приближение малых чисел Маха................ 4.8.2. Решение эллиптических уравнений............... 4.8.3. Расчет тепловых потоков методом Монте-Карло........ 4.9. Выводы................................. 5. Горение облаков углеводородных топлив при двухфазных выбросах в атмосферу 5.1. Образование облаков аэрозолей при выбросах в атмосферу...... 5.2. Основные уравнения........................... 5.4. Начальные и граничные условия.................... 5.5. Определяющие параметры при двухфазных истечениях........ 5.6. Эволюция двухфазного облака без зажигания............. 5.7. Огненный шар при зажигании двухфазного выброса.......... 5.8. Масштабные эффекты при горении двухфазных выбросов....... Введение Весьма широкий круг физических явлений природного и техногенного происхождения может быть охарактеризован как выброс инородного вещества в окружающую атмосферу. Явления, которые можно отнести к выбросам, весьма различны по своему масштабу, типам источника, фазовому составу и протекающим химическим процессам. При всем их разнообразии объединяющую роль играет возникновение в относительно однородной окружающей среде локализованной области с отличающимися от внешних свойствами, что определяет дальнейшую эволюцию, характер и степень взаимодействия с окружающей средой, а зачастую — и опасность выброса.

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

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

Образование, эволюция и горение топливного облака при выбросе горючего газа в атмосферу — сложный процесс, включающий целый ряд явлений: турбулентное смешение выброшенного вещества с атмосферным воздухом, приводящее к образованию горючей смеси; воспламенение от источника зажигания, диффузионное горение или горение предварительно перемешанных реагентов, протекающее в турбулентном режиме; тепловое излучение. Еще более широкий спектр физических процессов характерен для двухфазных выбросов, когда в атмосфере образуется облако, содержащее смесь паров и мелкодисперсных капель горючего вещества. Многообразие физико-химических явлений, сопровождающих образование, эволюцию и горение газовых и двухфазных выбросов, приводит к тому, что изучение этого класса течений возможно только с применением междисциплинарного подхода, совмещающего экспериментальные исследования и достижения нескольких теоретических дисциплин — гидродинамики, газовой динамики, теории горения и взрыва, механики многофазных сред, вычислительных методов.

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

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

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

Диссертация состоит из введения, пяти глав, заключения и списка литературы.

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

В главе 1 рассмотрено движение турбулентного осесимметричного термика в атмосфере с переменной по высоте плотностью. Модель построена на основе системы уравнений «глубокой» конвекции, в которой пренебрегается динамической сжимаемостью газа, но принята во внимание гидростатическая сжимаемость атмосферы.

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

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

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

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

Введение Глава 4 посвящена численному моделированию образования, горения и теплового воздействия огненного шара, возникающего при зажигании выброса углеводородного газа (метан, пропан) вблизи поверхности земли. Осесимметричное нестационарное течение описывается системой осредненных по Фавру уравнений сохранения совместно с (k ) моделью турбулентности, моделью турбулентного горения и глобально-кинетической схемой образования и выгорания частиц сажи. Проведены расчеты огненных шаров малого масштаба (с массой топлива от 2 до 20 г) без учета процессов радиационного теплопереноса, а также излучающих огненных шаров с массой топлива в диапазоне от 1 г до 1000 кг. Выявлена роль масштабных эффектов, получена зависимость безразмерного времени выгорания огненного шара от числа Фруда, хорошо согласующаяся с экспериментом. Методом Монте-Карло рассчитаны падающие на поверхность потоки теплового излучения, а также термические дозы, полученные за время горения облака. Рассчитанная доля энергии, переходящей в излучение, хорошо согласуется с литературными данными.

В главе 5 исследованы двухфазные выбросы при разгерметизации резервуара со сжиженным газом. Нестационарные течения смесей пара и капель горючего описываются на основе эйлерово-лагранжевого подхода с учетом обмена массой, теплом и импульсом между несущей (газовой) и дисперсной фазами. В качестве топлива использован пропан, масса выброса варьируется в диапазоне от 1 г до 1000 кг, температура хранения в резервуаре составляет от 268 до 351 K. Исследованы негорящие выбросы различного масштаба, выделены основные режимы испарения капель в облаке. Получена температурная, концентрационная и радиационная структура огненного шара на всех стадиях эволюции, от момента зажигания до полного выгорания топлива. Рассчитанные распределения температуры и продуктов горения сопоставлены с контурами видимого огненного шара, полученными в экспериментах Хасегавы и Сато. Отмечено хорошее соответствие размеров, формы и высоты подъема огненного шара в соответствующие моменты времени. Предложено единое описание безразмерного времени выгорания однофазных и двухфазных огненных шаров.

Введение На защиту выносятся следующие основные положения:

• численное и аналитическое исследование автомодельного движения плавучего облака в атмосфере с переменной плотностью • численное исследование подъема плавучих облаков (термиков) в стратифицированной сжимаемой атмосфере и вертикального переноса такими термиками газовых и дисперсных примесей • разработка физической теории выбросов конечной продолжительности, количественных критериев классификации режимов истечения и сгорания топлива • исследование формирования, эволюции и горения огненных шаров газовых и распыленных жидких топлив в атмосфере, установление внутренней структуры, радиационных характеристик, динамики развития горящих облаков • количественный анализ воздействия выбросов на окружающую среду, связанного с присутствием газовых и дисперсных примесей, теплового воздействия горящих выбросов • изучение роли масштабных явлений при горении газовых и двухфазных облаков, получение единых зависимостей для характеристик огненных шаров в широком диапазоне типов выбросов и параметров топлив Исследования выполнялись автором в течение 1986–99 гг в Институте проблем механики РАН согласно планам научно-исследовательских работ института и в рамках Государственной научно-технической программы «Безопасность» (пост. ГКНТ №1011 от 01.07.91 г), а также во время научной работы в Университете Центрального Ланкашира (Престон, Великобритания). Автор приносит глубокую благодарность профессору Г. М. Махвиладзе за постоянное внимание к работе, всемерную помощь и поддержку на всех этапах ее выполнения. Автор признателен д.ф.-м.н. О. И. Мелихову и сотрудникам лаборатории термогазодинамики ИПМ РАН за помощь и ценные обсуждения.

Введение Результаты, вошедшие в диссертацию, докладывались на следующих конференциях и симпозиумах:

• VIII Всесоюзном Симпозиуме по горению (Ташкент, 1986), VI Всесоюзном совещании по теоретическим и прикладным аспектам турбулентных течений (Таллинн, 1989), • Всесоюзной конференции по методам математического моделирования в задачах охраны природной среды и экологии (Новосибирск, 1991), • Конференции IMACS по математическому моделированию и прикладной математике (Москва, 1990), • Первой азиатской конференции по исследованиям и технологии борьбы с пожарами (Хэфэй, Китай, 1992), • 1-й и 2-й Международной конференции по химической промышленности и окружающей среде (Жерона, Испания, 1993; Альгедо, Сардиния, 1996), • Международном семинаре по проблемам пожаровзрывобезопасности (Престон, Англия, 1994), • Международном симпозиуме по крупномасштабным опасностям в промышленности и на море (Манчестер, Великобритания, 1995), • Конференции Британской информационной группы по исследованию взрывов (Абериствис, Уэльс, 1995), • 26 Международном симпозиуме по горению (Неаполь, Италия, 1996), • 5-м и 6-м Международном симпозиуме по научным основам пожаробезопасности (Мельбурн, Австралия, 1997; Пуатье, Франция, 1999), • 3 Международном семинаре по пожаровзрывобезопасности (Виндермер, Англия, 2000).

Введение Материалы диссертации неоднократно обсуждались на научных семинарах Института проблем механики РАН и научных семинарах Центра по исследованиям пожаров и взрывов Университета Центрального Ланкашира (Англия). Основные результаты, вошедшие в диссертацию, опубликованы в работах:

1) Г. М. Махвиладзе, О. И. Мелихов, С. Е. Якуш. Численное исследование подъема термика с частицами в стратифицированной атмосфере. VIII Вcесоюзный симпозиум по горению и взрыву. Кинетика и горение. Черноголовка, 1986, с. 16–20.

2) Г. М. Махвиладзе, О. И. Мелихов, С. Е. Якуш. Турбулентный осесимметричный термик в неоднородной сжимаемой атмосфере. Численное моделирование. Препринт № 303. ИПМ АН СССР, 1987, 67 c.

3) Г. М. Махвиладзе, О. И. Мелихов, С. Е. Якуш. О численном моделировании подъема турбулентного термика в неоднородной сжимаемой атмосфере. Изв. АН СССР, МЖГ, 1989, № 1, c. 72–80.

4) Г. М. Махвиладзе, О. И. Мелихов, С. Е. Якуш. Подъем турбулентного осесимметричного термика в неоднородной сжимаемой атмосфере. ПМТФ, 1989, № 1, с. 62–68.

5) Г. М. Махвиладзе, С. Е. Якуш. Эволюция запыленного термика и вынос аэрозольных частиц в верхние слои атмосферы. Препринт № 368, ИПМ АН СССР, 1989, 44 c.

6) Г. М. Махвиладзе, С. Е. Якуш. Автомодельный режим подъема термика в среде с переменной плотностью. В сб.: Турбулентные течения и техника эксперимента. Таллинн, 1989, с. 127–130.

7) Г. М. Махвиладзе, С. Е. Якуш. Перенос дисперсной примеси в атмосфере всплывающим термиком. Изв. АН СССР, МЖГ, 1990, № 1, с. 123–130.

8) Г. М. Махвиладзе, С. Е. Якуш. О влиянии дисперсной примеси на подъем запыленного теплового термика. ПМТФ, 1990, № 5, с. 69–77.

Введение 9) Г. М. Махвиладзе, С. Е. Якуш. Автомодельный осесимметричный термик в среде с переменной плотностью. Изв. АН СССР, МЖГ, 1991, № 4, с. 45–53.

10) G. M. Makhviladze, I. K. Selezniova, S. E. Yakush. Gaseous and particulate clouds in the unbounded atmosphere. In: A. A. Samarskii and M. P. Sapagovas (Eds.), Mathematical Modelling and Applied Mathematics, Elsevier, North Holland, 1992, pp. 279–287.

11) S. E. Yakush. Pollution of the atmosphere by re and explosion thermals. In: F. Weicheng and F. Zhuman (Eds.), First Asian conference on Fire Science and Technology, Int.

Academic Publishers, 1992, pp. 439–444.

12) S. E. Yakush. Atmospheric pollution by accidents. In: J. Casal (Ed.), Chemical Industry and Environment, Spain, Girona, 2-4 June 1993, v. 1, pp. 127–136.

13) G. M. Makhviladze, J. P. Roberts, S. E. Yakush, D. Davis. Criterion for the formation of a “cloud-like” release upon depressurisation of a gas vessel. IChemE Symp. Series, 1995, v. 139, pp. 97–112.

14) G. M. Makhviladze, J. P. Roberts, S. E. Yakush. Modelling of atmospheric pollution by explosions. Environmental Software, 1995, v. 10, No. 2, pp. 117–127.

15) G. M. Makhviladze, J. P. Roberts, S. E. Yakush. Turbulent buoyant thermal in a densitystratied atmosphere. Int. J. Heat Mass Transfer, 1996, v. 39, No. 7, pp. 1453–1462.

16) G. M. Makhviladze, J. P. Roberts, S. E. Yakush. Burning regimes for the nite-duration releases of fuel gases. Twenty Sixth Int. Symp. on Combustion, Naples, Italy, Jul – Aug 2, 1996. The Combustion Institute, Pittsburgh, PA, USA, 1996, v. 1, pp. 1549– 17) G. M. Makhviladze, J. P. Roberts, S. E. Yakush. On accidental fuel releases into the atmosphere. In: N. Piccinini and R. Delorenzo (Eds.) Chemical Industry and Environment II, Alghero, Sardinia, Italy, 1996. Politecnico di Torino, 1996, v. 3, pp. 873–880.

Введение 18) Г. М. Махвиладзе, Дж. П. Робертс, С. Е. Якуш. Образование и горение газовых облаков при аварийных выбросах в атмосферу. Физика горения и взрыва, 1997, Т. 33, 19) G. M. Makhviladze, J. P. Roberts, S. E. Yakush. Modelling the reballs from methane releases. In: Y. Hasemi (Ed.) 5th Int Symp. on Fire Safety Science, Melbourne, Australia.

IAFSS, 1997, pp. 213–224.

20) G. M. Makhviladze, J. P. Roberts, S. E. Yakush. Numerical modelling of reballs from vertical releases of fuel gases. Comb. Sci. and Techn., 1998, v. 132, No. 1-6, pp. 199– 21) Г. М. Махвиладзе, Дж. П. Робертс, С. Е. Якуш. Огненный шар при горении выбросов углеводородного топлива: I. Структура и динамика подъема. Физика горения и взрыва, 1999, Т. 35, №. 3, с. 7–19.

22) Г. М. Махвиладзе, Дж. П. Робертс, С. Е. Якуш. Огненный шар при горении выбросов углеводородного топлива: II. Тепловое излучение. Физика горения и взрыва, 1999, Т. 35, №. 4, с. 12–23.

23) G. M. Makhviladze, J. P. Roberts, S. E. Yakush. Combustion of two-phase hydrocarbon fuel clouds released into the atmosphere. Combustion and Flame, 1999, v. 118, No. 4, pp. 583–605.

24) G. M. Makhviladze, J. P. Roberts, S. E. Yakush. Modelling and scaling of reballs from single- and two-phase hydrocarbon releases. In: 6th Int Symp. on Fire Safety Science, Poitiers, France. IAFSS, 1999.

Глава Автомодельный турбулентный термик в атмосфере с переменной плотностью 1.1. Основные закономерности эволюции плавучих облаков и методы их исследования Экспериментальные и теоретические исследования свободно-конвективных течений, развивающихся в открытой атмосфере при возникновении объема газа с отличной от окружающей плотностью и обладающего за счет этого ненулевой плавучестью, насчитывают уже несколько десятилетий. Важнейшее свойство конвективных течений состоит в том, что на достаточном удалении от источника выброса течение «забывает» начальные условия и выходит на автомодельный режим. На этой стадии характеристики течения зависят от малого числа параметров, что позволяет плодотворно применять методы подобия и размерностей, получая законы изменения концентраций, температуры, скорости газа. В качестве предельных случаев свободно-конвективных течений различают плавучие струи, возникающие при длительном действии источника, и плавучие облака (термики). Отличие плотности газа от окружающей может быть обусловлено различными температурами или молекулярными массами газа и воздуха.

Предельные законы свободно-восходящих конвективных потоков от непрерывно действующего источника были установлены в работе [1]. Характеристики плаГлава 1. Автомодельный термик вучих струй в дальнейшем изучались на основе моделей с осреднением по сечению струи [2–4]. Широкий класс автомодельных решений, описывающих гидродинамическую, тепловую и концентрационную структуру турбулентных течений в свободной атмосфере при концентрационно-тепловом и динамическом воздействии получен в работах [5, 6]. В настоящее время изучение характеристик плавучих струй в неподвижной атмосфере и в сносящем потоке составляет самостоятельное направление гидродинамики и теории конвекции, что вызвано широким распространением этого класса течений во многих прикладных задачах. Подробный обзор основных достижений в изучении плавучих струй можно найти в [7–14].

Рассмотрим более подробно основные методы исследования и результаты по плавучим облакам (термикам), являющимся предметом изучения этой и последующей глав. При этом в настоящем разделе внимание сосредоточено на тех теоретических работах, где численно и аналитически решаются уравнения несжимаемой жидкости в приближении Буссинеска. Обзор имеющихся работ, в которых использовались уравнения, учитывающие сжимаемость среды, представлен в начале Главы 2.

Под термиком обычно понимают облако нагретого газа, возникающее в результате выделения определенного количества тепла. С точки зрения терминологии будем отличать термики от огненных шаров, рассмотренных в Главах 3–5: термик характеризуется быстрым выделением тепла, так что в дальнейшем его эволюция протекает без дополнительного энерговклада. В то же время под огненным шаром будем понимать горящее облако, в котором выделение тепла при горении топлива происходит в течение времени, сравнимого или большего времени развития конвекции. Наряду с высокой температурой, плавучесть термика может быть обусловлена и малой молекулярной массой выброшенного в атмосферу вещества, поэтому в общем случае часто говорят о концентрационно-тепловых термиках.

Основные характеристики плавучих облаков на автомодельной стадии движения были установлены теоретически и экспериментально еще в 50–60-х гг [2,3,15]. После того, как термик достаточно удаляется от источника, в задаче исчезает характерГлава 1. Автомодельный термик ный масштаб длины и единственным определяющим параметром становится полный запас плавучести термика B0, равный объемному интегралу от ускорения архимедовой силы, который в нейтрально стратифицированной атмосфере является интегралом движения. Радиус термика на этой стадии растет пропорционально высоте верхней кромки облака r = C zt, а вертикальная скорость (в осесимметричном случае) падает обратно пропорционально высоте v B0 zt. Поскольку скорость подъема равна производной координаты верхней кромки термика по времени, отсюда следует и закон подъема термика zt B0 t1/2 [2, 7, 8]. Условие сохранения интеграла от избыточной температуры по объему термика, выражающее закон сохранения тепловой энергии облака, приводит к соотношению r 3 = const, что дает закон спадания избыточной температуры в термике со временем t3/2.

Корневой закон движения осесимметричного термика был установлен экспериментально в опытах [15], где термики с отрицательной плавучестью создавались при выливании тяжелого соляного раствора в пресную воду. Отмечалось, что типичной формой термика является сфероид, несколько сплющенный в вертикальном направлении. Были подтверждены зависимости zt t1/2, v t1/2, а также найден тангенс угла расширения термика C = 0,25. Несколько большие значения угла раскрытия термика C = 0,28 были зафиксированы в экспериментах по подъему объема легкой жидкости в устойчиво стратифицированном соляном растворе [3]. Значения коэффициента пропорциональности в корневом законе зависимости координаты верхней кромки термика в широком диапазоне параметров плавучих облаков установлены в работе [16,17]. Обработка собственных экспериментальных данных, полученных при зажигании зарядов пороха с массой 5–20 г в лабораторных условиях и с массой 3– 27 кг в полевых условиях, а также литературных данных по динамике подъема термиков ядерных взрывов [18], позволила установить единую зависимость координаты верхней кромки термика:

где zt отсчитывается от виртуального источника. Эксперименты с концентрационныГлава 1. Автомодельный термик ми термиками малого масштаба показали также, что развитый турбулентный режим течения в термике реализуется при ламинарных числах Рэлея Ral > 3 · 1010, при этом автомодельная координата верхней кромки (равная коэффициенту пропорциональности в (1.1)) не зависит от числа Рэлея. Тангенс угла расширения термика на автомодельном участке движения в экспериментах [16, 17] равен C = 0,18–0,21.

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

Гипотеза о вовлечении использовалась, например, в работах [3, 19–23], эти модели отличаются лишь записью сил сопротивления, учетом присоединенной массы, предположениями о малости или существенности перепада плотностей газа в термике и окружающей атмосфере. В работе [19] наряду с силами плавучести учитывается циркуляция вихревого кольца и действие силы Жуковского, горизонтальная компонента которой растягивает вихревое кольцо в радиальном направлении, а вертикальная тормозит подъем. На автомодельном участке подъема модель дала угол расширения облака 0,32–0,34, что превышает экспериментальное значение 0,25 [15].

Модели с осредненными по объему параметрами использовались и для описания начальной стадии движения термика, на которой неподвижный газ приходит в вихревое движение. Так, в работе [24] построены решения, описывающие как начальную, так и автомодельную стадии процесса при нейтральной стратификации атмоГлава 1. Автомодельный термик сферы. На начальной стадии получены законы роста радиуса термика и координаты его верхней кромки, имеющие вид r = r0 + c1 t2, zt = zt,0 + c2 t2, где индексом 0 обозначены начальные значения, а c1, c2 — параметры, зависящие от начального перепада плотностей. На больших временах полученное решение выходит на асимптотику r t1/2, zt t1/2, причем зависимость от c1 и c2 при этом исчезает. Дальнейшее развитие эта теория получила в [25], где были учтены параметры турбулентности и стратификация атмосферы. Для определения параметра турбулентной диссипации использовались данные лабораторных экспериментов [26], в которых специальными неподвижными и вращающимися зондами измерялись характеристики турбулентного поля температуры. Отмечен резкий передний фронт термика и плавное спадание температуры в следе облака. К аналогичным выводам привели более поздние эксперименты по измерению внутренней структуры термиков больших масштабов [27].

Недостатком полученных в [25] решений является чрезмерно малый угол расширения облака C 0,014. Несоответствие между этой величиной и экспериментальным значением 0,25 было устранено в работе [28] путем учета присоединенной массы. Формирование вихревого течения на начальной стадии эволюции термика рассмотрено в [29], где на основе линеаризованного уравнения завихренности в случае невязкой несжимаемой жидкости получены соотношения для характерных времен, определяемых плавучестью, завихренностью и начальной скоростью. Процесс отрыва ламинарного термика от горизонтальной поверхности исследовался экспериментально в [30,31]. Формирование вихревого тороидального течения, изменение формы термика со временем и нарастание завихренности на начальной стадии подъема плавучего облака изучены теоретически и экспериментально в работах [32–34].

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

Высота зависания и размеры зависшего термика определены в работах [3,19,28]. Для стандартной атмосферы с температурным градиентом 6,5 C/км высота зависания Глава 1. Автомодельный термик H,[м] связана с выделившейся энергией Q0,[Дж] соотношением H = 1,87Q0 [3].

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

Аналитические решения, дающие пространственные распределения скорости, температуры и концентрации в термике, получены при большой вязкости среды в работах [4, 35–37]. Решение, как правило, строится путем разложения в ряд по малому параметру — числу Рэлея. В работе [36] методом последовательных приближений изучено распространение «ядра» плавучести и тороидального вихря. Отмечается, что вихревое кольцо распространяется сначала в боковом направлении, и лишь потом происходит формирование автомодельного всплывающего вихря. В работах [37–39] исследовано влияния сил плавучести на диффузию вихревого кольца. Решение уравнений движения и переноса тепла, полученное разложением по малому параметру (числу Рейнольдса), в качестве предельного случая содержит автомодельную стадию подъема всплывающего теплового термика.

Теория турбулентного термика при больших числах Рэлея была предложена в [40] и развита далее в [16, 41]. В этих работах уравнения несжимаемой жидкости решались в приближении вертикального погранслоя. Предполагалось подобие распределений ускорения сил плавучести и вертикальной скорости, для аппроксимации радиальных распределений искомых функций использовались гауссовы профили. Построена газодинамическая картина течения вблизи автомодельного термика и найдены координаты характерных точек термика в зависимости от величины коэффиГлава 1. Автомодельный термик циента турбулентности. Для стратифицированной атмосферы в [16, 41] рассмотрены колебания термика при зависании в устойчиво стратифицированной атмосфере, либо его разгон в неустойчиво стратифицированной среде. Полученные решения использовались для решения прикладных задач о выносе пассивно переносимой примеси в стратосферу, оценке тепловыделения при вулканических извержениях эксплозивного типа, генерации инфразвука всплывающим термиком, детонации концентрационного термика [16, 41, 42].

Первые попытки численного исследования конвективных элементов были предприняты в конце пятидесятых – начале шестидесятых годов [43–45]. Решались уравнения несжимаемой жидкости в приближении Буссинеска, зачастую без учета вязкости среды (роль эффективной вязкости играла схемная вязкость). Несовершенство численных методов и ограниченность вычислительных ресурсов позволили получить решения лишь на весьма грубых сетках и для ограниченного интервала времен.

Так, в [43] удалось на основе уравнений невязкой несжимаемой жидкости просчитать лишь начальную стадию эволюции на протяжении первых минут, дальнейший счет был невозможен из-за развития численной неустойчивости. Введение членов, описывающих турбулентную диффузию и вязкость, не сказалось существенным образом на получаемом решении и не смогло устранить указанные вычислительные трудности. Осесимметричный термик при малом начальном перепаде температур (порядка 1–3 ) моделировался в [45], расчеты проведены до времен порядка 30 мин. Отмечено наличие автомодельной стадии движения, на которой законы затухания максимальной температуры t3/2 и роста координаты верхней кромки t1/2 находятся в соответствии с результатами анализа размерностей и экспериментами [3, 15]. В то же время, термик в расчетах получался вытянутым по вертикали, тогда как эксперименты показывают, что при движении он должен сплющиваться [15].

Одна из проблем, с которой сталкивается расчет эволюции термиков, состоит в том, что геометрические размеры области, охваченной конвективным течением, увеличиваются со временем, поэтому на больших временах начинает сказываться влиГлава 1. Автомодельный термик яние границ расчетной области. В работе [46] была применена расширяющаяся со временем расчетная область, причем закон расширения выбирался таким, чтобы в новых переменных автомодельное решение было стационарным. Автомодельные решения были получены в широком диапазоне изменения коэффициентов переноса, при этом форма термиков на автомодельной стадии оказалась ближе к наблюдаемой в экспериментах, чем в работе [45]. Расширяющиеся по автомодельному закону сетки затем широко применялись для исследования автомодельной стадии подъема термика [16, 47–50].

Развитие осесимметричного облака в устойчиво стратифицированной атмосфере с тонким неустойчивым слоем численно моделировалось в [51]. Расчеты проводились при постоянных эффективных коэффициентах переноса, составляющих 50 м2 /с, выбор которых не обосновывался. На стадии зависания выше термика отмечалось возникновение зоны с обратной циркуляцией, что согласуется с данными наблюдений за природными облаками [52]. В работе [53] рассчитан подъем облака и конденсация содержащегося в нем водяного пара над постоянно действующим источником. В отличие от большинства выполненных в то время работ, использовались переменные (пропорциональные градиенту скорости) коэффициенты турбулентного переноса, а непосредственно над источником задавалась эффективная вязкость 300 м2 /с. В настоящее время численные методы широко применяются в метеорологии для моделирования облаков и других крупномасштабных течений [54].

Весьма подробное численное исследование различных стадий эволюции осесимметричного турбулентного термика в стратифицированной атмосфере выполнено в работах [16, 50]. Использовались уравнения несжимаемой жидкости в приближении Буссинеска, стратификация атмосферы учитывалась в уравнении притока тепла (и эквивалентном ему уравнении для ускорения силы плавучести). Для замыкания осредненных по Рейнольдсу уравнений введены постоянные эффективные коэффициенты турбулентной вязкости и температуропроводности, пропорциональные корню квадратному из полного запаса плавучести облака B0, причем величина коГлава 1. Автомодельный термик эффициента пропорциональности находилась из условия согласования автомодельной координаты верхней кромки с наблюдаемой в экспериментах величиной 4, (см. (1.1)). По результатам численных расчетов была построена единая зависимость автомодельной координаты верхней кромки от числа Рэлея и найдено, что экспериментальное значение автомодельной координаты достигается при числе Рэлея, равном 520. Найденные таким образом коэффициенты переноса, обеспечивающие правильный закон подъема на автомодельной стадии, использовались затем при расчете начальной стадии и стадии зависания. Рассмотрен целый ряд физических явлений, связанных с эволюцией турбулентных термиков — выход термика и вынос пассивной примеси в стратосферу, оценка энерговыделения вулканических извержений, образование окислов азота в высокотемпературном термике и их вынос в верхние слои атмосферы [55].

Таким образом, в настоящее время установлены основные закономерности и количественные характеристики эволюции плавучих облаков в стратифицированной атмосфере. Широко используется приближение Буссинеска, согласно которому плотность среды считается постоянной во всех членах уравнений, кроме архимедовой силы. При описании атмосферных течений в большинстве случаев учитывается лишь температурная стратификация атмосферы путем введения соответствующих источниковых членов в уравнение притока тепла. Предположение о постоянстве плотности атмосферы справедливо лишь в том случае, когда вертикальный размер области, охваченной конвективным течением, достаточно мал, так что весовая (изотермическая) сжимаемость газа не играет заметной роли. Это имеет место, например, в жидкостях, а также при «слабой» атмосферной конвекции, когда вертикальный масштаб течения намного меньше характерного масштаба изменения плотности атмосферы с высотой. Если же термик достаточно мощный, он может проникать на высоты порядка нескольких километров. При таких вертикальных масштабах течения эффекты сжимаемости атмосферы весьма существенны, поскольку плотность воздуха падает в несколько раз по сравнению с плотностью у земной поверхности. В таком Глава 1. Автомодельный термик случае в задаче имеется масштаб длины, связанный с вертикальным распределением плотности атмосферы, этот масштаб сохраняется и после того, как «забываются»

начальные условия. Поэтому при переменной плотности атмосферы автомодельного решения, строго говоря, не существует, а динамика подъема термика не может быть получена лишь из анализа размерностей. Тем не менее, корневой закон подъема (1.1) выполняется на определенном интервале времени даже в том случае, когда сжимаемость среды существенна, например, при подъеме термиков ядерных взрывов, высота зависания которых порядка 5–20 км, что сравнимо с вертикальным масштабом изменения плотности атмосферы.

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

1.2. Математическая модель термика в сжимаемой атмосфере Рассмотрим эволюцию осесимметричного термика, возникающего в стратифицированной атмосфере при мгновенном выделении некоторого количества тепла Q0.

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

Глава 1. Автомодельный термик Эволюция осесимметричного термика описывается в цилиндрической системе координат (r, z), начало которой совпадает с виртуальным источником, ось z направлена вертикально вверх, а ось r — в радиальном направлении. Предположим отклонения температуры и давления p от их значений в невозмущенной атмосфере малыми по сравнению с Ta и Pa. При этом плотность газа можно считать равной атмосферной плотности a везде за исключением членов, описывающих силу плавучести в уравнении для вертикальной проекции импульса. Сила плавучести, пропорциональная разности атмосферной и локальной плотностей, связана с отклонениями давления и температуры.

В линейном приближении уравнение состояния имеет вид где P = (/P )T = Pa, T = (/T )P = Ta. Поскольку рассматриваются только медленные (существенно дозвуковые) конвективные течения, ниже пренебрегается динамической сжимаемостью среды по сравнению с ее тепловым расширением, что эквивалентно |P p| |T |. Тогда выталкивающая сила, действующая на единицу объема, равна g(a ) = ga T. Подобное разложение часто применяется при исследовании глубокой конвекции [54]. Основное отличие этого подхода от классического приближения Буссинеска состоит в том, что линеаризация уравнения состояния (1.2) производится относительно локального значения невозмущенной плотности a (z), а не относительно некоторого фиксированного для всей области течения значения плотности. Это позволяет принять во внимание изменение атмосферной плотности с высотой, вызванное весовой сжимаемостью газа.

В рамках принятых предположений задача о расчете эволюции термика сводится к решению следующей системы уравнений:

a ur a v Глава 1. Автомодельный термик где u, v — радиальная и вертикальная компоненты скорости, µ и — эффективные коэффициенты динамической вязкости и теплопроводности, CP — теплоемкость газа, считающаяся постоянной, — отношение теплоемкостей при постоянном давлении и объеме, Rg = R /m — газовая постоянная (R = 8,31 Дж/моль·К — универсальная газовая постоянная, m — молекулярная масса газа), J — параметр стратификации, пропорциональный разности между локальным градиентом температуры и сухоадиабатическим градиентом g/CP. Эквивалентная форма представления параметра стратификации через градиент плотности атмосферы (второе выражение в правой части уравнения (1.7)) получено с использованием уравнения состояния совершенного газа Pa = a Rg Ta и соотношения CP = Rg /( 1).

Эффективные турбулентные коэффициенты вязкости µ и теплопроводности полагаются постоянными. Такая модель турбулентности хотя и является простейшей, позволяет тем не менее в расчетах на основе уравнений несжимаемой среды воспроизвести наблюдающийся в экспериментах корневой закон подъема термика. Возможность использования постоянных коэффициентов переноса при расчете свободно-конвективных течений обусловлена тем, что на автомодельной стадии характерный масштаб длины облака растет пропорционально B0 t1/2, а характерная скорость падает как B0 t1/2. Предполагая масштаб и скорость турбулентных пульсаций пропорциональными линейным размерам и средней скорости движения облака, получим, что произведение этих величин, дающее характерную кинематическую вязкость, не зависит от времени. Аналогичная ситуация справедлива, например, в осесимметричных турбулентных струях [12]. Подобное предположение использоваГлава 1. Автомодельный термик лось и в предыдущих работах, где были найдены аналитические решения, описывающие структуру термика [16, 40, 41].

Граничными условиями для системы уравнений (1.3)–(1.6) служат условия симметрии течения на оси r = 0 (u = 0, v/r = 0, /r = 0), а также условия отсутствия возмущений на бесконечности r 2 + z 2 (u = v = = 0).

Проведем анализ размерностей. Искомые функции u, v, зависят от пространственных координат r, z и времени t, а также от коэффициентов переноса µ,, параметров, характеризующих начальное состояние термика (его радиуса R 0 и тепловой энергии Q0 ) и параметров, задающих стратификацию атмосферы.

Как видно из вертикальных распределений невозмущенной температуры и плотности (1.7), стратификация атмосферы по температуре и по плотности характеризуется двумя линейными масштабами La = Rg T0 /g и LJ = g/J. Поэтому распределения температуры и плотности в невозмущенной атмосфере записываются в виде где индексом 0 обозначены параметры на уровне виртуального источника, и T — известные функции. Примером среды с переменной по высоте плотностью является изотермическая атмосфера, для которой из уравнения гидростатического равновесия следует a (z) = 0 exp(gz/Rg Ta ), так что (z) = exp(z/La ).

В дальнейшем примем, что температура меняется с высотой гораздо медленнее, чем плотность. Так, для стандартной атмосферы (J = 1,2 · 104 с2) зависимость Ta (z) близка к линейной, а a (z) — к экспоненциальной. При подъеме на высоту 10 км абсолютная температура уменьшается на 30%, тогда как плотность падает в 3–4 раза. Эти оценки позволяют пренебречь изменением коэффициента теплового расширения T = Ta с высотой и считать T T0, при этом вертикальный градиент температуры по прежнему учитывается в членах, описывающих температурную стратификацию среды. В то же время плотностная стратификация атмосферы обусловлена весовой сжимаемостью газа, следовательно, ее учет необходим даже в случае безразличной температурной стратификации (J = 0).

Глава 1. Автомодельный термик Интегрирование уравнения энергии (1.6) по пространству показывает, что полная тепловая энергия, запасенная в термике, сохраняется постоянной на начальной и автомодельной стадиях движения термика, пока влияние стратификации атмосферы незначительно:

Для целей обезразмеривания удобно наряду с тепловой энергии термика Q0 ввести полный запас плавучести термика B0, определив его как где = gT — ускорение силы плавучести. Определение B0 в (1.10) отличается от используемого в приближении Буссинеска [36, 41] тем, что в подынтегральное выражение входит функция (z), описывающая распределение плотности среды по высоте. Кроме того, введем числа Грасгофа и Прандтля, основанные на турбулентных значениях коэффициентов переноса:

Таким образом, искомые функции u, v и зависят от времени t, координат r, z, величин, характеризующих состояние атмосферы, свойства газа 0, T0, J,, CP, g, коэффициентов переноса µ, и параметров, описывающих начальное состояние термика — B0, R0. Воспользовавшись -теоремой, запишем эту зависимость в виде связи между безразмерными комплексами где n = 1/2, если f = u, v и n = 3/2 если f =, величина B0 /gR0 характеризует начальный перегрев в термике. Здесь и далее тильдой обозначены безразмерные величины.

Глава 1. Автомодельный термик Выражение (1.12) содержит три характерных времени ti, ts и td. Время ti определяет продолжительность начальной стадии, на которой формируется вихревой тор.

На этой стадии существенную роль играет параметр B0 /gR0, характеризующий начальное распределение плавучести по объему термика. Как показано в [16], продолжительность начальной стадии составляет ti 3R0 /B0, при t ti решение перестает зависеть от этого параметра. Величина ts = |J|1/2 определяет характерное время, через которое начинает сказываться температурная стратификация атмосферы: если J > 0 (устойчивая стратификация), то при t > ts происходит зависание термика, в противном случае — его разгон. Наконец, величина td представляет собой характерное время, за которое термик поднимается на высоту, где плотность окружающей среды существенно отличается от плотности на уровне виртуального источника 0, так что при t > td атмосфера не может считаться однородной. Время td может быть оценено как время, за которое верхняя кромка плавучего облака проходит расстояние, равное характерному вертикальному масштабу изменения атмосферной плотности La (см. уравнение (1.8)). Для термика, поднимающегося в несжимаемой атмосфере, как показано в [6, 41], на автомодельном участке подъема координата верхней кромки термика растет со временем как Следовательно, изменение окружающей плотности с высотой становится существенным по прошествии времени Из трех рассмотренных выше характерных времен величина ts зависит только от параметров атмосферы, тогда как два других времени ti и td зависят от параметров термика. В зависимости от состояния атмосферы и начальных параметров термика возможны различные соотношения между указанными временными масштабами, но, как правило, ti ts, td. Если запас плавучести термика относительно мал, так что td ts, то в течение всего периода подъема облака по закону zt t1/2 плотГлава 1. Автомодельный термик ность окружающей среды меняется незначительно и, следовательно, можно учитывать только температурную стратификацию атмосферы. Именно этот случай ранее подробно рассматривался в литературе на основе приближения Буссинеска. Можно сказать, что в [16, 40, 41] построено решение, отвечающее промежуточной асимптотике ti t ts td. Отметим, что если J > 0, то термик зависает на определенной высоте и решение остается справедливым и при t > ts. Если же J 0, термик поднимается неограниченно высоко и решение, полученное в предположении о несжимаемости среды, справедливо только при t td.

В настоящей работе рассматривается другой предельный случай, когда ti td ts. Типичная зависимость координаты верхней кромки термика, всплывающего в устойчиво стратифицированной атмосфере, схематически представлена на рис. 1.1, где отмечены основные стадии процесса: временной интервал 0 t ti соответствует начальной стадии (1), интервал ti t ts отвечает корневому закону подъема zt t1/2, причем на участке ti t td (2) справедливо автомодельное решение, описывающее подъем термика в несжимаемой атмосфере, а на участке td t ts (3) становится существенной плотностная стратификация среды. Наконец, при t ts (4) проявляется температурная стратификация атмосферы, приводящая (в случае устойчивой стратификации среды) к зависанию термика.

Представленное ниже приближенное решение, отвечающее промежуточной асимптотике ti t ts, при t td (стадия 2 на рис. 1.1) совпадает с решением для несжимаемой среды [16, 40, 41], а на стадии 3 описывает поведение термика в условиях существенного влияния плотностной стратификации атмосферы. При выводе определяющих уравнений в силу условия t ts опустим последний член в правой части уравнения притока тепла (1.6), содержащий параметр стратификации Для того, чтобы описать структуру термика на стадии, отвечающей корневому закону подъема, введем зависящие от времени масштабы длины L, скорости U и Глава 1. Автомодельный термик Рис. 1.1. Основные стадии подъема термика в атмосфере с переменной плотностью: 1 — начальная стадия, 2 — стадия, на которой атмосфера может считаться несжимаемой, 3 — стадия, на которой становится существенной плотностная стратификация среды, 4 — стадия зависания (при устойчивой стратификации атмосферы).

Глава 1. Автомодельный термик температуры, определив их как после чего сделаем замену переменных Такая (либо эквивалентная) замена переменных использовалась ранее при исследовании автомодельных термиков в несжимаемой среде [16, 37].

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

1.3. Численное моделирование термика в экспоненциальной изотермической атмосфере Рассмотрим подъем турбулентного осесимметричного термика в изотермической (Ta = const = T0 ) атмосфере, для которой вертикальное распределение плотности является экспоненциальным с характерным линейным масштабом La = Rg T0 /g, так что (z) = exp(z/La ).

Введем завихренность и функцию тока, определив их как Для исключения давления из уравнений движения применим операцию ротора к уравнениям (1.4), (1.5). При переходе к безразмерным переменным (1.16) система уравнений (1.3)–(1.6) сводится к виду Глава 1. Автомодельный термик где безразмерная завихренность и функция тока введены как = L /U, = /UL2, а распределение плотности атмосферы представлено в виде = exp( Gr1/4z ). Распределение плотности атмосферы зависит от безразмерного времени, поскольку время входит в определение характерного масштаба длины L (см. (1.15)).

Граничные условия на оси симметрии имеют вид = = / r = 0, на бесконечности = = = 0. Интегральное соотношение, выражающее постоянство тепловой энергии термика, в новых переменных приобретает вид Начальные условия при = 0 соответствуют автомодельному термику, всплывающему в несжимаемой атмосфере с постоянной плотностью. Система уравнений (1.17)–(1.19) решалась численно в прямоугольной области 0 r 4, 2 z на сетке, содержащей 50 100 узлов, интегрирование по времени производилось с шагом = 0,05. На каждом временном шаге распределения искомых величин находились путем глобальных итераций, осуществляемых методом последовательной верхней релаксации с прогонкой в горизонтальном и вертикальном направлениях.

Чтобы оценить влияние внешних границ расчетной области на получаемые решения, результаты расчетов сравнивались с данными, полученными в областях большеГлава 1. Автомодельный термик го размера 0 r 6, 3 z 9 и 0 r 8, 4 z 14. Изменение расчетной области отражалось в основном на положении внешних линий тока, при этом как структура термика, так и его положение (в частности, координата верхней кромки) менялись крайне незначительно. Так, разность координат верхней кромки zt, полученных в основной и промежуточной областях, составляла всего 1%, а в основной и наибольшей областях — 2%. Максимальные значения избыточной температуры отличались от полученных на основной сетке на 2% при решении в промежуточной и 2,6% при решении в наибольшей областях, максимальные значения функции тока отличались в пределах 2,5 и 5%. Таким образом, размеры расчетной области были достаточными чтобы свести к минимуму влияние свободных границ.

Рассмотрим теперь результаты, полученные при числе Грасгофа Gr = 400 и числе Прандтля Pr = 1. При таких параметрах автомодельная координата верхней кромки термика совпадает с наблюдаемой в экспериментах [17]. На рис. 1.2 представлена рассчитанная структура термика в момент времени = 0, соответствующий подъему автомодельного термика в несжимаемой атмосфере. Сплошные изолинии, соответствующие безразмерной избыточной температуры, построены на уровнях 0,1, 0,2,... 0,9 от максимальной температуры max = 4,62. Штриховыми линиями представлены изолинии функции тока (0,1, 0,2,... 0,9 от максимального значения max = 0,36). В правой части рисунка показано вертикальное распределение плотности атмосферы (). Автомодельная координата верхней кромки облака, опредеz ленная как высота, на которой избыточная температура на оси симметрии составляет 10% от максимального значения max, равна zt = 4,26, что близко к экспериментальному диапазону 4,3 4,4 установленному в [17] (см. уравнение(1.1)).

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

Структура плавучего облака и вертикальное распределение плотности атмосферы при = 1 показаны на рис. 1.3 (max = 4,82, max = 0,30, zt = 4,34). Аналогичные распределения при = 2 приведены на рис. 1.4 (max = 4,59, max = 0,26, zt = 4,53). К этому моменту плотность атмосферы на уровне верхней кромки облака Глава 1. Автомодельный термик в 10 раз меньше плотности на уровне виртуального источника, а расстояние, пройденное термиком, равно приблизительно 2,5La.

Сравнение рис. 1.2–1.4 показывает, что структура плавучего облака испытывает с течением времени определенные изменения: относительная ширина термика возрастает, верхняя кромка облака становится менее резкой. Вместе с тем, наиболее важным свойством полученных решений является то, что автомодельная координата верхней кромки термика zt и максимальная избыточная температура max меняются весьма слабо. Чтобы показать это более наглядно, на рис. 1.5 представлены временные зависимости автомодельной координаты верхней кромки термика zt ( ) (кривая 1) и максимальной избыточной температуры max ( ) (кривая 2), отнесенных к своим начальным значениям zt = 4,26, max = 4,62 соответственно. Штриховой линией нанесена зависимость от времени плотности атмосферы на уровне верхней кромки облака (t ( )). Можно заключить, что при изменении плотности атмосферы на порядок величины положение и параметры термика в автомодельных координатах остаются практически неизменными. Это означает, что в физических координатах термик поднимается в соответствии с корневым законом zt t1/2, а максимальная температура в облаке падает со временем как max t3/2, т. е., в соответствии с законами, которые дает анализ размерностей в случае автомодельного движения термика в несжимаемой атмосфере.

Глава 1. Автомодельный термик Рис. 1.2. Структура термика в автомодельных координатах при = 0. Справа приведено вертикальное распределение плотности атмосферы.

Глава 1. Автомодельный термик Рис. 1.3. Структура термика в автомодельных координатах при = 1.

Глава 1. Автомодельный термик Рис. 1.4. Структура термика в автомодельных координатах при = 2.

Глава 1. Автомодельный термик Рис. 1.5. Изменение автомодельной координаты верхней кромки термика zt (кривая 1) и максимальной избыточной температуры в термике (кривая 2) с ростом времени.

Функции отнесены к своим значениям при = 0 z = 4,26 и max = 4,62). Штриховая кривая — плотность атмосферы на уровне верхней кромки ( zt ).

Глава 1. Автомодельный термик 1.4. Приближенное аналитическое решение При больших числах Грасгофа и при Pr = 1 возможно найти аналитическое решение, описывающее термик, поднимающийся в атмосфере с переменной плотностью на стадии, отвечающей корневому закону роста координаты верхней кромки. Следуя работам [16, 40], будем использовать приближение вертикального погранслоя и пренебрежем градиентом давления p/z в уравнении для вертикальной компоненты количества движения (1.5). Это позволяет решать только уравнения неразрывности, вертикального импульса и температуры. Как и ранее, выполним переход к автомодельным координатам (1.16). Приведенные в предыдущем разделе результаты численного моделирования показывают, что в автомодельных координатах структура термика меняется со временем весьма медленно, поэтому опустим также нестационарные члены и будем искать решение, зависящее от только параметрически. В рамках данных упрощений уравнения неразрывности, импульса и энергии принимают вид Будем искать решение уравнений (1.21)–(1.23), предполагая профили вертикальной компоненты скорости v и избыточной температуры подобными и описывающимися функцией Глава 1. Автомодельный термик где профиль () описывает вертикальные распределения скорости и избыточной температуры в термике. В отличие от [16], под аргументом экспоненты введена пока неизвестная функция a(). Подстановка профилей (1.24) в уравнения (1.22), (1.23) показывает, что подобие возможно лишь при Pr = 1, это значение числа Прандтля, близкое к наблюдаемым в экспериментах [56], используется в дальнейшем.

Интегрирование уравнения (1.24) дает функцию тока = Gr1/4 1 exp(a2 ), после чего дифференцированием функции тока находится радиальная компонента скорости u = ()1 / z. После подстановки профилей (1.24) оба уравнения (1.22) и (1.23) приобретают одинаковый вид и могут быть записаны в виде уравнения для функции () (штрихом обозначается дифференцирование по z ):

Можно видеть, что если выбрать неизвестную функцию в горизонтальном профиле вертикальной скорости и избыточной температуры (1.24) в виде a() = Gr1/2 ()/4, второй и третий члены в фигурных скобках обращаются в нуль. Такой вид функции a() и будем использовать в дальнейшем. Осредним затем оставшиеся коэффициz енты данного уравнения, домножив каждый член на r d и проинтегрировав в предеr лах от нуля до бесконечности. Кроме того, введем новую независимую переменную = Gr1/4z. В результате получим следующее приближенное дифференциальное уравнение для функции (здесь и далее штрихом обозначается производная по ):

где аргумент функции (z/H), описывающей вертикальное распределение плотности атмосферы в новых переменных имеет вид ( ). Интегральное условие (1.9), или Глава 1. Автомодельный термик эквивалентное ему условие (1.20), приобретает вид Общее решение уравнения (1.25), удовлетворяющее граничным условиям = при = ±, имеет вид где Константа интегрирования 0 в уравнении (1.25) зависит от двух параметров — числа Грасгофа Gr и безразмерного времени. Значение 0 можно найти, подставляя решение (1.27) в интегральное условие (1.26).

Чтобы проанализировать решение (1.27) при произвольном распределении плотности атмосферы, проинтегрируем его один раз с учетом граничных условий = 0 при = ±. В результате получим Изоклины нуля этого уравнения имеют вид 1 = 0 и 2 = 2 + 4Gr1(1), точку их пересечения обозначим. Мы предполагаем, что функция () имеет единственный максимум в точке = max, в которой максимальны вертикальная скорость v и избыточная температура. Достаточное условие единственности точки максимума состоит в том, чтобы при < изоклина 2 проходила ниже 1 (т. е., 2 + 4Gr1(1) < 0), тогда как при > функция 2 возрастала бы монотонно (т. е., 2 + 4Gr1(1) > 0). Отметим, что этим условиям удовлетворяет широкий класс функций, в частности, они справедливы для экспоненциальной атмосферы, использованной выше при численном моделировании термика в среде с переменной плотностью.

Глава 1. Автомодельный термик Координата точки максимума max может быть найдена из условия пересечения функции () с изоклиной нуля 2, в этом случае выражение, стоящее в уравнении (1.29) в скобках обращается в нуль. Используя определение функции (см. уравнение (1.28)), представим условие максимума в виде = 4Gr1, или, окончательно, Разложим теперь функцию () в асимптотический ряд, воспользовавшись интегрированием по частям:

Судя по аналитическому решению, найденному для случая постоянной плотности атмосферы в [16, 40], а также по результатам численного моделирования, представленным выше, можно ожидать, что максимум функции () достигается при max = O(1). Вблизи максимума = Gr · O(1) — см. (1.28) поэтому при больших числах Грасгофа можно ограничиться лишь двумя первыми членами разложения, опуская члены порядка Gr4 или выше. В результате условие максимума (1.30) приобретает вид Интегрирование уравнения (1.26) для вычисления константы 0, входящей в решение (1.27), может быть произведено аналитически только в простейшем случае постоянной плотности атмосферы ( = 1). Соответствующие решения были построены и проанализированы в работах [16, 40], где показано, что функция () возрастает с высотой, достигая максимума в точке с координатой (в принятых здесь обозначениях) Глава 1. Автомодельный термик после чего резко падает практически до нуля. Такое поведение функции () позволяет аппроксимировать ее разрывной функцией = 2 при 0 1, = 0 при При произвольном распределении плотности атмосферы вычислить интеграл (1.26) и получить явные выражения для координаты максимума в общем случае не удается. Ниже рассмотрен частный случай распределения плотности атмосферы, в котором возможно получение явных окончательных формул. Именно, рассматривается функция распределения плотности атмосферы Эта функция удовлетворяет условию единственности максимума решения (), обсуждавшемуся выше. Она может рассматриваться в качестве модельного распределения, задающего при 0 спадание плотности атмосферы с высотой. Хотя уменьшение плотности при < 0 нефизично, эта область оказывает слабое влияние на решение, поскольку она охватывает лишь незначительную часть термика.

Рассмотрим для начала решения уравнения (1.25), полученные путем численного интегрирования. На рис. 1.6 представлены интегральные кривые, полученные при = 0, 1, 2, 3 и 4 (справа), а также соответствующие распределения плотности атмосферы () в те же моменты времени (слева). Отметим снова, что решение при = в точности соответствует подъему термика в несжимаемой среде, либо на той стадии, когда весовая сжимаемость атмосферы еще не проявляется. С продвижением термика во все более разреженные слои атмосферы, как видно из рис. 1.6, автомодельная координата верхней кромки облака остается практически неизменной (все кривые, соответствующие разным моментам времени пересекаются в одной точке).

Точка максимума функции () слегка сдвигается вниз, тогда как горизонтальный размер термика увеличивается как (max)1/2. Эти результаты находятся в хорошем соответствии с представленными выше результатами численного моделирования.

Наконец, рассмотрим результаты, полученные аналитически. Подставив распределение плотности (1.33) в уравнение (1.28), получим в явном виде вспомогательные Глава 1. Автомодельный термик Рис. 1.6. Изменение структуры термика, поднимающегося в атмосфере с плотностью, убывающей обратно пропорционально квадрату высоты. Профили () показаны при = 0, 1, 2, 3 и 4. В левой части рисунка приведены распределения плотности атмосферы в те же моменты времени.

Глава 1. Автомодельный термик функции и :

() = 1 + где B(p, q) — бета-функция, Bx (p, q) — неполная бета-функция [57], аргументом которой служит значение () = (1 + ( )2)1.

Чтобы найти константу 0, представим интеграл (1.26) в виде двух членов 1 и Первый из двух вспомогательных интегралов может быть вычислен точно, он выражается через бета-функцию B(p, q) как Чтобы вычислить интеграл 2, заменим функцию в подинтегральном выражении некоторой аппроксимирующей функцией. Будем искать асимптотику решения () при < max используя пилотную функцию () = 2(1+c 2), где c — константа, подлежащая определению. Подстановка этой функции в уравнение (1.29) показывает, что асимптотическое равенство левой и правой частей возможно только при c = 6/Gr. Поэтому аппроксимируем () разрывной кусочно линейной функцией где координата точки максимума max = (1 + c 2 )1/2 определяется интегральным соотношением (1.26). Тогда второй вспомогательный интеграл оценивается как Глава 1. Автомодельный термик Суммируя полученные значения интегралов 1 и 2, получим из интегрального соотношения (1.26) следующее выражение для константы 0 :

Подставляя это значение в условие на точку максимума (1.31), можно найти координату max. При этом удобно представить координату точки максимума в виде max = max max, где max — координата точки максимума при = 0 (т. е., в случае однородной по плотности атмосферы), а max — ее изменение со временем, которое может быть интерпретировано как «размытие» верхней кромки термика при продвижении в разреженные слои атмосферы. Вычисления дают следующую зависимость:

Функция max ( ) представлена на рис. 1.7 сплошной кривой, здесь же точками нанесены результаты, полученные при численном решении уравнения (1.25). Можно видеть, что полученное приближенное решение дает разность координат точки максимума с точностью не хуже 25–30%. Еще более высокой точности можно достичь если при аппроксимации функции () в уравнении (1.34) соответствующим образом подобрать наклон прямой, т. е., коэффициент c. Наилучшее совпадение с результатами численных расчетов достигается при c 8,8/Gr, соответствующая зависимость размытия верхней кромки от времени max ( ) отличается от (1.35) лишь коэффициентом во втором члене. Эта зависимость представлена на рис. 1.7 штриховой линией, отличие от результатов численных расчетов не превышает 5%.

Глава 1. Автомодельный термик Рис. 1.7. Изменение координаты точки максимума max ( ). Точки соответствуют численному решению, сплошной линией показано аналитическое решение (1.35), штриховая кривая отвечает решению (1.36).

Глава 1. Автомодельный термик 1.5. Выводы В настоящей главе проведено исследование подъема турбулентного осесимметричного термика в атмосфере с зависящей от высоты плотностью. В основе используемой модели лежит система уравнений «глубокой» конвекции, в которой пренебрегается динамической сжимаемостью газа, но принята во внимание гидростатическая сжимаемость атмосферы. Численные и аналитические решения получены для диапазона параметров, в котором изменения плотности среды с высотой проявляются во время подъема облака по корневому закону, т. е., до начала стадии зависания.

Численное интегрирование системы уравнений в переменных вихрь-функция тока проведено для изотермической атмосферы с экспоненциально спадающей плотностью. Расчеты проведены на интервале времени, в течение которого термик поднимается на высоты, где плотность атмосферы на порядок величины меньше плотности на высоте образования облака. Результаты расчетов показали, что при проникновении термика в более разреженные слои атмосферы уменьшение окружающей плотности вызывает дополнительное размытие верхней кромки и увеличение горизонтального размера термика. При этом, однако, автомодельная координата верхней кромки термика и максимальное значение безразмерной избыточной температуры меняются весьма слабо. Следовательно, законы роста координаты кромки облака zt t1/ и затухания максимальной избыточной температуры max t3/2, установленные с помощью анализа размерностей для несжимаемой среды, с хорошей точностью выполняются и для термика в сжимаемой атмосфере.

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

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

Глава Подъем термика и перенос дисперсных примесей в атмосфере 2.1. Математическое моделирование термиков на основе уравнений сжимаемого газа. Прикладные задачи В последние годы наблюдается быстрый прогресс в развитии методов численного моделирования течений жидкостей и газов применительно к задачам газовой динамики, гидродинамики, конвекции, механики многофазных и химически реагирующих сред [58–70]. Численное моделирование на основе одним из основных методов исследования плавучих облаков. Решение многомерных уравнений гидродинамики позволяет получать исследовать течения в тех диапазонах параметров, где неприменимы аналитические решения, рассмотренные в Главе 1. Приближение Буссинеска, традиционно используемое для описания течений несжимаемой жидкости в задачах естественной конвекции [64], в применении к течениям газа обладает весьма существенными ограничениями, прежде всего связанными с условием малости отклонений плотности газа от окружающей. Эти ограничения не выполняются, например, на начальной стадии эволюции высокотемпературного термика, а также для течений с большим вертикальным масштабом, когда существенную роль играет гидростатическая сжимаемость атмосферы. По мере развития численных методов и увеличения вычислительных мощностей все большее распространение получает моделирование Глава 2. Подъем термика и перенос примесей на основе полных уравнений динамики сжимаемого вязкого теплопроводного газа, позволяющих изучать течения со значительными перепадами плотности [66].

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

Образование вихревого кольца при подъеме осесимметричного термика изучалось численно в работах [72, 73], при этом использовалась область с подвижными (расширяющимися со временем) границами. Отмечалось, что сферический вначале объем нагретого газа через 16 с сворачивается в вихревое кольцо, а примерно через 30–40 с максимум температуры отходит от оси на периферию. Для некоторых начальных параметров удалось получить хорошее согласие с результатами работы [19].

Структура и динамика подъема приповерхностного термика, касающегося в начальный момент подстилающей поверхности, изучалась в аналогичной постановке в работе [74]. Использовались постоянные турбулентные коэффициенты переноса. Были исследованы два типа граничных условий на подстилающей поверхности — условия «прилипания» (равенство нулю тангенциальной компоненты скорости) и условие свободного проскальзывания (нулевой градиент тангенциальной скорости по нормали к поверхности). Качественное сравнение результатов расчетов с экспериментальными данными [30, 31, 34, 75] показало, что форма термика оказывается ближе к наблюдаемой в опытах при использовании условий прилипания.

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

Глава 2. Подъем термика и перенос примесей В работах [77, 78] численные расчеты проведены для различных соотношений между размерами термика и характерной высотой неоднородной атмосферы, а также при переменных коэффициентах обмена. Отмечено, что в зависимости от масштаба термика возможны качественно разные структуры вихревого течения. Рассмотрено также образование вихревого кольца при различной начальной геометрии термика.

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

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

Переменные коэффициенты турбулентного переноса, вычисляемые согласно k модели турбулентности, а также алгебраической модели турбулентности [82], использовались в расчетах крупномасштабных турбулентных термиков в [83], а для расчета приповерхностных термиков — в [84]. Исследование зависимости решений от начальных значений турбулентной вязкости и скорости ее диссипации показало сильную чувствительность модели к этим параметрам.

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

Так, в работах [85–87] рассмотрен дрейф термиков в сдвиговом поле ветра. Подъем и взаимодействие пары термиков в атмосфере рассмотрены в [88], где исследован процесс совместного всплывания двух соосных и находящихся первоначально на одинаковой высоте термиков. В первом случае показано, что верхний термик всплывает практически как уединенный, по закону t1/2, тогда как нижний движется как t2, что является следствием воздействия движения среды, наведенного первым термиком.

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

Для описания медленных конвективных течений с большим вертикальным масштабом в работе [91] использована гипозвуковая модель [66], хорошо зарекомендовавшая себя ранее при исследовании течений с малыми числами Маха. Чтобы учесть переменность плотности атмосферы с высотой, модель [66] была обобщена на случай произвольных значений параметра гидростатической сжимаемости. С точки зрения полноты описания течений эта модель занимает промежуточное место между полными уравнениями Навье-Стокса и классическим приближением Буссинеска. На основе модели в работе [91] исследовано всплытие в атмосфере системы термиков, выявлены эффекты взаимодействия плавучих облаков, рассмотрено изменение стратификации атмосферы, возникающее при подъеме бесконечной системы термиков.

Наряду с изучением свойств отдельных термиков, важное значение имеет решение прикладных задач, составляющей частью которых являются плавучие облака. В связи с экологическими проблемами большой интерес представляет изучение переноса всплывающим термиком мелкодисперсных частиц. Массовое загрязнение Глава 2. Подъем термика и перенос примесей верхних слоев атмосферы оптически активными частицами, экранирующими солнечное излучение, может привести к климатическим изменениям глобального масштаба [92–96], поэтому математическое моделирование подобных явлений представляется весьма актуальным.

Численное моделирование переноса всплывающей дисперсной примеси всплывающим термиком было проведено в работе [97], где вихревое кольцо описывалось на основе предложенного в [19] подхода, а примесь считалась пассивной (не оказывающей воздействия на движение газа) и описывалась на основе лагранжева подхода.

Результаты расчетов позволили, по крайней мере качественно, воспроизвести процесс формирования под облаком пылевой колонки, образованной выпадающими из термика частицами.

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

В работах [16, 50], наряду с упоминавшимся выше моделированием отдельных стадий эволюции термика на основе уравнений несжимаемой жидкости, исследован перенос всплывающим облаком пассивной примеси. Определена доля примеси, попадающая в стратосферу в результате подъема термика, а также построено распределение линейной плотности примеси по высоте в зависшем облаке.

В [102] путем численного решения уравнений сжимаемого газа изучен вынос горячим термиком аэрозольной примеси, находящейся в тонком цилиндрическом приземном слое. Отличительной чертой является использование модели «активной»

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

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

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

Глава 2. Подъем термика и перенос примесей 2.2. Уравнения движения термика с дисперсной примесью Пусть в начальный момент t = 0 в результате быстрого выделения тепла Q0 над плоской горизонтальной поверхностью образовалось неподвижное сферическое облако горячего газа радиусом R0, содержащее распределенные по объему мелкодисперсные частицы с общей массой Mp. Под действием сил плавучести термик начинает всплывать, увлекая дисперсную примесь, после зависания облака в атмосфере формируется распределение частиц, параметры которого важны для оценки экологических последствий загрязнения атмосферы.

Эволюция запыленного теплового термика изучается в предположениях, позволяющих рассматривать дисперсную двухфазную систему как совокупность двух взаимодействующих взаимопроникающих сплошных сред и использовать соответствующее математическое описание [104–106]. Такой подход правомерен если средние размеры частиц и расстояния между ними намного больше характерных молекулярно-кинетических размеров, но существенно меньше характерных размеров неоднородностей макроскопических параметров несущей среды. Кроме того, пренебрегается столкновениями частиц, их дроблением и испарением, объемная доля частиц и отношение истинных плотностей газа и примеси считаются малыми.

Будем рассматривать достаточно мелкие частицы, которые могут эффективно захватываться всплывающим облаком и увлекаться на большие высоты. Именно, частицы считаются мелкими настолько, что время их скоростной и температурной релаксации намного меньше времени развития конвекции (R0 /g)1/2, где g — ускорение силы тяжести. В рамках стоксова закона сопротивления время скоростной релаксакции одиночной частицы порядка d2 0 /18µl, где d — диаметр частиp цы, 0 — истинная плотность частицы, µl — молекулярная вязкость газа. Время скоростной релаксации оказывается сравнимым со временем развития конвекции при d2 18µl (R0/g)1/2 /0. Для крупномасштабных облаков (R0 103 м) и при µl = 2 · 105 кг/м2 с, 0 = 103 кг/м3 из этого условия получим диаметр частиц d 1 мм. Аналогичную оценку можно получить, сравнивая скорость витания одиГлава 2. Подъем термика и перенос примесей ночной частицы с характерной конвективной скоростью (R0 g)1/2. Для частиц, диаметр которых существенно меньше полученной выше величины (d 10 100 мкм) можно принять гипотезу мгновенной релаксации и считать скорости и температуры фаз равными [105, 106].

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

Здесь индекс p относится к дисперсной фазе, индекс m — к двухфазной среде в целом, P,, T — давление, плотность и температура, Rg = R/m — газовая постоянная, R = 8,31 Дж/моль·К — универсальная газовая постоянная, m — молекулярная масса газа, CV, Cp — теплоемкости газа при постоянном объеме и примеси, Yp — массовая доля дисперсной примеси в двухфазной смеси, µ, — постоянные эффективные коэффициенты турбулентной вязкости и теплопроводности, D — коэффициент турбулентной диффузии (предполагается, что D = const). Как показано в Главе 1, использование постоянных по пространству коэффициентов турбулентной вязкости и теплопроводности, несмотря на свою простоту, позволяет в расчетах воспроизвести важную особенность движения термика в атмосфере с переменной плотностью — наличие стадии подъема, на которой координата верхней кромки облака растет по корневому закону (пропорционально t1/2 ). Коэффициент пропорциональности в законе роста зависит от величины турбулентной вязкости, поэтому, выбирая надлежащим образом эффективные коэффициенты переноса, можно добиться в Глава 2. Подъем термика и перенос примесей расчетах правильного закона подъема термика. Более подробно процедура выбора коэффициентов переноса описана в последующих разделах.

Уравнение притока тепла (2.3) показывает, что удельная объемная теплоемкость смеси увеличивается по сравнению с чисто газовой (добавляется слагаемое p Cp ).

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

Распределение температуры в атмосфере описывается моделью международной стандартной атмосферы [96, 107], распределения невозмущенной плотности и давления находятся из уравнений состояния и гидростатического равновесия:

Рассматривается двухслойная модель атмосферы, в которой параметр стратификации J на высоте тропопаузы HT скачком меняет свое значение от J1 = 1,2 · 104 с2, характерного для тропосферы (0 z HT ), до величины J2 = 4,8 · 104 с2, типичной для стратосферы (z > HT ). Высота тропопаузы в зависимости от времени года составляет 10–16 км.

Начальные и граничные условия задаются следующим образом:

t=0:

Глава 2. Подъем термика и перенос примесей Введем безразмерные переменные, принимая в качестве масштабов давления, плотности и температуры их значения в невозмущенной атмосфере на уровне подстилающей поверхности: P = P0 = Pa (0), = 0 = a (0), T = T0 = Ta (0).

За пространственный масштаб можно выбрать либо начальный радиус термика R 0, либо характерный вертикальный масштаб атмосферы Rg T0 /g (как это сделано в Главе 1), что более предпочтительно, поскольку в дальнейшем варьируются только параметры термика, а состояние атмосферы остается неизменным. Определенное неудобство состоит в том, что для реальных термиков R0 Rg T0 /g и начальный безразмерный радиус облака оказывается много меньше единицы. Поэтому ниже вводится фиксированный масштаб длины L порядка начального размера термика R0, а все окончательные результаты приводятся к виду, не зависящему от конкретного значения L. Масштабы скорости и времени определяются как U = (L g)1/2, t = (L /g)1/2. При введении фиксированного масштаба длины L общее число определяющих параметров оказывается на единицу большим, чем в том случае, когда за линейный масштаб принимается характерный вертикальный масштаб изменения атмосферной плотности Rg T0 /g (в качестве «дополнительного» параметра выступает отношение L к Rg T0 /g). Однако в силу того, что ниже состояние атмосферы (а с ним и величина Rg T0 /g) фиксировано, данное отношение также не варьируется и появления дополнительного параметра не происходит.

В безразмерных переменных, обозначаемых тильдой, система уравнений (2.1)– (2.3) с начальными и граничными условиями (2.6)–(2.9) принимает вид Глава 2. Подъем термика и перенос примесей где p = Cp /CV — отношение теплоемкостей дисперсной и газовой фаз, g = (0, 1) — безразмерный вектор силы тяжести. В уравнения (2.10)–(2.17) входят безразмерные числа Рейнольдса Re = 0 L U/µ, Прандтля Pr = /CP µ, Шмидта Sc = /µ, и Маха Ma2 = U /Rg T0. Вертикальные распределения параметров невозмущенной атмосферы (2.5) в безразмерных переменных принимают вид где J = JL /g — безразмерный параметр стратификации атмосферы. Расчеты проводились при фиксированных параметрах что соответствует реальным параметрам атмосферы при следующих значениях масштабов:

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