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

Рассмотрим задачу Коши для уравнения вида

в которой скорость переноса v может быть функцией х. Для уравнения (6.1) можно предложить множество разностных схем, различающихся порядком аппроксимации, способом представления производных и т.д. Остановимся вначале на явных разностных схемах, в которых каждое уравнение системы содержит лишь одну неизвестную величин}", что позволяет последовательно вычислить значения решения на новом временном слое.

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

имеет очень жесткое условие устойчивости (т 2 vh) и нс используется для практических алгоритмов. Разностные схемы


являются условно устойчивыми. Для обеспечения их устойчивости необходимо, во-первых, выполнение условия Куранта Фридрихса - Леви (КФЛ):

а во-вторых, использование разностей навстречу потоку, т.е. применение схемы (6.3) при V > 0 и (6.4) при v 0.

Явная схема с разностями навстречу потоку. Если мы будем избирательно применять две предыдущие схемы, а именно, при v > > 0 схему (6.3), а при v

будет безразлична к направлению скорости и устойчива при условии v /h ^ 1. Нетрудно заметить, что односторонние разности в этой схеме берутся навстречу потоку (говорят, что схема обладает свойством mpanenopmuenoemu). Схем}" такого вида называют противопоточной или схемой с разностями навстречу потоку.

В случае уравнения с постоянным значением скорости переноса проблем с конструированием противопоточной разностной схемы нет. Выбирается соответствующая знаку скорости переноса разность, которая используется во всех узлах расчетной области. Условие (6.5) накладывает ограничение на соотношение шагов расчетной сетки. Обычно при заданном шаге по пространству из соотношения (6.5) определяют допустимый временной шаг т h/v.

Но если скорость переноса является функцией координаты (или времени), то выбор вида разностной аппроксимации необходимо осуществлять на основе анализа знака скорости переноса, например применяя условный оператор. Кроме гот, при переменной скорости переноса v = v(x) условие устойчивости нужно проверить для всех узлов сетки и из этого множества значений временного шага выбрать минимальный: т min,; h/vj.

В работе Куранта с соавторами (1952) был предложен интересный метод конструирования противопоточной схемы, в котором не использовался условный оператор. Важно отметить, что это не просто формальный прием, а подход, содержащий глубокие идеи, на основе которых можно сопоставлять и находить соответствие между противопоточными (несимметричными) и симметричными разностными схемами. К этому близка идея расщепления операторов разностных схем.

Представим скорость переноса в виде суммы ее положительной и отрицательной составляющих:

Это позволит представить оператор переноса в виде суммы двух операторов:

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

Если мы теперь в правой части (6.7) проведем элементарные преобразования и выделим симметричную разностную производную, то эта схема представится в виде

Можно сделать вывод, что нротивопоточная разностная схема (6.7) эквивалентна симметричной (6.2), в которую введена диссипативная добавка, обеспечивающая условную устойчивость схемы.

Схема Лакса. Эта схема была введена в практику вычислений на заре развития вычислительной газодинамики. II хотя упоминания о схеме такого типа встречались в работах разных авторов, общественное мнение связывает ее с именем американского математика Лакса (Lax, P.D.), опубликовавшего в 50-е годы серию работ по различным аспектам теории разностных схем. Применительно к уравнению переноса (6.1) эта схема имеет вид

Особенность схемы состоит в том, что для обеспечения ее устойчивости в аппроксимации производной но времени значение сеточной функции в узле (г, п) заменяется на полусумму значений в соседних узлах того же временного слоя. Эта операция обеспечивает при центральной аппроксимации пространственной производной условную устойчивость разностной схемы (при выполнении условия Куранта - Фридрихса - Леви v /h ^ 1).

Хотя здесь производная по х представлена со вторым порядком аппроксимации, схема из-за специфического представления вре- меннбй производной обладает значительной диссипацией. Это хорошо видно из первого дифференциального приближения:

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

где через а обозначено число Куранта. Из дифференциального приближения можно определить многие свойства этой схемы:

  • - схема становится недиссипативной при числе Куранта, равном единице;
  • - схема не чувствительна к направлению потока;

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

При уменьшении шага по времени диссипативные свойства схемы растут.

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

Схемы второго порядка. Разностные схемы, рассмотренные ранее, были схемам первого порядка (по пространственной или временной переменной). При построении схем второго порядка необходимо обеспечить повышенный порядок аппроксимации как по пространственной, гак и но временной переменой. Рассмотрим несколько схем такого типа.

Схема «чехарда». Схема второго порядка как по пространственной переменной, так и по времени простейшего типа может быть представлена в виде

Эту схему называют схемой с перешагиванием, но больше она известна под названием «чехарда» (leap-frog scheme). Схема является трехслойиой и строит решение по двум предыдущим временным слоям. Поэтому при ее применении возникают проблемы с началом вычислений, которое должно осуществляться каким либо другим методом.

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

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

которую будем рассматривать совместно с исходным уравнением (6.1) Это уравнение будем использовать для того, чтобы заменить в разложении временные производные пространственными. Это возможно, так как первая производная но времени выражается непосредственно из (6.1): du/dt = -vdu/dx. Вторая производная также легко находится из следующей цепочки соотношений:

Заметим, что это представление является точным лишь при постоянной скорости переноса: v = const. В противном случае оно носит приближенный характер, однако, если скорость переноса v(x) достаточно гладкая функция, его можно использовать для преобразований разностных соотношений, которые носят локальный характер.

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

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

называемую схемой Лакса Вендроффа. Эта схема была введена в практику вычислений вместе с рядом других в серии работ, опубликованных Лаксом и Всндроффом в 1960-1964 гг.

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

На первом полушаге вычислим промежуточное значение решения по простой схеме Лакса первого порядка. Этому промежуточному значению припишем верхний индекс п + 1/2 и будем иметь в виду, что используется также половинный шаг по времени. Применяя эту схему, получим значения решения на промежуточном временном слое: t = t n+l / 2 . При этом отметим, что из-за применения схемы Лакса, в которой на нижнем слое отсутствует центральный узел, решение воспроизводится на промежуточном слое также в системе полуцелых точек.

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


Второй полушаг состоит в вычислении решения на новом временном слое п + 1 на основе схемы с центральными разностями как по пространству, так и по времени - схемы «крест». Для вычисления пространственных производных используются значения решения на промежуточном слое в системе полуцелых точек, само решение восстанавливается в той же самой системе точек, в которой оно было определено к началу временного шага:

Соотношения (6.12) и (6.13) вместе определяют двухшаговую схему Лакса - Веидроффа. На первом ее этапе обеспечивается выполнение условий устойчивости. Этот этап называют иногда предиктором. Второй этап обеспечивает выполнение требуемой точности, и его называют корректором. Методы предиктора-корректора часто используются в вычислительной математике, при этом этап корректора может включать в себя итерационный блок.

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

Монотонность решения в схемах второго порядка. Последний член в правой части (6.11) имеет вид, отличный от вида диссипативных членов схем первого порядка (6.8) и (6.10). В данном случае он обеспечивает подавление ошибки, связанной с первым порядком аппроксимации производной по времени. Таким образом, данная схема является схемой второго порядка как но времени, гак и но пространственной переменной. Ее первое дифференциальное приближение уже не будет содержать диссипативный член, но в нем будет присутствовать дисперсионная составляющая с третьей производной, которая является причиной фазовых ошибок схемы. Можно ожидать, что данная схема будет слабо размазывать решение, но в области его резкого изменения могут появляться нефизические осцилляции, вызванные дисперсией.

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

С.К. Годуновым была установлена теорема о монотонности, занимающая одно из центральных мест в теории разностных схем. Согласно этой теореме, для линейного уравнения вида (6.1) не существует монотонных схем с порядком выше первого.

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

Схема Мак-Кормака. Это также двухшаговая схема второго порядка, безразличная к направлению потока. Ее более удобно продемонстрировать на консервативной форме уравнения переноса:

Схема состоит из двух последовательно выполняемых шагов:


На первом этапе (6.15) находят предварительное значение решения щ в узлах сетки на основе односторонней разностной схемы. По найденному таким образом решению вычисляют предварительные значения потоков / г. Далее, на основе односторонних схем, имеющих противоположное направление (6.16), определяется решение на следующем временном слое.

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

Размер: px

Начинать показ со страницы:

Транскрипт

2 МИНИСТЕРСТВО ОБРАЗОВАНИЯ И НАУКИ РФ НОВОСИБИРСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ Механико-математический факультет Кафедра математического моделирования Г. С. Хакимзянов, С. Г. Черный МЕТОДЫ ВЫЧИСЛЕНИЙ Часть 4. Численные методы решения задач для уравнений гиперболического типа Учебное пособие Новосибирск 014

3 ББК В.193 УДК Х 16 Рецензент канд. физ.-мат. наук А. С. Лебедев Издание подготовлено в рамках реализации Программы развития государственного образовательного учреждения высшего профессионального образования «Новосибирский государственный университет» на годы. Х 16 Хакимзянов, Г. С. Методы вычислений: В 4 ч. : учеб. пособие / Г. С. Хакимзянов, С. Г. Черный; Новосиб. гос. ун-т. Новосибирск: РИЦ НГУ, 014. Ч. 4: Численные методы решения задач для уравнений гиперболического типа. 07 с. ISBN Учебное пособие соответствует программе курса лекций «Методы вычислений», который читается на механико-математическом факультете НГУ. В его четвертой части излагаются основы численных методов решения начально-краевых задач для уравнений гиперболического типа, формулируются задачи для семинарских занятий, приводятся образцы контрольных работ и заданий для практических занятий на ЭВМ. Пособие предназначено для студентов и преподавателей математических специальностей высших учебных заведений. ISBN ББК В.193 УДК c Новосибирский государственный университет, 014 c Г. С. Хакимзянов, С. Г. Черный, 014

4 ОГЛАВЛЕНИЕ Предисловие Схемы для линейного уравнения переноса Свойство монотонности разностных схем Построение монотонных схем на основе метода дифференциального приближения Схемы для нелинейного уравнения переноса Схемы на адаптивной сетке для уравнения переноса Разностные схемы для уравнения колебаний струны Разностные схемы для гиперболической системы уравнений с постоянными коэффициентами Разностные схемы для системы нелинейных уравнений мелкой воды Разностные схемы для задач газовой динамики Контрольная работа по теме «Исследование разностных схем для уравнения переноса» Задания для лабораторной работы Ответы, указания, решения Библиографический список

5 Предисловие В четвертой части пособия изложены основы численных методов решения начально-краевых задач для уравнений гиперболического типа, сформулированы задачи по этой теме для семинарских занятий, приведены задания для практических занятий на ЭВМ и пример контрольной работы. Теоретические вопросы изложены достаточно кратко. Для более глубокого изучения рассматриваемых вопросов мы рекомендуем обратиться к учебнику С. К. Годунова и В. С. Рябенького , а также к книгам Г. И. Марчука , А. А. Самарского , А. А. Самарского и А. В. Гулина , А. А. Самарского и Е. С. Николаева , Б. Л. Рождественского и Н. Н. Яненко и учебным пособиям, изданным в НГУ . На лекциях рассматриваются теоретические вопросы, связанные с исследованием только конечно-разностных схем. В качестве примеров рассмотрены схемы для линейного уравнения переноса, нелинейного скалярного уравнения первого порядка, уравнения второго порядка, описывающего колебания струны, линейной системы уравнений первого порядка, системы нелинейных уравнений мелкой воды и уравнений газовой динамики. Каждый параграф сопровождается задачами, которые необходимо решить на семинарских занятиях. Многие задачи снабжены указаниями и подробными решениями. Дополнительные материалы для семинарских занятий можно найти в задачниках . В пособии приведены примеры заданий для практических занятий в компьютерных классах, даны рекомендации по выполнению заданий, обсуждаются вопросы, связанные с разработкой программ и представлением результатов. Дополнительные задания можно взять из методических пособий . Четвертая часть пособия имеет самостоятельную сквозную нумерацию параграфов и рисунков и самостоятельный библиографический список. Внутри параграфов для формул и утверждений (лемм и теорем) использована двухиндексная нумерация, например 4.. Ссылки на формулы, леммы, теоремы из предыдущих трех частей пособия даются добавлением спереди к их номеру цифры 1, или 3. Например, вместо «по формуле (4.) из пособия » мы пишем «по формуле (1.4.)», вместо «по теореме 8.3 из пособия » «по теореме.8.3». Авторы выражают глубокую признательность рецензенту Александру Степановичу Лебедеву за ценные советы и критические замечания, которые способствовали улучшению этого учебного пособия. 4

6 1. Схемы для линейного уравнения переноса 1.1. Некоторые сведения из теории гиперболических систем. Рассмотрим задачу Коши для линейной системы дифференциальных уравнений первого порядка u t + A u = f(x, t), < x <, 0 < t T, x u(x, 0) = u 0 (x), < x <. (1.1) Здесь u = (u 1,..., u m) T m-мерная вектор-функция переменных x, t, A вещественная m m матрица с элементами a i (x, t). Определение. Систему уравнений (1.1) будем называть гиперболической в некоторой области переменных (x, t), если в каждой точке этой области собственные значения λ 1, λ,..., λ m матрицы A вещественны и различны. Определение. Интегральная кривая x = x k (t) обыкновенного дифференциального уравнения dx dt = λ k(x, t) (1.) называется k-ой характеристикой системы уравнений (1.1). Предполагается, что элементы матрицы A обладают гладкостью, достаточной для того, чтобы через каждую точку плоскости (x, t) проходила единственная характеристика, отвечающая собственному значению λ k. Характеристики, проведенные через точку (x, t) (t > 0) в сторону убывания времени t, пересекут ось Ox в m различных точках. Упорядочим собственные значения гиперболической системы (1.1) (λ 1 (x, t) < λ (x, t) <... < λ m (x, t)) и через обозначим отрезок оси Ox, ограниченный точками пересечения этой оси с m-ой и первой характеристиками. Определение. Областью зависимости точки (x, t) для системы уравнений (1.1) называется множество точек верхней полуплоскости, ограниченное крайними характеристиками x = x m (t), x = x 1 (t) и отрезком . Область зависимости точки (x, t) изображена на рис. 1, а. Решение u системы (1.1) в точке (x, t) будет зависеть только от значений u 0 (x) на 5

7 отрезке . Следовательно, если начальные данные вне отрезка поменять на другие, то решение в точке (x, t) не изменится. Определение. Областью влияния точки (x 0, 0) называется множество точек (x, t) верхней полуплоскости, ограниченное крайними характеристиками системы (1.1), выходящими из (x 0, 0), т. е. характеристиками, соответствующими собственным значениям λ 1 и λ m. Область влияния точки (x 0, 0) показана на рис. 1, б. Если начальные данные изменить лишь в точке (x 0, 0), то решение гиперболической системы изменится только в точках (x, t), принадлежащих области влияния точки (x 0, 0). Предположим теперь, что нам вместо задачи Коши (1.1) нужно решить начально-краевую задачу на отрезке . Тогда в дополнение к начальным условиям необходимо задать краевые условия. Количество краевых условий на каждой из границ определяется количеством входящих внутрь области характеристик. Например, если через левую границу x = 0 внутрь области входит m 0 характеристик, т. е. m 0 собственных значений λ k положительны при x = 0, то на этой границе надо задать m 0 краевых условий. Если на границе x = l количество отрицательных собственных значений равно m l и, следовательно, ровно m l характеристик входит в область через правую границу, то на этой границе необходимо задать m l краевых условий. Поскольку собственные значения зависят от времени, то количество краевых условий на каждой из границ может меняться со временем. t dx dt = m λ m (x,t) dx dt = λ 1 t dx dt =λ 1 dx dt = m λ x l а x r x (x 0,0) б x Рис. 1. Характеристики системы уравнений (1.1), ограничивающие области зависимости точки (x, t) (а) и влияния точки (x 0, 0) (б) 6

8 Рассмотрим теперь однородную гиперболическую систему уравнений (1.1) с постоянными коэффициентами. Для постоянной матрицы A ее собственные векторы и собственные значения являются постоянными, т. е. не зависят от x и t. Пусть l k k-й левый собственный вектор матрицы A, отвечающий ее собственному значению λ k: l k A = λ k l k (k = 1,..., m). Умножим систему (1.1) слева на вектор l k: или где l k u t + l ka u x = 0. Это уравнение можно записать в следующем виде: l k u t + λ k l k u x s k t + λ s k k x = 0, = 0, (1.3) s k = l k u, k = 1,... m. (1.4) Решение s k (x, t) уравнения (1.3) переносится вдоль характеристики без изменения и потому вычисляется при t > 0 по начальному значению s k в точке пересечения k-ой характеристики с осью Ox: s k (x, t) = s k (x λ k t, 0). (1.5) Функции s k называются инвариантами Римана. 1.. Линейная модель мелкой воды. Простейшей математической моделью, в рамках которой можно описывать движение жидкости с поверхностными волнами, является линейная модель мелкой воды: η t + u 0 = 0, (1.6) x u t + g η = 0, (1.7) x η(x, 0) = η 0 (x), u(x, 0) = u 0 (x), (1.8) где η(x, t) возвышение поверхности жидкости над невозмущенным уровнем (см. рис.), u(x, t) скорость жидкости, η 0 (x) и u 0 (x) возвышение и скорость в начальный момент времени t = 0, 0 = const глубина бассейна, g = const ускорение свободного падения. 7

9 Систему уравнений (1.6), (1.7) можно записать в виде однородной системы (1.1) с матрицей A и вектором решения u: A = (0 0 g 0) (η, u = u). (1.9) Матрица A имеет два различных действительных собственных значения λ 1 = c 0, λ = c 0 = g 0, (1.10) поэтому система уравнений (1.6), (1.7) имеет гиперболический тип. Уравнения характеристик (1.) принимают такой вид: dx dt = c 0, dx dt = c 0, (1.11) поэтому характеристики являются прямыми линиями. Характеристики, проходящие через точку (x, t), t > 0, пересекают ось Ox в точках x l и x r, где x l = x c 0 t, x r = x + c 0 t. (1.1) Левые собственные векторы матрицы A, соответствующие собственным значениям (1.10), задаются формулами l 1 = (c 0, 0), l = (c 0, 0). (1.13) y 0 η y= (x,t) l x y=- 0 Рис.. Обозначения в задаче о распространении и трансформации волн в бассейне с вертикальными стенками Согласно (1.4) связь между инвариантами Римана r = s 1, s = s и исходными зависимыми переменными задается формулами r = c 0 η 0 u, s = c 0 η + 0 u, (1.14) 8

10 откуда η = r + s c 0, u = s r 0. (1.15) Из формулы (1.5) с учетом равенств (1.14) получаем формулы для решения задачи Коши в инвариантах r(x, t) = r(x λ 1 t, 0) = r(x + c 0 t, 0) = c 0 η 0 (x r) 0 u 0 (x r), (1.16) s(x, t) = s(x λ t, 0) = s(x c 0 t, 0) = c 0 η 0 (x l) + 0 u 0 (x l). (1.17) И наконец, используя соотношения (1.15), получаем точное решение задачи Коши (1.6), (1.7), (1.8) η(x, t) = η 0(x l) + η 0 (x r) + 0 u0(x l) u 0 (x r), c 0 u(x, t) = u 0(x l) + u 0 (x r) + c 0 η0(x l) η 0 (x r). 0 (1.18) При решении рассматриваемой начально-краевой задачи необходимо поставить по одному условию на каждом из концов отрезка . Будем, например, считать, что стенки бассейна являются непроницаемыми для жидкости, что означает равенство нулю скорости жидкости на этих стенках: u(0, t) = u(l, t) = 0. (1.19) Приведем в окончательном виде математическую формулировку задачи о движении жидкости с поверхностными волнами в ограниченном бассейне: найти непрерывное в замкнутой области D = решение η(x, t), u(x, t) следующей начально-краевой задачи η t + u 0 x = 0, u t + g η = 0, 0 < x < l, 0 < t T, x u(0, t) = u(l, t) = 0, 0 t T, η(x, 0) = η 0 (x), u(x, 0) = u 0 (x), 0 x l. (1.0) 1.3. Линейное уравнение переноса. Итак, если матрица A однородной гиперболической системы уравнений (1.1) постоянна, то такую систему можно свести к системе уравнений в инвариантах Римана, 9

11 при этом уравнения для инвариантов Римана не зависят друг от друга и каждое из них имеет вид u t + au x = 0, a = const. (1.1) Это уравнение является простейшим гиперболическим уравнением и называется линейным уравнением переноса. На этом уравнении можно изучать свойства разностных схем, применяемых для решения гиперболических систем уравнений. Рассмотрим для линейного уравнения переноса (1.1) задачу Коши u t + au x = 0, < x <, 0 < t T, u(x, 0) = u 0 (x), < x <. (1.) Характеристика x = x(t) уравнения (1.1) определяется уравнением dx dt = a, (1.3) т. е. является прямой с наклоном a к оси Ot. Следовательно, точное решение задачи Коши определяется по формуле u(x, t) = u 0 (x at). (1.4) График точного решения в момент времени t получается переносом графика начальной функции на величину at (в положительном направлении оси Ox, если a > 0 и наоборот). Для уравнения переноса с постоянным коэффициентом a легко выписать точное решение и для начально-краевой задачи. Пусть, например, a = const > 0. Тогда корректной будет следующая начальнокраевая задача u t + au x = 0, 0 < x l, 0 < t T, u(0, t) = µ 0 (t), 0 t T, u(x, 0) = u 0 (x), 0 x l, u 0 (0) = µ 0 (0). (1.5) Легко проверить, что если u 0 (x) и µ 0 (t) дифференцируемые функции, то решение задачи (1.5) определяется формулой u(x, t) = { u0 (x at) при t x/a, µ 0 (t x/a) при t x/a. (1.6) 1.4. Явная противопоточная схема. Перейдем теперь к изучению конечно-разностных схем решения линейного уравнения переноса. 10

12 Начнем с явной схемы с направленными против потока разностями (противопоточная схема) для начально-краевой задачи u t + au x = f(x, t), 0 < x l, 0 < t T, a = const > 0, u(0, t) = µ 0 (t), 0 t T, u(x, 0) = u 0 (x), 0 x l, u 0 (0) = µ 0 (0). (1.7) Всюду далее будем рассматривать только равномерные сетки, покрывающие замкнутую область D = . Построим следующую разностную схему u n + a un un 1 = f n, = 1,..., N, u n 0 = µ n 0, n = 0,..., M, u 0 = u 0(x), = 0,..., N, (1.8) аппроксимирующую задачу (1.7) с порядком O(+). Как и ранее, эту схему можно записать в операторном виде L u = f. Название противопоточная схема связано с тем, что если мы рассматриваем уравнение переноса как модельное уравнение для системы уравнений, описывающих течение жидкости или газа, и отождествляем коэффициент a со скоростью жидкости, то при положительной скорости, т. е. при a > 0, в схеме берутся левые разностные производные, использующие узел x 1, расположенный против «потока» (расположенный вверх по потоку). Введем равномерные нормы в пространстве сеточных функций U и пространстве правых частей F: где f F (= max u U max n u n C = max 0 N un, = max n un C, (1.9)) µn 0, (u 0) C, max f n n C, (1.30) f n C = max 1 N f n равномерные нормы на слое t = t n. С помощью принципа максимума можно доказать следующее утверждение. Теорема 1.1. Выполнение условия a 1 (1.31) 11

13 достаточно для устойчивости противопоточной схемы (1.8) в равномерной норме. Д о к а з а т е л ь с т в о. Пусть x узел сетки с номером 1 N. Перепишем разностное уравнение схемы в этом узле = (1 r)u n + ru n 1 + f n, где r = a/. Из условия теоремы следует, что 1 r 0, поэтому будет справедливой следующая оценка (1 r) u n +r u n 1 + f n (1 r) u n C +r u n C + f n C u n C + max m f m C. В граничном узле имеем следующую оценку 0 = µ n+1 0 max m µm 0. Следовательно, максимальное из левых частей этих неравенств не может превзойти максимального из двух чисел в правых частях этих неравенств: (C max max m) µm 0, u n C + max f m m C, а это и есть принцип максимума. Получили, что при условии (1.31) схема (1.8) удовлетворяет принципу максимума. Поэтому (см. теорему 3.1.1) она будет устойчивой в равномерной норме по начальным данным, краевым условиям и по правой части. Это же условие (1.31) является и необходимым условием устойчивости схемы (1.8), что следует из спектрального признака устойчивости Неймана. Докажем это. Возьмем гармонику u n = λ n e iφ (1.3) и подставим ее в однородное разностное уравнение. В результате для множителя перехода получим уравнение Следовательно, λ = 1 r (1 e iφ) = 1 r(1 cos φ) ir sin φ. λ = 1 r(1 cos φ) + r (1 cos φ) + r sin φ = 1

14 = 1 r(1 cos φ) [ r(1 cos φ) r(1 + cos φ)] = 1 r(1 cos φ)(1 r). Пусть в схеме (1.8) шаги и связаны законом предельного перехода r = a = const. (1.33) Тогда собственные числа λ (φ) не зависят от, поэтому необходимое условие устойчивости Неймана сводится к требованию или λ (φ) 1, φ R. (1.34) r(1 cos φ)(1 r) 0, φ R. (1.35) Очевидно, что это неравенство эквивалентно при a > 0 условию (1.31). Итак, условие (1.31) при a > 0 является необходимым и достаточным условием устойчивости противопоточной схемы в равномерной норме. Отметим, что при a < 0 схема (1.8) абсолютно неустойчива, поскольку в этом случае нарушается неравенство (1.34) (см. задачу 1.1). Какую же схему следует использовать при a < 0, когда поток распространяется справа налево? Отметим, что в этом случае корректной будет такая начально-краевая задача u t + au x = f(x, t), 0 x < l, 0 < t T, a = const < 0, u(l, t) = µ l (t), 0 t T, u(x, 0) = u 0 (x), 0 x l, u 0 (l) = µ l (l). (1.36) Для этой задачи возьмем следующую противопоточную схему u n + a un +1 un = f n, = 0,..., N 1, u n N = µn l, n = 0,..., M, u 0 = u 0(x), = 0,..., N, (1.37) которая аппроксимирует дифференциальную задачу (1.36) с порядком O(+). Используя принцип максимума и спектральный признак Неймана, можно показать, что схема (1.37) при a < 0 будет устойчива при выполнении условия a 1. С другой стороны, при a > 0 схема (1.37) будет абсолютно неустойчивой (см. задачу 1.). 13

15 Таким образом, мы построили две условно устойчивые явные схемы с направленными против потока разностями для уравнения переноса с постоянным коэффициентом a u n u n + a un un 1 + a un +1 un Они устойчивы при выполнении неравенства = f n, если a > 0, = f n, если a < 0. (1.38) a 1. (1.39) Во внутренних узлах сетки противопоточную схему (1.38) можно записать в виде одного уравнения u n + a + a u n un 1 + a a u n +1 un = f n. (1.40) Аналогично выглядит явная противопоточная схема и в случае знакопеременного коэффициента a(x, t). Например, если на границах отрезка выполнены условия a(0, t) > 0, a(l, t) < 0, 0 t T, то получим такую противопоточную схему где u n + a + un un 1 + a un +1 un = f n, = 1,..., N 1, u n 0 = µ n 0, u n N = µ n l, n = 0,..., M, (1.41) u 0 = u 0 (x), = 0,..., N, a + = an + a n, a = an a n, (1.4) которая аппроксимирует с порядком O(+) начально-краевую задачу u t + a(x, t)u x = f(x, t), 0 < x < l, 0 < t T, u(0, t) = µ 0 (t), u(l, t) = µ l (t), 0 t T, u(x, 0) = u 0 (x), 0 x l. (1.43) 14

16 С помощью принципа максимума можно доказать (см. задачу 1.10), что для устойчивости противопоточной схемы (1.41) с переменным коэффициентом a(x, t) достаточно выполнения условия max a(x, t) 1. (1.44) x,t 1.5. Схема Лакса. Далее для простоты изложения будем рассматривать начально-краевую задачу (1.7) с однородным уравнением переноса u t + au x = 0. (1.45) В схеме Лакса разностное уравнение, аппроксимирующее уравнение переноса (1.45), записывается так 0, 5 (u n +1 +) un 1 + a un +1 un 1 = 0, = 1,..., N 1. (1.46) Для локальной погрешности аппроксимации имеем выражение ψ n, = u tt u xx +..., поэтому при = O() схема Лакса не будет аппроксимировать уравнение переноса, а при законе предельного перехода r = a = const (1.47) будет аппроксимировать с порядком O(+). Таким образом, аппроксимация имеет место лишь при определенной связи между шагами и, т. е. схема Лакса принадлежит к классу условно аппроксимирующих схем. Для множителя перехода получаем формулу λ(φ) = cos φ ir sin φ. Следовательно, при законе предельного перехода (1.47) необходимое условие устойчивости схемы Лакса заключается в выполнении неравенства r 1, т. е. a 1. (1.48) 15

17 1.6. Схема Лакса Вендроффа. Разностные уравнения этой схемы выглядят так u +1/ 0, 5 (u n +1 +) un + a un +1 un = 0, / u n + a u +1/ (1.49) u 1/ = 0. Схема Лакса Вендроффа относится к семейству двухшаговых схем. В этой схеме сначала в полуцелых узлах x +1/ = x +/ по схеме Лакса вычисляются вспомогательные величины u +1/, относящиеся к моменту времени t n + /. Затем, на втором шаге, вычисляются значения искомой сеточной функции на (n + 1)-м слое по времени. Для исследования аппроксимации и устойчивости двухшаговых схем предварительно производится исключение из схемы вспомогательных величин u. В результате исключения получим одношаговую схему Лакса Вендроффа u n + a un +1 un 1 = a un +1 un + un 1, (1.50) которая, как легко проверить, аппроксимирует уравнение переноса (1.45) со вторым порядком по и. Для множителя перехода имеем такое выражение λ = 1 ir sin φ r sin φ. Поэтому необходимое условие устойчивости λ 1 будет равносильно выполнению неравенства (1 r sin φ) + r sin φ 1, или 1 4r sin φ + 4r4 sin 4 φ + 4r sin φ (1 sin φ) 1. Последнее неравенство равносильно условию r 1. Таким образом, необходимое условие устойчивости схемы Лакса Вендроффа совпадает с необходимым условием (1.48) устойчивости схемы Лакса Диссипация и дисперсия. Наряду с уравнением переноса u t + au x = 0, a = const (1.51) 16

18 рассмотрим еще два уравнения u t + au x = µu xx, µ = const > 0, (1.5) u t + au x + νu xxx = 0, ν = const. (1.53) Пусть начальная функция в задаче Коши для этих уравнений представляется в виде ряда Фурье u(x, 0) = m b m e imx. (1.54) Будем искать решение каждого из этих уравнений методом разделения переменных u(x, t) = b m λ t e imx = b m u m (x, t), (1.55) m m где u m (x, t) гармоника с волновым числом m u m (x, t) = λ t e imx, (1.56) λ подлежит определению. Действительная и мнимая часть гармоники представляют собой m-волны, длина l которых связана с волновым числом формулой l = π m. (1.57) Так как уравнения (1.51) (1.53) линейны, то поведение каждой из гармоник можно рассматривать независимо. Подставляя гармонику с волновым числом m в уравнение переноса (1.51), получаем или ln(λ) + aim = 0 λ = e aim. Следовательно, если гармоника (1.56) является решением уравнения переноса, то она имеет вид Обозначая ξ = x at, получаем u m (x, t) = e im(x at). (1.58) u m (x, t) = e imξ = u m (ξ, 0). (1.59) 17

19 Таким образом, в любой момент времени t > 0 гармоника u m получается сдвигом начальной гармоники на величину at. Следовательно, уравнение переноса описывает движение m-волн, которые независимо от их длины распространяются с постоянной скоростью v m = a без искажения своей формы. Легко проверить, что гармоника (1.56) будет решением второго уравнения (1.5), если ln(λ) + aim = µm или λ = e aim e µm, т. е. гармоника в этом случае имеет вид u m (x, t) = e µmt e im(x at). Следовательно, для всех гармоник происходит затухание амплитуды волн (диссипация волн). Поскольку m = π/l, то короткие волны затухают быстрее длинных. Скорость v m распространения волн не зависит от длины волн и равна по-прежнему a. За диссипацию волн отвечает член µu xx со второй производной от решения. Наконец, подстановка гармоники в уравнение (1.53) дает ln(λ) + aim + ν(im) 3 = 0 или откуда получаем, что λ = e im(a νm), u m (x, t) = e im(x (a νm)t). Таким образом, третье уравнение описывает движение волны без изменения ее амплитуды (без диссипации). Но скорость ее распространения зависит от длины волны v m = a νm. (1.60) Из этой формулы видно, что волны разной длины распространяются с различными скоростями (волны диспергируют). Более значительным изменениям подвергается скорость распространения коротковолновых возмущений (большие m). За дисперсию волн отвечает член νu xxx с третьей производной от решения. 18

20 Рассмотрев поведение отдельных гармоник, мы теперь сможем предсказать качественное поведение решения (1.55) задачи Коши для этих уравнений. Пусть, например, начальная функция u(x, 0) имеет вид ступеньки { 1, x 0, u(x, 0) = (1.61) 0, x > 0 и a > 0. Разложение такой функции в ряд Фурье (1.54) будет содержать весь набор гармоник. Решение задачи Коши для уравнения переноса (1.51) представляется в таком виде u(x, t) = m b m e im(x at) = m b m e imξ = u(ξ, 0), (1.6) т. е. решением задачи будет движущийся со скоростью a начальный профиль. Решение u(x, t) = m b m e µmt e im(x at) = m b m e µmt e imξ (1.63) задачи Коши для уравнения (1.5) с диссипативным членом, в котором короткие волны сильно затухают, будет иметь вид размазанной ступеньки. Наконец, решение u(x, t) = m b m e im(x (a νm)t) (1.64) задачи Коши для уравнения (1.53), в котором волны разной длины движутся с разными скоростями, имеет немонотонный, осциллирующий характер. Согласно формуле (1.60) при ν > 0 волны малой длины будут иметь скорость м еньшую, чем волны большой длины, а при ν < 0 наоборот. Поэтому осцилляции будут отставать от основного решения (описываемого первыми гармониками) при ν > 0 и, соответственно, перемещаться впереди при ν < Дифференциальное приближение разностной схемы. Вернемся к численному решению задачи Коши для уравнения переноса (1.51). В качестве начального профиля возьмем ступеньку { 1, x x0, u(x, 0) = (1.65) 0, x > x 0 19

21 и проведем расчет по явной противопоточной схеме u n + a un un 1 = 0, a = const > 0. (1.66) В результате получим решение в виде размазанной ступеньки (рис. 3), т. е. решение будет качественно таким же как и решение уравнения (1.5) с диссипативным членом. В чем дело? Ведь мы хотели решить уравнение переноса, в котором диссипативного члена нет. Дело в том, что мы искали численно решение не уравнения переноса, а решение разностной схемы. Таким образом, свойства решений аппроксимируемого дифференциального уравнения и аппроксимирующего разностного уравнения могут не совпадать. Как же в таком случае предсказать свойства решения разностного уравнения? y x 30 Рис. 3. Графики точного решения (штриховые линии) и численного решения (сплошные линии), полученного с помощью противопоточной схемы (1.66) в моменты времени t = 1 (1); t = 8 (); t = 15 (3). a = 1; x 0 = 10; a/ = 0, 5 Это можно сделать с помощью метода дифференциального приближения, с которым мы сейчас кратко познакомимся. Суть этого метода заключается в замене исходного разностного уравнения специальным дифференциальным уравнением, которое обладает всеми свойствами исследуемого разностного уравнения. Поэтому вместо исследования разностного уравнения исследуют это дифференциальное уравнение, что во многих случаях сделать гораздо проще. Получение дифференциального уравнения, соответствующего разностному уравнению, начинается с записи этого разностного уравнения в виде так называемой теоретической разностной схемы, в которой разностные операторы действуют в том же функциональном пространстве, что и аппроксимируемые ими дифференциальные операторы. Например, разностное уравнение (1.66) записывается в виде следующей теоретической разностной 0

22 схемы u(x, t +) u(x, t) u(x, t) u(x, t) + a = 0. (1.67) Решением такой схемы является функция u(x, t) непрерывных аргументов x и t в то время, как решением уравнения (1.66) является сеточная функция u, определенная только в узлах сетки. Пусть достаточно гладкая функция u(x, t) является решением теоретической разностной схемы (1.67). Подставим ее в эту схему и выразим u(x, t +) и u(x, t) через значения функции и ее производных в точке (x, t) по формуле Тейлора. В результате получим дифференциальное уравнение, эквивалентное разностной схеме (1.67) u t + au x + u tt + 6 u ttt a u xx + a 6 u xxx +... = 0. (1.68) Определение. Дифференциальное уравнение бесконечного порядка (1.68), полученное после разложения по формуле Тейлора решения u(x, t) теоретической разностной схемы (1.67), называется дифференциальным представлением разностной схемы (1.66). Некоторые свойства разностной схемы можно исследовать уже с помощью этого дифференциального представления, но для наших целей будет удобнее использовать другую форму дифференциального представления, получающуюся в результате исключения из (1.68) всех производных по времени кроме той, которая входит в аппроксимируемое уравнение (1.51), т. е. кроме u t. Покажем, например, как исключить производные по времени в членах порядка и. Для этого перепишем уравнение (1.68) с учетом слагаемых до порядка O() и O() u t + au x + u tt + 6 u ttt a u xx + a 6 u xxx = O() (1.69) и найдем с помощью полученного уравнения производную u t: u t = au x u tt 6 u ttt + a u xx a 6 u xxx + O() (1.70) Эту производную подставим в слагаемые уравнения (1.69), содержащие производные (u t) t и (u t) tt. Учитывая порядок малости коэффициентов при второй и третьей производных по времени, получаем, что в (u t) t 1

23 достаточно подставить производную (1.70), вычисленную с точностью O(+): u t = au x u tt + a u xx + O(+), (1.71) а в (u t) tt с точностью O(+): u t = au x + O(+). (1.7) В результате такой подстановки уравнение (1.69) примет следующий вид: u t + au x + (au x u tt + a) u xx + t 6 (au x) tt = = a u xx a 6 u xxx + O(), или u t + au x a u tx 4 u ttt + a 4 u txx a 6 u ttx = = a u xx a 6 u xxx + O(). (1.73) Выполнив подстановки в уравнение (1.69), далее аналогичные действия предпринимаем с уравнением (1.73). Теперь надо подставить производную u t, определенную из уравнения (1.73), в четыре слагаемых этого же уравнения: u t + au x a (au x + a u tx + a u xx) x 4 (au x) tt + + a 4 (au x) xx a 6 (au x) tx = a u xx a 6 u xxx + O(). После приведения подобных получим уравнение u t + au x a 1 u txx + a 4 u ttx = = a (a) (1 r) u xx + a u xxx + O(), 6 (1.74) в котором, в отличие от (1.69), нет вторых производных по времени. Оставшиеся в (1.74) смешанные производные u txx и u ttx вычислим на основе равенства (1.7): u txx = au xxx + O(+), u ttx = a u xxx + O(+). (1.75)

24 Следовательно, дифференциальное представление (1.74) примет вид u t + au x = a (1 r)u xx a 6 (r 3r + 1)u xxx + O(). (1.76) Таким образом, мы избавились от производных по времени при степенях и. Но производные по t пока остались при более старших степенях в правой части O(). Если продолжить описанную процедуру дальше, то в представлении (1.68) можно убрать производные по времени до сколь угодно высокого порядка. В результате получим дифференциальное представление схемы в виде или u t + au x = a (1 r)u xx + a 6 (1 r)(r 1)u xxx +... (1.77) u t + au x = k= c k k u x k. (1.78) Определение. Уравнение бесконечного порядка (1.78) называется П-формой дифференциального представления разностной схемы. Пусть разностная схема имеет порядки аппроксимации γ 1 и γ по и соответственно. Определение. Дифференциальное уравнение, полученное из П-формы дифференциального представления отбрасыванием членов порядка O(γ1+1, γ+1) и более высокого, называется первым дифференциальным приближением (п. д. п.) разностной схемы. Для противопоточной схемы (1.66) п. д. п. является дифференциальным уравнением второго порядка u t + au x = µu xx, µ = a (1 r), (1.79) которое, как видим, совпадает с уравнением (1.5) с диссипативным членом. Таким образом, при r 1 наша схема неявно вводит в аппроксимируемое уравнение переноса вязкость (диссипацию), которую называют аппроксимационной или схемной вязкостью. Наличие аппроксимационной вязкости и приводит к размазыванию начальной ступеньки. Определение. Свойство разностной схемы, обусловленное наличием в ее п. д. п. производных четного порядка, называется численной диссипацией. 3

25 П-форма дифференциального представления разностной схемы Лакса Вендроффа имеет вид u t + au x = a 6 (1 r)u xxx a3 8 r(1 r)u xxxx +..., а п. д. п. u t + au x + νu xxx = 0, ν = a 6 (1 r) (1.80) совпадает с уравнением (1.53) с дисперсионным членом. Следовательно, при r 1 схема Лакса Вендроффа неявно вводит в аппроксимируемое уравнение переноса дисперсию, поэтому решение разностной схемы может осциллировать (рис. 4). y Рис. 4. Графики точного решения (штриховые линии) и численного решения (сплошные линии), полученного с помощью схемы Лакса Вендроффа в моменты времени t = 1 (1); t = 8 (); t = 15 (3). a = 1; x 0 = 10; a/ = 0, 5 Определение. Свойство разностной схемы, обусловленное наличием в ее п. д. п. производных нечетного порядка, называется численной дисперсией. Подведем итог наших рассуждений. Для задач с плавно меняющимся решением, вклад в которое высокочастотных гармоник невелик, точность схемы Лакса Вендроффа выше точности противопоточной схемы. Если мы решаем численно задачу, в которой решение имеет резко меняющийся монотонный профиль, то применение противопоточной схемы первого порядка даст монотонный неосциллирующий профиль, но сильно сглаженный. Это результат действия численной диссипации. Схема Лакса Вендроффа, обладающая численной дисперсией, может дать немонотонные профили численного решения в окрестности разрыва или резкого изменения решения, искаженные нефизичными осцилляциями. x 4

26 З А Д А Ч И 1.1. Показать, что при a < 0 схема (1.8) абсолютно неустойчива. 1.. С помощью спектрального метода Неймана показать, что явная схема для уравнения (1.1) u n + a un +1 un = 0, n = 0,..., M 1, = 0, ±1, ±,... (1.81) при a > 0 абсолютно неустойчива С помощью спектрального метода Неймана вывести необходимое условие устойчивости трехслойной схемы «leap-frog» (схема с перешагиванием, схема «чехарда») для уравнения (1.1) u n 1 + a un +1 un 1 = 0, n = 1,..., M 1, = 0, ±1, ±,..., (1.8) если закон предельного перехода задан в виде (1.33) Определить порядок аппроксимации явной схемы с центральной разностью u n + a un +1 un 1 = 0, (1.83) построенной для уравнения переноса (1.1). С помощью спектрального метода Неймана исследовать устойчивость этой схемы, если закон предельного перехода задан в виде a = const. (1.84) 1.5. Определить порядок аппроксимации мажорантной схемы u n + a un +1 un 1 = a un +1 un + un 1, (1.85) построенной для уравнения переноса (1.1). С помощью спектрального метода Неймана исследовать устойчивость этой схемы, если закон предельного перехода задан в виде (1.84). 5

27 1.6. Определить порядок аппроксимации схемы Мак-Кормака u un + a un +1 un = 0, 0, 5 (u +) un / + a u u 1 = 0, (1.86) построенной для уравнения переноса (1.1). С помощью спектрального метода Неймана исследовать устойчивость этой схемы, если закон предельного перехода задан в виде (1.84) Определить порядок аппроксимации противопоточной схемы с весами u n + σa un (1 σ)a un un 1 = 0, (1.87) построенной для уравнения переноса (1.1) с коэффициентом a > 0. С помощью спектрального метода Неймана вывести необходимое условие устойчивости схемы (1.87), если закон предельного перехода задан в виде (1.84) Используя принцип максимума, исследовать устойчивость в равномерной норме неявной противопоточной схемы u n + a un+1 1 = f n+1, = 1,..., N, u n 0 = µ n 0, n = 0,..., M, u 0 = u 0(x), = 0,..., N, (1.88) построенной для задачи (1.7) Используя принцип максимума, найти достаточное условие устойчивости в равномерной норме противопоточной схемы с весами u n + σa un (1 σ)a un un 1 = f n+1/, u n 0 = µ n 0, n = 0,..., M, u 0 = u 0(x), = 0,..., N, (1.89) построенной для задачи (1.7). Здесь 0 σ 1. 6

28 1.10. Используя принцип максимума, доказать, что выполнение условия (1.44) достаточно для устойчивости противопоточной схемы (1.41) с переменным коэффициентом a(x, t) Получить п. д. п. (1.80) схемы Лакса Вендроффа Найти п. д. п. неявной схемы u n + a un+1 1 = 0, (1.90) построенной для уравнения переноса (1.1) с коэффициентом a > 0. Дать качественное объяснение поведения решения разностной схемы при t > 0, если в начальный момент времени t = 0 задана ступенька (1.61).. Свойство монотонности разностных схем.1. Одно из основных требований, предъявляемых к разностным схемам, состоит в том, что решение разностного уравнения должно передавать особенности поведения решения аппроксимируемого дифференциального уравнения. Рассмотрим, например, задачу Коши для линейного уравнения переноса u t + au x = 0, a = const > 0, < x <, t > 0, (.1) u(x, 0) = u 0 (x). (.) Если u 0 (x) неубывающая (невозрастающая) функция переменной x, то при любом фиксированном t > 0 решение u(x, t) задачи (.1), (.) также будет неубывающей (невозрастающей) функцией переменной x. Это следует из того, что в любой момент времени решение задается формулой u(x, t) = u 0 (x at). (.3) Естественно потребовать, чтобы и решение разностной схемы, аппроксимирующей задачу (.1), (.), тоже обладало аналогичным свойством. Но оказывается, что многие разностные схемы нарушают монотонность численного решения: вместо ожидаемых монотонных профилей получаются решения, содержащие нефизичные осцилляции (рис. 4). Причиной их возникновения является численная дисперсия разностных 7

29 схем, рассмотренная в предыдущем параграфе. В настоящем параграфе мы приведем условия, при выполнении которых разностная схема будет сохранять монотонность численного решения. Рассмотрим произвольную явную разностную схему = α b α u n +α, (.4) где α целое число, α = α 1, α 1 + 1,..., α, узлы x +α определяют шаблон схемы. Определение. Разностная схема (.4) называется схемой, сохраняющей монотонность численного решения (монотонной схемой), если она любую монотонную функцию u n переводит в монотонную на (n + 1)-м временн ом слое функцию, причем с тем же направлением роста. Пример.1. Аппроксимируем уравнение (.1) на равномерной сетке противопоточной схемой u n + a un un 1 = 0. (.5) Эта схема имеет первый порядок аппроксимации по и. Пусть сеточная функция u n на n-ом временн ом слое является монотонной, например, монотонно возрастающей функцией, т. е. u n un 1 для произвольного. В этом случае при выполнении условия устойчивости схемы, имеющего вид aæ 1, где æ = /, получим 1 = (u n aæ(u n u n 1)) (u n 1 aæ(u n 1 u n)) = = (1 aæ) (u n u n 1) + aæ(u n 1 u n) 0. Итак, решение монотонно возрастает и на (n + 1)-ом слое. Таким образом, противопоточная схема (обладающая диссипацией при aæ < 1) является схемой, сохраняющей монотонность. Пример.. Покажем, что схема Лакса Вендроффа (1.49) (не обладающая диссипацией при aæ < 1) не сохраняет монотонность численного решения. Пусть начальная функция для уравнения (.1) имеет вид (1.61) { 1, при x 0, u 0 (x) = 0, при x > 0. 8

30 Следовательно, начальная сеточная функция { u 0 1, при 0, = u 0 (x) = 0, при > 0 является монотонно убывающей. Перепишем рассматриваемую схему в виде одношаговой схемы (1.50), а затем в виде схемы (.4) с коэффициентами b 1 = a æ + aæ, b 0 = 1 a æ, b 1 = a æ aæ. (.6) Тогда нетрудно убедиться, что на первом слое по времени имеет место равенство 1, при 1, u 1 b = 1 + b 0, при = 0, b 1, при = 1, 0, при. При aæ < 1 схема устойчива, но b 1 + b 0 > 1, т. е. сеточная функция u 1 не является монотонно убывающей. Монотонность схем для уравнений с постоянными коэффициентами можно исследовать, пользуясь следующей теоремой . Теорема.1. Для того чтобы разностная схема (.4) с постоянными коэффициентами b α сохраняла монотонность, необходимо и достаточно выполнение при всех α условий b α 0. (.7) Д о к а з а т е л ь с т в о. Необходимость. Предположим, что схема (.4) сохраняет монотонность, однако существует отрицательный коэффициент b α0 < 0. Возьмем монотонно возрастающую функцию u n = { 0, < α0, 1, α 0. (.8) Тогда u0 n+1 1 = b α u n α b α u n 1+α = α α = b α b α = b α0 < 0, α α 0 α α

31 т. е. функция не является монотонно возрастающей, и, следовательно, схема (.4) не сохраняет монотонность, что противоречит исходному предположению. Полученное противоречие доказывает, что все коэффициенты b α неотрицательны. Достаточность. Пусть b α 0 и u n монотонная функция, например, монотонно возрастающая функция. Тогда 1 = α b α u n +α α b α u n 1+α = = α b α (u n +α u n 1+α) 0, т. е. также монотонно возрастающая функция. Таким образом, схема (.4) сохраняет монотонность. Возвратимся вновь к примерам.1 и., причем теперь не будем предполагать, что a > 0. Противопоточная схема для уравнения (.1) при произвольном знаке коэффициента a выглядит так: где u n Перепишем ее в виде (.4) + a + un un 1 a + = a + a + a un +1 un, a = a a. = 0, (.9) где = b 1 u n 1 + b 0 u n + b 1 u n +1, (.10) b 1 = æa +, b 0 = 1 æ a, b 1 = æa. При выполнении условия устойчивости a æ 1 (.11) все эти коэффициенты неотрицательны. Кроме того, они являются постоянными, поэтому, согласно теореме.1, противопоточная схема (.9) сохраняет монотонность решения при условии (.11). Схема Лакса Вендроффа устойчива при том же условии (.11), что и противопоточная схема, и ее можно записать в виде (.10) с коэффициентами (.6), откуда видно, что при условии a æ < 1 один из 30

32 коэффициентов b 1 или b 1 отрицателен. Согласно теореме.1, отсюда следует, что схема Лакса Вендроффа, имеющая второй порядок аппроксимации по и, не сохраняет монотонность численного решения. Но, возможно, существуют другие схемы второго порядка аппроксимации, которые обладают свойством монотонности. Оказывается, что таких схем нет. В работе показано, что для линейного уравнения переноса (.1) невозможно построить монотонную схему с постоянными коэффициентами второго порядка аппроксимации... Рассмотрим теперь схему (.4) с переменными коэффициентами b α. Будет ли для таких схем условие (.7) неотрицательности коэффициентов достаточным для сохранения монотонности численного решения? Оказывается, нет. Приведем соответствующий пример. Пример.3. Пусть решается задача Коши для уравнения u t + a(x)u x = 0, (.1) где a(x) строго возрастающая положительная ограниченная функция: 0 < a(x) < 1 и a > 0. Возьмем для решения этой задачи схему с переменными коэффициентами 0, 5 (u n +1 +) un 1 + a u n +1 un 1 = 0, (.13) где a = a(x), x узел равномерной сетки. Выписанная схема является аналогом схемы Лакса (1.46), которая сохраняет монотонность численного решения (см. задачу.1). Будем считать, что для любого выполнено условие æa < 1, (.14) гарантирующее устойчивость схемы (.13) в равномерной норме по начальным данным: C u 0 C. (.15) Запишем схему (.13) в виде (.4): = b 1, u n 1 + b 1, u n +1, (.16) где b 1, = 1 + æa, b 1, = 1 æa, 31

33 при этом коэффициенты b α снабжены дополнительным индексом, поскольку они являются переменными коэффициентами и изменяются при переходе от одного узла к другому. В силу условия (.14) оба коэффициента положительны, однако схема (.13) не сохраняет монотонность численного решения. В самом деле, взяв монотонно возрастающую функцию { u n 0, < 0, = 1, 0, убеждаемся, что на (n + 1)-м слое по времени имеет место равенство = 0, при < 1, b 1, 1, при = 1, b 1,0, при = 0, 1, при 1. Но b 1, 1 > b 1,0, поэтому сеточная функция возрастающей. не является монотонно Приведенный пример показывает, что для схем с переменными коэффициентами должны использоваться другие признаки монотонности, нежели признак (.7), указанный в теореме.1. Теорема.. Пусть коэффициенты разностной схемы = b 1, u n 1 + b 0, u n + b 1, u n +1 (.17) удовлетворяют в каждом узле x условию Тогда выполнение при всех условий b 1, + b 0, + b 1, = 1. (.18) b ±1, 0, b 1, + b 1, 1 1 (.19) необходимо и достаточно для того, чтобы схема (.17) с переменными коэффициентами сохраняла монотонность численного решения. Д о к а з а т е л ь с т в о. Запишем схему (.17) с переменными коэффициентами, удовлетворяющими условию (.18), в следующем виде: = u n b 1, (u n u n 1) + b1, (u n +1 u n). (.0) 3

34 Тогда +1 = un +1 b 1,+1 (u n +1 u n) + b1,+1 (u n + u n +1). Следовательно, +1 un+1 = (u n +1 u n) (1 b 1,+1 b 1,) + (+ b 1,+1 u n + u n (+1) + b 1, u n u n) (.1) 1. Необходимость. Пусть схема (.17) монотонна. Докажем, что ее коэффициенты удовлетворяют неравенствам (.19). Предположим, что это не так и какие-то из условий (.19) не выполняются в некотором узле x 0, например b 1,0 < 0. Положим Из (.1) тогда следует u n = { 0, если < 0, 1, если un+1 0 = b 1,0 < 0, т. е. функция не является монотонно возрастающей, что противоречит исходному предположению о монотонности схемы (.17). Аналогично проверяются и остальные неравенства в (.19). Достаточность. Пусть в каждом узле x коэффициенты схемы (.17) удовлетворяют неравенствам (.19) и функция u n является монотонной, например монотонно возрастающей. Тогда из равенства (.1) следует, что функция также будет монотонно возрастающей функцией. Теорема доказана. Нетрудно проверить, что коэффициенты схемы (.16) не удовлетворяют второму из условий (.19) теоремы.3, поэтому эта схема не является схемой, сохраняющей монотонность численного решения. Дадим другую формулировку теоремы.. Теорема.3. Для того чтобы конечно-разностная схема u n + C 1/ un x, 1/ C+ +1/ un x,+1/ = 0, (.) сохраняла монотонность численного решения, необходимо и достаточно выполнение при всех условий где æ = /, u n x,+1/ = un +1 un C ± +1/ 0, C +1/ + C+ +1/ 1 æ, (.3). 33

35 Д о к а з а т е л ь с т в о. Схему (.) можно переписать в виде (.17), при этом b 1, = æc 1/, b 1, = æc + +1/, b 0, = 1 æc 1/ æc+ +1/. Тогда для коэффициентов b α выполняется равенство (.18), а условия (.19) эквивалентны условиям (.3). Замечание. В работе доказано, что выполнение неравенств (.3) достаточно для того, чтобы схема (.) была TVD-схемой (Total Variation Diminising Sceme), т. е. схемой, решение u n которой при любом n 0 удовлетворяет условию невозрастания полной вариации TV () TV (u n), (.4) где под полной вариацией сеточной функции u n понимается величина TV (u n) = = u n +1 u n. (.5) В настоящее время TVD-схемы и их разнообразные модификации применяются при решении многих задач с разрывными решениями. Причина столь большой популярности этих методов заключается в том, что они дают неосциллирующие профили решения, высокую разрешимость в области разрывов и сохраняют высокую точность в областях гладкости решения. Современные TVD-схемы высокого порядка аппроксимации основаны на тех или иных способах восстановления (реконструкции) значений функций на границах ячеек по их значениям в центрах соседних ячеек. При этом шаблон схемы является переменным и зависит от поведения численного решения. Алгоритмы реконструкции основываются на использовании специальных ограничителей потоков , которые строятся так, чтобы схема с ограничителями обладала TVDсвойством (.4)..3. Монотонизация схемы Лакса Вендроффа. Если начальная функция при t = 0 задана в виде ступеньки, то на следующих слоях по времени мы будем получать по схеме Лакса Вендроффа ступеньку, искаженную осцилляциями (см. рис. 4). Но оказывается, что схему Лакса Вендроффа можно модифицировать так, чтобы она стала обладать 34

36 TVD-свойством (.4), а значит, согласно теореме.3, стала бы схемой, сохраняющей монотонность численного решения. Однако коэффициенты модифицированной схемы уже не будут постоянными, они могут зависеть от решения на n-м слое, т. е. модифицированная схема будет нелинейной. Рассмотрим уравнение переноса (.1) в случае a = const > 0. Схема Лакса Вендроффа (1.50) может быть переписана так или u n +a un x,+1/ + un x, 1/ a () u n x,+1/ un x, 1/ = 0, (.6) или u n + au n x, 1/ + a (1 aæ) un x,+1/ un x, 1/ u n = 0, (.7) + au n x, = a (1 aæ) un xx,. (.8) П. д. п. (1.79) противопоточной схемы содержит в правой части диссипативный член 0, 5a (1 aæ) u xx, а в представлении (.8) такой же диссипативный член в разностной форме имеет противоположный знак. Таким образом, схема Лакса Вендроффа представляется в виде монотонной схемы с направленной против потока разностью, дополненной так называемым антидиффузионным членом, который устраняет диссипативный член в п. д. п. противопоточной схемы, превращая ее в схему Лакса Вендроффа. Уменьшая антидиффузионный член в местах возможного появления осцилляций, можно попытаться предотвратить их. Регулировать антидиффузионный член в схеме Лакса Вендроффа (.7) будем с помощью функции-ограничителя Φ(ξ) некоторого аргумента ξ: u n +au n x, 1/ + a (1 aæ) ((Φu n x) +1/ (Φu n x) 1/) = 0. (.9) Если Φ 0, то имеем монотонную противопоточную схему первого порядка аппроксимации. Если же Φ 1, то получаем схему Лакса Вендроффа второго порядка аппроксимации на гладких решениях, но осциллирующую на разрывных решениях. 35

37 В разностной схеме (.9) Φ +1/ = Φ(ξ +1/). В качестве дискретного аргумента ξ +1/ выберем величину u n x, 1/ ξ +1/ = u n при u n x,+1/ 0, x,+1/ (.30) 1 при u n x,+1/ = 0. На осциллирующем решении отношение u n x, 1/ /un x,+1/ становится отрицательным, поэтому при ξ +1/ < 0 полагаем, что Φ +1/ = 0. Далее будем считать, что функция Φ = Φ(ξ) непрерывного аргумента ξ также принимает нулевые значения при ξ < 0. Более того, предполагая, что функция-ограничитель является непрерывной, полагаем, что Φ(ξ) 0 при всех ξ 0. Далее рассмотрим случай, когда ξ +1/ > 0. Будем подбирать функцию-ограничитель таким образом, чтобы схема удовлетворяла TVDусловию (.3) и сохраняла второй порядок аппроксимации на гладких решениях. Для этого преобразуем модифицированную схему Лакса Вендроффа (.9) к виду (.): или u n + au n x, 1/ + a (1 aæ) ((Φ ξ u n [ + a aæ ((Φ) ξ) +1/ +1/ Φ 1/) u n x, 1/ = 0, Φ 1/)] u n x, 1/ = 0. Таким образом, коэффициенты схемы (.9), записанной в виде (.), определяются по формулам [ C + +1/ = 0, C 1/ = a aæ ((Φ))] ξ Φ 1/. +1/ Согласно теореме.3, условие 0 C 1/ 1 æ (.31) будет гарантировать, что схема Лакса Вендроффа с введенной в нее функцией-ограничителем будет сохранять монотонность численного решения. Далее мы предполагаем, что условие устойчивости схемы Лак- 36

38 са Вендроффа выполнено, т. е. aæ 1. Тогда для того чтобы неравенства (.31) были справедливыми для всех aæ 1, необходимо и достаточно выполнения неравенств (Φ) ξ +1/ Φ 1/, а для этого достаточно потребовать выполнения для всех следующих неравенств: (Φ) 0, 0 Φ +1/. ξ +1/ Область в плоскости переменных Φ и ξ, в которой выполняются эти неравенства, изображена на рис. 5, а. Если график функции Φ = Φ(ξ) лежит в этой области, то модифицированная схема (.9) будет сохранять монотонность численного решения. Φ Φ= Φ Φ = Φ = ξ Φ = ξ Φ=ξ 1 1 Φ = а ξ б ξ Рис. 5. а в заштрихованной области модифицированная схема Лакса Вендроффа (.9) является TVD-схемой; б в области с двойной штриховкой модифицированная схема Лакса Вендроффа является TVD-схемой второго порядка аппроксимации Итак, далее будем считать, что Φ(ξ) = 0 при ξ 0, 0 Φ(ξ) min(, ξ) при ξ > 0. (.3) Исследуем теперь порядок аппроксимации модифицированной схемы, предполагая, что непрерывная функция Φ = Φ(ξ) удовлетворяет 37

39 следующим дополнительным ограничениям: Φ(ξ 1) Φ(ξ) L ξ 1 ξ, ξ 1, ξ, (.33) Φ(1) = 1, (.34) т. е. потребуем, чтобы функция Φ = Φ(ξ) удовлетворяла условию Липшица с некоторой постоянной L > 0 и график этой функции проходил через точку (1, 1). Перепишем модифицированную схему Лакса Вендроффа (.9) в виде исходной схемы Лакса Вендроффа (.7) с добавочным членом где u n + au n x, 1/ + a (1 aæ) (u n x,+1/ un x, 1/) + + a (1 aæ) Rn = 0, (.35) R n = (Φ +1/ 1) u n x,+1/ (Φ 1/ 1) u n x, 1/. (.36) Пусть u = u(x, t) достаточно гладкое решение задачи Коши (.1), (.). Подставим это решение в выражение (.36), сохранив в нем все прежние обозначения, но учитывая, что теперь u n x,+1/ = u(x +1, t n) u(x, t n). (.37) Очевидно, что если на n-м слое по времени функция u(x, t n) является линейной, u(x, t n) = Bx + C, то R n 0. Используя условия (.33), (.34), нетрудно проверить, что для квадратичной функции u(x, t n) = Ax + Bx + C (A 0) равенство R n = O() имеет место для всех узлов произвольного числового промежутка (α, β), не содержащего точку экстремума x = B/A. В общем случае справедливо следующее утверждение. Лемма.1. Пусть выполнены условия (.33), (.34) и достаточно гладкое решение задачи Коши (.1), (.) удовлетворяет на некотором числовом отрезке [α, β] условию u x (x, t n) 0 x [α, β]. (.38) Тогда R n = O() x (α, β). (.39) 38


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

ФЕДЕРАЛЬНОЕ АГЕНТСТВО ПО ОБРАЗОВАНИЮ НОВОСИБИРСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ Механико-математический факультет Г. С. Хакимзянов, С. Г. Черный МЕТОДЫ ВЫЧИСЛЕНИЙ Часть 3. Численные методы решения задач

Теория устойчивости разностных схем 1 Устойчивость решения задачи Коши по начальным данным и правой части Пусть B банахово (то есть полное нормированное) пространство функций, заданных в некоторой области

Основные понятия теории разностных схем. Примеры построения разностных схем для начально-краевых задач. Большое количество задач физики и техники приводит к краевым либо начальнокраевым задачам для линейных

Дифференциальные уравнения. 1999. Т.35. 6. С.784-792. УДК 517.957 ОДНОЗНАЧНАЯ РАЗРЕШИМОСТЬ КРАЕВЫХ ЗАДАЧ ДЛЯ ЭЛЛИПТИЧЕСКИХ УРАВНЕНИЙ С НЕЛИНЕЙНОСТЯМИ Ю. В. Жерновый 1. Введение. Постановка задачи. Наиболее

Разностная аппроксимация начально-краевой задачи для уравнения колебаний. Явная (схема «крест») и неявная разностные схемы. Рассмотрим несколько вариантов разностной аппроксимации линейного уравнения колебаний:

Глава IV. Первые интегралы систем ОДУ 1. Первые интегралы автономных систем обыкновенных дифференциальных уравнений В этом параграфе будем рассматривать автономные системы вида f x = f 1 x, f n x C 1

Разностная аппроксимация начально-краевой задачи для уравнения теплопроводности. Понятие явной и неявной схемы. 1 Разностная аппроксимация уравнения теплопроводности Рассмотрим различные варианты разностной

Теория устойчивости разностных схем 1 Операторно-разностные схемы 1.1 Введение Пусть B банахово (то есть полное нормированное) пространство функций, заданных в некоторой области G R m, и пусть u(t) абстрактная

Уравнения переноса. Схемы «бегущего» счета Рассмотрим ряд наиболее часто используемых разностных схем, аппроксимирующих начально-краевые задачи для линейного уравнения переноса: u t + c(x, t) u x = f(x,

Скалько Юрий Иванович Цыбулин Иван Шевченко Александр Волновое уравнение второго порядка Волновое уравнение в форме уравнения второго порядка записаывается как 2 u t 2 = c2 2 u x 2 + f Дополним уравнение

МЕТОДЫ ВЫЧИСЛЕНИЙ Лекторы: проф. Б. И. Квасов, проф. Г. С. Хакимзянов 5 6 семестры 1. Математические модели и вычислительный эксперимент. Классификация уравнений математической физики. Примеры корректных

Разностные схемы для уравнения колебаний в многомерном случае Для многомерных уравнений колебаний можно составить аналог схемы «крест» и неявной схемы. При этом явная схема «крест» так же, как и в одномерном

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

Уравнения В алгебре рассматривают два вида равенств тождества и уравнения Тождество это равенство которое выполняется при всех допустимых) значениях входящих в него букв Для тождества используют знаки

Простейшие способы исследования разностных схем на устойчивость Напомним, что разностная схема L h y h = ϕ h (x), x ω h, l h y h = χ h (x), x γ h, аппроксимирующая краевую или начально-краевую задачу Lu

ГЛАВА УСТОЙЧИВОСТЬ ЛИНЕЙНЫХ СИСТЕМ В этой главе исследуется устойчивость самого простого класса дифференциальных систем линейных систем В частности, устанавливается, что для линейных систем с постоянными

Сибирский математический журнал Январь февраль, 2001. Том 42, 1 УДК 517.929 МЕТОДЫ ИССЛЕДОВАНИЯ УСТОЙЧИВОСТИ СИСТЕМ С ЛИНЕЙНЫМ ЗАПАЗДЫВАНИЕМ Б. Г. Гребенщиков Аннотация: Изложены методы исследования асимптотической

Глава 1 Дифференциальные уравнения 1.1 Понятие о дифференциальном уравнении 1.1.1 Задачи, приводящие к дифференциальным уравнениям. В классической физике каждой физической величине ставится в соответствие

ЛЕКЦИИ 8 9 Теорема Хилле Иосиды S 3. Определение и элементарные свойства максимальных монотонных операторов Всюду на протяжении этих двух лекций символом H обозначено гильбертово пространство со скалярным

Федеральное агентство по образованию Федеральное государственное образовательное учреждение высшего профессионального образования ЮЖНЫЙ ФЕДЕРАЛЬНЫЙ УНИВЕРСИТЕТ Р. М. Гаврилова, Г. С. Костецкая Методические

ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ Московский государственный технический университет имени Н.Э. Баумана Факультет «Фундаментальные науки» Кафедра «Математическое моделирование» À.Í. Êàíàòíèêîâ, À.Ï. Êðèùåíêî

Модуль Тема Функциональные последовательности и ряды Свойства равномерной сходимости последовательностей и рядов Степенные ряды Лекция Определения функциональных последовательностей и рядов Равномерно

Дифференциальные уравнения первого порядка разрешенные относительно производной Теорема существования и единственности решения В общем случае дифференциальное уравнение первого порядка имеет вид F ()

ГЛАВА: Метод конечных разностей. Лекция 5: Устойчивость разностных схем (10 слайдов, 6 рисунков) Слайд 1: Классификация РС по типам устойчивости. По типам устойчивости выделяют следующие РС: абсолютно

МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ ГРАЖДАНСКОЙ АВИАЦИИ В.М. Любимов, Е.А. Жукова, В.А. Ухова, Ю.А. Шуринов М А Т Е М А Т И К А Р Я Д Ы ПОСОБИЕ по изучению дисциплины и контрольные задания

Лекция 9 Линеаризация диффе6ренциальных уравнений Линейные дифференциальные уравнения высших порядков Однородные уравнения свойства их решений Свойства решений неоднородных уравнений Определение 9 Линейным

ФЕДЕРАЛЬНОЕ АГЕНТСТВО ПО ОБРАЗОВАНИЮ МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ СТРОИТЕЛЬНЫЙ УНИВЕРСИТЕТ Институт Фундаментального Образования Факультет Общенаучных Кафедр - ФОК Колебания бесконечной струны. Формула Даламбера.

Лекция 3 Теорема существования и единственности решения скалярного уравнения Постановка задачи Основной результат Рассмотрим задачу Коши d f () d =, () = Функция f (,) задана в области G плоскости (,

Методы построения разностных схем Однородные схемы для уравнения второго порядка с переменными коэффициентами Под однородными разностными схемами понимаются разностные схемы, вид которых не зависит ни

ВАРИАЦИЯ И ЭКСТРЕМУМ ФУНКЦИОНАЛА А. Н. Мягкий Интегральные уравнения и вариационное исчисление Лекция Пусть задан функционал V = V , y(x) M E. Зафиксируем функцию y (x) M. Тогда любую другую функцию

Экономичные разностные схемы для многомерных задач математической физики. Схема переменных направлений для начально-краевой задачи для уравнения теплопроводности в прямоугольнике. Как уже было показано

ПОНЯТИЕ ПРОИЗВОДНОЙ ФУНКЦИИ Пусть имеем функцию определенную на множестве X и пусть точка X - внутренняя точка те точка для которой существует окрестность X Возьмем любую точку и обозначим через называется

Уравнения гиперболического типа. Колебания бесконечной и полубесконечной струны. Метод Фурье Метод Фурье Стоячие волны 4 Лекция 4.1 Уравнения гиперболического типа. Колебания бесконечной и полубесконечной

Предварительные сведения теории разностных схем 1 Формулы суммирования по частям и разностные формулы Грина для сеточных функций Получим ряд соотношений, которые в дальнейшем будем использовать при исследовании

Лекция 8 4 Задача Штурма-Лиувилля Рассмотрим начально-краевую задачу для дифференциального уравнения в частных производных второго порядка описывающего малые поперечные колебания струны Струна рассматривается

II ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ Дифференциальные уравнения первого порядка Определение Соотношения, в которых неизвестные переменные и их функции находятся под знаком производной или дифференциала, называются

Лекция 5 5 Теорема существования и единственности решения задачи Коши для нормальной системы ОДУ Постановка задачи Задача Коши для нормальной системы ОДУ x = f (, x), () состоит в отыскании решения x =

Составитель ВПБелкин 1 Лекция 1 Функция нескольких переменных 1 Основные понятия Зависимость = f (1, n) переменной от переменных 1, n называется функцией n аргументов 1, n В дальнейшем будем рассматривать

Глава 4. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ И ИХ СИСТЕМ В этой главе рассматриваются основные численные методы решения задачи Коши для обыкновенных дифференциальных уравнений

Методы решения сеточных уравнений 1 Прямые и итерационные методы В результате разностной аппроксимации краевых задач математической физики получаются СЛАУ, матрицы которых обладают следующими свойствами:

Описание вычислительных моделейequatio Capter Sectio.. Разностные схемы для уравнений параболического типа Рассмотрим вначале простейшее уравнение теплопроводности: u, t uxx cost (.) Рис.. Введем в области

МИНИСТЕРСТВО ОБРАЗОВАНИЯ И НАУКИ РОССИЙСКОЙ ФЕДЕРАЦИИ Государственное образовательное учреждение высшего профессионального образования «НАЦИОНАЛЬНЫЙ ИССЛЕДОВАТЕЛЬСКИЙ ТОМСКИЙ ПОЛИТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ»

Задача Коши для волнового уравнения. Формула Даламбера 37, 438, I, II, 385, 439, 445, 37, III, IV, 37, 446.. 37 Найти общее решение уравнения u tt a u xx..) Шаг. Находим замену переменных Способ через

МЕТОДИЧЕСКИЕ УКАЗАНИЯ К РАСЧЕТНЫМ ЗАДАНИЯМ ПО КУРСУ ВЫСШЕЙ МАТЕМАТИКИ «ОБЫКНОВЕННЫЕ ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ РЯДЫ ДВОЙНЫЕ ИНТЕГРАЛЫ» ЧАСТЬ Ш ТЕМА РЯДЫ Оглавление Ряды Числовые ряды Сходимость и расходимость

Кафедра математики и информатики Элементы высшей математики Учебно-методический комплекс для студентов СПО, обучающихся с применением дистанционных технологий Модуль Дифференциальное исчисление Составитель:

ГЛАВА. УСТОЙЧИВОСТЬ ЛИНЕЙНЫХ СИСТЕМ 8 степень со знаком +, из полученного следует, что () π возрастает от до π. Итак, слагаемые ϕ i() и k () +, т. е. вектор (i) ϕ монотонно ϕ монотонно возрастают при

ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ ÌÃÒÓ Московский государственный технический университет имени Н.Э. Баумана Факультет «Фундаментальные науки» Кафедра «Математическое моделирование» À.Í. Êàíàòíèêîâ,

Глава 6 Основы теории устойчивости Лекция Постановка задачи Основные понятия Ранее было показано, что решение задачи Коши для нормальной системы ОДУ = f, () непрерывно зависит от начальных условий при

Глава 9. Численные методы. Лекция 4. Разностный метод Эйлера решения задачи Коши для дифференциальных уравнений.. Дифференциальная и разностная задачи Эйлера. Определение. Дифференциальной задачей Эйлера

ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ Общие понятия Дифференциальные уравнения имеют многочисленные и самые разнообразные приложения в механике физике астрономии технике и в других разделах высшей математики (например

Линейные и нелинейные уравнения физики Уравнение Лапласа в полярной системе координат. Старший преподаватель кафедры ВММФ Левченко Евгений Анатольевич 518 Глава 5. Уравнения эллиптического типа 25.2. Разделение

Лекция 3 Устойчивость равновесия и движения системы При рассмотрении установившихся движений уравнения возмущенного движения запишем в виде d dt A Y где вектор-столбец квадратная матрица постоянных коэффициентов

Числовые ряды Числовая последовательность Опр Числовой последовательностью называют числовую ф-цию, определенную на множестве натуральных чисел х - общий член последовательности х =, х =, х =, х =,

5 Степенные ряды 5 Степенные ряды: определение, область сходимости Функциональный ряд вида (a + a) + a () + K + a () + K a) (, (5) где, a, a, K, a,k некоторые числа, называют степенным рядом Числа

Министерство образования Российской Федерации МАТИ - РОССИЙСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНОЛОГИЧЕСКИЙ УНИВЕРСИТЕТ им К Э ЦИОЛКОВСКОГО Кафедра Высшая математика В В Горбацевич К Ю Осипенко Уравнения с частными

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

Обобщение на случай уравнения Хопфа схемы П.Лакса имеет вид

Здесь, очевидно, используется дивергентный вид уравнения (3.6).

Упражнения . Рассмотрим схему Лакса - Вендроффа для уравнения Хопфа. Пусть начальные условия для задачи Коши поставлены следующим образом: u(x, 0) = ch - 2 (x) . Тогда уравнение Хопфа имеет первый интеграл : . Проверить, что приведенная выше схема является консервативной , т.е. в ней на сеточном уровне автоматически выполняется тот же закон сохранения.

Построить аналогичную схему с использованием характеристической формы записи уравнения Хопфа (3.9). Будет ли она консервативной?

Схема условно устойчива при выполнении условия Куранта (точнее, обобщения условия Куранта)

Здесь и ниже, как и ранее в (3.7), f = 0, 5u 2 . При этом предполагается, что течение достаточно гладкое, момент градиентной катастрофы еще не наступил, в решении нет ударных волн и других разрывов.

Схема Куранта - Изаксона - Риса . Обобщение схем КИР на квазилинейный случай (при использовании дивергентной формы записи уравнений) очевидно.

Схема устойчива при выполнении условия Куранта

Обобщение схемы Лакса - Вендроффа (схема предиктор - корректор). Для квазилинейных уравнений (а также линейных уравнений с переменными коэффициентами, неоднородных уравнений и т.п.) схема Лакса - Вендроффа становится более сложной. Для ее построения необходимо ввести так называемые полуцелые точки (точки с дробными индексами). На первом этапе (предиктора) значения в полуцелых точках вычисляются по приведенной выше схеме - обобщению на квазилинейный случай схемы Лакса:

на втором этапе (корректор) используется схема "чехарда" (трехслойная схема на крестообразном шаблоне, которая не входит в семейство (3.8)):

Схема Лакса - Вендроффа принадлежит к так называемым центральным схемам. Ее шаблон симметричен. На первом этапе рассчитываются значения сеточной функции в полуцелых точках шаблона на промежуточном слое (t m - 1/2 , x m - 1/2) , (t n + 1/2 , x m + 1/2) , на втором этапе вычисляется решение на верхнем слое в точке (t n + 1 , x m) . Схема устойчива при выполнении условия Куранта.

Аналогично строятся схемы Лакса - Вендроффа для линейных неоднородных уравнений.

Нецентральная схема Мак - Кормака (предиктор - корректор).

Как и приведенная выше схема Лакса - Вендроффа, схема МакКормака состоит из двух этапов. Рассмотрим построение схемы МакКормака для однородного уравнения (3.7). Первый этап (предиктор) имеет вид

т.е. используется схема "явный правый уголок". Второй этап - корректор:

Таким образом, расчет на первом этапе по схеме "правый уголок", на втором - "левый уголок".

Другая схема Мак - Кормака имеет вид

Такие разностные схемы называют нецентральными . К их преимуществам относят отсутствие полуцелых индексов, более простую постановку граничных условий. В линейном случае схемы Мак - Кормака совпадают со схемой Лакса - Вендроффа. Схемы имеют второй порядок аппроксимации по обеим переменным, схемы устойчивы при выполнении условия Куранта.

Схема Русанова (центральная схема третьего порядка точности).

Для построения схемы Русанова вводятся не только полуцелые точки, но и два слоя промежуточных точек с дробными индексами. Первый этап схемы Русанова (переход к слою 1/3) имеет вид

ее второй этап есть схема "чехарда"

а третий этап

На первом этапе производится расчет по схеме Лакса, на втором - по схеме "крест" ("чехарда"). Последнее слагаемое третьего этапа вводится для обеспечения устойчивости схемы (член, пропорциональный разностной аппроксимации 4 - й производной).

Схема является условно устойчивой при выполнении условия Куранта и условия .

Нецентральная схема Уорминга - Кутлера - Ломакса 3 - го порядка точности.

Первый этап:

Второй этап:

Третий этап:

Последний член добавляется для устойчивости схемы, которая является условно устойчивой при выполнении условий Куранта.