Моделирование отклика верхнего слои океана на атмосферное воздействие

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

Эта модель была реализована для моделирования эволюции характеристик верхнего слоя в Северной Атлантике (Мошонкин, Дианский, 1992). В работах (Дианский и др., 1999; Глазунов и др., 2001) она использовалась для средних широт Северного полушария в Атлантическом и Тихом океанах как океанический блок совместной модели для изучения атмосферного отклика на аномалии ТПО в средних широтах. Более ранние версии модели тестировались в условиях высоких широт и тропиков. В работе (Залесный, Мошонкин, 1988) с использованием данных океанской станции погоды «Bravo» (56,5° с.ш.) показано, что модель дает удовлетворительные результаты в Субарктике (инверсия температуры в термоклине), когда глубокая конвекция развивается над слоем мощной термической инверсии, гидростатическая устойчивость в котором сохраняется стратификацией солености. В работе (Мошонкин, Харендупракаш, 1991) модель хорошо воспроизводила характеристики ВПС в Бенгальском заливе, когда принципиален учет потока плавучести, связанного с распреснением с поверхности речными водами.

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

Используем для апробации модели ежесуточные данные наблюдений на океанской станции погоды (ОСП) «Charlie» («С»), которая находится в 50- 100 км к северу от субполярного мегафронта и Северо-Атлантического течения (52,75° с.ш., 35,5° з.д.) в полосе западного переноса в атмосфере, где повторяемость циклонов максимальна (так называемая зона «шторм-треков»). В Сборнике климатолого-статистических данных по ОСП «С» (разделы 1 и 2, 1984) приведены ежечасные актинометрические и метеорологические наблюдения, профили температуры и солености на каждые три часа за период с 01.01.1976 г. по 31.12.1980 г. Все эти данные с помощью осреднения сводились к суточному шагу по времени (1827 сут) (Мошонкин, Дианский, 1994). Пропусков в данных после 15.11.1976 г. не было. Горизонтальные градиенты температуры воспроизводились по данным объективного анализа поля ТПО, который проводился раз в пять суток в ГМЦ РФ. В середине февраля 1978 г. происходило вторжение в район ОСП субарктического мегафронта почти месячной продолжительности (Мошонкин, 1982). Поэтому для моделирования выбран период в 445 суг с 15.11.1976 г. но 02.02.1978 г. В качестве начальных условий использовались профили температуры и солености на 15.11.1976 г.

На рис. 6.6 показаны потоки Qo, Q, и скорость ветра на ОСП «С», поступающие на вход в модель в качестве граничных условий. Кроме температуры и солености для сравнения с модельными результатами привлечена реализация ежесуточной толщины ВПС. Проведен тщательный подбор параметров модели (см. прил. 6.2). На рис. 6.7 показаны реализации в 445 сут с 15.11.1976 г. по 02.02.1978 г. натурных и модельных температур воды и положения нижней границы ВПС при оптимальных параметрах. Среднеквадратичная ошибка модельного прогноза температуры ВПС равна 0,29°С. Экстремальные ошибки связаны с горизонтальной неоднородностью характеристик ВПС в данном районе и адвекцией, которые не учтены в настоящем модельном эксперименте. Интересно отметить, что для модельной ошибки при воспроизведении сезонного цикла на станции «Рара», где адвекция считается незначительной, получена среднеквадратичная ошибка 0,26°С (Alexander, Penland, 1996).

Многолетнее воспроизведение процессов в верхнем слое океана одной лишь локальной моделью затруднительно из-за незнания адвекции. К сожалению, мы не можем непосредственно и с необходимой точностью восстановить полную адвекцию тепла в этот период. Тем не менее по данным наблюдений реконструируются средние горизонтальные градиенты температуры в ВПС путем билинейной интерполяции пентадных данных объективного анализа ГМЦ СССР в пятиградусной окрестности станции «С» (Moshonkin, Diansky, 1995). Поэтому непосредственно для модели может быть рассчитана адвекция аномальными горизонтальными дрейфовыми течениями в поле средней температуры океана (см. прил. 6.2).

Приведем пример многолетнего использования балансовой модели для воспроизведения аномалий температуры воды в ДСО для разных временных масштабов. Для этой цели применялась схема моделирования, согласно которой внутрисезонные аномалии в верхнем слое океана Т'_т„Л^г->^+) Д^дя момента 1_/+[ -t_i + Д/ вычисляются по схеме, показанной на рис. 6.8. В модель поступает суммарный профиль Г(д,С) = Т'„юА²^,)+ 7,л(²^,-) • На вход модели подаются текущие значения аномалий бюджета тепла поверхности океана () , и аномалий дрейфовой адвекции с(_а^ (см. прил. 6.2), а также суммарная динамическая скорость в воде у поверхности океана. В результате рассчитывается новый профиль Г(д,/.₊₁) и вычисляется искомый профиль аномалий в деятельном слое на момент /,_+|: Т'пшА², *_м) = Дд,6_+!) - Тоь(²,Ь) • ^Затем формируется профиль для следующего расчетного шага: Г(г,) = Т'_тоА², *₁₊,) + ^?,₊₁) • Обратим

здесь внимание на различие суммарных профилей Т„ и Т. Интересно отметить, что все другие рассмотренные схемы были неустойчивы.

Рис. 6.6. а - потоки тепла Оо (1) и Ог (2) и б - модуль скорости ветра на 27 м на ОСП «С» с 15.11.1976 г. по 02.02.1978 г., шаг - 1 сут

Рис. 6.7. Временной ход температуры (а) и нижней границы В11С (б) по данным наблюдений (I) и модельному прогнозу (2) ( = 8 , с = 35 , =8,2 = 0,052 м ¹ , Я = 300 м).

Станция погоды «С», с 15.11.1976 г. по 02.02.1978 г., шаг- 1 сут

На рис. 6.9, а показаны межгодовые аномалии ТПО. В схеме эксперимента (рис. 6.8) Т_оЬ(0д) составляется из этих аномалий ТПО и климатического хода. На рис. 6.9, 6 демонстрируется способность модели воспроизводить вну гри- сезонные аномалии Т' за пятилетний период наблюдений на станции погоды

«С» (1827 сут с 01.01.1976 г. по 31.12.1980 г.). Здесь был использован полный наблюденный профиль солености при расчете. Рис. 6.9, 6 показывает хорошее согласие с наблюдениями, что говорит о том, что в средних широтах внутрисе- зонные аномалии температуры поверхности океана в основном формируются квазилокальным атмосферным воздействием на океан. Сравнение межгодовых и сезонных аномалий ТПО показывает, что они имеют приблизительно одинаковую амплитуду.

Рис. 6.8. Схема многолетнего моделирования вну грисезонных аномалий характеристик верхнего слоя океана

На рис. 6.9, в показан спектр внугрисезонных аномалий ТПО (рассчитан через корреляционную функцию с использованием спектрального окна Парзена (Дженкинс, Ватте, 1971)). Спектральные максимумы в интервале частот менее 0,12 сут ¹ значимы с 95%-ной вероятностью. Особо выделяются пики на частотах 0,021 сут ¹ (период - 48 сут) и 0,05 суг ¹ (20 сут). Остальные экстремумы (периоды 12,8; 10,4; 9,1; 7,8 и 6,7 сут) находятся в пределах продолжительности 2-4 естественных атмосферных синоптических периодов. Скорость уменьшения спектральной плотности, в целом, обратно пропорциональна квадрату частоты. Это вполне согласуется с положением об красношумном отклике океана (см., например, (Hasselmann, 1976; Frankignoul, Hasselmann, 1977; Frankignoul, 1985; Newman, Story, 1985; Пигербарг, 1989)).

Рис. 6.9. Временной ход межгодовых (а) и вну грисезонных (6) аномалий ТОО по данным наблюдений (сплошная линия) на станции погоды «С» за период с 1.01.1976 г. но 31.12.1980 г. (1827 сут). На 1рафике (о) показаны воспроизведенные аномалии ТПО с помощью локальной модели верхнего слоя океана (пунктир) с учетом реального среднего годового хода температуры и солености в слое 0-300 м. Спектральная плотность внутрисезонных аномалий ТПО (в)

Представленная здесь модель верхнего слоя использовалась также в работах (Мошонкин, Дианский, 1992, 1993) для моделирования внугрисезонных аномалий ТПО для зимнего сезона на акватории Северной Атлантики (20- 62,5° с.ш.) по данным, полученным в ходе Первого глобального эксперимента ПИГАП (ПГЭП), проведенного с декабря 1978 г. по декабрь 1979 г. По данным ПГЭП были рассчитаны в указанной акватории все компоненты потока тепла с суточным шагом. Так же с суточным шагом по данным ПГЭП были рассчитаны компоненты напряжения фения ветра. Изменениями солености за счет осадков и испарения пренебрегалось. Компоненты аномального переноса рассчитываются через аномалии напряжения ветрового трения по нестационарной модели Экмана (Гилл, 1986). При численной реализации адвективные компоненты правой части (6.18) из прил. 6.2 вычислялись с предыдущего шага по времени, с использованием метода направленных разностей при нулевых потоках на боковых границах. Полные потоки аномального экмановского дрейфа рассчитывались в предположении постоянства правой части внутри шага интегрирования путем точного решения уравнений экмановской модели. Инерционные колебания при этом исключались на каждом временном шаге, поскольку его длина составляла одни сутки. Эксперименты показали, что при гаком подходе нестационарная модель Экмана близка к стационарной, поскольку главная изменчивость в полных потоках дрейфовой скорости определяется изменчивостью поля вегра. Адвекция аномалий температуры Т’ средним течением вычислялась на основе данных о климатических течениях, рассчитанных для зимнего сезона методом адаптации в работе (Демин и др., 1991). Составляющие средней климатической скорости и и V были получены с помощью осреднения исходных

массивов горизонтальных скоростей в верхнем слое океана (/? =75 м). Толщина этого слоя близка в среднем по акватории Северной Атлантики к величине, соответствующей зимним условиям (Мошонкин, Дианский, 1992).

На рис. 6.10 показаны временные коэффициенты Фурье ЗУЭ-мод в изменчивости модельных и натурных Т', по данным наблюдений ПГЭП для периода с 1 декабря 1978 г. но 2 марта 1979 г. в области 20-62,5° с.ш. Здесь модельные аномалии температуры поверхности океана были рассчитаны с использованием процедуры коррекции потока тепла на границе «атмосфера - океан» (Мошонкин, Дианский, 1992). Коэффициент корреляции между приведенными кривыми равен 83%, что свидетельствует о хорошей согласованности модельных и натурных Т'.

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

Рис. 6.10. Временной ход коэффициентов Фурье для первых мод сингулярного разложения между натурными и модельными температурными аномалиями, полученными но данным Г1ГЭГ1 за период с 01.12.1978 г. по 02.03.1979 г. (в безразмерных единицах)

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

Кроме того, уравнение для кинетической энергии турбулентности в этой модели, в силу необходимости интегрировать его в пределах ВПС, записывается в упрощенном виде. Поэтому сейчас в моделях циркуляции используется ряд параметризаций коэффициентов турбулентности, в которых коэффициенты зависят от локальных характеристик стратификации и вертикального сдвига скорости течения (см., например, (Монин, Обухов, 1953; Pacanovsky, Philander, 1981 и др.)).

Такие параметризации можно рассматривать как некие асимптотические представления для уравнения КЭТ.

Рис. 6.11. Пример средних за год модуля скорости среднего течения (слева) и модуля скорости

турбулентных пульсаций (справа, 'Jlk ) для Северной Атлантики, см/с, координаты модельные; слой 0-30 м; 1957 г.

Система уравнений модели циркуляции не описывает турбулентные процессы («подсеточная» физика), поэтому необходимо либо параметризовать коэффициенты обмена, либо добавить дополнительные уравнения, воспроизводящие эффекты турбулентности. Заметим также, что в моделях циркуляции океана расчет потоков тепла на границе атмосферы и океана идет, как правило, с использованием модельной температуры поверхности океана, а она принципиально зависит от турбулентного перемешивания. Таким образом, задание граничных условий также во многом определяется типом используемой модели турбулентности. Современное состояние моделирования турбулентных процессов в моделях циркуляции, использующихся для реконструкции изменчивости климата или средних климатических состояний, дано, например, в работах (Burchard et al., 1999; Warner et al., 2005; Мошонкин и др., 2014). Наиболее популярной моделью турбулентности, используемой сейчас в моделях циркуляции, является модель Мелора-Ямады уровня 2,5 (см. замечания о ней в (Warner et al., 2005)). Она используется в широко применяемой модели РОМ. Это положение зафиксировано на сайтах таких моделей, как MOM, ОРА, РОМ. Отмечается, что использование такого рода модели турбулентности дает наилучшие результаты, но очень существенно увеличивает время расчетов. Поэтому актуально как применение модели турбулентности достаточно полной физической формулировки (далее - МТ), так и использование такой численной схемы, которая позволила бы построить наиболее экономичный расчетный алгоритм.

Такой алгоритм оказывается возможным построить на основе методов многокомпонентного расщепления по физическим процессам и геометрическим координатам, представленных, например, в (Zalesny et al., 2010). Важно, что при этом численная схема МТ аналогична схемам модели циркуляции, для подюиочения в которую она предназначена. Это делает блок турбулентности адаптивным к модели циркуляции и позволяет использовать в МТ программные модули модели циркуляции (например, модули переноса и диффузии). Описание физической постановки и численного алгоритма МТ дано в прил. 6.3. Данная МТ присоединена к модели крупномасштабной циркуляции океана в ^-системе координат с приведенной глубиной сг = (д — ?)/ (Н -?), где г - обычная вертикальная координата; Н - глубина места; д - уровень океана. Модель циркуляции основывается на полных «примитивных» уравнениях для горизонтальных компонентов скорости течения, гидростатики, неразрывности, эволюции тепла и соли и уравнении состояния (см. физическую формулировку и описание методов решения, например, в (Мозйопкт е! а1., 2011)). В качестве новой системы координат вместо обычной географической выбирается двухполюсная ортогональная система с целью избежать численного шума у Северного полюса. Полюсы новой системы расположены на географическом экваторе за пределами расчетной области в точках с координатами 120° з.д. и 60° в.д.

Расчеты проведены для области Северной Атлантики (севернее 30° ю.ш.) с морями, Северного Ледовитого океана (СЛО) и Берингова моря. В Тихом океане жидкая граница проходит по проливам Алеутских островов. Модель имеет пространственное разрешение (1/4)°. Топография дна океана берется из массива данных ЕТ0Р05. Модельная глубина ограничена минимальной величиной 5 м. Шаг интегрирования по времени один час. По вертикали задаются 40 ст-уровней с учащающимся к поверхности океана разрешением в верхнем слое, что позволяет лучше описывать наиболее интенсивные здесь турбулентные процессы. На жидких границах от поверхности до дна задаются ежемесячные климатические значения температуры и солености по данным наблюдений. Учет стока основных рек, влияющих на распределение солености, осуществляется подобно выпадению осадков вблизи устьев рек. Коэффициенты крупномасштабной боковой диффузии для температуры и солености берутся одинаковыми (как функции глубины и широты), с максимальным значением 2 10⁶ см²/с на физическом экваторе и плавным уменьшением в два раза с глубиной и к высоким широтам.

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

Для изучения качества воспроизведения климатической изменчивости океана численные эксперименты проведены как с МТ, так и с параметризацией (Pacanovsky, Philander, 1981) на 60 лет для 1948-2007 гг. Для задания граничных условий (напряжение трения ветра, потоки тепла и влаги на каждый шаг модели циркуляции по времени, равный часу) использованы атмосферные характеристики из последней модификации массива данных CORE (см. сайт). Для коротковолновой радиации учитывается ее проникающая способность. Результаты расчетов, проведенных по отмеченной акватории Атлантики и Арктики, в районе упомянутого в предыдущих параграфах корабля погоды «Charlie» (57,75° с.ш. и 35,5° з.д.), расположенного в зоне частого прохождения атмосферных циклонов (так называемая зона «шторм-треков»), сравнивались с ежедневными без пропусков наблюдениями по температуре и солености в течение полного годового цикла (слой 0-2 км, пример 1979 г.) (Мошонкин и др., 2014). Показано, что МТ количественно и качественно улучшает воспроизведение термохалинной структуры океана по сравнению с асимптотическими параметризациями (Монин, Обухов, 1953; Pacanovsky, Philander, 1981).

Рис. 6.12. Изменения в поле скорости течений при замене параметризации (Pacanovsky, Philander, 1981) в модели циркуляции (слева) на модель турбулентности (Мошонкин и др., 2014) (справа). Гольфстрим, средние в слое 0-100 м за десять лет (1998 2007 гг.). Координаты модельные, ордината -44 соответствует мысу Гаттерас. Изолинии 50 и 100 см/с дополнительно маркируют величины модуля скорости

Временные масштабы турбулентности менялись от секунд, в слоях высоких градиентов плотности, до 8-9 мин в хорошо перемешанных слоях, существенно варьируя по пространству. Частота диссипации КЭТ в слое 0-50 м, где турбулентность максимальна, наиболее часто менялась в диапазоне от 0,002 до 0,060 Гц. На рис. 6.11 показано типичное пространственное распределение модуля скорости среднего течения и модуля скорости турбулентных пульсаций в Северной Атлантике для приповерхностного слоя океана. В климатическом плане основные источники КЭТ приурочены к наиболее интенсивным течениям. В осенне-зимний период выделяется устойчивая область максимальных значений КЭТ у побережья Северной Америки, обусловленная как прямым ветровым воздействием, так и конвекцией из-за охлаждения океана холодной воздушной массой, выносимой с берега на океан. Проведены эксперименты с локальным и полным вариантами МТ. Сравнение результатов показало очень существенную роль адвекции КЭТ течениями. На большей части акватории величины турбулентных пульсаций yflk равны 1-8 см/с, что в ряде районов сопоставимо с величиной средней скорости течения.

Рассмотренная МТ позволяет более адекватно, чем асимптотические параметризации, воспроизводить особенности пространственного распределения температуры и солености в океане на больших масштабах времени. Это столь существенно меняет распределение плотности воды, что заметно меняются даже крупномасштабные особенности циркуляции. В качестве примера на рис. 6.12 показана средняя за 1998-2007 гг. «модельная» циркуляция в районе Гольфстрима. Видим, что при использовании МТ растет область высоких скоростей, а поворот Гольфстрима от берега воспроизводится на 150-200 км ближе к мысу Гаттерас (климатическая средняя «точка отрыва» Гольфстрима), тогда как модель циркуляции с параметризацией (Pacanovsky, Philander, 1981) «переносит» этот поворот существенно на север, что не соответствует действительности. МТ показала принципиально лучшее качество воспроизведения термоха- линной и динамической структуры океана, что позволяет получить новые не только количественные, но и качественные результаты при воспроизведении изменчивости климата океана. В частности, использование МТ позволило качественно улучшить воспроизведение таких важных климатических особенностей циркуляции, как поворот Гольфстрима у мыса Гаттерас в Атлантике и вихрь Бофорта в Арктике. Заметим также, что использование МТ приближает поле уровня океана к наблюдаемым величинам и конфигурациям.