Нелинейные инварианты дискретной системы уравнений динамики моря в квазистатическом приближении Текст научной статьи по специальности «Физика»

ТЕРМОГИДРОДИНАМИКА ОКЕАНА И АТМОСФЕРЫ

Научная статья

УДК 551.465 ЕБ№ ЖБииМ

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

С. Г. Демышев н

Морской гидрофизический институт РАН, Севастополь, Россия н demyshev@gmail.com

Поступила в редакцию 16.01.2023; одобрена после рецензирования 03.02.2023; принята к публикации 05.07.2023.

Аннотация

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

Методы и результаты. Используется метод неопределенных коэффициентов, при котором вводятся новые неизвестные, что позволяет удовлетворять дополнительным условиям. Получены схемы, которые обеспечивают сохранение одновременно температуры в первой и в К-й (К> 1) степени и солености в первой и Ь-й (Ь > \) степени. Найденные аппроксимации температуры и солености на гранях бокса при полиномиальной зависимости плотности от температуры и солености приводят к дивергентному виду уравнения адвекции плотности. Такой вид обеспечивает выполнение закона сохранения и полной энергии, и суммы кинетической и динамической потенциальной энергии в дискретной постановке. На основе анализа циркуляции в Черном море в 2011 г. показано, что при увеличении степени инвариантов имеют место следующие эффекты: обостряются градиенты во фронтальных зонах в поле температуры, усиливаются процессы подъема более соленых вод в центре моря и опускания более пресных по его периферии, уменьшается плотность мелкомасштабных особенностей в поле вертикальной скорости. Выводы. Получена дифференциально-разностная схема модели динамики моря в квазистатическом приближении, которая обладает рядом нелинейных инвариантов, соответствующих непрерывной задаче. Результаты проведенного расчета черноморской циркуляции для реальных условий 2011 г. показали, что наличие инвариантов степени больше двух позволяет уточнить особенности циркуляции на малых масштабах.

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

Благодарности: работа выполнена при финансовой поддержке гранта Российского научного фонда 23-27-00141.

Для цитирования: Демышев С. Г. Нелинейные инварианты дискретной системы уравнений ди-намики моря в квазистатическом приближении II Морской гидрофизический журнал. 2023. Т. 39, № 5. С. 557-583. ЕБЫ .Г№8ииМ.

© Демышев С. Г., 2023

Original article

Nonlinear Invariants of a Discrete System of the Sea Dynamics Equations in a Quasi-Static Approximation

S. G. Demyshev и

Marine Hydrophysical Institute, Russian Academy of Sciences, Sevastopol, Russia H demyshev@gmail.com

Abstract

Purpose. The study is purposed at obtaining the approximations providing the presence of discrete nonlinear invariants for the difference system of the sea dynamics equations in the absence of external forces, friction and diffusion, and at analyzing the features of the resulting schemes at the example of calculating the Black Sea circulation for 2011.

Methods and Results. The method of undetermined coefficients at which the new unknowns are introduced is applied, that makes it possible to satisfy the additional conditions. The schemes providing simultaneous preservation of temperature in the first and the K-th (K > 1) degrees and salinity in the first and the L-th (L > 1) degrees, were obtained. The approximations of temperature and salinity found on the box faces with a polynomial dependence of density on temperature and salinity lead to a divergent form of the density advection equation. This form provides fulfilling the law of conservation both of the total energy and the sum of kinetic and dynamic potential energy in a discrete formulation. Based on the analysis of circulation in the Black Sea in 2011, it is shown that at increase of the degree of invariants, the following effects take place: the gradients in the temperature field in the frontal zones as well as the processes of the saltier water upwelling in the sea center and the fresher water downwelling along its periphery are intensified, and the intensity of small-scale features in the vertical velocity field decreases.

Conclusions. A discrete dynamical model in a quasi-static approximation was obtained. It has a number of nonlinear invariants corresponding to the continuous problem. The results of calculating the Black Sea circulation for real conditions in 2011 showed that presence of the degree invariants exceeding two made it possible to specify the circulation features on small scales.

Keywords: numerical modeling, kinetic energy, potential energy, discrete conservation laws, Black Sea, cyclonic circulation, anticyclonic eddies, sea dynamics

Acknowledgments: The work was carried out with financial support of the Russian Science Foundation grant 23-27-00141.

For citation: Demyshev, S.G., 2023. Nonlinear Invariants of a Discrete System of the Sea Dynamics Equations in a Quasi-Static Approximation. Physical Oceanography, 30(5), pp. 523-548.

1. Введение

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

В работе [2] построен разностный аналог теоремы Нетер для одномерного уравнения Эйлера, который вследствие перехода в пространство сеточных функций существенно отличается от оригинала дополнительными условиями. Построение соответствующих разностных схем тесно связано со свойством их инвариантности, определение которого было дано в работе [3]. Свойство инвариантности относительно групп преобразований является необходимым условием применения теоремы Нетер для получения законов сохранения. Инвариантные разностные схемы построены для обыкновенных дифференциальных уравнений второго порядка [4]. В работе [5] получены разностные аналоги для линейных и нелинейных одномерных волновых уравнений с условием симметрии и законов сохранения. Для одномерных уравнений мелкой воды в лагранжевых координатах получена новая инвариантная конечно-разностная схема. Она обладает локальными законами сохранения энергии, массы, центра масс и импульса [6]. В приложении к задачам вычислительной геофизики закон сохранения полной энергии используется как норма дискретного решения, что повышает его устойчивость. В работах [7-9] на основе энергетически устойчивых схем исследуются свойства численного решения уравнений На-вье - Стокса. На тестовых задачах показывается точность и работоспособность представленных алгоритмов. Ряд работ посвящены схеме КАБАРЕ, которая вводилась как новая трехслойная явная разностная схема с пространственным расщеплением временной производной [10]. Она тщательно на высоком уровне проработана для одномерного уравнения переноса без диссипации [11, 12], что позволяет аккуратно воспроизводить свойства решения соответствующего гиперболического уравнения.

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

В качестве примера использования энергетического анализа можно указать на следующие работы. В [13] представлены оценки пространственного распределения и трансформации вихревой энергии в модельном океаническом бассейне. Авторы показывают, что локальная генерация и диссипация вихревой энергии компенсируется так называемыми нелокализованными потоками из других частей бассейна, особенно по периферии вихревых структур. В [14] на основе анализа распределения средней кинетической энергии вихрей, рассчитанной по данным спутниковой альтиметрии, проведены оценки пространственного распределения источников и стоков вихревой энергии для реально наблюдаемых мезомасштабных структур. Влияние горизонтального разрешения моделей на описание каскадов энергии исследовалось в [15]. Авторы доказывают, что при расчете бюджетов кинетической (КЭ) и потенциальной (ПЭ) энергии прямые и обратные каскады энергии точнее описываются в модели, разрешающей структуры с размерами, меньшими бароклинного радиуса

деформации Россби. Модельные оценки сезонной изменчивости вихревой кинетической энергии, сделанные в [16], свидетельствуют, что сезонное увеличение теплового воздействия атмосферы вызывает интенсификацию мезомас-штабной изменчивости в поле скорости у поверхности океана. Анализ энергетического цикла Лоренца для всего Мирового океана по данным реанализа выполнен в работе [17], где получены среднемноголетние оценки компонентов энергии и представлены основные различия между энергетическими циклами океана и атмосферы. Методика оценок вихревой энергии и механизмов ее трансформации, предложенная Лоренцом, широко применяется для региональных исследований циркуляции. Так, например, район течения Куросио рассмотрен в [18], вихревая изменчивость в Красном море изучается в [19], подробный анализ механизмов образования мезомасштабных вихрей в Охотском море представлен в [20].

Для Черного моря в ряде исследований выполнены оценки кинетической энергии вихрей и течений на основе данных наблюдений [21, 22]. Исследование годовой и междесятилетней изменчивости доступной потенциальной энергии по данным массива наблюдений за 1910-1998 гг. представлено в [23]. Энергетические балансы полузамкнутых морей (в том числе и Черного моря) рассчитаны в работе [24], где показано, что работа силы плавучести является одним из основных факторов формирования мезомасштабных особенностей динамики в полузамкнутых морях.

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

2. Законы сохранения в дифференциальной постановке

Запишем уравнения модели идеальной жидкости в приближении Бус-синеска и несжимаемости морской воды:

+ Г)у + ™и2=-£<;х-1(Р'+ Е)х, (1)

Ро

V, + + /)и + = -^^у - - (Р'+ Е)у, (2)

^ Ро

и, + V,, + = О, (3)

г

Р = ^Ро^ + ^ | РФ = ^Ро^ + Л (4)

о

(иТ) х + (VI) (^Т) г= О, (5)

!\+(и!1х + (и?), + = О, (6)

Р =С (Т, !!). (7)

При г = 0 = ; на дне при г = Щх, у) и =v = ^ = 0, Тг = = 0. (8)

На боковых стенках: для меридиональных -

и = Vx = 0, Тх = !х = о, (9а)

для зональных участков границы -

V = иу = 0, Ту = !у = 0, (96)

Начальные условия: при ( = (о

(Т, 5) = (ПП и=и°, v=v0, С = (Ю)

Здесь: ^ = vx - u , E =

1 u2 + V2

В системе (1)-(10) и^,м> - компоненты вектора скорости, направленные по осям х, у, г соответственно; Т, Р - температура, соленость и давление морской воды;/- параметр Кориолиса; С, - приведенный уровень моря; g -ускорение свободного падения.

В выражении (7) 0(Т, 5) - полином произвольной степени от температуры и солености. В общем случае уравнение состояния включает давление, но во многих моделях используется соотношение (7), один из вариантов которого рекомендуется Межправительственной океанографической комиссией 1. В соответствии с ним схемы сохранения строились с учетом полиномиальной зависимости плотности от температуры и солености. Для дальнейшего упрощения записи выкладок полагаем, что ро =1 гр/см3.

1 The International Thermodynamic Equation of Seawater - 2010: Calculation and use of thermo-

dynamic properties. UNESCO, 2010. 196 p. (Intergovernmental Oceanographic Commission Manuals and Guides ; No. 56).

MOPCKOËrHflPOOH3H^ECKHË®YPHA^ tom39 №5 2023 561

2. 1. Уравнения адвекции температуры, солености, плотности

Введем обозначение

1

ф>я =| фdz, <ф = — Ц

1 н н

ф>к = — | Ц §dxdydz, V = 111 dxdydz,

<

V

о а г о а г

где О2 - горизонтальная поверхность на уровне 2.

В силу выполнения условия (3) мы можем уравнения (5) и (6) записать в недивергентном виде:

Т + птх + УТУ + м>Тг = 0, (11)

+ иБх + уБу + м>Бг = 0. (12)

Из выражений (3), (11), (12) и краевых условий (8)-(9) следует, что для любых К > 1, Ь > 1

ТК + иТК + уТуК + '^Т2К = 0, (13)

^ + и£,ьх + + м>Бь2 = 0. (14)

Возвращаясь к дивергентной форме записи выражений (13), (14) и интегрируя по области, имеем первые два инварианта (К > 1, Ь > 1)

< ТТ >V +< >п° = 0, < ^ >V +< >п° = 0. (15)

Заметим, что уравнение (15) соответствует

< тК ^ + < <тК >°°=V | 1 Я ТК)= о,

- ? ^ н

< # >7 + < <Х >°° = V I IА (^ )dxdydz = 0.

У О - ? п,

Пусть О - полином от Т и 5", тогда из выражений (11) и (12) получаем уравнение адвекции плотности

Р т +{и р )х + (■ур )у + (м> р )г = 0. (16)

Дивергентный вид выражения (16) обеспечивает сохранение интеграла от

1 д н

плотности, аналогично уравнению (15)---| ^рdxdydz = 0.

V ^

2. 2. Уравнения для скорости изменения кинетической и потенциальной энергии

Из уравнений (1)-(4) получаем уравнение бюджета кинетической энергии Е, + [к(^ + Р'+ Е)]х + [V+ Р'+ Е)]у + [м+ Р'+ Е)\ = ^р м. (17)

Из уравнения (17) при интегрировании по бассейну следует

< Е + < сЕ + ё(52X >°° =< ОТ . (18)

Из вида уравнений (17), (18) можно увидеть, что

1 ^ н

< Е, +< с,Ео >П0 =--1 \\Edxdydz = 0.

У ^ Qz

Полагаем, что потенциальная энергия имеет вид П = — g z р . Тогда уравнение скорости изменения потенциальной энергии записывается следующим образом:

ТТ, + (иП) х + (уП) , + (мП) г = - ^рм. (19)

Из уравнения (19) получаем III инвариант:

2 1 д н

< Е< + П, + < с Е + g(с2), >°° = - - | Ц (Е + Ъ)dxdydz = 0. (20)

2 V ОТ _ ? Qz

2. 3. Уравнения для скорости изменения динамической потенциальной энергии

Важной физической характеристикой динамики моря и атмосферы является доступная потенциальная энергия, которую ввел Э. Лоренц в работе [27]. Суммарная энергия так же, как общая, при условии адиабатичности процессов сохраняется или является инвариантом. В силу нелинейности уравнения состояния возникает неоднозначность в определении доступной потенциальной энергии (в работе [27] предполагалась линейная зависимость плотности от температуры). Поэтому в последнее время в ряде работ для анализа результатов численных экспериментов используется аномалия потенциальной энергии - динамическая потенциальная энергия [28].

В непрерывном случае полагаем, что

1 1

р(х, у, z, г) = р*(х, у, z, г) + р5 (£), где р* = -1 —Цр(х, у, z, t)dxdy

Т т ^ ^ у

Я, (21)

х - время интегрирования. В отличие от определения, данного в работе [29], в уравнении (21) предполагается независимость р* от времени. Введем обозначения: Вре = — g z р*.

С учетом уравнения неразрывности, подставляя разложение (21) в выражение (16) получаем

Р * + ир I + Ф У> + ™ РI + ™ Р= 0- (22)

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

Бре + (иБре + {уБре) + {м^Бре ^ - (gzwp5 ) + gp'w = 0. (23)

В силу отсутствия потоков через стенки бассейна, учитывая выражение (23), получаем IV инвариант

<Е + <;А + Цре >У +<g(£■)t = ^|]{{(£ + вре= 0. (24)

2 V дti< £

3. Законы сохранения в дискретной модели

Аппроксимируем бассейн с неровным дном боксами, центрам которых соответствуют целочисленные значения индексов /, ], к (/ = /1,..., j = ] 1,...,]ы, к = 1,..., К/,}), граням - / + 1/2,у + 1/2, к + 1/2. Горизонтальные размеры боксов (кх, ку) постоянные, по вертикали используется неравномерная аппроксимация (кк2 = 1/2 - 1/2, ккк = zk+1 - ^ ). Разностные операторы имеют вид (для/ к - аналогично)

фу = 4+/2,у,к + 11/2у 5 ,= 4+/2у - Ф/-1/2,у,к ф = ^

X

(25)

{ф}°к =7^, {ФГ = VЕ К^Лк, V = £ Кк]кхку.

11к г,] У г, у к=1 /у к=1

На горизонтах zk рассчитываются температура, соленость и горизонтальные компоненты скорости, на горизонтах Zk+l/2 - вертикальная скорость [25]; 0.к - площадь поверхности на горизонте zk. Распределение переменных ука-занонарис. 1.

Рассмотрим уравнения, непрерывные по времени. Конечно-разностная система уравнений модели со вторым порядком точности по пространственным переменным (с точностью до неравномерного шага) неоднократно выписывалась [25].

Рис. 1. Распределение переменных в боксе (г, j, к) F i g. 1. Distribution of variables in the box (i,j, к)

3. 1. Инварианты дифференциально-разностных уравнений адвекции температуры, солености и плотности

Запишем уравнение состояния в точке (i,j, к) в виде полинома

Рил = ЦТи*') = = > (26)

п=о ш=о

где п > 0, т > 0, апт - постоянные. Заметим, что N М необязательно

должны быть целыми.

Для сохранения интеграла по пространству от плотности при условии адиабатичности и отсутствия внешних источников необходимо получить аппроксимацию нелинейных слагаемых таким образом, чтобы наряду сТи5 сохранялись ТК и Б1 соответственно:

V f i ii(TK}^=0' V I J Я()dxdydz=(27)

1 <L

V 81

-

В результате интеграл по объему от каждого члена ряда (26) в силу положительности температуры и солености ограничен.

Уравнение неразрывности в дискретной форме в точке {¡, ], к) с учетом обозначений (25) имеет вид

§ А у л + § у^; у ,к + § zW, у л = 0. (28)

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

Введем в точках (/+1/2,у, к) новые неизвестные Тм12ук,8м12ук,рм12ук, вточках

(1,у+\/2,к)-Тш1иД.^,ри+1/2>к, в точках (/',у; к+ 1/2)- ^лМ/2Д.¿^р,.Лк+1/2. Запишем уравнения адвекции температуры и солености в точке (/, у, к)

с1Т

+ (и, уТ ук ) + §у (V. мТ, м ) + § z ( wh уЛТш уЛ ) = 0, (29)

+ 5Х (Ыи Д. м) + 5 у (V, Л Дм) + 5 z К ьЛук) = о. (30)

Далее все рассуждения и выкладки, проведенные для температуры, для солености будут идентичны. Заметим, что Т1+} л,Т, з у+112Л,Ти} л+1/2

и 1,2,у,к, 5и}+1/2,к, 5г,1/2 покане определены.

Учитывая дискретное уравнение неразрывности (28), перепишем уравнение (29) в виде

-Т

Л +

Ч+1/2,у,к И г+1/2,у,к

+

1/2,у,к Т/,у,к ) иг-1/2,у,к -1/2,у,к Т/,у,к )^ ^х +

^¡у+1/2,к у+1/2,к _ Т у ,к ) _ Ц, у -1/2,к у -1/2,к _ Т, у ,к ) ^ \ +

(Ти уМУ2 - Ту ,к )- W,j,k-1/2 (Ти у,к-!/2 - Т, у,) ](# )"' = 0. (31)

Wг,jк+1/2 V Тг,j,к+1/2

Аналогичным образом поступаем с уравнением (30). Пусть функционал 2, например, от температуры и его производная записываются следующим образом:

2 _ 2 (Т \ -2' - у-к - 2' 2' - -2'-У 'к

2, у к ~ 2г, у к , у щ Л _ 2, у ,к ' гДе 2, у щ ~ 1ГГ •

М М -Тик

Умножим уравнение (31) на 2г у к И результат преобразуем с учетом уравнения неразрывности (28). Получаем уравнение для 21 у к

Л

^ ^^^^^^у+1/2,к [2г,у ,к + й,у ,к ,у +1/2,к _ Тг,у ,к _ V¿,у-1/2,к [2г,у ,к _ Й,у ,к ,у ,к _ Тг,у-1/2,к ^у ^ + у,к+1/2 ,к ^ ^^^^^^,к ^у,к+1/2 " Ту ,к )] " у ,к-1/2 [_<2,у,к - 2и,к ^¡у,к - ,к-1/2 )]} ) = 0-

(32.1)

Предположим, что - к удовлетворяет уравнению адвекции, которое запишем в виде

dQ к 1

+1/2,],к а,] ,к +1/2, ]к - а -1/2,] ,к Qi,],к -1/2,]к

dt

+ш +,} +1/2,к ад,] ,к ' г,] +1/2,к " (й,,-1/2,к " 0,- ,к ,, -1/2,к

+[(й,-к+1/2 " а,у,к )М,-к+1/2 " ,к-1/2 " а,у,к К,М-1/2 ](Аг Г +

, ¡,к , у,к+^,к+,к)=

(32.2)

Тогда из уравнений (32.1) и (32.2) получаем выражение, например, для Т+1/2,(аналогичнодля ^у+^,^-л+1/2)

Т

(0+1,-, кТ+1, - ,к 0+1, - ,к) (й, - ,кТг, - ,к а, - ,к )

i+1/2,-,к

^+1,-,к ],к

Полагаем, что а, у к = Т у к . Тогда температура (соленость аналогично) аппроксимируется на граняхбокса (/, к) следующимобразом:

Т

К

1+1/2, ],к

К-1

К

К

К

К

Т+1,],к ~Т/,],к тК-1 _ ТК-1

V '+1.,/.к -,к У

Т

К

К-1

К

' 1!, /'+1/2,к

К

тК

т/, / к+1/2

К -1 К

К

К

Т^К ГТ1К \

т/, / к+1 _ т/, / к ТК-1 _ тК-1 ^ /к+1 / к ]

К

Т I,/+1,к _ Т I,/,к

тК-1 _ тК-1

V './'+1.к /^Л У

(33)

если Qi, /,к = ^,к Б

X-1

/ ОХ

,/, к /+1/2,/,к

ь Л

0ь

, / к+1/2

X -1

(

Ь

ь

оЬ _ оЬ Д

"Л-, / к+1 /к

оЬ-1 _ оЬ-1

V , ] ,к+1 , / к ]

V Ми^ ¡,/к

оь

' I, -+1/2,к

Ь-1

/ оЬ

ь

ь \

+1,к к \ •+1,к •,] к у

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

Т

К-1

(

I, ],к +1/2

К

Т/ +

Т

'к-1

I, ],к +1

Ф(Т, -к+1,Т, -' к ) у

С"

' ],к+1/2

Ь-1

(

ь

Б шк +

к+1

1,к+1'°, 1,к )

•}к' у

(34)

где

К—Т ь-2

^(т Т ТК~п~2тп о о ^ — ^ оь~т~2от

^ ^^ ] к) ¿^ j,k+\1i,/к> т Г/, ] к) /к+г^!,-к'

п=0

т=0

Таким образом, из аппроксимаций выражений (33) и (34) следует, что при интегрировании по всей области сохраняются T, S и TK, SL (K>2,L> 2), то есть выполняются следующие законы сохранения:

dT K dt

de,

', j tk

dt

i j ,1/2

= 0,

dSLjk dt

dÇ

— SL

dt

i j ,1/2

= o.

Эти выражения соответствуют интегралам (27) (I и II инварианты). Для того чтобы обеспечить закон сохранения полной энергии, необходимо адекватным образом описать работу силы плавучести в уравнениях для скорости изменения кинетической и потенциальной энергии. Для этого необходимо из уравнений адвекции тепла и соли (29), (30) с учетом выполнения соотношений (33) и (34) получить уравнение для плотности, которое должно иметь дивергентный вид.

Пусть уравнение адвекции плотности записывается в виде, в котором

Р ¿+1/2, у,к ,Ри у+1/2, к > Р и у, к+1/2 неизвестные:

dj dt

+ 5 х fa j Jk Pi, ) + 5 v (v, j Jk Pi, ) + 5 * (i,kPi, j Jk ) = 0-

(35)

Плотность в точке (, , у, к) имеет вид р } к = Ог} к (Т,ул, Sh} к ). Полагаем, что выполняется

dp

i, j,k

dt

(Y dTi j k / y dS^ j k

g'< j* Jk 's

'S dt

где

/ v dG{ j k lG,jk )T ~

dT

i, j ,k

(Gi.jk l

dG.

i, j,k

dS

i, j ,k

Умножаем уравнения (29) на (0,,у,к)т и (30) на (Ог,у,к)8. Проводя соответствующие преобразования, учитывая уравнение неразрывности (28) и требования дивергентной формы уравнения (35), получаем соотношение, например, для

Р,+1/2,у,к (для Р,,1/2,к ,Р,,у,к+1/2 - аналогично)

Рг+1/2, j,k

T

^ (Gi+1 ,j,k) )t + (Gi, j,k) )t (Gi+1, j,k )t Ti+\, j,k + G j,k )t Ti, j,k

г+1/2, j,k

+

S,

i+1/2, j,k

(Gi+l,\kOs +(Gi,i,k))s 1 (Gi+l,/,k )sSi+\,j,k +(Gi,i,k )sSi,

+

GM,j,k +G,,j,k

(36)

Дивергентный вид уравнения адвекции плотности, как следствие уравнений (29), (30), обеспечивается аппроксимацией плотности на гранях бокса в виде формулы (36).

РI+1/2,у,к, РI,у+1/2,*' РI,у,к+1/2 в общем случае должны иметь вид

К-1Т*._к - ТК

N М

п = ±_I +1--■к '-'к У Уа пТП-1 Ст-1

К; +1/2,-,к т^-1 грК-\ ^¿.^ п,т-1'" г+1/2,-к /

К Т^, - Т

г+ 1/2, - к /+1/2, -, к

+

! +1, -, к -1/, -, к п=1 т=1 уЬ

i+1,- к ,- к

ь-1 БЬ^-к -БЬ

к V 1 V 1 мчТ^п-! от-1

Ь оЬ-Г ^ ^ап-\,ттТ/+1/2, - ,к 0/+1/2, -, к

Ь +1, - к ~ , - к п=1 т=1

^ ^ (ап,т-1пТ +1/2,-,кБ/+1/2,-к +

п=1 т=1

+ап-1,ттТ/п1/2, - ,к0!+1/2, - ,к ) + ,т ^п 1/2, - ,кБ! +1/2, - ,к > (37а)

К-1 т/К. -.к - Т

N М

п = __1 1+1,з,к 1,з,к Vуа 0т_1

Кг,-+1/2,к т^ грК-1 грК-1 п,т-1'" г,-+1/2,к^/,-+1/2,к

К Т/ +1, -, к _ Т/, -, к п=1 т=1

Ь — 1 ,■+1 к ~

У

+

N М

Л',] +1,к ,]к ^ м^Т^п-1 ат-1 У

-1 -\,ттТ/ ,]+У2кБ/,

Т ^Ь' ' п-1,т /, ]+1/2,к /, ]+1/2,к

Ь °г, 1 +1,к _ , 1 к п=1 т=1

N М

2 2 (ап,т-\пТ/п +1/2,кБ/,}^И2,к + п=1 т=1

———т—у М _п-Й—у

+ап- 1,ттТ/,4-1/2,к°г ,-+1/2,к ) + ^ ^ ап,гя Т/+1/2,к0/,-+1/2,к •

(376)

п=1 т=1

К ТК N М

V — ЛТК — Т

п = К_1 1 ■ ]-к+1 1 ■ -^ пТп-1 от-1 '

Кг, - ,к+1/2 ^ грК-1 грК-1 п.иг-г" /, - ,к+ 1/2° г, - к+1/2

К Т/, - ,к+1 _ Т/, - к п=1 т=1

ь -1 - - NM

+

/, /■ ,к+1 /, /■ ,к

фп-1 от-1

Ь оЬ-1 оЬ-1 ^^ап-1,т тТ1, - ,к+1/2 ^ ,к+1/2

Ь ,-,к+1 _ ,-,к п=1 т=1

^ ^ (ап,т-\пТп,к+1/2,-к+1/2 +

п=1 т=1

+ап-1,т ГпТ"-^к+1/2 , к+1/2 ) + Т П-, к+1/2 , к+1/2 • (37в)

п=1 т=1

Из выражений (33), (34) следует, что при интегрировании по всей области сохраняются Тг ^ к и Т^к, °г}к и . Заметим, что соотношения (37)

в определенном смысле не зависят от схемы сохранения для температуры и солености. А именно: в формулах (33) и (37а)-(37в) могут стоять разные значения N М, К, Ь . Например, аппроксимация уравнения адвекции температуры

обеспечивает сохранение Ть (выражение (33)), выражения (37) приводят к дивергентному виду уравнения адвекции плотности, а уравнение состояния имеет наивысшую степень от Т, равную N.

Таким образом, аппроксимации (37а)-(37в) обеспечивают дивергентный вид уравнения адвекции плотности (35) и выполняется

фг .к 1К Г ]По 1 д н гг

' ' РI 3 1/2 ( = О, ЧТО соответствует--J JJ (p)dxdydz = О

3к [ I ,3

\ I dt " ) V Ы _ ? ^

(III инвариант).

Заметим, что, в отличие от схемы КАБАРЕ [11, 12], аппроксимации температуры, солености и плотности обеспечивают выполнение соответствующих законов сохранения и имеют «нелинейный» вид (30), (37а)-(37в). В соотношениях (33) и (37) на гранях бокса участвуют значения Т, Б только в двух точках.

Анализ вывода соотношений (36) показывает, что полученный результат легко обобщается на общий случай зависимости О от г (г = 1,..., Я, где Я - целое число) функций К :

О'к=О . (38)

Пусть О - дифференцируемый функционал по каждой К

dG

О„,

1,3, к

К dFr.k

и соответствующие производные ограничены. Полагаем, что каждая К в точке (1,' к) удовлетворяет уравнению адвекции

- + 8 и . К-г к + § V.. К-г к + 5 ж . К-г к = 0, (39)

X 1,3,к 1,3,к У 1,3,к 1,з,к z 1,],к 1,],к ' V '

где г =1, ..., Я. Для каждого уравнения (39) (г = 1, ..., Я) проведем соответствующие преобразования, получим

_+ (рг _гг Ъ И -(Кг -¥г Ли И +

^ ^ 1+1/2,3,к 1 1,3,к У 1+1/2,3,к х У1 1-1/2,к 1 1,к Зи1-У2,к х

+(р[1+т,к - К},к )у,,,+т.кК - ^з-тл - Х,}-т,кК + (40)

+Щк+т)-1 - 32У1 + +¥г.ь (8 и + 5 у + 5ж ,к) = 0.

г,],к V X У z г,],к

dGi .к , dFlr.k

Выполняется следующее соотношение: —-= ^ (О-,3к)рГ ^ . Тогда

после преобразования уравнения (40) с учетом определения (38), умножения его на ОрГ и суммирования полученного результата имеем

й0..к , г г 1 г г х ^ + ^ (01,у,к [(р+1/2, у,к -р,у,к )и,+\12,у,кКх ~ (р-1/2,у,к - р,у,к )и,-\12,у,кКх +

V г,у+1/2,к 1 г,у,к * г,у+1/2,к'у V г,у-1/2,к 1 г,у,)!у1,у-1/2,к'у

у1 - у2 -ру х^пК )-1 +

к (5.к + 5Дк + 5ук)] - о. (41)

Предположим, что для функционала О выполняется уравнение адвекции. С учетом выражения (41) получаем рекуррентное соотношение для О на гранях (/ + 1/2,], к) бокса (г, у', к) (для 0..+1Пк и О к+У2 - аналогично)

О = V Г Р' + 0г,у,к (0г+1,у,к )р'Рг+\,у,к + (01,у,к )р' у,к +

^Гг+\12,у,к 2-й |_г+У2,],к ~ ~ ~ _|

'

+ 1+1,у,к 1,у,к 2 '

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

йО

—^ + 8 и О .к + Ъ\. О к+ § ™ О к = 0

и тем самым обеспечивает равенство нулю интеграла по объему от О. Полученное свойство в точности соответствует непрерывной постановке. Заметим, что разностное соотношение (42) не зависит от конкретного вида

Р+1/2ук, Ру+У2к, Рук+1/2- ^ качестве примера приведем уравнение состояния

из работы 1 в виде

р к =У У а Т ".Л"1 .к. (43)

^ у к />/ > пт г, у к г, у к ^ '

п=0 т=0

В соответствии с уравнением (43) требуем, чтобы при адиабатическом приближении и отсутствии внешних сил интегралы по объему от Ту к и к

не менялись по времени. Тогда выражение (37а) с учетом выражения (43) записывается в виде (для у к аналогично)

. т ■» , т3 Т + Т2 Т2 + Т Т3 + Т4 ъ Г -х -х

_ _ 1+1. у .к^11+1. у У ¡. 1 .к^1'+1. у .к11. у .к^11+1. у .к111 ¡. 1 .к ЧГ |а пТ Х , а пТ С Х,

Н+1/2, у ,к - тЗ , Т2 Т | Т Т2 | Т3 ¿_1\ап,йП1 /+1/2,у ,к """"пД /+1/2,7,к°/, у,к I

Э Т г+1.7 ,к + 1.7 ,кТ ¡, у ,к + Т г+1.7 ,кТ ¡, у ,к + Т ¡,] ,к п=1 V У

С +

+ у 0 у Та»-ПТ,+1/2,ап,0пТ+1/2,+ап-1,1Т+1/2,у1/2,+ ап,1Т1+1/2,у,кС1+У2,

3. 2. Дискретное уравнение сохранение полной энергии

В отсутствие треиия и внешних сил дискретное уравнение для КЭ имеет вид [25]

dE.

3—

dt

+5,

+ 5.

и13 ,к I Е13 ,к ^ 8 ^13 ^ (р, 3 ,к ^ + 5

+ 5,

3 л I Е•,.3 ,к 8 ^'3 ^р,.3 л )

* \8?X,} + (Рг3,к )

г, 3 ,к >

ж

1,3 ,к

(Е ), , + ""к (Е ),

+

г3 —

= 8

ж, М Р,

и—

В уравнении (44) использованы следующие обозначения:

(44)

и ■ ■ 1 ,1 и • • 1

( Е } = г'3 'к+1 1 ■3 — ( Е ^

\^и31,3—+1/2 0 5 31,3—+1/2

^, 3—+{^1,3—

Е = (ми У + (уи У

ил 2 •

Уравнение для плотности в отсутствие диффузии и в адиабатическом приближении имеет вид (35), где р. 3 —, р+1/2>3—, рз3+У2—, р,3_!/2 удовлетворяют соотношениям (37а)-(37в).

Тогда уравнение для скорости изменения потенциальной энергии

Пг 3 к = — gzкрг 3 к записывается следующим образом [25]:

dП

3—

dt

5 х (ч 3П г, 3 ,к ) + 5 У (у1 _ 3 Ц, 3 — ) + 5z 3 — Ц_ 3 — ) = - я

Ж Р ,ЛА

И:

(45)

В уравнении (40) Пг+у =-gzk рг+1/2>Л—, Цу+У2— =-gzk р^1/2>—, Пг,3,к +1/2 _ _8 к +1/2 Рг,3,к +1/2 •

Отметим, что разностный аналог работы силы плавучести в уравнении (45) идентичен слагаемому в уравнении (44) для Е13 к в случае дивергентного

вида уравнения для плотности (35).

Проинтегрируем уравнения (44) и (45) в разностном смысле по пространству, тогда дискретный аналог полной энергии удовлетворяет уравнению (IV инвариант)

dt к dt

(ии,1) . </ С

dt

+ -

«Л

К

«Л

= 0. (46)

Закон сохранения полной энергии (46) - аналог непрерывного инварианта (20) - выполняется при условии точного соответствия аппроксимации работы силы плавучести (правая часть уравнений (44) и (45)) в обоих уравнениях.

3. 3. Дискретное уравнение для скорости изменения динамической потенциальной энергии

Полагаем, что (Qk+1/2 - площадь поверхности на горизонте z k+1/2)

1 i2 1 ( \ Pk+i/2 = Т"— - Zp ■,j,k+i/2hhy h, где ht - шаг по времени.

M 2 tl ) tl "k+1/2 V '' j )

Тогда предположим:

* s * t * s

fi, j,k = fi, j ,k + fit 5 fi+1/2, j ,k = fi+1/2, j, k + fit ' fi, j+1/2,k = fi, j+1/2,k + fit >

fi,j,k +1/2 _ Рг, j,k+1/2 + fit+1/2' (^7)

Учитывая уравнение неразрывности (28), получаем

Pi'J,k- + 5 (V ,.,,О*,.,,) + 5„ (v.,.,,o.\,,U 5, (ж ,.,,о.*,,,) +

dt

Z (5,k + ,k + 5zWi, j,k ) + W,j,k Z5zfit = 0.

(48)

Пусть Ор к = — ^^ р*.. к . Тогда следствием обозначений (47) и уравнения

(48) является уравнение для динамической потенциальной энергии (к) в адиабатическом приближении, при отсутствии диффузии и внешних источ-

ников.

dDpe

г, jk

dt = -g-

+

Sx ()+ 5y (Vu**Pjk)+ Sz {whhka-k)] + (gzk A fts) =

KWjjPijk

h

(49)

В уравнении (49) использованы обозначения

aip 1/2, j,k SZ k fi+1/2, j,k >ai, j+1/2,k gz k fi, j+1/2,k' ai, j,k+1/2 SZ k+1/2 fi, j,k+1/2"

Тогда выполняется следующее свойство:

d£.

gzk-

dt

dt

+ <!

-y

dt

dt

j

dt

> = 0. (50)

Аппроксимации (47) обеспечивают выполнение уравнения (49) и, следовательно, (50) (V инвариант, который соответствует его непрерывному аналогу (24)). Заметим, что закон сохранения (50) можно получить и в более общем

случае - зависимости средней плотности по времени рS = рS (z, t).

4. Численные эксперименты

Рассмотрим примеры расчета циркуляции Черного моря для 2011 г. [29] с использованием некоторых полученных аппроксимаций.

В численных экспериментах задан равномерный шаг по горизонтальным координатам 1,6 км, по вертикали использовалось 27 горизонтов со сгущением в верхнем слое моря. Сток причерноморских рек учитывался в соответствии с [30] и составлял около 340 км, из них на реки северо-западной части моря (Дунай, Днестр, Днепр, Южный Буг) приходится примерно 78 %, нареки Кавказа (Риони и более мелкие реки) - 13 % и на реки Турции (Ешиль-Ирмак, Кы-зыл-Ирмак, Сакарья) - 5 %. Реки Кавказа в силу их многочисленности были представлены вдоль кавказского побережья в виде трех источников. Соленость в устьях рек принималась равной нулю. Температура воды в устьях рек, кроме рек Турции, задавалась из работы [31]. Температура рек Турции задавалась равной температуре прибрежных вод моря.

Предполагалось, что в верхнебосфорском течении температура и соленость та же, что и в море. В нижнебосфорском потоке температура принималась равной 16 °Си соленость -35 %о.

Для задания атмосферного воздействия использовались данные SKIRON за 2011 г., вертикальное перемешивание описано на основе теории Меллора -Ямады [31]. Начальные условия в этом расчете соответствовали 1 января 2011 г.

Было проведено три численных расчета. Они различались разностными схемами аппроксимации уравнений адвекции температуры и солености. В первом (I) использовалась «традиционная» схема, которая обеспечивала сохранение T,, j,k, Tu,k, Su j,k > S?,j ,k > во втором (II) - Th j k, Tfj k, Sh j л, Slj л и в третьем (III) - Tj k, TVjk' Sj jk, S)]-k ■ Выбор обусловлен следующими соображениями.

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

Начальные условия для второго и третьего расчетов соответствуют 1 февраля 2011 г. по результатам первого эксперимента.

Время интегрирования для второго и третьего расчетов составило 20 сут. Все рисунки приведены на 21 февраля 2011 г.

Рис. 2 демонстрирует, что бассейн был охвачен обширным циклоническим круговоротом с двумя центрами в западной и восточной его частях. Между берегом и Основным Черноморским течением возникали и эволюционировали мезомасштабные вихри размером от нескольких километров до десятков. Имеет место качественное соответствие в полях уровня (рис. 2), а количественные различия в структуре уровня между тремя расчетами незначительные и наблюдаются в районах интенсивной временной изменчивости Основного

Рис. 2. Приведенный уровень моря в I (а), II (b), III (с) экспериментах F i g. 2. Reduced sea level in the experiments I (a), II (b) and III (с)

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

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

т, °с

28 30 32 34 36 38 40 ° в.д.

ъ

т,° с

c

Рис. 3. Температура на горизонте Зм в I (а), II (b), III (с) экспериментах F i g. 3. Temperature at the 3-m horizon in the experiments I (a), II (b) and III (c)

Наибольшая разница в полях температуры наблюдается в верхнем слое моря в прибрежных областях (рис. 4). Выделяются области на северо-западном шельфе, зона аппвелинга около берегов Крыма и вдольбереговая полоса у Анатолийского и Кавказского побережий. Максимальная разница между значениями температуры в I и II, I и III экспериментах (рис. 4) по модулю составляет 2 °С для обоих вариантов и средняя - 0,14 °С.

а г, °с

b

Рис. 4. Разница между значениями температуры на горизонте 3 м в I и II (а), I и III (b) расчетах F i g. 4. Difference between the temperature values at the 3-m horizon in the calculations I and II (a) and land III (b)

Поле солености имеет структуру, которая соответствует обширному циклоническому круговороту, охватывающему в этот период глубоководную часть моря, и содержит области менее соленой воды во вдольбереговой зоне моря. Наблюдается четко выраженная тенденция - при использовании схем с инвариантами большей степени (II и III расчеты - см. рис. 5,b,c) усиливается эффект подъема глубинных вод в центральной части моря и опускания менее соленых вод по периферии циклонического круговорота. Об этом свидетельствует структура изолинии 18,35 %о, которая служит условным маркером. Ее относительно гладкая структура демонстрирует, что повышение порядка инварианта не приводит к усилению вычислительного шума в модели, который проявляется в виде мелкомасштабной компоненты.

Рис. 5. Соленость на горизонте 30 м в I (а), II (b), III (c) экспериментах. Изолиния 18,35 %о выделена красным цветом

F i g. 5. Salinity at the 30-m horizon in the experiments I (a), II (b) and III (c). The 18.35 %o isoline is highlighted in red

В поле солености наибольшая разница имеет место в области втока рек и в прибосфорском районе. В центральной части моря наблюдаются ее отрицательные значения, что соответствует превышению солености во II и III экспериментах по сравнению с I расчетом, и положительные значения на свале глубин северо-западного шельфа и вдоль береговой области северной и западной частей моря (рис. 6). Экстремальные значения составляют приблизительно 1и2%о для обоих вариантов расчетов.

30 32 34 36 38 40 В.Д.

а

as,% о

28 30

b

Рис. 6. Разница между значениями солености на глубине 30 м в 1и11 (а), I и III (Ь) расчетах F i g. 6. Difference between the salinity values at the 30-m depth in the calculations I and II (a) and I and III (b)

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

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

Рис. 7. Вертикальная скорость в I (а)и III (Ь) расчетах (западная часть моря) на глубине 400 м F i g. 7. Vertical velocity at the 400-m depth in the calculations I (a) and III (Ь) (western part of the

Аналогичная картина наблюдается и в восточной части моря (рис. 8). Крупномасштабные черты в структуре вертикальной скорости во всех трех расчетах примерно одни и те же, но интенсивность мелкомасштабной изменчивости ниже во II и III расчетах.

Рис. 8. Вертикальная скорость в I (а) и III (Ь) расчетах (восточная область моря) на глубине 400 м

F i g. 8. Vertical velocity at the 400-m depth in the calculations I (a) and III (Ь) (eastern part of the 5. Заключение

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

TK (K > 1), SL (L > 1), р, (E + TT), (E + Dpe).

Легко увидеть, что выражения (31.1), (31.2) и (32) могут быть обобщены на случай произвольного вида дифференцируемого функционала Qt. k , зависящего от Ti -к ,и при условии выполнения уравнения неразрывности. При использовании схем, обладающих инвариантами TKjk, SLj k, где KrnL могут быть

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

При выводе разностного уравнения адвекции плотности от функционала Gt j k требуется его дифференцируемость от Тг j k и j k . Нетрудно обобщить

полученные результаты на случай, когда Gt j k зависит от R функций (формула

(42)). Предполагая, что каждая переменная функционала удовлетворяет уравнению адвекции типа (35) в поле скорости, удовлетворяющей выражению (28), получаемуравнениедля Gijk в дивергентом виде.

Привлекательной особенностью полученных аппроксимаций является возможность независимо от друг друга использовать схемы, обеспечивающие сохранение TK k, SLj k при K> 2, L > 2, и полной энергии при уравнении состояния, имеющем вид полинома от температуры и солености степени N > 1, M> 1.

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

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

1. Noether E. Invariante Variationsprobleme II Nachrichten von der Gesellschaft der Wissenschaften zuGöttingen. 1918. Vol. 2. P. 235-257.

2. Дородницын В. А. Конечно-разностный аналог теоремы Нетер II Доклады Академии наук. 1993. Т. 328, № 6. С. 678-682.

3. Самарский А. А., Мажукин В. И., Матус П. П. Инвариантные разностные схемы для дифференциальных уравнений с преобразованием независимых переменных II Доклады Академии наук. 1997. Т. 352, № 5. С. 602-605.

4. Дородницын В. А., Капцов Е. И. Инвариантные разностные схемы для обыкновенных дифференциальных уравнений второго порядка, обладающих симметриями. М. : ИПМ им. М. В. Келдыша, 2014. 42 с. (Препринты ИПМ им. М.В. Келдыша ; № 16). URL: http://library.keldysh.ru/preprint.asp?id=2014-16 (дата обращения: 30.09.2023).

5. Cheviakov А. F., Dorodnitsyn V. А., Kaptsov Е. I. Invariant Conservation Law-Preserving Discretizations of Linear and Nonlinear Wave Equations II Journal of Mathematical Physics. 2020. Vol. 61, iss. 8. 081504. doi:10.1063/5.0004372

6. Капцов Е. И. Численная реализация инвариантной схемы для одномерных уравнений мелкой воды в лагранжевых координатах. М. : ИПМ им. М. В. Келдыша, 2019. 28 с. (Препринты ИПМ им. М. В. Келдыша; № 108). doi:10.20948/prepr-2019-108

7. An energy-stable scheme for incompressible Navier-Stokes equations with periodically updated coefficient matrix/L. Lin [etal.] II Journal of Computational Physics. 2020. Vol. 418. 109624. doi:l 0.1016/j.jcp.2020.109624

8. Chen H., Sun S., Zhang T. Energy stability analysis of some fully discrete numerical schemes for incompressible Navier-Stokes equations on staggered grids II Journal of Scientific Computing. 2018. Vol. 75, iss. 1. P. 427-456. https://doi.org/10.1007/sl0915-017-0543-3

9. Lin L., Yang Zh., Dong S. Numerical approximation of incompressible Navier-Stokes equations based on an auxiliary energy variable II Journal of Computational Physics. 2019. Vol. 388. P. 122. doi:10.1016/j.jcp.2019.03.012

10. Головизнин В. М., Самарский А. А. Разностная аппроксимация конвективного переноса с пространственным расщеплением временной производной II Математическое моделирование. 1998. Т. 10, № 1. С. 86-100.

11. Афанасьев Н. А., Головизнин В. М., Соловьев А. В. Схема КАБАРЕ с улучшенными дисперсионными свойствами для систем линейных дифференциальных уравнений гиперболического типа II Вычислительные методы и программирование. 2021. Т. 22, вып. 1. С. 67-76. https://doi.org/10.26089/NumMet.v22rl05

12. Validation of the low dissipation computational algorithm CABARET-MFSH for multilayer hydrostatic flows with a free surface on the lock-release experiments / V. M. Goloviznin [et al.] II Journal of Computational Physics. 2022. Vol. 463. 111239. https://doi.Org/10.1016/j.jcp.2022.111239

13. Grooms I., Nadeau L.-P., Smith K. S. Mesoscale eddy energy locality in an idealized ocean model II Journal of Physical Oceanography. 2013. Vol. 43, iss. 9. P. 1911-1923. https://doi.Org/10.1175/JPO-D-13-036.l

14. Review on mechanical energy of ocean mesoscale eddies and associated energy sources and sinks /X. Shang [et al.] II Journal of Tropical Oceanography. 2013. Vol. 32, iss. 2. P. 24-36. doi:10.3969/j.issn,1009-5470.2013.02.003

15. Kjellsson J., Zanna L. The Impact of horizontal resolution on energy transfers in global ocean models II Fluids. 2017. Vol. 2, iss. 3. 45. https://doi.org/10.3390/fluids2030045

16. Seasonal variability of eddy kinetic energy in a global high-resolution ocean model / J. K. Rieck [et al.] II Geophysical Research Letters. 2015. Vol. 42, iss. 21. P. 9379-9386. doi:10.1002/2015GL066152

17. An estimate of the Lorenz energy cycle for the World Ocean based on the 1/10° STORM/NCEP simulation / J.-S. von Storch [et al.] II Journal ofPhysical Oceanography. 2012. Vol. 42, iss. 12. P. 2185-2205. https://doi.Org/10.1175/JPO-D-12-079.l

18. Yang Y., San LiangX. On the seasonal eddy variability in the Kuroshio extension II Journal of Physical Oceanography. 2018. Vol. 48, iss. 8. P. 1675-1689. doi:10.1175/JPO-D-18-0058.1

19. The eddy kinetic energy budget in the Red Sea / P. Zhan [et al.] II Journal of Geophysical Research: Oceans. 2016. Vol. 121, iss. 7. P. 4732-4747. doi:10.1002/2015JC011589

20. Stepanov D. V. Mesoscale eddies and baroclinic instability over the eastern Sakhalin shelf of the Sea of Okhotsk: a model-based analysis// Ocean Dynamics. 2018. Vol. 68, iss. 10. P. 13531370. doi:10.1007/sl0236-018-1192-2

21. Кубряков А. А., Станичный С. В. Синоптические вихри в Черном море по данным спутниковой альтиметрии II Океанология. 2015. Т. 55, № 1. С. 65-77. EDN TGWBIX. doi:10.7868/S0030157415010104

22. Menna M., Poulain P.-M. Geostrophic currents and kinetic energies in the Black Sea estimated from merged drifter and satellite altimetry data II Ocean Science. 2014. Vol. 10, iss. 2. P. 155165. https://doi.org/10.5194/os-10-155-2014

23. Суворов А. М., Шокурова И. Г. Годовая и междесятилетняя изменчивость доступной потенциальной энергии в Черном море II Морской гидрофизический журнал. 2004. № 2. С. 29-41.

24. Cessi P., Pinardi N., Lyubartsev ^.Energetics of semienclosed basins with two-layer flows at the strait II Journal of Physical Oceanography. 2014. Vol. 44, iss. 3. P. 967-979. https://doi.Org/10.1175/JPO-D-13-0129.l

25. Демышев С. Г. Энергетика климатической циркуляции Черного моря. Ч. I. Дискретные уравнения скорости изменения кинетической и потенциальной энергий II Метеорология и гидрология. 2004. № 9. С. 65-80. EDNPGCNXF.

26. Демышев С. Г. Аппроксимация силы плавучести в численной модели бароклинных течений океана II Известия РАН. Физика атмосферы и океана. 1998. Т. 34, № 3. С. 404-412.

27. Lorenz E. N. Available potential energy and the maintenance of the general circulation II Tellus. 1955. Vol. 7, iss. 2. P. 157-167. https://doi.org/10.3402/tellusa.v7i2.8796

28. Roquet F. Dynamical potential energy: a new approach to ocean energetics II Journal of Physical Oceanography. 2013. Vol. 43, iss. 2. P. 457-476. https://doi.Org/10.1175/JPO-D-12-098.l

29. Demyshev S. G., Dymova O. A. Numerical analysis of the Black Sea currents and mesoscale eddies in 2006 and 2011 II Ocean Dynamics. 2018. Vol. 68, iss. 10. P. 1335-1352. doi:10.1007/sl0236-018-1200-6

30. Гидрометеорология и гидрохимия морей СССР. Т. 4. Черное море. Вып. 1. Гидрометеорологические условия. СПб. : Гидрометеоиздат, 1991. 428 с.

31. Mellor G. L., Yamada T. Development of a turbulence closure model for geophysical fluid problems II Reviews of Geophysics. 1982. Vol. 20, iss. 4. P. 851-875. https://doi.org/10.1029/RG020i004p00851

Об авторе:

Демышев Сергей Германович, зав. отделом теории волн, главный научный сотрудник,

Морской гидрофизический институт РАН (299011,Россия, г. Севастополь, ул. Капитанская, 2),

доктор физико-математических наук, Scopus Author ID: 6603919865, SPIN-код: 1848-2350, Isti-

naResearcherID (IRID): 17369115, ResearcherID: C-1729-2016, ORCID ID: 0000-0002-54052282, demyshev@gmail.com