Дискуссионный математический форумМатематический форум
Математический форум Math Help Planet

Обсуждение и решение задач по математике, физике, химии, экономике

Теоретический раздел
Часовой пояс: UTC + 3 часа [ Летнее время ]
MathHelpPlanet.com RSS-лента Математического форума

Часовой пояс: UTC + 3 часа [ Летнее время ]


Разностные схемы решения уравнений в частных производных 1-го порядка

Разностные схемы решения уравнений
в частных производных 1-го порядка


Рассмотрим проблему конструирования разностных схем для решения уравнения (8.4):


A(x,t)\cdot u_{x}+ B(x,t)\cdot u_{t}+ C(x,t)\cdot u= G(x,t),

где u(x,t) — искомая функция; A(x,t),\,B(x,t),\, C(x,t) — коэффициенты, G(x,t) — правая часть, B(x,t)\ne0. Это уравнение относится к эволюционным. Основная идея состоит в том, что после замены дифференциального уравнения его конечно-разностной аппроксимацией получаются формулы, явно или неявно выражающие значения решения для одного момента времени через значения решения в предыдущий момент времени. Таким образом, если известно решение в начальный момент времени, можно шаг за шагом найти решение для всех последующих моментов. Чтобы проиллюстрировать общий подход, рассмотрим типичную задачу для уравнения переноса.


Пример 8.2. Построить разностные схемы для задачи u_t+u_x=0,\,-\infty<x<+\infty,~ 0<t<T, u(x,0)=\psi(x),\,-\infty<x<+\infty, где \psi(x) — заданная функция.


Решается задача Коши, в которой A(x,t)=1,~ B(x,t)=1,~ C(x,t)=0,~ G(x,t)=0. Для анализа различных разностных схем, используемых для нахождения приближенного решения, найдем сначала аналитическое решение поставленной задачи. Оно может быть найдено различными методами, в частности методом характеристик, который опирается на следующий факт: начальное условие в некоторой точке x при t=0 переносится в плоскости Oxt вдоль линии, называемой характеристикой. Рассмотрим один из вариантов этого метода.


1. Решаем систему двух обыкновенных дифференциальных уравнений (уравнений характеристик), где s,\,\mu — параметры, 0<s<\infty,\,-\infty<\mu<+\infty\colon


\left\{\begin{aligned}&\frac{dx}{ds}=A(x,t), & & x(0)=\mu,\\ &\frac{dt}{ds}= B(x,t), & & t(0)=0.\end{aligned}\right.
(8.33)

Поскольку A(x,t)=1,~ B(x,t)=1, то система имеет вид


\left\{\! \begin{aligned}&\frac{dx}{ds}=1, \quad x(0)=\mu;\\ &\frac{dt}{ds}=1,\quad t(0)=0.\end{aligned}\right.

Отсюда x(s)=s+C_1,~ t(s)=s+C_2 и x(0)=C_1=\mu,~ t(0)=C_2=0. Следовательно, характеристики описываются уравнениями x(s)=s+\mu,~ t=s.


Заметим, что в общем случае результатом п.1 будут формулы связи переменных x,\,t и s,\,\mu\colon x=x(s,\mu),~ t=t(s,\mu), где параметр s изменяется вдоль характеристик, а параметр \mu изменяется вдоль прямой линии t=0.


Исключая s, получаем x=t+\mu,\,-\infty<\mu<+\infty, т.е. каждому значению ц соответствует прямая линия (рис. 8.4).


График линейной функции в зависимости от параметра

2. Используя координаты s,\,\mu и то, что с учетом (8.33) справедливо \frac{du}{ds}= u_x\,\frac{dx}{ds}+ u_t\,\frac{dt}{ds}= A(x,t)u_x+ B(x,t)u_t, сведем исходное уравнение в частных производных к обыкновенному дифференциальному уравнению


\frac{du}{ds}= C \bigl(x(s,\mu),t(s,\mu)\bigr)u= G \bigl(x(s,\mu),t(s,\mu)\bigr),\quad u(0)= \psi(\mu).
(8.34)

Для рассматриваемого примера уравнение (8.34) имеет вид \frac{du}{ds}=0,~ u(0)= \psi(\mu). Отсюда u(s,\mu)=C, u(0,\mu)= C=\psi(\mu), следовательно, u(s,\mu)= \psi(\mu).


3. Найдем решение u(x,t) исходной задачи, выражая s и \mu через x и t (см. п.1). Так как s=t,~ \mu=x-t, то u(x,t)= \psi(x-t). На каждой характеристике значение решения постоянно, определяется значением параметра \mu и равно \psi(\mu). Например, значение решения в точке B то же, что и в точке A, и равно y(0) (рис. 8.4).


Заметим, что если функция \psi(\mu) имеет вид \psi(\mu)= \begin{cases} \psi_0(x)\ne0,& x<0,\\ 0,& x \geqslant 0,\end{cases} решение задачи в точках, расположенных на характеристике, проходящей через точки A и B, и правее этой характеристики, равно нулю, а в точках, расположенных левее этой характеристики, не равно нулю.


Перейдем к конструированию различных явных и неявных разностных схем. В качестве сетки D_h возьмем совокупность точек (x_{i},t^n), где x_{i}= ih,~ t^n=n\tau, i=0,\pm1,\pm2,\ldots; n=0,1,\ldots,N;~ T=N\tau. При построении разностных схем производную u_x можно аппроксимировать по формулам (8.20)–(8.22), а производную u_t — по формулам (8.27)–(8.29). Сформируем некоторые из них.


Для получения первой схемы применим (8.20), (8.27):


\begin{gathered}\frac{\widehat{u}_{i}^{n+1}-\widehat{u}_{i}^{n}}{\tau}+ \frac{\widehat{u}_{i+1}^{n}-\widehat{u}_{i}^{n}}{h}=0,\quad i=0,\pm1,\pm2,\ldots; n=0,1,2,\ldots,N-1;\\[2pt] \widehat{u}_{i}^{0}= \psi(x_{i}),\quad i=0,\pm1,\pm2,\ldots\end{gathered}

Полученную схему здесь и далее можно представить в форме \mathcal{L}_h\widehat{u}_{h}= f_h, где


\begin{aligned}& \mathcal{L}_h= \begin{cases}\frac{\widehat{u}_{i}^{n+1}-\widehat{u}_{i}^{n}}{\tau}+ \frac{\widehat{u}_{i+1}^{n}-\widehat{u}_{i}^{n}}{h}\,,& i=0,\pm1,\pm2,\ldots; n=0,1,2,\ldots,N-1;\\ \widehat{u}_{i}^0,& i=0,\pm1,\pm2,\ldots;\end{cases}\\[4pt] f_h= \begin{cases}0,& i=0,\pm1,\pm2,\ldots; n=0,1,2,\ldots,N-1;\\ \psi(x_{i}),& i=0,\pm1,\pm2,\ldots \end{cases} \end{aligned}

Решение на (n+1)-м слое запишется в виде


\widehat{u}_{i}^{n+1}= \widehat{u}_{i}^{n}\! \left(1+ \frac{\tau}{h}\right)-\widehat{u}_{i+1}^{n} \frac{\tau}{h},\quad i=0,\pm1,\pm2,\ldots; n=0,1,2,\ldots,N-1.
(8.35)

Узлы, участвующие в расчетах, образуют трехточечный шаблон H_3(x_{i},t^{n})= \bigl\{(x_{i},t^{n}), (x_{i+1},t^{n}), (x_{i},t^{n+1})\bigr\}, изображенный на рис. 8.5, а. Проанализируем эту схему, используя выводы, полученные при применении метода характеристик.


Если функция \psi(x) имеет вид \psi(x)= \begin{cases}\psi_0(x)\ne0,& x<0,\\ 0,& x \geqslant 0,\end{cases}, то \widehat{u}_{i}^{0}=0,~ i=0,1,2. Тогда по формуле (8.35) при любом соотношении между \tau и h получаем \widehat{u}_{i}^{0}=0,~ i=0,1,2;~ n=1,2,\ldots, т.е. в первом квадранте решение уравнения тождественно равно нулю. Это противоречит точному решению (см. рис. 8.4), следовательно, сходимость решения разностной задачи к решению дифференциальной задачи отсутствует. Таким образом, схема (8.35) имеет первый порядок аппроксимации по \tau и по h и является абсолютно неустойчивой. Следовательно, она не может быть использована на практике и приведена здесь из методических соображений.


Узлы трехточечного шаблона

Для получения второй схемы применим (8.21) и (8.27), т.е. в отличие от предыдущего случая производная их аппроксимируется левой разностью:


\begin{gathered}\frac{\widehat{u}_{i}^{n+1}-\widehat{u}_{i}^{n}}{\tau}+ \frac{\widehat{u}_{i+1}^{n}-\widehat{u}_{i}^{n}}{h}=0,\quad i=0,\pm1,\pm2,\ldots; n=0,1,2,\ldots,N-1;\\[2pt] \widehat{u}_{i}^{0}= \psi(x_{i}),\quad i=0,\pm1,\pm2,\ldots\end{gathered}

Решение на (n+1)-м слое запишется в виде


\widehat{u}_{i}^{n+1}= \widehat{u}_{i}^{n}\! \left(1-\frac{\tau}{h}\right)+ \widehat{u}_{i+1}^{n} \frac{\tau}{h},\quad i=0,\pm1,\pm2,\ldots; n=0,1,2,\ldots,N-1.
(8.36)

Узлы, участвующие в расчетах, образуют трехточечный шаблон H_3(x_{i},t^n)= \bigl\{ (x_{i-1},t^n), (x_{i},t^{n}), (x_{i+1},t^{n+1})\bigr\}, изображенный на рис. 8.5,б. Схема имеет первый порядок аппроксимации по \tau и по h и является условно устойчивой (является устойчивой при \frac{\tau}{h} \leqslant 1).

Правая часть формулы (8.36) является формулой линейной интерполяции. Она определяет значение функции, отстоящей от точки (x_{i},t^{n}) на расстояние, равное т. Это значение переносится по характеристике в точку (x_{i},t^{n+1}). Поведение характеристик в данной задаче изображено на рис. 8.6. Проанализируем более подробно свойства этой схемы для трех случаев соотношения шагов h и \tau сетки D_h.


При \frac{\tau}{h}=1 разностная схема в данной задаче дает точное решение, так как \widehat{u}_{i}^{n+1}= \widehat{u}_{i-1}^{n}, т.е. решение на следующем временном слое, передающееся вдоль характеристик (помечены на рисунке пунктирной линией), остается постоянным (рис. 8.6, а).


Разностные схемы для трёх случаях

При \frac{\tau}{h}<1 и \frac{\tau}{h}>1 по характеристике переносится значение, соответствующее точке, отмеченной крестиком (рис. 8.6, б и в). При \frac{\tau}{h}<1 разностная схема дает решение, соответствующее истинному с точностью, которой обладает линейная интерполяция. При \frac{\tau}{h}>1 можно получить решение, существенно отличающееся от точного. Например, пусть начальное условие задано в виде \psi(x)= \begin{cases}\psi_0(x)\ne0,& x<0,\\ 0,& x \geqslant 0,\end{cases}. Тогда \widehat{u}_{i}^{0}=0,~ i=0,1,2,\ldots, и, следовательно, \widehat{u}_{i}^{1}=0 при i=1. Точное решение в этом узле определяется значением \psi_0(h-\mu) (см. п.3 метода характеристик), поэтому погрешность может быть практически любой.


Для получения третьей схемы применим (8.22) и (8.27):


\begin{gathered}\frac{\widehat{u}_{i}^{n+1}-\widehat{u}_{i}^{n}}{\tau}+ \frac{\widehat{u}_{i+1}^{n}-\widehat{u}_{i-1}^{n}}{h}=0,\quad i=0,\pm1,\pm2,\ldots;\quad n=0,1,2,\ldots,N-1;\\[2pt] \widehat{u}_{i}^{0}=\psi(x_{i}),\quad i=0,\pm1,\pm2,\ldots \end{gathered}
(8.37)

Узлы, участвующие в расчетах, образуют четырехточечный шаблон H_4(x_{i},t^{n})= \bigl\{(x_{i-1},t^{n}), (x_{i},t^{n}), (x_{i+1},t^{n}), (x_{i},t^{n+1})\bigr\}, изображенный на рис. 8.5,в. Схема является абсолютно неустойчивой [35] и имеет первый порядок аппроксимации по \tau и второй — по h.


Для получения четвертой схемы (схемы Лакса) применим (8.22) и модификацию (8.27), где вместо u_{i}^{n} берется среднее \frac{u_{i+1}^{n}+ u_{i-1}^{n}}{2}\colon


\begin{gathered} \frac{\widehat{u}_{i}^{n+1}-\dfrac{\widehat{u}_{i+1}^{n}+ \widehat{u}_{i-1}^{n}}{2}}{\tau}+ \frac{\widehat{u}_{i+1}^{n}-\widehat{u}_{i-1}^{n}}{2h}=0,\quad i=0,\pm1,\pm2,\ldots;\quad n=0,1,2,\ldots,N-1;\\[2pt] \widehat{u}_{i}^{0}=\psi(x_{i}),\quad i=0,\pm1,\pm2,\ldots \end{gathered}
(8.38)

Узлы, участвующие в расчетах, образуют трехточечный шаблон H_3(x_{i},t^{n})= \bigl\{(x_{i-1},t^{n}), (x_{i+1},t^{n}), (x_{i},t^{n+1})\bigr\}, изображенный на рис. 8.5,г. Схема является условно устойчивой (является устойчивой при \frac{\tau}{h} \leqslant 1). Заметим, что простейшая модификация схемы (8.37), а именно использование среднего значения, привела к тому, что разностная схема из абсолютно неустойчивой превратилась в условно устойчивую. Формулы (8.35)–(8.38) явно выражают решение на (n+1)-м временном слое через решение на n-м слое. Поэтому соответствующие им схемы называются явными двухслойными.


Для получения пятой схемы применим формулу (8.20), записанную на (n+1)-м временном слое, и (8.27):


\begin{gathered}\frac{\widehat{u}_{i}^{n+1}-\widehat{u}_{i}^{n}}{\tau}+ \frac{\widehat{u}_{i+1}^{n+1}-\widehat{u}_{i}^{n+1}}{h}=0,\quad i=0,\pm1,\pm2,\ldots;\quad n=0,1,2,\ldots,N-1;\\[2pt] \widehat{u}_{i}^{0}=\psi(x_{i}),\quad i=0,\pm1,\pm2,\ldots \end{gathered}
(8.39)

Эта схема относится к неявным, использует трехточечный шаблон H_3(x_{i},t^{n})= \bigl\{(x_{i},t^{n}), (x_{i},t^{n+1}), (x_{i+1},t^{n})\bigr\} (рис. 8.5,д) и является абсолютно устойчивой. Схема (8.39) имеет первый порядок аппроксимации по \tau и по h.


Для получения шестой схемы применим (8.21), записанную на (n+1)-м временном слое, и (8.27):


\begin{gathered}\frac{\widehat{u}_{i}^{n+1}-\widehat{u}_{i}^{n}}{\tau}+ \frac{\widehat{u}_{i}^{n+1}-\widehat{u}_{i-1}^{n+1}}{h}=0,\quad i=0,\pm1,\pm2,\ldots;\quad n=0,1,2,\ldots,N-1;\\[2pt] \widehat{u}_{i}^{0}=\psi(x_{i}),\quad i=0,\pm1,\pm2,\ldots \end{gathered}
(8.40)

Эта схема также относится к неявным, использует трехточечный шаблон H_3(x_{i},t^{n})= \bigl\{(x_{i-1},t^{n+1}), (x_{i},t^{n+1}), (x_{i},t^{n})\bigr\} (рис. 8.5,е) и является абсолютно устойчивой. Схема имеет первый порядок аппроксимации по \tau и по h.


Для получения седьмой схемы будем использовать средние значения выражений, определяемых формулой (8.20), записанной на (n+1)-м и n-м временных слоях, и формулой (8.27), записанной при x=x_{i} и x=x_{i+1}\colon


\begin{gathered}\frac{\dfrac{\widehat{u}_{i}^{n+1}-\widehat{u}_{i}^{n}}{\tau}+ \dfrac{\widehat{u}_{i+1}^{n+1}-\widehat{u}_{i+1}^{n}}{\tau}}{2}+ \frac{\dfrac{\widehat{u}_{i+1}^{n}-\widehat{u}_{i}^{n}}{h}+ \dfrac{\widehat{u}_{i+1}^{n+1}-\widehat{u}_{i}^{n+1}}{h}}{2}=0,\quad i=0,\pm1,\pm2,\ldots;\quad n=0,1,2,\ldots,N-1;\\[2pt] \widehat{u}_{i}^{0}=\psi(x_{i}),\quad i=0,\pm1,\pm2,\ldots \end{gathered}

Эта схема относится к неявным, использует четырехточечный шаблон H_4(x_{i},t^{n})= \bigl\{(x_{i},t^{n}), (x_{i},t^{n+1}), (x_{i+1},t^{n}), (x_{i+1}, t^{n+1})\bigr\} (рис. 8.5,ж) и является абсолютно устойчивой. Схема имеет второй порядок аппроксимации по \tau и по h. Второй порядок аппроксимации обеспечивается путем применения оператора осреднения аппроксимационных операторов на (n+1)-м и л-м временных слоях.


Пятая, шестая и седьмая схемы относятся к неявным двухслойным.

Математический форум (помощь с решением задач, обсуждение вопросов по математике).
Кнопка "Поделиться"

Часовой пояс: UTC + 3 часа [ Летнее время ]


Яндекс.Метрика

Copyright © 2010-2019 MathHelpPlanet.com. All rights reserved