dt дх

где f = 0.52. Ранее к линейному аналогу уравнения (10.8) в форме (9.2) применялась схема Лакса - Вендроффа (9.16). Ввиду нелинейности функции F замена дй/дР на эквивалентную производную будет значительно сложнее, чем это было в. п. 9.1.3. Так, из уравнения (10.8) следует

дЧ dF д г dF\ dF . дй . dF

~ dxKdt)dt~dt~ дх

dt дх dt дх \ dt J dt dt dx

где А = дР/дй. В случае уравнения Бюргерса А = й. Использование вышеприведенных соотношений дает

Если ввести подходящую форму дискретизации для д{АдР/дх)/дх, то при применении схемы Лакса - Вендроффа к уравнению (10.8) получим

+ {-шУ - - iMr-Fix)], (10.10)

(10.5) на оператор LxF, который в случае положительного и имеет форму

Ff + l-Pf-l (/-2-3F/-1 + 3 --/+i)

8 случае отрицательного и выражение (10.6) заменяется на

,4). Pfi-Pj-x , (M-3F/ + 3,-F,.,)

Формулы (10.6) и (10.7) имеют ошибки аппроксимации порядка 0(Ajc2) независимо от выбора за исключением случая

9 = 0.5, когда эти ошибки имеют порядок ©(Ал:). Для рассматриваемой в п. 10.1.4 задачи о распространении ударной волны явная четырехточечная схема со сдвигом вверх по потоку, соответствующая формулам (10.6) и (10.7), включена как составная часть в программу BURG (рис. 10.4) для случая ME = 3,

Уравнения, определяющие неустановившиеся одномерные течения невязкой жидкости с ударными волнами (§ 10.2 и 14.2), по своей структуре аналогичны консервативной форме невязкого уравнения Бюргерса, т. е.

+ ? = 0, (10.8)

где в применении к уравнению Бюргерса Л/+1/2 = Uj+1/2 == = 0.5{Uf + Uf+i) и т. д. Схема (10.10) обладает ошибкой ап-лроксимации порядка 0{М, Ах) и устойчива при

\UmaxAt/Ax\ 1.0.

С учетом явного появления якобиана А после полушагов J- 1/2 и / + 1/2 схема (10.10) является в вычислительном плане менее экономичной, чем эквивалентный алгоритм (9.16) для линейного уравнения конвекции. Более экономичный метод вычисления получается в результате замены схемы (10.10) на эквивалентный ей двухэтапный алгоритм

;+i/2 = 0.5( ;? + ,) - o.SiFl+i-Fj), (10.11) Г = / - (FUm - FU,2). (10.12)

Формула (10.11) дает возможность получить промежуточное решение t/y i/2 обычно на временном слое (n-f 1/2)Д, позволяющее вычислить F*i2 и т. д., подставляя результат в правую часть формулы (10.12). Двухэтапная схема (10.11), (10.12) эквивалентна схеме (9.16), если F = u. Для упрощения программирования значения F при / -f 1/2 вводятся в память на место F{j), Это сделано в программе BURG (рис. 10.4) для варианта МЕ = 2.

В историческом плане [Richtmyer, Morten, 1967] двухэтапный алгоритм (10.11), (10.12) оказался весьма успешным для предсказания поведения невязкого сжимаемого потока. Попытки распространить его на изучение потока с вязкостью, сделанные, например, в работе [Thommen, 1966], были не столь удачными из-за формального снижения точности до порядка О (At, Ах). В итоге результаты для стационарного состояния были достаточно точными, но шаг по времени подвергался ограничению, связанному с условием устойчивости Куранта - Фридрихса - Леви. При решении нестационарных задач шаг по времени необходимо было подвергнуть еще более суровому ограничению, чтобы добиться приемлемой точности. Обобщение алгоритма (10.11), (10.12), сделанное Томменом и ставящее своей целью дискретизацию уравнения (10.3), можно представить в форме

Mi+ri2 = oMul + иГ) - 0.5 (FUi - Fl) +

+ 0.5s [0.5 ( ; , - 2н? + И/\,) + 0.5 (и; - 2ы + / +2)]. (ЮЛЗ)

где S = vAt/Ax. Эта схема включена в программу BURG (ME = 2). Из рассмотрения формулы (10.13) видно, что расчет члена с вязкостью охватывает четыре узловые точки; это вызывает необходимость задания добавочных граничных условий для вычисления г/у+г. Схема (10.13), (10.14) приводит к устойчивым решениям, если

Л/(Л2Д/ + 2г)<Лх2. (10.15)

Как указывается в книге [Peyret, Taylor, 1983], явный и практически приемлемый критерий получения устойчивого решения имеет вид АхУ {2v + \А\ Ах).

Совокупность формул (10.13), (10.14) можно интерпретировать в качестве одного из членов семейства Sp, у введенного в работе [Lerat, Peyret, 1975]. Это семейство методов обсуждается в книге [Peyret, Taylor, 1983]. Родственное семейство Sp для расчетов по невязкому уравнению Бюргерса обсуждается в п. 14.2.2. Параметр у определяет пропорцию распределения при вычислении члена с вязкостью, который в формуле (10.13) центрируется в узловых точках / и /+ 1. В случае схемы Томмена а = р = у = 0.5.

Общая классификация разностных схем для решения невязкого уравнения Бюргерса (10.2), включающая в себя и семейство Sp, дается в работе [Yanenko et al., 1983]. Эта классификация охватывает также и неявные схемы, причем относится, в частности, к уравнениям Эйлера для одномерного-неустановившегося течения (10.40), (10.41). Яненко и его соавторы описывают и дают ссылки на работы в отношении многих схем советского происхождения, малоизвестных западным авторам.

10.1,3, Неявные схемы

Применение неявных схем к нелинейным уравнениям типа уравнения Бюргерса не является столь же непосредственным как в случае линейных уравнений (гл. 9). Вариант неявной схемы Кранка - Николсона в применении к уравнению (10.3) принимает вид

Х = - 0.5L. (Ff + Fr) + 0.5vL.. {и } + (Ю. 16)

где Аи+ = и+-и, L = i-U О, 1)/2Ах и L = (l, -2, 1)/Ах Чтобы воспользоваться чрезвычайно эффективным алгоритмом Томаса (п. 6.2.2), необходимо свести соотношение (10.16) к трехдиагональной линейной системе уравнений относительно