Объединенная численная модель течений, волнения и транспорта наносов озера Донузлав

Математическое моделирование

морских систем

УДК 551.465

В.В. Фомин, В.А. Иванов

Объединенная численная модель течений, волнения и транспорта наносов озера Донузлав

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

Озеро Донузлав является акваторией активного хозяйственного использования. Анализ связанных с этим экологических последствий требует специальных комплексных исследований течений, волнения и транспорта наносов в этом районе. До настоящего времени известна лишь одна работа, посвященная расчетам циркуляции в озере [1]. Динамические процессы в озере в значительной мере определяются ветровой циркуляцией вод и индуцированными поверхностным волнением осциллирующими течениями. Раздельное моделирование этих процессов имеет ограниченных научный интерес и не отвечает современным требованиям комплексного гидрометеорологического обеспечения хозяйственной деятельности в данном регионе. Одним из возможных путей повышения адекватности моделирования динамики озера является учет взаимодействия между течениями и волнами. Поскольку ветровые течения и поверхностное волнение имеют несоизмеримые пространственно-временные масштабы, прямое численное моделирование этих процессов в рамках одной системы гидродинамических уравнений при современном уровне развития вычислительной техники является практически нереализуемой задачей. Альтернативным подходом может быть совместное использование моделей циркуляции вод и волновых спектральных моделей с одинаковым пространственно-временным разрешением в комбинации с полуэмпирическими теориями пограничных слоев. В последние несколько лет такие объединенные модели (<coupling models) начинают находить применение в задачах многокомпонентного регионального моделирования [2-4].

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

© В.В. Фомин, В.А. Иванов, 2006

мерной численной модели течений [5,6] и модифицированном варианте спектральной волновой модели SWAN (Simulating Waves Nearshore) [7].

1. Описание модели

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

1.1. Поля течений

Для расчета полей течений используется гидродинамический блок, базирующийся на примитивных уравнениях движения вязкой несжимаемой жидкости в приближении Буссинеска и гидростатики. Следуя [8], мгновенные характеристики полей представляются в виде суммы компонент: осредненной (<р ), волновой (фу,) и компоненты, описывающей турбулентные флуктуации

(?'):

Ф = р + <р„+<р'9

где черта означает осреднение по характерному для волновой компоненты периоду времени. Предполагается также, что ср' = 0, а волновая компонента

не коррелируется с турбулентными флуктуациями: ср„(р' - 0. После подстановки (1) в уравнения Навье - Стокса и осреднения слагаемых получаем систему уравнений для осредненного движения (знак осреднения и термин «ос-редненный» для простоты записи и изложения далее опускаются):

du r 1 др д , ч д , . д / л 1

--/v +--- = — (тхх) + — (г ) + — (тxz) н

dt р дх дх ду dz р

dv 1 др д , ч д , ч д . ч 1 , .

йг р ду дх } ду дг р

др ди ду ди>

— + gp = 0 , — +— +-= 0.

дг дх ду дг

Здесь и,у,\у - компоненты скорости течений по осям х,у,г соответственно (ось г направлена вертикально вверх); d/dt -ид!дх + уд!ду + м>д!дг - субстанциональная производная; р- давление; р - средняя плотность воды; g - ускорение силы тяжести; / - параметр Кориолиса. Следует отметить, что упрощенный вид уравнения баланса импульса по вертикали (гидростатическое приближение) справедлив вне узкого придонного пограничного слоя.

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

г>у=2

ду

Тхх ~ { , тху - тух - Ам

ду ди — +—

дх ду

;

ТХ2 - ' ^

Коэффициент вертикальной турбулентной вязкости Км в (6) рассчитывается по дифференциальной модели Меллора - Ямады [9]. Коэффициент горизонтальной турбулентной вязкости Ам вычисляется с помощью формулы Смагоринского [10].

Система уравнений (2) - (4) отличается от обычно используемой при моделировании циркуляции наличием дополнительных слагаемых: м + | дру„и„ [ др\у„и„

х дх ду дг ' ^

м | др„ + ру„у„ | др\у„у„

у дх ду дг

представляющих дополнительный поток импульса в жидкости, создаваемый волновой компонентой движения [11]. Эти величины определяются из решения волновой задачи (пункт 1.2).

Допустим, что бассейн ограничен свободной поверхностью г/(х,у,г) и

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

Граничные условия по вертикали имеют вид

дг/ дт] д^ „ д\3 та ™ = + + > при г =

дг дх о у дг р

дН дН Ж дМ ть __

* = -и-—V—, Км- = -± при 2 — —Н , (9) дг дх д у дг р

где и = (м, у) - вектор горизонтальной скорости течений.

Касательные напряжения на поверхности раздела воздух - вода та вычисляются по скорости ветра и10 на высоте Юм, коэффициенту трения Сй и плотности воздуха ра на основе зависимости:

та=Р„С, |и10|и10.

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

Коэффициент Сй считается зависящим от скорости ветра и характеристик волнения (пункт 1.3), а Сь определяется с учетом взаимодействия волн и течений в придонном пограничном слое (пункт 1.4).

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

К задаче (2) - (11) применяется а -преобразование вертикальной координаты:

* =

где D = H + r] - динамическая глубина. Далее преобразованные уравнения представляются в виде двух систем. Первая система уравнений описывает баротропную компоненту решения (осредненные по глубине горизонтальные компоненты скорости течений и уровень), вторая - отклонения горизонтальных компонент скорости от соответствующих осредненных компонент и вертикальную скорость. Полученные системы решаются численно на прямоугольной сетке С с шагами Ах, Ау по горизонтали. Для более детального разрешения пограничных слоев по вертикальной координате модель допускает использование неравномерного шага по вертикали. Аппроксимации операторов адвекции проводятся с помощью монотонной схемы Superbee, построенной на основе принципа уменьшения общей вариации решения TVD (Total Variation Diminishing) [12]. Как известно, применение этого класса схем дает возможность удовлетворительно моделировать гидродинамические поля при резких изменениях их пространственных градиентов. Для интегрирования по времени в модели применяется метод предиктор-корректор. Предиктором служит центрально разностная схема, а корректором - схема Кранка - Ни-колсона. Причем система для баротропной компоненты решается с более мелким шагом по времени, величина которого определяется условием Куранта для длинных волн. Следует отметить, что применение метода предиктор-корректор дает возможность практически повысить устойчивость алгоритма и проводить расчеты с достаточно большим шагом по времени по сравнению с другими методами интегрирования, хотя и требует дополнительных вычислительных ресурсов.

1.2. Поверхностное волнение

Для описания полей волнения используется концепция частотно-углового спектра Е {cr,6,x,yj), где а,в - частотная и угловая координаты; x,y,t - горизонтальные координаты и время. Соответствующее эволюционное уравнение для определения Е имеет вид закона сохранения энергии [7]:

Iw+ffcAO+lM+Trk")*^)^. <13)

ot дх ду дсо дО со

Здесь N -Е/со - плотность волнового действия; величины сх,су,са),св -скорости переноса вдоль соответствующих пространственных и частотно-угловых координат. Выражения для этих скоростей определяются в предположении медленно меняющихся характеристик топографии дна и течений на основе соотношений геометрической оптики [7].

В данной работе для численного решения уравнения (13) применяется модифицированный вариант модели SWAN [7], в которой функция источника определяется в виде

S = Sin+Snl + Sm+Sbf + Sdib,

где слагаемые имеют следующий смысл: Sin - источник генерации волн ветром; Snl - нелинейные взаимодействия спектральных гармоник; Swc - диссипация энергии вследствие обрушения волновых гребней (whitecapping); Sbj ~ диссипация энергии, обусловленная донным трением; Sdip - обрушение волн

на критических глубинах (depth-induced breaking). В оригинальном варианте SWAN источник генерации волн Sin определяется в предположении, что коэффициент трения Cd зависит только от скорости ветра U10. В модифицированном варианте используется обобщение параметризации учитывающее его зависимость не только от скорости ветра, но и от характеристик волнения (пункт 1.3). Это приводит к тому, что в волновой и циркуляционной задачах касательные напряжения ветра вычисляются единообразно.

Решение задачи (13),(14) дает возможность получить оценки ряда спектральных характеристик волнения, в частности высоту значительных волн (significant wave height) hs и среднеквадратическую высоту волн hr:

А, =4^, Лг=л&5/2,

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

2 /г(г)2

О л»,

а также оценить волновые компоненты импульса в уравнениях гидродинамического блока. Следуя [4], будем задавать в (2),(3) средние по глубине значения Мх, Му, что для условий мелкой воды является оправданным допущением. С учетом сказанного имеем следующие выражения:

D

дх ду

дх ' ду

£>

Здесь Я у - симметричный тензор волновых радиационных напряжений, компоненты которого имеют вид [11]:

+ 5уу=8р((пвт2е + п-\/2)Е),

^ху =$ух = Вр(пьтвсо$0Е).

В (18) величина и^^/с^ , где сё,ср - групповая и фазовая скорости,

определяемые из дисперсионного соотношения для линейных гравитационных волн.

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

/55// 0233-7584. Мор. гидрофиз. жури., 2006, №2 47

чения. Для определения придонного трения Sbf используется эквивалентный масштаб шероховатости дна (пункт 1.4).

1.3. Коэффициент поверхностного трения

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

с [ 1,14 х10~3 при |U10|< 10м/с,

J [(0,49 + 0,065|U10|)xl0"3 при |ü10|> 10м/с.

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

Cd=K2/ln2

где л: = 0,4 - постоянная Кармана, z = 10m, а z0- параметр шероховатости морской поверхности, характеризующий гидродинамические свойства подстилающей поверхности и определяющий режим турбулентного перемешивания на границе раздела атмосфера - море. Как видно из (10) и (20), даже при постоянной скорости ветра шероховатость поверхности моря может приводить к пространственной неоднородности касательных напряжений ветра.

Определение реальных значений z0 является предметом интенсивных экспериментальных исследований (см. обзор в [15]). Эти исследования показывают, что отношение z0/hr хорошо коррелируется с крутизной // и возрастом волн cw/|U10|, где ст - фазовая скорость, соответствующая максимуму в волновом спектре. Предложен ряд формул для определения z0, использующих один или оба из указанных параметров. Эти формулы можно представить степенной зависимостью вида:

z0=ahr

г \'т ст

/Л

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

1.4. Коэффициент придонного трения

В предположении логарифмической зависимости профиля скорости течений вблизи дна для коэффициента придонного трения Сь имеем выражение, подобное (21):

Cb=K2/ln2

Здесь хъ - параметр шероховатости донной поверхности; 2 - расстояние от дна до точки, в которой определяется коэффициент трения. Если характеристики волнения не учитываются, то

гь=кь!30,

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

При наличии волновых движений, в соответствии с теорией Гранта -Мадсена [16], модифицированной к условиям мелкой воды, параметр шероховатости определяется из следующих соотношений:

1Р

•KJ 30,

к be - кь

24

^ Л £/„. к.

С/,

♦Г

и.

Здесь к.

Ьс

U. = , U.c = f u.w = Jfj2Uw.

- так называемый кажущийся масштаб шероховатости (apparent bottom roughness); Uф - скорость трения при взаимодействии волн и течений; Umc - компонента скорости трения, обусловленная только течениями; Uw-амплитуда скорости придонных волновых течений; fw - коэффициент волнового сопротивления;^ =Uw/g) ; со = 2п/Тр\ Тр- период волн, соответствующий максимуму энергии в волновом спектре.

Алгоритм вычисления коэффициента трения следующий. Вначале оценивается Uw с помощью линейной теории волн:

Uw=ao)/sh

где a-hsl 2 - амплитуда волны, к - волновое число, определяемое по дисперсионному соотношению со . Далее находится /„. по эмпирическому соотношению:

1

+ log

10

1

= -0,08 +log

ю

После этого коэффициент трения определяется из (22),(24),(25) методом итераций. В качестве начального приближения выбирается Сь с предыдущего временного шага. В соответствии с [16] величина 2 в (22) считается заданной и равной расстоянию от дна до ближайшей расчетной точки по вертикали.

1.5. Транспорт наносов

Изменения профиля донной поверхности Н(х,у,/) за счет процессов эрозии (размыва) и аккумуляции (намыва) дна определяются на основе решения уравнения деформации вида [17]:

дН_ dt

1

дх ду

где - компоненты вектора потока (расхода) наносов q \ р,, а - плот-

ность и пористость донных осадков. Суммарный поток наносов складывается из потоков донных (ць) и взвешенных ^Л.)наносов. Для расчета этих компонент используется подход, основанный на методике Ван-Рийна [17, 18]. Согласно этой методике,

Здесь те- вектор эффективных придонных тангенциальных напряжений [17];

мерный диаметр частиц; V- коэффициент молекулярной вязкости морской воды; \]а - осредненный по глубине вектор горизонтальной скорости течений; В-динамическая глубина; F- расчетный параметр, учитывающий вид вертикального распределения взвешенных наносов, определяется согласно [18]; С/7 - объемная концентрация взвешенных наносов у дна, вычисляемая с учетом скорости оседания частиц. В отличие от [17], в настоящей работе при вычислении тс вместо Иа используется \]с. Физически это более оправданно, поскольку из того, что \}а = 0, не следует равенства нулю скорости течений в придонном слое, а именно эта скорость вносит основной вклад в придонные касательные напряжения.

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

Скорости течений и частотно-угловой спектр волнения рассчитывались на сетке с горизонтальным разрешением Ах =Ду=200м. Шаг интегрирования по времени Д/ = 5 мин. По сг-координате использовалось 11 уровней с постоянным шагом. Значения баротропных компонент скорости течений, уровня и характеристики волнения (высота, направление и период) на жидких границах расчетной области определялись линейной интерполяцией из решения вспомогательных задач для всего Черного моря на сетке размером 5 км. Отклонения горизонтальных компонент скорости от баротропных компонент на жидких границах рассчитывались с помощью условий излучения. Мини-

тс - критическое тангенциальное напряжение; Л /]13 - безраз-

2. Реализация модели и выбор параметров

мальное значение коэффициента вертикальной турбулентной вязкости равнялось 10"4 м2/с. Для расчета коэффициента горизонтальной турбулентной вязкости безразмерная константа Смагоринского принималась равной 0,2.

При расчете волнового спектра применялась частотно-угловая сетка с разрешением 15° по угловой координате и с 26 неравномерно расположенными узлами по частоте. При этом минимальная частота волн сох =0,05Гц. Параметр шероховатости морской поверхности определялся с учетом характеристик волнения по формуле (21) при следующих значениях констант [15]: а = 8072 ; т = 0 ; п = 4,5 . Крутизна волн оценивалась отношением кг! Ьт, где Ьт — длина волн, соответствующая максимуму энергии в волновом спектре. Выбор этих значений коэффициентов связан с тем, что для условий мелкой воды шероховатость морской поверхности в большей степени зависит от крутизны волн, а не от их возраста. Максимальное значение высот волн перед обрушением определялось из неравенства йтах < 0,83Н . Донное трение

вычислялось с учетом эквивалентного масштаба шероховатости кЬс.

Расчет потоков наносов и деформаций дна проводился при рх= 2650 кг/м'; 5 = 1,6; £ = 0,4. В соответствии с результатами натурных наблюдений за транспортом наносов в районе, примыкающем к озеру [19], средний диаметр донных наносов <а?л = 0,3 мм (песок средней зернистости).

3. Анализ численных экспериментов

Озеро Донузлав - сильно вытянутый и сужающийся в северо-восточном направлении бассейн с изрезанной береговой чертой, соединенный с морем узким проливом шириной ~ 600 м. В озере имеется несколько котловин (-20 м) и глубоководный канал вдоль его оси (рис. 1, а). Указанные особенности морфометрии в значительной степени определяют характер ветровой циркуляции и волнения в озере. Сосредоточим внимание на исследовании динамики южной части озера, поскольку именно здесь динамические процессы проявляются наиболее отчетливо. Более подробно рассматриваемая область показана на рис. 1,6.

Выполним моделирование установившейся штормовой ситуации, вызванной постоянным ветром. Как показывает анализ данных ветрового режима (обобщенная информация из банка данных Морского отделения Укр-НИГМИ за период 45 лет любезно предоставлена Л.Н. Репетиным), в рассматриваемом районе продолжительные штормовые ситуации возникают преимущественно в осенний и зимний период. Причем наибольшую повторяемость имеют штормовые ветры северо-восточного (1,19%), восточного (1,23 %) и южного (0,54 %) направлений со скоростями ~ 15 м/с.

Численные эксперименты проводились для указанных направлений ветра при |и10|=15 м/с по следующему сценарию. В начальный момент времени

(/ = 0) вода в озере находится в состоянии покоя. При / > 0 на водную поверхность начинает действовать однородный по пространству ветер, который через 3 ч достигает своего предельного значения и далее остается неизменным.

Как показал анализ предварительных расчетов, режим циркуляции вод в озере можно считать практически установившимся через 24-30 ч после начала действия ветра, а режим волнения - через 8 - 10 ч.

J_i_I_i_I_!_I__I_L

КМ

Рис. 1. Рельеф дна (м) озера Донузлав 52 ISSN 0233-7584. Мор. гидрофиз. журн., 2006, № 2

3.1. Течения и уровень

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

Р и с. 2. Распределение полных потоков (м2/с) для северо-восточного (я), восточного (б) и южного (<?) направлений ветра

по два локальных экстремума (рис.2, а)%. Циклонический круговорот примыкает к северному берегу, а антициклонический - к южному. Между ними располагается узкая интенсивная полоса переноса водных масс, ориентированная вдоль оси озера противоположно ветру. Когда ветер имеет восточное направление (рис.2, б), интегральная циркуляция во многом подобна предыдущему случаю. Главное отличие заключается в том, что северный локальный максимум внутри циклонического круговорота имеет заметно больший размер, а южный максимум в этом круговороте смещен на юг к северной косе. При южном ветре картина существенно отличается от двух предыдущих случаев. Характер циркуляции меняется на противоположный (рис.2, в). У северного берега теперь образуется антицйклонический круговорот, а у южного - циклонический. Оба круговорота содержат по два максимума. Потоки воды между круговоротами направлены теперь на юго-запад. Следует отметить, что во всех трех рассмотренных случаях за пределами озера на взморье наблюдается интегральный перенос водных масс вдоль берега в западном направлении, что определяется результатами моделирования циркуляции для всего моря.

Анализ трехмерной структуры полей течений выявил следующие характерные особенности. Наиболее интенсивные течения со скоростями до 40 -60 см/с на поверхности проявляются вдоль берегов, где глубины не превышают 1 - 3 м (рис.3, а). Направление этих течений не меняется с увеличением глубины и в основном совпадает с генеральным направлением ветра. Лишь с наветренной стороны мысов имеют место существенные отклонения от направления ветра. Вдали от берегов над глубоководными участками акватории отмечается некоторое ослабление поверхностных течений. С глубиной течения перестраиваются и в них начинают проявляться круговороты, качественно повторяющие вихревые структуры,' обнаруженные в интегральной циркуляции (рис. 3, б) (поле получено путем интерполяции с а -координатной сетки на плоскость ^ = -4 м ). При этом вдоль оси канала возникают компенсационные течения, направления которых противоположны направлению ветра. Характерно, что направления течений в придонном слое (<т = -.1) соответствуют направлениям полных потоков (рис. 3, в).

Во всех рассмотренных случаях в проливе, соединяющем озеро с морем, наблюдается двухслойная структура течений. При северо-восточном и восточном ветре до глубины 2 м течения направлены в сторону моря: на горизонте 3 м они слабы и Имеют разнонаправленный характер; начиная с глубины 4 м, они направлены внутрь озера. При южном ветре картина обратная: до 2 м течения направлены в озеро, а ниже 3 м - в сторону моря.

Из-за относительно малых размеров озера и наличия пролива распределение уровня в озере в значительной мере связано с тенденциями уровня всего Черного моря. Как показывают расчеты для всего Черноморского бассейна, ветер восточного и северо-восточного направлений формирует в районе, непосредственно примыкающем к озеру, сгоны величиной 2-5 см. При этой ситуации часть воды из озера вытесняется через пролив в открытое море, что приводит к глобальному понижению уровня на всей или почти на всей его акватории. При северо-восточном ветре максимальные сгоны (~ 30 с.м) возникают в самой северной оконечности озера. Нагоны же образуются лишь в

мелководных участках, непосредственно примыкающих к основаниям кос, которые отделяют озеро от моря. Здесь повышения уровня вблизи берега достигают 10-15 см. При восточном ветре сгонам (до 20 - 25 см) подвержена вся акватория за исключением юго-западной оконечности озера, где уровень повышается до 25 см. Напротив, южный ветер формирует нагон по всей

Р и с. 3. Распределение вектора горизонтальных течений при а = О (а), z = -5 м (б) и сх = -1 (в) для северо-восточного ветра

3.2. Волнение

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

ISSN 0233-7584. Мор. гидрофиз. журн., 2006, № 2 55

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

Р и с. 4. Высота значительных волн (м) и среднее направление волнения для северо-восточного (,а), восточного (б) и южного (в) ветра

На рис.4 показано пространственное распределение высот значительных волн (м) и среднее направление волнения в южной части озера для трех направлений ветра. Изолинии высот даны через 0,1 м. Как видно, при северовосточном ветре поле высот волн имеет вытянутую вдоль оси озера конфигурацию (рис.4, а). Зона максимального волнения (-0,8 м) расположена в центре акватории над глубоководной котловиной. В мелководных участках озера развитие волнения ограничено влиянием диссипации волновой энергии и обру-

шением на критических глубинах. Высоты волн здесь не превышают 0,5-0,6 м. Направления ветра и волнения в основном совпадают. Однако при выходе на мелководье волны за счет влияния рефракции отклоняются на 20-30° относительно генерального направления ветра. При восточном ветре, за счет уменьшения эффективных разгонов, интенсивность волнения на 10% слабее, чем при северо-восточном. Область максимумов высот волн смещена на запад, в сторону наветренного берега (рис.4, б). Отчетливо проявляются эффекты рефракции у мыса в районе с. Озеровка. При южном ветре режим волнения существенно иной (рис.4, в). В этом случае возникает упоминаемый выше дополнительный фактор. Как видно, из открытой части моря в сторону берега приходят штормовые волны высотой - 2,2 м, которые, распространяясь через пролив и трансформируясь в нем, образуют на выходе в озеро локальную область максимумов высот со значениями 1,3-1,5 м, ориентированную в северном направлении. Таким образом, несмотря на то что наиболее эффективные условия для генерации волнения создаются при северо-восточном ветре, при южном ветре высоты волн на входе в озеро могут существенно превышать высоты волн, генерируе-

км

Р и с. 5. Амплитуда скорости придонных волновых течений (см/с) для северо-восточного (а) и южного (б) ветра

На рис, 5 представлены результаты расчетов по формуле (26) поля амплитуд придонных волновых течений (см/с) для северо-восточного и южного ветра. Видно, что при северо-восточном ветре (рис.5, а) выделяются две зоны интенсификации скоростей 30 см/с), которые расположены в южной части бассейна с глубинами 2-3 м. В котловине на входе в озеро и в глубоководном канале скорости не превышают 10 см/с. При южном ветре (рис,5, б) в зоне интенсивного волнения, проникающего из открытого моря, у северного бере га формируется область максимумов со скоростями 40-60 см/с» При этом в проливе волновые течения достигают 40 см/с.

3.3. Транспорт наносов

При штормовых условиях волнение и течения в озере создают в придонном слое касательные напряжения, достаточные для взвешивания и переноса донных осадков, что в конечном итоге может приводить к заметным деформациям дна. Поскольку для стационарных течений и волнения значения П не зависят от времени, результат интегрирования уравнения (28) дает следую» щее выражение для величин деформаций дна: Д# = ПД/0, где Д/0 - время действия устойчивой синоптической ситуации. Функцию П в (28) можно интерпретировать как скорость деформации дна в единицу времени. Значения Л > 0 соответствуют процессам аккумуляция, а Г2 < 0 - процессам эрозии.

На рис. 6 показаны пространственные распределения потоков наносов (кг-м'^с"1) для северо-восточного (рис. 6, а) и южного (рис. 6, б) ветра. Тем» ным цветом выделены зоны интенсификации. Как видно, направление и величина потоков неравномерны по пространству. При северо-восточном ветре проявляются две области усиления потоков со значениями 0,2 - 0,3 кг -м 1 -с л. Первая область расположена у мыса, примыкающего к с. Озеровка. С восточной стороны мыса потоки направлены на юго-запад. После огибания мыса направление потоков наносов меняется на северо-западное. Вторая область находится у южного побережья со стороны п. Мирный. Здесь прослеживается перенос в западном направлении. На глубинах, превышающих 8 м, потоки наносов практически отсутствуют. При южном ветре интенсивность потоков в целом увеличивается примерно в 2 раза. Более отчетливо проявляется зона максимумов (0,7-1,0 кг-м^с"1) у северного берега. При этом генеральное направление движения наносов меняется на противоположное. Характерно, что положение зон интенсификации транспорта наносов совпадает с положением областей максимумов амплитуд придонных волновых течений (рис. 5). При этом направление потоков в целом соответствует направлению придонной циркуляции вод.

Распределения деформаций дна (см) в южной части озера за период 1 сут приведены на рис. 7 для северо-восточного (рис. 7, а) и южного (рис. 7, б) направлений ветра. Видны чередующиеся области эрозии и аккумуляции донного материала. Такое чередование является следствием закона сохранения общей массы наносов, описываемого уравнением (28). За сутки макси-мальные величины эрозии и аккумуляции достигают соответственно следующих значений: для северо-восточного ветра - 4,7 и 3,4 см; для южного ветра-10,7 и 8,1 см.

Р и с. 6. Направление и интенсивность потоков наносов для северо-восточного (а) и южного (б) ветра

Р и с. 7. Распределение деформаций дна (см) для северо-восточного-(а) и южного (б) ветра

Как показывает сравнительный анализ, учет в формуле (39) только одной компоненты, соответствующей донным наносам, приводит к уменьшению максимальных деформаций дна примерно в 5 раз по сравнению с вариантом учета обеих составляющих потока наносов. Отсюда следует вывод, что для рассмотренных случаев определяющий вклад в формирование донных деформаций в озере вносит циркуляция вод, а основная роль волн сводится к взвешиванию донных осадков.

3.4. Влияние волнения на течения и уровень

Представляет интерес оценить количественные различия между результатами объединенного и раздельного моделирования полей циркуляции и волнения. Поскольку этот вопрос требует отдельного рассмотрения, в рамках данной работы ограничимся оценками, касающимися влияния волнения на поля течений и уровень моря. Были выполнены численные эксперименты по моделированию полей течений и уровня без учета волнения, но с той же параметризацией турбулентной вязкости. Параметр шероховатости донной поверхности в (22) задавался выражением (23), а коэффициент поверхностного трения определялся по формуле (19). Как оказалось, его величина по акватории озера в среднем на 45% больше средней величины коэффициента поверхностного трения, рассчитанного с учетом волнения по формулам (20), (21). Аналогичное сопоставление для коэффициента придонного трения показало, что средняя величина этого коэффициента при учете взаимодействия волн и течений в придонном пограничном слое примерно в 2 раза больше, чем без его учета. Таким образом, объединенная модель дает, с одной стороны, больший приток энергии от ветра, но, с другой, - большую диссипацию энергии в придонном пограничном слое. Кроме того, объединенная модель дает более интенсивное поле течений. Например, при северо-восточном ветре максимумы скоростей в приповерхностном слое увелйчиваются на 12%, при южном ветре - на 40%.

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

Р и с. 8. Распределение функции тока (м7с) для южного ветра: а - без учета волнения; б - с учетом волнения; в- аномалии функции тока

Следующим интересным эффектом, который воспроизводит объединенная модель, является волновой подъем {set-up) и волновое опускание {set-down) среднего уровня моря, которые обусловлены наличием дополнительного потока импульса в жидкости, создаваемого пространственной неоднородностью поля волновых радиационных напряжений (18). Следует отметить, что эти явления имеют теоретическое обоснование и также проявляются в результатах лабораторных и натурных исследований [11,20]. Для выделения этого эффекта удобно рассматривать величины относительных аномалий уровня 5г\ = 100(77,. - T]n)/rjc, где ?]с, г/п - уровни, рассчитанные с учетом и без учета волнения. Пример распределения Sr/ для южного ветра показан на рис. 9. Как видно, с внешней стороны кос образуются узкие полосы отрицательных аномалий уровня (понижение на 10-13 %). Это есть проявление волнового опускания уровня при подходе волн по нормали к берегу в его классическом виде. В самом озере картина становится более сложной. Здесь также проявляются области отрицательных аномалий. Причем их положение хорошо коррелируется с положением зон интенсификации ветрового волнения. В пределах зоны наиболее интенсивного волнения прослеживается локальная область максимума отрицательных аномалий. У основания южной косы формируется область положительных аномалий (10-40%), соответствующая волновому подъему уровня. В целом такое пространственное распределение аномалий подтверждает известный теоретический результат о том, что интенсивные волны при выходе на мелководье понижают средний уровень моря.

3 4 5 6 7 8 9 10 11 12

КМ

Р и с. 9. Относительный вклад (%) волнения в изменения среднего уровня моря при южном ветре

Выводы

Разработана объединенная численная модель течений, волнения и транспорта наносов озера Донузлав. Модель основана на а -координатном гидродинамическом блоке, спектральной волновой модели и полуэмпирических

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

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

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

Пространственное распределение потоков наносов характеризуется наличием локальных максимумов, положение которых совпадает с положением максимумов амплитуд придонных волновых течений. Направление перемещения наносов в целом соответствует направлению придонной циркуляции вод. Определяющий вклад в формирование донных деформаций вносит циркуляция, а основная роль волн сводится к взвешиванию донного материала. При скорости ветра 15 м/с максимумы донных деформаций, рассчитанные на период 1 сут, для песка средней зернистости достигают ± 9 см.

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

СПИСОК ЛИТЕРАТУРЫ

1. Андросович А.И., Иванов В.А., Михайлова Э.Н., Шапиро Н.Б. Моделирование ветровых течений в озере Донузлав // Морской гидрофизический журнал - 1996.- № 2 - С. 15-26.

2. Mastenbroek С, Burgers G., Janssen P. The dynamical coupling of a wave model and a storm surge model through the atmospheric boundary layer // J. Phys. Oceanogr- 1993 - 23 -P. 1856-1866.

3. Davies A.M., Lawrence J. Examining of influence of wind and wind wave turbulence on tidal currents, using a tree-dimensional hydrodynamic model including wave-current interaction // Ibid.- 1994.-24.- P. 2441-2459.

4. Xie L., Wu K., Pietrafesa L., Zhang C. A numerical study of wave-current interaction through surface and bottom stress: wind-driven circulation in the South Atlantic Bight under uniform winds//J. Geophys. Res.-2001.-C8.-P. 16841-16855.

5. Фомин B.B. Численная модель циркуляции вод Азовского моря. - Научные труды Укр-НИГМИ.- 2002.- Вып. 249.- С. 246-255.

6. Фомин В.В., Алексеев Д.В., Иванча Е.В., Черкесов Л.В. Исследование гидродинамических полей в северо-западной части Черного моря, возникающих при прохождении циклона // Доклады Национальной академии наук Украины - 2004 - № 5 - С. 137-142.

64

ISSN 0233-7584. Мор. гидрофиз. журн., 2006, № 2

7. Booij N., Ris R., Holthuijsen L A third-generation wave model for coastal regions. Model description and validation // J. Geophys. Res.- 1999.- ¡04 (C4)~ P. 7649-7666.

8. Mellor G.L The three-dimensional current and wave equations // J. Phys. Oceanography-

2003,™ 23, m 9.- P. 1978-Î989,,

9. Mellor G.L., Yamada T. Development of a turbulence closure model for geophysical fluid problems // Rev. Geophys. Space Phys,- 1982.- 20.- P. 851-875.

10.. Smagorinsky 1 General circulation experiments with primitive equations. I. The basic experiment // Mon. Weath. Rev.- 1963.- 91,- P. 99-164.

11. Филлипс О. M. Динамика верхнего слоя океана. - JL: Гидрометеоиздат, 1980 - 319 с.

Î2 Sweby Р. К. High resolution schemes using flux limitera for hyperbolic conservation laws // J. Numcr. Ail - 1984.- 21.- P. 995 1011.

13. Large W.G., Pond S. Open ocean momentum fluxes in moderate to strong winds // i. Phys. Oceariogr.- 1981.~Ц.~P. 324-336.

14. Китайгородский СЛ. Физика взаимодействия атмосферы и океана - Л.: Гидрометеоиздат,

1970,,- 284 с.

15. Taylor Р.К., УеHand M J. The dependence of sea surface roughness on the height and steepness of the waves //1. Phys. Oceanog.- 2001.- Ц (2).- P. 572-590.

16. Grant W.D., Madsen O.S. Combined wave and current interaction with a rough bottom // J. Geophys. Res.™ 1979.-84 - P 1797-1808.

17. Блатов A.C, Иванов В Л, Гидрология и гидродинамика шельфовой зоны Черного моря -Киев: Наук, думка, 1992.- 242 с.

18. Van Rijn LC Sediment transport; part Î: bed load transport // J. Hydraul. Eng.- 1984.- 10.-P, 1431-1456.

19. Анциферов C.M., Косьян P.Д. Взвешенные наносы в верхней части шельфа.- М.: Наука,

1986.-224 с.

20. Ле Блон П., МайсекЛ. Волны в океане. Т.2.- М.: Мир, 1980.- 363 с.

Морское отделение Украинского научно- исследовательского Материал поступил

гидрометеорологического института, в редакцию 21.12.04

Севастополь После доработки 01.02.05 Морской гидрофизический институт HAH Украины, Севастополь

ABSTRACT Numerical model of the Donuzlav lake dynamics permitting to calculate simultaneously currents, sea level, waves and deposit transport is submitted. The model is based on a hydrodynamical block and a spectral wave model. For the characteristic storm situations, the features of the integrated water circulation and three-dimensional structure of currents, and also wind waves fields are investigated; the estimates of values of deposit fluxes and bottom deformations are obtained. Presence of intensive vertical structures in the current field induced by the bottom topography is revealed. Essential wave intensification in the southern part of the lake occurring during penetration of storm waves through the passage is found. It is shown that consideration of the waves in the modeling results in qualitative changes,of the circulation structure in the lake and also in formation of the well-pronounced areas of wave lifting and sea level lowering.