Является ли Кранк-Николсон устойчивой схемой дискретизации для уравнения реакция-диффузия-адвекция (конвекция)?

Я не очень знаком с общими схемами дискретизации для PDE. Я знаю, что Кранк-Николсон является популярной схемой для дискретизации уравнения диффузии. Это также хороший выбор для адвекции?

Я заинтересован в решении уравнения реакции-диффузии-адвекции ,

$\frac{\partial u}{\partial t} + \nabla \cdot \left( \boldsymbol{v} u - D\nabla u \right) = f$

где - коэффициент диффузии вещества а - скорость.$D$$u$$\boldsymbol{v}$

Для моего конкретного приложения уравнение может быть записано в виде,

$\frac{\partial u}{\partial t} = \underbrace{D\frac{\partial^2 u}{\partial x^2}}_{\textrm{Diffusion}} + \underbrace{\boldsymbol{v}\frac{\partial u}{\partial x}}_{\textrm{Advection (convection)}} + \underbrace{f(x,t)}_{\textrm{Reaction}}$

Вот схема Кранка-Николсона, которую я применил,

$\frac{u_{j}^{n+1} - u_{j}^{n}}{\Delta t} = D \left[ \frac{1 - \beta}{(\Delta x)^2} \left( u_{j-1}^{n} - 2u_{j}^{n} + u_{j+1}^{n} \right) + \frac{\beta}{(\Delta x)^2} \left( u_{j-1}^{n+1} - 2u_{j}^{n+1} + u_{j+1}^{n+1} \right) \right] + \boldsymbol{v} \left[ \frac{1-\alpha}{2\Delta x} \left( u_{j+1}^{n} - u_{j-1}^{n} \right) + \frac{\alpha}{2\Delta x} \left( u_{j+1}^{n+1} - u_{j-1}^{n+1} \right) \right] + f(x,t)$

Обратите внимание на термины и . Это позволяет схеме перемещаться между:$\alpha$$\beta$

• $\beta=\alpha=1/2$ Crank-Niscolson,
• $\beta=\alpha=1$ полностью неявно
• $\beta=\alpha=0$ полностью явно

Значения могут быть разными, что позволяет диффузионному члену быть Crank-Nicolson, а адвекционным термином - что-то еще. Какой самый стабильный подход, что бы вы порекомендовали?

Это хорошо сформулированный вопрос и очень полезная вещь для понимания. Коррок прав, когда отсылает вас к анализу фон Неймана и книге Левека. Я могу добавить немного больше к этому. Я хотел бы написать подробный ответ, но на данный момент у меня есть время только на короткий:

С вы получаете метод, который является абсолютно стабильным для произвольно больших размеров шагов, а также с точностью до второго порядка. Однако метод не является L- стабильным, поэтому очень высокие частоты не будут затухать, что является нефизическим.$\alpha=\beta=1/2$

С вы получаете метод, который также безусловно стабилен, но только с точностью до 1-го порядка. Этот метод очень диссипативный. Это L- стабильный.$\alpha=\beta=1$

Если вы берете , ваш метод можно понимать как применение аддитивного метода Рунге-Кутты к полудискретизации с центрированной разностью. Анализ стабильности и точности для таких методов значительно сложнее. Очень хорошая статья о таких методах здесь .$\alpha\ne\beta$

Какой подход рекомендовать, сильно зависит от величины , вида исходных данных, с которыми вы имеете дело, и точности, которую вы ищете. Если допустима очень низкая точность, то - очень надежный подход. Если умеренный или большой, то проблема диффузионно-доминирующая и очень жесткая; обычно даст хорошие результаты. Если очень мало, то для конвективных членов может быть выгодно использовать явный метод и разматывание более высокого порядка.$D$$\alpha=\beta=1$$D$$\alpha=\beta=1/2$$D$

Вообще говоря, вы захотите использовать неявный метод для параболических уравнений (диффузионная часть) - явные схемы для параболических уравнений в частных производных должны иметь очень короткий временной шаг, чтобы быть стабильным. И наоборот, для гиперболической части (адвекции) вам понадобится явный метод, так как он дешевле и не нарушает симметрию линейной системы, которую вы должны решить, используя неявную схему для диффузии. В этом случае вы хотите избежать центрированных различий, таких как и переключиться на односторонние различия по соображениям стабильности.$(u_{j+1}-u_{j-1})/2\Delta t$$(u_j-u_{j-1})/\Delta t$

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

Основная идея состоит в том, чтобы предположить, что в вашем дискретном приближении можно записать сумму плоских волн , где - волновое число, а - частота. Вы вбиваете плоскую волну в свое приближение к PDE и молитесь, чтобы она не взорвалась. Мы можем переписать плоскую волну как и мы хотим убедиться, что .$e^{i(k j\Delta x-\omega n\Delta t)}$$k$$\omega$$\xi^n e^{ikj\Delta x}$$|\xi|\le 1$

В качестве иллюстрации рассмотрим обыкновенное уравнение диффузии с полностью неявным дифференцированием:

$\frac{u^{n+1}_j-u^n_j}{\Delta t} = D\frac{u^{n+1}_{j-1}-2u^{n+1}_j+u^{n+1}_{j+1}}{\Delta x^2}$

Если подставить в плоскую волну, а затем разделить на и , мы получим уравнение$\xi^n$$e^{ikj\Delta x}$

$\frac{\xi-1}{\Delta t} = D\frac{e^{-ik\Delta x}-2+e^{ik\Delta x}}{\Delta x^2}\xi$

Очистите это немного сейчас, и мы получим:

$\xi = \frac{1}{1+\frac{2D\Delta t}{\Delta x^2}\left(1-\cos k\Delta x\right)}$ .

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

$\frac{u_j^{n+1}-u_j^n}{\Delta t} = v\frac{u_{j-1}^n-u_{j+1}^n}{2\Delta x}$

и посмотреть , что вы получите. (На этот раз у него будет мнимая часть.) Вы обнаружите, что , что печально. Отсюда мое предупреждение, что вы не используете его. Если вы можете сделать это, то у вас не должно быть особых проблем с поиском стабильной схемы для полного уравнения адвекции-диффузии.$\xi$$|\xi|^2 > 1$

Тем не менее, я бы использовал полностью неявную схему для диффузионной части. Измените разность в адвективной части на если и если и выберите временной шаг так, чтобы , (Это условие Куранта-Фридрихса-Леви .) Оно является точным только для первого порядка, поэтому вы можете посмотреть схемы дискретизации более высокого порядка, если это вас касается.$u_j-u_{j-1}$$v > 0$$u_j-u_{j+1}$$v < 0$$V\Delta t/\Delta x \le 1$