Информационное обеспечение безопасности мореплавания и экологической безопасности в зонах лимитированных глубин в системе «Река – море» средствами математического моделирования Текст научной статьи по специальности «Физика»

ИНФОРМАЦИОННОЕ ОБЕСПЕЧЕНИЕ БЕЗОПАСНОСТИ МОРЕПЛАВАНИЯ И ЭКОЛОГИЧЕСКОЙ БЕЗОПАСНОСТИ В ЗОНАХ ЛИМИТИРОВАННЫХ ГЛУБИН В СИСТЕМЕ «РЕКА - МОРЕ» СРЕДСТВАМИ МАТЕМАТИЧЕСКОГО МОДЕЛИРОВАНИЯ

В.Н. Молчанов, кандидат физико-математических наук, доцент. Институт прикладной экологии, Санкт-Петербург

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

Ключевые слова: экологическая безопасность, загрязненность, аварийные выбросы, движение водотока

THE INFORMATIONAL SUPPORT OF THE NAVIGATION SECURITY AND THE ECOLOGICAL SECURITY IN ZONES OF LIMITED DEPTHS IN THE SYSTEM «RIVER - SEA»BY MATHEMATICAL MODELING

V.N. Molchanov. Saint-Petersburg institute of Applied Ecology

In this article an author offers new equations of motion of viscous liquid under the action of the gravitational potential in the hydrostatic approximation. These mathematical models allow us to calculate the extreme safety level of the surface coastal water areas of the sea, bays and rivers (high tides, low tides), and water pollution from ships with accidental emissions.

Key words: ecological security, pollution, accidental emissions, motion of viscous liquid

Для получения уравнений наклоненного движения - движения, порождённого действием ускорения силы тяжести, пространственной размерности «два» (X,Y,t) и «один» (X,t) в качестве исходных использовались уравнения Навье-Стокса в приближении Буссинеска.

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

двумерное (плоское) представление уравнений движения на законных основаниях сохраняет все дифференциальные операторы примитивных уравнений.

На рис. 1 представлено схематически взаимодействие статических и динамических компонент наклоненного движения жидкости.

Рис. 1. Схема наклоненного движения:

Р - гидростатическое давление; g - ускорение силы тяжести; отметки уровенной поверхности

неустановившегося движения вязкой жидкости от переменных расходов воды (или от действия дополнительных к расходам внешних сил: ветер, прилив, градиенты плотности); h - отметка дна; в-углы наклона водной поверхности установившегося движения водотока

Исходные уравнения движения можно представить в виде:

1 дР

L(u)=Fx , Fx= - -—+ gix, ix = sin p, (1)

p дх

1 дР

L(v)=Fy , Fy = _!др., (2)

P дУ

1 дР

0=Fz , Fz = - -—+ giz, iz = cos p. (3)

p д2

Плоскость координат X0Y расположена на водной поверхности установившегося (от заданного расхода и градиентов плотности) наклоненного движения, и горизонтальная составляющая силы гравитационного потенциала действует только в направлении оси 0X; ноль оси 0Z (рис. 1) находится на динамически подвижной уровенной поверхности установившегося (от заданных по начальным данным расхода жидкости и поля плотности) движения, высотная отметка уровенной поверхности Е может иметь разные знаки; отрицательное направление оси 0Z совпадает с положительным направлением действия гидростатического давления; отметка дна является переменной функцией z: ze [е|, z* ],

z*=-h(x,y,t). Координата t может быть сохранениа, если формализовать заносимость водотока.

Необходимо заметить, что ось 0X ориентирована всегда по направлению стрежневого потока водотока, а ось 0Y соответственно по нормали к ней. При этом, направление стрежневого потока и проекции гравитационного потенциала на ось 0X совпадают. Таким образом, пространственная извилистость русла (для математического моделирования) водотока пропадает, а различие участков русла реки (водотока) достигается заданием (начальными данными) соответствующих каждому отрезку русла углов наклона уровенной поверхности установившегося (от заданного расхода и градиентов плотности) движения водотока. В плоском (X,Y,t) представлении движения жидкости в водотоке в любой точке поперечного сечения водотока определён вектор скорости V(х, y, t) движения жидкости и,

таким образом, движение может уклоняться от стрежневого потока. Уклонение определяется соотношением глубин водотока, а также направлением, модулем и продолжительностью действия на него внешних сил. В уравнении гидростатики (3) содержится проектор вектора ускорения силы тяжести g на ось 0Z - cos в. Предложенный на рис. 1 треугольник сил демонстрируется в школьной физике, однако в гидромеханике гидростатического приближения он предлагается впервые. Давление P движущейся под действием ускорения силы тяжести жидкости меняется в зависимости от угла наклона уровенной поверхности в . P = gpH cos в. При в = 0 (водоём) P достигает максимума, а при в = 90 0 (свободное падение воды) P=0 и достигает минимума. Если P е (0°,90°), то оно «обслуживает наклонённое движение» - водоток. В одномерном представлении движения жидкости (X,t), например лотке, меняющем углы наклона подложки и его поперечное сечение, второй оси Y просто нет. Если адресуем действие гравитационного потенциала в двумерном (X,Y,t) представлении только оси X, то, предположив постоянство всех суммируемых функций в уравнении (1) и (2) по оси Y, мы получим точный и понятийно корректный перевод двумерного представление в одномерное для уравнеий движения вязкой жидкости под действием гравитационного потенциала. Обратно: в одноменрном представлении движение по определению не зависит от извилистости водотока - нет второй координаты.

Правая часть наиболее динамически полного уравнения движения для водотока после интегрирования уравнения (1) по вертикали z с использованием уравнения (3) будет иметь вид:

т^* тт пД 1 ф 31 . н SB ...

Fx =- g H cos Р(— + ) + gH sin Р(— -£ +1). (4)

2 р dx дх 2 дх

Процедура вывода правой части уравнения движения для водотока изложена в работе [1].

Первый член выражения (4) «побуждает движение» в водотоке и (при в=0) в водоеме, второй - только в водотоке. Если u2=gH, и - групповая скорость длинной волны, то первый член осуществляет (волновым механизмом переноса) действие (на левую часть уравнения движения) градиента плотности и градиента уровенной поверхности, а второй член переносит действие гравитационного потенциала «в чистом виде».

Заметим, что для водотока в двумерном представлении

р =-н (н - ^3),

у 2 р дх дх

Эта правая часть не содержит угла наклона водной поверхности в и аналогична записи для водоёма. Что отвечает принятой (1)-(3) аксиоматике. В уравнении (2) нет проектора гравитационного потенциала.

Полную правую часть для проекции на ось 0Х обобщенного уравнения движения водотока и водоёма размерности (X, У, 1) можно записать в виде:

F х =-gH cos в

H 1 dp д% H дв „

---^ + tge(--- +1)

2 p дх дх 2 дх

-alV + (тх -тxd)cosв •

Здесь при соотношениях : а = 1,р = 0 - водоем; а = 0,р^ 0 - водоток; IV - эффект уклонения движения (Кориолиса) от вращения Земли; l = 2raz sin ф - параметр Кориолиса: Ю z - угловая скорость вращения Земли; sin ф - широта места; V = vH, v - составляющая трехмерного вектора скорости движения вязкой жидкости по оси 0Y; тх - касательное напряжение трения ветра о водную поверхность; тxD - трение о дно.

Обобщённые двумерные уравнения, начальные и краевые условия для водотока и водоёма будут иметь вид:

дt дx

= -gH еоэ в

ш - K АГ \

H дx У H

д дy

иV _ к_д_ГЦ

H ду \ Н

-а IV =

Н 1 др дк Н дв ---^ + tgв(--- +1)

2 р дх дх 2 дх

(т х-ги\и\ / Н2 )еов в

(5)

ЗУ +_д_ дt дх

- к Л V

Н дх У Н

ду

V. - к АГ V

Н ду УН

+ а =

= - gH

Н1 др+дк

2 р ду ду

+т

-^ц / н 2

(6)

дН дЬ + дV = 0 ; дt дх ду

(7)

дС, д —- + —

дt дх

С, и - к АГС!

Н дх У Н

ду

с, V - к д(С

Н ду У Н

= 0;

(8)

Р = Р(^, Т).

(9)

Начальные условия:

к,и V |,=0 =к ,и " V ";

С I = С*- ;

, к=о ^ 1 '

в ук=о=в;; К л,=0=К у

(10) (11) (12)

Краевые условия:

Ст\ г, = С т ;

ип 1г, = 0;

(Ь)п |Г2 = 0, Ц |Г2 = и* (п, V;

(13)

(14)

(15)

+

д

+

д

+

Ц

п 1Г,

= и*

(16)

(С )П

= 0:

(17)

С,-

(С)' = 0,

ип 1г 2 < о

(С )п = 0, ип Г > 0'

(18)

где и,У= (и;у)Н; и,\ - эйлеровы скорости размерности (х,у,2Д); Ь^- = И у - заданная матрица глубин; К, К = (к; к ) Н - заданные коэффициенты турбулентного обмена количеством движения и веществом, температурой; р^ - заданный угол наклона водной поверхности

^х,

Ж

- касательное

водотока; g - ускорение силы тяжести; р - плотность воды; тху ху

напряжение трения ветра о водную поверхность; у - эмпирический коэффициент -энергетическая добавка или вычет, который в общем случае может содержать и проектор компонентов приземного вектора ветра на оси координат; х - коэффициент придонного трения; С, = с,Н , , = 1,2,3; с1 = Б - соленость; с2 = Т - температура воды; с3 = с - пассивная

примесь; Г1 - твердая граница контура (берег); Г 2 - жидкая граница контура; Г 3 - источники количества движения и примесей; п - внешняя нормаль контура;

* - заданные значения; т - источники примесей; п, V - амплитуды и фазы прилива.

Запись одномерного уравнения движения для водотока с точностью до начертания букв аналогична двумерному представлению (5):

+_д_

д дх

02 - К д (о

ш

дх I ш

= -gш еоэ Р

Н 1 др Н др „ч — +-----tgв(--- +1)

дх 2 р дх 2 дх

+ (т х-т хД )еов р, (19)

2

где 0 = иНВ = иш - расход воды, В - ширина водотока, ш = НВ - площадь поперечного сечения водотока.

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

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

- для «главной оси 0Х» уравнеия имеют одинаковый вид;

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

- переход от двумерного представления к одномерному и обратно (при постановках сопряжённых краевых задач) осуществляется простыми процедурами замены исходных параметров уравнений и краевых условий;

- моделирование динамической системы «водоток-водоём-водоток» в любых композициях их пространственных размерностей («два» и «один») осуществляется простыми процедурами.

Алгоритмическое сопряжение водотока и водоема осуществляется по следующему правилу: как только заданная максимальная текущая координата водотока (супремум) превосходит последующую, и она (текущая координата) становится наименьшей (инфимумом) координатой водоема, в уравнениях (5) и (6) для всех последующих (возрастающих) координат осуществляется операция: Р = 0 и а = 1. Обратное сопряжение -

«водоем-водоток» очевидно (а = 0, Р; = Рр). Координаты принадлежности водоток, водоем

заданы. Подобные процедуры позволяют без труда производить расчеты в сложной динамической системе «река-море».

О краевых условиях

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

Для одномерного представления для количества движения, например:

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

Можно заметить, что плоское представление эстуарийной зоны водотока, где происходит смешение пресных и солёных вод водотока и водоёма, формирующих зону динамического разрыва течений по вертикали, вполне пригодно для вычисления уровенной поверхности их объединения. Это связано с тем, что гидростатическое давление жидкости P = gpH, водоём или P = gpH cos в, водоток является непрерывной по вертикали функцией, в отличие от течений. Тогда, задавшись средневзвешенными по вертикали значениями плотности воды в качестве начальных данных, мы можем получить вполне правдоподобные уровенные поверхности эстуарийной зоны водотока или (и) «губы моря». А это очень важно для получения оценок безопасности мореплавания в зонах лимитированных глубин, произведённых средствами математического моделирования.

На сетке Карского моря для выявления вклада плотности воды в динамику уровенной поверхности моря был проведён такой численный эксперимент: задавалась начальными данными матрица средневзвешенной по вертикали солёности воды. На южной «жидкой границе» с Обской губой она бралась 7 0/оо, а на северной, в районе Новой Земли - 33 0/оо. Считалась задача «стационирование по начальным данным». Превышение уровенной поверхности на южной границе области решения составляло 80 см. То есть горизонтальный пространственный градиент солёности в 1 0/00 «порождал» повышение уровня в 3 см. Это совпадает по порядку величин с оценками, производившимися десятилетиями раньше по «динамическому методу» Н.Н. Зубова, заимствованому от норвежского метеоролога Бьёркнеса: геострофическая модель движения атмосферы. Здесь геострофическая компонента движения является корректной частью общей модели движения вязкой жидкости.

На сетке Обской губы проводились комплексные численные эксперименты на «плоской модели движения» вязкой жидкости (5)-(18). Брался ветер 10 м/с и считались на ЭВМ БЭСМ-6 численные решения (5)-(18) на временном интервале 2 суток и полного стационирования. По отчётам Амдерминской обсерваторией УГМС: репрезентативным ветром (средневзвешенно осреднённым на интервале действия) для штормовых явлений на акватории Обской губы является ветер силою 10 м/с, продолжительностью действия 2 суток. В работе [3] приведены многочисленные результаты рсчётов уровенных поверхностей, скоростей потоков и загрязнений на сетке Обской губы. Сопоставление максимальных сгонно-нагоных колебаний уровенной поверхности в южной части Обской губы (наиболее чувствительной по глубине для судоходства) показало точное совпадение расчётов с наблюдённой Амдерминской обсерваторией УГМС с натурой.

Хорошие результаты моделирования сгонно-нагонных колебаний уровенных

At = k, k = 0,1,2...

Q *(t, L), At = k +1

поверхностей водоёмов дают основание для предложения методики «обеспечение безопасности мореплавания в зоне действия лимитированных глубин» на трассе Северного Морского Пути (СМП).

Предлагается:

- построить сеточные области (матрицы, вектор-столбцы) численного решения (5)-(18) по глубинам устьев рек, губ, заливов, прибрежных акваторий морей в двумерном и одномерном представлениях необходимых для СМП водоёмов и водотоков;

- адаптировать предложенные автором алгоритмы и программы к выбранным водоёмам и водотокам;

- провести отладочные расчёты на ЭВМ по выбранным объектам;

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

- построить диагноз-прогностические альбомы матриц уровенных поверхностей для выбранных объектов трассы СМП;

- для северо-западного участка трассы СМП (Баренцево, Белое моря) «достроить» алгоритмы и программы приливо-отливной компаненты колебаний уровенного режима;

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

На рис. 2 показаны результаты расчётов уровенных поверхностей на сетке Обской

губы.

показано направление действия ветра силою 10 м/с. В разрывах изолиний - уровень в м. Заштрихованные квадратики в С-З части губы - области «полного сгона» - нулевого значения динамической глубины (Ь+ ^),

где исходная глубина 0.5 м

Рис. 2. Результаты расчётов уровенных поверхностей

Ветер - 10 м/с, продолжительность действия 2 суток. Хорошо видна неочевидность ожиданий: С-З ветер продуцирует одинаковый по модулю нагон в южной части губы с В ветром, при этом, в северной части губы, С-З - 0.8м нагон, а при В ветре - аналогичный по модулю (-0,8 м) сгон. Поворот ветра на 45 0 (С-В) не меняет отметки уровня в северной части губы (сгон: -0,8 м), но в 2 раза изменяет отметки уровня в южной части губы. Эти неожиданности показывают объективность, достоинства математического моделирования, как инструмента принятия решений по безопасности мореплавания в зонах лимитированных глубин. Более подробное описание представлено В.Н. Молчановым в работе «Трилогия» [4].

Литература

1. Молчанов В. Н. Новые уравнения наклонённого движения вязкой жидкости размерности (Х,У,1) и (ХД), постановки краевых задач // Естественные и антропогенные аэрозоли: сб. трудов III Междунар. конф. СПб., 2001. С. 519-528.

2. Молчанов В.Н. Гидродинамическая модель циркуляции в водоеме произвольной формы с произвольным набором загрязняющих источников, имеющих диффузионный характер: материалы Всесоюз. симпозиума «Океанографические аспекты охраны вод от химических загрязнений». М.: АН СССР, 1975. С. 88-93.

3. Молчанов В.Н. Двумерные модели динамики потоков и загрязнения устьевых взморьев: дис. ... канд. физ.-матем. наук. Л., 1979. 17 с.

4. Молчанов В.Н. Трилогия. Избр. СПб.: Стратегия будущего, 2011. 229 с.