Hybrid numerical-analytical method for solving problems of salt ion transport in membrane systems with axial symmetry

Full Text

Введение

Гидродинамические системы с осевой симметрией в настоящее время широко используются и интенсивно исследуются [1–8]. Так, например, в работе [9] найдено точное решение в рамках уравнений Эйлера закрученных осесимметричных стационарных течений идеальной несжимаемой жидкости, а в работе [10] в рамках уравнений Навье–Стокса рассмотрены нестационарные осесимметричные течения однородной вязкой несжимаемой жидкости, в которых осевая и окружная скорости зависят только от радиуса и от времени, а радиальная скорость равна нулю. Произведено расщепление краевых задач для рассматриваемого типа течений на две задачи, каждая из которых содержит две неизвестные функции (давление и одна из компонент скорости). Актуальным является исследование аналогичных задач, но с дополнительными процессами переноса вещества за счет конвекции, диффузии, электромиграции и т.д. Изучению свойств переноса ионов соли в системах с вращающейся жидкостью, в том числе методом вращающегося мембранного диска, посвящены ряд работ [1–8]. Эта работа является продолжением работы [3], где приведена базовая модель, основанная только на основных законах сохранения и естественных граничных условиях и не использующая никаких гипотез для упрощения. Базовая модель служит для проверки исследования основных свойств и вывода упрощенных моделей на основе этих свойств. Достоинствами базовой модели являются общность и отсутствие подгоночных параметров, к недостаткам можно отнести сложность расчета при высоких скоростях вращения мембранного диска и скачке потенциала, а также реальных концентрациях раствора. Это связано с большими градиентами концентрации и потенциала вблизи катионнообменной мембраны (КОМ). Для преодоления указанных ограничений предлагается разбиение области решения на две подобласти: первая — основная часть области (ОЧО) вдали от катионообменной мембраны; вторая — область возрастания катионов (ОВК) вблизи КОМ, границей которой служит точка локального минимума концентрации катионов около КОМ. Для каждой из этих областей из базовой модели выводятся упрощенные модельные задачи, удобные как для численного решения (в ОЧО), так и для аналитического (в ОВК). Аналитическое решение зависит от того, больше или меньше ток, протекающий через систему, по сравнению с предельным током [11]. В допредельном режиме ОВК совпадает с квазиравновесной областью пространственного заряда (ОПЗ), однако в сверхпредельном режиме ОВК должна содержать наряду с квазиравновесной ОПЗ еще и промежуточную область (слой), необходимую для сращивания решений. Решения в этих подобластях сращиваются, что позволяет искать решение при большой угловой скорости, начальной концентрации и скачке потенциала. Модель без области возрастания катионов назовем упрощенной моделью без ОВК. Аналогичная модель переноса ионов соли для плоского канала впервые была предложена, по-видимому, в работе [12], однако в этой работе не было исследовано решение в ОВК и поэтому не было произведено сращивание.

1. Формулировка базовой модели и модели без ОВК в задачах с осевой симметрией

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

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

1.1. Область исследования

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

Рис. 1. Исследуемое сечение области $\Omega$ и ее границы: 1 — глубина раствора, где выполняется условие электронейтральности; 2 — ось симметрии; 3 — катионнообменная мембрана; 4 — открытая граница
[Figure 1. Investigated cross-section of the region $\Omega$ and its boundaries: 1—depth of the solution where the electroneutrality condition is satisfied; 2—axis of symmetry; 3—cation-exchange membrane; 4—open boundary]

1.2. Система уравнений

В области $\Omega$ для описания переноса ионов соли используются связанная система уравнений Навье–Стокса с объемной электрической силой $\vec{f}$ и Нернста–Планка–Пуассона, которые в цилиндрической системе координат с учетом осевой симметрии принимают следующий вид:
\[ \begin{equation*}
\begin{array}{c}
\rho \frac{\partial u}{\partial t}+\rho \Bigl( u\frac{\partial u}{\partial r}-\frac{v^2}{r}+w\frac{\partial u}{\partial z} \Bigr)+\frac{\partial p}{\partial r}=
\eta \Bigl[ \frac{1}{r}\frac{\partial }{\partial r}\Bigl( r\frac{\partial u}{\partial r} \Bigr)-\frac{u}{r^2}+\frac{{\partial^2}u}{\partial z^2} \Bigr]+f_r,
\\
\rho \frac{\partial v}{\partial t}+\rho \Bigl( u\frac{\partial v}{\partial r}-\frac{uv}{r}+w\frac{\partial v}{\partial z} \Bigr)=
\eta \Bigl[ \frac{1}{r}\frac{\partial }{\partial r}\Bigl( r\frac{\partial v}{\partial r} \Bigr)-\frac{v}{r^2}+\frac{{\partial^2}v}{\partial z^2} \Bigr]+f_\varphi ,
\\
\rho \frac{\partial w}{\partial t}+\rho \Bigl( u\frac{\partial w}{\partial r}+w\frac{\partial w}{\partial z} \Bigr)+\frac{\partial p}{\partial z}=\eta \Bigl[ \frac{1}{r}\frac{\partial }{\partial r}\Bigl( r\frac{\partial w}{\partial r} \Bigr)+\frac{\partial^2 w}{\partial z^2} \Bigr]+ f_z;
\end{array}
\end{equation*} \]
\[ \begin{multline}\label{eq:eq1}
\vec{N}_i= \Bigl( \frac{F}{RT}{z_i}{D_i}{C_i}{E_r}-{D_i}\frac{\partial {C_i}}{\partial r}+{C_i}u \Bigr)
\vec{e}_r+ {C_i} v \vec{e} _{\varphi }+ {}
\\
{}+\Bigl( \frac{F}{RT}{z_i}{D_i}{C_i}{E_z}-{D_i}\frac{\partial {C_i}}{\partial z}+{C_i}w \Bigr) \vec{e} _z,
\quad i=1, 2;
\end{multline} \tag{1} \]
\[ \begin{equation}
\frac{\partial {C_i}}{\partial t}=-\frac{1}{r}\frac{\partial }{\partial r} \bigl( r N_{i, r}\bigr)-\frac{\partial N_{i, z} }{\partial z},\quad i=1, 2;
\end{equation} \tag{2} \]
\[ \begin{equation}
\frac{\partial^2 \Phi }{\partial r^2}+\frac{\partial^2 \Phi }{\partial z^2}=
-\frac{1}{r}\frac{\partial \Phi }{\partial r}-\frac{F}{\varepsilon }\bigl( z_ 1 C_1 + z_ 2 C_2\bigr),
\end{equation} \tag{3} \]
\[ \begin{equation}
\vec{I}= F( z_1 \vec{N}_1 +z_2 \vec{N}_{2}),
\end{equation} \tag{4} \]
где
\[ \begin{equation*}
\begin{array}{l}
f_r =\varepsilon \Bigl( \frac{1}{r}\frac{\partial }{\partial r}\Bigl( r\frac{\partial \Phi }{\partial r} \Bigr)+
\frac{\partial^2 \Phi }{\partial z^2} \Bigr) \frac{\partial \Phi }{\partial r}=
\varepsilon \Bigl( \frac{\partial^2\Phi }{\partial r^2}+\frac{\partial^2 \Phi }{\partial z^2} \Bigr)
\frac{\partial \Phi }{\partial r}+\varepsilon \frac{1}{r}
\Bigl( \frac{\partial \Phi }{\partial r} \Bigr)^{2}, \\
f_\varphi =0, \\
f_z =\varepsilon \Bigl( \frac{1}{r}\frac{\partial }{\partial r}\Bigl( r\frac{\partial \Phi }{\partial r} \Bigr)+\frac{{\partial^2}\Phi }{\partial z^2} \Bigr)\frac{\partial \Phi }{\partial z}=
\varepsilon \Bigl( \frac{{\partial^2}\Phi }{\partial r^2}+\frac{{\partial^2}\Phi }{\partial z^2} \Bigr)\frac{\partial \Phi }{\partial z}+\varepsilon \frac{1}{r}\frac{\partial \Phi }{\partial r} \frac{\partial \Phi }{\partial z};
\end{array}
\end{equation*} \]
$u$, $v$, $w$ — радиальная, азимутальная и аксиальная компоненты скорости течения раствора; $\vec{N}_1$, $\vec{N}_2$, $\vec{C}_1$, $\vec{C}_2$ — потоки и концентрации катионов и анионов в растворе соответственно; $z_1$, $z_2$ — зарядовые числа катионов и анионов; $\vec{I}$ — плотность тока; $D_{1}$, $D_{2}$ — коэффициенты диффузии катионов и анионов соответственно; $\Phi$ — потенциал электрического поля; $\vec{E}=-\nabla \Phi$ — напряженность электрического поля; $\varepsilon$ — диэлектрическая проницаемость электролита; $F$ — постоянная Фарадея; $R$ — газовая постоянная; $T$ — абсолютная температура; $t$ — время; $\rho$ — плотность; $\eta$ — динамическая вязкость; $p$ — давление.

1.3. Краевые условия для базовой модели

Граница 1 является входом для раствора; для концентраций выполняется условие электронейтральности; задается постоянная начальная концентрация $C_{0}$; для скорости ставится условие отсутствия нормального напряжения $\bigl(\nabla \vec{u}+ (\nabla \vec{u})^{\top}\bigr)\vec{n}=0$; задается нулевое давление; граница 1 считается также эквипотенциальной (анодом) — $\Phi (t,r,0)= d_{\phi }$.

Граница 2 является осью симметрии.

Граница 3 соответствует вращающейся идеально селективной КОМ и является выходом для катионов, концентрация которых постоянна и определяется емкостью мембраны — $C_1 (t,r,H)=C_{km}$. КОМ предполагается идеально селективной и поэтому для анионов задается условие отсутствия потока — $-\vec{n}\cdot \vec{N}_ 2=0$; для потенциала $\Phi (t,r,H)=0$ задается азимутальная скорость в виде $v(t,r,H)=\omega r$, а остальные компоненты скорости принимаются равными нулю.

Граница 4 является выходом для раствора; ставится условие выноса ионов только конвективным потоком — $\vec{N}_{i}=-u\cdot {C_i}$, $i=1, 2$; для потенциала задается условие $-\vec{n}\cdot \bigl(r\frac{\partial \Phi }{\partial r}, \frac{\partial \Phi }{\partial z}\bigr)^{\top}=0$; для скорости предполагается условие отсутствия нормального напряжения — $\bigl(\nabla \vec{u}+(\nabla \vec{u})^{\top}\bigr)\vec{n}=0$.

В начальный момент времени концентрация ионов соли предполагается постоянной — $C_1 (0,r,z)=C_2 (0,r,z) = C_0$, а раствор — неподвижным.

1.4. Краевые условия для модели без ОВК

Все начальные и граничные условия остаются такими же, как в базовой модели, за исключением условия на границе ОВК, для которой задается концентрации катионов — $-\vec{n}\nabla C_1 \bigr|_{z={z_m}}=0$, где $z_{m}$ — координата минимума концентрации, зависящая от угловой скорости вращения диска и начальной концентрации. Поскольку толщина области возрастания катионов (ОВК) мала, т.е. $z_m \approx H$, это условие может быть заменено на условие $-\vec{n}\nabla C_1 \bigr|_{z=H}=0$. Таким образом, модель без ОВК будет определена во всей ячейке. Ниже показано, что ее решение близко к решению базовой модели в ОЧО, но, естественно, их решения отличаются в ОВК.

1.5. Сравнительный анализ численных результатов исследования базовой модели и модели без ОВК

Были проведены исследования для разнообразных исходных данных (концентрация, угловая скорость и т.д.), соответствующих сверхпредельному режиму, когда возникают электроконвективные вихри. В качестве примера ниже приведены результаты вычислений при следующих условиях:¹ $C_0 =0.01$ моль/м$^3$, $\omega = \pi /2$, $d _{\phi }=1.5$ В, радиус вращающегося мембранного диска $r=1$ мм, высота ячейки — 2 мм.

На рис. 2 приведены линии тока жидкости в базовой модели и для модели без ОВК, которые показывают образование электроконвективного вихря, причем линии тока раствора практически полностью совпадают.

Рис. 2. Линии тока жидкости при $t=1000$ c в базовой модели (a), в базовой модели около катионнообменной мембраны (b), в модели без области возрастания катионов (c), в модели без области возрастания катионов около катионнообменной мембраны (d)
[Figure 2. Fluid streamlines at $t = 1000$ s in the base model (a), in the base model near the cation-exchange membrane (b), in the model without the cation increase region (c), in the model without the cation increase region near the cation-exchange membrane (d)]

Концентрации, потенциал и т.д. зависят от $t$, $r$ и $z$. Для того чтобы сделать обозримой их зависимость от $z$, рассмотрим графики при фиксированных $r$ начиная со значения $r=0.1$ мм и до значения $r=0.9$ мм с шагом 0.1 мм. На всех приведенных ниже рисунках по оси абсцисс показаны координаты высоты ячейки \(z\), причем \(z=0\) соответствует глубине раствора, где выполняется условие электронейтральности, а \(z=2\) мм — условной границе раствор/КОМ. По оси ординат показаны требуемые величины (концентрации, потенциал, компоненты скорости и т.д.).

На рис. 3, a показаны профили концентрации катионов $C_{1}$ в базовой модели, по которым можно сделать вывод, что концентрация сначала почти постоянна, а затем линейно убывает, что соответствует поведению в области электронейтральности, далее образуется расширенная ОПЗ и концентрация достигает своего минимума, а потом нелинейно возрастает в квазиравновесной ОПЗ (рис. 3, b).

На рис. 3 показаны графики концентрации катионов в базовой модели и модели без ОВК, из которых можно сделать вывод, что они совпадают с большой точностью в области электронейтральности и расширенной ОПЗ, однако около КОМ отличаются. В базовой модели концентрации достигают минимума в некоторой точке $z_{m}$, расположенной вблизи КОМ, и далее увеличиваются до граничного значения $C_{1,m}$, в то время как в модели без ОВК концентрации продолжают монотонно убывать, причем в области $(z_{m}, 1]$ практически постоянны, что объясняется разницей в граничном условии для катионов около КОМ.

Графики концентрации анионов $C_{2}$, показанные на рис. 4 для базовой модели и в модели без ОВК, ведут себя примерно одинаково, включая область около КОМ.

Рис. 3. Концентрация $C_{1}$ при $t=1000$ c и разных сечениях по $r$ в базовой модели (a), в базовой модели вблизи катионнообменной мембраны (b), в модели без области возрастания катионов (c), в модели без области возрастания катионов вблизи катионнообменной мембраны (d)
[Figure 3. Concentration $C_{1}$ at $t = 1000$ s and different cross sections for $r$ in the base model (a), in the base model near the cation-exchange membrane (b), in the model without the cation increase region (c), in the model without the cation increase region near the cation-exchange membrane (d)]

Рис. 4. Концентрация $C_{2}$ при $t=1000$ c и разных сечениях по $r$ в базовой модели (a), в базовой модели вблизи катионнообменной мембраны (b), в модели без области возрастания катионов (c), в модели без области возрастания катионов вблизи катионнообменной мембраны (d)
[Figure 4. Concentration $C_{2}$ at $t = 1000$ s and different cross sections for $r$ in the base model (a), in the base model near the cation-exchange membrane (b), in the model without the cation increase region (c), in the model without the cation increase region near the cation-exchange membrane (d)]

Рис. 5. Разность концентраций $C_1-C_2$ при $t=1000$ c и разных сечениях по $r$ в базовой модели (a), в базовой модели без учета квазиравновесной области пространственного заряда (b), в модели без области возрастания катионов (c)
[Figure 5. Concentration difference $C_1-C_2$ at $t = 1000$ s and different cross sections for $r$ in the base model (a), in the base model without taking into account the quasi-equilibrium region of the space charge (b), in the model without the cation increase region (c)]

Для анализа распределения пространственного заряда рассмотрим разность концентраций $C_{1}-C_{2}$, которая представляет собой плотность распределения пространственного заряда, нормированную на число Фарадея $F$. Из анализа графиков разности концентраций $C_{1}-C_{2}$ (рис. 5, b, c) в базовой модели и модели без ОВК следует их совпадение с высокой точностью в области электронейтральности и расширенной ОПЗ, однако, как показано на рис. 5, a, c базовой модели дополнительно формируется квазиравновесная ОПЗ.

На рис. 6 показано визуальное совпадение с большой точностью графиков потенциала базовой модели и модели без ОВК при небольших $C_{k,m}$ — порядка $C_{0}$. При этом график разницы потенциалов этих двух моделей (рис. 7) показывает небольшие отличия (максимум которых достигается при $z=0.1$ мм и в сечениях от $r=0.3$ мм до $0.5$ мм). Это связано с тем, что при данных условиях скачок потенциала в ОВК мал.

Рис. 6. Графики потенциала при $t=1000$ c и разных сечениях по $r$ в базовой модели (a), в модели без области возрастания катионов (b)
[Figure 6. Plots of the potential at $t = 1000$ s and different cross sections for $r$ in the base model (a), in the model without the cation increase region (b)]

Рис. 7. Графики разности потенциалов между базовой моделью и моделью без области возрастания катионов при $t=1000$ c и разных сечениях по $r$
[Figure 7. Plots of the potential difference between the base model and the model without the cation increase region at $t = 1000$ s and different cross sections for $r$]

Если $C_{k,m}$ значительно больше $C_{0}$, например, в тысячу и более раз, что вполне реально, то тогда разница будет более существенной и она компенсируется скачком потенциала в ОВК, рассчитанным аналитически. 

По графикам радиальной скорости (см. рис. 8) базовой модели и модели без ОВК можно сделать вывод об их высокой степени совпадения.

На рис. 9 и рис. 10 приведены графики азимутальной и аксиальной компонент скорости для базовой модели и модели без ОВК соответственно, которые показывают их практически полное совпадение.

Рис. 8. Графики радиальной скорости при $t=1000$ c и разных сечениях по $r$ в базовой модели (a), в модели без области возрастания катионов (b)
[Figure 8. Plots of the radial velocity at $t = 1000$ s and different cross sections for $r$ in the base model (a), in the model without the cation increase region (b)]

Рис. 9. Графики азимутальной скорости при $t=1000$ c и разных сечениях по $r$ в базовой модели (a), в модели без области возрастания катионов (b)
[Figure 9. Plots of the azimuthal velocity at $t = 1000$ s and different cross sections for $r$ in the base model (a), in the model without the cation increase region (b)]

Рис. 10. Графики аксиальной скорости при $t=1000$ c и разных сечениях по $r$ в базовой модели (a), в модели без области возрастания катионов (b)
[Figure 10. Plots of the axial velocity at $t = 1000$ s and different cross sections for $r$ in the base model (a), in the model without the cation increase region (b)]

Рис. 11. Графики аксиальной скорости при $t=1000$ c и разных сечениях по $r$ в базовой модели (a), в модели без области возрастания катионов (b)
[Figure 11. Plots of the axial velocity at $t = 1000$ s and different cross sections for $r$ in the base model (a), in the model without the cation increase region (b)]

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

2. Аналитическое решение краевой задачи в ОВК

Аналитическое решение в ОВК зависит от величины скачка потенциала, а именно допредельного или сверхпредельного. В допредельном режиме ОВК совпадает с квазиравновесной ОПЗ, а в сверхпредельном режиме эти области очень близки. Из численного исследования базовой модели следует, что толщина ОВК не зависит от радиальной координаты \(r\), радиальная (\(u\)) и аксиальная (\(w\)) компоненты скорости близки к нулю, а азимутальная (\(v\)) компонента скорости почти постоянна в ОВК. Поэтому в уравнении (1) останется только третья компонента:
\[ \begin{equation*}
N_{i , z}=\frac{F}{RT}{z_i}{D_i}{C_i}{E}_{z}-{D_i}\frac{\partial C_i}{\partial z},\quad i=1, 2.
\end{equation*} \]

Кроме того, толщина ОВК также практически не зависит и от времени $t$, т.е. наблюдается квазистационарный режим. В результате левая часть уравнения (2) равна нулю, а в правой части в силу независимости от $r$ равна нулю производная потока по $z$, кроме того, в уравнениях (3) и (4) также остаются составляющие, зависящие только от $z$. Таким образом, для нахождения неизвестных функций в ОВК получается следующая система уравнений:
\[ \begin{equation}
N _{i , z}=\frac{F}{RT}{z_i}{D_i}{C_i}{E}_{z}-{D_i}\frac{\partial C_i}{\partial z},\quad i=1, 2;
\end{equation} \tag{5} \]
\[ \begin{equation}
\frac{\partial N_{i , z}}{\partial z} =0,\quad i=1, 2;
\end{equation} \tag{6} \]
\[ \begin{equation}
\frac{\partial^2 \Phi }{\partial z^2} =-\frac{F}{\varepsilon } ( z_1 C_1 + z_ 2 C_ 2 );
\end{equation} \tag{7} \]
\[ \begin{equation}
I_z =F( z_1 N_{1, z}+ z_2 N_{2 , z}).
\end{equation} \tag{8} \]

Для асимптотического решения системы (5)–(8) выполним переход в области ОВК ($z_m\leqslant z\leqslant H$) к безразмерному виду с использованием характерных величин: $z^{( u )}={z}/{H}$; $i^{( u)}= C_i/C_0$; $\vec{N}_i ^{( u )}=\vec{N}_{i}/N_0$; $D_i^{( u )}= D_i/D_0$; $I_z ^{(u)}= I_z/I_0$; $\Phi^{( u )}=\Phi/\Phi_0$; $E_z^{(u)}= E_z/E_0$; $\varepsilon ^{( u)}= \varepsilon / b_0$. Здесь $C_0$ — начальная концентрация раствора; $H$ — высота ячейки; $D_0$ — коэффициент диффузии электролита; $\Phi _0 = R T_0/F$ — тепловой потенциал; $I_0 = F D_0 C_0 / H$ — характерный ток; $N_0 = D_0 C_0 / H$ — характерный поток ионов; $E_0 = \Phi _0 / H$ — отношение теплового потенциала к характерной ширине; $b_0 = H^2 F C_0 / \Phi_0$ — характерная величина, имеющая размерность электрической постоянной.

Параметр $\varepsilon ^{(u)}$ является малым параметром, имеющим смысл отношения квадрата толщины области пространственного заряда к квадрату высоты ячейки. Другими словами, это отношение (квадратов) толщин заряженной и незаряженной областей раствора.

Параметр $\varepsilon ^{(u)}$ может быть записан в виде $\varepsilon ^{(u)} \approx 0.58\cdot 10^{-12} / C_0 $.

При начальной концентрации $C_0=0.01$ моль/м$^3$, использованной для численного расчета, получаем $\varepsilon ^{(u)}\approx 5.8\cdot 10^{-11}$. Поэтому в дальнейшем $\varepsilon ^{(u)}$ считается малым параметром.

Система уравнений (5)–(8) в безразмерной форме, например, для раствора NaCl, с учетом идеальной селективности КОМ принимает следующий вид (индекс «$u$» для простоты записи опущен):
\[ \begin{equation}
\frac{d C_1}{dz} =C_1 E_z -I_z ,
\quad
\frac{d C_2}{dz} =- C_2 E_z ,
\quad
\varepsilon \frac{d E_z}{dz} = C_1 - C_2,
\end{equation} \tag{9} \]
где $E ( z, \varepsilon )=- d \Phi / dz $ — напряженность электрического поля, а $\Phi( z, \varepsilon )$ — его потенциал; $C_1 ( z,\varepsilon )$, $ C_ 2 ( z,\varepsilon )$, $I_z$ — искомые концентрации катионов и анионов соответственно; $I_z$ — ток, равный потоку катионов; $\varepsilon >0$ — малый параметр.

В качестве краевых условий в безразмерном виде при $z=z_m/H$ ставятся условия сращивания с упрощенной моделью, а при $z=1$ — условия
\[ \begin{equation}
C_1\bigr|_{z=1} = C_{1,M},
\quad
\Bigl[ - C_2 \frac{\partial \Phi}{\partial z} + \frac{\partial C_2}{\partial z} \Bigr]_{z=1}=0,
\quad
\Phi\bigr|_{z=1}=0.
\end{equation} \tag{10} \]

2.1. Решение в допредельном случае

Система уравнений (9) представляет собой сингулярно-возмущенную задачу. Проведя замену
$<br />\xi = ({z-1})/\sqrt{\varepsilon }, \; \;<br />E_z = \widetilde{E} ( \xi ,\varepsilon )/ \sqrt{\varepsilon }, \; \;<br />C_1 ( z,\varepsilon )= \widetilde{C}_1 ( \xi ,\varepsilon ), \; \;<br />C_2 ( z,\varepsilon )= \widetilde{C}_2 ( \xi ,\varepsilon )<br />$
можно убедиться, что в допредельном режиме ОВК совпадает с квазиравновесной ОПЗ.

Действительно, эти замены при малых $\varepsilon $ в начальном приближении сводят систему уравнений к классической системе уравнений Больцмана–Дебая, которая и описывает квазиравновесную ОПЗ [11]
\[ \begin{equation}
\frac{d \widetilde{C}_1}{d\xi }=\widetilde{C}_{1} \widetilde{E},
\quad
\frac{d \widetilde{C}_2}{d\xi } =-\widetilde{C}_{2} \widetilde{E},
\quad
\frac{d \widetilde{E}}{d\xi }= \widetilde{C}_1 -\widetilde{C}_{2}
\end{equation} \tag{11} \]
с соответствующими краевыми условиями, которая имеет аналитическое решение:
\[ \begin{equation}
E ( z, \varepsilon )=\frac{ \sqrt{-\alpha } }{\sqrt{\varepsilon }}
\frac{4\sqrt{\beta } \exp\bigl(\sqrt{-\alpha }({z-1})/{\sqrt{\varepsilon }}\bigr)}
{1-\beta \exp\bigl(\sqrt{-4\alpha }({z-1})/{\sqrt{\varepsilon }}\bigr) },
\end{equation} \tag{12} \]
\[ \begin{equation}
C_1( z, \varepsilon )=\frac{\varepsilon}{2} \frac{dE}{dz}+\frac{\varepsilon E^2}{4}-\frac{\alpha }{2},
\end{equation} \tag{13} \]
\[ \begin{equation}
C_2 ( z, \varepsilon )=-\frac{\varepsilon }{2}\frac{dE}{dz}+\frac{\varepsilon E^2 }{4}-\frac{\alpha}{2},
\end{equation} \tag{14} \]
где
$<br />\alpha = - \bigl( C_1 ( z_m , \varepsilon )+ C_2 (z_m , \varepsilon ) \bigr)\approx<br />-\bigl( C_1 ( 1, \varepsilon )+ C_2 ( 1, \varepsilon ) \bigr)<0,<br />$
$\beta $ — некоторое положительное число, которое определяется из условия
$<br />\widetilde{E} ( 0, \varepsilon )=\sqrt{2 (C_{1,m}+\alpha )},<br />$
выведенного из граничных условий (10) с использованием первого интеграла системы уравнений (11) и эквивалентного условию $C_1 ( 1, \varepsilon )=C_{1,m}$.

При допредельном режиме нет расширенной ОПЗ, поэтому $E ( z_m -0, \varepsilon )$ — значение численного решения (решение в ОЧО в точке $z_m $) должно быть конечным и, соответственно, $E( z_m +0, \varepsilon )$ (значение аналитического решения в ОВК) тоже должно быть конечным. Это условие выполняется, если $z_m =1-k\sqrt{\varepsilon } | \ln \varepsilon |$, тогда получаем $E(z_m +0, \varepsilon ) = 4\sqrt{\beta }\sqrt{-\alpha }$, при $\varepsilon \to 0$, если взять $k= {1}/({2\sqrt{-\alpha }})$. Равенства $E(z_{m}-0,\varepsilon )=E ( z_{m}+0, \varepsilon )$ можно добиться с использованием уже следующего первого приближения.

Сравнение точек $z_m $, полученных при допредельном режиме в результате численного решения базовой модели, где $z_m$ вычисляется как точка, в которой достигается минимум концентрации, и аналитического решения в ОВК показывает их совпадение с высокой точностью: в базовой модели $z_m \approx 1.976$ мм, а в модели ОВК $z_m \approx 1.999$ мм.

Концентрации сращиваются в силу выбора $\alpha$ и формул (13), (14).

2.2. Решение в сверхпредельном случае

В сверхпредельном режиме наряду с квазиравновесным погранслоем возникает и расширенная ОПЗ [11], где концентрации ионов малы (рис. 11) (имеют порядок $O ( \sqrt{\varepsilon } )$) и поэтому процедура сращивания в допредельном режиме не может быть использована и должна быть модифицирована. Возьмем сечение $r=0.5$ мм и построим графики концентрации катионов и анионов, а также функцию $C_ 1-C_2$ (см. рис. 11). Из этого рисунка видно, что в сверхпредельном состоянии ОВК состоит из двух частей: квазиравновесной ОПЗ и промежуточной области, которая служит для сращивания решения в квазиравновесной ОПЗ и в расширенной ОПЗ, являющейся частью ОЧО. Причем на правой границе промежуточной области сращивается напряженность электрического поля, а на левой границе — концентрация катионов. Точка $z_{m}$, где производная концентрации равна нулю, близка к 1. При допредельном режиме $z_m$ и $z_k$ совпадают, при сверхпредельном режиме $z_m <z_k$ (см. рис. 11).

При асимптотическом решении краевой задачи для системы (9) в области $( z_m , 1 ]$ в сверхпредельном режиме примем естественное для этого режима предположение об отсутствии в этой области анионов, то есть $C_ 2 ( z, \varepsilon )\equiv 0$, $z\in ( z_m , 1 ]$.

Тогда получаем систему уравнений
\[ \begin{equation*}
\frac{d C_1}{dz}= C_ 1 E_z - I_z ,
\quad
\varepsilon \frac{d E_z}{dz} =C_1,
\end{equation*} \]
которая имеет аналитическое решение в общем виде:
\[ \begin{equation}
E (z, \varepsilon)=\frac{1}{\sqrt{\varepsilon }}
\frac{4b\sqrt{\beta } \exp\bigl( b (z-1)/\sqrt{\varepsilon } \bigr)}
{1-\beta \exp\bigl( 2b(z-1)/\sqrt{\varepsilon } \bigr)},
\end{equation} \tag{15} \]
\[ \begin{equation}
C_1 ( z, \varepsilon )=
\frac{4b^2 \sqrt{\beta } \exp\bigl( b ( z-1)/ \sqrt{\varepsilon } \bigr)
\bigl( 1+\beta \exp\bigl( 2b (z-1)/ \sqrt{\varepsilon } \bigr) \bigr)}
{ \bigl( 1-\beta\exp \bigl( 2b(z-1)/ \sqrt{\varepsilon } \bigr) \bigr)^2 },
\end{equation} \tag{16} \]
где $b>0$.

2.3. Сращивание решений и определение констант

Возьмем
\[ \begin{equation*}
z_k =1- k_1 \sqrt{\varepsilon }+\cdots \to 1,\quad
\varepsilon \to 0+.
\end{equation*} \]

Для сращивания напряженности $E ( z, \varepsilon )$ в точке $z_{k}$ применим соотношение
\[ \begin{equation*}
E ( z_k +0, \varepsilon )=E ( z_k -0, \varepsilon ) ,
\end{equation*} \]
где в $E ( z_k +0, \varepsilon )$ используется решение в интервале $ (z_k , 1]$, а в $E ( z_ k -0, \varepsilon )$ — решение в расширенной ОПЗ $E ( z, \varepsilon )= \sqrt{2 ( I_z z- I )} / \sqrt{\varepsilon }$, которое продолжается в интервал $(z_m , z_ k )$. Таким образом получаем, что
\[ \begin{equation*}
\frac{1}{\sqrt{\varepsilon }}
\frac{4b\sqrt{\beta } \exp(- k_1 b ) }
{1-\beta \exp( -2k_1 b) }
=\frac{\sqrt{2 ( I_z - I)}}{\sqrt{\varepsilon }}.
\end{equation*} \]
Отсюда
\[ \begin{equation*}
k_1= -\frac{1}{b} \ln
\Bigl( \frac{-4b\sqrt{\beta }+\sqrt{8\beta ( 2 b ^ 2 + I_ z - I )}}
{2\beta \sqrt{2 ( I _ z - I )}} \Bigr).
\end{equation*} \]

Положим
\[ \begin{equation*}
z_ m =1- k_ 2 \sqrt{\varepsilon } | \ln\varepsilon |+\cdots \to 1,\quad \varepsilon \to 0+.
\end{equation*} \]

Условие сращивания концентрации $C_ 1 ( z, \varepsilon )$ в точке $z_{m}$ имеет вид
\[ \begin{equation*}
C_1 ( z_m -0, \varepsilon )= C_1 (z_m +0,\varepsilon ),
\end{equation*} \]
где в $C_1 (z_m -0, \varepsilon )$ используется решение в расширенной ОПЗ $C _ 1 ( z, \varepsilon ) = I_z \sqrt{\varepsilon }/ \sqrt{2 ( I_z z- I)}$, а в ${{C}_{1}}\left( {{z}_{m}}+0,\varepsilon \right)$ — решение из интервала $( z_k , 1]$, продолженное в интервал $( z_ m , z_ k )$. Откуда получим, что
\[ \begin{equation*}
k_2 =\frac{1}{2b}, \quad
\beta =\Bigl( \frac{ I_ z } {4 b ^ 2 \sqrt{2 ( I _ z - I )}} \Bigr)^ 2 .
\end{equation*} \]

Для нахождения $b$ используем условие $ C_ 1 ( 1,\varepsilon)= C_{1,m}$, откуда
\[ \begin{equation*}
b=\frac{1}{2} \sqrt{\frac{C_{1,m} (1-\beta )^2} {\sqrt{\beta } ( 1+\beta )}}.
\end{equation*} \]

3. Алгоритм гибридного численно-аналитического решения

1. Численно решаем модель без ОВК и находим $C_ 1 ( 1 , \varepsilon )$, $C_2 ( 1, \varepsilon )$.
2. Находим численно скачок потенциала для упрощенной модели в ОЧО. Находим скачок потенциала для базовой, используя соотношение
   \[ \begin{multline*}
   \Phi _0 = \int_0^1 E (z,\varepsilon )dz = {} \\
   {}=
   \int_0^{z_m} E(z,\varepsilon )dz+ \int_{z_m}^{z_ k} E (z,\varepsilon )dz+
   \int_{z_k}^1 E(z,\varepsilon )dz= {} \\
   {}= \int_0^{z_m} E(z,\varepsilon )dz + \int_{z_m}^{z_k} E(z,\varepsilon )dz-
   \int_{z_k}^1 \frac{d C_2 }{C_ 2}= {}
   \\
   {}=
   -\ln \frac{C_2(1,\varepsilon )}{C_2 (z_k ,\varepsilon )}+
   \int_0^{z_m} E (z,\varepsilon )dz+
   \int_{ z_ m }^{ z_ k } E (z,\varepsilon )dz .
   \end{multline*} \]
   С учетом того, что $z _k \approx z _m \approx 1$, получаем
   \[ \begin{equation*}
   \Phi_0 \approx -\ln \frac{C_2 (1,\varepsilon )}{ C_ 2 ( z_k , \varepsilon )}+
   \int_0^1 E (z,\varepsilon )dz \quad
   \text{или} \quad \Phi _0\approx \Phi_{\text{ОВК}} + \Phi_{\text{ОЧО}}.
   \end{equation*} \]
   Здесь первое слагаемое $\Phi_{\text{ОВК}}$ — скачок потенциала в области возрастания катионов (около КОМ), а второе слагаемое $ \Phi_{\text{ОЧО}}$ — скачок потенциала, численно рассчитанный в ОЧО.
3. Находим аналитическое решение в ОВК по формулам (12)–(14) при допредельном токовом режиме и по формулам (15), (16) при сверхпредельном.
4. Используя шаги 1) и 3), получаем решение базовой задачи.

Заключение

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

Конкурирующие интересы. Заявляем, что в отношении авторства и публикации этой статьи конфликта интересов не имеем.
Авторский вклад и ответственность. Е.В. Казаковцева — разработка численно-аналитического метода и программы численных экспериментов; компьютерное и математическое моделирование; интерпретация полученных результатов; работа с черновиком и переработанным вариантом рукописи. А.В. Коваленко — компьютерное моделирование; обработка и анализ экспериментальных данных; подготовка первичного варианта рукописи. А.В. Письменский — обработка и анализ экспериментальных данных, интерпретация полученных результатов. М.Х. Уртенов — идея исследования; математическое моделирование; разработка численно-аналитического метода; интерпретация полученных результатов; подготовка первичного варианта рукописи; работа с черновиком и переработанным вариантом рукописи. Авторы несут полную ответственность за предоставление окончательной рукописи в печать. Окончательная версия рукописи была одобрена всеми авторами.
Финансирование. Исследование выполнено за счет гранта Российского научного фонда № 24–19–00648, https://rscf.ru/project/24-19-00648/.
Благодарность. Авторы благодарны рецензентам за тщательное прочтение статьи, ценные предложения и комментарии.

---

¹Эти данные соответствуют конструкции электрохимической ячейки с вращающимся мембранным диском, приведенной в [2].

About the authors

Ekaterina V. Kazakovtseva

Kuban State University

Author for correspondence.
Email: vivkaterina@mail.ru
ORCID iD: 0009-0003-0040-0880
SPIN-code: 4895-4042
https://www.mathnet.ru/person208186

Senior Lecturer; Dept. of Data Analysis and Artificial Intelligence

Russian Federation, 350040, Krasnodar, Stavropolskaya st., 149

Anna V. Kovalenko

Kuban State University

Email: savanna-05@mail.ru
ORCID iD: 0000-0002-3991-3953
SPIN-code: 3693-4813
Scopus Author ID: 55328224000
http://www.mathnet.ru/person112835

Dr. Techn. Sci., Associate Professor; Head of Department; Dept. of Data Analysis and Artificial Intelligence

Russian Federation, 350040, Krasnodar, Stavropolskaya st., 149

Alexander V. Pismenskiy

Kuban State University

Email: archer812@mail.ru
ORCID iD: 0000-0003-4046-2229
SPIN-code: 9932-7747
Scopus Author ID: 13004856800
https://www.mathnet.ru/person208187

Cand. Phys. & Math. Sci.; Head of Department; Dept. of Applied Mathematics

Russian Federation, 350040, Krasnodar, Stavropolskaya st., 149

Makhamet Kh. Urtenov

Kuban State University

Email: urtenovmax@mail.ru
ORCID iD: 0000-0002-0252-6247
SPIN-code: 7189-0748
Scopus Author ID: 6603363090
http://www.mathnet.ru/person119069

Dr. Phys. & Math. Sci., Professor; Dept. of Applied Mathematics

Russian Federation, 350040, Krasnodar, Stavropolskaya st., 149

References

Supplementary files