«ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ВЗАИМОДЕЙСТВИЯ АТМОСФЕРЫ С ВОДОЕМАМИ СУШИ ...»

Московский Государственный Университет

имени М. В. Ломоносова

Научно-исследовательский вычислительный центр

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

УДК 551.5

Степаненко Виктор Михайлович

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ

ВЗАИМОДЕЙСТВИЯ АТМОСФЕРЫ С ВОДОЕМАМИ СУШИ

Специальность 25.00.30 – метеорология, климатология и агрометеорология Диссертация на соискание ученой степени кандидата физико-математических наук

Научный руководитель:

член-корреспондент РАН, доктор физико-математических наук В. Н. Лыкосов Москва – 2007 г.

Посвящается светлой памяти Михаила Арамаисовича Петросянца

ОГЛАВЛЕНИЕ

Введение………………………………………………………………………………… Глава 1. Обзор современного состояния проблемы………………………………… Глава 2. Термодинамическая модель водоема и ее верификация………………….. 2.1. Физическая постановка задачи и описание модели…………………………………….. 2.2. Численная реализация модели……………………………………………………………. 2.3. Верификация модели……………………………………………………………………… Глава 3. Параметризация водоемов в мезомасштабной атмосферной модели и ее верификация

3.1. Описание мезомасштабной модели……………………………………………………… 3.2. Верификация мезомасштабной модели и параметризации водоемов………………… Глава 4. Агрегирование турбулентных потоков над гидрологически неоднородной сушей………………………………………………………………………………….. 4.1. Метод эффективных параметров……………………………………………………….. 4.2. Мозаичный метод………………………………………………………………………... 4.3. Другие методы агрегирования турбулентных потоков………………………………... 4.4. Методика верификации мозаичного метода агрегирования………………………….. 4.5. Ошибки мозаичного метода в дневное время суток…………………………………... 4.6. Об одном эффекте нелинейности формул аэродинамического метода……………… 4.7. Ошибки мозаичного метода в ночное время суток……………………………………. 4.8. Выводы……………………..…………………………………………………….. Заключение…………………………………………………………………………… Список литературы……………………………………………………………………

ВВЕДЕНИЕ

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

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

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

Водные объекты в ряде регионов занимают значительную часть площади (рис. 1) и оказывают большое влияние на региональный климат.

Рис 1. Доля площади суши, занятая водными объектами, болотами и переувлажненными территориями (Cogley, 1991) Очевидно, что для этих регионов исключительно актуальна проблема реакции озер, водохранилищ и их экосистем на современные и будущие изменения климата (Эдельштейн, 2005). Она приобретает особенное значение, поскольку как раз в регионах с повышенной плотностью гидрологической сети в высоких широтах в XX в. наблюдалось хорошо выраженное потепление климата. Сказанное относится к северу территории США, Канаде, Скандинавским странам, Карелии, Западной Сибири (рис. 1). Например, в Карелии (Климат Карелии: изменчивость и влияние на водные объекты и водосборы, 2004) среднегодовая приземная температура за 100 лет наблюдений (1890 – 1990 гг.) выросла на 0.6 °С, что привело к заметным последствиям в гидрологическом режиме:

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

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

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

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

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

согласно последнему, под гидрологической неоднородностью понимается неоднородное распределение по территории характеристик гидрологического режима, в первую очередь – величин стока.

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

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

Для достижения этой цели решались следующие задачи:

• определение наиболее актуальных направлений исследований в области • построение и верификация с использованием данных наблюдений одномерной модели термодинамики водоема;

• включение одномерной модели в мезомасштабную атмосферную модель в качестве параметризации водоемов и верификация совместной модели с привлечением доступных данных наблюдений;

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

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

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

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

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

Автор благодарит коллектив кафедры метеорологии и климатологии географического факультета МГУ и, в особенности, А. В. Кислова, Н. Ф. Вельтищева, И.

В. Тросникова и Б. А. Семенченко за конструктивные замечания по существу работы.

В заключение автор выражает признательность своим родителям Е. О. Бароновой и М. М. Степаненко, которые своим примером способствовали выбору научной деятельности в качестве профессии.

Представленные в данной работе исследования выполнены при финансовой поддержке РФФИ (гранты 04-05-64898 и 07-05-00200), ИНТАС (гранты 00-189 и 01-2132) и Совета ученых географического факультета МГУ (грант молодых ученых 2006 г.).

ГЛАВА 1.

ОБЗОР СОВРЕМЕННОГО СОСТОЯНИЯ ПРОБЛЕМЫ

Во введении обоснована необходимость использования моделей внутренних водоемов суши в задачах прогноза погоды и климата. Ниже мы рассмотрим основные результаты, полученные научным сообществом к настоящему времени в области численного моделирования водоемов и использования моделей водоемов в задачах взаимодействия атмосферы с сушей. Следует отметить, что литература по этим направлениям накопилась довольно обширная, поэтому настоящий обзор не претендует на полноту, и ссылки на конкретные работы часто приводятся в качестве примера, в то время как существует множество работ других авторов близкой тематики и методологии исследования. Актуальность рассматриваемой в настоящей работе проблемы можно метеорологического общества начиная с 1990 г. вышло около 110 статей, в названии которых присутствует слово «озеро», причем подавляющее их количество посвящено различным аспектам взаимодействия озер с атмосферой.

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

Например, для Великих Американских озер разработана система оперативного прогноза волнения, температуры и других параметров состояния водоема, причем в качестве граничных условий используются данные расчетов по мезомасштабной модели ММ (Kelley et al., 1998).

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

Трехмерные вихреразрешающие модели водоемов находят применение также при совместном моделировании пограничного слоя атмосферы и перемешанного слоя водоема (океана) (Пушистов и Шлычков, 2001). Расчеты по совместной модели оказываются вычислительно весьма затратными и могут быть эффективно осуществлены только при параллельной реализации программного кода (Glazunov and Lykosov, 2003). Совместное вихреразрешающее моделирование является полезным инструментом изучения сложных обменных процессов, происходящих в обеих средах; в частности, его результаты могут быть использованы для построения параметризаций нелокального турбулентного обмена.

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

Одномерное приближение основано на эмпирическом факте, что интенсивность горизонтального теплообмена в водоемах, как правило, на несколько порядков меньше интенсивности вертикального теплообмена. Как показано в многочисленных работах, одномерные модели демонстрируют достаточно хорошее совпадение с натурными наблюдениями (Mironov et al., 2006; Goudsmit, 2002; Гранкина, 2006). Тем не менее, не для всех водоемов одномерное приближение можно считать оправданным: как показано в некоторых работах (Гречушникова, 2001), в водохранилищах тепловой эффект горизонтальной адвекции может быть сопоставим с эффектом вертикального теплообмена. Из этого следует, в модели существенно проточного водоема необходимо вводить по крайней мере одно горизонтальное измерение, параллельное ведущему течению.

используемых в задачах взаимодействия водоемов с атмосферой. Большинство из них основано на решении уравнения теплопроводности конечно-разностными методами.

Пожалуй, первые работы, в которых в блоках взаимодействия с подстилающей поверхностью атмосферных моделей стали использоваться одномерные модели водоемов, появились в США в начале 90-х годов XX в. (Bates et al., 1993; Hostetler et al., 1993)2. Они были посвящены одной из наиболее практически значимых проблем взаимодействия озер с атмосферой – моделированию зимних штормов над крупными незамерзающими озерами.

Одной из ключевых проблем построения одномерных моделей водоема является параметризация коэффициента турбулентности. Как известно, турбулентность в озерах развивается под влиянием динамических факторов (трение ветра, ветровое волнение, сейши и т.д.) и эффектов стратификации (в общем случае, температурной и соленостной) (Алексеевский, 2006; Wuest and Lorke, 2003), которые желательно адекватно учесть при вычислении этого коэффициента. К настоящему времени широкое применение получили диагностические формулы для коэффициента турбулентности (Rodi, 1993; Engelund, 1978;

Simoes, 1998), которые учитывают влияние стратификации и сдвига скорости течений через число Ричардсона. В частности, они используются в схемах тепло- и влагообмена на Возможно, первый опыт параметризации озер в атмосферных моделях описан в работе (Pitman, 1991), в которой наличие водоемов в ячейке атмосферной модели учитывалось через повышенную влажность почвы.

поверхности суши, предназначенных для использования в моделях общей циркуляции атмосферы (Technical description of community land model, 2004), поскольку позволяют значительно экономить время вычислений. В то же время, в моделях, работающих в исследовательском режиме, все чаще используется «k-» параметризация турбулентности (Монин и Яглом, 1965; Mellor and Yamada, 1974), которая является более содержательной с физической точки зрения, чем диагностические схемы, но значительно уступает последним по вычислительной эффективности. Эта параметризация основана на решении двух прогностических уравнений относительно кинетический энергии турбулентности E (ТКЕ) и скорости и диссипации, а коэффициент турбулентности рассчитывается по значениям этих двух величин по известной формуле Колмогорова. Параметризация «k-»

критикуется рядом исследователей, в частности, за то, что при выводе уравнения для использованы недостаточно обоснованные гипотезы. Кроме того, как показывают численные эксперименты (см. главу 2), эта параметризация весьма чувствительна к выбору входящих в нее полуэмпирических безразмерных коэффициентов. Тем не менее, на настоящем этапе «k-» параметризация, по всей видимости, еще далеко не исчерпала потенциала развития: в частности, предложены ее модификации, учитывающие производство ТКЕ за счет сейш (Goudsmit et al., 2002).

исследователей предлагает так называемую «k-» параметризацию (Umlauf et al., 2003), в которой уравнение для скорости диссипации ТКЕ заменяется уравнением для – частоты турбулентности. Так или иначе, вопрос о степени физической адекватности той или иной параметризации турбулентности, по всей видимости, в настоящее время далек от удовлетворительного разрешения. Значительного продвижения в этом вопросе можно ожидать от сравнительных численных экспериментов по моделированию различных типов турбулентных течений, с привлечением всех предлагаемых параметризаций. В частности, эта идея реализуется авторами известной одномерной модели пограничного слоя океана (водоема) GOTM (Umlauf et al., 2006).

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

В задачах оперативного прогноза погоды особенное значение приобретает вычислительная эффективность всех блоков прогностической модели. Поэтому в рамках моделей краткосрочного прогноза весьма перспективным представляется использование полуэмпирических моделей водоемов. Одна из таких моделей создана в немецкой службе погоды DWD (Mironov et al., 2006). Уравнения этой модели получены аналитическим интегрированием уравнения теплопроводности с привлечением эмпирического вида вертикального профиля температуры, и представляют собой обыкновенные дифференциальные уравнения относительно нескольких термодинамических переменных – температуры слоя перемешивания, температуры дна водоема и др. Оказалось, что построенная таким образом модель достаточно адекватно воспроизводит термический режим реальных озер, и в этом отношении почти не уступает более сложным дифференциальным моделям (Dutra et al., 2006); при этом время вычислений в этой модели более чем на порядок меньше, чем для дифференциальных моделей.

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

Например, двухмерные модели использовались для воспроизведения подледной конвекции в весеннее время года (Pushistov and Ievlev, 2000).

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

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

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

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

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

Вследствие разницы потока явного тепла между сушей и водной поверхностью в береговой зоне развиваются бризы (Бурман, 1965), имеющие, как правило, выраженный суточный ход направления ветра: днем ветер дует с водоема на сушу в направлении, близком к нормали к береговой линии, ночью – в обратную сторону. Бризы в значительной степени определяют климат прибрежных районов, в частности, уменьшая дневные максимумы температуры над сушей (Вельтищев, 1989). Скорости ветра в бризах тропических широт могут достигать 10 м/с и более. Кроме того, давно отмечено, что бризовая циркуляция существенно влияет на перенос атмосферных аэрозолей и газовых примесей (Бурман, 1965). Этот факт имеет большое значение для мегаполисов, расположенных на берегах крупных озер, например, для Чикаго (оз. Мичиган), Детройта и Кливленда (оз. Эри). Результаты мезомасштабного численного моделирования (Eastman et al., 1995) показали, что наблюдаемые превышения ПДК некоторых загрязняющих веществ в воздухе этих городов могут быть в значительной мере обусловлены переносом примеси к источнику возвратным течением бризовой ячейки.

Зимой в высоких широтах подавляющее количество озер замерзает и покрывается снегом, как и окружающая суша. Вследствие этого, свойства подстилающей поверхности (альбедо, теплоемкость, теплопроводность) водоема мало отличаются от тех же свойств суши, и свойственный ему в летнее время внутренний пограничный слой, по всей видимости, зимой практически не выражен. Поэтому зимой бризы практически не наблюдаются. Важное исключение составляют крупные незамерзающие озера, над которыми нередко возникает сильно неустойчивая стратификация при вторжении арктических воздушных масс. При этом развивается интенсивная конвекция, сопровождаемая выпадением большого количества снега. Более того, если имеются благоприятные синоптические условия, то над такими озерами может образоваться мезоциклон, в котором гряды конвективной облачности вытягиваются к центру низкого давления, приобретая спиралевидную форму, легко распознаваемую на спутниковых снимках (Forbes and Merrit, 1984). Перечисленные выше явления подробно изучены для густонаселенного района Великих Американских озер, где они производят значимый экономический эффект. За последнее десятилетие в мировой литературе вышло, по меньшей мере, несколько десятков статей, посвященных разным аспектам зимних конвективных явлений над этими озерами. Первые работы, в которых расчет этих конвективных явлений производился посредством мезомасштабной атмосферной модели, а термодинамика озер – на основе одномерной модели, появились в начале 90-х гг. ХХ в.

(Bates et al., 1993). В них было показано, что наличие Великих Американских озер в центре континента увеличивает количество зимних осадков и повышает приземные зимние температуры в этом регионе. Аналогичные оценки влияния водоемов на региональный климат проводились впоследствии и для других регионов, например, для района озер штата Невада (Hostetler et al., 1993) и для центральной Канады (Mackay, 2006).

Оценке влияния водоемов суши на глобальный климат посвящена работа (Bonan, 1995). В ней при воспроизведении современного климата моделью общей циркуляции атмосферы в качестве параметризации водоемов использована одномерная термодинамическая модель. Наиболее значимый отклик современного климата на присутствие озер оказался приурочен к крупнейшим озерным системам мира: Великим Американским и Канадским озерам, озерам Африки, карело-финской системе и озерам Западной Сибири. Средняя приземная температура летом в этих районах оказалась на 1- °С ниже, поток скрытого тепла на 10-45 Вт/м2 больше, а поток явного тепла на 10-30 Вт/м меньше, чем в контрольном эксперименте (без озер). Зимний отклик метеорологического режима для большинства озерных регионов оказался незначимым, за исключением Великих Американских озер, где после включения озерной параметризации выросло количество осадков на 0.4 мм/сут, и Карелии, где средняя зимняя приземная температура выросла на 10 °С. Очевидно, что вывод о столь высокой чувствительности зимних температур в Карелии к наличию озер требует дополнительной проверки в аналогичных экспериментах с региональными климатическими моделями.

Особый интерес представляет оценка отклика региональных гидрологических систем на изменения климата в будущем; этот вопрос неоднократно обсуждался в литературе. Чтобы произвести эту оценку, привлекаются, как правило, данные расчетов моделей общей циркуляции атмосферы, осуществленных по одному или нескольким сценариям МГЭИК (Межправительственная группа экспертов по изменению климата), которые выступают в роли входных данных для некоторой гидрологической модели, или интерпретируются непосредственно. Например, в работе (Кондратьев и Бовыкин, 2003) данные сценарных расчетов глобального климата XXI в. задавались на вход в модели системы «водоем-водосбор», а в работе (Арпе и др., 2000) годовой сток с водосбора р.

Нева рассчитывался как разность модельных осадков и испарения. Метододогически похожая методика оценки региональных гидрологических последствий глобального потепления с использованием простой диагностической модели испарения с суши представлена в работе (Климат Карелии: изменчивость и влияние на водные объекты и водосборы, 2004). Во всех перечисленных работах отмечается, что в XXI в. на северозападе Европейской территории России следует ожидать роста годового количества осадков; при этом отклик речного стока и уровней озер будут определяться как годовыми осадками, так и их внутригодовой изменчивостью, в частности, количеством твердых осадков, а также ростом испарения с водосборов и реакцией биоты. В частности, ожидается снижение уровня некоторых озер Карелии в теплое время года.

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

Итак, гидрологическая интерпретация климатических расчетов, базирующаяся на диагностических гидрологических моделях различного уровня сложности, позволила получить богатый багаж количественных оценок и качественных выводов относительно будущих изменений гидрологической системы суши, а также вариации ее параметров в прошлом. Однако методология этой интерпретации имеет недостаток, который может оказаться в некоторых случаях существенным: в климатических расчетах, как правило, не учитывается в достаточно полной мере обратное влияние изменений в гидрологической системе на климат. Таким образом, гидрологическая система оказывается в этих расчетах изолированной от климатической системы. В то же время, как показано в работе (Bonan, 1995) (см. выше), водоемы оказывают существенное влияние на региональный климат, и «приближение изоляции» может привести к значительным ошибкам в оценке термического и влажностного режима атмосферы. Эта «изоляция» может быть преодолена развитием гидрологических блоков в моделях общей циркуляции атмосферы. Можно надеяться, что это позволит существенно уточнить современные представления о будущих изменениях в озерно-речных системах планеты.

Отметим, что гидрологическая интерпретация расчетов по атмосферным моделям используется не только в численных экспериментах на климатическом масштабе времени, также в расчетах на меньших масштабах, в частности, в оперативных прогнозах паводков (Rivin and Heise, 2006).

Важным обстоятельством при включении параметризации озер в атмосферную модель становится то, что подавляющее большинство озер суши при современном пространственном разрешении этих моделей оказываются подсеточными объектами. Для учета подсеточной неоднородности суши в общем случае предложены различные методы (Koster and Suarez, 1992), получившие название методов агрегирования. Наиболее физически обоснованным методом можно считать так называемый мозаичный метод (Avissar and Pielke, 1989). Имеется большое количество работ, в которых оценивается точность методов агрегирования, но в них преимущественно рассматривается неоднородность растительности и влажности почвы. В то же время, насколько нам известно, для случая плотной гидрологической сети аналогичные оценки практически не проводились. Оценка применимости мозаичного метода для случая гидрологически неоднородной поверхности суши обсуждается в главе 4, там же приводится более детальный обзор литературы по этой тематике.

Анализ литературных источников и Интернет-ресурсов показал, что в настоящее время в мире существует несколько одномерных моделей водоема, предназначенных для использования в атмосферных моделях. Это блок термодинамики водоема в модели Community Land Model, версия 3.0 (Technical description of community land model, 2004), параметризации водоемов в климатической модели Института метеорологии общества Макса Планка (последняя версия 5.0) (Tsuang et al., 2001), региональной модели Hirlam (Ljungemyr et al., 1996), Канадской региональной модели (McFarlane and Flato, 2000) и модель Flake (Mironov, 2006). К настоящему времени накоплено уже достаточное количество расчетов, осуществленных с этими моделями как в автономном режиме (когда атмосферное воздействие считается известным), так и в режиме совмещения с атмосферной моделью; так что, по-видимому, уже назрела необходимость широкого сравнения этих моделей между собой, аналогично тому, как это делалось в рамках проектов PILPS (для моделей почвы) и SNOWMIP (для моделей снежного покрова).

В заключение отметим, что перспективным приложением одномерных моделей водоемов является их включение в системы усвоения данных на суше. Системы усвоения данных уже достаточно давно и успешно применяются в системах оперативного прогноза погоды (Фролов и др., 2000). Несколько позднее концепция усвоения данных вошла в практику моделирования циркуляции океана (Саркисян, 2000). В последние годы стали развиваться первые системы усвоения данных наблюдений на суше, например, LDAS (Land data assimilation system, NASA/NOAA/NCEP) (Rodell et al., 2004) и ELDAS (European land data assimilation system), также появились работы, посвященные различным аспектам создания подобных систем (Ettema and Viterbo, 2004). Очевидно, что реалистичное восстановление пространственного распределения термодинамических характеристик суши, являющееся основной задачей подобных систем, невозможно без применения физически адекватной параметризации водных объектов.

ГЛАВА 2.

ТЕРМОДИНАМИЧЕСКАЯ МОДЕЛЬ ВОДОЕМА

И ЕЕ ВЕРИФИКАЦИЯ

2.1. Физическая постановка задачи и описание модели Как показано в главе 1, построение моделей водных объектов суши востребовано в задачах взаимодействия этих объектов с атмосферой. При этом следует отметить, что в Сибири и Северной части Европейской территории России, и во многих других регионах с плотной гидрологической сетью водоемы занимают значительно большую часть суши, чем водотоки, поэтому в первом приближении наличием водотоков можно пренебречь.

Существенной особенностью мелких водоёмов (Malm et al., 1997; Pavlov, 1995) является то, что горизонтальная неоднородность распределения термодинамических параметров в них в большинстве случаев очень мала. Это позволяет в первом приближении пренебречь горизонтальным теплообменом и рассматривать перенос тепла только в вертикальном направлении. Кроме того, как показывают многочисленные теоретические и экспериментальные исследования, одномерное приближение также применимо к расположенному под водоёмом грунту, и к формирующимся в зимний период снежному и ледяному покровам.

2.1.1. Теплоперенос в теле водоёма. В основу описания термодинамического режима водоёма в настоящей модели положено одномерное уравнение теплопроводности.

Предполагается, что начало направленной вниз вертикальной координаты z совмещено с уровнем свободной поверхности водоёма, так что область, в которой ищется решение, представляет собой отрезок [0, h], где h=h(t) – характерная (средняя) глубина водоёма, а t – время. Уравнение теплопроводности с учетом поглощения солнечной радиации в теле водоема имеет вид:

Здесь с – теплоемкость воды, – её плотность, – коэффициент турбулентной теплопроводности, ограниченный снизу значением молекулярного коэффициента теплопроводности, Т – температура, S – поток солнечной радиации в толще водоёма, Bw = r E – скорость приращения слоя воды на верхней границе (водный баланс на поверхности водоема), r – интенсивность осадков, E – скорость испарения с поверхности водоема. Второе слагаемое в правой части (2.1) отражает эффект сдвига начала координаты z за счет осадков и испарения на свободной поверхности водоема. В z– системе координат область решения меняется во времени, что приводит к необходимости соответствующим образом адаптировать конечно-разностную сетку в процессе численного интегрирования. В этой связи удобно перейти от исходной вертикальной координаты z к новой вертикальной координате =, которая меняется в пределах [0, 1] и позволяет ввести фиксированную сетку. В переменных (, t ) уравнение (2.1) принимает вид:

Для расчета потока солнечной радиации в толще водоёма применяется широко используемая в различных исследованиях экспоненциальная зависимость где S (0) - суммарная солнечная радиация на поверхности водоема, – альбедо водной поверхности, е – коэффициент экстинкции, принимающий значения от 0.5 до 1 м-1, а e – доля солнечной радиации, поглощенной тонкой пленкой на поверхности водоема (Deas et al., 2000) (в настоящей версии принята равной 0.4).

Турбулентный теплообмен является основным механизмом вертикального обмена тепла и массой в водоёме. Его параметризация в предлагаемой модели осуществлена в нескольких вариантах, приведенных ниже.

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

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

причём температура T в данной формуле выражена в градусах Цельсия.

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

где s - коэффициент теплопроводности на поверхности водоема, b - некоторый коэффициент молекулярной теплопроводности, V – скорость ветра, V0 – скорость ветра, при которой коэффициент теплопроводности достигает максимального значения max.

Очевидно, что описанная параметризация не лишена недостатков: отсутствует, в частности, теоретическое обоснование задания временного интервала, через который производится мгновенное перемешивание, а также выбора параметров формулы для. Другими словами, эти параметры могут сильно различаться в разных расчёта географических регионах, а также между отдельными водоемами. Главное же достоинство рассмотренной параметризации – простота её алгоритмической реализации. В литературе (см., например, (Simoes, 1998)) предлагаются более сложные полуэмпирические параметризации турбулентной теплопроводности, некоторые из них также используются в предлагаемой модели и описаны в п. 2.1.5.

2.1.3. Параметризация турбулентности на основе уравнения для кинетической энергии турбулентности с использованием пути смешения. В рамках данного подхода используется формула:

кинетическая энергия турбулентности (черта сверху – знак осреднения по Рейнольдсу, величины со штрихами – отклонения от среднего значения), Cet = 0.072 - безразмерная константа. Кинетическая энергия турбулентности рассчитывается с помощью следующего уравнения (Монин и Яглом, 1965):

где k = l E - коэффициент турбулентности, E = Ce E ( E = 1 - безразмерная константа), Ce = 0.09, слагаемое представляет собой суммарное производство кинетической энергии турбулентности за счет сдвига скорости и за счет эффекта плотностной стратификации:

а - скорость диссипации кинетической энергии турбулентности, u и v – широтная и меридиональная составляющие скорости течения, которые вычисляются согласно формулам, описанным в п. 2.1.6.

Скорость диссипации энергии турбулентности рассчитывается по формуле Колмогорова:

Данная параметризация, получившая название «E-l» параметризации, обладает, безусловно, более богатым физическим содержанием, чем «эмпирическая». В то же время, использование в ней пути смешения l, относительно которого нет надёжных экспериментальных данных, вносит некоторую неопределенность в описание турбулентности.

2.1.4. Параметризация турбулентности на основе уравнений для кинетической энергии турбулентности и скорости ее диссипации. В рамках данной, получившей широкое признание и распространение (Лыкосов, 1992), параметризации (для краткого ее обозначения используются два названия - «k-» или «E-» параметризация) коэффициент турбулентной теплопроводности определяется следующим образом:

а коэффициент турбулентности k вычисляется по формуле:

Величина E рассчитывается с помощью уравнения (2.6), а скорость ее диссипации находится из уравнения:

где = Сe ( = 1.3 - безразмерная константа), а C1 - функция числа Рейнольдса Re :

где = 1.007 *106 м2/с – молекулярная вязкость воды, безразмерная константа C0 принята равной 1.9.

Заметим, что во многих моделях океанических и атмосферных пограничных слоев последнее слагаемое в уравнении (2.11) записывается в виде где C1 и C2 - различные константы. Это приводит к нарушению естественного условия баланса производства и диссипации ТКЕ P =, которое должно иметь место при стационарной и однородной по пространству турбулентности. В самом деле, применяя условия однородности и стационарности к уравнению (2.11) со слагаемым (2.13), получаем:

Очевидно, что используемое в настоящей модели уравнение (2.11) свободно от этого недостатка.

Определению безразмерных коэффициентов Сe, Cet, C0, E,, и др., которые входят в «E-» параметризацию, посвящено большое количество работ. Их значения определялись как в аналитических исследованиях решений (2.6, 2.11), так и при сравнении результатов моделирования с данными наблюдений. В наших численных экспериментах показано, что численное решение (2.6, 2.11) чувствительно к значениям безразмерных коэффициентов, что можно считать недостатком данной параметризации. В частности, использование «атмосферных» значений этих коэффициентов согласно работе (K. G. Rao et al., 1996), которые незначительно отличаются от «водных», приводит к тому, что модель занижает реальную интенсивность конвективного перемешивания в случае неустойчивой плотностной стратификации. Вследствие этого, осенняя конвекция в водоемах умеренных и высоких широт не воспроизводится, и в течение этого периода в водоеме сохраняется неустойчивая стратификация.

Граничное условие для энергии турбулентности при = 0 (на верхней границе) имеет следующий вид:

где k we - безразмерный коэффициент, учитывающий генерацию турбулентности за счет обрушения поверхностных волн (в ряде работ, и в настоящей модели принимается равным 100), - модуль напряжения трения (формулы расчета напряжения трения приведены в пп. 2.1.6 и 2.1.12). На нижней границе используется условие (2.14), но в левой части со знаком “+” и при kwe = 1. Граничные условия для скорости диссипации на верхней и нижней границах также задаются в форме потоков:

где = 0.38 - постоянная Кармана (на нижней границе в левой части (2.15) используется знак «+»).

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

параметризация.

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

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

1) Параметризация Никурадзе (Rodi, 1993) 2) Параболическая параметризация (Engelund, 1978) 3) Параметризация ре-нормализационной группы (Simoes, 1998) В этих формулах использованы следующие обозначения: lm - путь смешения, z = h z вертикальная координата, направленная вверх, = 0.38 - постоянная Кармана, C = 1.5 и C1 = 100 - безразмерные постоянные, u* - скорость трения в поверхностном слое водоема, ( x) = max(0, x), Ri - число Ричардсона, определяемое по формуле:

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

2.1.6. Расчет динамики водоема. Для вычисления зональной u и меридиональной v компонент скорости течения в описываемой модели используются уравнения типа диффузионными слагаемыми):

При этом формулируются следующие граничные условия: на свободной поверхности = 0 (граница водоем-атмосфера) поток импульса считается непрерывным:

где в правых частях записаны, соответственно, зональная и меридиональная параметризации турбулентных потоков в приземном слое (см. п. 2.1.12). Аналогичные граничные условия используются и на границах раздела Г водной среды с твердой поверхностью (вода-лёд, вода-грунт), только потоки импульса задаются в этом случае по формулам Шези (Чеботарев, 1975):

где C z - коэффициент Шези, определяемый шероховатостью поверхности. В этих формулах знак “+” используется в случае поверхности вода-грунт, а “–“ в случае поверхности вода-лед.

2.1.7. Водный баланс водоема. Для определения глубины водоема h записывается уравнение водного баланса водоёма:

где r – интенсивность осадков, Es – скорость испарения с поверхности, Rs – поверхностный сток, Rb – водообмен тела водоёма с нижележащим грунтом.

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

Отметим, что формула (2.23) строго говоря, справедлива только в том случае, если площадь зеркала водоема слабо зависит от его уровня.

2.1.8. Теплообмен в слое льда. Если в процессе вычислений оказывается, что температура поверхности водоёма опускается ниже 0 °С, то образуется первичный слой льда, толщина которого в модели принята равной 1 см. В этом слое также решается уравнение (2.1а), но при этом полагается, что суммарная солнечная радиация (за вычетом отраженной) полностью поглощается на поверхности льда и вниз не проникает. Это означает, что в уравнении (2.1а) не учитывается четвертый («радиационный») член в правой части. Кроме того, в качестве коэффициента теплопроводности используется значение молекулярного коэффициента для льда. На границе раздела лёд-вода в качестве граничного условия задаётся температура фазового перехода 0 °С. Таяние льда на границе с атмосферой происходит при температуре фазового перехода и его скорость определяется тепловым балансом на этой границе. При этом считается, что стаявшая вода мгновенно добавляется к водному слою.

Отсутствие в модели объемного поглощения солнечной радиации в ледовом покрове с проникновением части радиации под него исключает, в частности, возможность развития подледной конвекции в весенний период. Как показывают данные наблюдений (Jonas et al., 2003) и моделирования (Mironov et al., 2002; Pushistov and Ievlev, 2000) подледная конвекция является важным механизмом вертикального теплообмена и выравнивания профилей гидрофизических характеристик в водоеме. Поскольку моделирование подледной конвекции не входит в задачи настоящей работы, соответствующая модификация уравнения теплопроводности в слое льда не используется.

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

2.1.9. Тепловлагоперенос в снежном покрове. В зимнее время на замерзший водоём могут выпадать твердые осадки, образуя снежный покров, который в модели характеризуется распределением по вертикали двух основных параметров – температуры и удельного содержания жидкой влаги. Их эволюция рассматривается в координатах (z, t) и описывается следующей системой уравнений (Володина и др., 2000):

Здесь L – удельная теплота плавления, Ffr – скорость замерзания, W – удельное содержание жидкой влаги, – инфильтрационный поток жидкой влаги в снежном покрове. Поток в разностном выражении выглядит так:

где hg – гидравлическая проводимость, (в модели hg = 0.01 м/с), П – пористость снега, Whc – константа, характеризующая водоудерживающую способность снега (в модели Whc=0.04), z – шаг сетки модели по вертикальной координате. Кроме процессов, описываемых системой (2.24), в модели также используется параметризация процесса гравитационного оседания (уплотнения) снежного покрова во времени. На нижней границе снежного покрова (поверхность раздела «снег-лёд») температура и поток тепла предполагаются непрерывными. Для определения температуры поверхности снежного покрова используется уравнение теплового баланса, описанное в пункте 2.1.11.

2.1.10. Тепловлагоперенос в слое грунта под водоемом. В основу описания процессов тепловлагопереноса в грунте под водоёмом положена модель, представленная в работе (Володин и Лыкосов, 1998), в которой состояние почвы характеризуется температурой, содержанием жидкой, твердой и газообразной влаги. Поскольку под телом водоёма грунт должен быть насыщенным жидкой (при промерзании – твердой) влагой, то содержанием водяного пара в нем можно пренебречь. Если, кроме того, пренебречь влагопроводностью грунта за счет градиента температуры3, то соответствующая система уравнений принимает вид:

Здесь W – коэффициент влагопроводности, а I - удельное содержание льда. Как видно из этой системы, в почве рассматриваются процессы диффузии тепла и влаги, инфильтрация жидкой влаги, а также процессы промерзания/таяния воды. Коэффициенты, определяющие интенсивность этих процессов, зависят от переменных состояния грунта – T, W, I.

Это часто принимаемое допущение, основанное на малости коэффициента термовлагопроводности.

Гравитационный поток влаги в почве определяется по формуле, приведённой в работе (Clapp and Hornberger, 1978):

где soil,max - максимальное значение гравитационного потока, Wmax - максимально возможное содержание жидкой влаги (оно достигается при заполнении ею всего объема пор, свободных ото льда), bsoil – безразмерный показатель (нижним индексом max обозначено максимальное значение соответствующей величины, а индексом soil – величина, зависящая от типа почвы). Суммарная теплоёмкость почвы складывается из теплоёмкостей сухой почвы, воды и льда:

а коэффициенты диффузии влаги и тепла вычисляются в соответствии с формулами:

где величина Pf определяется через потенциал влажности (измеряется в м):

На границе с водоёмом задается непрерывность температуры и теплового потока, а также поток жидкой влаги, определяемый степенью насыщенности водой верхних горизонтов грунта. На нижней границе слоя грунта (в описываемых ниже экспериментах ее глубина задавалась от 20 до 100 м под дном водоёма) потоки тепла и влаги задаются равными нулю.

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

где S – суммарный поток солнечной радиации, Ea – поток встречного длинноволнового излучения атмосферы, Es – собственное излучение поверхности, Hs и LEs – потоки явного и скрытого тепла, соответственно, – альбедо поверхности. В переходные сезоны, когда на незамерзший водоем могут выпадать осадки в твердом виде, или случается дождь в энергетический вклад этих процессов. Альбедо водной поверхности зависит от высоты Солнца и вычисляется по формуле:

где h0 - высота Солнца.

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

1. Аэродинамический метод с коэффициентами, рассчитываемыми согласно теории подобия Монина-Обухова (Монин и Обухов, 1954). Соответствующие формулы имеют вид:

где s - напряжение трения, cp – теплоёмкость воздуха при постоянном давлении, а – плотность воздуха, CH, CE и CM– коэффициенты обмена для температуры влажности воздуха и импульса, соответственно, и q – потенциальная температура и удельная влажность на некотором уровне в приземном слое, s и qs – те же величины на поверхности земли, V – скорость ветра на некотором уровне в приземном слое.

Коэффициенты обмена CH, CE и CM задаются по полуэмпирическим зависимостям, предложенным в работах (Paulson, 1970; Businger, 1971; Beljaars et al., 1991).

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

2. Вторая параметризация (Louis, 1979) также построена на основе теории подобия Монина-Обухова. Однако в отличие от предыдущей, коэффициенты обмена определяются через число Ричардсона таким образом, что формулы (2.33) приобретают явный вид.

Вследствие этого параметризация становится более экономичной в вычислительном отношении. Подпрограммы данной параметризации заимствованы из почвенной модели ISBA (Mahfouf et al., 1995). В главе 2 данная параметризация будет описана более подробно.

2.2. Численная реализация модели Для решения уравнений модели используется разностная схема, основанная на методе расщепления по физическим процессам с неявным представлением диффузионных членов и с явным представлением слагаемых, ответственных за другие процессы.

Аппроксимация пространственных производных выполнена на сдвинутой сетке, изображенной на рис. 2.1.

Рис. 2.1. Конечно-разностная сетка модели водоема по координате Значения турбулентной кинетической энергии и скорости диссипации турбулентной рассчитываются в дробных узлах сетки, что позволяет аппроксимировать диффузионные слагаемые в уравнениях переноса тепла и импульса со вторым порядком точности. Кроме того, благодаря использованию сдвинутой сетки, слагаемые производства ТКЕ за счет сдвига скорости и эффекта стратификации в (2.7) также аппроксимируются конечноразностными аналогами второго порядка. Температура и компоненты скоростей течений рассчитываются в целых узлах сетки.

Шаг интегрирования по времени в разных экспериментах составлял от 5 до водоема.

2.2.1. Аппроксимация уравнения теплопроводности. Конечно-разностный аналог уравнения (2.1а) в отсутствие ледостава имеет следующий вид:

Здесь использована традиционная индексация узлов разностной сетки: верхний индекс обозначает номер шага по времени, нижний – номер узла по вертикальной координате.

Система линейных уравнений (2.34) имеет трехдиагональную матрицу и решается методом прогонки.

2.2.2. Аппроксимация и метод решения уравнения теплового баланса. Для вычисления температуры водоема на границе с атмосферой используется уравнение Здесь также приняты традиционные обозначения: индекс s означает величину на поверхности, ( s + 1) - величину, относящуюся к соседнему к поверхности узлу сетки, а – на уровне в атмосфере, на котором задаются метеорологические величины, звездочками * указаны величины, которые известны из данных натурных наблюдений или вычисляются в атмосферной модели.

В уравнении (2.35) неизвестной величиной является только Ts j +1, поскольку при заданном Ts j +1 правая часть (2.35) находится после решения уравнения теплопроводности (2.34) с температурой на границе Ts j +1. Трактуя (2.35) как уравнение относительно Ts j +1, его можно решить известными приближенными методами: методом касательных, методом хорд, и др. (Бахвалов, 2006). В рассматриваемой модели уравнение (2.35) решается в два этапа. На первом этапе происходит сканирование возможного диапазона температур, которое заканчивается нахождением отрезка температур, на концах которого невязка уравнения (2.35) имеет разный знак. На втором этапе внутри этого отрезка ищется корень Ts j +1 методом хорд. Как показали численные эксперименты, описанный итерационный процесс всегда сходится с высокой точностью (невязка уравнения (2.35) не превышает 0. Вт/м2) и достаточно быстро (не более 10-15 итераций).

2.2.3. Аппроксимация уравнений для турбулентной кинетической энергии и скорости её диссипации. Для аппроксимации уравнений турбулентной кинетической энергии и скорости диссипации используется схема Кранка – Николсон (Марчук, 1978), обеспечивающая второй порядок аппроксимации по времени. Согласно этой схеме уравнения (2.6, 2.11) заменяются следующими конечно-разностными аналогами:

Здесь использованы следующие обозначения разностных операторов:

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

После вычисления Ei ++1 2 и ij+11/ 2 по формулам (2.9) и (2.10) вычисляется коэффициент турбулентной теплопроводности ij+1 / 2.

2.2.4. Аппроксимация уравнений движения. При аппроксимации уравнений движения (2.20) также используется схема Кранка – Николсон:

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

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

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

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

интегрирования 30 лет);

2) Озеро Сырдах (Якутия, 70 км от Якутска, резко-континентальный климат, длительность интегрирования 20 лет);

3) Озеро Монте-Ново (Португалия, субтропический климат, длительность интегрирования 13 лет);

4) Озеро вблизи Тикси (устье р. Лены, арктический климат, длительность 5) Озеро Коссенблаттер (Германия, умеренный климат западных побережий, длительность интегрирования 14 дней);

6) Озеро Вендюрское (Карелия, умеренно-континентальный климат, длительность интегрирования 5 дней).

В экспериментах 1, 2, 3 воспроизводилась сезонная и межгодовая изменчивость водного и термодинамического режима озер, а в экспериментах 4, 5, 6 – суточная и синоптическая изменчивость.

2.3.1. Источники данных наблюдений в приземном слое атмосферы. Как следует из описания модели (п. 2.1), для ее интегрирования в качестве граничных условий необходимо задавать значения метеорологических величин в приземном слое атмосферы.

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

• данные стандартных метеорологических измерений на синоптических станциях;

(http://cdiac.esd.ornl.gov/ftp/ndp048) и содержащего ряды многолетних наблюдений на 225 станциях бывшего СССР. В настоящей работе данные этого архива использовались при моделировании водоема вблизи Колпашево (Томская область, правобережье р. Обь) и оз. Сырдах;

• данные специальных полевых экспериментов; эти данные использованы в экспериментах с озерами вблизи Тикси, оз. Вендюрское и оз. Коссенблаттер;

• данные реанализа; привлечены при моделировании оз. Монте-Ново.

Более подробная информация о подготовке метеорологических данных приводится ниже.

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

Для расчёта потока суммарной солнечной радиации на горизонтальную площадку S использовалась формула Кондратьева (Матвеев, 2000):

где 0 – эмпирическая функция высоты Солнца, – оптическая толщина атмосферы для интегрального потока, принимаемая равной 0.105, n – балл облачности в долях единицы, csh = 0.5607 – эмпирический коэффициент, h0 – высота солнца. Приход солнечной радиации на горизонтальную площадку на верхней границе атмосферы определяется известной формулой:

где S 0 - солнечная постоянная, – широта, – склонение Солнца, – часовой угол.

Встречное длинноволновое излучение атмосферы Ea рассчитывалось как функция температуры и влажности на высоте 2 м и балла облачности (Idso, 1981):

В этих выражениях T2 - температура воздуха на уровне 2 м, e2 - парциальное давление водяного пара на уровне 2 м, - постоянная Стефана - Больцмана. Собственное излучение поверхности задается известной формулой Стефана – Больцмана.

2.3.3. Воспроизведение температуры снега в Колпашево (Томская область). В Томской области находится множество озер небольшой глубины. В нашем распоряжении не имелось данных измерений на каком-нибудь конкретном озере, с которыми можно было бы сравнивать результаты моделирования. Поэтому объектом моделирования было выбрано некоторое абстрактное озеро с характерными для множества озер Томской области характеристиками. Глубина этого модельного озера составила 2 м. В качестве входных данных для модели использованы ряды измерений на метеорологической станции Колпашево за период с 1936 по 1966 гг. Начальные условия (вертикальное распределение температуры и других переменных состояния в начальный момент времени) в данном эксперименты не играли важной роли, поскольку достаточно быстро (по сравнению с полным периодом интегрирования) «забывались» моделью. Верификация модели произведена по параметру, измеряемому на метеорологической станции, – температуре поверхности снежного покрова. Сравнение результатов моделирования и данных наблюдений в Колпашево основывалось на предположении, что температура поверхности снега над замерзшим водоёмом, рассчитываемая моделью, близка к температуре поверхности снега на суше, которая измерялась в Колпашево (точность измерений – 0.5 °С). Это предположение представляется естественным, поскольку в условиях резко континентального климата ледовый покров озер достигает значительной толщины (1-1.5 м), и, следовательно, влияние водоема на тепловой режим вышележащего снежного покрова невелико. Результаты этого сравнения для двух первых месяцев 1961 г.

приведены на рис. 2.2. Как видно, рассчитанные кривые и кривые наблюдений хорошо коррелируют. В то же время, имеет место систематическое занижение моделью среднемесячного значения температуры на ~2 °С. Это, связано, по всей видимости, с несовершенством параметризации потоков тепла на подстилающей поверхности при устойчивой стратификации, которая характерна в зимних условиях, особенно, в районах континентального климата.

Рис.2.2. Температура поверхности снега в Колпашево по данным моделирования и измерений (январь – февраль 1961 г.) экспериментов. Озеро Сырдах – сравнительно большое, находящееся в цепи таких же, по размерам и глубине, озёр, соединяющихся только в периоды редких для этих территорий многоводий ручьями типа «травяных речек». Оно вытянуто в направлении с СЗ на ЮВ и занимает большую часть площади озерной котловины. Размеры озера составляют: в длину около 2 км, в ширину до 1 км, с площадью зеркала воды около 2 км2. Средняя глубина воды в озере составляет величину 4.5 м, максимальная достигает 12 м. Ледостав начинается в первой половине октября, сход ледяного покрова приходится на конец мая.

Озеро расположено над мерзлыми породами, мощность которых составляет 280 – 320 м, и под ним существует сквозной талик (Павлов и Тишин, 1981).

Моделирование водного и теплового режимов оз. Сырдах и нижележащих грунтов производилось за период 1965 – 1984 гг., причем в качестве входных данных использовались данные метеорологических наблюдений на ближайшей станции (г.

Якутск). В качестве начальных вертикальных профилей всех переменных состояния задавались их типичные климатические распределения для месяца начала интегрирования (декабря). Результаты моделирования сопоставлялись с данными натурных измерений (Павлов и Тишин, 1981), которые проводились на самом озере в 1976 – 1977 гг. Следует отметить, что метеостанция «Якутск» находится на расстоянии около 70 км от оз. Сырдах, поэтому атмосферные условия в описываемом эксперименте заданы с некоторой, вообще говоря, неизвестной ошибкой. Очевидно, что указанное обстоятельство является одним из источников расхождения расчетных данных с наблюдениями (см. ниже).

2.3.5. Моделирование сезонного промерзания и испарения на озере Сырдах.

Эволюция водного, ледяного и снежного покровов в данном озере по данным моделирования представлена на рис. 2.3. По результатам сравнения модельных и натурных данных (Павлов и Тишин, 1981; Pavlov, 1995) можно сделать следующие выводы:

- максимальное промерзание модельного озера происходит в малоснежные зимы, а минимальное промерзание – в многоснежные зимы, что соответствует закономерности, наблюдаемой в природе;

- характерные глубины промерзания (толщина слоя льда) в модельном водоеме находятся в пределах 0.7 – 1.5 м, что совпадает с реальными значениями для озёр Центральной Якутии;

- ледостав по результатам моделирования начинается в начале октября, а заканчивается в конце мая, что также согласуется с данными наблюдений;

- модельное испарение с поверхности озера за теплый период составляет в среднем 400 мм, что близко к наблюденному значению (450 мм).

2.3.6. Моделирование талика под озером Сырдах. На рис. 2.4 в виде распределения термоизоплет представлены результаты моделирования теплового режима талика и мерзлого грунта под озером (нулевая термоизоплета нанесена жирным пунктиром). Как видно из рисунка, талик устойчиво существует под озером в течение лет (1965 – 1984 гг.). Его нижняя граница колеблется от 1.2 до 2.0 м под дном озера. В то же время, согласно данным измерений (Тишин, 1979) глубина талика под оз. Сырдах достигает десятков метров, а в центральной части озера, по всей видимости, превышает 100 м. Такое несоответствие модельных и эмпирических данных объясняется тем, что образование глубокого талика - процесс, существенно превышающий по длительности период интегрирования модели (20 лет), так что нижняя граница модельного талика просто «не успевает» достигнуть глубины наблюдаемой границы.

Рис. 2.3. Промерзание/таяние и снегонакопление в оз. Сырдах по данным моделирования Рис. 2.4. Распределение температуры, °С, в грунте под оз. Сырдах по результатам моделирования. Глубина (вертикальная координата) выражена в значениях log10 ( z + 1), где z – глубина, отсчитываемая от дна озера вниз 2.3.7. Моделирование теплового режима водных масс оз. Сырдах. По данным моделирования были получены среднемесячные вертикальные профили температуры в оз.

Сырдах за период сентябрь 1976 – август 1977 гг. и затем сравнены с данными измерений за тот же период (Павлов и Тишин, 1981). Эти результаты представлены на рис. 2.5.

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

моделирования. В ходе интегрирования модели в фиксированный срок (22 ч.) в разные дни июля 1977 г. записывались вертикальные распределения турбулентных характеристик E,, k. Пример такого распределения показан в табл. 2.1. В ней также приведены типичные значения этих характеристик для водоемов, согласно обзорной работе (Wuest and Lorke, 2003).

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

(Карелия). Измерения на оз. Вендюрское (Карелия, 62°10' с.ш., 33°10' в.д.) производились в период 18-22 июля 2003 г в рамках проекта INTAS-01-2132, и в настоящее время они доступны на сайте http://nwpi.krc.karelia.ru/intas-01-2132/. Измерялись температура поверхности и дна озера, глубина слоя перемешивания, а также метеорологические величины в приводном слое, включая поток солнечной радиации и радиационный баланс.

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

Рис. 2.6. Температура поверхности оз. Вендюрское по данным измерений и по результатам моделирования (18-22 июля 2003 г.) На рис. 2.6 отражены два модельных результата: один получен с использованием процедуры инициализации начального профиля температуры типа «раскрутка» (перед началом основного периода интегрирования модель интегрируется с повторяющимися несколько раз граничными условиями первого дня), а другой – без нее. Как видно, процедура «раскрутки» позволяет заметно улучшить качество воспроизведения моделью температуры озера в первый день измерений.

На основании проведенных экспериментов можно заключить, что модель в целом удачно воспроизводит суточную и межсуточную изменчивость температуры (коэффициент корреляции между модельным и экспериментальным временными рядами составляет 0.8), мгновенные ошибки не превосходят 1.5 °С. Одной из возможных причин ошибок модели может быть наличие эффектов трехмерной динамики водоема (сгоннонагонных явлений, сейш, и др.), которые в одномерной модели не учитываются.

2.3.10. Моделирование термического режима озера Коссенблаттер (Германия).

Озеро Коссенблаттер находится в Германии, в 50 км к юго-востоку от Берлина. Средняя глубина озера составляет около 2 м. Используемые в описываемом эксперименте данные измерений на этом озере были получены в ходе эксперимента LITFASS-98 (Beyrich, 2000).

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

Задачей настоящего численного эксперимента, наряду с оценкой качества воспроизведения термического режима рассматриваемого озера, была оценка эффекта мелководья на турбулентные потоки явного и скрытого тепла. Эффект мелководья заключается в том, что при небольших глубинах водоема волны на поверхности водоема становятся круче и чаще разрушаются, приводя к турбулизации приповерхностного слоя воздуха. Последнее, в свою очередь, вызывает увеличение турбулентных потоков. Для учета этого эффекта Г.Н. Паниным была предложена соответствующая модификация формул теории подобия Монина-Обухова (Panin et al., 1996; Панин и др., 2006):

где H TSW и ETSW - поток явного тепла и испарение с поверхности на мелководье, H T и ET поток явного тепла и испарение с поверхности глубокого водоема, рассчитываемые по классическому варианту теории подобия, hw - средняя высота волн, H - глубина водоема.

Поскольку измерения волнения на оз. Коссенблаттер не производились, высота волн рассчитывалась по эмпирической формуле следуя (Давидан и др., 1985).

В численных экспериментах по воспроизведению термического режима оз.

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

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

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

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

Применение поправки на мелководье к потокам явного и скрытого тепла согласно проиллюстрировано на рис. 2.8 и рис. 2.9. На этих рисунках, в целях наглядности, приведены результаты расчета для 50-часового временного отрезка, выделенного из всего периода интегрирования. Как видно, поправка в целом является небольшой (среднее значение для потока явного тепла + 0.4 Вт/м2, для потока скрытого тепла + 7.3 Вт/м2), причем при преобладании положительных значений поправки встречаются и отрицательные (диапазон значений поправки для потока явного тепла от -8.1 до + 21. Вт/м2, для поток скрытого тепла - от – 4.7 до 63.8 Вт/м2). Относительная величина поправки для обоих потоков тепла составляет величину порядка нескольких процентов.

Максимальные значения поправки на мелководье соответствуют случаям с сильным ветром, и, соответственно, максимальным волнением на поверхности озера. Это можно видеть из сопоставления рис. 2.8 и рис. 2.9 с рис. 2.10, на котором изображен временной ход скорости ветра по данным наблюдений.

Несмотря на то, что поправка Г.Н. Панина вносит небольшое изменение в величину турбулентных потоков, ее использование при интегрировании модели на период более нескольких суток приводит к существенному изменению рассчитанной температуры поверхности озера. Как видно на рис. 2.11, к концу периода интегрирования модели температура, рассчитанная с учетом поправки (2.41, 2.42), оказалась ниже температуры в контрольном эксперименте на величину порядка 1 °С.

Температура, С Рис. 2.7. Временной ход температуры поверхности оз. Коссенблаттер (8 – 21 июня, 1998) по данным наблюдений и моделирования Поток явного тепла, Вт/м Рис. 2.8. Временной ход потока явного тепла над оз. Коссенблаттер, рассчитанный с учетом и без учета эффекта мелководий Поток скрытого тепла, Вт/м Рис. 2.9. Временной ход потока скрытого тепла над оз. Коссенблаттер, рассчитанный с учетом и без учета эффекта мелководий Скорость ветра, м/с Рис. 2.10. Временной ход скорости ветра над оз. Коссенблаттер по данным измерений Температура, С Рис. 2.11. Временной ход температуры поверхности оз. Коссенблаттер, рассчитанный с учетом и без учета эффекта мелководий (8-21 июня 1998 г.) термического режима озера Монте-Ново (Португалия). В данном эксперименте (Dutra et al., 2006) оценивалась способность модели воспроизводить временной ход температуры поверхности озера, а также вертикального профиля температуры в различные сезоны года.

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

Озеро Монте-Ново находится в Португалии (38° 42’ с.ш., 7°30’ з.д.). В качестве входных метеорологических данных в численных экспериментах были использованы данные реанализа ERA-40 (Uppala et al., 2005) на высоте 10 м. Перед тем, как использовать данные реанализа, была произведена их валидация путем сравнения с измерениями на метеорологической станции Портел, расположенной поблизости от рассматриваемого озера. Сравнивались среднемесячный суточный ход температуры и скорости ветра по данным реанализа и по данным измерений для различных сезонов года.

Оказалось, что среднемесячный суточный ход температуры воспроизводится реанализом достаточно успешно – ошибки находятся в пределах 1-2 °С. В то же время суточный ход скорости ветра, полученный по данным реанализа, существенно отличается от реального (рис. 2.12) Рис. 2.12. Суточный ход скорости ветра, осредненный для различных месяцев года, по данным измерений (пунктирная линия) и реанализа ERA-40 (сплошная черная линия), осреднение произведено за период 1990-2002 гг То, что в данных ERA-40 практически не находит отражение суточный ход скорости ветра, может быть объяснено отсутствием в модели, производящей реанализ, эффекта местных ветров. Местные ветры, развивающиеся за счет термической неоднородности подстилающей поверхности, обусловливают, как правило, максимум скорости ветра в дневные часы, что и наблюдается на рис. 2.12, за исключением июля.

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

Поскольку ветер является критическим фактором интенсивности турбулентных потоков над озером, скорость ветра по данным ERA-40 была скорректирована так, чтобы она точнее отражала реальный суточный ход.

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

На рис. 2.13 изображен временной ход температуры поверхности оз. Монте-Ново за рассматриваемый период по результатам моделирования и данным наблюдений. Как видно, модель достаточно адекватно описывает годовой термический цикл озера в условиях субтропического климата.

Рис. 2.13. Временной ход температуры поверхности озера Монте-Ново по данным измерений (ромбики) и моделирования (сплошная линия) На рис. 2.14 изображены результаты серии экспериментов по воспроизведению вертикального профиля температуры в озере в различные сезоны года. В экспериментах варьировалась глубина озера: задавалось 17 различных глубин от 6 до 14 м. Как видим, наиболее устойчивая стратификация, связанная с интенсивным радиационным прогревом верхних слоев, по результатам расчетов, имеет место в июле. Менее устойчивая стратификация имеет место в переходные сезоны года. Интересно отметить, что в октябре, в интервале глубин 6-8 м наблюдается повышенное перемешивание по вертикали, что, возможно связано с объединением поверхностного и придонного пограничных слоев. В январе наблюдается нейтральная стратификация, что согласуется с данными наблюдений.

Рис. 2.14. Результаты моделирования вертикального профиля температуры в оз.

Монте-Ново в различные сезоны года (осреднение за период 1990-2002 гг.) 2.3.12. Моделирование термического режима мелкого озера вблизи Тикси.

Данные измерений на этом озере были заимствованы их архива проекта GAME-Siberia (Ohata, 1997). Озеро является самым мелким из тех, с которыми производились численные эксперименты в настоящей работе: максимальное значение глубины за период наблюдений (июль 1998 – сентябрь 1999), составило 70 см. Площадь озера – около 1.5 га.

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

Для получения начальных данных в численных экспериментах использовалась процедура инициализации типа «раскрутки».

На рис. 2.15 приведено сравнение рассчитанной температуры поверхности озера и измеренной в ходе упомянутого проекта. Как видно, модель успешно воспроизводит все основные особенности реального временного хода температуры озера, обусловленного как суточной, так и синоптической изменчивостью состояния атмосферы. Следует отметить, что такое хорошее совпадение модельных и натурных данных было достигнуто благодаря учету теплового воздействия вечной мерзлоты на термический режим озера и использованию поправки на мелководье (2.41, 2.42). Оба фактора имеют большое значение для термического режима данного озера из-за его небольшой глубины. При пренебрежении любым из двух факторов модель существенно завышает реальную температуру озера.

Аккуратный учет теплового воздействия вечной мерзлоты на термодинамику водоема был достигнут благодаря более реалистичному заданию начального профиля температуры в слое грунта под водоемом. На нижней границе слоя грунта задавалась температура -5 °С, на верхней границе, совпадающей с дном озера 10 °С. Интерполяция между этими двумя значениями производилась посредством экспоненциальной функции на глубинах 0-8 м под дном водоема, ниже температура предполагалась постоянной и равной значению на нижней границе (-5 °С).

Температура, С Рис. 2.15. Временной ход температуры поверхности мелкого озера в районе Тикси по результатам моделирования и данным измерений (июль 1998 г.).

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

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

В частности, оказалось, что для успешного моделирования мелких озер необходимо учитывать эффект обрушения волн на мелководье (здесь подтверждены результаты, полученные ранее Г. Н. Паниным (Панин и др., 2006)), а также достаточно точно задавать начальное распределение температуры в слое грунта под водоемом. Именно учет этих двух факторов позволил реалистично воспроизвести термический режим мелкого (глубиной 0.7 м) озера вблизи Тикси. В то же время, для более глубоких водоемов (например, для оз. Вендюрского, глубиной 8 м), чувствительность рассчитываемой температуры к учету эффекта мелководий оказалась пренебрежимо малой.

Хорошее согласие с наблюдениями, которое было достигнуто при моделировании температуры снежного покрова в Колпашево, во-первых, продемонстрировало физическую адекватность используемого в модели блока снежного покрова (Володина и др., 2000), а во-вторых, подтвердило предположение о том, что наличие водоема слабо влияет на температуру снежного покрова.

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

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

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

ГЛАВА 3.

ПАРАМЕТРИЗАЦИЯ ВОДОЕМОВ В МЕЗОМАСШТАБНОЙ

АТМОСФЕРНОЙ МОДЕЛИ И ЕЕ ВЕРИФИКАЦИЯ