УДК 551.465
Н.Б. Шапиро
Численная модель субмаринной разгрузки пресных вод в карстовой полости
Обсуждаются результаты численных экспериментов по воспроизведению структуры гидрологических полей в карстовой полости (район м. Айя, Южный берег Крыма), наблюдаемой во время второй специализированной экспедиции в сентябре 2008 г. Для описания процесса субмаринной разгрузки пресных вод используется трехмерная гидростатическая модель с заданными априори границами слоев (а- и/или 2-уровни) в приближении «твердой крышки» и с процедурой «конвективного приспособления». Установлено, что только хорошее разрешение по вертикали позволяет количественно правильно описать тонкий верхний распресненный слой моря. Показано, что для адекватного представления о влиянии рельефа дна на верхний слой моря желательно использовать смешанную сетку по вертикали, сочетающую а- и 2-координаты. Представлены результаты численных экспериментов, позволяющие судить о влиянии местоположения, количества и мощности источников пресной воды, расположенных на боковой границе полости.
Ключевые слова: численное моделирование, гидростатическая модель, субмаринная разгрузка пресных вод, карстовая полость.
Данная работа представляет собой продолжение исследования по численному моделированию процесса субмаринной разгрузки в карстовой полости, расположенной в районе м. Айя на Южном берегу Крыма. Результаты моделирования представлены в [1, 2], причем в [2] приводятся и обсуждаются данные наблюдений, полученные во время двух специальных экспедиций. В [1] была предпринята попытка воспроизвести ситуацию, наблюдаемую во время первой экспедиции в сентябре 2007 г.
В настоящей работе представлены результаты численных экспериментов по воспроизведению структуры гидрологических полей в той же карстовой полости во время второй экспедиции, которая проводилась в сентябре 2008 г., через год после первой. Так же четко в приповерхностном слое наблюдались распресненные воды - с минимальной соленостью 10%о, что на две единицы меньше, чем в 2007 г. В глубинных слоях соленость, как и прежде, составляла 18,2%о. В верхнем слое моря температура воды была 22,3 - 23,3°С, с глубиной она увеличивалась и достигала 24,3°С. В целом температура воды была значительно выше, чем в 2007 г., когда из-за сгона вод она менялась от 14,7 -15°С на поверхности моря до 14,2°С на дне.
Структура течений в полости качественно была такой же, как во время первой экспедиции, а именно: в тонком верхнем слое толщиной ~ 0,2 м вода вытекала из полости (со скоростью ~ 4 см/с), глубже она, наоборот, втекала в нее.
По данным о скорости течения и солености на разрезе - входе в полость - на основе стационарного баланса массы и соли был оценен суммар-
© Н.Б. Шапиро, 2011
66
0233-7584. Мор. гидрофиз. журн., 2011, № 5
ный дебит Q0 источников пресной воды, разгружающихся в полости. Он получился равным 6000 м3/сут » 0,06 м3/с, т. е. примерно в два раза больше, чем в 2007 г., что, вероятно, связано с большей величиной осадков, выпавших над Крымом в 2008 г.
Для моделирования ситуации в карстовой полости, как и в предыдущей работе [1], решается следующая задача. Пусть в начальный момент времени в области, включающей карстовую полость, имеет место покой и задаются меняющиеся только по вертикали поля температуры и солености. Затем через определенные трещины в боковых границах карстовой полости начинает просачиваться пресная вода с нулевой соленостью, температурой 15°С (в данном случае более низкой, чем температура воды в море) и заданными априори расходами. Пресная вода вследствие конвективной неустойчивости практически мгновенно всплывает к поверхности моря, перемешиваясь при этом с соленой морской водой, и в виде распресненного слоя выносится в открытое море. Из моря глубинными течениями в полость вносится соленая морская вода. В результате в полости формируется трехмерная (по существу двухслойная) структура полей течений, температуры и солености.
Численные эксперименты проводились в рамках той же модели, что и в работе [1], детали численной схемы подробно описаны в [3]. Это N-слойная гидростатическая модель с заданными априори границами слоев с использованием приближения «твердой крышки» и процедуры «конвективного приспособления» с сохранением запасов тепла и соли при появлении неустойчивой стратификации по плотности (подъем пресных вод к поверхности моря явно не описывается). Положение границ слоев Zk, k = 0, 1, ..., N, в работе [1] задается формулой
Zk(x, y) = ак H(x, y), (1)
где 0 £ sk £ 1; s0= 0; sN = 1; Z0 = 0 - невозмущенная поверхность моря; ZN = H - дно моря.
Уравнения модели получаются интегрированием по вертикали в пределах слоев исходной системы уравнений с учетом граничных условий на поверхности моря, дне и внутренних границах слоев. При задании Zk по формуле (1) получаем такие же уравнения, как и при аппроксимации дифференциальных уравнений в о-координатах (x, y, о = z/H), что позволяет в данном случае считать используемую модель дискретной версией модели в о-коор-динатах.
Благодаря приближению «твердой крышки» интегральное уравнение неразрывности имеет вид
Ux + Vy = 0, (2)
где U, V- компоненты полного потока вдоль осей X, Yсоответственно.
Уравнение (2) позволяет ввести интегральную функцию тока Y согласно формулам
U = - Yy, V = Yx. (3)
На поверхности моря потоки массы и соли полагаются равными нулю. Действие ветра не учитывается (в полость ветер практически не проникает,
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5 67
а над рассматриваемой открытой частью моря его действием пренебрегаем). Поток тепла на поверхности моря пропорционален разности температур воды T\ и воздуха Т0:
ОГ = 1 (Т - Т\). (4)
На дне принимаются условия обтекания с трением и отсутствия потоков тепла и соли. На внутренних границах слоев используются условия непрерывности потоков импульса, тепла и соли, а также гипотезы замыкания для этих потоков, обладающие свойством транспортивности [3].
На твердых боковых границах ставятся условия прилипания для горизонтальных компонент скорости течения и условия равенства нулю нормальных составляющих потоков тепла и соли. В трещинах на боковых границах полости, а именно в местах, помеченных на рис. \ точками А, В, С, В, в слоях к\ £ к £ к2 предполагается поступление (по нормали к границе) пресных вод с соленостью 5"* = 0%о, температурой Т = \5°С и известными расходами ОА, ОВ, ОС, Ов. При этом предполагается, что в пределах слоев [к\, к2] скорость течения с глубиной не меняется.
А
У, м: 45т 40т 35т 30т 25т 20т
15т
10т
5-
5 10
15 ' 20 X, м
10 ' '15 ' 20 X, м
а
б
Р и с. 1. Рельеф дна (Н, м) - а, положение станций и разрезов, которые используются при анализе результатов расчета, - б
68
0233-7584. Мор. гидрофиз. журн., 2011, № 5
5
На границе полости задаются расходы воды, т. е. распределение полного потока и, фактически, интегральная функция тока. Поскольку используется приближение «твердой крышки», интегральный расход воды через открытую границу известен, в каждый момент времени он равен суммарному расходу воды, поступающей из подземных источников. Для того чтобы определить распределение этого расхода вдоль открытой границы, используется условие свободного протекания для полного потока. Это условие следует из предположения, что вода течет по нормали к границе. Например, на «южной» границе при Y = 0
U = 0, Э V/ Э y = 0 (5)
или для интегральной функции тока
Э Y /Э n = 0, (6)
где Э /Э n - производная по нормали к границе.
Для скорости течения на открытой границе также ставится условие свободного протекания, но для бароклинных компонент скорости (отклонений от баротропных, средних по глубине, компонент), например, на «южной» границе при Y = 0, оно имеет вид
щ = 0, Э (Vk - V/H)/Э y = 0, (7)
где uk, vk - компоненты скорости течения в k-м слое вдоль осей X, Y.
Для температуры и солености на открытой границе ставятся условия транспортивности для суммарных (адвективных плюс диффузионных) нормальных потоков тепла и соли, например, на границе Y = 0
- Khk(Tk)y + Vktk = v+ (tk)l + v- tk,
- khk(Sk)y + Vksk = vk (sk)L + v- sk, (8)
где Vk = vkhk - компонента потока воды в k-м слое, причем Vk ° vk + vk , vk = = max (0, Vk) > 0 - втекающий, vk = min(0, Vk) £ 0 - вытекающий из области поток; Tk, Sk - температура и соленость воды в k-м слое. Видно, что в k-м слое в рассматриваемую область втекает вода с фоновой температурой и соленостью (Tk)L, (Sk)L, а вытекает - со своей температурой и соленостью.
Численные эксперименты. Как показано в работе [1], модель качественно правильно описывает трехмерную структуру формирующихся в карстовой полости полей течений, температуры и солености. Однако имело место количественное несоответствие между рассчитанным и наблюдаемым распределениями солености на поверхности моря. Поверхностная соленость, как известно, является наиболее показательной характеристикой процесса субмаринной разгрузки подземных вод. Рассчитанная в [1] соленость оказалась выше наблюдаемой как внутри, так и вне карстовой полости. Кроме того, представлялась неестественной рассчитанная структура течений у открытой границы области, где образовывались круговороты, проявляющиеся
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
69
в интегральной циркуляции и в глубинных течениях. Для подавления этих круговоротов у границы специально вводилась «буферная зона» с завышенными значениями коэффициента рэлеевского (пропорционального скорости течения) трения. Как показали специальные численные эксперименты, если вместо условия свободного протекания на границе просто задавать распределение полного потока (равномерное вдоль границы), то круговороты у открытой границы вообще не образуются.
Это показывает, что образование данных круговоротов является артефактом, оно связано с неправильной аппроксимацией условия свободного протекания (5), (6). Фактически условие аппроксимировалось с первым порядком точности, хотя, как и уравнения, должно аппроксимироваться со вторым порядком точности. Далее все численные эксперименты проводились с более точной аппроксимацией граничного условия (детали аппроксимации приводятся в приложении). Отметим, что уточнение аппроксимации позволило избавиться от нереальных особенностей в поле течений и необходимости вводить «буферную зону».
Итак, перейдем к моделированию ситуации в карстовой полости, имеющей место в период второй экспедиции. На рис. 1 приведен рельеф дна в рассматриваемой области, показано положение источников пресной воды, станций и разрезов, которые используются при демонстрации результатов расчета. Отметим, что дебит субмаринной разгрузки считается известным, однако количество источников пресной воды, их положение, величина расходов поступающей воды неизвестно.
Вначале был проведен расчет при тех же значениях параметров, что и в работе [1], а именно в рамках 8-слойной модели на равномерной по ö-коор-динате сетке (ok = k/N) с коэффициентами горизонтальной (A=102 см2/с, К = 0) и вертикальной (Az = 50 см2/с, Kz = 0,1 см2/с) вязкости и диффузии соответственно, l = 0. Шаги сетки по горизонтали А х = А y = 0,5 м. Суммарный дебит подземных источников Q0 полагался равным 6000 м3/сут (или »0,06 м3/с), расходы источников QA = 500 м3/сут, QB = 500 м3/сут, QC = = 2500 м3/сут, Qd = 2500 м3/сут. Предполагалось, что подземные воды втекают в три внутренних слоя (k1 = 4, k2 = 6) с температурой 15°С и соленостью 0%о. Как уже указывалось, в начальный момент времени скорость течений равна нулю, температура и соленость в полости и на открытой границе меняются только по вертикали, причем соленость увеличивается с глубиной по линейному закону от 17,9 до 18,2%о, а температура, учитывая условия второй экспедиции, наоборот, уменьшается от 24 до 23°С. На открытой границе фоновое распределение температуры и солености сохраняется в течение всего расчета. Расчет (с шагом по времени А t = 0,72 с) проводился на 10 ч, когда все поля выходили на стационарный режим.
На рис. 2 приведены результаты, полученные к концу расчета. Показана интегральная функция тока, соленость и температура в поверхностном слое, векторы средней по глубине скорости течения и скорости течений в поверхностном и придонном слоях. Отметим, что структура полей течений и солености близка к структуре полей, полученных при моделировании результатов первой экспедиции (с уточненной аппроксимацией граничного условия). Поля температуры, естественно, отличаются, поскольку первая экспедиция про-70 ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
водилась во время сгона, когда температура воды в море была значительно ниже температуры втекающей через трещины пресной воды. Видно, что у открытой границы области нет неестественных особенностей в поле течений, которые описаны в [1]. В то же время, как и в [1], четко проявляются особенности, причем во всех гидрологических полях, связанные с поднятием дна около входа в карстовую полость.
Р и с. 2. Интегральная функция тока ¥ м3/с, соленость S%o, температура T°C в поверхностном слое (вверху), векторы скорости среднего по глубине течения в поверхностном (k = 1) и придонном (k = 8) слоях (внизу). Приведены максимальные значения скорости в см/с; N = 8
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5 71
Подчеркнем, что при выбранных значениях параметров, как и при воспроизведении результатов первой экспедиции, поверхностная соленость получается значительно выше наблюдаемой. Такое несоответствие может быть связано как с особенностями модели, так и с неточным заданием внешних факторов, например местоположения источников пресной воды и их расходов.
Как показали численные эксперименты, результаты расчетов существенно зависят от вертикального разрешения модели - во-первых, от количества слоев и, во-вторых, от вида вертикального разрешения, а именно использования о- или Z-координаты. Кроме того, оказалось, что поверхностная соленость в полости существенно зависит от интенсивности источников пресной воды Qa, Qb, расположенных в вершинной части полости. Задание в [1] относительно больших расходов воды у источников QC, QD, расположенных около входа в полость, основывалось на субъективных представлениях участников экспедиции.
В связи с этим рассмотрим и сравним результаты численных экспериментов, проведенных в рамках 8- и 40-слойной моделей с различным разрешением по вертикали: с равномерной по о-координате сеткой (сетка о) и в случае, когда границы слоев горизонтальные и все слои, кроме придонного, имеют одинаковую толщину (сетка Z). В 8-слойной модели толщина 7 верхних слоев равна 40 см, в 40-слойной модели толщина 39 верхних слоев — 7 см, т. е. все границы слоев располагаются над поднятием дна у входа в полость. Будем считать, что пресная вода поступает только из одного источника, расположенного в вершине полости (QA = 0,06м3/с, QB = QC = QD = 0), причем во все слои (k1 = 1, k2 = N). Последнее предположение не является принципиальным, а сделано только для адекватного сопоставления результатов, полученных при различном числе слоев. Подчеркнем, что в модели [3] предусмотрено использование любого разрешения по вертикали путем произвольного задания топографии границ слоев Zk(x, y).
На рис. 3 приведены вертикальные профили солености и температуры на нескольких станциях, расположенных внутри (ст. 1, 4, 10) и вне (ст. 7) полости. Видно, что увеличение числа слоев (от 8 до 40) приводит к существенному уменьшению поверхностной солености и толщины распресненного слоя, что приближает результаты расчета к данным наблюдений. Отметим, что этому способствовало также (даже в 8-слойной модели) задание достаточно мощного источника в вершинной части полости (с расходом, равным суммарному дебиту источников в предыдущем варианте). Видно, что в 40-слойной модели вертикальная структура профилей солености мало чувствительна к виду вертикального разрешения и близка к структуре профилей, полученной в 8-слойной модели с сеткой Z.
Перечисленные выводы можно отнести также к структуре вертикальных профилей температуры. Подчеркнем, что теперь, в отличие от первой экспедиции, проводившейся во время сгона вод [1], в поверхностном слое находятся воды с пониженной температурой. Это естественно, так как пресные воды по сравнению с морскими имеют более низкую температуру (15°C). Поднимаясь к поверхности, они приводят не только к уменьшению солености, но и к охлаждению поверхностных вод.
72 ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
8 11 14 17 8
9 12 15 18 8
л
и £ 2
ел
1-4
3
4
0 1
л к „
£ 2 Е?
1-4
з Н
Ст. 1
18 20 22 24 Т
4
0 1 2
3 -
4 -
5
6
7 ^
Ст. 4
19 21 23 Т
1 -2 -
3 -
4 -
5 -
6 -7 -
9 12 15 18 8 10 13 16 8 0 1
2 -
3
4
5
6
7
8 -9 :
2 -
3 -
4 -
0 1 2
3 -
4 -
Ст. 10
19 21 23 Т
19 21 23 Т
0 1 : 2 :
3 :
4 :
5 :
6 :
7 :
8 : 9 :
Р и с. 3. Вертикальные профили солености и температуры на указанных станциях в случае действия только одного источника пресной воды, расположенного в вершинной части полости (вл = 0,06 м3/с, вв = вс = QD = 0); сплошные кривые - при N = 40, штриховые - при N = 8, жирные - на сетке а, тонкие - на сетке 2
На рис. 4 для тех же вариантов N = 40, N = 8, сетки а и 2) показаны распределения солености и температуры в верхнем слое моря (к = 1). Увеличение числа слоев приводит к уменьшению солености непосредственно в полости и в открытой части моря, причем на обеих сетках. Видно, что распределения солености, полученные на разных сетках, значительно отличаются друг от друга. На сетке а в поле солености в верхнем слое проявляется влияние поднятия дна, особенно при грубом разрешении N = 8). В то же время на сетке 2 такого влияния практически нет, по крайней мере при хорошем вертикальном разрешении N = 40). Это означает, что на самом деле рельеф дна практически не влияет на распределение поверхностной солености. Влияние, проявляющееся при использовании сетки а, связано с переменной глубиной слоя. Поэтому к интерпретации результатов, рассчитанных при использовании сетки а, нужно относиться с осторожностью. Отметим, что указанные особенности присущи и распределению поверхностной температуры.
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
73
числе слоев (N = 40, N = 8) и различных сетках
На рис. 5 приведены распределения интегральной функции тока Y и векторы средней по глубине скорости течения (U/H, V/H), рассчитанные для тех же вариантов, что и распределения солености и температуры. Сразу же отметим, что в интегральной циркуляции и в поле средних по глубине течений хорошо виден так называемый СЭБИР — совместный эффект бароклин-ности и рельефа дна, проявляющийся в появлении круговоротов над поднятием дна. Видно, что при увеличении числа слоев и использовании сетки Z этот эффект несколько сглаживается.
74 ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
N=40 Сетка О N=8 N=40 Сетка Ъ N=8
и/н
Р и с. 5. Распределения интегральной функции тока ¥ м3/с и векторы средней по глубине скорости течения при различном числе слоев (^ = 40, N = 8) и различных сетках. Приведены максимальные значения скорости
На рис. 6 приведены распределения солености и температуры во внутренних слоях, рассчитанные в рамках 40-слойной модели с сеткой а. Видно, что на фоне монотонного изменения солености и температуры вдоль полости над поднятием дна четко выделяется зона пониженной солености и пониженной температуры. В то же время при расчете на сетке как видно из рис. 7, в распределениях полей солености и температуры в приведенных приповерхностных слоях поднятие дна практически не проявляется. Кроме того,
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
75
не приводя рисунков, отметим, что в более глубоких слоях соленость и температура, исключая небольшую окрестность около источника пресной воды, практически не меняются по горизонтали. Это указывает на то, что рельеф дна не влияет на поля температуры и солености не только на поверхности моря, но, по существу, вообще в верхней его части, по крайней мере над поднятием дна. Особенности полей солености и температуры над поднятием дна, полученные в модели с ^-координатами, обусловлены, как уже указывалось, не влиянием рельефа дна, а переменностью глубины границы слоев.
S(k=2, Z=0.05H) S(k=4, Z=0.1H) S(k=8, Z=0.2H) S(k=16, Z=0.4H)
T(k=2, Z=0.05H) T(k=4, Z=0.1H) T(k=8, Z=0.2H) T(k=16, Z=0.4H)
Р и с. 6. Распределения S%o (вверху) и T°C (внизу) во внутренних слоях, рассчитанные в модели с сеткой a (N = 40)
76 ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
S(k=2, г=14см) S(k=5, г=35см)
S(k=10, г=70см) S(k=20, г=1.4м)
м 4540 353025 20 15 10 5
5 10 15 20 м 5 10 15 20 м 5 10 15 20 м 5 10 15 20 м
T(k=2, г=14см) T(k=5, г=35см) T(k=10, г=70см) T(k=20, г=1.4м)
м 45 40 35 30 25 20 15 10 5
5 10 15 20 м 5 10 15 20 м 5 10 15 20 м 5 10 15 20 м
Р и с. 7. Распределения 5%о (вверху) и Т°С (внизу) во внутренних слоях, рассчитанные в модели с сеткой 2 N = 40)
Заметим, что интерполяция полей солености и температуры, рассчитанных на сетке о, на постоянные горизонты не спасает положения. Это хорошо видно на рис. 8, где приведены распределения солености на горизонтах, полученные с учетом, что соленость не меняется по вертикали в пределах слоев (ступенчатая интерполяция) или меняется линейно от слоя к слою (линейная интерполяция). В последнем случае считается, что соленость определяется в середине слоев. В интерполированных полях опять проявляется влияние поднятия дна, только теперь, что естественно, это зона повышенной солености. ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5 77
Поскольку интерполированные с сетки а поля значительно отличаются от полей на сетке для получения полей именно на горизонтах необходимо, по-видимому, с самого начала проводить расчеты в модели с 2-координатами.
Ступенчатая интерполяция
8(к=2, Ъ=14см) 8(к=5, Ъ=35см) 8(к=10, Ъ=70см) 8(к=20, Ъ=1.4м)
м
45
40
35
30Н
25Н
20т
15т
10т 5т
5 10 15 20 м
5 10 15 20 м 5 10 15 20 м линейная интерполяция
5 10 15 20 м
Р и с. 8. Распределения S%o на отдельных горизонтах, рассчитанные в модели с сеткой а (Ы = 40), после ступенчатой (вверху) и линейной (внизу) интерполяции
78
ISSN 0233-7584. Мор. гидрофиз. журн, 2011, № 5
На рис. 9 показаны векторы скорости течений в различных слоях, полученные в модели с сеткой о, на рис. 10 — то же в модели с сеткой 2. Отметим, что и в первом, и во втором случаях видны особенности в поле течений, связанные с влиянием рельефа дна, которые, однако, проявляются различным образом в зависимости от вида сетки. В целом в верхних слоях вода вытекает из полости в открытое море, а в глубинных слоях, наоборот, втекает.
и(к=1, г=0.025И) и(к=4, г=0.1Н)
и(к=8, г=0.2И) и(к=16, г=0.4И)
4540353025 20 15 10 5
15 20 м
10 15 20 м
10 15 20 м 5 10 15 20 м
и(к=20, г=0.5И)
и(к=30, г=0.75И) и(к=35, г=0.875И) и(к=40, г=И)
м : 45т 40т 35т 30т 25т
20т
15т 10т 5т
5 10 15 20 м 5 10 15 20 м 5 10 15 20 м 5 10 15 20 м
Р и с. 9. Векторы скорости течения в различных слоях, рассчитанные в модели с сеткой о (V = 40). Приведены максимальные значения скорости в см/с
м
0233-7584. Мор. гидрофиз. журн., 2011, № 5
79
и(к=1, Ъ=7см) и(к=5, Ъ=35см) и(к=10, Ъ=70см) и(к=20, Ъ=1.4м)
и(к=25, Ъ=1.75м) и(к=30, Ъ=2.1м) и(к=35, Ъ=2.45м) и(к=40, 2.73м</<Н)
Р и с. 10. Векторы скорости течения в различных слоях, рассчитанные в модели с сеткой 2 N = 40). Приведены максимальные значения скорости в см/с
Отметим, что при выбранных значениях параметров расход вытекающей из карстовой полости воды получается равным = -0,24 м3/с, а расход втекающей в грот морской воды в2 = 0,18 м3/с, причем и в2 заметно превышают расход поступающей подземной пресной воды в0 = -0,06 м3/с.
На рис. 11 показаны распределения солености, температуры и компоненты скорости течения V, направленной вдоль полости, на двух разрезах -
80
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
продольном (вдоль полости) и поперечном (вход в полость), полученные в модели с сеткой Z. Хорошо видна во всей полости практически двухслойная структура всех полей с четко выраженными скачками солености и температуры, расположенными на достаточно малой глубине. В модели с сеткой о получается практически такая же картина. Отметим, что по сравнению с расчетами с более грубым вертикальным разрешением (Ы = 8) полученная термохалинная структура значительно лучше соответствует данным наблюдений.
а б
Р и с. 11. Распределения солености, температуры и компоненты скорости течения V, см/с, на разрезах вдоль полости (а) и на входе в нее (б). Штриховыми линиями показаны изохалина 12%о, изотерма 21°С и изотахи V > 0
Для того чтобы представить себе влияние местоположения источников и их мощности, на рис. 12 показаны распределения солености в поверхностном слое, полученные в рамках 40-слойной модели с сеткой Z. Видно, что один источник пресной воды, расположенный вблизи входа в грот, не приводит к появлению распресненной воды внутри полости. В то же время наличие одного источника в вершинной части полости приводит к существенному понижению поверхностной солености и в полости, и в прилегающей области моря. Более того, значения солености получаются близкими к измеренным при меньших (по сравнению с основным случаем, когда QA = 0,06 м3/с) расходах этого источника (при QA = 0,03 м3/с и QA = 0,02 м3/с). Наличие несколь-
0233-7584. Мор. гидрофиз. журн., 2011, № 5
81
ких источников пресной воды, как можно судить из рис. 12, влияет на распределение поверхностной солености (особенно в море вне полости) и, возможно, приближает результаты к реальности. Однако чтобы делать определенные выводы о местоположении и мощности источников пресной воды, а также о суммарном дебите источников пресной воды, необходимы дополнительные расчеты и самое главное - более детальные натурные измерения.
8(к=1)
м 45403530252015 10 5
5 10 15 20 м 5 10 15 20 м а б
5 10 15 20 м 5 10 в
15 20 м г
45403530252015105
5 10 15 20 м 5 10 15 20 м 5 10 15 20 м
д е ж
5 10 15 20 м
з
м -
Р и с. 12. Распределения S% в поверхностном слое моря (к = 1), рассчитанные в модели с сеткой 2 N = 40), при различной интенсивности источников пресной воды: QA = Qв = Qc = 0, & = 0,06 м3/с - а; QA = Qв = Qc = QD = 0,015 м3/с - б; QA = QD = 0,03 м3/с, Qв = Qc = = 0 - в; QA = Qв = 0,005 м3/с, Qc = QD = 0,025 м3/с - г; QA = 0,03 м3/с, Qв = Qc = QD = = 0 - д; QA = 0,02 м3/с, Qв = Qc = QD = 0 - е; QA = QD = 0,01 м3/с, Qв = Qc = 0 - ж; QA = 0,02 м3/с, Qв = Qc = 0, QD = 0,01 м3/с - з
82
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5
В заключение сделаем вывод, что модель, несмотря на приближение гидростатики, может качественно и даже количественно правильно описывать трехмерную структуру полей температуры, солености и скорости течений в карстовой полости, формирующихся в результате субмаринной разгрузки подземных вод. При этом необходимо использовать достаточно хорошее разрешение по вертикали, позволяющее правильно описать структуру тонкого верхнего распресненного слоя воды. Кроме того, для адекватного представления о пространственной структуре термохалинных полей и о влиянии рельефа дна на распределение температуры и солености желательно использовать смешанную сетку, сочетающую о- и Z-координаты.
Приложение. Рассмотрим аппроксимацию условия (5), (6) на «южной» границе У = 0, а именно на линии у = -1/2, учитывая, что в численной модели используется сетка В (по терминологии Аракавы). На этой сетке интегральная функция тока у определяется в центрах боксов, а компоненты полного
потока и горизонтальной скорости течения, например Ц+1/2, у+1/2, Ума, у+1/2, и+1/2, у+1/2, к, Уг+1/2, у+1/2, к, - в узлах боксов. Граница области проходит через узлы боксов.
Условие Э¥ / Эу = 0 при У = 0 в работе [1] расписывалось на двух уровнях:
¥ = ¥
г',-1 ',1 '
а интегральная функция тока рассчитывалась (методом верхней релаксации) во внутренних точках области с номерами у > 1.
Учитывая, что компоненты полного потока связаны с функцией тока формулами, имеющими второй порядок точности,
иг'+1/2, у+1/2 = +1, у ¥г, у+1 ¥г'+1, у+1 )/(2Ду),
^'+1/2, у+1/2 = ( ¥г'+1, у+1 + ¥
+1, у ', у+1
] )/(2Дх),
граничное условие (5) необходимо расписывать на трех уровнях, а именно:
¥ = ¥ ¥ = ¥
',1 ',2 ',-1 ',1
Теперь интегральная функция тока рассчитывается во внутренних точках области с номерами у > 2. В результате получаем аппроксимацию условий при У = 0:
и = 0, Э V/ Э у = 0
со вторым порядком точности, которые согласуются с интегральным уравнением неразрывности (2). Кстати, условия для полного потока на твердой границе также расписываются на трех уровнях.
Условия для бароклинных компонент скорости течения (7), как и в работе [1], записываются в виде
Мг+1/2, -1/2, к = 0, У;+1/2, -1/2, к — ^'+1/2, -1/2/Н+1/2, -1/2 = ^'+1/2, 1/2, к — ^'+1/2, 1/2/Н+1/2, 1/2.
Граничные условия (8) для температуры и солености, которые определяются в центрах боксов, аппроксимируются так же, как и в работе [1].
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5 83
СПИСОК ЛИТЕРАТУРЫ
1. Шапиро Н.Б. Моделирование трехмерной структуры гидрологических полей в карстовой полости под действием субмаринной разгрузки подземных вод // Морской гидрофизический журнал.- 2010. - № 1. - С. 46 - 62.
2. Иванов В.А., Прусов А.В., Рябцев Ю.Н., Шапиро Н.Б. Физические механизмы смешения морских вод с водами субмаринной разгрузки // Современные проблемы океанологии. -Севастополь: МГИ НАН Украины, 2009. - 90 с.
3. Михайлова Э.Н., Шапиро Н.Б. Опыт воспроизведения пространственно-временной изменчивости термохалинных полей в Севастопольской бухте // Морской гидрофизический журнал.- 2008. - № 5. - С. 23 - 39.
Морской гидрофизический институт НАН Украины, Материал поступил
Севастополь в редакцию 11.05.10
После доработки 24.05.10
АНОТАЦ1Я Обговорюються результати чисельних експерименив з вщтворення структури пдролопчних полiв у карстовш порожниш (район мису Айя, Швденний берег Криму), спо-стережуванш тд час друго! спецiалiзованоl експедицп у вересш 2008 р. Для опису процесу субмаринного розвантаження прюних вод використовуеться триви]шрна пдростатична модель iз заданими апрiорi границями шарiв (а- та/або 2^вш) в наближенш «твердо! кришки» i з процедурою «конвективного пристосування». Встановлено, що лише хороше просторове роздiлення по вертикалi дозволяе кiлькiсно правильно описати тонкий верхнш розпрiснений шар моря. Показано, що для адекватного уявлення про вплив рельефу дна на верхнш шар моря бажано використовувати змшану Ытку по вертикал^ поеднуючу а- i 2-координати. Представленi результати чисельних експерименив, якi дозволяють судити про вплив мюцеположення, кiлькостi та потужност джерел прюно! води, розташованих на бiчнiй границi порожнини.
Ключовi слова: чисельне моделювання, пдростатична модель, субмаринне розвантаження прюних вод, карстова порожнина.
ABSTRACT Results of numerical experiments aimed at reconstructing the structure of hydrophysical fields in a karst cavity (region of cape Aiya, the Crimean Southern coast) observed during the second specialized expedition in September, 2008 are discussed. To describe the process of sub-marine freshwater discharge, the three-dimensional hydrostatic model with a priori prescribed boundaries of the layers (a- and/or Z-levels) is used in approximation of a «rigid lid» including the procedure of «convective adjustment». It is found that only good vertical resolution permits to carry out quantitatively correct description of the upper thin sea layer with low salinity. It is shown that in order to obtain adequate notion on the relief effect upon the upper sea layer, the mixed vertical grid combining a- and Z- coordinates should be advisably applied. Presented are the results of numerical experiments permitting to define the effect of location, quantity and capacity of freshwater sources located on the cavity lateral boundary.
Keywords: numerical modeling, hydrostatic model, sub-marine freshwater discharge, karst cavity.
84
ISSN 0233-7584. Мор. гидрофиз. журн., 2011, № 5