Научная статья на тему 'Схема расщепления для задач фильтрации газа с химическими реакциями'

Схема расщепления для задач фильтрации газа с химическими реакциями Текст научной статьи по специальности «Математика»

CC BY
115
52
i Надоели баннеры? Вы всегда можете отключить рекламу.
Область наук

Аннотация научной статьи по математике, автор научной работы — Снытников В. Н., Юрченко Е. М.

Работа выполнена при частичной поддержке NASA, грант JURRISS № NAG5-8717. Исследуется возможность применения схем расщепления для решения задач динамики реакционноспособной газовой среды. Реализована одна из схем расщепления Г.И. Марчука второго порядка точности по времени с различными вариантами разностных схем на дробных шагах. Сравнение численных результатов, полученных по схемам расщепления, с аналитическими нелинейными решениями показало первый порядок сходимости численного решения к точному. Установлено, что схемы расщепления адекватно передают профиль волны химической реакции, обеспечивают допустимую погрешность при расчете температуры и, таким образом, могут рассматриваться как основа при создании более сложных алгоритмов решения задач физико-химической газодинамики.

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Splitting scheme for the problems of gas filtration with chemical reactions

The applicability of application of splitting schemes to the solution of problems on dynamics of reactive gaseous medium is investigated. One of the splitting schemes suggested by G. I. Marchuk has been implemented with various versions of difference schemes on fractional steps. The comparison of numerical results obtained by splitting schemes with analytical nonlinear solution showed the first order of convergence of the numerical solution to the exact solution. It has been found that the splitting schemes adequately represent the wave profile of a chemical reaction with an admissible error of temperature calculation and therefore can be considered a reliable basis in working out more complicated solution algorithms for the problems of physical and chemical gasdynamics.

Текст научной работы на тему «Схема расщепления для задач фильтрации газа с химическими реакциями»

Вычислительные технологии

Том 6, № 5, 2001

СХЕМА РАСЩЕПЛЕНИЯ ДЛЯ ЗАДАЧ ФИЛЬТРАЦИИ ГАЗА С ХИМИЧЕСКИМИ

РЕАКЦИЯМИ*

В. Н. Снытников Институт катализа им. Г. К. Борескова СО РАН, Новосибирск, Россия

Е.М. ЮРЧЕНКО

Институт вычислительных технологий СО РАН, Новосибирск, Россия

e-mail: snyt@catalysis.nsk.su

The applicability of application of splitting schemes to the solution of problems on dynamics of reactive gaseous medium is investigated. One of the splitting schemes suggested by G.I. Marchuk has been implemented with various versions of difference schemes on fractional steps. The comparison of numerical results obtained by splitting schemes with analytical nonlinear solution showed the first order of convergence of the numerical solution to the exact solution. It has been found that the splitting schemes adequately represent the wave profile of a chemical reaction with an admissible error of temperature calculation and therefore can be considered a reliable basis in working out more complicated solution algorithms for the problems of physical and chemical gasdynamics.

Введение

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

Особенностью рассматриваемых моделей является одновременное наличие в них конвективного и диффузионного переноса, изменение концентраций реагентов в химических реакциях вместе с выделением или поглощением энергии. Каждый из процессов в отдельности изучен достаточно подробно. Так, процессы переноса описываются гиперболическими уравнениями ([1], см. приведенную в работе библиографию), а процессы распространения тепла (диффузии вещества) — уравнениями параболического типа. Для более сложных математических моделей была предложена идея построения численных методов так, чтобы на каждом шаге рассматривался один хорошо исследованный процесс [2]. Такой подход известен как метод дробных шагов с расщеплением по физическим процессам

* Работа выполнена при частичной поддержке NASA, грант JURRISS №NAG5-8717.

© В.Н. Снытников, Е.М. Юрченко, 2001.

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

Методы расщепления были развиты, в частности, в работах [2, 3]. Известны работы по исследованию схем расщепления для нелинейных задач газовой динамики [4]. В то же время для задач горения возможности таких схем изучены сравнительно слабо. Более того, применимость метода расщепления к решению подобных задач часто ставится под сомнение из-за экспоненциально быстрого изменения функциями своих значений при малых изменениях, например, температуры в скоростях химических реакций. Обычно задачи динамики газа в реакционноспособной среде решаются численно с применением всевозможных простых схем (например, [5]) даже без рассмотрения вопросов корректности применения и устойчивости схем. Трудность в исследовании свойств схем для таких задач заключается в сильной нелинейности уравнений. Кроме того, математическая модель реакции — диффузии — может иметь решения волнового типа (волны горения или детонации) , в то время как отдельные операторы (диффузии или источников) такими свойствами не обладают. К дополнительным проблемам в исследовании разностных методов следует отнести известность лишь узкого круга аналитических решений для нетривиальных моделей теории горения, на которых можно отрабатывать численные схемы.

В представленной работе исследовалась возможность применения схем расщепления на примере схемы Г. И. Марчука [2] для решения задач теории горения. В первую очередь схема расщепления была использована для решения неоднородного уравнения с переносом, диффузией и нелинейным экспоненциальным слагаемым типа приближения Франк-Каменецкого для скорости химической реакции. Для этого уравнения известно нетривиальное аналитическое решение волнового типа, с которым можно сравнивать численные результаты. На втором этапе изучалась квазигомогенная модель фильтрации реагента через пористый слой катализатора с одной необратимой химической реакцией. Для этой модели ранее было найдено автоволновое решение [6], которое использовано в работе для оценки точности численных схем. В работе также приведены примеры расчета нестационарных задач с использованием созданных алгоритмов.

1. Схема расщепления, варианты реализации

Будем искать численное решение начально-краевой задачи в области П = [0, то] х [xl, xr], ограниченной по времени интервалом (0, tmax), на регулярной сетке ит х Uh, где ит = {ti = гт, г = 0,...,M, т = tmajM}, Uh = {х3 = jh, j = 0,...,N, h = (Xr - XL)/N}. Схема расщепления [2] выглядит следующим образом:

' (E + 2Л1У-2/3 = (E - 2Лг)^г-1, (E + 2Л2)^г-1/3 = (E - Т-Л2)>Г2/3,

рг+1/3 = рг-1/3 + 2 т f, (1)

(E + 2Л2^г+2/3 = (E - 2Л2^г+1/3, (E + 2Л1)^г+1 = (E - 2Л1)^г+2/3-

Эта схема представляет собой абсолютно устойчивый разностный аналог второго порядка аппроксимации по т эволюционной задачи

—и

— + Аи = / (¿), и = д при г = 0

в предположении, что А — положительно полуопределеннный дифференциальный оператор, т.е. (Ар, и) > 0. Оператор А представляется в виде суммы операторов более простой структуры:

A = £ Aa, (2)

а=1

п

где Аа = Аа(г) > 0. Предполагаем, что Л-7 = ^ Ла, где Ла = Аа(^). В случае уравнений

а=1

фильтрации оператор Л состоит из двух слагаемых: Л = Л1 + Л2, которые аппроксимируют первую и вторую производные по пространству соответственно. Таким образом, на первом шаге считается перенос, на втором — диффузия, а на третьем шаге, относительно которого схема симметрична, рассчитывается вклад стоков или источников энергии, который описывается функцией /.

В схеме расщепления численное решение конвективного переноса было реализовано тремя способами: направленной разностью, центральной разностью и с помощью ЕСТ-ал-

горитма. Центральный дробный шаг, на котором решается нелинейное уравнение =

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

Результаты, полученные по схеме расщепления, сравнивались с простой неявной схемой:

- и,*1- Л - 2^+1 + ^+1

и-^ + Vй , = Б^1-% + ^ + /](и), (3)

т п П2 ^

где / = /(и*). Такая схема имеет порядок аппроксимации 0(т + П).

2. Решение начально-краевой задачи для одномерного неоднородного нелинейного уравнения

Рассматривается одномерное неоднородное нестационарное уравнение

-и -и ,-2и .

Ж + 1'ТХ. = кзХ? + / (и), (4)

где /(и) = а + ЬвХи. Это уравнение обладает известным аналитическим решением типа простой волны и = и(х — иг), где и = V ± ([7]):

2

u(t, x) = — — ln л

±ß + cexp ^±jjVkx — T ßVkv^j ^

Приведем уравнение (4) к безразмерному виду преобразованием

t = tot', x = Lx', u = uou'.

Таблица 1 Погрешности решений, полученных разными методами

Ки = •т/Л = 1

т = 0.025 Л = 0.05 т = 0.0125 Л = 0.025 т = 0.0025 Л = 0.005

2.9 ■ 10-2 1.5 ■ 10-2 2.9 ■ 10-3

9.7 ■ 10-3 4.5 ■ 10-3 9.8 ■ 10-4

1.0 ■ 10-4 2.6 ■ 10-5 1.2 ■ 10-6

3.9 ■ 10-3 1.9 ■ 10-3 3.7 ■ 10-4

8.4 ■ 10-3 4.2 ■ 10-3 8.3 ■ 10-4

8.8 ■ 10-5 2.2 ■ 10-5 8.9 ■ 10-7

8.4 ■ 10-3 4.2 ■ 10-3 8.3 ■ 10-4

Схема

Ки = ит/Л = 0.1

т = 0.025 Л = 0.05

т = 0.0125 Л = 0.025

т = 0.0025 Л = 0.005

1.2

10

-2

5.9

10

-3

1.2 10

-3

9.7

10

-3

4.9

10

-3

9.8 • 10

-4

4.0

10

-5

1.0

10

-5

4.01 10

-7

5.3

10

-3

2.6

10

-3

5.3 10

-4

8.4

10

-4

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

4.2

10

-4

8.3 10

-5

6.9

10

-5

1.7

10

-5

6.9 10

-7

8.4

10

-4

4.2

10

-4

8.3 10

-5

Для этого выберем масштабы скорости V, длины Ь и времени Ь0 так, чтобы Ь0 = ЬД>, а также и0 = — 2/Л. Переходя в уравнении к новым переменным и опуская штрихи, получаем

ди ди 2и „ - ) дЬ дх дх2 '

где /(и) = —9 + 7в-2и. Решение уравнения принимает вид

(6)

и(Ь, х) = 1п

+ Сехр (±^Шх — Ь(9 ±Щ

(7)

где Б = к/¿V, 9 = Ь0^2, 7 = 9в2, параметры 9 и 7 только положительны. Таким образом, уравнение (6) зависит от трех параметров, в отличие от исходного, зависящего от пяти.

В табл. 1 представлены значения максимальных отклонений численного решения начально-краевой задачи для уравнения (6) от аналитического. Параметры сетки включены в таблицу, остальные параметры задачи выбраны следующими: Б = 1, 7 =1, 9 = 1, х € [0, 6], Ь € [0, 4]. Цифрами обозначены схемы:

1. Простой разностный аналог.

2. Схема расщепления (направленная разность, точное решение).

3. Схема расщепления (центральная разность, точное решение).

4. Схема расщепления (направленная разность, предыдущий слой).

5. Схема расщепления (центральная разность, предыдущий слой).

6. Схема расщепления (ЕСТ-алгоритм, точное решение).

7. Схема расщепления (ЕСТ-алгоритм, предыдущий слой).

Как видно из представленных результатов, схема расщепления при использовании центральной разности для аппроксимации конвективного переноса одновременно с заданием аналитического решения для нелинейного уравнения д^/дЬ = /(<р) (в таблице строка под номером 3) обеспечивает порядок 0(т2 + к2), что позволяет сделать вывод о правильности реализации алгоритма. Действительно, если поставить в правую нелинейную часть известное решение, то уравнение перестанет быть нелинейным. Мы получим именно такое

1

2

3

4

5

6

7

уравнение, для которого схема расщепления строилась как схема второго порядка аппроксимации по времени. Центральная разность позволяет иметь согласованный с порядком по времени порядок по пространству 0(П2). Эта же схема во всех остальных случаях обеспечивает первый порядок. Исключение составляет алгоритм, в котором для аппроксимации переноса применялся ЕСТ-алгоритм (для нелинейного уравнения использовалось известное решение) — он дает не только порядок 0(т2 + П2), но и лучший результат, т.е. наименьшее отклонение от известного решения. Простая неявная схема при вычислениях показала отклонения порядка 0(т + П), что и следовало ожидать.

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

Несмотря на то, что схема расщепления дает в общей сложности первый порядок аппроксимации и сходимости, ее использование для данной задачи оказалось более эффективным. Она позволяет применять для получения заданной точности более крупную сетку по сравнению с простой неявной схемой. Например, максимальное отклонение тах 8vh ~ 4.2- 10-3 с помощью простой неявной схемы можно получить на 240 пространственных точках и 3200 временных слоях, а по схеме расщепления (в таблице строка под номером 7) та же точность получается при т = 0.0125, П = 0.025, что составляет 320 временных слоев и 240 пространственных точек.

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

Дальнейшее исследование свойств схемы расщепления производилось на следующей системе уравнений:

-М -М ,

-Т -Т - 2Т

С8Рв(1 — е)+ есдрдV— = + (^К^, N), (8)

( Е

К(Т, N) = Мк ехр( — т

где р3,о3 и Л — плотность, теплоемкость и теплопроводность твердой фазы; V — скорость фильтрации; рд и сд — плотность и теплоемкость газа; ( — тепловой эффект реакции, скорость которой К(Т, N).

Эта система уравнений представляет квазигомогенную модель фильтрации газа в пористой среде с химической реакцией [6]. Для заданной пористости е она описывает изменение концентрации N реагента с его низким содержанием в газовой фазе вместе с изменением температуры среды Т в предположении интенсивного теплообмена между газом и твердой фазой.

Система (8) обладает автомодельным решением [6]. Оно получено при решении полной краевой задачи

¿У 1 — У (

¿т = * т-у ехр1— Т,

(9)

Т ^ 0, У ^ 0, т 1, У 1

обозначениях

N T E

Y = 1 T =—, Z = —,

N0 Tn Tn

в

S= (lz£) 2 = (1 - z)csPs ß = QNo

v (ev)2 QNo\ y - l) ' Y £cgpg ' TnCgpg'

Решение зависит от модифицированного числа Зельдовича Z и числа Sv, имеющего достаточно общий смысл для гетерогенных сред с химическими реакциями; Sv — это соотношение между потоком тепла, переносимого за счет теплопроводности на масштабе химической реакции, и потоком энергии, связанной с химической реакцией с поправочным коэффициентом порядка единицы, зависящим от теплофизических свойств двухфазной среды.

График автомодельного решения — зависимость степени превращения от температуры — представлен на рис. 1 (Z = 5, ln Sv = 6.659). На рис. 2 приведены соответствующие графики зависимости температуры и концентрации от автоволновой переменной £ = x — ut.

Приведем систему (8) к безразмерным величинам:

t = tot'. x = Lx', T = TnT'' N = NoN'. to = L/v.

Опуская штрихи после замены, получим систему с меньшим количеством параметров:

dN + д^ = щ dt dx '

dT dT

dt Vt dx

д2 T Dt— + OW,

W = N exp ( — Z

dx2

(10)

в обозначениях

Z = E/Tn, vt = eCgPg ,, Csps (1 — e)

kL \ QkNoL П = —, DT =---——, O

ev' CsPs(1 — e)LV CsPs(1 — e)T,nv'

0 0.2 0.4 0.6 0.8 1

Рис. 1. Автомодельное решение.

N Т 1.0 0.8 0.6 0.4 0.2 0.0 -0.2

-16 -14 -12 -10 -8 -6 -4 -2"

Рис. 2. Концентрация N(£) и температура Т(£).

;

; N /

; ! / » /

; / в # / /

; / /

; Т

"мм мм мм мм 1 1 1 1 1111 мм

Для системы (10) схема расщепления выглядит следующим образом:

' (Е + N1-2/3 = (Е — Nг-1, (е + 2VTЛ^ Г-2/3 = (Е — 2^1) Г-1, (Е — 2БтЛ^ Тг-1/3 = (е + 2БтЛ^ Г-2/3, N *+1/3 = N*-1/3 — 2тп • Ц Т1), Тг+1/3 = Т1-1/3 + 2тв • Ц1^ *,Т*),

(Е — 2 БтЛ^ Тг+2/3 = (Е + 2 Г+1/3,

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

(Е + 2vтЛг) Т+1 = [Е — 2ит^) Г+2/3, (Е + 2Л1 ^^ = (Е — 2л1)Nг+2/3,

где Л1, Л2 — операторы, аппроксимирующие -/-х и -2/-х2 соответственно. Для сравнения взята простая неявная схема

(11)

щ+1 — щ Щ+1 — л . .

-1-1 + ——:—— = —пЦ (N^тл,

т .п

Т^+1_т^ Т^+1_т^+1

2 2 + Vт - 2-1

т

п

Б

т

Т^+1 _ 2Т^+1 + т^+1

^--2 + 2-1 + ец (Щ}т*).

п2

(12)

4. Результаты расчетов

Для значений скорости волны и = 5 • 10-4 и и = — 1, 01 • 10-3 проведено две серии расчетов с начальными автомодельными профилями N и Т. Обе волны, хотя и направлены в разные стороны, обладают очень маленькой скоростью, т. е. решение задачи практически стационарное. В табл. 2 приведены результаты расчетов. Цифрами в ней обозначены: 1 — схема расщепления с направленной разностью, 2 — схема расщепления с ЕСТ-алгоритмом, 3 — простая неявная схема. Видно, что все схемы показали первый порядок сходимости, а точность вычислений существенно зависела от способа расчета конвективного переноса. Так, при использовании направленной разности погрешность расчета концентрации составила

около 2 % на самой крупной сетке, а температуры — 5 • 10-4 %. Схема с ЕСТ-алгоритмом дает примерно одинаковую погрешность для обеих функций: концентрации — 0.1 %, температуры — 0.04 %. При и = -1.01 • 10-3 погрешности несколько больше. Например, на самой крупной сетке для направленной разности концентрация — 2 %, температура — 0.2 %. Для ЕСТ-алгоритма — 2.5 и 1.3% для концентрации и температуры соответственно. Неявная схема показала такие же результаты, как и схема расщепления с направленной разностью. Таким образом, можно сказать, что для малых скоростей волны схемы расщепления применимы и дают результаты не хуже, чем простая неявная схема.

Таблица 2

Отклонения max 5

Схема т = 0.01, px = 150 т = 0.005, px = 300 т = 0.0025, px = 600

N T N T N T

и = 5 ■ 10-4

1 2, 09 ■ 10-2 4,54■10-6 1, 08 ■ 10-2 2, 33 ■ 10-6 5, 53■10-3 1, 3 ■ 10-6

2 1, 94 ■ 10-3 4, 27■10-4 9, 23 ■ 10-4 1,13 ■ 10-4 5,66■10-4 3, 63 ■ 10-5

3 2, 09 ■ 10-2 4,54■10-6 1, 08 ■ 10-2 2, 33 ■ 10-6 5, 53■10-3 1, 3 ■ 10-6

и = -1, 01 ■ 10-3

1 2, 73 ■ 10-2 2,11■10-3 1, 39 ■ 10-2 1,11 ■10-3 7, 03 ■ 10-3 5, 67 ■ 10-4

2 2,47■10-2 1,26■10-2 9, 97 ■ 10-3 3, 62 ■ 10-3 3,58■10-3 9, 30 ■ 10-4

3 2, 73 ■ 10-2 2,11■10-3 1, 39 ■ 10-2 1,11 ■10-3 7, 03 ■ 10-3 5, 67 ■ 10-4

Для больших значений скоростей ш = ±1 все схемы показали удовлетворительные результаты, за исключением схемы с FCT-алгоритмом. Применимость последней схемы при больших скоростях требует дополнительного исследования.

Интересен случай задания для системы (10) произвольных начальных данных. Будет ли решение иметь форму волны или нет? Если да, то какая будет скорость волны? Чтобы проверить это, зададим начальные условия в виде

N (»•«) = N0 + ИТГ • Т (0.x) = Ti + To + J--/+ 1 + AT«-*'»1. (13)

где Ti, Ni — значения температуры и концентрации в момент времени t = 0 на левом краю, а То + Ti, N0 — на правом краю; а, ß, A, r, X0, Xi, X2 — константы, определяющие вид начального распределения температуры и концентрации. Расчеты проведены по схеме расщепления Г. И. Марчука с направленными разностями при следующих значениях параметров: Z = 5, п = 6.7■ 103, в = 6■ 103, а = 6, ß =15, A = 1.5, r = 0.3, N0 = 0, T0 = 0.5, Ni = 1, Ti = 5 ■ 10-4. Остальные параметры отличались для различных вариантов расчетов и равнялись: vT = 1, DT = 0.001, X0 = 3.7, Xi = 5, X2 = 4.8 для варианта 1 (рис. 3); vT = 0.01, DT = 0.1, X0 = 11.7, Xi = 13, X2 = 12.8 для варианта 2 (рис. 4); вариант 3 (рис. 5) отличается от варианта 2 значением коэффициента теплопроводности Dt = 1.

На рис. 3 — 5 приведены профили концентрации (сплошные линии) и температуры (пунктирные линии) в моменты времени ti = 1, t2 = 2, t3 = 3, t4 = 4, t5 = 5, для всех трех вариантов. Жирными линиями показаны профили, соответствующие начальному

Рис. 3. Вариант расчета 1.

N Т

-

- Ж | V* ^хФ- : \

- 1 1 А : • «1 : • »» • 1&

- ¿5 и ¿3 1 1 1 '| 1 1 1 г »' \ Щ : 1 / 1 IV-1- 11

х

10

15

20

3

2

1

0

5

Рис. 4. Вариант расчета 2.

моменту времени. При рассмотрении рис. 3 видно, что в близкие к начальному моменты времени £ ^ 1 в зоне высокой температуры протекает реакция с выделением тепла. Количество реагента в этой зоне, несмотря на его приток с левой границы, уменьшается, фронт концентрации смещается влево. В то же время конвективный перенос тепла уводит вправо зону высоких температур. Начиная с некоторого момента времени температура на фронте становится близкой к автоволновому значению. С этого момента фронт перемещается вправо. Поскольку высокотемпературная зона смещается вправо со скоростью, большей скорости фронта, падая по амплитуде за счет теплопроводности, значение температуры на фронте волны уменьшается. Скорость волны при этом медленно меняется от своего значения ^ ^^ 1. Полученный к моменту времени £ = 5 режим еще далек от

-

- п 11 1 1 1 1 1 1 1 1 1 1 —1'\ 1

1 ~ $ - ч 1 - - - 1 1 1 1 1 1 1 1 1 / г _ 11

0 5 10 15 20

Рис. 5. Вариант расчета 3.

установления. Снижение скорости конвективного переноса тепла с одновременным повышением теплопроводности меняет скорость волны на отрицательное значение (рис. 4, 5). При этом более низкая теплопроводность в варианте 2 ведет к более высокому значению температуры Tf ^ 1.5 за фронтом волны по сравнению с вариантом 3 (Tf ~ 1.15). На момент времени ¿3 волна имеет установившиеся профиль и другие характеристики. Скорость распространения волны ш совпадает с теоретической скоростью автомодельной волны, которая вычисляется по формуле ш = (пУуTf — в)/(цTf — в) [6]. Таким образом, при задании произвольных начальных данных решение в одной локальной области принимает форму волны, профиль которой зависит, в частности, от коэффициента теплопроводности Ду.

Сравнение результатов расчетов, полученных с помощью схемы расщепления и простой неявной схемы, показало, что при уменьшении числа Куранта разница между численными решениями сокращается. В частности, при параметрах сетки т =1 • 10-4, рх = 120 эта разница не превышала 0.24%. Такой результат согласуется с теорией и позволяет сделать вывод о том, что получаемые решения сходятся к решению задачи.

Заключение

Расчеты, проведенные по схемам расщепления для решения нестационарных задач филь-тации газа с химическими реакциями, показали, что:

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

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

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

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

Список литературы

[1] Карамышев В. Б. Монотонные схемы и их приложения в газовой динамике. Новосибирск, 1994.

[2] ЯнЕнко Н. Н. Метод дробных шагов решения многомерных задач математической физики. Новосибирск: Наука, 1967.

[3] Марчук Г. И. Методы расщепления для решения нестационарных задач // ЖВМ и МФ. 1995. Т. 35, №6.

[4] КовЕня В. М., ЯнЕнко Н. Н. Метод расщепления в задачах газовой динамики. Новосибирск: Наука, 1981.

[5] Новикова И. И., Писаренко В. Н., Жукова Т. Б. Интенсификация гетерогенно-каталитического процесса гидрирования ацетиленовых углеводородов в изопреновой фракции // ТОХТ. 2000. Т. 34, №1. С. 87-91.

[6] Снытников В. Н. Уединенные бегущие волны в гетерогенной среде с химическими реакциями // Физика горения и взрыва. 1998. Т. 34, №3. С. 29-36.

[7] Зайцев В.Ф., Полянин А. Д. Дифференциальные уравнения с частными производными: Справочник. М.: Междунар. программа, 1996.

Поступила в редакцию 10 апреля 2001 г., в переработанном виде —18 июня 2001 г.

i Надоели баннеры? Вы всегда можете отключить рекламу.