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

Двухслойные разностные схемы. Начнем рассмотрение численных методов решения уравнений и систем гиперболического типа с простейшего случая одномерного уравнения переноса. Оно имеет вид (10.3):

Для определенности пока положим с > 0.

Корректная постановка начальной и краевой задач требует корректного задания начального и граничного условий. Они имеют вид.

Решением задачи Коши для уравнения (10.3) является бегущая волна:

где с — скорость переноса, а функция ЧДх — ct) определяется из начальных или граничных условий. Характеристики уравнения имеют вид x-ct = = const и при постоянной скорости переноса являются прямыми линиями. Решение однородного уравнения (10.3), как отмечено выше, остается постоянным вдоль характеристики, поэтому начальные и граничные условия переносятся вдоль этих линий. В случае неоднородного уравнения вдоль характеристики оно превращается в обыкновенное дифференциальное уравнение.

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

Приведем вид двухпараметрического семейства двухслойных разностных схем 1-го и 2-го порядков аппроксимации¹:

где введены следующие обозначения:

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

Схема Лакса (трехточечная схема) получается при у = 0, q = 1. Ее порядок аппроксимации 0(т + h^[1][2]т). Здесь и далее порядок аппроксимации приводится для модельного уравнения переноса с постоянными коэффициентами. В случае переменных коэффициентов в схему надо внести необходимые изменения. При этом естественно линейное уравнение рассматривать как частный случай квазилинейного (см. далее).

Схема является условно устойчивой, т. е. она устойчива при выполнении условия Куранта, а = cx/h < 1. Отметим здесь, что величина, а играет определяющую роль при исследовании разностных схем на аппроксимацию и устойчивость. Она называется числом Куранта. Исследование разностной схемы на устойчивость для линейного эволюционного уравнения с постоянными коэффициентами можно провести с использованием спектрального признака (Неймана)^[2].

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

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

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

Схема Куранта — Изаксона — Риса (КИР) получается при у = О,.

q = asignc. Ее порядок аппроксимации 0(т + К). Схема КИР условно устойчива, т. е. она также устойчива при выполнении условия Куранта a = = ст/h < 1. Приведем разностные уравнения для схемы КИР во внутренних точках расчетной области:

Схему КИР в соответствии с формой шаблона иногда называют схемой «левый уголок» (для с > 0) или схемой «правый уголок» (для с < 0).

Рассмотренные схемы КИР очевидным способом обобщаются на случай квазилинейного уравнения.

Схема КИР тесно связана с численными методами характеристик. Дадим краткое описание идеи таких методов.

Две последние полученные схемы (при разных знаках скорости переноса) можно интерпретировать следующим образом. Построим характеристику, проходящую через узел (t_n+{, х_т), значение в котором нам необходимо определить, и пересекающую слой t_n в точке х? = х_т — ст. Для определенности считаем, что скорость переноса с положительна.

Проведя линейную интерполяцию между узлами х_т__х и х_т на нижнем слое по времени, получим.

Далее перенесем вдоль характеристики значение и_п(х?) без изменения на верхний слой t_n+{, т. е. положим м" ⁺¹ = и_п(х'). Последнее значение естественно считать приближенным решением однородного уравнения переноса. Тогда.

или, возвращаясь от числа Куранта к сеточным параметрам,.

т.е. мы пришли другим способом к уже известной схеме «левый уголок», устойчивой, если a = cz/h 0. При a > 1 точка пересечения характеристики, выходящей из узла (?_я+1, х_т), с п-м слоем, но времени расположена левее узла (?_я, х_т__{). Таким образом, для отыскания решения г/"⁺¹ используется уже не интерполяционный, а экстраполяционный процесс, который оказывается неустойчивым.

Неустойчивость схемы «правый уголок» при с > 0 также очевидна Для доказательства этого можно использовать либо спектральный признак, либо условие Куранта, Фридрихса и Леви^[4]. Аналогичные рассуждения можно провести и для случая с < 0 и схемы «правый уголок».

Неустойчивая четырехточечная схема получается при у = 0, q = 1, ее порядок аппроксимации 0(т + И^[4]). Сеточные уравнения для разностной схемы имеют вид.

Схема Лакса — Вендроффа возникает при у = 0, q = а. Порядок аппроксимации схемы Лакса — Вендроффа есть 0(т² + /г²). Схема устойчива при выполнении условия Куранта, а = cx/h < 1.

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

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

которое получается путем дифференцирования исходного уравнения (10.3) сначала по времени t, затем по координате х и вычитания одно из другого получившихся равенств — следствий исходного дифференциального уравнения.

Далее, заменяя вторую производную во втором слагаемом в правой части с точностью до 0(/г²), получим новую разностную схему, аппроксимирующую исходное дифференциальное уравнение с точностью 0(т² + /г²). Сеточные уравнения для схемы Лакса — Вендроффа во внутренних узлах расчетных сеток есть.

Неявная шеститочечная схема возникает при q = 0; Если у = 0,5, ее порядок аппроксимации 0(т² + Л²), если у = 1 — 0(т + И^[4]).

Неявная нецентральная схема. При q = asignc и у = 0,5 порядок аппроксимации — 0(х² + h), при у = 1 — 0(т + /?).

Спектральный признак устойчивости для эволюционных уравнений.

Для широкого класса эволюционных (зависящих от времени) задач исследование устойчивости можно осуществить с помощью спектрального признака, который в случае разностной задачи с постоянными коэффициентами состоит в следующем. Заменяем правую часть разностного уравнения нулем, краевую задачу — задачей Коши, функцию ф", — гармоникой е'^шт и ищем решение для задач с одной пространственной переменной в виде.

где о — произвольное число, 0 < со < 2л.

Для устойчивости разностной схемы необходимо, чтобы спектр X = л (со) лежал в круге |Л| < 1 + ст, где с не зависит от шага интегрирования по времени т.

Обоснование спектрального принципа следующее. Рассматривается линейное уравнение эволюционного типа.

где А — дифференциальный оператор с постоянными коэффициентами. Для этого уравнения решается задача Коши: при? = 0 задается и₀(х) и ищется решение «(?, х),? > 0.

Пусть ср (?, k) — преобразование Фурье функции м (?, х), а (/е) — частотная характеристика оператора А. Беря преобразование Фурье от правой и левой частей уравнения, получим уравнение.

решение которого есть.

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

где F_t — оператор (типа свертки) с частотной характеристикой р (?) = е"(*Х Задача Коши для рассматриваемого уравнения будет поставлена корректно в том и только в том случае, когда р (?) ограничено при всех? > 0.

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

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

коэффициентов см. ниже. Коэффициенты q_m зависят от шагов разностной сетки т и h.

Пусть, кроме того,.

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

Запись разностной схемы в приведенном виде задает разностный оператор послойного перехода

Тогда мы имеем.

Из этого равенства следует ограниченность разностного оператора G. Но для функции и" (х) можно записать.

Поскольку || О || < ||G'||^W, оператор, переводящий щ (рс) в и^п_} ограничен. Однако для устойчивости разностной задачи необходимо, чтобы метод выдерживал счет со сколь угодно малыми шагами сетки, т. е. при т, /г —> 0.

При т —" 0 величина п = Г/т стремится к бесконечности. Если при этом норма Gⁿ неограниченно возрастает, то сколь угодно малые ошибки в задании и₀(х) могут привести к сколь угодно большим возмущениям функции uⁿ(x)_f т. е. возникает вычислительная неустойчивость.

Устойчивость задачи Коши означает, что непрерывная зависимость и^п{х) от и₀(х) равномерна по т, h.

Считаем, что фиксирована некоторая функция h = h (т), такая что выполняется следующее условие:

Тогда q_m зависят лишь от т, a G = G_r Для рассматриваемого уравнения переноса такая связь очевидна, если считать число Куранта фиксированным.

Положим п = Т/т.

Разностное уравнение устойчиво, если || G" || < М, где М не зависит от п (но может зависеть от 7).

Достаточное условие устойчивости (Неймана) есть ||G_T|| < 1 + ст, где с не зависит от сеточных параметров. Покажем, что в этом случае разностная схема удовлетворяет определению устойчивой схемы.

Действительно, тогда.

Более сильное требование ||G_T|| < 1 иногда называют условием строгой устойчивости.

Строго устойчивые схемы могут существовать лишь для дифференциальных уравнений, решения которых подчиняются принципу максимума:

В пространстве L² вычисление нормы оператора основано на равенстве.

где X — частотная характеристика оператора. Она ищется следующим образом. Подставляем решение в специальном виде и^п(х) = е^{, кх} в запись разностной схемы, тогда и^п+1(х) = Х (к)е'^кх. После преобразований получаем.

Такая форма записи удобна при исследовании на устойчивость разностных схем в пространстве неопределенных коэффициентов. Обозначив hk = ф, получаем следующее условие устойчивости, основанное на вычислении спектра оператора послойного перехода л (ф):

для всех значений ф, 0 < ф < 2тт (условие Неймана), или.

Подробнее о спектральном признаке устойчивости см. работу В. С. Рябенького^[7].

• [1] Шокип Ю. И., Яненко Н. Н. Метод дифференциального приближения. Применениек газовой динамике. Новоеибирск: Наука, 1985.
• [2] Подробнее см. ниже и в литературе, например: Рябенький В. С. Метод разностныхпотенциалов для некоторых задач механики сплошной среды.
• [3] Подробнее см. ниже и в литературе, например: Рябенький В. С. Метод разностныхпотенциалов для некоторых задач механики сплошной среды.
• [4] См. работу: Рябенький В. С. Метод разностных потенциалов для некоторых задач меха ники сплошной среды.
• [5] См. работу: Рябенький В. С. Метод разностных потенциалов для некоторых задач меха ники сплошной среды.
• [6] См. работу: Рябенький В. С. Метод разностных потенциалов для некоторых задач меха ники сплошной среды.
• [7] Рябенький В. С. Метод разностных потенциалов для некоторых задач механики сплошной среды.