Вычислительные технологии
Том 7, № 2, 2002
ПОДОБИЕ СТРУКТУР КОНВЕКТИВНЫХ ПОГРАНИЧНЫХ СЛОЕВ АТМОСФЕРЫ
И ВОДОЕМА: РЕЗУЛЬТАТЫ ЧИСЛЕННЫХ ЭКСПЕРИМЕНТОВ С ВИХРЕРАЗРЕШАЮЩИМИ МОДЕЛЯМИ *
В. А. Шлычков, П.Ю. Пушистов Институт водных и экологических проблем СО РАН Новосибирск, Россия e-mail: slav@ad-sbras.nsc.ru
The possibility of the simultaneous development of turbulent penetrative convection in the atmosphere and the lake during the night-time cooling of a near-surface water layer is shown. The numerical investigation of the hydrodynamics equations has been made on the basis of large-eddy simulation approach. The results of calculations of convective processes in both natural media obtained by usage of observed data are presented. The identity of the spatial structure and the energy transformations in the boundary layers of the atmosphere and the lake is shown.
Учитывая экспериментально установленные факты подобия в структурах пограничных слоев водоема и атмосферы [1, 2], а также на основе результатов анализа составляющих теплового баланса в приводном слое воздуха [2, 3] можно предположить высокую степень вероятности возникновения и развития проникающей турбулентной конвекции в нижней части атмосферного пограничного слоя (АПС) над водоемом в ночной период. Явление сопряженной конвекции в АПС и верхнем перемешанном слое водоема (ВПС) относится к числу малоизученных. Вместе с тем имеется обширная литература, описывающая результаты наблюдений и моделирования проникающей конвекции в ВПС при ночном охлаждении поверхности озер и океанов [3 - 5]. То же относится к явлению формирования конвективного АПС в случае адвекции холодной воздушной массы с суши на относительно теплую поверхность водоема [6, 7].
В силу сложности явления сопряженной конвекции возможности инструментальной идентификации ограничены; по этой же причине практически отсутствуют и данные наблюдений. В настоящей работе предпринята попытка численного моделирования ансамблей проникающей конвекции в ВПС и АПС при ночном выхолаживании с параметрическим учетом горизонтального переноса холодного воздуха с береговой зоны на теплый водоем (ночной бриз). Для этих целей используется базовое предположение о двойной
* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований, грант №99-05-64735.
© В. А. Шлычков, П. Ю. Пушистов, 2002.
структуре поля турбулентности, причем крупновихревая составляющая спектра воспроизводится явно, а мелкомасштабная параметризуется. В западной литературе такой метод получил название LES-метода (Large Eddy Simulation) [6, 8].
1. Уравнения модели атмосферы
Согласно [9] в конвективном АПС можно выделить три основных типа разномасштабных процессов: 1) упорядоченные средние движения; 2) движение когерентных структур ("крупные вихри") с характерным вертикальным размером, сравнимым с толщиной пограничного слоя; 3) изотропная мелкомасштабная турбулентность.
Введем прямоугольную декартову систему координат (х,у,г), в которой ось г направлена вертикально вверх. Представим искомую вектор-функцию ф = (и^,и,в,п), где — компоненты вектора скорости вдоль осей х, у, г; в — потенциальная температура; п — аналог давления, в виде суммы
ф =Ф + фф. (1)
Здесь функция Ф(г,Т) = (и,У, 0, в, П) описывает детерминированное осредненное течение, а ф' = (и', V1, и', в', п') — квазиупорядоченное поле турбулентных крупных вихрей конвективной природы.
Определим Ьх,Ьу как горизонтальные размеры области, в которой развивается проникающая конвекция в виде ансамбля спонтанно образующихся термиков, и примем на боковых границах области предположение о периодичности процессов вдоль х, у. Подставляя представление (1) в уравнения мезомасштабной гидротермодинамики атмосферы [10], проведем их осреднение в горизонтальной плоскости по правилам
ф7 =0, ф = Ф, (2)
Ьх Ьу - 1 1 [ [
где f =--(f) dx dy. В результате получим систему уравнений, описывающую
Lx ^у } }
о о
горизонтально однородное течение в АПС:
ди „„ „ ч д ди д_
дУ" д дУ д_
Ж = "'(и - и") + ТгкТ7 - Шт>'
дв . , д ,, дв д -¡г- ,„,
Ж + А(К'© = ТгКтЖ - Тгви- (3)
где ис, Ус — составляющие геострофического ветра, обусловленные барическим градиентом в свободной атмосфере; I — параметр Кориолиса; К — коэффициент вертикального турбулентного обмена подсеточного масштаба; Кт = К/Рг, Рг — число Прандтля в АПС; Adv© — функция, описывающая в параметрической форме адвекцию температуры. Здесь и ниже штрихи у конвективных отклонений опущены.
Систему уравнений для описания конвективных процессов в АПС получим путем покомпонентного вычитания (3) из соответствующих уравнений исходной системы:
du dU дп , ^ д тлди д — + w— = - — + lv + DXy u + — К— + — dt дг дх дг дг дг
1V дУ дп д тдг д — + т— = — ---1и + ПХу V + — К— + —
1Ь 'дг ду дг 'дг дг
1т дж л „ ^ д тгдт
-17 = + Х9 + °хут + 7Г Ктг>
аЬ дг дг дг
1в д0 п л д тдв д —
— + т— = Вху е + — Кт— + — те, аЬ дг дг дг дг
ди дv дт ^ (4
дх ду дг '
< д д д д
где — = — + (и+и)——+ (У+V)——+ т—--оператор индивидуальной производной; Вху =
аЬ дЬ дх ду дг
д д д д
—Кх——-Ку—--оператор горизонтального турбулентного обмена с коэффициентами
дх дх ду ду
Кх, Ку; Л — параметр плавучести.
Турбулентное замыкание подсеточных масштабов проводится на основе простой модели Смагоринского для коэффициентов горизонтального обмена Кх, Ку и путем решения уравнений полуэмпирической (к — е)-теории турбулентности для расчета К. Основные соотношения модели приведены в [11]. Значения Ьх и Ьу определялись по результатам численных экспериментов с последовательным увеличением области интегрирования по горизонтали вплоть до выполнения требования инвариантности средних характеристик по отношению к множеству реализаций конвективного ансамбля [12].
2. Уравнения модели водоема
Разделение процессов в водоеме на средние и обусловленные конвекцией проведем таким же путем, как это было сделано выше для АПС. Используя представление (1) применительно к вертикальному столбу воды в срединной квазиоднородной части водоема, получим систему уравнений для среднего горизонтально однородного течения в ВПС, имеющего характерный вертикальный масштаб порядка 10 м [13]:
ди ~ д у ди д
ттт = IV + — К—--— ит ,
дЬ дг дг дг
дУ ~ д у дУ д
-тгт = —и + т-К---— гуту ,
дЬ дг дг дг
дТ = АКт дТ — д ту (5)
дЬ дг дг дг '
где Т — осредненная по сечению столба воды площадью Ьх х Ьу температура в водоеме; Т — конвективные отклонения температуры (здесь и далее символы с тильдой соответствуют
обозначениям АПС, но относятся к водной среде). Оператор осреднения / аналогичен по своей структуре оператору /, введенному в (2), с размерами области осреднения Ьх, Ьу. Система уравнений для описания проникающей конвекции в водоеме имеет вид
Ни ди 1 др у д у ди д
— + т— = —-7- + IV + Вху и + — К— + — 1Ь дг у дх дг дг дг
(V дУ 1 др п д ~ дЬ д
— + т— = -----1и + ПХуV + — К— + — ьт ,
№ дг р ду дг 'дг 'дг
(т 1 др , а п , п ~ I д п дт
кг = - ~р дХ+9вт Т+Дхут+
— + т—= б Т + —к— + д
а дг ху дг дг дг
ди дЬ дт о (6)
дх ду дг '
где р — возмущение давления; р — средняя плотность воды; 9 — ускорение силы тяжести; вт — коэффициент термического расширения воды.
3. Краевые и начальные условия. Метод решения
Определим вертикальную структуру области следующим образом: Н — верхняя граница АПС, к — толщина приводного слоя воздуха, для которого выполняется предположение о постоянстве турбулентных потоков, г = 0 — уровень поверхности раздела вода — воздух, г = Н — нижняя граница области в водоеме.
Для осредненных уравнений АПС (3) и ВПС (5) поставим следующие условия:
и = V = V
о,
О©
дг
7я при г = Н;
(7)
ди дУ д©
К— = с«|и|и, К— = си|и|У, -рсрК— = Qт при г = к; (8)
рК
.д_и
дг
рК
ди
дг
и = У
рК
0,
.дУ_ дг
рК
дУ
дТ
дг
дг
7я при г
Т = Т0 при
0;
Н,
(9) (10)
где 7н = 3 °С/км — коэффициент стратификации стандартной атмосферы; си — коэффициент сопротивления; Qт — поток тепла в атмосферу от водоема; Т0 — температура поверхности вод; Пн — устойчивая температурная стратификация водоема на глубине. В данной постановке мы ограничились априорным заданием величин Qт и Т0. Такой подход позволяет выяснить возможность развития конвективной неустойчивости на основе реальных значений Qт, Т0, взятых из наблюдений.
Сформулируем краевые условия для систем (4), (6). Формально эти условия должны быть согласованы с (7) - (10) и должны дополнять их. Можно, однако, считать, что конвективные структуры в АПС интенсивно развиваются выше приводного слоя, т. е. при г > к. Это полностью соответствует теоретическому и экспериментальному представлениям о характере процессов вблизи подстилающей поверхности [3]. Согласно этому запишем
дф + сдф дЬ дг
0 при г = Н;
и = V = т = 0, и = V = т = 0,
9 = 90(£,х,у) при г = к; Т = Т0(£,х,у), при г = 0;
(11)
(12) (13)
2
Ж - СI =0 при * = Н, ^ (14)
где в0, Т0 — случайные возмущения температуры малой амплитуды; С, С — фазовые скорости гравитационных волн. Условия (12), (13) позволяют раздельно провести численную реализацию двух сложных трехмерных систем (4) и (6), не задаваясь вопросами склейки решений при * = 0 (тем более, что области определения решений в горизонтальной плоскости, вообще говоря, не совпадают при различных размерах Ьх,у,Ьху). Реализация этих условий оказалась возможной благодаря исходному выделению средних течений и отклонений в (1).
В качестве начальных условий для (3) и (5) задавались следующие
Ф = Ф0, Ф = Ф0 при Ь = Ь0, (15)
где Ф0, Ф0 — стационарные решения систем (3), (5) при отсутствии конвекции и А^© = 0.
Переходя к изложению метода решения, заметим, что вихреразрешающее моделирование требует значительных вычислительных ресурсов: так, типичное значение временного шага в ЬЕБ-моделях, обсуждаемых в [6], составляет 1 с. Для принятого здесь сеточного разрешения 128 х 128 х 60 узлов задача описания суточной моды с указанным значением ДЬ является весьма трудоемкой. В данной численной модели шаг по времени взят равным 15 с и выбран из соображений аппроксимации физических процессов. Возможность увеличения шага достигнута за счет применения неявных конечно-разностных схем на основе метода расщепления [14]. Кроме того, принятое согласно (1) раздельное описание средних полей и конвективных отклонений допускает постановку радиационных краевых условий (11), (14) на крупновихревой компонент решения, которые обеспечивают свободный выход из области быстрых гравитационных волн. Перечисленные особенности делают данную модель доступной для расчетов на персональных компьютерах класса РеП;шт-3.
Геометрическая область решения представляет собой параллелепипед, в котором введена регулярная сетка. Дискретизация уравнений проводится в терминах исходных переменных "скорость — давление" с использованием основных [г,], к) и полуцелых 1 1.1
г + 2,] + 2 ,к + узлов. Сеточная структура полей позволяет аппроксимировать уравнение неразрывности со вторым порядком точности в пределах одной элементарной ячейки и сконструировать схему с сохранением вторых моментов. Конечно-разностные аналоги исходной системы получены из эквивалентной консервативной формы уравнений.
Для интегрирования по времени использован неявный метод расщепления в варианте "предиктор — корректор", модифицированный в целях повышения устойчивости алгоритма. Этот метод оказался предпочтительнее более простого метода покомпонентного расщепления по соображениям, связанным со спецификой развития проникающей конвекции.
Обозначим через А оператор, содержащий слагаемые адвективного переноса и турбулентного обмена, а в оператор В включим линейные части уравнений, традиционно учитываемые на этапе динамического согласования. Согласно методу "предиктор — корректор" задача на элементарном временном интервале решается в два этапа: вначале определяются вспомогательные поля ф1+1/4, ф1+1/2 последовательным обращением операторов А и В на основе чисто неявных схем, а итоговое решение получается уточнением вспомогательного решения по явной схеме с полным оператором [А + В).
Однако применение явных схем нежелательно при описании гидродинамически неустойчивых течений, которые на определенных этапах развития могут экспоненциально возрастать со временем. Это соображение подводит к идее объединения этапов коррекции
и согласования в один блок, сохраняющий неявный способ реализации. Таким образом, приходим к схеме, которая является обобщением метода "предиктор — корректор" за счет исключения явного этапа в расчетном цикле. Схема тестировалась на некоторых известных решениях типа рэлей — бенаровской конвекции.
4. Результаты расчетов
В качестве исходной информации для модели АПС задавался поток явного тепла, полученный в ходе натурного эксперимента на оз. Красное летом 1984 г. [3]. Согласно данным [3], в 20 ч местного времени Qт был близок нулю, далее его значение увеличивалось, достигая максимума 25 Вт/м2 в 02 ч, после чего начинало постепенно уменьшаться. В соответствии с синоптической ситуацией в период наблюдений (штиль или слабый ветер) зададим и а = У а = 0. На интервале с 20 ч до 04 ч утра происходит выхолаживание приводного слоя воздуха со средней скоростью 1.5 °С/ч. Это приводит к формированию неустойчивой стратификации в АПС, которая обусловливает развитие конвективных движений в виде всплывающих термиков и вертикальных струй.
На рис. 1, а представлен рассчитанный по полной модели вертикальный профиль осред-ненной потенциальной температуры воздуха (кривая 1). Распределение оказалось типичным для конвективного АПС, а именно: отчетливо выделяется слой перемешивания толщиной Hi ~ 500 м со слабоустойчивой стратификацией и противоградиентным потоком тепла; выше формируется характерный инверсионный перегиб профиля В; в нижней части области расположен тонкий (несколько десятков метров) слой со сверхадиабатическим градиентом Вг.
При моделировании конвекции в водной среде суточный ход поверхностной температуры Т0(1) задавался из наблюдений, отвечающих тому же временному периоду, что и для атмосферы. Температура воды изменялась от 18 °С в вечернее время до 15.8 °С в утренние часы, т. е. происходило ночное выхолаживание поверхности водоема. Охлаждение поверхности приводит к формированию плотностной неустойчивости верхнего слоя,
Рис. 1. Вертикальное распределение потенциальной температуры АПС (а) и температуры ВПС (кривые разнесены по горизонтальной оси, Т\г=о = 16.8°С) при £ = 03 ч (б): кривые 1 соответствуют полной задаче с Ла = Уа = 0, кривые 2 — с Ла = 5 м/с, пунктир 3 — градиентно-диффузионная модель с Л а = 0.
которая разрешается в виде ансамбля тонущих переохлажденных термиков, и происходит их замещение теплыми массами воды с нижних горизонтов. Устойчивая в начальный момент стратификация разрушается с образованием слоя перемешивания в профиле средней температуры. Качественно эволюция средней температуры ВПС близка к зеркально отраженной структуре поля ©(¿) в АПС (рис. 1, б, кривая 1).
Точечные кривые 2 на рис. 1 соответствуют варианту с Ца = 5 м/с (слабый ветер в АПС и вызванное им дрейфовое течение в ВПС). Увеличение толщины перемешанного слоя, особенно заметное в водоеме (кривая 2 на рис. 1, б), обусловлено вынужденной конвекцией динамического происхождения. Определенный интерес представляют результаты расчета с чисто диффузионной К-моделью, те. когда в системах (3), (5) конвективные потоки не учитываются (в'Ш =Тй = 0). Из рис. 1 (пунктирные линии 3) видно, что модель без учета свободноконвективного обмена неудовлетворительно воспроизводит характерные черты распределения температуры по вертикали: в обеих средах наклоны Вг, Т отрицательны, и, следовательно, полного разрушения неустойчивости не происходит. Данный вывод полностью согласуется с результатами [15].
Стадия максимального развития конвекции приходится на 02-04 ч местного времени, причем максимальные вертикальные скорости в АПС достигают 1.5 м/с, а экстремум й в ВПС составил 1.8 см/с. Представление о пространственной структуре поля крупновихревой турбулентности дает рис. 2, на котором изображены контуры плавучих струй в виде изолиний восходящих (а) и нисходящих (б) токов в АПС и ВПС соответственно. Видим, что в АПС конвективный ансамбль представляет собой систему нерегулярно расположенных разномасштабных когерентных структур. Их горизонтальные масштабы укладываются в диапазон реально наблюдаемых в атмосфере конвективных образований.
Рассматривая поле й в водной среде (рис. 2, б), отметим его своеобразную "сетчатую" макроструктуру. В отличие от атмосферы (поле а), где термики изолированы друг от друга, в воде образуются нитевидные конгломераты погружающихся конвективных элементов. Это связано со сравнительно высокой интенсивностью конвекции в ВПС и обусловлено коагуляцией близко расположенных структур. У всплывающих объемов воды подобных свойств самоорганизации не наблюдается. С увеличением глубины интенсивность вертикальных токов ослабевает, связность структур исчезает и конгломераты распадаются на ряд отдельных термиков. Средние горизонтальные масштабы конвективных элементов составляют 630 м в АПС и 7.2 м в ВПС.
Рис. 2. Горизонтальная структура поля вертикальной скорости в АПС (а) и ВПС (б).
а
б
ж, м
О
2
4
ж, км
Рис. 3. Сечение в плоскости (хг) полей отклонений температуры в АПС (а) и ВПС (б).
Рис. 3 иллюстрирует вертикальное сечение полей в, Т в тот же момент времени, что и рис. 2 (с целью улучшения разрешения рисунка по горизонтали показан лишь фрагмент области). Для доказательства подобия структуры течений в АПС и ВПС вертикальная координата на рис. 3 нормирована значениями толщины слоев перемешивания И|И.|, которые определялись по профилям конвективного потока тепла. Над крупными термиками в АПС (поле а) формируются холодные "шапки" (г/И. ~ 1), за счет которых реализуется механизм вовлечения. Приведенная структура конвективных полей сходна с полученной в [7], где моделируется конвекция в морском пограничном слое. Структура поля Т в воде (рис. 3, б) имеет аналогичный характер, за тем исключением, что направленные в глубину вертикальные струи холоднее окружающей воды, а отклонения температуры в слое термоклина положительны. Из сопоставления полей, приведенных на рис. 2, 3, видно, что характерные масштабы термических структур воды и воздуха отличаются на два порядка.
Отметим ряд общих для АПС и ВПС закономерностей, присущих турбулентно-конвективным явлениям в этих геофизических средах. Основным источником энергии конвективных движений является поток тепла, подводимый (отводимый) к поверхности раздела. Он обеспечивает накопление доступной потенциальной энергии в нижних слоях АПС и верхних слоях столба воды. Положительная плавучесть обусловливает вертикальное смещение частиц от исходного состояния, т. е. генерируется ад-компонент кинетической энергии. Ускоренное движение частицы по вертикали приводит в действие механизм ее замещения, в результате чего возникает течение с горизонтальной конвергенцией и пополняются и, V-составляющие кинетической энергии. Результаты численного моделирования энергетики систем (4) и (6) подтверждают отмеченные закономерности и показывают однотипность энергетических преобразований в конвективных пограничных слоях водоема и атмосферы. Вопросы энергетического подобия этих природных сред в терминах скоростей диссипации обсуждались в [1, 16].
5. Заключение
Сделан важный шаг на пути построения объединенной модели АПС и ВПС, способной описать крайне сложный механизм энергомассообмена между водоемом, прилегающей сушей и атмосферой. Представленные модели конвекции в атмосфере и водоеме позволяют
вполне адекватно описывать как поля средних течений, так и тонкую пространственно-временную структуру турбулентных потоков с когерентными структурами, возникающих при развитой проникающей конвекции над относительно теплым водоемом в период ночного охлаждения.
Обнаружено качественное подобие пространственных структур полей средней температуры и конвективных пульсаций в атмосфере и водоеме. Подобие обусловлено сходными механизмами вертикального переноса энергии и выражается в идентичности строения конвективных слоев перемешивания, наличия инверсии в АПС и термоклина в ВПС, тонких слоев с неустойчивой стратификацией вблизи поверхности раздела. Сопоставление конвективных полей показывает их близость по свойству генерации системы термиков с собственным внутренним масштабом, определяемым физическими характеристиками среды. Характерные временные масштабы АПС и ВПС примерно одинаковы, а пространственные различаются на два порядка.
Авторы выражают благодарность чл.-корр. РАН В. Н. Лыкосову и д. ф.-м. н. В. М. Маль-бахову за полезные дискуссии.
Список литературы
[1] Brubaker J. M. Similarity structure in convective boundary layer of a lake // Nature. 1987. Vol. 330, No. 6150. P. 742-745.
[2] Imberger J., Patterson J. C. Physical limnology // Advances in Appl. Mech. 1989. Vol. 27. P. 303-475.
[3] Панин Г. Н. Тепло- и массообмен между водоемом и атмосферой в естественных условиях. М.: Наука, 1985. 206 с.
[4] Бояринов П. М., Петров М. П. Процессы формирования термического режима глубоких пресноводных водоемов. Л.: Наука, 1991. 176 с.
[5] Shau T. J., Greeg M. C. Turbulence in a oceanic convective mixed layer //J. Phys. Oceanogr. 1986. Vol. 16. P. 1777-1798.
[6] Neustadt F. T.M., Mason P. J., Moeng C.-H., Schumann U. Large-Eddy Simulation of the Convective Boundary Layer: a Comparison of four Computer Codes // Turbulent Shear Flow. N.Y., Berlin: Springer-Verlag, 1991. P. 343-367.
[7] Walter B. A. The Evolution of organized mesoscale convection in the marine planetary boundary layer during cold air outbreaks // Third Conf. on Mesoscale Processes. 1987. Vancouver. B.C., Canada. P. 202-203.
[8] Stevens B., Lenschow D. H.. Observations, experiments, and Large Eddy simulation // Bull. Am. Met. Soc. 2001. Vol. 82, No. 2. P. 283-294.
[9] КононЕнко С.М., Малъбахов В.М., Пушистов П.Ю. Численная модель ансамбля нерегулярной мезомасштабной конвекции в планетарном пограничном слое Экмана // Тр. ЗСРНИГМИ. 1978. Вып. 41. C. 91-100.
[10] ГУТМАН Л. Н. Введение в нелинейную теорию мезометеорологических процессов. Л.: Гидрометеоиздат, 1969. 295 с.
[11] Shlychkov V.A., PuSHISTOV P. Yu. Application of Eddy-Resolving Models for Penetrating Turbulent Convection in the Atmosphere and Deep Lakes // Bull. Nov. Comp. Center, Num. Model. in Atmosph., etc. ICM MG Publ., Novosibirsk: 2000. Vol. 6. P. 79-88.
[12] ВАН Мигем Ж. Энергетика атмосферы. Л.: Гидрометеоиздат, 1977. 327 с.
[13] Пушистов П. Ю. Когерентные структуры в глубоких водоемах умеренных широт (наблюдения, оценка возможности математического моделирования) // Фундаментальные проблемы воды и водных ресурсов на рубеже третьего тысячелетия: Тр. Междунар. науч. конф. Томск, 2000. С. 431-435.
[14] Ковеня В. М., ТАРНАВСКИй Г. А., Черный С. Г. Применение метода расщепления в задачах газовой динамики. Новосибирск: Наука, Сиб. отд-ние, 1990. 248 с.
[15] Bolding K., Stips A., BURCHARD H., Eifler W. A comparison between observed and modelled dissipation rates under convective conditions // Physical Processes in Natural Waters: Proc. 5th Workshop. Irkutsk, 2000. P. 80-83.
[16] КИТАЙГОРОДСКИЙ С. А. Влияние опрокидывания ветровых волн на локальное взаимодействие между атмосферой и океаном // Изв. РАН. Физика атмосферы океана. 1997. Т. 33, №6. C. 828-836.
Поступила в редакцию 2 октября 2001 г., в переработанном виде — 5 декабря 2001 г.