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

4.3.1. Системы типа «реакция – диффузия»

Рассмотрим $s$ компонент, реагирующих между собой в ходе $R$ независимых реакций при постоянной температуре. Изменение концентрации во времени и пространстве может быть описано системой локальных уравнений баланса массы для $s$ выбранных компонент
\[
\frac{\partial c_{i}}{\partial t}+\operatorname{div} \mathbf{j}_{i}=f_{i} .
\]

Здесь $c_{i}$-молярная концентрация компоненты $i$, $f_{i}$ описывает возникновение $i$-й компоненты в результате $R$ независимых реакций (при этом $f_{i}=\sum_{j=1}^{R} v_{i j} r_{j}$, где $v_{i j}$-стехиометрический коэффициент $i$-й компоненты в $j$-й реакции и $r_{j}$-скорость $j$-й реакции), $\mathbf{j}_{i}$ – молярная плотность потока компоненты $i$. Если учитывать только диффузионные и конвективные составляющие потока, то имеет место соотношение
\[
\mathbf{j}_{i}=\mathbf{j}_{d i}+\mathbf{v} c_{i},
\]

где $\mathbf{v}$-вектор локальной мгновенной скорости смеси; в дальнейшем мы будем полагать $\mathbf{v}=\mathbf{0}$.

Диффузионный поток $\mathbf{j}_{\mathrm{d} i}$ линейно зависит от градиентов концентраций:
\[
\mathbf{j}_{\mathrm{d} i}=-\sum_{k=1}^{s} D_{i k} \operatorname{grad} c_{k}
\]

где $s$ – число компонент. Комбинируя формулы (4.3.1) и (4.3.3), получаем
\[
\frac{\partial c_{i}}{\partial t}=\operatorname{div}\left(\sum_{k=1}^{s} D_{i k} \operatorname{grad} c_{k}\right)+f_{i}\left(c_{1}, \ldots, c_{s}\right), \quad i=1,2, \ldots, s .
\]

Здесь $D_{i k}$ – коэффициенты диффузии и взаимной диффузии. Если положить в дальнейшем $D_{i i}=$ const и $D_{i j}=0$ при $i
eq j$

Рис. 4.9. Одномерная двухкомпонентная система типа «реакция – диффузия».
и рассматривать только две компоненты реакции, считая при этом задачу одномерной, то соответствующие уравнения баланса можно записать в виде (см. рис. 4.9)
\[
\begin{array}{l}
\frac{\partial c_{1}}{\partial t}=\frac{D_{1}}{L^{2}}\left(\frac{\partial^{2} c_{1}}{\partial z^{2}}\right)+f_{1}\left(c_{1}, c_{2}\right), \\
\frac{\partial c_{2}}{\partial t}=\frac{D_{2}}{L^{2}}\left(\frac{\partial^{2} c_{2}}{\partial z^{2}}\right)+f_{2}\left(c_{1}, c_{2}\right) .
\end{array}
\]

Здесь $L$ – размер системы, $z \in(0,1)$ – безразмерная координата и $j_{1}=-D_{1}\left(\partial c_{1} / \partial z\right), j_{2}=-D_{2}\left(\partial c_{2} / \partial z\right)$.

Анализ поведения систем типа «реакция-диффузия» для случая двух компонент представляет собой достаточно общую задачу. В дальнейшем мы будем использовать следующие обозначения: $x=c_{1}, y=c_{2}, D_{x}=D_{1}, D_{y}=D_{2}, f=f_{1}, g=f_{2}$; для искомых функций $x(z, t)$ и $y(z, t)$ имеем
\[
\begin{array}{l}
\frac{\partial x}{\partial t}=\frac{D_{x}}{L^{2}} \frac{\partial^{2} x}{\partial z^{2}}+f(x, y), \\
\frac{\partial y}{\partial t}=\frac{D_{y}}{L^{2}} \frac{\partial^{2} y}{\partial z^{2}}+g(x, y) .
\end{array}
\]

Выбор начальных и граничных условий для системы (4.3.7) зависит от конкретной физической ситуации. В случае задания на границе рассматриваемой пространственной области постоянных значений концентраций мы будем говорить о граничных условиях 1-го рода (ГУ1), или условиях Дирихле:
\[
\begin{array}{lll}
z=0: & x(0, t)=\bar{x}_{0}, & y(0, t)=\bar{y}_{0}, \\
z=1: & x(1, t)=\bar{x}_{1}, \quad \bar{y}(1, t)=\bar{y}_{1} .
\end{array}
\]

Важный частный случай ГУ1:
\[
\bar{x}_{0}=\bar{x}_{1}=\bar{x}, \quad \bar{y}_{0}=\bar{y}_{1}=\bar{y},
\]

где $\bar{x}$ и $\bar{y}$ представляют собой решение уравнений
\[
f(\bar{x}, \bar{y})=0, \quad g(\bar{x}, \bar{y})=0 .
\]

При таких граничных условиях есть очевидное стационарное решение уравнений (4.3.7), однородное по пространству:
\[
x(z) \equiv \bar{x}, \quad y(z) \equiv \bar{y} .
\]

Другим часто встречающимся типом граничных условий являются условия 2 -го рода (ГУ2, или условия Неймана)
\[
z=0,1: \frac{\partial x}{\partial z}=\frac{\partial y}{\partial z}=0 .
\]

Указанные условия характеризуют непроницаемость границ области для компонент $x$ и $y$. Тривиальное стационарное решение (4.3.11), очевидно, удовлетворяет и ГУ2.

Граничные условия 3-го рода (ГУЗ), которые описывают частичную проницаемость границ системы для компонент $x$ и $y$, имеют вид
\[
\begin{array}{l}
z=0: \alpha_{x 0} \frac{\partial x}{\partial z}+\beta_{x 0} x=\gamma_{x 0}, \quad \alpha_{y 0} \frac{\partial y}{\partial z}+\beta_{y 0} y=\gamma_{y 0}, \\
z=1: \alpha_{x 1} \frac{\partial x}{\partial z}+\beta_{x 1} x=\gamma_{x 1}, \quad \alpha_{y 1} \frac{\partial y}{\partial z}+\beta_{y 1} y=\gamma_{y 1} .
\end{array}
\]

Для всех трех типов граничных условий мы ввели здесь сокращения ГУ1, ГУ2 и ГУЗ, которые в дальнейшем (в данной главе и в гл. 6) будут часто использоваться для упрощения записи. Заметим, что система может иметь на своей левой и правой границах граничные условия различных типов. Так, например, если нас интересует симметричное относительно центра промежутка решение для ГУ1 (4.3.9), то мы можем рассмотреть это решение на половинном промежутке, т. е. для $z \in[0,1 / 2]$, причем в точке $z=0$ мы задаем ГУ1, а в точке $z=1 / 2$ ГУ 2 .

Начальные условия для системы (4.3.7) имеют вид
\[
t=0: x(z, 0)=x_{0}(z), \quad y(z, 0)=y_{0}(z)
\]

и описывают начальное распределение концентраций (концентрационные профили).

Уравнения (4.3.7) с приведенными выше начальными и граничными условиями, кроме изотермической системы типа «реакция – диффузия», могут описывать также, к примеру, некоторые задачи экологии [4.42].

Рассмотрим теперь асимптотическое поведение уравнений (4.3.7) с граничными условиями типа ГУ2 для случая очень больших интенсивностей массопереноса [4.43]. При этом $\partial^{2} x / \partial z^{2} \rightarrow 0$ и $\partial^{2} y / \partial z^{2} \rightarrow 0$, в результате чего вместо уравнений (4.3.7) мы получаем систему, в которой пространственные градиенты концентраций отсутствуют:
\[
\frac{d x}{d t}=f(x, y), \quad \frac{d y}{d t}=g(x, y) .
\]

Уравнения (4.3.15) описывают систему с идеальным перемешиванием (систему с сосредоточенными параметрами).

Положим теперь один из коэффициентов диффузии равным нулю. Это можно сделать в тех случаях, когда величины коэффициентов диффузии существенно различаются между собой, либо какая-нибудь компонента системы связана (неподвижна). Тогда найденное решение упрощенной задачи может служить аппроксимацией решения исходной задачи. Например, если $D_{x} \gg D_{y}$, то построенная таким образом аппроксимационная модель принимает вид
\[
\begin{array}{l}
\frac{\partial x}{\partial t}=\frac{D_{x}}{L^{2}} \frac{\partial^{2} x}{\partial z^{2}}+f(x, y), \\
\frac{\partial y}{\partial t}=g(x, y) .
\end{array}
\]

Система (4.3.7) в случае задания ГУ2 обладает еще одной интересной особенностью. Зная решение $x(z), y(z)$ на промежутке $z \in[0,1]$ при заданном $L$, мы можем с помощью «сложения профилей» построить решение $\tilde{x}(\tilde{z}), \tilde{y}(\tilde{z})$ при $\tilde{L}=m L$, где $m$-натуральное число, воспользовавшись для этого следующим способом:
\[
\tilde{z} \in\left[\frac{k}{m}, \frac{k+1}{m}\right]: \tilde{x}\left(\tilde{z}=\frac{k+z}{m}\right)=x(\varphi(z)), \quad 0 \leqslant k \leqslant m-1,
\]

где функция $\varphi(z)$ определена при $z \in[0,1]$ как
\[
\begin{array}{ll}
\varphi(z)=z, & \text { если } k \text { четное, } \\
\varphi(z)=1-z, & \text { если } k \text { нечетное. }
\end{array}
\]

Аналогичный подход – «сложение решений» – можно использовать и для периодических решений. Для того чтобы лучше уяснить себе смысл операции «сложения решений», читателю рекомендуется изобразить этот процесс графически для $m=2$.

4.3.1.1. Задача 11. Система «реакция – диффузия» для кинетики типа «брюсселятор»
Если использовать модельную кинетику типа «брюсселятор» (см. задачу 7), то функции $f$ и $g$ в уравнениях (4.3.7) принимают вид
\[
\begin{array}{l}
f(x, y)=A-(B+1) x+x^{2} y, \\
g(x, y)=B x-x^{2} y .
\end{array}
\]

Тривиальное стационарное решение (4.3.11) в этом случае записывается как
\[
x(z) \equiv A, \quad y(z) \equiv B / A .
\]

Роль параметров в этой задаче играют $D_{x}, D_{y}, L, A, B$.

4.3.1.2. Задача 12. Система «реакция-диффузия», случай SH-кинетики
Если рассматривать распределенную систему с кинетикой, описываемой SH-моделью (см. задачу 4), то функции $f$ и $g$ в уравнениях (4.3.7) представляются в виде
\[
\begin{array}{l}
f(x, y)=\alpha\left(v_{0}+x^{\gamma}\right) /\left(1+x^{\gamma}\right)-x(1+y), \\
g(x, y)=x(\beta+y)-\delta x .
\end{array}
\]

Решение системы уравнений $f=g=0$, в отличие от модели типа «брюсселятор», здесь может быть найдено только численно (при заданных значениях параметров $\alpha, \delta, v_{0}>0, \beta, \gamma>$ $>1$ ). При этом в определенном диапазоне изменения параметров можно получить несколько решений. Как и выше, задачу 12 можно рассматривать с граничными условиями всех трех типов. Отметим, что данная задача имеет восемь параметров: $D_{x}, D_{y}$, $L, \alpha, \delta, v_{0}, \beta, \gamma$.

4.3.1.3. Задача 13. Система «реакция-диффузия» в случае модели Майнхардта
Одной из наиболее известных моделей морфогенеза в настоящее время является модель системы «реакция-диффузия» с кинетикой типа активатор-ингибитор, предложенная Майнхардтом [4.44-4.47]. Эта модель описывает пространственную дифференциацию ткани в процессе эмбриогенеза при почти симметричных начальных условиях.

Модель описывается системой двух уравнений «реакциядиффузия» типа (4.3.7). Функции $f$ и $g$ в данном случае имеют івид
\[
\begin{array}{l}
f(x, y)=\rho_{0} \rho+c \rho x^{2} / y-\mu x, \\
g(x, y)=c^{\prime} \rho^{\prime} x^{2}-v y .
\end{array}
\]

Величины $\rho$ и $\rho^{\prime}$ характеризуют здесь плотность источников для веществ, действующих как активатор и ингибитор, а $\rho_{0}, c, \mu$, $c^{\prime}, v$ – положительные параметры. У этой системы существует единственное тривиальное решение $\bar{x}, \bar{y}$ :
\[
\begin{array}{l}
\bar{x}=\rho_{0} \rho / \mu+c \rho v / c^{\prime} \rho^{\prime} \mu, \\
\bar{y}=c^{\prime} \rho^{\prime} \bar{x}^{2} /
u .
\end{array}
\]

Для задачи 13 можно использовать граничные условия всех трех типов – ГУ1, ГУ2 или ГУЗ. В данной задаче имеется 10 параметров: $D_{x}, D_{y}, L, \rho, \rho^{\prime}, \rho_{0}, c, \mu, c^{\prime}, v$.

4.3.2. Задача 14. Трубчатый неизотермический реактор с аксиальным перемешиванием
Модель трубчатого реактора с аксиальным перемешиванием используется как для гомогенных, так и для гетерогенных (каталитических) реакторов ([4.48], [4.49]). В последнем случае полезна упрощенная, псевдогомогенная модель, основанная на предположении, что гетерогенную систему «катализатор-реакционная смесь» можно заменить гомогенной средой с некими эффективными характеристиками. Такая псевдогомогенная модель и формулируется в задаче 14.

Рассмотрим трубчатый реактор с теплопередачей через стенку (см. рис. 4.10). Будем считать, что реакционная смесь полностью перемешивается в радиальном направлении; тем самым мы будем рассматривать только продольные градиенты концентраций компонент и температуры. Далее, предположим, что плотность потока компонент в продольном направлении описывается соотношением вида $\sim-D_{\mathrm{e}}(\partial c / \partial l)$, аналогичным

закону Фика, где $D_{\mathrm{e}}$-эффективный коэффициент диффузии, $c$-концентрация компоненты.

Аналогично будем предполагать, что плотность потока тепла в продольном направлении задается соотношением вида – $k_{\mathrm{e}}(\partial T / \partial l)$, где $k_{\mathrm{e}}$ – эффективный коэффициент теплопроводности.

Рис. 4.10. Трубчатый реактор с аксиальным переносом тепла и массы; граничные условия типа Данквертса.
Предположим, далее, что в реакторе протекает реакция первого порядка, описываемая выражением для скорости реакции вида $k_{\infty} c \exp (-E / R T)$ ( $E$ – энергия активации, $R$ – газовая постоянная) и с тепловым эффектом (энтальпией реакции) $-\Delta H_{\mathrm{r}}$. Мы будем считать, что плотность теплового потока через стенку трубки определяется выражением $4 U\left(T-T_{c}\right) / d$, где $U$-соответствующий коэффициент теплопередачи, $d$-диаметр трубки и $T_{c}$ – температура теплообменной (например, охлаждающей) среды вне трубки.

Мы будем предполагать, что в реакторе имеется катализатор. Пусть плотность жидкости равна $\rho_{\mathrm{f}}$, а ее теплоемкость на единицу объема при постоянном давлении постоянна и равна $C_{p}$; плотность же катализатора равна $\rho_{\mathrm{s}}$, а его теплоемкость. на единицу объема также постоянна и равна $C_{p s}$. Тогда уравнения баланса массы и энергии можно представить в виде
\[
\begin{array}{c}
\varepsilon_{\mathrm{p}} \frac{\partial c}{\partial \tau}=D_{\mathrm{e}} \frac{\partial^{2} c}{\partial l^{2}}-v \frac{\partial c}{\partial l}-k_{\infty}\left(1-\varepsilon_{\mathrm{p}}\right) c \exp (-E / R T) \\
{\left[\varepsilon_{\mathrm{p}} \rho_{\mathrm{f}} C_{p \mathrm{f}}+\left(1-\varepsilon_{\mathrm{p}}\right) \rho_{\mathrm{s}} C_{p \mathrm{~s}}\right] \frac{\partial T}{\partial \tau}=k_{\mathrm{e}} \frac{\partial^{2} T}{\partial l^{2}}-\rho_{\mathrm{f}} C_{p \mathrm{f}} v \frac{\partial T}{\partial l}-} \\
-\frac{4 U}{d}\left(T-T_{c}\right)+\left(1-\varepsilon_{\mathrm{p}}\right)\left(-\Delta H_{\mathrm{r}}\right) k_{\infty} c \exp (-E / R T) .
\end{array}
\]

Здесь $\tau$-время, $\varepsilon_{\mathrm{p}}$ – доля объема, занятая жидкостью, ш $v$ – (постоянная) скорость жидкости. Начальные условия следующие:
\[
\tau=0: c(l, 0)=c^{0}(l), \quad T(l, 0)=T^{0}(l) .
\]

Пусть теперь на границе области выполняются граничные условия типа Данквертса (рис. 4.10) (см. [4.8])
\[
\begin{aligned}
\tau>0, \quad l=0:-D_{\mathrm{e}} \frac{\partial c}{\partial l} & =v\left(c_{0}-c\right), \\
-k_{\mathrm{e}} \frac{\partial T}{\partial l} & =v \rho_{\mathrm{f}} C_{p \mathrm{f}}\left(T_{0}-T\right), \\
\tau>0, \quad l=L: \frac{\partial c}{\partial l} & =\frac{\partial T}{\partial l}=0 .
\end{aligned}
\]

Здесь $c_{0}$ и $T_{0}$-концентрация и температура до входа в реактор. Введем следующие безразмерные переменные и параметры:
\[
\begin{array}{l}
z=l / L, \quad t=\tau v / L \varepsilon_{\mathrm{p}}, \quad y=1-c / c_{0}, \quad \Theta=\frac{T-T_{0}}{T_{0}} \cdot \frac{E}{R T_{0}}, \\
\mathrm{Pe}_{\mathrm{M}}=\frac{v L}{D_{e}}, \quad \mathrm{Pe}_{\mathrm{H}}=\frac{v \rho_{\mathrm{f}} C_{p \mathrm{f}} L}{k_{\mathrm{e}}}, \quad \gamma=\frac{E}{R T_{0}}, \\
B=\gamma \frac{\left(-\Delta H_{\mathrm{r}}\right) c_{0}}{\rho_{\mathrm{f}} C_{p \mathrm{f}} T_{0}}, \quad \mathrm{Da}=\frac{L\left(1-\varepsilon_{\mathrm{p}}\right)}{v} k_{\infty} \exp (-\gamma), \\
L e=\frac{\varepsilon_{\mathrm{p}} \rho_{\mathrm{f}} C_{p \mathrm{f}}+\left(1-\varepsilon_{\mathrm{p}}\right) \rho_{\mathrm{s}} C_{p \mathrm{~s}}}{\rho_{\mathrm{f}} C_{p \mathrm{f}} \varepsilon_{\mathrm{p}}}, \quad \beta=\frac{4 U L}{d 0 \rho_{\mathrm{f}} C_{p \mathrm{f}}}, \\
\Theta_{c}=\frac{T_{c}-T_{0}}{T_{0}} \cdot \frac{E}{R T_{0}} \text {. } \\
\end{array}
\]

Здесь $y$ обозначает конверсию, $\Theta$ – безразмерную температуру, $z$-безразмерную координату, $t$-безразмерное время, Ре Ре $_{\text {H }}$ – числа Пекле для массы и тепла, $\gamma$ – безразмерную энергию активации, $B$-безразмерное адиабатическое повышение температуры, Dа – число Дамкёлера, $\beta$ – безразмерный коэффициент теплопередачи, $\Theta_{c}$ – безразмерную температуру охлаждающей среды и Le – число Льюиса.
В переменных (P14-6) мы получаем систему уравнений
\[
\frac{\partial y}{\partial t}=\frac{1}{\mathrm{Pe}_{M}} \frac{\partial^{2} y}{\partial z^{2}}-\frac{\partial y}{\partial z}+\mathrm{Da}(1-y) \exp \frac{\Theta}{1+\Theta / \gamma},
\]

Le $\frac{\partial \Theta}{\partial t}=\frac{1}{\mathrm{Pe}_{\mathrm{H}}} \frac{\partial^{2} \Theta}{\partial \boldsymbol{z}^{2}}-\frac{\partial \Theta}{\partial \boldsymbol{z}}+B \mathrm{Da}(1-y) \exp \frac{\Theta}{1+\Theta / \gamma}-\beta\left(\Theta-\Theta_{c}\right)$

с граничными условиями
\[
\begin{array}{l}
z=0: \mathrm{Pe}_{\mathrm{M}} y=\frac{\partial y}{\partial z}, \quad \mathrm{Pe}_{\mathrm{H}} \Theta=\frac{\partial \Theta}{\partial z}, \\
z=1: \frac{\partial y}{\partial z}=\frac{\partial \Theta}{\partial z}=0 .
\end{array}
\]

Данная задача имеет восемь параметров: $\mathrm{Pe}_{\mathrm{M}}, \mathrm{Pe}_{\mathrm{H}}, \gamma, B, \mathrm{Da}, \beta$, $\Theta_{c}$ и Le. Заметим, что уравнения (P14-7), (Р14-8) отличаются от уравнений (4.3.7), описывающих систему типа «реакциядиффузия», конвективными слагаемыми, а именно слагаемыми $\partial y / \partial z$ и $\partial \Theta / \partial z$. При этом граничные условия (Р14-9) принадлежат к типу ГУ3, а условия (Р14-10) – к типу ГУ 2.

Если аксиальное перемешивание осуществляется в слабой степени (эффективные диффузия и теплопроводность очень малы, $\left.\mathrm{P}_{\mathrm{H}} \rightarrow \infty, \mathrm{P}_{\mathrm{M}} \rightarrow \infty\right)$, то вместо уравнений (P14-7), (P14-8) мы получаем гиперболические уравнения первого порядка, которые описывают трубчатый реактор идеального вытеснения:
\[
\begin{array}{c}
\frac{\partial y}{\partial t}+\frac{\partial y}{\partial z}=\mathrm{Da}(1-y) \exp \frac{\Theta}{1+\Theta / \gamma} \\
\operatorname{Le} \frac{\partial \Theta}{\partial t}+\frac{\partial \Theta}{\partial z}=B \mathrm{Da}(1-y) \exp \frac{\Theta}{1+\Theta / \gamma}-\beta\left(\Theta-\Theta_{c}\right) .
\end{array}
\]

При этом условия (Р14-9) и (P14-10) заменяются условиями
\[
z=0: y=0, \quad \Theta=0 .
\]

Соотношения (P14-11)-(P14-13) используются иногда для нахождения приближенного решения исходной задачи при больших значениях $\mathrm{Pe}_{\mathrm{H}}, \mathrm{Pe}_{\mathrm{M}}$.

4.3.3. Задача 15. Трубчатый неадиабатический реактор
с аксиальным перемешиванием (двухфазная модель)
В случаях, когда температура и концентрация в объеме жидкости и на поверхности частиц катализатора резко отличаются, приходится рассматривать двухфазную модель. В этой модели, кроме переменных $T$ и $c$, характеризующих жидкость, вводятся температура поверхности катализатора $T^{*}$ и концентрация на его поверхности $c^{*}$.

Здесь мы рассмотрим вариант этой модели, предложенный Лью и Амундсоном [4.50]. Помимо предположений, введенных при описании задачи 14 , будем считать, что реакция протекает на поверхности катализатора. При составлении соответствующих уравнений баланса мы будем учитывать конвективный пере-

нос, межфазовый тепло- и массообмен, теплоотдачу стенок реактора. Потоки компонент реакции и поток тепла с внешней поверхности катализатора в жидкую реакционную смесь будем описывать с помощью коэффициентов массоотдачи и теплоотдачи, считая потоки пропорциональными разности концентраций $c^{*}-c$ или разности температур $T^{*}-T$. При этом мы не учитываем перенос тепла и массы внутри частиц катализатора. Площадь внешней поверхности частиц катализатора, отнесенную к единице объема, обозначим через $a$.

Если рассматривать реакцию первого порядка, то для случая стационарного режима уравнения баланса вещества и энтальпии можно представить в форме [4.51]
\[
\begin{array}{c}
D_{\mathrm{e}} \frac{d^{2} c}{d l^{2}}-v \frac{c^{*} c}{d l}+k_{c} a\left(c^{*}-c\right)=0, \\
k_{\mathrm{e}} \frac{d^{2} T}{d l^{2}}-v C_{p} \frac{d T}{c^{\prime} l}+h a\left(T^{*}-T\right)-\frac{4 U}{d}\left(T-T_{c}\right)=0 .
\end{array}
\]

Уравнения баланса массы и энтальпии на внешней поверхности катализатора принимают вид
\[
\begin{array}{c}
k_{\mathrm{e}} a\left(c^{*}-c\right)+k_{\infty} \exp \left(-\frac{E}{R T^{*}}\right) c^{*}=0, \\
h a\left(T^{*}-T\right)-k_{\infty} \exp \left(-\frac{E}{R T^{*}}\right) c^{*}\left(-\Delta H_{\mathrm{r}}\right)=0 .
\end{array}
\]

Граничные условия на входе и выходе задаются соотношениями (P14-4) и (P14-5). Вводя безразмерные переменные по. формулам (Р14-6) и, кроме того, полагая
\[
J_{M}=\frac{k_{\mathrm{e}} a L}{v}, \quad J_{\mathrm{H}}=\frac{h a L}{v C_{p}}, \quad \omega=\frac{c_{0}-c^{*}}{c_{0}}, \quad \theta=\frac{T^{*}-T_{0}}{T_{0}} \gamma,
\]

получим:
\[
\begin{array}{c}
\frac{1}{\mathrm{Pe}_{M}} \frac{d^{2} y}{d z^{2}}-\frac{d y}{d z}+J_{M}(\omega-y)=0 \\
\frac{1}{\mathrm{Pe}_{\mathrm{H}}} \frac{d^{2} \Theta}{d z^{2}}-\frac{d \Theta}{d z}+J_{\mathrm{H}}(\theta-\Theta)-\beta\left(\Theta-\Theta_{c}\right)=0, \\
J_{M}(\omega-y)-\mathrm{Da}(1-\omega) \exp \frac{\theta}{1+\theta / \gamma}=0 \\
J_{\mathrm{H}}(\theta-\Theta)-B \mathrm{Da}(1-\omega) \exp \frac{\theta}{1+\theta / \gamma}=0 .
\end{array}
\]

Граничные условия в безразмерном виде задаются соотношениями (Р14-9) и (P14-10).

Комбинируя (Р15-9) и (Р15-8), имеем
\[
\theta=\Theta+B \frac{J_{M}}{J_{\mathrm{H}}}(\omega-y) .
\]

Воспользовавшись (P15-10), исключим из соотношения (Р15-8) величину $\theta$ :
\[
g(\omega)=J_{M}(\omega-y)-\mathrm{Da}(1-\omega) \exp \left[\frac{\Theta+B \frac{J_{M}}{J_{\mathrm{H}}}(\omega-y)}{1+\left[\Theta+B \frac{J_{M}}{J_{\mathrm{H}}}(\omega-y)\right] / \gamma}\right]=0 .
\]
(P15-11)
Таким образом, задача описывается системой двух дифференциальных уравнений (P15-6) и (P15-7) относительно неизвестных функций $y(z)$ и $\Theta(z)$ (краевая задача) и одним нелинейным (алгебраическим) уравнением (P15-11) для определения величины $\omega$. Функция $\theta(z)$ при этом выражается через $y(z)$ и $\omega$ по формуле (P15-10). В общем случае данная задача имеет 9 параметров: $\mathrm{Pe}_{\mathrm{M}}, \mathrm{P}_{\mathrm{H}}, J_{\mathrm{M}}, J_{\mathrm{H}}, B, \mathrm{Da}, \gamma, \beta$ и $\Theta_{c}$ и представляет собой одну из проблем, в которых речь идет о совместном решении системы дифференциальных и алгебраических уравнений (см. также задачу (4.3.16) для случая установившегося режима).

4.3.4. Задача 16. Неизотермическая модель внутренней диффузии в частице пористого катализатора
Большинство каталитических реакций протекает на пористых частицах катализатора. Процессы тепло- и массообмена внутри пористой частицы и на ее поверхности часто существенным образом влияют на результирующую скорость реакции. Поэтому анализ указанных процессов имеет большое практическое значение. Математические модели, описывающие взаимодействие процессов тепло- и массообмена с реакцией на пористом катализаторе, обычно рассматриваются для частиц трех геометрических форм: бесконечная пластина или одна пора, проходящая через частицу $(a=0)$; бесконечный цилиндр $(a=1)$; и, наконец, шар ( $a=2$ ). В дальнейшем параметр $a$ будет определять собой форму частицы. Введем следующие предположения [4.11]:

1. Сложные процессы переноса внутри пористой структуры можно описать с помощью постоянных эффективных коэффициентов диффузии $D_{\mathrm{e}}$ и теплопроводности $\lambda_{\mathrm{e}}$.
2. Можно рассматривать одиночную частицу катализатора с коэффициентами формы $a=0,1,2$.
3. Концентрация исследуемой компоненты остается постоянной на всей внешней поверхности частицы; мы будем обозначать эту концентрацию как $c_{\mathrm{s}}$. Аналогично, соответствующую постоянную температуру мы обозначим $T_{\mathrm{s}}$. Концентрация и температура в ядре текущей жидкости ( $c_{\mathrm{B}}$ и $T_{\mathrm{B}}$ ) также считаются постоянными. Тепло- и массообмен на внешней поверхности частицы может быть описан с помощью постоянных коэффициентов теплоотдачи $\alpha$ и массоотдачи $k_{c}$. При указанных предположениях уравнения баланса массы и энтальпии принимают вид
\[
\begin{array}{c}
\frac{\partial c}{\partial \tau}=D_{\mathrm{e}}
abla^{2} c-r(c, T), \\
C_{p} \frac{\partial T}{\partial \tau}=\lambda_{\mathrm{e}}
abla^{2} T+\left(-\Delta H_{\mathrm{r}}\right) \cdot r(c, T) .
\end{array}
\]

Здесь $C_{p}$-теплоемкость на единицу объема псевдогомогенной реагирующей среды «катализатор-реакционная смесь» и $\left(-\Delta H_{\mathrm{r}}\right)$ – теплота реакции.

Рассмотрим теперь каталитическую реакцию с соотношением для скорости реакции, задаваемым в форме
\[
r=k_{\infty} c \exp (-E / R T) .
\]

Если рассматривать введенные выше формы частиц катализатора, то уравнения баланса (P16-1) и (P16-2) можно переписать в виде
\[
\begin{array}{c}
\frac{\partial c}{\partial \tau}=D_{\mathrm{e}}\left(\frac{\partial^{2} c}{\partial x^{2}}+\frac{a}{x} \frac{\partial c}{\partial x}\right)-k_{\infty} c \exp (-E / R T), \\
C_{p} \frac{\partial T}{\partial \tau}=\lambda_{\mathrm{e}}\left(\frac{\partial^{2} T}{\partial x^{2}}+\frac{a}{x} \frac{\partial T}{\partial x}\right)+\left(-\Delta H_{\mathrm{r}}\right) k_{\infty} c \exp (-E / R T) .
\end{array}
\]

Граничные условия, описывающие перенос тепла и массы на внешней поверхности частицы, имеют вид ( $R$ может обозначать, например, радиус частицы)
\[
\begin{aligned}
x=R: \quad k_{c}\left(c_{\mathrm{B}}-c\right) & =D_{\mathrm{e}} \frac{\partial c}{\partial x} . \\
\alpha\left(T_{\mathrm{B}}-T\right) & =\lambda_{\mathrm{e}} \frac{\partial T}{\partial x},
\end{aligned}
\]

При этом соответствующие профили концентрации и температуры считаются симметричными относительно центра частицы $(x=0)$ и, следовательно,
\[
x=0: \frac{\partial c}{\partial x}=\frac{\partial T}{\partial x}=0 .
\]

Применяя уравнения (P16-3), (Р16-4) к случаю установившегося состояния, можно получить соотношение между концентрацией и температурой внутри частицы
\[
T-T_{\mathrm{s}}=\left(-\Delta H_{\mathrm{r}}\right) \frac{D_{\mathrm{e}}}{\lambda_{\mathrm{e}}}\left(c_{\mathrm{s}}-c\right),
\]

которое принято называть формулой Пратера ( $T_{\mathrm{s}}$ и $\boldsymbol{c}_{\mathrm{s}}$-соответственно температура и концентрация на поверхности частицы).
Введем безразмерные переменные
\[
\begin{array}{c}
r=x / R, \quad y=c / c_{\mathrm{B}}, \quad \Theta=\frac{E}{R T_{\mathrm{B}}^{2}}\left(T-T_{\mathrm{B}}\right), \quad \gamma=\frac{E}{R T_{\mathrm{B}}}, \\
\beta=\frac{\left(-\Delta H_{\mathrm{r}}\right) c_{\mathrm{B}} D_{\mathrm{e}}}{\lambda_{\mathrm{e}} T_{\mathrm{B}}}, \quad \Phi=R \sqrt{\frac{k_{\infty}}{D_{\mathrm{e}}}}, \quad \mathrm{Lw}=\frac{\lambda_{\mathrm{e}}}{C_{p} D_{\mathrm{e}}}, \quad(\mathrm{P} 16-9) \\
\mathrm{Sh}=\frac{k_{\mathrm{e}} R}{D_{\mathrm{e}}}, \quad \mathrm{Nu}=\frac{\alpha R}{\lambda_{\mathrm{e}}}, \quad t=\frac{\tau \lambda_{\mathrm{e}}}{C_{p} R^{2}} .
\end{array}
\]

В этих переменных уравнения баланса в безразмерной форме перепишутся так:
\[
\begin{array}{c}
\mathrm{Lw} \frac{\partial y}{\partial t}=\frac{\partial^{2} y}{\partial r^{2}}+\frac{a}{r} \frac{\partial y}{\partial r}-\Phi^{2} y \exp \left(\frac{\Theta}{1+\Theta / \gamma}\right), \\
\frac{\partial \Theta}{\partial t}=\frac{\partial^{2} \Theta}{\partial r^{2}}+\frac{a}{r} \frac{\partial \Theta}{\partial r}+\gamma \beta \Phi^{2} y \exp \left(\frac{\Theta}{1+\Theta / \gamma}\right) .
\end{array}
\]

Граничные условия при этом принимают вид
\[
\begin{array}{r}
t>0, r=1: 1=y+\frac{1}{\mathrm{Sh}} \frac{\partial y}{\partial r}, \\
0=\Theta+\frac{1}{\mathrm{Nu}} \frac{\partial \Theta}{\partial r}, \\
t>0, r=0: \frac{\partial y}{\partial r}=\frac{\partial \Theta}{\partial r}=0 .
\end{array}
\]

Если интенсивность тепло- и массообмена на внешней поверхности частицы оказывается высокой ( $\mathrm{Sh} \rightarrow \infty, \mathrm{Nu} \rightarrow \infty$, $c_{\mathrm{s}}=c_{\mathrm{B}}, T_{\mathrm{s}}=T_{\mathrm{B}}$ ), то граничные условия (P16-13) переходят в условия
\[
y(1, t)=1, \quad \Theta(1, t)=0 .
\]

Формула Пратера (Р16-8) в этом случае приобретает вид
\[
\Theta=\gamma \beta(1-y) .
\]

Задача 16 имеет в общем случае 7 параметров: $\gamma, \beta, Ф, a, L w$, $\mathrm{Sh}, \mathrm{Nu}$.

4.3.5. Задача 17. Уравнения Навье – Стокса для случая течения жидкости между двумя бесконечными соосными вращающимися дисками

Рассмотрим установившееся, вращательное и осесимметричное течение вязкой несжимаемой жидкости. Пусть область, в которой перемещается жидкость, ограничена двумя бесконечными вращающимися дисками. Соответствующие уравнения Навье-Стокса в цилиндрических координатах $(r, \Theta, z)$ можно записатьв виде [4.52]
\[
\begin{array}{c}
u \frac{\partial u}{\partial r}+w \frac{\partial u}{\partial z}-\frac{v^{2}}{r}=-\frac{1}{\rho} \frac{\partial p}{\partial r}+v\left(
abla^{2} u-\frac{u}{r^{2}}\right), \\
u \frac{\partial v}{\partial r}+w \frac{\partial v}{\partial z}+\frac{u v}{r}=v\left(
abla^{2} v-\frac{v}{r^{2}}\right), \\
u \frac{\partial w}{\partial r}+w \frac{\partial w}{\partial z}=-\frac{1}{\rho} \frac{\partial p}{\partial z}+v
abla^{2} w .
\end{array}
\]

Здесь $(u, v, w)$-составляющие скорости по координатам ( $r, \Theta$, $z$ ), $p$-давление, $\rho$-плотность и $v$ – кинематическая вязкость. Oператор $
abla^{2}$ в цилиндрических координатах имеет вид
\[

abla^{2}=\frac{\partial^{2}}{\partial r^{2}}+\frac{1}{r} \frac{\partial}{\partial r}+\frac{\partial^{2}}{\partial z^{2}} .
\]

К уравнениям (P17-1)-(P17-3) добавляется уравнение неразрывности
\[
\frac{1}{r} \frac{\partial(r u)}{\partial r}+\frac{\partial w}{\partial z}=0 .
\]

Мы будем считать, что один из дисков располагается в плоскости $z=0$ и вращается с угловой скоростью $\Omega$. Другой диск лежит в плоскости $z=d$ и вращается с угловой скоростью $S \Omega$, где $S \in[-1,1]$.

Граничные условия для уравнений (Р17-1)-(Р17-4) («условия прилипания») записываются в виде
\[
\begin{array}{l}
z=0: u=0, \quad v=\Omega r, \quad w=0, \\
z=d: u=0, \quad v=S \Omega r, \quad w=0 . \\
\end{array}
\]

Введем предположения, приводящие задачу (Р17-1)(Р17-5) к краевой задаче для системы обыкновенных дифференциальных уравнений.

Будем искать решение, у которого осевая составляющая скорости зависит только от осевой координаты, т. е.
\[
w \equiv w(z)
\]

Далее предположим, что при $r=0$ составляющая скорости $u$ имеет конечное значение. Тогда, интегрируя уравнение (P17-4) по переменной $r$, мы получаем ( $A(z)$ – постоянная интегрирования)
\[
r u+A(z)=-\frac{r^{2}}{2} w^{\prime}(z) .
\]

Подставив $r=0$, получаем $A(z) \equiv 0$ и, следовательно,
\[
u=-\frac{r}{2} w^{\prime}(z) \text {. }
\]

Из уравнения (Р17-3) имеем
\[
w \cdot w^{\prime}(z)=-\frac{1}{\rho} \frac{\partial p}{\partial z}+v w^{\prime \prime}(z),
\]

после чего, интегрируя по $z$, находим ( $B(r)$ – постоянная интегрирования)
\[
\frac{\omega^{2}}{2}=-\frac{p}{\rho}+v w^{\prime}(z)+B(r) .
\]

Дифференцирование соотношения (P17-8) по $r$ дает
\[
\frac{1}{\rho} \frac{\partial p}{\partial r}=\frac{d B(r)}{d r} .
\]

Из уравнения (Р17-1) после подстановки формул (Р17-7) и (Р17-9) получаем
\[
\frac{d B(r)}{d r}-\frac{v^{2}}{r}=r \varphi(z) .
\]

Интегрирование (P17-10) с учетом граничных условий (P17-5) дает
\[
B(r)=\frac{1}{2} r^{2} \Omega^{2}+\frac{1}{2} r^{2} q+c,
\]

где $q$ и $c$ – некоторые постоянные величины. Сравнивая теперь формулы (Р17-9) и (Р17-10), имеем
\[
v=\operatorname{rg}(z)
\]

а комбинируя формулы (Р17-8) и (Р17-11), находим
T. e.
\[
\frac{p}{\rho}=\frac{w^{2}}{z}-v w^{\prime}(z)+\frac{r^{2}}{2}\left(\Omega^{2}+q\right)+c,
\]
\[
\frac{p}{\rho}=P(z)+\frac{r^{2}}{z}\left(\Omega^{2}+q\right) .
\]

Соотношения (P17-7), (Р17-12) и (P17-14) определяют вид соответствующих функций для составляющих скорости и для давления. Принимая во внимание эти соотношения, воспользуемся подстановкой
\[
\begin{array}{c}
\xi=z / d \\
u=r \Omega F(\xi), \\
v=r \Omega G(\xi), \\
w=(v \Omega)^{1 / 2} H(\xi), \\
\frac{p}{\rho}=v \Omega P(\xi)+\frac{1}{2} k \Omega^{2} r^{2} .
\end{array}
\]
:Здесь $k$ – не определенная пока постоянная.
Обозначим $\mathrm{Re}=\Omega d^{2} / v$ и преобразуем уравнения (P17-1), (P17-2) и (P17-4) с помощью формул (P17-15a), .., (P17-15e) ж виду
\[
\begin{array}{l}
F^{\prime \prime}=\sqrt{\operatorname{Re}} H F^{\prime}+\operatorname{Re}\left(F^{2}-G^{2}+k\right), \\
G^{\prime \prime}=2 \operatorname{Re} F G+\sqrt{\operatorname{Re}} H G^{\prime}, \\
H^{\prime}=-2 \sqrt{\operatorname{Re}} F .
\end{array}
\]

Функции $F, G$ и $H$ являются безразмерными. Функцию $p$, выраженную с помощью функций $F$ и $H$, мы можем найти из уравнения (P17-3). Используя формулы (Р17-15), граничные условия (Р17-5) представим в безразмерном виде
\[
\begin{array}{l}
\xi=0: F=H=0, \quad G=1, \\
\xi=1: F=H=0, \quad G=S . \\
\end{array}
\]
Если значения параметров $\operatorname{Re}$ и $S$ заданы, то система уравнений (P17-16)-(P17-18) (с соответствующими граничными условиями) представляет собой краевую задачу пятого порядка с одной неизвестной постоянной $k$. Соотношения (Р17-19), (P17-20) определяют шесть граничных условий. Тем самым задача поставлена полностью.