Научная статья на тему 'Регуляризация метода прямых в условиях машинной точности с примерами в гидродинамике со свободной поверхностью'

Регуляризация метода прямых в условиях машинной точности с примерами в гидродинамике со свободной поверхностью Текст научной статьи по специальности «Математика»

CC BY
64
25
i Надоели баннеры? Вы всегда можете отключить рекламу.
Область наук
Ключевые слова
ГИДРОДИНАМИКА СО СВОБОДНОЙ ПОВЕРХНОСТЬЮ / ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ / РЕГУЛЯРИЗАЦИЯ / HYDRODYNAMICS WITH FREE SURFACE / NUMERIC SIMULATIONS / REGULARIZATION

Аннотация научной статьи по математике, автор научной работы — Шамин Р. В.

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

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Regularization of the method of lines with a high (computer) precisionwith examples on hydrodynamics of ideal fluid with a free boundary

We consider a realization of the lines method for numerical solution of evolution equations under the conditions of a high (computer defined) precision. We develop an effective algorithm of the lines method, which allows us to perform numerical calculations even in the case when the standard method turns to be unstable. We consider the examples of a successful application of the proposed approach for solutions of important hydrodynamics problems.

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

Вычислительные технологии

Том 13, № 5, 2008

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

Р. В. Шлмин

Институт океанологии им. П.П. Ширшова РАН, Москва, Россия

e-mail: roman@shamin.ru

We consider a realization of the lines method for numerical solution of evolution equations under the conditions of a high (computer defined) precision. We develop an effective algorithm of the lines method, which allows us to perform numerical calculations even in the case when the standard method turns to be unstable. We consider the examples of a successful application of the proposed approach for solutions of important hydrodynamics problems.

Введение

Погрешности округлений при вычислениях — принципиальная проблема вычислительной математики. Одним из возможных способов преодоления трудностей вычислений в условиях машинной точности является применение методов интервального анализа (см., например, [1-3]). Однако использование интервального исчисления требует специального программирования численных методов, снижает скорость вычислений и требует увеличения объема оперативной памяти. В настоящей работе предложен простой и эффективный способ “борьбы” с особенностями машинной точности. Этот способ (см. разд. 2) — один из вариантов спектрального фильтра. В то же время предлагаемый метод можно трактовать как регуляризацию некорректных задач по методу квазирешений. По сравнению с другими спектральными фильтрами этот метод более эффективный в рассматриваемых задачах. Разумеется, предложенный подход может быть эффективным лишь в некоторых вычислительных задачах. В настоящей работе рассматривается применение численных методов для получения приближенных решений систем уравнений, описывающих нестационарное течение идеальной жидкости со свободной поверхностью, где с помощью нашего подхода удалось получить важные результаты в задачах океанологии и математической гидродинамики.

При описании нелинейной динамики поверхности идеальной жидкости мы используем уравнения, полученные в работах [4, 5]. Эти уравнения эквивалентны системе уравнений Эйлера со свободной поверхностью. Уравнения, исследуемые в разд. 3, представляют собой системы нелинейных интегродифференциальных уравнений в частных производных, разрешенных относительно производных по времени. Мы используем методы прямых для сведения уравнений в частных производных к системам обыкновенных

* Работа выполнена при финансовой поддержке Президентской программы “Ведущие научные школы РФ” (грант № НШ-7550.2006.2) и международного гранта INTAS Ref. Nr 05-1000008-8014.

© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2008.

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

Мы рассматриваем следующие задачи гидродинамики со свободной поверхностью: течение идеальной жидкости со свободной поверхностью с конечной глубиной и начальную стадию неустойчивости Рэлея—Тейлора с исследованием границы между тяжелой идеальной жидкостью и вакуумом в условии, когда жидкость находится над вакуумом. В начальный момент свободной поверхности придается возмущение и исследуется развитие неустойчивости Рэлея—Тейлора. Эта модель — классический пример неустойчивости (см., например, [6, 7]).

Необходимо отметить, что предлагаемый метод не в полной мере подавляет численную неустойчивость в задачах гидродинамики со свободной поверхностью, а лишь “продлевает жизнь” численной схеме. Описываемые подходы с успехом применялись в исследованиях в Институте океанологии им. П.П. Ширшова РАН и в Институте теоретической физики им. Л.Д. Ландау РАН.

1. Постановка проблемы

Пусть Н — сепарабельное гильбертово пространство. Рассмотрим в Н, вообще говоря, нелинейный оператор А, определенный на подпространстве Рс Н. Предположим, что в Н можно выбрать ортонормированный базис } с Р.

Через Ск([0,Т]; Н), 0 < Т < то обозначим пространство к раз непрерывно дифференцируемых функций на [0,Т] со значениями в Н с нормой

Определение 1. Функция и € С*([0,Т]; Н) такая, что и(Ь) € Р при Ь € [0; Т], и удовлетворяющая (1) называется решением задачи (1).

Рассмотрим проекционный метод сведения к системе обыкновенных дифференциальных уравнений для задачи (1). Для конечного N введем проектор Р^ : Н ^ Н по формуле

«Исхода) = тах(ІК*)ІІн + Ии(к)С011я )•

и ІЄ[0,Т ]

Будем рассматривать абстрактную задачу Коши:

и'(і) = А[и(і)], і Є [0,Т],

(1)

и(0) = ф, ф Є Н•

При фиксированном N будем рассматривать функцию им(і) = ^^=1 (і)^к, где

а1^ Є С 1[0,Т] являются решением задачи Коши

(им)'(і) = Рм(А[и(і)]), і Є [0,Т],

(2)

им (0) = Рм ф.

Функция им полностью определяется вектор-функцией ам(і) = (ам(і),...,аМ(і))Т• Перепишем задачу (2) в виде системы ^го порядка обыкновенных дифференциальных уравнений относительно ам. Обозначим ^(ам) = (^\(ам), • •. , Рм(ам))т, где

" N

ЕаМ (і)^

І=1

^к) • Будем рассматривать следующую систему уравнен

Fk К (*)) = ^А

ний, эквивалентную (2):

(ам)'(£) = ^(ам(Ь)), Ь € [0,Т],

(3)

ам (0) = фм,

где фм = (фь... ,фм)Т, Фк = (Ф,^к)н.

Правую часть в уравнении (3) можно трактовать как функцию ^ ^ .

Приведем схему Рунге—Кутты четвертого порядка для системы (3). Через т обозначим шаг по переменной Ь. Соответственно, Ьт = тт, т = 0,... , М. Предположим, что на рассматриваемом интервале существует единственное решение задачи (3). Численные решения получим по следующей формуле:

(Ьт) = (Ьт-1) + ^ (к1 + 2к2 + 2кз + к4),

6

т-

к1 = ^ (ам (іт—1)),

к2 = ^(ам ^іт—1) + 2 к1 ^

к3 = ^(аМ ^т-1) + 2к2^

(4)

Ьт-1) + ^ к2

к4 = ^ (ам (Ьт-1) + ткз), т = 1,..., М.

В реальности при вычислениях на ЭВМ по схеме (4), как правило, возникают ошибки в силу ограниченной машинной точности. Основной источник ошибок — вычисление функции ^. Элементы вектора ам трактуются как коэффициенты ряда Фурье, поэтому погрешности в вычислении этой функции могут повлечь катастрофические последствия для проведения расчетов. Типична ситуация, когда по модулю эти коэффициенты быстро стремятся к нулю, а при вычислениях приходится иметь дело с числами, по модулю значительно меньшими гарантированных значащих цифр. Вычисляя по схеме (4), мы наблюдаем значения аи, которые связаны с точными значениями следующим соотношением:

|ам(*т) - ам(*т)| = С, т =1,...,М.

Конкретные значения ам(Ьт) зависят от реализации счета по схеме (4). Скажем, что последовательность (Ьт) вычислительно неустойчива, если числовая последовательность С не стремится к нулю при N ^ то. В этом случае относительная погрешность вычисления быстро стремится к бесконечности при N ^ то.

2. Регуляризация в случае вычислительной неустойчивости

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

Пусть выбрано число д > 0. Введем функцию Я, : ^ по формуле

Я,(Ж1,Ж2, . . . ,жм) = (г,(Ж1),Г,(Х2), . . . ,Г,(жМ)),

где

, . { х, |х| > д

г«(х) = 1 0, |х| < д

Величины будем рассчитывать по формуле

(Ьт) = Л, (ам (Ьт-1) + т/6(к1 + 2к2 + 2кз + к4)). (5)

Величины к1, к2, к3, к4 вычисляются согласно схеме (4). Выбор параметра д > 0 осуществляется из эмпирических соображений. Имеет смысл выбирать д на порядок больше, чем максимальное из значений ^т.

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

Покажем, что алгоритм, задаваемый формулой (5), является одной из реализаций метода квазирешений в теории регуляризации некорректных задач. Сделаем следующие предложения о разрешимости задач (1), (3).

Условие 1. Пусть существует Т > 0 такое, что:

1) при Ь € [0,Т] существует единственное решение задачи (1);

2) решение задачи (1) принадлежит пространству С2([0,Т]; Р1), где множество Р1 компактно в Р;

3) для достаточно больших N при Ь € [0,Т] существует единственное решение задачи (2), и это решение принадлежит пространству С2([0,Т];Р1).

Методы, позволяющие проверить выполнение условия 1, получены в работах [8, 9]. Введем функционал невязки для задачи (1) по формуле

3(и) = ||и' — А[и]||С([0,Т];Я) + Ни(0) — фНя.

Очевидно, что функционал 3 : С 1([0,Т]; Р) ^ К является непрерывным. Функция и € С 1([0,Т]; Р) является решением задачи (1) тогда и только тогда, когда 3(и) = 0.

В условии 1 нахождение приближенных решений задачи (1) может быть сведено к нахождению последовательности и1, минимизирующей функционал 3 на множестве С2([0,Т]; Р1). Строить минимизирующую последовательность можно с помощью решений задачи (2).

Теорема 1. Пусть выполнено условие 1, тогда имеем

Иш I(им) = 0,

N ^ж

где ^ суть решения задачи (2). Более того, ||^ — и*||С 1([о,т];Р) ^ 0, при N ^ то, где и* — решение задачи (1).

Доказательство. Подставим приближенные решения uN в функционал невязки

1 ^ ) — А ] ||С([0,Т];Н) + 11^(0) — ^НН =

Ж

= Н (и^ — ^ А [и^ — (1 — ^) А [и^ 11С([0,Т];Н) + ^ |2 =

fc=N +1

Ж

= Н(1 — РМ)А [и^ IIС ([0,Т ];Н) + X/ ^ |2'

к=М +1

Второе слагаемое в последней сумме стремится к нулю при N ^ то. Рассмотрим величину

Ж

Н(1 — ) А [иМ] 11С ([0,Т];Н) = 4шш^ X! |(А И , ^ )Н |2' (6)

*е[ ’ ] к=М+1

Поскольку при всех £ € [0,Т] имеем им(£) € ©1? то в силу компактности Р1 в Р и непрерывности оператора А : Р ^ Н получаем, что множество А[им(£)] ограничено в Н равномерно по N и £ € [0,Т]. Поэтому величина (6) также стремится к нулю при N ^ то.

В силу компактности вложения пространства С2([0, Т]; Р1) в С 1([0, Т]; Р) из последовательности им можно извлечь подпоследовательность иМр, сходящуюся в С 1([0, Т]; Р). Легко видеть, что пределом этой подпоследовательности является единственное решение задачи (1), функция и*. Действительно, из непрерывности функционала I следует, что I( Иш иМ;) = Иш I (им0 = 0. На самом деле и сама последовательность им сходится к и*. Пусть Ишвир ||uN — и* ||С1([0,т];р) = 1!ш Ци^ — и* ||С1([0,т];р). В силу компакт-

N^1^ ’ ’ Nq ^ж

ности в С 1([0,Т]; Р) подпоследовательности и^ можно считать, что эта подпоследовательность имеет пределом V*. Но тогда V* тоже будет решением задачи (1), а значит, имеем V* = и*. В итоге получаем, что Ишвир ||uN — и*||с 1([0,т];р) = 0. Следовательно,

N ^ж

имеем Иш ||uN — и* IIС 1([0Т];р) =0. □

N^ж и ’ " '

Теорема 1 является обоснованием известного в теории регуляризации некорректных задач метода квазирешений. В этой теореме рассматривается последовательность решений задач (2), однако вместо функций uN, которые на практике недостижимы, можно рассматривать различные приближения в пространстве С2([0,Т]; Р1). Эти приближения можно получать с помощью метода Рунге—Кутты (4) и подходящего метода интерполяции по переменной £.

Пусть фиксировано достаточно большое N. Через uNM обозначим последовательность, полученную по схеме (4) по М точкам и последующим интерполированием с использованием стандартных кубических сплайнов по узлам {£т}. Соответственно, через и;^м, д > 0 обозначим аналогичную последовательность, полученную по схеме (5).

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Теорема 2. Пусть выполнено условие 1. Предположим, что uNM, и^™ € С2([0,Т]; Р1) и Иш ||uNM — uN||с2([0,т];р1) = 0 при фиксированном N. Тогда для любого

M —— <Ж ’ ’

£ > 0 существуют N, М такие, что I (и^М < £, также существуют, возможно, другие, N, М и д > 0 такие, что I (и^М < £.

Доказательство. В силу теоремы 1 существует такой номер N0, что I (и^ < е/2 для всех N > N0. Введем обозначение ДNM = uNM — uN. Из предположений теоремы следует, что ДNM € С2([0,Т]; Р1). Заметим, что по построению uNM(0) = uN. Имеем

^ + | (и‘,и) — А|и*'и ' ' 1,2

ІН.

З (и"°М) = 3 (им° + ДМ°М) = II (им°)7 - А [и*° + ДМ°М] + (ДМ°М)' ІІС([0)Т];Я) +

+ ||и^° (О)-,01ІЯ < 1 (и^°) -А [и^° + Д^°М ] ||С([0,Т];Я) + 11 (ДМоМ) ||с([0,Т];Я) + ||и^° (0)-

Поскольку Ііт ||ДМ°М|С2([0 Т];р1) = 0 и в силу непрерывности оператора А : Р ^ Н

М

существует такое М0, что при всех М > М0 имеем

|(ДМ°М ) ([0>Т];Я) < Є/4,

|3 (ик°) - I (и"°)' - А К° + Д"°М] ||С(|0.Г];Н)1 < £/4.

Следовательно, для всех N > N0 и М > М0 имеем

3 (иММ) < £.

Пусть теперь N1 и Мі таковы, что 3 (и^1 Мі) < є/2. Обозначим Д

и*'1"*1) < £/2. Обозначим Д^М = и^М - м^М.

По построению функций и и в силу свойств кубических сплайнов имеем

Иш ЦД^^||с2([0Т];р1 = 0. Оценим значение функционала I (uNlMl). В силу непрерыв-

д—0

ности функционала I : С 1([0,Т]; Р) ^ К существует такое д0 > 0, что для всех д < д0 имеем ( ) ( ) ( ) ( )

II (и^1) — I (и^1) | = II (и^1) — I (uNlMl + Д^1) | < £/2' Следовательно, при всех N > ^, М > М1, д < д0 выполнено

I (иГ) < £'

Теорема 2 дает возможность строить минимизирующую последовательность для функционала I с помощью приближений и^™, которые можно строить конструктивно. А теорема 1 гарантирует, что минимизирующая последовательность сходится к решению задачи (1). В конкретных приложениях условия теорем 1 и 2 часто могут быть легко проверены. В частности, в задачах гидродинамики со свободной поверхностью, которые рассматриваются в следующем разделе, условия этих теорем выполнены [9, 10].

3. Примеры в гидродинамике

3.1. Поверхностные волны идеальной жидкости

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

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

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

Считая движение жидкости потенциальным, имеем

!(ж,у, Ь) = УФ(ж,у,Ь),

где !(ж,у,Ь) — двумерное поле скоростей, Ф(ж,у,Ь) — потенциал. Из условия несжимаемости жидкости = 0 следует, что потенциал скоростей удовлетворяет уравнению Лапласа:

ДФ(ж, у, Ь) = 0.

С этим уравнением связываются следующие граничные и начальные условия:

(П + ФжПж — Фу )|у=п(ж ,*) = 0

пк=о = По (ж),

Ф|*=0 = Фо(ж,у).

Здесь д — ускорение поля тяжести.

Удобно ввести величину Ф(ж, Ь) = Ф(ж, п(ж, £),£), которая является значением потенциала на свободной поверхности (см. [4]). В.Е. Захаровым было установлено, что переменные п и Ф — это канонически сопряженные величины, т. е.

— Л < у < п(ж, Ь),

— то < ж < то, Ь > 0.

дЬ £Ф, 5 Ф = £Н

где гамильтониан Н совпадает с полной энергией жидкости Н = Т + и,

ж п(х ,*)

Следуя работе [4], совершим конформное отображение области, занимаемой жидкостью в плоскости (ж,у), в полупространство в переменных (и,!'):

— то < и < то, —Л < V < 0.

После преобразования поверхность п(х,і) может быть представлена в параметрическом виде

у = у(и, і), х = и + Х(и,і).

Переменные Х(и,і) и у(и,і) связаны соотношением

У = Я[х],

где оператор Я — интегральный оператор вида

Я[/](и) =

/(и)

2Л, 7 8ІпЬ(п/2Л,(и' — и))

— те

Обратный к Я оператор Т имеет вид

ї-і [ /(и)

[/](и) 1 — е—п/(^(и—“'))

^и;.

Как показано в работе [4], переменные у(и,£) и Ф(и,£) полностью описывают движение жидкости и подчиняются следующей системе интегродифференциальных уравнений, разрешенных относительно производных по времени:

Я[Ф«]

+ х,

Я[Ф„]

Ф. = —

Фи — (я[Фи])2

27

— Т

Я[Фи]

7

Ф« — 0У,

(7)

(8)

где 3 = хи + уи — якобиан отображения.

Уравнения (7)-(8) будем рассматривать с периодическими граничными условиями. Пусть и £ (0, 2п), а функции у, Ф представлены в виде

у(и,«)= £ Ук(*)е‘к“, Ф(и,і)= £ Фк(*}е‘к“.

(9)

,=-те ,= -те

Операторы Я и Т в фурье-представлении имеют простой вид:

у, = Я, ж,, Я, = г 1апЬ кН,

Хй = Т, у,, Т, = -г еоШ кН.

Операция дифференцирования по переменной и в фурье-представлении имеет обычный вид:

Д, = гк.

Для вычисления операторов Я, Т, Д будем применять быстрое преобразование Фурье.

Численный эксперимент 1. Рассмотрим стоячие волны идеальной жидкости с конечной глубиной. Выберем следующие параметры: глубина Н = 6.0, ускорение свободного падения д = 1.0, шаг по времени т = 0.001. В качестве начальных условий возьмем следующие функции:

у(0, и) = 0.01 сое и, Ф(0,и) = 0.

те

те

те

Рис. 1. Счет по схеме (4), і = 20.0 Рис. 2. Счет по схеме (5), і = 100.0

Аппроксимируем ряды (9) конечными суммами с N = 512 слагаемыми и будем применять численные схемы, описанные выше. Рассмотрим спектры решений в логарифмическом масштабе: по оси абсцисс отложим номер гармоники к, а по оси ординат — ^ |ук | и ^ |Ф& |. Приведем результаты численного эксперимента, которому был присвоен номер ZF-2-3 в системе учета вычислительных экспериментов www.calcs.ru. В этом эксперименте расчет проводился по схеме (4) — на рис. 1 показан спектр решения при і = 20.0. Гармоники с номерами от 0 до 25 характеризуют собственно решение, горизонтальная часть спектра с гармониками от 26 до 360 обусловлена конечностью машинной точности при вычислениях. Растущая часть спектра с гармониками от 360 до 511 — результат вычислительной неустойчивости схемы (4) в условиях машинной точности.

В численном эксперименте ZF-2-4 применялась модифицированная схема метода прямых (5), с выбранным д = 10-12, с тем же начальным условием. На рис. 2 показан спектр решения при і = 100.0. В отличие от случая на рис. 1 здесь мы не наблюдаем вычислительной неустойчивости.

3.2. Неустойчивость Рэлея—Тейлора

Неустойчивое течение Рэлея—Тейлора будем рассчитывать, используя уравнения Дьяченко, полученные в [5]. Эти уравнения — модификация уравнений (7)-(8) для случая бесконечной глубины.

Образуем пару комплексных функций г(т, і) = х(т, і) + іу(т, і) и на действительной оси П(и,і) = Ф(и,і) + іН[Ф(и,і)], где т = и + іи. Введем новые переменные Я(т,і) и V(т,і) по следующим формулам:

Я(т,і) = —, V (т,і) = і .

Функции Я и V аналитичны в нижней полуплоскости и удовлетворяют следующим условиям:

Я(т,і) ^ 1, |т| ^ то, 1тт < 0,

V(т,і) ^ 0, |т| ^ то, 1тт < 0.

Как показано в работе [5], функции Я и V удовлетворяют следующей системе интег-родифференциальных уравнений:

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Я* = *(иЯад — Цад Я), (10)

V* = г(Ц^ — Вад Я) + д(Я — 1), (11)

где и = Р^Я + VЯ), В = P(VV), Р = -(I + гИ) — интегральный оператор. Эти

уравнения называются уравнениями Дьяченко. Разрешимость уравнений (10)—(11) рассматривалась в [8, 9]. Сходимость численных методов была установлена в работе [10].

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

Численный эксперимент 2. Для проведения численных расчетов используем метод прямых, приняв N = 512. Выберем следующие параметры: ускорение свободного падения д = —10.0, шаг по времени: т = 0.0001. В качестве начальных условий возьмем следующие функции:

Я(0, ад) = 1 + 0.01 ехр(—гад), V (0,ад) = 0.

Используем результаты численного эксперимента И.Т-3-5, где расчет проводился по схеме (4). На рис. 3 показан спектр решения при £ = 0.3. Первые гармоники от 0 до 10 характеризуют решение, а возрастающие по модулю гармоники (от 10) выражают факт неустойчивости. В данном примере мы сталкиваемся с неустойчивостью самих уравнений. Возникающие погрешности вычислений быстро разрушают численный счет.

В численном эксперименте И.Т-3-6 проверим нашу модифицированную схему (5) в этом случае с параметром q = 10-12. На рис. 4 дан спектр решения при £ = 1.55. Здесь мы не наблюдаем какой-либо вычислительной неустойчивости. Разумеется, дальнейший счет по схеме (5) требует увеличения размерности N и приводит к разрушению гладкого решения, но это следствие неустойчивости течения Рэлея—Тейлора. На рис. 5 показан профиль свободной поверхности течения при £ = 1.55.

Рис. 3. Счет по схеме (4), Ь = 0.3

Рис. 4. Счет по схеме (5), Ь = 1.55

Рис. 5. Профиль свободной поверхности при t = 1.55

Автор выражает глубокую благодарность академику РАН В.Е. Захарову и доктору физ.-мат. наук А.И. Дьяченко за постановку задачи и постоянное внимание к работе, а также доктору физ.-мат. наук Л.Б. Чубарову за ценные замечания и полезное обсуждение.

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

[1] Шокин Ю.И. Интервальный анализ. Новосибирск: Наука, 1981. 112 с.

[2] Калмыков С.А., Шокин Ю.И., Юлдашев З.Х. Методы интервального анализа. Новосибирск: Наука, 1986. 224 с.

[3] Алефельд Г., Херцвергер Ю. Введение в интервальные вычисления. М.: Мир, 1987. 360 с.

[4] Дьяченко А.И., Захаров В.Е., Кузнецов Е.А. Нелинейная динамика свободной поверхности идеальной жидкости // Физика плазмы. 1999. Т. 22, № 10. С. 916-928.

[5] Zakharov V.E., Dyachenko A.I., Vasilyev O.A. New method for numerical simulation of a nonstationary potential flow of incompressible fluid with a free surface // European J. of Mechanics B/Fluids. 2002. Vol. 21. P. 283-291.

[6] Taylor G. The instability of liquid surface when accelerated in a direction perpendicular to their planes. I // Proc. Roy. Soc. 1950. T. 201. Ser. A, N 1065. P. 192-196.

[7] ИнногАмов Н.А. Турбулентная стадия тейлоровской неустойчивости // Письма в ЖТФ. 1978. Т. 4, вып. 12. C. 743-747.

[8] Шамин Р.В. О существовании гладких решений уравнений Дьяченко, описывающих неустановившиеся течения идеальной жидкости со свободной поверхностью // Докл. РАН. 2006. Т. 406. № 5. C. 112-113.

[9] Шамин Р.В. К вопросу об оценке времени существования решений системы Коши— Ковалевской с примерами в гидродинамике со свободной поверхностью // Современная математика. Фундаментальные направления. 2007. Т. 21. C. 133—148.

[10] Шамин Р.В. Об одном численном методе в задаче о движении идеальной жидкости со свободной поверхностью // Сиб. ЖВМ. 2006. Т. 9, № 4. C. 325-340.

[11] НАлимов В.И. Нестационарные вихревые волны // Сиб. математ. журнал. 1996. Т. 37, № 6. C. 1356-1366.

[12] Craig W., Sulem C. Numerical simulation of gravity waves // J. Comput. Phys. 1993. Vol. 108. P. 73-83.

Поступила в редакцию 28 июня 2007 г.

i Надоели баннеры? Вы всегда можете отключить рекламу.