Эффективный метод решения модифицированной задачи Орра - Зоммерфельда для анализа неустойчивости течений в Арктическом бассейне Текст научной статьи по специальности «Математика»

УДК: 551.465, 519.624 MSC2010: 34L16

ЭФФЕКТИВНЫЙ МЕТОД РЕШЕНИЯ МОДИФИЦИРОВАННОЙ

ЗАДАЧИ ОРРА — ЗОММЕРФЕЛЬДА ДЛЯ АНАЛИЗА НЕУСТОЙЧИВОСТИ ТЕЧЕНИЙ В АРКТИЧЕСКОМ БАССЕЙНЕ

© С. Л. Скороходов, Н. П. Кузьмина

ФИЦ ИУ РАН, ИО РАН ул. Вавилова, д. 40, Нлхимовский проспект, д. 36 Москва, Российская Федерация e-mail: [email protected], [email protected]

Accurate method of solving modified Orr-Sommerfeld problem for

analysis OF the oceanic currents instability in the arctic basin.

Skorokhodov S. L., Kuzmina N. P.

Abstract. The efficient numerical method for solving the modified Orr — Sommerfeld problem has been elaborated. The govering equation defines the operator pencil of polynomial type with spectral parameter "c " entering the equation and boundary condition. This equation describes long-wave stable and unstable perturbations of geostrophic currents with linear vertical shear. The model includes vertical density diffusion; it is used to study generation of large scale intrusions in the Arctic Ocean.

The method for evaluation of eigenfunctions and eigenvalues is based on using power expansions at boundary and internal points of the layer and smooth matching of that expansions. The equation for Wronskian of 4 basic solutions, W(c) = 0, enables us to compute the discrete spectrum of the problem. The numerical analysis of the spectrum for odd eigenfunctions shows, that the first eigenvalue ci corresponds to unstable currents, and all the rest eigenvalues correspond to stable currents.

The asymptotics of eigenvalue ci for small value of parameter of the problem was studied and the main term of the asymptotics was found. Numerical data for eigenvalues cn have been obtained for a wide range of the Peclet number (modified Reinolds number), R e (0,106), with high accuracy 10-20.

The results obtained confirm and supplement early published analytical considerations justifying the conclusion that geostrophic current can be unstable due to vertical diffusion of density.

Keywords: spectral problem, Orr — Sommerfeld equation, eigenvalue, eigenfunction, unstable flow

Введение

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

I + иМ£) % - Ц-М| = , (!)

где р(х,у, г,Ь) — нормированное на отсчетную плотность отклонение давления от среднего состояния, ось х направлена на восток, ось у — на север, ось г — вверх, Ь — время, и (г) — зональная компонента скорости геострофического течения, К — коэффициент вертикальной диффузии. Течение рассматривается в слое г Е [—Но, Но] толщиной 2Н0 с боковыми границами у = 0 и у = Ь. Возмущения плотности д, зональной и и меридиональной V компонент скорости выражаются через возмущение давления в виде

1 др 1 др 1 др д д дг , и / ду , V / дх , где д — ускорение свободного падения, / = 2 П вт^ — параметр Кориолиса.

Краевыми условиями для уравнения (1) являются отсутствие протекания на верхней и нижней границах слоя,

I + и (г) I )дР — и-(г) | = кЦ? , г = ±Н , (3)

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

1 д? 1 д?

0. (4)

1 dp f dx

1 др у=о / дх у=Ь

Используя для зональной скорости и (г) распределение, соответствующее течению Пуазейля и (г) = в(Н^ — г2), вводя безразмерные переменные

- х - У г г . вНо2

х = ь, у = Ь, * = Но, Ь = ,

где вН,2 — характерный масштаб скорости течения, перепишем задачу (1)-(4) относительно безразмерных переменных, опуская крышки:

I + (1 — 9? + 2дх = К д4?, (5)

дЬ дх дг2 дх К дг4

д . (л 2ч 9\9р . 0 дР 1 д 3Р ,,

дЪ + (1 - + 2*д* = ' г = ±1' (6)

= 0, (7)

У=1

др др дх у=о дх

где Я = вИ^ / (КЬ) — число Пекле (аналог числа Рейнольдса).

Решение р(х,у, г,г), с учетом условий (7), будем искать в виде бегущей волны вдоль координаты х:

р(х, у, г, г) = Ф(г) егк(х-с1) яп(пу), (8)

где к — произвольное волновое число, а константа с — неизвестная скорость волны. Подстановка (8) в уравнения (5) и (6) приводит к задаче для искомой функции Ф(г):

(1 - г 2 - с) Ф ''(г) + 2Ф(г) = -^Ф(/у>(*), (9)

гкЯ

(1 - г2 - с) Ф '(г) + 2 г Ф(г) = -1- Ф "'(г), г = ±1. (10)

гкЯ

Уравнение (9) имеет четвертый порядок и помимо двух краевых условий (10) требует еще двух граничных условий. Различные варианты таких условий рассмотрены в [4]; в данной работе мы остановимся на равенстве нулю потоков плотности (плавучести) на верхней г =1 и нижней г = -1 границах слоя:

Ф ''(г) = 0 , г = ±1. (11)

Задача (9)—(11) является спектральной задачей для поиска собственных функций (далее СФ) Фп(г) и соответствующих им собственных значений (далее СЗ) сп.

Уравнение (9) соответствует модифицированному уравнению Орра-Зоммер-фельда (см. [5], [6], [7]) и определяет операторный пучок полиномиального типа (см. [8], [9]), причем спектральный параметр с входит также и в краевое условие (10).

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

1. Метод решения

Вычисление СФ и соответствующих СЗ основано на методе, построенном в [10] для расчета СЗ задачи Орра-Зоммерфельда.

Для нахождения нечетных СФ и соответствующих СЗ сп при любых значениях параметра кЯ представим регулярную функцию Ф(г) в виде двух разложений в точках г = -1 и г = 0;ее нечетное продолжение на отрезок [0,1] даст искомое решение

при z Е [-1, 1].

Ф(г) = ^ а„(1 + г)п , Ф(г) = ^ Ьпг2п+1. (12)

п=0 п=0

Оба разложения (12) сходятся при всех |г| < то. Подставляя эти ряды в уравнение (9), получаем рекуррентные соотношения для коэффициентов ап и Ьп:

¿кК [2п(п + 1)ап+1 + (2 — п(п — 1))ап — с(п + 1)(п + 2)ап+2]

bn+2

ап+4 (п + 1)(п + 2)(п + 3)(п + 4)

гкК[(п + 1)(2п + 3)(1 — с)Ьп+1 + (1 — п(2п + 1))Ьп] (п + 1)(2п + 3)(2п + 4)(2п + 5) Для учета краевого условия (11) в точке г = —1 положим

а2 = 0,

а для удовлетворения условия (10) используем соотношение

n

0,1,

«3

6

(cai + 2ao).

(13)

(14)

(15)

Теперь построим 2 линейно-независимых решения Ф1(г) и Ф2(г), тождественно удовлетворяющих краевым условиям (10) и (11). Для этого зададим коэффициенты

(1) (1) (2) (2) a0 , ai и ад , a такими

o1

(1)

0, а(11) = 1; а02) = 1, а^ = 0,

а все последующие а^ и ап2) вычислим по (14), (15) и (13) для ап.

Далее построим 2 независимых решения Ф3(г) и Ф4(г), задав коэффициенты Ь,3",

Т (3) ,(4) ,(4)

Ь1 и Ь0 , Ь1 такими

Ь03) = 0, Ь13) = 1; Ь04) = 1, Ь14) =0, (17)

а все последующие Ы3 и Ь^ вычислим по (13) для Ьп.

В окрестности точки г = — 1 общее решение уравнения (9) с краевыми условиями (10), (11) представимо в виде комбинации

Ф{-1}(г) = ^Ф^г) + 4Ф2(г) , (18)

а в окрестности г = 0 общее нечетное решение — в виде комбинации

Ф{о}(г)= ^(г) + ^(г), (19)

где константы ... , Е С — произвольные, а Ф1(г),... , Ф4(г) — разложения (12) с соответствующими коэффициентами.

(2)

(2)

(16)

(3)

Для гладкой сшивки общих решений (18) и (19) выберем произвольную точку г* Е (-1, 0) и потребуем в ней совпадения функций и ее первых трех производных:

^Ф^*) + ^Ф2т)(г*) = 4Ф3т)(г*) + 4Ф4т)Ы, т = 0,1, 2, 3. (20)

Все последующие производные решений Ф{-1}(г) и Ф{0}(г) в точке г* совпадут в силу уравнения (9) четвертого порядка.

Решая систему (20) относительно коэффициентов ¿1,^2,^3,^4, находим искомую функцию Ф(г) в виде (18) и (19), удобном в окрестности точек г = -1 и г = 0, соответственно.

Для нетривиального решения линейной системы (20) необходимо равенство нулю определителя Вронского Ш(Ф1, Ф2, Ф3, Ф4; с, г* ) :

Ш(Фь Ф2, Фз, Ф4; с; г*) =0 , (21)

причем зависимость функций Фр (г) от спектрального параметра с здесь указана в форме Ш(...; с;...). Уравнение (21) является основным для поиска СЗ с, его решение будем строить с помощью итерационного метода Ньютона,

с(п+1) = с(п) - Ш(...; с(п); г*) [дШ1"; с— 1 -1 , „ = 0,1,... (22)

* дс

Начальная итерация с(0) выбиралась с помощью метода непрерывного продолжения по параметру кЯ задачи (9)—(11) и с использованием обобщенного принципа аргумента:

N = Л- / ШИ ¿с. £ с. =- /сЩИ ¿с. (23)

2пг / Ш(...; с) ' ^ р 2пг / Ш(...; с) ' 1 ;

Б Р=1 Б

где N — число нулей функции Ш(...; с) внутри некоторого контура Д, а ср — ее нули внутри Д. Отметим здесь, что полюсов функция Ш(Ф1, Ф2, Ф3, Ф4; с; г*) не имеет.

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

Разработанный метод позволил вычислить дискретный спектр СЗ ст задачи (9)-(11) для нечетных СФ в широком диапазоне параметров кЯ Е (0,106), причем относительная погрешность расчетов не превышала 10-20. Точность расчетов контролировалась с помощью увеличения мантиссы вычислений, изменения длины рядов Тейлора в разложениях (12) и вариации точки сшивки г* Е (-1, 0).

2. Численные результаты

Приведем данные расчетов СЗ для нечетных СФ при некоторых параметрах кК.

На Рис. 1 представлены в комплексной плоскости "с" первые 33 СЗ сп, п = 1,..., 33, для кК = 104. Последующие СЗ имеют 1т(сп) < —1, причем И,е(сп) ^ 2 при п ^ то.

Особую значимость имеют СЗ, для которых 1т(сп) > 0. Из представления (8) решения р(х,у,г,Ь) легко видно, что искомые решения с этими СЗ оказываются неустойчивыми с увеличением времени Ь > 0 и приводящими к росту малых возмущений. На Рис. 1 таким СЗ является с1 = 1.00707106781 + г 0.00707106781.

На Рис. 2 представлены первые 11 СЗ сп, п = 1,..., 11, для кК = 103. Последующие СЗ имеют 1т(сп) < —1 и Ие(сп) ^ | при п ^ то. Неустойчивой СФ здесь соответствует с1 = 1.02236017527 + г 0.02236103467.

Проведенный численный анализ выявил важное свойство поведения первого СЗ с1 при уменьшении параметра кК, чему соответствует увеличение коэффициента вертикальной диффузии в исследуемом слое. Приведем некоторые значения с1 = с1(кК) при кК е (0.01,1).

kR Re(ci)

1 1.33191395963352088677

0.5 1.33297527533996627539

0.2 1.33327589778426305712

0.05 1.33332974197241993087

0.01 1.33333318967469796631

Im(c1)

0.01678348554345884397 0.00844694134830805473 0.00338504503318573967 0.00084654211204427524 0.00016931201943078500

Эти данные показывают, что с1 ^ | при кК ^ 0. Докажем это.

3. Асимптотика c1 при kR ^ 0

Перепишем задачу (9)-(11) в более удобной форме

Ф""^) = ikR[(1 - z2 - с)Ф"(г) + 2Ф(г) ] , (24)

Ф"^) = ikR[(1 - z2 - c)Ф/(z)+2zФ(z)] , Ф//(z) = 0 , z = ±1. (25)

Представим нечетное решение Ф^) и СЗ c1? стремящееся при kR ^ 0 к конечному значению c10), в виде асимптотического разложения по степеням малого параметра kR:

Ф^; z) = ^0(z) + kR^1(z) + ..., c1(kR) = c10) + kRc(11) + ... (26)

Подстановка (26) в (24), (25) приводит к цепочке краевых задач для нечетных функций <^п(г), п = 0,1,...

0. Для п = 0 получаем:

) = 0,

^>"(±1) = 0, ^'(±1) = 0. Нечетным решением (27) является

^o(z) = AlZ

с произвольной константой A1.

(27.1)

(27.2)

(28)

1. Для n = 1 имеем:

^i'"(z) = i (1 - z2 - c1O))^0/(z) + 2^o(z)

^1"(±1) = i i-c1O)^0(±1) ± 2po(±1)l , 0i(±1) = 0.

(29.1)

(29.2)

С учетом вида ^0(z) из (28), находим нечетное решение (z) из уравнения (29.1):

<^(z) = + Bz + B3z3 ,

где B1 и B3 — искомые константы. Краевые условия (29.2) приводят к системе

iAi + 6Вз = iAi [-c10) + 2 , ^ + 6B

нетривиальное решение которой существует лишь при

iA1

,(0)

4 3

(30)

Запишем окончательный вид полученной функции ^i(z):

^Н2^z5 + Biz - ^z3 ,

где константа произвольна. Из представления (26) и соотношения (30) следует равенство

4

lim cT(kR) = - , fcfi^o 3

что и завершает доказательство, а результаты численных расчетов предыдущего раздела это иллюстрируют.

Остается отметить, что все остальные СЗ cn(kR), n = 2, 3,..., при kR ^ 0 стремятся к то, приближаясь к прямой Re(cn) = 2/3, а Im(cn) ^ — то.

0

1

Рис. 1. Собственные значения cn при kR = 104.

Заключение

На основе высокоточного численного метода построено решение модифицированной спектральной задачи Орра-Зоммерфельда с граничными условиями, типичными для океана. Исследуемое уравнение определяет операторный пучок полиномиального типа, причем спектральный параметр "с" входит также и в граничное условие. Модельное уравнение описывает длинноволновые устойчивые и неустойчивые возмущения геострофического течения с линейным вертикальным сдвигом скорости с учетом вертикальной диффузии плотности и применяется для исследования образования крупномасштабных интрузий в Арктическом бассейне.

Метод нахождения собственных функций и собственных значений " с" основан на использовании разложений решения в граничной и внутренней точках слоя и последующей гладкой сшивке этих разложений. Условие равенства нулю вронскиана независимых решений позволяет вычислить дискретный спектр "с" задачи. Численно показано, что первое из собственных значений ci соответствует неустойчивому течению, а все остальные — устойчивому. Получен главный член асимптотики с1 при

1 1 1 ЕЮЕМ УАШЕЗ 1 в 1 1 зоштюыз

- к = 1, Р = 10^ -

- о о

- О О О о о -

О

- О -

- О , о 1 1 1

О 0.2 0.4 0 6 0.0 1 1.2 1.4

Рис. 2. Собственные значения сп при кЯ = 103.

стремлении к нулю параметра задачи. Численные расчеты проведены для широкого диапазона числа Пекле (модифицированное число Рейнольдса), Я Е (0,106), и выполнены с точностью 20 верных дес. знач. цифр; они полностью соответствуют проведенному асимптотическому анализу С1.

Полученные результаты подтверждают и дополняют аналитические рассмотрения в [3], [4] и дают основание для вывода о том, что геострофическое течение в ограниченном по вертикали океаническом слое может быть неустойчиво из-за влияния диффузии плотности, причем неустойчивые моды не являются неустойчивостью критического слоя.

Работа выполнена при финансовой поддержке РФФИ (проекты 16 — 01 — 00781 и 15 — 05 — 01479).

Описок литературы

1. Pedlosky, J. (1992) Geophysical Fluid Dynamics. New York: Springer.

2. Cushman-Roisin, B. (1994) Introduction to the Geophysical Fluid Dynamics. New Jersey: Prentice Hall, Englewood Cliffs.

3. Кузьмина, Н. П. Об одной гипотезе образования крупномасштабных интрузий в Арктическом бассейне / Н. П. Кузьмина // Фундаментальная и прикладная гидрофизика / Т. 9. — № 2, 2016. — C. 15-26.

KUZMINA, N. P. (2016) A possibility of large scale intrusions generation in the Arctic Ocean under stable-stable stratification: An analytical consideration. Ocean Sci. Discussion, doi: 10.5194/os. 15.

4. Kuzmina, N. P. A possibility of large scale intrusions generation in the Arctic Ocean under stable-stable stratification: An analytical consideration // Ocean Sci. Discussion, doi: 10.5194/os-2016-15, 2016.

5. Drazin, R. G., Reid, W. H. (1981) Hydrodynamic Stability. Cambridge: Cambridge Univ. Press.

6. Reddy, S. C., Schmid, P. J., Heningson, D. S. Pseudospectra of the Orr-Sommerfeld operator // SIAM J. Appl. Math. — 1993. — V. 53. — No 1. p. 15-47.

7. Шкаликов, А. А. Как определить оператор Орра-Зоммерфельда? / А. А. Шкаликов // Фундаментальная и прикладная гидрофизика / Т. 53. — № 4, 1998. — C. 36-43.

SHKALIKOV, A. A. (1998) How to define the Orr-Sommerfeld opetator?. Vestnik Moskovskogo Universiteta. Seriya 1. Matematika. Mekhanika. 53 (4). p. 36-43.

8. Гохберг, И. Ц. Введение в теорию линейных несамосопряженных операторов / И. Ц. Гохберг, М. Г. Крейн. — M.: Наука, 1965. — 448 c.

GOHBERG, I. and KREIN, M. (1965) Introduction to the theory of linear nonselfajoint operator. Moscow: Nauka. 448 p.

9. Копачевский, Н. Д. Спектральная теория операторных пучков. Курс лекций / Н. Д. Копачевский. — Симферополь: ООО "ФОРМА", 2009. — 128 c.

KOPACHEVSKII, N. D. (2009) Spectral theory of operator pencils. Lectures on Mathematics. Simferopol: OOO "FORMA". 128 p.

10. Скороходов, С. Л. Численный анализ спектра задачи Орра-Зоммерфельда / С. Л. Скороходов // Журнал вычисл. матем. и матем. физ / Т. 47. — № 10, 2007. — C. 1672-1691.

SKOROKHODOV, S. L. (2007) Numerical analysis of the spectrum of the Orr-Sommerfeld problem. Comput. Math. and Math. Phys. 47 (10). p. 1603-1621.