Конечно-разностная аппроксимация уравнения потенциальной завихренности для стратифицированной несжимаемой жидкости и пример его использования при расчете циркуляции Черного моря. Часть I. Дифференциально-разностное уравнение потенциальной завихренности идеальной жидкости Текст научной статьи по специальности «Физика»

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

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

УДК 551.465 EDN: BCKKBN

Конечно-разностная аппроксимация уравнения потенциальной завихренности для стратифицированной несжимаемой жидкости и пример его использования при расчете циркуляции Черного моря

Часть I. Дифференциально-разностное уравнение потенциальной завихренности идеальной жидкости

С. Г. Демышев

Морской гидрофизический институт РАН, Севастополь, Россия н [email protected]

Поступила в редакцию 09.06.2023; одобрена после рецензирования 25.07.2023; принята к публикации 18.01.2023.

Аннотация

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

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

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

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

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

Для цитирования: Демышев С. Г. Конечно-разностная аппроксимация уравнения потенциальной завихренности для стратифицированной несжимаемой жидкости и пример его использования при расчете циркуляции Черного моря. Часть I. Дифференциально-разностное уравнение потенциальной завихренности идеальной жидкости // Морской гидрофизический журнал. 2024. Т. 40, № 2. С. 165-179. EDN BCKKBN.

Original article

Finite-Difference Approximation of the Potential Vorticity Equation for a Stratified Incompressible Fluid and an Example of its Application when Modeling the Black Sea Circulation

Part I. Finite-Difference Equation of Potential Vorticity of Ideal Fluid

S. G. Demyshev

Marine Hydrophysical Institute, Russian Academy of Sciences, Sevastopol, Russia H [email protected]

Abstract

Purpose. The study is purposed at deriving the discrete equations of absolute vorticity and potential vorticity for a three-dimensional stratified incompressible fluid as an exact consequence of the finite-difference equations of sea dynamics in the field of a potential mass force in the adiabatic approximation provided that viscosity and diffusion are absent. The purpose also consists in analyzing both the features of two-dimensional projections of the absolute vorticty equation onto the coordinate planes and the three-dimensional equation of potential vorticity.

Methods and Results. In order to determine the discrete analogues of absolute and potential vorticity, introduced is the additional grid on which the finite-difference equations for the components both of absolute vorticty and potential vorticity are written down. Two-dimensional analogues of the three-dimensional equation of absolute vorticity on the planes (x, y), (y, z) and (x, z) are obtained; they possess the feature of preserving vorticity, energy and enstrophy (square of vorticiy). A discrete equation for potential vorticity of a stratified incompressible fluid is derived from the finite-difference system of three-dimensional equations of sea dynamics in the adiabatic approximation at the absence of viscosity and diffusion.

Conclusions. In the case of a linear equation of state, the discrete equations of absolute vorticity and potential vorticity which are the exact consequence of finite-difference formulation are obtained. The equation of potential vorticity is of a divergent form, and two-dimensional analogues of the absolute vorticity equation on the planes (x, y), (y, z) and (x, z) have two quadratic invariants that provide preservation of the average wave number.

Keywords: discrete equation, dynamics of sea, kinetic energy, vortex, potential vorticity, Ertel invariant

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

For citation: Demyshev, S.G., 2024. Finite-Difference Approximation of the Potential Vorticity Equation for a Stratified Incompressible Fluid and an Example of its Application when Modeling the Black Sea Circulation. Part I. Finite-Difference Equation of Potential Vorticity of Ideal Fluid. Physical Oceanography, 31(2), pp. 149-160.

Введение

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

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

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

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

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

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

1 Scott J. R., Chang S.-C. A new flux conserving Newton's method scheme for the two-dimensional, steady Navier-Stokes equations. NASA, 1993. 50 p. (NASA Technical Memorandum ; 106160). URL: https://ntrs.nasa.gov/api/citations/19930019437/downloads/19930019437.pdf (date of access: 29.03.2024).

2 Scott J. R. A new flux-conserving numerical scheme for the steady, incompressible Navier-Stokes equations. NASA, 1994. 49 p. (NASA Technical Memorandum ; 106520). URL: https://ntrs.nasa.gov/api/citations/19940024968/downloads/19940024968.pdf (date of access: 29.03.2024). МОРЖОЙ ГИДРОФИЗИЧЕСКИЙ ЖУРНАЛ том 40 № 2 2024 167

при сохранении массы и импульса движения. Для 2- и 3-мерных уравнений Навье - Стокса в приближении несжимаемости исследуется их дискретизация в формулировке сохранения энергии - импульса - углового момента [7]. Показано, что при кососимметричной линеаризации уравнений Навье - Стокса сохраняется энергия, при ньютоновской - импульс и угловой момент. На основе численных расчетов делается вывод, что линеаризация с двумя шагами Ньютона на каждом временном шаге обеспечивает сохранение всех трех параметров на больших временных интервалах.

В работе [8] получено обобщение схемы Аракавы - Лэмба [9] при дискретизации высокого порядка, но с точностью до аппроксимации по времени. Установлено, что для периодических задач с использованием операторов суммирования по частям высокого порядка выполняются свойства симметрии, законы сохранения энергии и энстрофии.

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

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

Уравнение потенциальной завихренности стратифицированной несжимаемой жидкости в дифференциальной постановке

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

^ + = + (1)

а Ро Ро

уи = 0 , (2)

^ + с11У(Ш) = 0, (3)

7ГТ

— +div(^tJ) = 0, (4)

dt

P =f(T,S). (5)

Введены обозначения: U = (u,v,w) - компоненты вектора скорости течения по осям (х, у, z). направленным на восток, север и вертикально вниз соответственно; g = (0, 0, g) - ускорение свободного падения ; (Г, S, Р, р) - температура, соленость, давление и плотность морской воды; ро = 1 г/см3 (в дальнейшем полагаем давление и плотность нормированными на ро); / = = (0,0, fz) - параметр Кориолиса, где fz = 2ю 8Шф; ю - угловая скорость

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

В уравнении (1) введены абсолютный вихрь скорости и кинетическая энергия движения:

\ =

0w ОМ 0W , z = 5М +

J Л Л ' J Л Л ' J Л Л </ ' V У

ду dz dz дх дх ду

2 2 2

J7 и + v + w Г1\

E = Р 0---. (7)

При z = 0 w = -qt, при z = H(x, y) w = 0. (8)

На боковых стенках условия непротекания:

для меридиональных и =0, для зональных участков границы v =0. (9)

Для температуры и солености ставятся условия отсутствие потоков.

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

при t = to (Т, S) = (Р, <S°), u = u°,v = v°,w = w°. (10)

Применим к уравнению (1) операцию rot, тогда для с; получаем

^ + Vx(lxU) = Vx(gp). (11)

at

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

3 Кочин Н. Е., Кибель И. А., Розе Н. В. Теоретическая гидромеханика. 6-е изд., испр. и доп. М. : Физматгиз, 1963. 1 т. 583 с.

МОРСКОЙ ГИДРОФИЗИЧЕСКИЙ ЖУРНАЛ том 40 № 2 2024 169

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

В силу определения абсолютного вихря (6) имеем

VI = 0. (12)

Тогда уравнение (11) можно переписать в следующем виде:

5 + и(^)-1(Уи) = Ух(£р). (13)

Если рассматривается двумерное движение в плоскостях (х, у), (х, ¿), (у, ¿), то для каждой проекции можно указать следующие законы сохранения. В плоскости (х, ¿), например, уравнения (2) и (13) переписываются в виде

ди дw дР дР дР — + — = 0, — + и— + V— = 0. дх дъ д дх дг

Тогда выполняются законы сохранения (- функция тока):

£ 15п=1 I(5')2' п=!1')2*п=0

и свойство антисимметрии:

з (у ,5у) = - з (5у у).

В случае выполнения двух квадратичных законов сохранения можно показать 4 [12], что среднее волновое число не зависит от времени. Указанное свойство означает, что не может быть систематической передачи энергии к движениям с высокими волновыми числами. При невыполнении этого свойства в численной модели может происходить ложный (вычислительный) каскад энергии к минимально разрешаемым на сетке (двухшаговым) волнам. Групповая скорость двухшаговых волн равна нулю, поэтому происходит накапливание энергии на этих масштабах.

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

Из уравнений (3)-(5) следует уравнение адвекции плотности

|^ + (Иу(рй) = 0. (14)

В силу равенства (12) выполняется

V [и(У1)] - V р (УЦ)] - V [V X (£ р )] = 0. (15)

4 Дымников В. П. Вычислительные методы в геофизической гидродинамике. М. : Отдел вычислительной математики АН СССР, 1984. 148 с.

170 МОРСКОЙ ГИДРОФИЗИЧЕСКИЙ ЖУРНАЛ том 40 № 2 2024

Из уравнения адвекции плотности (14) получим три уравнения для градиентов плотности. Для этого уравнение (14) продифференцируем последовательно по х, у, z. Полученные уравнения умножим на ^х , ^у 2, уравнение (13) - на соответствующие производные от р . После простых преобразований, учитывая соотношения (14) и (15), получаем уравнение для PV [2]

—+ (С/У)ю = 0, оо = 1Ур, (16)

д1

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

д_ а

| ю йО = 0. (17)

Равенство (17) означает, что потенциальная завихренность ю сохраняется в каждой частице при ее движении. Напомним, что уравнение (16) характеризует основное свойство стратифицированной несжимаемой жидкости при выполнении следующих условий: потенциальности массовой силы (сила тяжести) и несжимаемости морской воды, отсутствия вязкости и диффузии плотности в уравнениях движения [4]. Полученный Г. Эртелем результат обладает двумя замечательными свойствами. Если вместо плотности мы рассмотрим любой лагранжевый инвариант (например, 9). то со = УЭ также будет инвариантом [4]. Второе свойство заключается в том, что для несжимаемой жидкости в силу вида уравнения (16) инвариантом является любой

дифференцируемый по ю оператор ^. Умножим (16) на -= ^. Тогда по-

йт

лучим

йд

= 0, — |^(ю) й О = 0, йО = со^^

д О О

Л

Дискретное уравнение абсолютного вихря скорости движения и законы

сохранения

Аппроксимируем бассейн с неровным дном боксами, центрам которых соответствуют целочисленные значения индексов i,у, k (г = г,...,^, = ^,...,, k = 1, ..., К/,у), граням - /'+1/2, у+1/2, k+1/2. Горизонтальные размеры боксов (К, К) постоянные, по вертикали используется неравномерная аппроксимация (^ = гк+У2 - г^д, К+У 2 = %+1 - ) [13,14].

Разностные операторы имеют вид (для у, k - аналогично)

фх Л+У2^ +Ф,-1,2,лк , 5 Ъшп,^-Ф,-1,2,лк , У2 2 + 52ф ,

тг,^ ,к ^ хтг,^ ,к у^ ' х,у\г,],к х*г^,к утг,] ,к'

(18)

О

{Ф}"к =7Г^Ф^А, (ФГ = -11 £Ф,,j,ЛКНу, V = Х I.

" к ,, j V ,, j к=1 ,, j к=1

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

/-1/2,7+1/2, £-1/2

i,j,k-m/

,7+1/2,^-1/2

ж"

А-1/2

/+1 /2j+1 /2, £- 1 /2

¿+1/2,7,£-1/2

ь.ьк. _

f^Fi/tr

j+l/2:

/+1/2.

/+1/2|/-1/2,£ £+1/2

^/+1/2,7+1/2, А

hk

/ ^/-1/2,7,£+1/2

/

J&5 j-mj-m, £+1/2

7-1/2,7+1/2,£+Ш /

fWij;Fm- 1 7-1/2,£+1/2_

) /+1/2,7+1/2,£+1/2

f+1/2,7-I/2,£+1/2

Рис. 1. Распределение переменных в боксе (/,7, В вершинах бокса, обозначенных звездочкой, определена (ю), на его ребрах - компоненты абсолютной завихренности скорости

S/,j+1/2,£ > S/+1/2,j,£+1/2' S/+1/2, j+1/2,£

F i g. 1. Distribution of variables in box (i, j, £). At the box vertices indicated by an asterisk, PF (ю) is determined, and at its edges - the components of absolute vorticity

5i, j+1/2,£ >^/+1/2, j,£+1/2 '^i+1/2, j+1/2,£

С учетом обозначений (6), (7) и (18) выпишем дифференциально-разностные уравнения модели (1)-(5) (дифференциальные по времени) [9, 14]:

du+

dt

" Ч *« ■ ], .„2, „, + [ ' ]

M+1/2, j, £

= ^x (Ei+1/2, j/£ + p+1/2, j,£ )-

dv,

i, j+1/2, £ dt

[u, 5 ' ] /, j+1/2, £ -[ X ] /, j+1/2, £ =-S „ ( E/,, j+1/2, £ + P,, j+,/2, £ ),

(19.1)

(19.2)

г, у,к+1/2

Ж

^ 1„к+1/2 + 1„к+1/2 = ""5 г (+1/2 + +1/2 ) + £Рг,Лк+1/2, (19.3)

-л, у,к+1/2

]г, у,к+1/2

5 *иг, у,к + 5 уУ,, у,к + 5 М, у ,к = 0,

йТ

и ,к + 5* (и^Тук) + 5у к) + 5г ) = 0,

й

йБ,

г,У,к

Ж

+5* (и, у, А у, к)+5 у (V, у, А у,к) + 5г К у А у, к) = 0,

рг, 1 ,к = аТи,,к + РБг-,,,к .

(20) (21)

(22) (23)

В соответствии с разностными операторами (18) компоненты вихря скорости (рис. 1) имеют вид

£*/+1/2,к+1/2 = 5у (М!г

у+1/2,к+1/2 "у V'Кг,у'+1/2,к+1/^ "г Ч^у+^к+Ш

) - 5г (V,,

Л

пу+1/2,],.к+1/2 = 5г (и г+1/2, у, к+1/2 ) - 5х (+1/2,у,к+1/2 X

(24)

4г

= 5* (Уг+1/2, у+1/2,к ) - 5 у (и

) +

г+1/2,у+1/2,к х г+1/2,у+1/2,к у г+1/2,у+1/2,к у+1/2

Из аппроксимации (24) следует, что в точках / + 1/2, у + 1/2, k + 1/2 выполняется (аналог (12))

5 * 4 * + 5 у 4у + 5 г 4г = 0. (25)

Полагаем, что слагаемые в квадратных скобках в левой части уравнений (19.1)-(19.3) записываются в виде

[V, 4г]

1

г+1/2, у,к

-у

-у

IV,

г+1/2,у,к 4г+1/2,у,к

Vi+1, у,к 4 г+1/2, у ,к + Vi, у,к 4 г -1/2, у, к

[™, 4^ ]г+1/2,у,к

2wi+1/2,j,к 4f+1/2,^■,k

^г+1, у ,к 4 Г+1/2, у ,к + ^г ,у,к 4 ^у-1/2, у,к

[и, 4']г,у+1/2,к

2и

г, у +1/2,к 4г, j+1/2,k

* иг, у+1,к 4 г, у +1/2,к + иг, у,к 4 г, у -1/2,к

(26)

[^, 4* ]г,у+1/2,к

^г^+Ш,. 4i,j+1/2,k

г Wi, у+1,к 4 *, у+1/2,к + ^г, у,к ^ -1/2,к

[и,4"]

_ 1

г, у ,к+1/2 3

2и 4г

2иг, у ,к+1.2 4г, у ,к+1/2

'* иг, у ,к+14 о,.+1/2 + иг, у,к 4у у ,к-1/2

[' • 4' ]

1

i,j,k+1/2 3

2ч

-у

\ j,k+141+1/2, у,к + \ у,к 41-1/2, у,к

i,j,k+1/2 4г, у,к+1/2

Конечно-разностное представление (26) обеспечивает в случае двумерного бездивергентного движения в плоскостях (х, у), (х, z), (у, z) выполнение двух квадратичных законов сохранения: кинетической энергии и энстрофии (квадрата вихря) - и свойства антисимметрии [9, 14].

Учитывая выражение (24) и соотношения (26), выпишем уравнения для компонентов абсолютного вихря скорости - для 4* в точке /', у +1/2, к +1/2 , для 4у в точке г +1/2, у +1/2, к и для 4г в точке г +1/2, у +1/2, к (аналог уравнения (13)):

й4 *

+ 5 у (К4 * ]) + 5 г ([м,4 * ]) - 5 у ([и, 4у ]) - 5 г ([и, 4г ]) = g5 у р-

ш й4у

+ 5* ([и, 4у ]) + 5 г ([^,4у ]) - 5* ([V,4 * ]) - 5 г ([V,4г ]) = -в5* р ж

+ 5* ([и, 4г ]) + d у ([v,xг ]) - 5* ([м,4 * ]) - 5 у ([м,4у ]) = 0.

ш

(27.1)

(27.2)

(27.3)

Проведем анализ системы уравнений (27.1)-(27.3) с точки зрения выполнения законов сохранения.

Рассмотрим, например, двумерное бездивергентное движение в плоскости (х, z). В этом случае уравнение неразрывности (20) преобразуется к виду

5 и. ., + 5 ж ., = 0

* г,у ,к г г,у,к

и позволяет ввести функцию тока:

иг+1/2,к = 5гУу+1/2,к , Мг,к +1/2 = 5*Уу+1/2 .

Из уравнения неразрывности (28) следует, что

у у __Г у

Чм+1/2, к+1/2 _ 5г V,+1/2 , к+1/2 + 5 * ^+1/2 , к+1/2

Уравнение вихря упрощается и преобразуется [9] к виду

(28)

Ж4у+1/2, к+1/2 1 V"! т , 3 *=1

ж

4у+1/2, к+1/2 ) _ 0,

г+1/2, к+1/2' ^г+1/2, к+1/2

(29)

где

^ = 5* (V у+1/2, к+1/2 )5у (4у+1/2, к+1/2 ) - 5У (V г+1/2, к+1/2 )5* (4у+1/2, к+1/2 ),

з2 = 5*

1/2, к+1/2 5 у (4у+1/2, к+1/2 )

- 5,

^у+1/2, к+1/2 5* (4у+1/2, к+1/2 )

Зз = 5 у

4у+1/2, к+1/2 5* (Уу+1/2, к+1/2 )

У V

Л (-- 5.

г+1/2, к+1/2 "Д т г+1/2, к+1/2 V У V

5г'+1/2, к +1/2 5у (Уу+1/2, к+1/2 )

(30)

*

*

Выражение (30) означает, что уравнение (29) обладает двумя квадратичными инвариантами:

( У ^ у 2 Л

— и ик -П (5*^у+1/2,к+1/2 ) + (5 г ¥у+1/2,к+1/2 )

, / 4+1/2,к+1/2П*Пг = 0, , / ,

г, к

=0,

У

ЖУ (4 у )2ККк = 0

, / , (4г+1/2,к+1/2 ) П*Пг = 0.

йГТ.

(31)

Аппроксимация нелинейных слагаемых в уравнении (29) обеспечивает свойство антисимметрии:

3(^+1/2, к+1/2, 4у+1/2, к+1/2 ) = -3(4у+1/2, к+1/2, ¥'+1/2 , к+1/2 ). (32)

Аналогичная ситуация имеет место и для движения в плоскостях (у, z) и (х, у):

ЖУ 4* ЬЬк = ЖУ

Л / 4у+1/2,к+1/2КуК = ^ //

(5 у V**+1/2 ,к+1/2 ) + (5 г V *+1/2,к+1/2 )

2

^уКг

(33)

= Ж / (4у+1/2,к+1/2 )2Ку1К = 0,

ЙГ у,к

3 (¥*+1/2 , к+1/2,4*+1/2, к+1/2 ) = - (4*+1/2 , к+1/2 ¥*+1/2, к+1/2 )

Ж / 4>+1/2,У+1/2К*Ку = Ж /

^(5 * ¥г+1/2, у+1/2 )2 + (5 у ¥г+1/2, у+1/2 )2^

(34)

=^ / (4г+1/,у+1/2 )2 к*кУ=0,

3(¥г+1/2, у+1/2, 4г+1/2, у+1/2 ) = -3(4г+1/2, у+1/2, ¥+1/2, у+1/2 ).

Наряду со свойством сохранения полной энергии [5] соотношения (31)-(34) обеспечивают в системе уравнений (1)-(5) с краевыми (8)-(9) и начальными (10) условиями наличие двух квадратичных дискретных инвариантов.

Дискретное уравнение потенциальной завихренности стратифицированной несжимаемой жидкости

Рассмотрим свойства системы (27.1)-(27.3) в случае трехмерного движения. Запишем эту систему в следующем виде:

с*

-*

^ 1 г+1/2, 7+1/2,к+1/2 = §5 у (рг+1/2 , у+1/2,к+1/2 ), (35.1)

' Ъг+1/2, у+1/2,к+1/2

+ 1 г+1/2, у+1/2,к+1/2 = 55у (рг+1/2, у+1/2,к+1/2

Гу у

■эг+1/2, у+1/2,к+1/2

*г

йГ

г+1/2, у+1/2,к+1/2 "§5* (рг+1/2 , у+1/2,к+1/2 ), (35.2)

*

31+1/2, у+1/2,к+1/2 уг п /о- оч

^ + 1 г+1/2, у+1/2,к+1/2 = 0, (35.3)

где

(36)

у+1/2,к+1/2 = 5 у ([^4* ])г,у+1/2,к + 5 г ([^4* ])г, j+1/2,k -

- 5 у ([и ,4> ])г, у+ 1/2,к - 5 г ([и ,4 г ])г, у + 1/2,к ,

1у+1/2,у,к+1/2 = 5* ([и,4у ])г+1/2,у,к+1/2 + 5г ^ ])г+1/2,у,.+1/2 -

- 5* ([V ,4 ' ])г+1/2 , у,к+1/2 5 г ([V ,4 2 ])г+1/2 , у ,к+1/2

1г+1/2,у +1/2,к = 5* ([и,4 ])г+1/2,у +1/2,к + 5у ([у, 4 ])г+1/2,у +1/2,к --5* ([М, 4* ])г+1/2,у+1/2,к - 5у ([М, 4^ ])г+1/2,у+1/2,к.

Нетрудно убедиться, что в соответствии с соотношением (25), уравнениями (27.1)-(27.3) в точке / + 1/2, у + 1/2, k + 1/2 и определением (36) получаем

5* (1г+1/2, у+1/2,к+1/2 ) + 5у (1у+1/2, у+1/2,к+1/2 ) + 5г(1г+1/2, у+1/2,к+1/2 ) = 0. (37)

Заметим, что равенство (37) выполняется потому, что плотность линейно зависит от температуры и солености. Полагаем, что из уравнений (21)-(23) следует дискретное уравнение для плотности, которое имеет вид

^ + 5* (и, у ,к рг, у.к ) + 5 у (V, у,. р, у,к ) + 5г К, у.к р, у.к ) = 0. (38)

Запишем уравнение (38) в точках (г, у +1/2, к +1 / 2), (г +1/2, у, к +1 / 2), (г +1/2, у +1/2, к) соответственно:

-уг

й рг, у+1/2,к+1/2 г. (-у\ - (-у\

^ + 5* (иг, j+1/2,k+1/2рг,у +1/2,к+1/2 ) + 5у (Vг■,у+1/2,к+1/2рг,у +1/2,к+1/2 ) + (39]^

-уг

+5г у+1/2,к+1/2рг,у+1/2,к+1/2 ) = 0,

-*г

^ + 5* (иг+1/2,у,к+1/2рг+1/2,у,к+1/2 ) + 5у (Vi+1/2,j,k+1/2рг+1/2,у,к +1/2 ) + ^9 2)

-*г

+5г (Мг+1/2,у,к+1/2рг+1/2,у,к+1/2 ) = 0

А-*У

й рг+1/2,у+1/2,к ~ --*у- - --*у-

■ + 5* (иг+1/2,у+1/2,крг+1/2,у+1/2,к ) + 5у (У;+1/2,у+1/2,крг+1/2,у+1/2,к ) +

^^ г + 1/2, у + 1/ 2,к I г + 1/2, у + 1/2,к ^ ух г + 1/2, у + 1/2,Л1 г + 1/2, у + 1/2,к ^ (39 3) -*у

+5г (н>г+1/2,у+1/2,крг+1/2,у+1/2,к ) = 0.

Продифференцируем в разностном смысле уравнение (39.1) по х, (39.2) по у, (39.3) по z. В результате получаем

г

г

й0 х (р,+1/2„/+1/2,к+1/2) - ( „х МП П --+ Ох (Я1+1/2,/+1/2,к+1/2 ) - 0 I40.1)

~Х2

й0у(р1+1/2,/+1/2,к+1/2) + О(Ку ) - 0 (402)

^ + Оу (Я1+1/2,/+1/2,к+1/2) — 0, (402)

й0 г (Р1+1/2,/+1/2,к+1/2) ,г / П1 \ л /"ЛПОЧ -2/-+ 0г (Я1+1/2,/+1/2,к+1/2 ) — 0, (40.3)

где обозначения Щ/+1/2,к+1/2> Я1+1/2,/,к+1/2, Я,+1/2,/+1/2,к очевидны.

Отметим свойства введенных функций в уравнениях (39.1)-(39.3) и (40.1)-(40.3). Они имеют вид

-х -у -г

~уг —хг —ху -хуг

Р1,/+1/2,к+1/2 — Рг +1/2,/,к+1/2 — РК+1/2,/+1/2,к — Рг+1/2,/+1/2,к+1/2 , (41)

-х -у -г -хуг

Яг,/+1/2,к+1/2 — Я1+1/2,/,к+1/2 — Я1+1/2,/+1/2,к — Ох (М1+1/2,/+1/2,к+1/2Р1+1/2,/+1/2,к+1/2 ) +

-хуг

+0у (У1+1/2,/+1/2,к+1/2Рг+1/2,/+1/2,к+1/2 ) + 0г +1/2,/+1/2,к+1/2Р1+1/2,/+1/2,к+1/2 )

Полагаем, что разностный аналог PV определяется в точке г +1/2, / +1/2, к +1/2 и записывается следующим образом: ъм12 .+1/2¿+1/2 —

^¿+1/2,/+1/2,к+1/2 0 х (Р 1+1/2,/+1/2,к+1/2 ) + ^¿+1/2,/+1/2,к+1/2 0 у (Р 1+1/2,/+1/2,к+1/2 ) + _^ (42)

+^1+1/2,/+1/2,к+1/2 0 г (Р 1+1/2,/+1/2,к+1/2 ).

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

Для того чтобы получить разностный аналог теоремы Эртеля (15), умножим уравнения (35.1)-(35.3) на

—уг —хг —ху

0х (Р 1+1/2,/+1/2,к+1/2), 0у (Р 1+1/2,/+1/2,к+1/2 ), 0г (р1+1/2,/+1/2,к+1/2 X а систему (40.1)-(40.3)

~2

на ^+1/2,/+1/2,к+1/2 ^+1/2,/+1/2,к+1/2 ^+1/2,/+1/2,к+1/2 соответственно. В результате в точке i + 1/2, у + 1/2, k + 1/2 получим

-х — уг -у —хг -г —ху -х -у -г

^ + Тх 5х (Ру) + Ту 5у (Р ) + Тг 5г (Р ) + 0хЯх + 5уЯу + 0гЯг — 0.

ш

Так как г удовлетворяют соотношению (25), Тх, Ту, Тг - (37), учи-

тывая равенства (41), в точке i + 1/2, у + 1/2, k + 1/2 получаем дифференциально-разностное уравнение потенциальной завихренности (дифференциальной по времени) стратифицированной несжимаемой жидкости в дивергентом виде

— + 5(Тхруг + Ях) + 5 (Турхг + уЯу) + 8,(Тгр^ + ^ гЯг) — 0. (43) йг

X

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

Вид нелинейных слагаемых в уравнении (43), строго говоря, не соответствует их дифференциальному аналогу (15). Интеграл по области от — ра-

dt

вен нулю при выполнении следующих соотношений на границах:

Y*j+1/2,k+1/2 = 0, Rlj+1/2,k+1/2 = 0 при i = h, i = iN , Tf+1/2, j,k+1/2 = 0 R+1/2, j,k+1/2 = 0 при j = j, j = Jm , (44)

Y+1/2,j+1/2,k = 0, R;+1/2,j+1/2,k = 0 при k = 1/2, k = K, J + 1/2 ,

Из соотношений (44) следует, что необходимы дополнительные краевые условия, которые отсутствуют в исходной постановке.

Заключение

Выписан дифференциально-разностный аналог абсолютного вихря скорости для системы уравнений идеальной жидкости без квазистатического приближения. Проекции этого уравнения на двумерные подпространства (x, y), (y, z), (x, z) сохраняют энергию, вихрь, энстрофию и обладают свойством антисимметрии. Для того чтобы получить известную схему Аракавы - Лэмба для уравнений мелкой воды, необходимо выписать отличную от (26) исходную разностную систему уравнений для горизонтальных скоростей.

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

СПИСОК ЛИТЕРАТУРЫ 1. Rossby C.-G. Planetary flow patterns in the atmosphere // Quarterly Journal of the Royal Meteorological Society. 1940. Vol. 66, iss. S1. P. 68-87. https://doi.org/10.1002/j.1477-870X.1940.tb00130.x

2. Ertel H. Ein neuer hydrodynamischer Wirbelsatz // Meteorologische Zeitschrift. 1942. Vol. 59, no. 9. S. 277-281.

3. Hoskins B. J., McIntyre M. E., Robertson A. W. On the use and significance of isentropic potential vorticity maps // Quarterly Journal of the Royal Meteorological Society. 1985. Vol. 111, iss. 470. P. 877-946. https://doi.org/10.1002/qj.49711147002

4. Жмур В. В., Новоселова Е. В., Белоненко Т. В. Потенциальная завихренность в океане: подходы Эртеля и Россби с оценками для Лофотенского вихря // Известия Российской академии наук. Физика атмосферы и океана. 2021. Т. 57, № 6. С. 721-732. EDN SRKASA. https://doi.org/10.31857/S0002351521050151

5. Капцов Е. И. Численная реализация инвариантной схемы для одномерных уравнений мелкой воды в лагранжевых координатах // Препринты ИПМ им. М. В. Келдыша. 2019. № 108. С. 1-28. EDN DHLHXH. https://doi.org/10.20948/prepr-2019-108

6. Bihlo A., Popovych R. O. Invariant discretization schemes for the shallow-water equations // SIAM Journal on Scientific Computing. 2012. Vol. 34, no. 6. P. B810-B839. https://doi.org/10.1137/120861187

7. Efficient discretizations for the EMAC formulation of the incompressible Navier-Stokes equations / S. Charnyi [et al.] // Applied Numerical Mathematics. 2019. Vol. 141. P. 220-233. https://doi.org/10.1016/j .apnum.2018.11.013

8. Sorgentone C., La Cognata C., Nordstrom J. A new high order energy and enstrophy conserving Arakawa-like Jacobian differential operator // Journal of Computational Physics. 2015. Vol. 301. P. 167-177. https://doi.org/10.1016/j.jcp.2015.08.028

9. Arakawa A., Lamb V. R. A potential enstrophy and energy conserving scheme for the shallow water equations // Monthly Weather Review. 1981. Vol. 109, iss. 1. P. 18-36. https://doi.org/10.1175/1520-0493(1981)109<0018:APEAEC>2.0.C0;2

10. Salmon R. A general method for conserving quantities related to potential vorticity in numerical models // Nonlinearity. 2005. Vol. 18, no. 5. P. R1-R16. https://doi.org/10.1088/0951-7715/18/5/R01

11. Sugibuchi Y., Matsuo T., Sato S. Constructing invariant-preserving numerical schemes based on Poisson and Nambu brackets // JSIAM Letters. 2018. Vol. 10. P. 53-56. https://doi.org/10.14495/jsiaml.10.53

12. Arakawa A. Computational design for long-term numerical integration of the equations of fluid motion: Two-dimensional incompressible flow. Part I // Journal of Computational Physics. 1966. Vol. 1, iss. 1. P. 119-143. https://doi.org/10.1016/0021-9991 (66)90015-5

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

14. Демышев С. Г. Численные эксперименты по сопоставлению двух конечно-разностных схем для уравнений движения в дискретной модели гидродинамики Черного моря // Морской гидрофизический журнал. 2005. № 5. C. 47-59. EDN YURAEH.

15. Демышев С. Г. Нелинейные инварианты дискретной системы уравнений динамики моря в квазистатическом приближении // Морской гидрофизический журнал. 2023. Т. 39, № 5. С. 557-583. EDN JWSUUM.

Об авторе:

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

ФГБУН ФИЦ МГИ (299011, Россия, г. Севастополь, ул. Капитанская, 2), доктор физико-

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

IstinaResearcherID (ШГО): 17369115, ResearcherID: C-1729-2016, ORCID ID: 0000-00025405-2282, [email protected]