Моделирование двумерных и трехмерных бароклинных возмущений в стратифицированных бассейнах методами полнонелинейной численной модели Текст научной статьи по специальности «Физика»

УДК 532.5

С.В. Семин

МОДЕЛИРОВАНИЕ ДВУМЕРНЫХ И ТРЕХМЕРНЫХ БАРОКЛИННЫХ ВОЗМУЩЕНИЙ В СТРАТИФИЦИРОВАННЫХ БАССЕЙНАХ МЕТОДАМИ ПОЛНОНЕЛИНЕЙНОЙ ЧИСЛЕННОЙ МОДЕЛИ

Рассмотрена полнонелинейная численная модель MITgcm. На ее основе проведено моделирование эволюции наклонного пикноклина, показана значимость учета рельефа дна при анализе бароклинных волновых процессов. Проведено моделирование локализованного импульсного возмущения, заданного на основе солито-на слабонелинейной теории волн. Численно продемонстрирован закон спадания амплитуда в результате цилиндрической расходимости солитоноподобного возмущения.

Ключевые слова: численное моделирование, внутренние волны, солитон, цилиндрическая расходимость.

Введение

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

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

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

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

© Семин С.В., 2011.

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

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

Другим вариантом может служить система ROMS (Régional Ocean Modeling System) [8,9,10], в основе которой лежат гидростатические примитивные уравнения для расчета океанических процессов. Стоит уточнить, что вообще говоря, эта система состоит из трех моделей: океанической, атмосферной и волновой. Это позволяет ей описывать эффекты от взаимодействия процессов в различных средах. Однако такое усложнение модели привело к сильному увеличению расчетного времени и заставило разработчиков упростить математическую модель. В результате гидростатическое приближение уравнений гидродинамики не позволяет точно описывать поведение внутренних волн, потому данная модель непригодна для достижения целей этой работы. Кроме того отсутствие полной и хорошо структурированной документации усложняет процесс исследований.

В результате наш выбор пал на глобальную модель циркуляции, созданную в Массачусетсом Технологическом Институте - MITgcm (Massachusetts Institute of Technology global circulation model) [11,12,13]. Океаническая ее часть основывается на уравнениях гидродинамики трехмерной несжимаемой жидкости в приближении Буссинеска с учетом ускорения Кориолиса. Однако в отличие от ROMS, MITgcm позволяет учитывать негидростатические процессы, то есть в данном случае эволюция вертикального потока напрямую зависит от диффузионных, вязкостных и адвективных свойств жидкости, что предоставляет возможности для более точного описание динамики внутренних волн.

В основе системы MITgcm лежит грамотно спроектированное вычислительное ядро. Оно позволяет запускать расчеты, как в однопроцессорном режиме, так и распределенно (в smp-режиме или при помощи кластера). И наконец, отличная документированность теоретических основ и архитектурных особенностей системы выгодно отличают MITgcm от других рассмотренных комплексов.

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

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

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

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

Основные уравнения

В основе океанической модели программного комплекса MITgcm лежит система уравнений Навье-Стокса с приближением Бассинеска (1)-(7).

^ +—(V + У ьрт +.

яг ро

Ф,

мы

У иРмы ) = ,

ыг

дг

= -Ь,

-мы

Я*. , дРмы = С

вг

дг

д*

У иУ +—= 0

дг

Ь = ¿(0,5, г),

яг " б0, ^ = .

яг 5

(1) (2)

(3)

(4)

(5)

(6)

(7)

где V - вектор горизонтальной составляющей скорости, р0 - характерная плотность, р5 -давление на поверхности жидкости, рыг - гидростатическое давление, - негидростатическое давление, * - вертикальная составляющая скорости, еш - коэффициент негидроста-тичности (равен единице, когда учитываются негидростатическое давление), Ь - плавучесть, 0 - температура, 5 - соленость.

В случае решения примитивных уравнений негидростатическое давление не учитывается, а вертикальная составляющая скорости рассчитывается как интеграл от дивергенции горизонтального вектора скорости. В дальнейших расчетах использовался негидростатический режим, когда давление представлено тремя слагаемыми в (1), а в правой части уравнений учитываются все слагаемые (8)-(12):

{иГ - {- 20у 81п(ф)+ 20г сов(ф)}+ ЛиУ2ы + А4У4и + ^

Г ] (8)

С = -У-У и-<

г

С = -У-Уу -

УГ

г

и 21ап(ф)I

- {- 20и б1П(ф)}+ ЛиУ2У + Л4У4У +

С* = - У -У* +

и2 + V2

+ 20и соз(ф) + А

д

дг

2 +

б0 =У-[^У©]-0У-У + ^

(9)

(10) (11)

<

г

г

QS = V[£VS]-SV-V + FS

(12)

где □ - угловая скорость вращения Земли, А - коэффициент горизонтальной вязкости, А -коэффициент вертикальной вязкости, А - бигармонический коэффициент вязкости, К - тензор диффузии (диагональная матрица коэффициента диффузии), ^ - внешнее воздействие

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

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

Эксперимент с отклоненным пикноклином

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

а

б

Рис. 1. Начальное распределение плотности:

а - эксперимент с ровным дном; б - параболическое дно

Начальное распределение плотности задавалось по формуле:

p(x, z) = p0dp - tanh

lz - ax - z„

pyQ

d

pyc

(13)

где р0 - характерное значение плотности для выбранной области, й- коэффициент аномалии плотности, г - вертикальная координата, х - горизонтальная координата, а - коэффициент наклона пикноклина, г - глубина залегания пикноклина, й - ширина пикноклина. Для последующих экспериментов коэффициенты имели следующие значения: р0 = 1000 кг/м3, йр = 0,005, а = 0,0002, грус = 4,7 м, ёрус = 1 м .

Расчетное время для каждого эксперимента составляло 50 часов с шагом по времени в 5 секунд. Размеры области задавались таким образом, чтобы объем жидкости в обоих случаях был одинаковым, таким образом, область с ровным дном задавалась 6667 х 10 м, а область с параболическим дном - 10000x10 м. Пространственная дискретизация для обоих экспериментов совпадала: 500x50 точек. Поэтому шаги по горизонтальной оси были разные: 13.3 м для случая ровного дна и 20 м - для параболического; шаги по вертикальной координате в обоих экспериментах одинаковы.

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

Рис. 2. Эволюция поля плотности для эксперимента с ровным дном

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

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

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

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

Рис. 3. Эволюция плотности для эксперимента с параболическим дном

Эксперимент с распространением локализованного начального возмущения

вида солитон

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

Распределение плотности (рис. 4) задавалось по формуле:

р(х. z) = р0 - dр - 1апЬ

(---Ф. --)-)

d

(14)

рус

где п(х, г) = Г(х) • Ф(z), ^(х) - горизонтальная структура - солитон уравнения Гарднера с известными параметрами:

р(х) =-^^, (15)

1 + ВеЬ (ух)

А = В - ^ , (16)

а1

У =

(в2 - l)x2

6а ^

(17)

где В - свободный параметр, определяющий форму солитона (этим параметром и отличались приведенные эксперименты), остальные параметры принимали значения: г =-30 м,

■5 1

dmc = 4 м, dp = 0,005, c = 1,8544 м/с, р = 731,5959 м3/с, а =-0,052069 -,

с

pyc "р

л-3 1

а =-2,4859 -10-

м • с

Density anomaly distribution for t=0 h kg/m3 Density anomaly distribution for t=0 h kg/m3

Рис. 4. Начальное распределение плотности:

а - солитон с параметром B = 0,1; б - B = 0,6

Начальное поле горизонтальной скорости определялось по формуле: йФ

и(х,г) = сЕ(х) —. (18)

йг

Вертикальная составляющая высчитывалась средствами модели MITgcm согласно выбранному режиму моделирования (гидростатический, квази-гидростатический, негидростатический).

Поскольку в обоих случаях возникает множество высокочастотных возмущений, то в качестве результатов моделирования будет нагляднее представить х-1 диаграммы распространения начального возмущений. На рис. 5 и рис. 6 изображены эти диаграммы для изо-пикн, залегающих на 28-, 30- и 32-метровой глубине.

В обоих случаях область имела размеры 5000x100 м, расчетное время составляло 6 часов с шагом по времени 0.5 с, пространственная дискретизация 500 х100 шагов. Таким образом, шаг по горизонтальной оси - 1 0 м, а по вертикальной оси - 1 м.

Рис. 5. X— диаграммы эволюции изопикнических поверхностей на глубине 28, 30 и 32 м для начального возмущения с параметром В = 0,1

Рис. 6. X— диаграммы эволюции изопикнических поверхностей на глубине 28, 30 и 32 м для начального возмущения с параметром В = 0,6

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

Эксперимент с цилиндрической расходимостью солитоноподобной волны

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

На рис. 6 представлен пикноклин, залегающий на 29-метровой глубине в случае распространения возмущения с параметром В = 0,1 и В = 0,6 .

Pycnocline (-29.0 m) with density anomaly -1.2 fort=0.0 min m Pycnocline (-29.0 m) with density anomaly -1.2 fort=0.0 min m

Рис. 7. Начальное состояние изопикнической поверхности, залегающего на 29 метровой глубине:

а - солитон с параметром В = 0,1; б - В = 0,6

Наличие третьего измерения сильно замедляет расчеты, поэтому решено было уменьшить дискретизацию расчетной области до следующих параметров: 150 х150х50 точек, сама же расчетная область в этом случая имела размеры: 10000x10000x100 м. Однако в результате такого изменения стало невозможно отслеживать положение центральной изопик-ны, где плотностная аномалия равна нулю, потому далее будут представлены рисунки для изопикн, залегающих на 29- и 31-метровой глубине.

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

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

Как видно из представленных рис. 8 и рис. 9, амплитуда солитона при распространении спадает как для случая с параметром В = 0,1, так и для В = 0,6 .

Рис. 8. R-t диаграммы эволюции изопикнических поверхностей на глубине 29 и для начального возмущения с параметром B = 0.1 в трехмерном случае

31 м

Рис. 9. R-t диаграммы эволюции изопикнических поверхностей на глубине 29 и 31 м для начального возмущения с параметром B = 0.6 в трехмерном случае

Рис. 10. Законы спадания амплитуды, полученные аналитически и численно

В работе [15] теоретически и численно проанализированы трансформация цилиндрических уединенных волн в рамках уравнения Гарднера, а так же определены зависимости

амплитуд солитонов от лучевой координаты для различных сочетаний внешних параметров. На рис. 10 эти зависимости представлены черными линиями, а полученные в рассмотренных экспериментах численные данные на основе полнонелинейной модели MITgcm - синими.

Как видно из графиков рис. 10 в начале расчетов происходит адаптация полей начального возмущения и формирование полнелинейной локализованной структуры. Для возмущения B = 0,6 из-за меньшей амплитуды начального поля на рис. 10 видно неплохое согласие результатов аналитических слабонелинейных и численных полнонелинейных расчетов, тогда как для возмущения большей амплитуды (с В = 0,1) различия более заметны.

Выводы

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

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

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

Библиографический список

1. Чернов, А.Г. Динамика внутренних гравитационных волн в трехслойной среде / А.Г. Чернов,

0.Е. Полухина, А.А. Куркин // Двенадцатая Всероссийская научная конференция студентов-физиков и молодых ученых («ВНКСФ - 12», Новосибирск): материалы конф. / Новосиб. гос. ун-т, 2006. С. 624.

2. Чернов, А.Г. Исследование динамики внутренних волн в трехслойной жидкости / А.Г. Чернов, О.Е. Полухина, А.А. Куркин // XI нижегородская сессия молодых ученых. Естественнонаучные дисциплины: материалы конф. - Н. Новгород: Изд-во Гладкова О.В., 2006. С. 11.

3. Chernov, A.G. Dynamics of two-dimensional long internal gravity waves in a three-layer fluid / A.G. Chernov, O.E. Poloukhina, A.A. Kurkin // Advanced Problems in Mechanics: вook of Abstracts. XXXIII Summer School - Conference. St. Petersburg. 2006. P. 26.

4. Чернов, А.Г. Нелинейная динамика длинных внутренних волн / А.Г.Чернов [и др.] // Тринадцатая Всероссийская научной конференция студентов-физиков и молодых учёных (ВНКФ-13, Ростов-на-Дону, Таганрог): материалы конф. / Екатеринбург-Ростов-на-Дону-Таганрог: Изд-во АСФ России. 2007. Т. 1. С 548-550.

5. Bezruk, I.V. Large amplitude long nonlinear internal gravity waves in stratified basins: models and dynamics / I.V. Bezruk [at al.] // EGU General Assembly, Vienna, Austria, April 15-20, 2007. Geophysical research abstracts. V. 9. EGU07-A-05326, European Geosciences Union 2007.

6. Chernov, A.G. Nonlinear dynamics of internal gravity waves in a three-layer fluid / A.G. Chernov,

1.V. Bezruk, O.E. Polukhina // EGU General Assembly, Vienna, Austria, April 15-20, 2007. Geophysical research abstracts. V. 9. EGU07-A-05321, European Geosciences Union 2007.

7. Езерский, А.Б. Динамика солитонов возбуждаемых в резонаторах на поверхности мелкой воды / А.Б. Езерский [и др.] // Прикладная нелинейная динамика. 2004. Т. 1-2. №1-2. С. 138-158.

8. Haidvogel, D.B. Numerical Ocean Circulation Modeling / D.B. Haidvogel, A. Beckmann // Numerical Ocean Circulation Modeling. 1999. Р. 318.

9. Song, Y.T. A general pressure gradient formulation for ocean models. Part I: Scheme design and diagnostic analysis // Monthly Weather Review (MWR), 1998. Р. 3213-3230.

10. Song, Y.T. A general pressure gradient formulation for ocean models. Part II: Energy, momentum and bottom torque consistency / Y.T. Song, D.G. Wright // Monthly Weather Review (MWR). 1998. Р.3231-3247.

11. Adcroft, A. Overview of the formulation and numerics of the MITgcm / A. Adcroft [at al.] // Proceedings of the ECMWF seminar series on Numerical Methods, Recent developments in numerical methods for atmosphere and ocean modeling. 2004. Р. 139-149.

12. Adcroft, A. Numerical algorithms for use in a dynamical model of the ocean / Ph.D dissertation, 1995.

13. Hill, C. Application of a parallel navier-stokes model to ocean circulation in parallel computational fluid dynamics / C. Hill, J. Marshal // Implementations and Results Using Parallel Computers, 1995. Р. 545-552.

14. Езерский, А.Б. Возбуждение солитонов внутренних волн в резонаторах / А.Б. Езерский, О.Е. Полухина, Е.Н. Пелиновский // Изв. АИН РФ. Прикладная математика и механика. 2004. Т. 6. С. 90-97.

15. Полухина, О.Е. Цилиндрическая расходимость уединенных волн в рамках обобщенного уравнения Гарднера / О.Е. Полухина, Н.М. Самарина // Изв. РАН. Физика атмосферы и океана, 2007. № 6 Т. 43. C. 818-825.

Дата поступления в редакцию 02. 07.2011

S.V. Semin

TWO AND THREE DIMENSIONAL BAROCLINIC PERTURBATIONS MODELING IN CLOSED STRATIFIED BASINS BY METHODS OF FULLY NONLINEAR NUMERICAL MODEL

In the present work fully nonlinear numerical model MITgcm considered. On its basis modeling of evolution tilting pycnocline spent, the importance of a relief of a bottom is proved at carrying out of similar experiments. Modeling of fully nonlinear soliton, set on a basis of weakly nonlinear theory of waves spent. The law of falling off amplitude as a result of solitons cylindrical divergence numerically proved.

Key words: numerical modeling, internal waves, soliton, cylindrical divergence.