Нелинейная монотонизация схемы К.И. Бабенко для численного решения квазилинейного уравнения переноса

Размещено на http://www.allbest.ru//

Размещено на http://www.allbest.ru//

РОССИЙСКАЯ АКАДЕМИЯ НАУК

Ордена Ленина Институт прикладной математики

им. М. В. Келдыша

Нелинейная монотонизация схемы К. И. Бабенко для численного решения квазилинейного уравнения переноса

Т. А. Александрикова

М. П. Галанин

Москва - 2003

Аннотация

Работа посвящена тестированию предложенного авторами варианта нелинейной монотонизации схемы К.И. Бабенко для численного решения квазилинейного уравнения переноса, а также его сравнению с другими известными конечно-разностными схемами. Для сравнения использовались явная и неявная схемы с левой разностью первого порядка аппроксимации, схема Lax-Wendroff'a второго порядка аппроксимации, а также монотонизованная схема “Кабаре”.

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

T.A. Alexandrikova, M.P. Galanin

Nonlinear monotonization of K. I. Babenko scheme for the numerical solution of the quasi - linear advection equation

Abstract

The paper is dedicated to testing the authors offered variant of nonlinear monotonised К.I. Babenko scheme for numerical solution of quasi - linear advection equation, as well as its comparison with other well - known finite - difference schemes. There were used for the comparison the implicit and explicit schemes with the left difference of first - order approximations, the Lax - Wendroff scheme of second order approximations, as well as monotonised "Cabaret" scheme.

The presented information about the errors of numerical solutions allows to compare a quality of the schemes. Benchmark analysis of the errors has shown that the offered scheme gives higher accuracy in the broad range of Courant numbers, as well as smaller "smudging" decision of shocks.

Введение

Уравнение переноса является представителем фундаментальных уравнений математической физики [1], широко используемым для описания движения сплошной среды. Многие задачи газовой динамики, магнитной гидродинамики, астрофизики связаны с решением уравнений такого типа. Учет движения элемента сплошной среды приводит к уравнению переноса квазилинейного типа. Основным следствием нелинейности этих уравнений является “опрокидывание” начального профиля за конечное время и возникновение ударных волн. Изучение закономерностей распространения ударных волн представляет особый интерес в исследовании течений сплошной среды. Наличие разрывов в решении представляет определенные трудности, при этом точность численного решения зависит от качества воспроизведения этих разрывов. На настоящий момент существует большое количество разностных схем для решения задач такого типа. Однако большинство из них дают “расплывчатое” решение, в том числе “размазанный” фронт ударных волн, а схемы повышенного порядка аппроксимации могут приводить к искажениям численного решения, в том числе осцилляциям, т.е. немонотонности решения [2]. Различными способами можно повысить качество численного решения. Один из них связан с использованием неоднородных конечно-разностных схем, в которых реализуется контроль за перемещением разрыва и изменение алгоритма вычислений в его окрестности. Однако основная проблема в использовании таких схем заключается в том, что для нелинейных уравнений гиперболического типа характерно появление разрывных решений даже при гладких начальных условиях, так что положение особых точек заранее неизвестно. К тому же число разрывов может меняться со временем и к каждому из них нужно приспосабливать алгоритм расчета отдельно. Поэтому особый интерес вызывает исследование однородных схем повышенного порядка аппроксимации, а также построение новых квазимонотонных схем на их основе. При этом существуют разнообразные способы монотонизации. В данной работе монотонизация схемы К. И. Бабенко проведена при помощи введения в разностное уравнение искусственной вязкости. Также в данной работе рассмотрена предложенная в [3] схема, монотонизованная с использованием алгоритма, основанного на знании области зависимости точного решения.

В настоящей работе представлен вывод предложенного авторами варианта нелинейной монотонизации схемы К.И. Бабенко для квазилинейного уравнения переноса. Также приведены результаты сравнительного анализа новой и некоторых других известных конечно-разностных схем (явной и неявной схем с левой разностью [2, 4], схемы Lax-Wendroff'a [4] и монотонизованной схемы “Кабаре” [3, 5, 6]). Сравнение численных результатов, полученных по рассмотренным схемам, проведено на системе тестов, аналогичной использованной в работах [4, 5, 7].

Решения получены для финитных начальных условий и различных чисел Куранта на пространственно-временном прямоугольнике достаточно больших размеров. Для определения ошибок численных решений использованы конечномерные анал оги норм пространств C, L1, L2. В работе представлены рисунки с точным и численными решениями для некоторых моментов времени, а также таблицы ошибок для тех же моментов времени и для всего временного промежутка. Такой способ представления результатов позволяет качественно и количественно оценивать характеристики схем.

Авторы выражают благодарность Елениной Т.Г. за внимание, советы и полезные обсуждения.

Работа выполнена при частичной финансовой поддержке Российского фонда фундаментальных исследований (проект РФФИ № 03-01-00461).

1. Постановка задачи

Будем рассматривать задачу Коши для квазилинейного уравнения переноса

, , (1.1)

с финитными начальными условиями u(x,0)= u0(x) вида:

1. “треугольник”,

2. “прямоугольник”,

3. “левый треугольник”,

4. “правый треугольник”,

5. “ступенька вниз”,

6. “ступенька вверх”.

Точное решение уравнения (1.1) имеет вид функционального соотношения

u(x,t)= f(x-ut),

где функция f определяется начальным возмущением при :

Поясним решение (1.1), полученное методом характеристик [8]. Заметим, что выражение является полной производной функции u(x,t):

вдоль характеристики G, которая задается уравнением

.

Из выражения (1.1) получаем . Это означает, что функция u(x,t) принимает постоянное значение

u(x,t) = u0()

на кривой G, уравнение которой имеет вид

x = + u0()t.

Приведем точные решения для указанных выше начальных условий:

“треугольник”:

“прямоугольник”:

“левый треугольник”:

“правый треугольник”:

“ступенька вниз”:

“ступенька вверх”:

Сравнение численного решения с точным будем проводить с использованием конечномерных аналогов норм в пространствах C, L1, L2 на =(-,+)[0,T]:

, , ;

и на R=(-,+) при t=ti :

, , .

2. Нелинейная монотонизация схемы К. И. Бабенко

Как известно [5], cхема К.И. Бабенко [9] имеет второй порядок аппроксимации по x и t, но является немонотонной. Cхема является неявной, однако это не мешает вычислять искомое решение на новом временном слое “бегущим счетом”. При числе Куранта, равном единице, для линейного уравнения переноса схема дает точное решение.

Для численного решения введем равномерную (для простоты) пространственно-временную сетку

где h, - шаги разностной сетки по x и t соответственно.

Здесь и далее будем использовать стандартные обозначения для сеточных величин: , где верхний индекс - номер временного слоя, нижний - номер узла по х. Решение на слое n считаем известным.

Перепишем уравнение (1.1) в виде

.

Схему Бабенко (“квадрат” [9]) для рассматриваемого уравнения можно записать следующим образом:

(2.1)

где , или

(2.2)

Пусть , тогда (2.2) перепишем в виде:

Или

(2.3)

где через =a/h обозначено число Куранта.

Заметим, что первые два слагаемых (2.3) дают обычную схему с левой разностью.

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

(2.4)

Заметим, что при ? 0 получим исходную схему (2.1), (2.3), а при ? 1 - схему с левой разностью.

Запишем (2.4) в виде

Где

. (2.5)

Предположим, что y 0 во всей пространственно-временной области, тогда схема (2.5) удовлетворяет принципу максимума [10], если выполнены неравенства

(2.6)

Отметим, что буквальное условие применимости принципа максимума отличается от второго неравенства (2.6): в нем фигурирует . Однако из второго неравенства (2.6) в предположении следует нужное условие.

Введем параметр и положим =(R,),

-1=(R-1,-1). Заметим, что (2.7)

Неравенства (2.6) выполнены, если функция =(R,) имеет вид [3]:

(2.8)

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

Параметр R* в данной работе подбирается экспериментально. Для численных расчетов, как и в [5], использовалось R* =1.2.

Уравнение (2.4) является нелинейным. Из (2.7) следует, что искомую величину можно вычислить по известному R:

(2.9)

Используя полученную зависимость =(R,), найдем соотношения для R. Распишем уравнение (2.4):

(2.10)

которое с учетом (2.7) приобретает вид

(2.11)

Обозначим правую часть (2.11) через b:

. (2.12)

Отметим, что (2.11) представляет собой нелинейное уравнение относительно R с известной правой частью. При этом левая часть (2.11) - монотонная функция R. Учитывая вид функции =(R, ), получим следующие расчетные формулы:

(2.13)

Замечания:

при значение b не определено, однако из (2.7) следует, что в этом случае R = 0.

из (2.13) видно, что R не определено при b = 0, но тогда из (2.12) и (2.10) непосредственно вычисляется искомое значение .

Таким образом, используя формулы (2.9), (2.12), (2.13), можно определить искомую величину . Отметим, что зависит от , поэтому численно величина определяется итерационным способом. В данной работе в качестве первого приближения выбирается равной , для которой определяется соответствующее значение г и по известному b находится значение R. Далее по формуле (2.9) определяется неизвестная величина , тем самым осуществляется переход на следующую итерацию при условии невыполнения критерия окончания итерационного процесса.

3. Другие разностные схемы

Для сравнения предложенной схемы в данной работе были рассмотрены некоторые широко известные конечно-разностные схемы, а именно: явная и неявная схемы с левой разностью, схема Lax-Wendroff'a, а также предложенный в [3] метод прыжкового переноса. В большинстве из них число Куранта выражается через решение не так, как в § 2. Однако вид его достаточно очевиден, так что мы сохраняем обозначение. Ниже представлена информация об основных свойствах этих схем. Более подробно о них можно узнать в указанной литературе.

Схема с левой разностью [2, 4]

Явная схема, имеющая первый порядок аппроксимации по пространству и по времени. Схема устойчива при числе Куранта 1, при этом в случае линейного уравнения переноса для = 1 дает точное решение.

Неявная схема с левой разностью [2, 4]

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

Схема Lax-Wendroff'a [4]

Явная схема, аппроксимирующая исходное уравнение со вторым порядком по времени и по пространству. Является устойчивой при 1, при =1 для линейного уравнения переноса дает точное решение.

Метод прыжкового переноса [3, 5, 6, 11]

В работе [3] предложена схема, построенная с помощью монотонизации схемы “Кабаре” [3], и на ее основе разработан алгоритм метода прыжкового переноса.

Трехслойная схема “Кабаре”, аппроксимирующая уравнение конвективного переноса со вторым порядком точности, является бездиссипативной и обладает улучшенными по сравнению с классическими линейными схемами дисперсионными свойствами [6]. Она устойчива при 1, дает точное решение в случае линейного уравнения переноса для = 1 и = 0.5, однако является немонотонной.

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

4. Организация вычислительного эксперимента

Аналогично работам [4, 5, 7] численное решение поставленной задачи (1.1) получено для всех указанных форм начальных профилей при следующих значениях параметров: l=520, l1=0, l2=20, T=1000, h=1. Число Куранта () ограничено сверху величиной m, которая принимала значения: 0.1, 0.25, 0.5, 0.9. Для неявной схемы с левой разностью также были проведены расчеты при m = 3. Таким образом, временной шаг принимал те же значения, что и число Куранта. Для определения точности полученного при этом решения вычислялись относительные ошибки с использованием норм

, ,

интегрально по временному промежутку 0 t T, а также

, ,

локально для момента времени tj.

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

5. Результаты численного исследования

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

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

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

Явная схема с левой разностью

Схема передает решение весьма грубо, разрывы решений “размазываются” на несколько шагов (от 4 до 10). Размазывание решения тем сильней, чем меньше . Решения типа волн разрежения размазываются на большее количество шагов по сравнению с ударными волнами. Этот дефект решения усиливается со временем. Для начальных профилей №№ 1, 2, 4 на тех интервалах времени, на которых происходит эволюция точного решения без изменения амплитуды, амплитуда численного решения непрерывно уменьшается. Поведение решения в эти моменты времени и в последующие, сопровождающиеся образованием и дальнейшим распространением ударных волн, практически ничем не отличается. Моменты образования ударных волн трудно прослеживаются. Скорость распространения ударных волн не зависит от и совпадает со скоростью распространения ударных волн в точном решении. монотонизация схема бабенко квазилинейный уравнение

Таблицы ошибок для явной схемы с левой разностью (локально при t = T)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9982	0,9369	0,9504	0,2737	0,0328	0,2941
L1	0,0633	0,0448	0,0539	0,0433	0,0093	0,0027
L2	0,2069	0,1468	0,1806	0,0641	0,0121	0,0253
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9979	0,9495	0,9558	0,2973	0,0306	0,2724
L1	0,0621	0,0448	0,0538	0,0424	0,0083	0,0024
L2	0,2048	0,1514	0,1829	0,0639	0,0116	0,0236
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9972	0,9672	0,9641	0,3471	0,0261	0,2318
L1	0,0601	0,0448	0,0538	0,0409	0,0068	0,0021
L2	0,2014	0,1591	0,1869	0,0656	0,0092	0,0205
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9957	0,9874	0,9753	0,4397	0,0138	0,5071
L1	0,0569	0,0447	0,0538	0,0389	0,0058	0,0022
L2	0,1958	0,1717	0,1931	0,0726	0,0064	0,0327

Таблицы ошибок для явной схемы с левой разностью (интегрально по t)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,521	0,7489	0,7146	0,6176	0,3427	0,7489
L1	0,0712	0,0499	0,0548	0,0712	0,0055	0,0066
L2	0,1937	0,1115	0,1502	0,1661	0,0112	0,0523
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5174	0,7263	0,7362	0,6083	0,3585	0,7263
L1	0,0687	0,0478	0,0543	0,0678	0,0051	0,0061
L2	0,1888	0,1084	0,1533	0,158	0,0104	0,0499
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,4939	0,6757	0,7752	0,5932	0,3906	0,6757
L1	0,0644	0,0445	0,0534	0,0620	0,0041	0,0051
L2	0,1803	0,1057	0,1591	0,1433	0,0089	0,0456
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,4879	0,8217	0,8049	0,5421	0,45	0,8217
L1	0,0575	0,0402	0,0522	0,0521	0,0034	0,0045
L2	0,1662	0,1173	0,1698	0,1161	0,0064	0,0483

Неявная схема с левой разностью

Решение, полученное при помощи этой схемы, повторяет точное еще хуже, чем в случае аналогичной явной схемы. Для этой схемы характерно еще более сильное размазывание разрывов (от 7 до 15 шагов по пространству), которое усиливается с увеличением числа Куранта. Как видно из рисунков 5.2.1-5.2.6, для численного решения характерно "выгибание" заднего фронта на начальных моментах времени, что усиливается при увеличении . С течением времени зависимость численного решения от числа Куранта ослабевает. Для начальных профилей №№ 1, 2, 3, 4 численное решение, полученное по явной и неявной схемам с левой разностью, эволюционирует в волну треугольного профиля с распространяющейся вправо ударной волной и уменьшающейся амплитудой, что соответствует точному решению. Однако максимум амплитуды численного решения отстает от максимума точного. На рисунках численному решению, наиболее близкому к точному, соответствует = 0.1.

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

Таблицы ошибок для неявной схемы с левой разностью (локально при t = T)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9986	0,9176	0,9428	0,3071	0,0054	0,5485
L1	0,0649	0,0449	0,0539	0,0446	0,0098	0,0017
L2	0,2095	0,1407	0,1774	0,0651	0,0086	0,0257
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9988	0,9012	0,9367	0,3326	0,0058	0,5517
L1	0,0661	0,0451	0,0539	0,0455	0,0106	0,0018
L2	0,2115	0,1362	0,1751	0,0667	0,0093	0,0262
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9991	0,8704	0,9257	0,3759	0,0065	0,5561
L1	0,0681	0,0452	0,0540	0,0473	0,0118	0,0019
L2	0,2147	0,1286	0,1711	0,0704	0,0103	0,0271
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9994	0,8123	0,9064	0,4459	0,04345	0,5615
L1	0,0712	0,0455	0,0542	0,0505	0,01763	0,0047
L2	0,2197	0,1167	0,1646	0,0784	0,02047	0,0415
г = 3	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	1	0,4158	0,7721	0,7639	0,0148	0,5748
L1	0,0873	0,0512	0,0556	0,0686	0,0232	0,0037
L2	0,2441	0,0722	0,1319	0,1319	0,0209	0,0352

Таблицы ошибок для неявной схемы с левой разностью (интегрально по t)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5255	0,7309	0,6944	0,6198	0,3247	0,7309
L1	0,0746	0,0528	0,0557	0,0758	0,0073	0,0034
L2	0,2	0,1163	0,1465	0,1763	0,0118	0,0368
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5264	0,6885	0,6851	0,6184	0,3128	0,6907
L1	0,0772	0,0551	0,0564	0,0791	0,0079	0,0035
L2	0,2047	0,1201	0,1441	0,1835	0,0125	0,0366
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5249	0,6374	0,6705	0,6193	0,2956	0,6424
L1	0,0813	0,0588	0,0576	0,0844	0,00894	0,0037
L2	0,2121	0,1269	0,1407	0,1948	0,01373	0,0364
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5178	0,6813	0,6119	0,6119	0,2519	0,6887
L1	0,0877	0,0652	0,0596	0,0929	0,0995	0,0091
L2	0,2237	0,1417	0,1363	0,2117	0,0171	0,0583
г = 3	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5213	0,6763	0,4811	0,6035	0,2153	0,5748
L1	0,1191	0,0986	0,0729	0,1313	0,0173	0,0069
L2	0,2718	0,2035	0,1353	0,2707	0,0233	0,0469

t=6 t=10



t=70 t=170

Рис. 5.2.1 Численное решение по неявной схеме с левой разностью для “треугольника” ( = 0,1; 0,5; 0,9)

t=18 t=40

t=90 t=160

Рис. 5.2.2 Численное решение по неявной схеме с левой разностью для “прямоугольника” ( = 0,1; 0,5; 0,9)

t=60 t=160

Рис.5.2.3 Численное решение по неявной схеме с левой разностью для “левого треугольника” ( = 0,1; 0,5; 0,9)



t=5 t=20

Рис. 5.2.4 Численное решение по неявной схеме с левой разностью для “правого треугольника” ( = 0,1; 0,5; 0,9)



t=20 t=80

Рис. 5.2.5 Численное решение по неявной схеме с левой разностью для “ступеньки вниз” ( = 0,1; 0,5; 0,9)



t=10 t=40

Рис. 5.2.6 Численное решение по неявной схеме с левой разностью для “ступеньки вверх” ( = 0,1; 0,5; 0,9)

3. Схема Lax-Wendroff'a

Для данной схемы характерно явление нормальной дисперсии. Образование ударных волн сопровождается появлением осцилляций за их фронтом. Для начальных профилей №№ 2, 4, 6 характерно возникновение осцилляций на левом фронте образующейся волны разрежения. Эти осцилляции прогрессируют со временем, что приводит к увеличению погрешности вычислений (особенно в норме С). Так, при = 0.5 для начального профиля № 2 в момент времени t = 10 относительная ошибка вычислений в норме С равнялась 0.17, а в момент t = T - 1.076.

Расплывание решения проявляется в меньшей мере, чем у решений, полученных по явной и неявной схемам с левой разностью. Разрывы решений передаются этой схемой лучше, размазывание разрывов происходит на 3-5 шагах.

Таблицы ошибок для схемы Lax-Wendroff'a (локально при t = T)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	1,4261	4,2852	1,4544	6,0711	2,2936	0,3209
L1	0,2426	6,3731	0,2661	12,842	1,4617	0,0056
L2	0,4473	5,3225	0,4876	9,2826	1,7353	0,0313
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,3402	1,4331	1,0226	2,0332	0,6961	0,2894
L1	0,0742	1,1533	0,0914	2,2578	0,1116	0,0029
L2	0,1313	1,2521	0,2484	2,0841	0,2461	0,0263
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,4511	1,0761	0,7487	0,8188	0,1495	0,3139
L1	0,0543	0,3498	0,0521	0,6471	0,0259	0,0026
L2	0,1291	0,4583	0,1769	0,7067	0,0479	0,0254
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	1	1,9891	1	2,7844	1,8783	1
L1	0,5184	2,6625	0,4483	5,2051	0,5182	0,2184
L2	0,7692	2,3194	0,7377	4,0324	0,6205	0,4668

Таблицы ошибок для схемы Lax-Wendroff'a (интегрально по t)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,7204	1,9162	0,8304	1,9224	2,2932	0,7472
L1	0,2043	2,0458	0,2549	4,0616	0,3489	0,0116
L2	0,3255	2,1123	0,4379	3,5791	0,7501	0,0527
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,7142	0,6417	0,7118	0,7588	0,6961	0,6417
L1	0,0754	0,5434	0,0949	1,0543	0,0325	0,0058
L2	0,1587	0,5963	0,2505	0,9867	0,0983	0,0384
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,7223	0,6093	0,7395	0,7727	0,3734	0,5303
L1	0,0603	0,2083	0,0612	0,3878	0,0148	0,0045
L2	0,1671	0,2523	0,1994	0,4165	0,0417	0,0324
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,8262	1	0,7562	0,8494	1,0426	1
L1	0,4144	1,7086	0,3655	2,8081	0,1244	0,2142
L2	0,6521	1,3794	0,6293	2,1012	0,2012	0,4619

4. Метод прыжкового переноса

По сравнению с указанными выше схемами решение, полученное с использованием метода прыжкового переноса, лучше передает точное. Основным достоинством метода является способность точно передавать разрывы типа ударных волн. Однако в силу специфики данного метода решение меняется скачкообразно, поэтому точное положение ударных волн можно наблюдать только через определенные интервалы времени. Так, для начального профиля № 5 при = 0.5 решение на каждом 4-ом временном слое скачком смещается вправо на один шаг по x, догоняя тем самым точное решение. Таким образом, ошибка вычислений в эти моменты времени нулевая. Однако увеличение числа Куранта приводит к появлению дефектов численного решения в виде “ступенек” и неверных скоростей движения ударных волн. Из рис. 5.4.3 видно, что для начального условия № 4 на волне сжатия, независимо от числа Куранта, образуются “ступеньки”, а рис. 5.4.2, 5.4.5 показывают, что при = 0.9 ударные волны в численном решении движутся с меньшими скоростями, чем ударные волны в точном решении. На рис. 5.4.1-5.4.5 численному решению, наиболее близкому к точному, соответствует = 0.5.

Таблицы ошибок для метода прыжкового переноса (локально при t = T)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	6,8416	1,0153	1,0355	6,6942	0,0019	1
L1	12,6125	0,0304	0,0571	12,6935	0,0016	0,0019
L2	8,1811	0,2135	0,2947	8,8841	0,0018	0,0442
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	6,7684	1	1	5,5631	0,0018	1
L1	11,7261	0,0304	0,0142	9,6542	0,0017	0,0019
L2	8,1673	0,2103	0,1446	7,8144	0,0018	0,0448
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	1	1	1	6,6931	0,0018	1
L1	0,2248	0,0986	0,1095	8,3942	0,0019	0,0019
L2	0,5629	0,3772	0,3989	6,2856	0,0024	0,0448
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	4,7944	3,5594	1	7,0426	0,2626	1
L1	4,6922	1,6669	0,1746	6,1921	0,1554	0,2458
L2	4,6991	2,1543	0,4995	4,6883	0,1954	0,4958

Таблицы ошибок для метода прыжкового переноса (интегрально по t)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	1	1	0,9581	1	0,5183	1
L1	6,2321	0,0232	0,0557	6,1564	0,0012	0,0019
L2	4,4215	0,1822	0,2894	4,6653	0,0034	0,0439
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	1	1	0,9644	1	0,5502	1
L1	5,5273	0,0276	0,0089	4,3492	0,0012	0,0019
L2	4,1864	0,1832	0,1213	3,8435	0,0035	0,0439
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	1	1	0,9754	1	0,6252	1
L1	0,2053	0,0850	0,0971	3,7451	0,0014	0,0019
L2	0,5182	0,3117	0,3521	3,0393	0,0039	0,0439
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,9292	1	0,9928	1	0,9	1
L1	2,0043	0,7447	0,1543	3,1267	0,0774	0,2357
L2	2,2779	1,0711	0,4394	2,7382	0,1552	0,4854

t=18 t=40

Рис. 5.4.1 Численное решение для метода прыжкового переноса для “прямоугольника” (=0,5; 0,9)

t=60 t=160

Рис. 5.4.2 Численное решение для метода прыжкового переноса для

“левого треугольника” (=0,5; 0,9)

t=5

Рис. 5.4.3 Численное решение для метода прыжкового переноса для

“правого треугольника” ( = 0,5; 0,9)

t=10 t=40

Рис. 5.4.4 Численное решение для метода прыжкового переноса для

“ступеньки вверх” ( = 0,5; 0,9)

t=20 t=80

Рис. 5.4.5 Численное решение для метода прыжкового переноса для

“ступеньки вниз” ( = 0,5; 0,9)

5. Монотонизованная схема “Бабенко”

Полученное по данной схеме решение наилучшим по сравнению с рассмотренными выше схемами образом передает точное во всем изученном диапазоне изменения чисел Куранта. У схемы не наблюдается явлений нормальной и аномальной дисперсии, расплывание решения незначительно - размазывание разрывов происходит на меньшее количество шагов по сравнению с другими рассмотренными схемами (1-3 шага). Наблюдается небольшое смещение вправо левого фронта образующейся волны разрежения (для начальных профилей №№ 2, 4, 6). Рисунки 5.5.1-5.5.6 показывают, что численные решения при различных числах Куранта практически не отличаются друг от друга, при этом численному решению, наиболее близкому к точному, соответствует = 0.1.

Таблицы ошибок для схемы К. И. Бабенко (локально при t = T)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,8768	1,0191	1,0136	1,0154	0,0079	0,0896
L1	0,0419	0,0615	0,0752	0,0728	0,0016	0,0007
L2	0,1316	0,2524	0,2572	0,2229	0,0023	0,0082
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,6804	1,0231	1,0217	1,0183	0,0076	0,0898
L1	0,0183	0,0564	0,0554	0,0569	0,0014	0,0007
L2	0,0994	0,2644	0,2592	0,2418	0,0020	0,0082
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5657	1,0246	1,0245	1,0229	0,0088	0,0952
L1	0,0128	0,0584	0,0553	0,0641	0,0013	0,0008
L2	0,0821	0,2733	0,2723	0,2639	0,0020	0,0087
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,1683	1,0254	1,0263	1,0216	0,0095	0,5001
L1	0,0056	0,0643	0,0581	0,0703	0,0017	0,0021
L2	0,0253	0,2815	0,2822	0,2775	0,0023	0,0221

Таблицы ошибок для схемы К. И. Бабенко (интегрально по t)

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5475	0,8818	0,8634	0,5623	0,4534	0,8818
L1	0,0208	0,0453	0,0559	0,0509	0,0011	0,0040
L2	0,0758	0,2016	0,2449	0,1786	0,0034	0,0493
г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5267	0,8727	0,8706	0,5495	0,4565	0,8413
L1	0,0132	0,0472	0,0501	0,0501	0,0013	0,0038
L2	0,0675	0,2142	0,2473	0,1931	0,0034	0,0477
г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,5183	0,8725	0,8748	0,5499	0,4628	0,7566
L1	0,0124	0,0505	0,0524	0,0593	0,0019	0,0034
L2	0,0657	0,2196	0,2564	0,2166	0,0035	0,0428
г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	0,4525	0,9127	0,8955	0,4811	0,4591	0,9127
L1	0,0115	0,0555	0,0549	0,0674	0,0012	0,0041
L2	0,0726	0,2312	0,2607	0,2312	0,0039	0,0345



t=6 t=10

t=70 t=170

Рис. 5.5.1 Численное решение по схеме К. И. Бабенко для “треугольника” ( = 0,1; 0,5; 0,9)

t=18 t=40

t=90 t=160

Рис. 5.5.2 Численное решение по схеме К. И. Бабенко для “прямоугольника” ( = 0,1; 0,5; 0,9)

t=60 t=160

Рис. 5.5.3 Численное решение по схеме К. И. Бабенко для “левого треугольника” ( = 0,1; 0,5; 0,9)



t=5 t=20

Рис. 5.5.4 Численное решение по схеме К. И. Бабенко для “правого треугольника” ( = 0,1; 0,5; 0,9)

t=20 t=80

Рис. 5.5.5 Численное решение по схеме К. И. Бабенко для “ступеньки вниз” ( = 0,1; 0,5; 0,9)

t=10 t=40

Рис. 5.5.6 Численное решение по схеме К. И. Бабенко для “ступеньки вверх” ( = 0,1; 0,5; 0,9)

6. Сравнение схем

На основе полученных численных результатов были сравнены рассмотренные схемы. Ниже приведены таблицы сравнения, содержащие названия тех схем, которые при данном числе Куранта дали наименьшую относительную ошибку в соответствующей норме. Если в ячейке таблицы указаны несколько схем, то ошибки вычислений, полученные для данных схем, либо совпадают, либо незначительно отличаются друг от друга. При этом используются обозначения: явная схема с левой разностью - ЯЛ, неявная схема с левой разностью - НЛ, схема Lax-Wendroff'a - LW, метод прыжкового переноса - ПП, схема К.И. Бабенко - Б.

Сравнение схем по локальным (при t = T) ошибкам

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	Б	НЛ	НЛ	ЯЛ	ПП	Б
L1	Б	ПП	НЛ, ЯЛ	НЛ, ЯЛ	Б	Б
L2	Б	НЛ, ЯЛ	НЛ	ЯЛ	ПП	Б

1) При г = 0.1 в нормах С и L2 для начальных условий №№ 1, 6 лучшие результаты показала предложенная в данной работе схема К.И. Бабенко, для условий №№ 2, 3 - неявная схема с левой разностью и для условия № 5 - метод прыжкового переноса. В норме L1 лучшей оказалась монотонизованная схема К.И. Бабенко в случае начальных условий №№ 1, 5, 6, метод прыжкового переноса - для начального профиля № 2 и схемы с направленными разностями - для условий №№ 3, 4.

г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	ЛВ	НЛ	НЛ	ЯЛ	ПП	Б
L1	Б	ПП	ПП	ЯЛ	Б	Б
L2	Б	НЛ	НЛ	ЯЛ	ПП	Б

2) При г = 0.25 в норме С для начального условия № 1 лучшие результаты получены с использованием схемы Lax-Wendroff'a, неявной схемы с левой разностью - для условий №№ 2, 3, метода прыжкового переноса - для условия № 5 и монотонизованной схемы К.И. Бабенко - для условия № 6. В норме L1 лучшей оказалась схема К.И. Бабенко для начальных условий №№ 1, 5, 6, метод прыжкового переноса - для условий №№ 2, 3 и явная схема с левой разностью - для условия № 4. В норме L2, как и при г = 0.1, для начальных условий №№ 1, 6 лучшие результаты показала схема К.И. Бабенко, для условий №№ 2, 3 - неявная схема с левой разностью, для условия № 4 - явная схема с левой разностью и для условия № 5 - метод прыжкового переноса.



г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	ЛВ	НЛ	ЛВ	ЯЛ	ПП	Б
L1	Б	ЯЛ, НЛ	ЛВ	ЯЛ	Б	Б
L2	Б	НЛ	НЛ, ЛВ	ЯЛ	Б, ПП	Б

3) При г = 0.5 в норме С для начальных условий №№ 1, 3 лучшей является схема Lax-Wendroff'a, неявная и явная схемы с левой разностью - для условий №№ 2 и 4 соответственно, метод прыжкового переноса - для условия № 5 и монотонизованная схема К.И. Бабенко - для условия № 6. В норме L1 лучшей оказалась схема К.И. Бабенко для начальных условий №№ 1, 5, 6, явная схема с левой разностью - для условия №№ 2, 4 и схема Lax-Wendroff'a - для условия № 3. В норме L2 для начальных условий №№ 1, 5, 6 лучшие результаты показала схема К.И. Бабенко, для условий №№ 2, 3 - неявная схема с левой разностью и для условия № 4 - явная схема с левой разностью.

г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	Б	НЛ	НЛ	ЯЛ, НЛ	Б	Б, ЯЛ
L1	Б	ЯЛ, НЛ	ЯЛ, НЛ	ЯЛ	Б	Б
L2	Б	ЯЛ, НЛ	НЛ	ЯЛ, НЛ	Б	Б

4) При г = 0.9 во всех нормах для начальных профилей №№ 1, 5, 6 минимальную ошибку дает схема К.И. Бабенко, а для условий №№ 2, 3, 4 - схемы с направленной разностью.

Сравнение схем по интегральным (по времени) ошибкам

г = 0.1	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	ЯЛ, НЛ	ЯЛ, НЛ	НЛ	Б	НЛ	НЛ
L1	Б	ПП	ЯЛ	Б	Б, ПП	ПП
L2	Б	ЯЛ, НЛ	НЛ	ЯЛ	Б, ПП	НЛ

1) При г = 0.1 в норме С для начального условия № 4 лучшие результаты показала монотонизованная схема К.И. Бабенко, для остальных начальных условий - неявная схема с левой разностью. В норме L1 лучшей оказалась схема К.И. Бабенко в случае начальных условий №№ 1, 4, 5, метод прыжкового переноса - для начальных профилей №№ 2, 5, 6 и схема с левой разностью - для условия № 3. В норме L2 минимальная ошибка получена с использованием схемы К.И. Бабенко для условий №№ 1, 5, явной и неявной схем с левой разностью - для условий №№ 2, 4 и № 3 соответственно.

г = 0.25	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	ЯЛ	ЛВ	НЛ	Б	НЛ	ЛВ
L1	Б	ПП	ПП	Б	Б	ПП
L2	Б	ЯЛ	ПП	ЯЛ	Б	НЛ

2) При г = 0.25 в норме С для начального условия № 4 лучшие результаты дает монотонизованная схема К.И. Бабенко, для начальных условий №№ 3, 5 - неявная схема с левой разностью, для условий №№ 2, 6 - схема Lax-Wendroff'a и для условия № 1 - явная схема с левой разностью. В норме L1 лучшей оказалась схема К.И. Бабенко в случае начальных условий №№ 1, 4, 5, метод прыжкового переноса - для остальных начальных профилей. В норме L2 минимальную ошибку дает схема К.И. Бабенко для условий №№ 1, 5, явная схема с левой разностью - для условий №№ 2, 4, неявная схема с левой разностью - для условия № 6 и метод прыжкового переноса - для условия № 3.

г = 0.5	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	ЯЛ	ЛВ	НЛ	Б	НЛ	ЛВ
L1	Б	ЯЛ	Б	Б	Б	ПП
L2	Б	ЯЛ	НЛ	ЯЛ	Б	ЛВ

3) При г = 0.5, как и при г = 0.25, в норме С для начального условия № 4 наименьшую ошибку дает монотонизованная схема К.И. Бабенко, для начальных условий №№ 3, 5 - неявная схема с левой разностью, для условий №№ 2, 6 - схема Lax-Wendroff'a и для условия № 1 - явная схема с левой разностью. В норме L1 лучшей оказалась схема К.И. Бабенко для начальных условий №№ 1, 3, 4, 5, метод прыжкового переноса - для условия № 6 и явная схема с левой разностью - для условий № 2. В норме L2 лучшие результаты показывает схема К.И. Бабенко для условий №№ 1, 5, явная схема с левой разностью - для условий №№ 2, 4, неявная схема с левой разностью - для условия № 3 и схема Lax-Wendroff'a - для условия № 6.

г = 0.9	треугольник	прямоуголь-ник	левый треугольник	правый треугольник	ступенька вверх	ступенька вниз
C	Б	НЛ	НЛ	Б	НЛ	НЛ
L1	Б	ЯЛ	Б, ЯЛ	ЯЛ	Б	Б
L2	Б	ЯЛ	НЛ	ЯЛ	Б	Б

4) При г = 0.9 в норме С для начальных условий №№ 1, 4 лучшая схема - схема К.И. Бабенко, для остальных начальных условий - неявная схема с левой разностью. В норме L1 лучшие результаты дает схема К.И. Бабенко для начальных условий №№ 1, 3, 5, 6, явная схема с левой разностью - для условий №№ 2, 4. В норме L2 лучшей оказалась схема К.И. Бабенко для условий №№ 1, 5, 6, явная схема с левой разностью - для условий №№ 2, 4, неявная схема с левой разностью - для условия № 3.

Заключение

Предложенная авторами в данной работе нелинейная монотонизованная схема К.И. Бабенко более высокое качество решения квазилинейного уравнения переноса во всем рассмотренном диапазоне чисел Куранта. Несомненным достоинством этой схемы является способность передавать разрывы решения с их размазыванием на наименьшее количество шагов по сравнению с другими рассмотренными схемами. Отметим также, что предложенный в [9] метод прыжкового переноса, несмотря на его способность точно воспроизводить численное решение, для некоторых начальных условий при определенных числах Куранта дает более низкое качество решения, что, в частности, связано с появлением ступенек, приводящих к сильному искажению решения.

Литература

Тихонов А. Н., Самарский А. А. Уравнения математической физики. - М.: Наука, 1972 г., 736 с.

Калиткин Н. Н. Численные методы. - М., Наука, Физматлит, 1978 г., 512 с.

Головизнин В. М., Карабасов С. А. Нелинейная коррекция схемы “Кабаре”. // Мат. Моделирование, № 12 (1998 г.). с. 107-123.

Галанин М. П., Еленина Т. Г. Сравнительный анализ разностных схем для линейного уравнения переноса. // Препринт ИПМ им. М. В. Келдыша РАН, 1998 г., № 52, 33 с.

Галанин М. П., Еленина Т. Г. Нелинейная монотонизация разностных схем для линейного уравнения переноса. // Препринт ИПМ им. М. В. Келдыша РАН, 1999 г., № 44, 33 с.

Головизнин В. М., Карабасов С. А. Метод прыжкового переноса для численного решения гиперболических уравнений. Точный алгоритм для моделирования конвекции на эйлеровых сетках. // Препринт ИБРАЭ № IBRAE-2000-04, Москва, 2002 г., 40 с.

Галанин М. П., Еленина Т. Г. Нелинейная монотонизация схемы К. И. Бабенко (“квадрат”) для уравнения переноса. // Препринт ИПМ им. М. В. Келдыша РАН, 2002 г., № 4, 26 с.

Рождественский Б. Л., Яненко Н.Н. Системы квазилинейных уравнений. - М.: Наука. Физматлит,1968 г.

Бабенко К. И., Воскресенский Г. П. Численный метод расчета пространственного обтекания тел сверхзвуковым потоком газа. // ЖВМиМФ, 1961 г., т. 1, № 6, с. 1051-1060.

Самарский А. А. Теория разностных схем. - М.: Наука, Физматлит, 1977 г., 656 с.

Федоренко Р. П. Введение в вычислительную физику. - М.: Изд-во МФТИ, 1994 г., 528 с..

Галанин М. П., Еленина Т. Г. Тестирование разностных схем для линейного уравнения переноса. // Препринт ИПМ им. М. В. Келдыша РАН, 1999 г., № 40, 42 с.

Дж. Уизем. Нелинейные волны. - М.: Мир, 1975 г.

Самарский А. А., Попов Ю. П. Разностные методы решения задач газовой динамики. - М.: Наука. Физматлит, 1980 г., 352 с.

Размещено на Allbest.ru