Математическое моделирование массопереноса в электромембранных системах в гальванодинамическом режиме с учетом электроконвекции и реакции диссоциации/рекомбинации молекул воды

Полный текст

1. Введение

Электромембранные системы (электродиализные аппараты, нано- и микрофлюидные устройства) применяются для обессоливания и деионизации воды, производства кислот и оснований, при преобразовании и производстве энергии и в других областях [1–3]. В основе принципов функционирования электромембранных систем лежит селективный перенос ионообменных мембран, то есть проницаемость для ионов одного знака (противоионов) и ограничение движения ионов другого знака (коионов). В структуре вольт-амперной кривой электромембранной системы выделяются следующие характерные участки [4, 5]:

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

Согласно современным представлениям основным механизмом сверхпредельного массопереноса для рассматриваемых в данной работе разбавленных растворов электролитов является электроконвекция (ЭК) [6, 7], то есть увлечение молекул жидкости ионами, образующими пространственный заряд у ионоселективной поверхности под действием электрической силы [6]. В режиме сверхпредельного тока на границе раздела раствор/мембрана формируется расширенная область пространственного заряда (ОПЗ), действие электрического поля на который генерирует значимое электроконвективное течение [6], которое интенсифицирует процесс массопереноса. Наряду с появлением ОПЗ при сверхпредельных токах фиксируется смещение показателя $\mathrm{pH}$ раствора, что объясняется диссоциацией молекул воды с образованием ионов водорода $\mathrm{H}^+$ и гидроксила $\mathrm{OH}^-$ [8]. Отклонение от равновесия в реакции $\mathrm{H}_2\mathrm{O} ↔^{k_d}_{k_r} \mathrm{H}^+ + \mathrm{OH}^-$ (где $k_d$ и $k_r$ — константы скорости диссоциации и рекомбинации воды соответственно) приводит к появлению новых носителей заряда (ионов $\mathrm{H}^+$ и $\mathrm{OH}^-$), что снижает эффективность обессоливания растворов [9], уменьшает толщину расширенной ОПЗ и величину заряда [10, 29], вызывает уменьшение интенсивности ЭК [12, 13]. Поэтому особое значение имеет математическое моделирование сверхпредельного переноса ионов с учетом ЭК и реакции диссоциации/рекомбинации молекул воды (РДРВ), которые оказывают противоположное действие на эффективность электромембранных процессов.

Одним из основных подходов математического описания переноса ионов в мембранных системах является использование уравнения Нернста–Планка, учитывающего диффузию и миграцию ионов, а также конвекцию раствора. Уравнение Нернста–Планка совместно с уравнением материального баланса описывают как стационарные, так и нестационарные процессы, учитывая, в частности, протекание химических реакций. На основе уравнений Нернста–Планка и условия электронейтральности построены модели переноса ионов соли и ионов водорода и гидроксила, образующихся в результате расщепления воды, в стационарном состоянии трехслойной системы (ионообменная мембрана и два смежных с ней слоя раствора электролита) Г. Гроссманом [14] и в неперемешиваемом слое электролита у ионообменной мембраны И. Рубинштейном [15].

И. Рубинштейн и Л. Штильман на основе системы уравнений Нернста–Планка и уравнения Пуассона для электрического потенциала построили одномерную модель переноса ионов $1\,{:}\,1$-валентного электролита в неперемешиваемом слое раствора у поверхности идеально селективной мембраны при постоянном скачке потенциала в стационарном состоянии [16]. Численное решение соответствующей краевой задачи модели показало существование расширенной ОПЗ у границы раствор электролита/мембрана, которая значительно больше, чем область равновесного двойного электрического слоя, при наложении достаточно высокого скачка потенциала [16]. В. И. Заболоцкий и соавт. [10] обобщили одномерную модель И. Рубинштейна и Л. Штильмана (то есть модель переноса ионов соли $1\,{:}\,1$-валентного электролита в неперемешиваемом слое раствора у поверхности идеально селективной мембраны) [16] для описания переноса ионов четырех видов: катионов и анионов соли, ионов водорода и гидроксила, образующихся в результате диссоциации молекул воды. Показано, что появление ионов водорода и гидроксила в обедненном диффузионном слое уменьшает величину пространственного заряда и толщину ОПЗ [10]. Помимо численного анализа математической модели переноса ионов в стационарном состоянии для анализа влияние РДРВ на перенос также получено асимптотическое решение [17, 18].

М. Х. Уртенов и соавт. [19] с использованием модели переноса ионов с учетом РДРВ в диффузионном слое у ионообменной мембраны в стационарном состоянии показали, что интенсивное некаталитическое расщепление воды происходит в расширенной ОПЗ и условие химического равновесия выполняется всюду в диффузионном слое, за исключением области рекомбинации и ОПЗ. Позже подобный подход (на основе системы уравнений Нернста–Планка–Пуассона) применен для описания переноса в одномерном сечении электродиализного канала, образованного между анионо- и катинообменной мембранами для стационарного [20] и нестационарного [21] случаев. Результаты численного моделирования работ [20, 21] показали, что РДРВ может создать ОПЗ в центральной части мембранного канала обессоливания.

В работах [22, 23] математическая модель стационарного переноса ионов соли с учетом РДРВ в диффузионном слое [19] была модифицирована с целью описания влияние температурных эффектов, связанных с РДРВ: система уравнений Нернста–Планка и Пуассона дополнена уравнением теплопроводности, а коэффициент скорости диссоциации определен как функция температуры.

В. В. Никоненко и соавт. в работе [24] получили численное и полуаналитическое решение задачи описания стационарного переноса ионов и ионов водорода и гидроксила, образующихся в результате диссоциации воды. Рассматриваемая в работе [24] модель основана на системе уравнений Нернста–Планка и Пуассона и учитывает отклонение от локальной электронейтральности в обедненном диффузионном слое. Аналитическое решение получено путем разделения обедненного диффузионного слоя на электронейтральную часть, расширенную область пространственного заряда и квазиравновесную зону у поверхности мембраны [24].

А. В. Коваленко в работе [25] для исследования влияния РДРВ и ЭК предложена двумерная модель на основе системы уравнений Нернста–Планка–Пуассона–Навье–Стокса для половины канала обессоливания у катионообменной мембраны. Было показано, что диссоциация воды приводит к увеличению значения скачка потенциала, при котором развивается ЭК [25]. В работе [21] область рассмотрения данной модели была расширена, она включает канал обессоливания от анионо- до катионообменной мембраны. Показано, что когда концентрации ионов водорода и гидроксила становятся сопоставимыми с концентрацией ионов соли, РДРВ может вызвать появление ОПЗ в центральной части канала обессоливания. В недавней работе [26] на основе модели переноса в электродиализном канале обессоливания исследовано совместное влияние РДРВ, ЭК и спейесеров на массоперенос.

В работах [8–26] заложены основы теоретического анализа влияния РДРВ на массоперенос в электромембранных системах. Работы [21, 25, 26] посвящены изучению массопереноса с учетом как ЭК, так и РДРВ. В них выполнено исследование только для потенциодинамического (когда задается скачок потенциала) режима электрического поля. При этом в исследованиях и практике функционирования электромембранных систем важное значение имеет и гальванодинамический режим (когда задается плотность протекающего через систему тока), так как для характеристики процесса массопереноса используются различные конкретные критические значения тока, такие, например, как предельный ток, экзальтационный ток, ток Харкаца и т.д.

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

2. Математическая модель

Рассмотрим покоящийся слой разбавленного раствора бинарного электролита у поверхности катионообменной мембраны, через который протекает постоянный ток плотностью $i$. Селективные свойства мембраны зададим граничными условиями. Плотность, температуру и диэлектрическую проницаемость раствора примем постоянными. Пусть $x$ и $y$ — нормальная и касательная координаты к поверхности мембраны соответственно; $H$ и $L$ — толщина и длина рассматриваемой области, рис. 1.

Рис. 1. Двумерное сечение слоя раствора $\mathrm{NaCl}$ у поверхности катионообменной мембраны (КОМ). Схематически показаны потоки ионов соли $\mathrm{Na}^+$, $\mathrm{Cl}^-$, а также водорода $\mathrm{H}^+$ и гидроксила $\mathrm{OH}^-$, образующихся в расширенной области пространственного заряда (ОПЗ) в сверхпредельном режиме
[Figure 1. Two-dimensional cross-section of a layer of $\mathrm{NaCl}$ electrolyte at the surface of a cation-exchange membrane (CEM). The flows of salt ions $\mathrm{Na}^+$, $\mathrm{Cl}^-$, as well as hydrogen $\mathrm{H}^+$ and hydroxyl $\mathrm{OH}^-$, generated in the extended space charge region (SCR) in the overlimiting mode]

Плотность потоков ионов соли (катионов $\mathrm{Na}^+$, $n=1$, и анионов $\mathrm{Cl}^-$, $n=2$) и продуктов диссоциации воды (ионов $\mathrm{H}^+$, $n=3$, и $\mathrm{OH}^-$, $n=4$) определяется уравнениями Нернста–Планка
\[ \begin{equation}
\vec{j_n}=-\frac{F}{RT}z_n D_n c_n \nabla \varphi - D_n \nabla c_n + c_n \vec{V}
\end{equation} \tag{1} \]
и материального баланса [27]
\[ \begin{equation}
\frac{\partial c_n}{\partial t}=-\operatorname{\mathrm{div}} \vec{j}_n +R_n,
\end{equation} \tag{2} \]
где $\vec{j}_n$, $c_n$, $D_n$, $z_n$ — плотность потока, концентрация, коэффициент диффузии и зарядовое число $n$-го иона соответственно; $\varphi$ — потенциал электрического поля; $\vec{V}$ — скорость течения раствора; $F$ — постоянная Фарадея; $R$ — газовая постоянная; $T$ — абсолютная температура; $R_n$ — скорость формирования $n$-го иона. Скорости формирования ионов соли равны нулю ($R_1=R_2=0$), так как они не участвуют в химической реакции. Скорость диссоциации молекул воды определяется величиной $k_d c_{\mathrm{H}_2\mathrm{O}}$ (где $c_{\mathrm{H}_2\mathrm{O}}$ — концентрация воды), скорость рекомбинации — $k_r c_3 c_4$. Результирующая скорость образования ионов $\mathrm{H}^+$ и $\mathrm{OH}^-$ равна разности скоростей диссоциации и рекомбинации: $R_3=R_4=k_d c_{\mathrm{H}_2\mathrm{O}}-k_r c_3 c_4=k_r(K_w-c_3 c_4)$, где $K_w$ — ионное произведение воды ($K_w=10^{-8} \text{моль}^2/\text{м}^6$).

Поле электрического потенциала с учетом пространственного заряда, обусловленного нарушением баланса катионов и анионов, определяется уравнением Пуассона
\[ \begin{equation}
\varepsilon_0 \varepsilon_r \Delta \varphi=-F\sum_{n=1}^4 z_n c_n,
\end{equation} \tag{3} \]
где $\varepsilon_0$ — электрическая постоянная, $\varepsilon_r$ — диэлектрическая проницаемость раствора, правая часть уравнения (3) представляет собой плотность пространственного заряда $\rho=F\sum_{n=1}^4 z_n c_n$.

Течение раствора электролита с учетом электрической силы, действующей на пространственный заряд, описывается уравнениями Навье–Стокса и неразрывности несжимаемой жидкости:
\[ \begin{equation}
\frac{\partial \vec{V}}{\partial t}+(\vec{V} \nabla) \vec{V}=-\frac{1}{\rho_0}\nabla P + \nu \Delta \vec{V}-\frac{1}{\rho_0}\vec{f_e}, \quad \operatorname{\mathrm{div}} \vec{V}=0,
\end{equation} \tag{4} \]
где $P$, $\rho_0$, $\nu$ — давление, плотность и кинематическая вязкость раствора электролита; $\vec{f}_e=\rho \nabla \varphi$ — плотность силы электрического поля, действующей на заряд плотностью $\rho$.

Величины $\varphi$, $\vec{V}$, $P$, $\vec{j}_n$, $c_n$, $n=1, 2, 3, 4$, в системе уравнений (1)–(4) являются неизвестными функциями координат $x$, $y$ и времени $t$.

В предлагаемую математическую модель включены следующие граничные условия: для скорости течения раствора задается условие прилипания на границе раздела раствор/мембрана ($x=H$) и нулевая нормальная скорость на остальных границах:
\[ \begin{equation}
\vec{V}(H,y,t)=0, \quad
\vec{V} \cdot \vec{n}(0,y,t)=0, \quad
\vec{V} \cdot \vec{n}(x,0,t)=0, \quad
\vec{V} \cdot \vec{n}(x,L,t)=0,
\end{equation} \tag{5} \]
где $\vec{n}$ — нормаль к границе.

На границе $x=0$ концентрации ионов соли обоих знаков фиксированы и равны начальной концентрации электролита; концентрации ионов $\mathrm{H}^+$ и $\mathrm{OH}^-$ равны $\sqrt{K_w}=10^{-4} \text{ моль}/\text{м}^3$ (что соответствует нейтральному раствору):
\[ \begin{equation}
c_1(0,y,t)=c_0, \quad
c_2(0,y,t)=c_0, \quad
c_3(0,y,t)=\sqrt{K_w},\quad
c_4(0,y,t)=\sqrt{K_w}.
\end{equation} \tag{6} \]

Перечислим граничные условия при $x = H$: концентрация катионов соли задается постоянной; условие непрерывности потока ионов на границе раздела раствор/мембрана для анионов соли; ионы $\mathrm{H}^+$, $\mathrm{OH}^-$, формирующиеся в результате диссоциации воды, беспрепятственно переносятся из зоны реакции через катионообменную мембрану [21]; рассматриваются некаталитическая диссоциация и плотность потока инжекция ионов гидроксила $\mathrm{OH}^{-}$ с поверхности мембран, $j_{4c}$ не учитывается:
\[ \begin{equation}
\begin{array}{c}
c_1(H,y,t)=N_c c_0, \\
\Bigl( -\dfrac{F}{RT} z_2 D_2 c_2 \dfrac{\partial \varphi}{\partial x}-D_2 \dfrac{\partial c_2}{\partial x} + c_2V_y \Bigr) (H,y,t)=\dfrac{T_{2c}}{F z_2}i_x(H,y,t), \\
\dfrac{\partial c_3}{\partial x} c_1(H,y,t)=0, \\
\Bigl( -\dfrac{F}{RT} z_4 D_4 c_4 \dfrac{\partial \varphi}{\partial x}-D_4 \dfrac{\partial c_4}{\partial x} + c_4V_y \Bigr) (H,y,t)=j_{4c}=0,
\end{array}
\end{equation} \tag{7} \]
где $(H, y, t)$ — аргументы соответствующих выражений.

На границах $y=0$ и $y=L$ приняты «мягкие» условия [28, 29]:
\[ \begin{equation}
\frac{\partial c_n}{\partial y} (x,0,t)=0, \quad
\frac{\partial c_n}{\partial y} (x,L,t)=0, \quad n=1, 2, 3, 4.
\end{equation} \tag{8} \]

При моделировании гальванодинамического режима с использованием уравнения Пуассона для потенциала одна из продольных границ системы, параллельных поверхности мембраны, предполагается эквипотенциальной (задается некоторое постоянное значение скачка потенциала, например ноль), а на другой границе задается граничное условие, связывающее производную потенциала и плотность тока. Зададим нулевой потенциал на границе $x=H$:
\[ \begin{equation}
\varphi (H,y,t)=0.
\end{equation} \tag{9} \]

Гальванодинамическое граничное условие на границе $x = 0$ выведем на основе следующего соотношения для плотности полного тока [30, 31]:
\[ \begin{equation}
\vec{i}_{tot}=F\sum_{n=1}^4 z_n \vec{j}_n-\varepsilon_0\varepsilon_r \frac{\partial(\nabla \varphi)}{\partial t},
\end{equation} \tag{10} \]
где первое слагаемое в правой части является током проводимости, а второе — током смещения.

Пусть граница рассматриваемого слоя раствора электролита с объемом ($x=0$) достаточно удалена от поверхности мембраны, у которой формируется ОПЗ и развивается ЭК при протекании тока. Тогда на этой границе плотность тока смещения пренебрежимо мала, а также можно принять равномерное распределение нормальной составляющей плотности тока, т. е. $i_x(0,y,t)=i$, где $i$ — плотность задаваемого тока.

С учетом данных упрощений из уравнения (10) можно записать следующее граничное условие при $x=0$:
\[ \begin{equation}
\frac{\partial \varphi}{\partial x} (0,y,t)=-\frac{RT}{F^2}
\Bigl( \frac{i+\sum_{n=1}^4 F z_n D_n \frac{\partial c_n}{\partial x}}{\sum_{n=1}^4 z^2_nD_nC_n} \Bigr) (0,y,t).
\end{equation} \tag{11} \]

На границах $y=0$ и $y=L$ для потенциала устанавливается «мягкое» граничное условие
\[ \begin{equation}
\frac{\partial \varphi}{\partial x} (x,0,t)=0,\quad
\frac{\partial \varphi}{\partial x} (x,L,t)=0.
\end{equation} \tag{12} \]

Начальные условия соответствуют предположению, что концентрации ионов соли равны начальной концентрации раствора электролита во всей рассматриваемой области (кроме границы $x=H$); концентрации ионов водорода и гидроксила равны $10^{-4} \text{моль}/\text{м}^3$, что соответствует нейтральному раствору:
\[ \begin{equation}
\begin{array}{l}
c_1(x,y,0)=c_2(x,y,0)=c_0, \\
c_3(x,y,0)=c_4(x,y,0)=\sqrt{K_w}, \\
\varphi(x,y,0)=0, \quad 0\leqslant x<H.
\end{array}
\end{equation} \tag{13} \]

Численное решение краевой задачи (1)–(9), (11)–(13) найдено методом конечных элементов с использованием пакета Comsol Multiphysics [28]. Уравнения для концентрации ионов (2) и уравнение Пуассона (3) построены на основе модуля общей формы уравнения в частных производных General Form PDE; уравнения Навье–Стокса и неразрывности (4) реализованы с помощью модуля ламинарного течения Laminar Flow [28]. Дискретизация по пространственным координатам полей концентрации и потенциала использует квадратичные интерполяционные функции Лагранжа; в модуле ламинарного течения используются элементы второго порядка для составляющих скорости и линейные элементы для поля давления [32]. Вычисления по времени реализуются с использованием сегрегированного узла с неявным методом выбора временного шага BDF (формула обратного дифференцирования) [32]. Каждая отдельная итерация включает выполнение двух отдельных шагов: первый шаг — расчет концентраций ионов и потенциала; на втором этапе рассчитываются скорость и давление. Шаг по времени определяется так, чтобы выполнялось требование для относительной погрешности (устанавливается $10^{-3}$). Кроме того, при решении краевых задач моделей устанавливалось ограничение на максимальный шаг по времени, чтобы погрешность расчета плотности тока на границах составляла менее 1 %. Все расчеты проводились с использованием процессора Intel(R) Core(TM) i9-10900K (10 ядер 3.70 ГГц), 64 ГБ оперативной памяти.

Влияние РДРВ на характеристики массопереноса зависит от соотношения концентрации ионов соли и продуктов диссоциации, а следовательно, и от исходной концентрации раствора электролита. При описании процессов массопереноса в мембранных системах плотность тока рассматривается относительно предельного значения, которое прямо пропорционально исходной концентрации электролита $c_0$ [33]. Поэтому для оценки влияния РДРВ на массоперенос выполнены расчеты при различных значениях концентрации $c_0$ для следующих значений плотности тока:
\[ \begin{equation*}
\begin{array}{ll}
i=10 \text{А}/\text{м}^2 &\text{для} ~c_0=10 \text{моль}/\text{м}^3; \\
i=1 \text{А}/\text{м}^2 &\text{для} ~c_0=1 \text{моль}/\text{м}^3; \\
i=0.1 \text{А}/\text{м}^2 &\text{для} ~c_0=0.1 \text{моль}/\text{м}^3.
\end{array}
\end{equation*} \]

Рассматривается раствор электролита $\mathrm{CuSO}_4$ при температуре $T=298$ K плотностью $\rho_0=1002$ кг/м$^3$ и вязкостью $\nu =0.89 \cdot 10^{-6}$ м$^2$/с; коэффициенты диффузии ионов задаются значениями $D_1=0.72 \cdot 10^{-9}$ м$^2$/c, $D_2 =1.065\cdot 10^{-9}$ м$^2$/c, $D_3=9.34 \cdot 10^{-9}$ м$^2$/c, $D_4=5.23 \cdot 10^{-9}$ м$^2$/c; числа переноса катионов в катионообменной мембране и растворе имеют значения $T_{1C} = 1$ и $t_1=0.403$ соответственно; зарядовые числа ионов следующие: $z_1=1$, $z_2= -1$, ${z_3=1}$, $z_4= -1$. Толщина и длина рассматриваемого слоя раствора электролита ${H=10^{-3}}$ м и $L=10^{-3}$ м соответственно.

Указанные параметры системы соответствуют экспериментальному исследованию, выполненному в работе [34], с той разницей, что длина рассматриваемой системы выбрана меньшей. Это позволило уменьшить вычислительную трудоемкость задачи без значительного искажения хронопотенциограммы, так как вынужденное течение раствора отсутствует и обессоливание раствора не зависит от продольной координаты. Для упрощения численного решения отношение концентрации катионов при $x=H$ к $c_0$ (используемое в граничном условии (7)) было принято равным $N_c=1$. Это значение меньше, чем в реальных системах [35]; однако в работе [36] показано, что значение $N_c$ существенно не влияет на решение задачи в расширенной ОПЗ.

Большие градиенты концентраций ионов и потенциала электрического поля вблизи границы раздела раствор/мембрана обуславливают значительную вычислительную сложность рассматриваемой задачи. Поэтому для упрощения численного решения задачи была построена неоднородная расчетная сетка, в которой 33 пограничных слоя с коэффициентом растяжения 1.16 и толщиной первого слоя 0.001 мкм размещены у правой границы. В остальной части рассматриваемой области расчетная сетка построена с использованием треугольной сетки (опция Free Triangular) [32]. Полученная сетка состоит из 163 911 элементов области и 2 811 граничных элементов. Структура расчетной сетки определена на основе следующего алгоритма:

• проводятся расчеты хронопотенциограммы $d \varphi_{k_1}(t)$ при некоторой плотности тока на интервале времени от $0$ до установления квазистационарного состояния $\tilde{t}$ для расчетной сетки, состоящей из $k_1$ элементов;
• затем сетка уточняется (и содержит теперь $k_2$ элементов), проводятся расчеты значения $d \varphi_{k_2}(t)$;
• затем определяется максимальная относительная погрешность расчета скачка потенциала по формуле $\gamma_1=\max\limits_{t \in [0,\tilde{t}]} ( |d \varphi_{k_1}-d\varphi_{k_2}|/d \varphi_{k_2} ) \cdot 100\%$.

Приведенные шаги повторяются до тех пор, пока не выполнено условие ${\gamma_{k_i}\leqslant 5\,\%}$.

3. Результаты

Чтобы оценить влияние ЭК и РДРВ на массоперенос, для каждого значения концентрации электролита выполнены следующие варианты расчетов:

• без учета ЭК и РДРВ;
• с учетом ЭК и без учета РДРВ;
• без учета ЭК с учетом РДРВ;
• с учетом ЭК и РДРВ.

На рис.  2, a показаны хронопотенциограммы, рассчитанные для указанных вариантов учета механизмов сверхпредельного тока, а также хронопотенциограмма, полученная экспериментально для катионообменной мембраны Neosepta CMX и раствора электролита $\mathrm{CuSO}_4$ концентрацией 10 моль/м$^3$ при плотности тока 10 А/м$^2$ и площади мембраны $3\times4$ мм$^2$ [34].

Рис. 2. Хронопотенциограммы $d_{\varphi}'(t)$ и средние значения толщины электроконвективного вихревого слоя $d_{ec}(t)$, рассчитанные для концентрации $c_0=10$ моль/м$^3$, при плотности тока $i=10$ А/м$^2$ (a), (b); $c_0=1$ моль/м$^3$, $i=1$ А/м$^2$ (c), (d); $c_0=0.1$ моль/м$^3$, $i =0.1$ А/м$^2$ (e), (f). Расчеты с учетом (сплошные линии) и без учета (пунктирные линии) ЭК, с учетом (синие лини) и без учета (красные линии) РДРВ. На рис. (a) также приведена экспериментальная хронопотенциограмма для катионообменной мембраны Neosepta CMX в растворе $\mathrm{CuSO}_4$ 10 моль/м$^3$ (черная линия) при $i=10$ А/м$^2$, полученная в работе [34]
[Figure 2. Chronopotentiograms $d_{\varphi}'(t)$ and average values of the thickness of the electroconvective vortex layer $d_{ec}(t)$, calculated for concentration $c_0=10$ mol/m$^3$ at current density $i=10$ A/m$^2$ (a), (b); $c_0=1$ mol/m$^3$, $i=1$ A/m$^2$ (c), (d); $c_0=0.1$ mol/m$^3$, $i=0.1$ A/m$^2$ (e), (f). Calculations with (solid lines) and without (dashed lines) considering EC, and with (blue lines) and without (red lines) considering DRRW. Fig. (a) also shows an experimental chronopotentiogram for the Neosepta CMX cation exchange membrane in the 10 mol/m$^3$ $\mathrm{CuSO}_4$ (black line) at $i=10$ A/m$^2$ from reference [34]

На хронопотенциограмме, рассчитанной для $c_0=10$ моль/м$^3$ с учетом ЭК и РДРВ, можно выделить участок начального быстрого омического роста скачка потенциала (при $t < t'=10^{-8}$ c). На рис. 2 этот участок отброшен и используется приведенный скачок потенциала $d_\varphi'(t)=d_\varphi(t)-d_\varphi(t')$, так как это исключает влияние омического сопротивления, которое является функцией расстояния между измерительными электродами, толщины мембраны и некоторых других параметров, которые трудно найти, тогда как они не имеют значения для поведения мембраны и не учитываются в модели.

Далее следует участок монотонного роста скачка потенциала, вызванный электродиффузионным обессоливанием раствора. Скорость роста скачка потенциала быстро возрастает после уменьшения концентрации на границе раздела раствор/мембрана. В момент времени $\tau \approx 71.15$ с касательная к профилям концентрации в электронейтральной области вблизи поверхности мембраны проходит через 0 при $x = H$.

После этого момента начинается формирование расширенных ОПЗ, а на хронопотенциограмме отмечается быстрое увеличение скачка потенциала. При $t \approx 71.5$ с развивается электроконвективное течение, что вызывает снижение скачка потенциала, за которым следует его рост с меньшей скоростью; рост скачка потенциала сопровождается увеличением интенсивности ЭК.

С течением времени скачок потенциала начинает колебаться относительно некоторого значения, система переходит в квазистационарное состояние.

В экспериментальных измерениях хронопотенциограмм ионообменных мембран [34] момент времени, предшествующий быстрому росту скачка потенциала, определяется как переходное время. Для аналитической оценки времени перехода используется уравнение
\[ \begin{equation}
\tau_S=\frac{\pi D}{4} \frac{c_0 F z_1}{T_{1C}-t_1} \frac{1}{i^2}.
\end{equation} \tag{14} \]

Уравнение (14) было получено Сандом (Henry J. S. Sand) на основе теоретического анализа переноса ионов в бесконечном диффузионном слое, в котором выполняется условие локальной электронейтральности [37]. Согласно теории Санда, переходное время определяется по моменту, когда концентрация противоионов у поверхности мембраны уменьшается до нуля, а скачок потенциала быстро возрастает. Для рассматриваемых параметров системы согласно соотношению (14) переходное время Санда $\tau_S \approx 70.55$ с. Переходное время, рассчитанное по предлагаемой модели, определено по моменту времени, когда касательная к концентрационным профилям (в произвольном нормальном сечении рассматриваемого слоя, например, в сечении при $y=0.5L$) в электронейтральной области вблизи поверхности мембраны проходит через 0 при $x=H$.

Для всех рассматриваемых значений концентрации раствора электролита (см. таблицу) значения переходного времени, рассчитанные без учета ЭК и РДРВ, совпадают друг с другом, а также с аналитической оценкой Санда (отличие составляет менее 0.03 %). Для концентрации раствора $c_0=10$ моль/м$^3$ ЭК и РДРВ слабо влияют на переходное время: для всех четырех типов расчетов различие значений переходного времени не превышает 1 %. Для концентрации раствора $c_0=1$ моль/м$^3$ ЭК увеличивает переходное время (0.2 с, то есть 0.3 %), а РДРВ увеличивает на 3.63 с ($\approx$ 5 %). Для концентрации раствора $c_0=0.1$ моль/м$^3$ в расчетах с учетом РДРВ система стабилизируется со временем без достижения предельного состояния.

Структура рассчитанных хронопотенциограмм качественно совпадает со структурой экспериментальной кривой [34], в которой выделяются следующие характерные участки (рис. 2): медленный рост скачка потенциала, предельное состояние быстрого роста скачка потенциала, переходная стадия развития ЭК и квазистационарное состояние.

Экспериментальное переходное время составляет $\tau_{exp} \approx 77$ с, что на 8.4 % больше значения, рассчитанного с использованием модели. Подобное поведение экспериментальной и расчетной кривых наблюдается при $t < \tau$. Переходная стадия развития ЭК отличается более существенно.

На рассчитанной хронопотенциограмме в момент времени $\tau$ наблюдается резкое снижение скачка потенциала, связанное с развитием ЭК. В экспериментальной ХП скорость роста скачка потенциала плавно уменьшается из-за появления дополнительных механизмов, основным из которых является ЭК [34]. По-видимому, это различие связано с геометрической и электрической неоднородностью поверхности, что влияет на развитие ЭК, в то время как в расчете предполагается, что поверхность мембраны является совершенно однородной.

Различие хронопотенциограмм, рассчитанных с учетом и без учета РДРВ (рис. 2), можно объяснить с использованием временных зависимостей отношения концентраций продуктов диссоциации и ионов соли. На рис. 3 приведены зависимости от времени отношений средних значений концентрации ионов водорода $c_{3 av}$ и гидроксила $c_{4 av}$ к соответствующим значениям ионов соли $c_{1 av}$ и $c_{2 av}$, которые рассчитаны с учетом ЭК и РДРВ по формуле
\[ \begin{equation*}
c_{n av}(t)=\frac{1}{HL} \int_0^H \int_0^L c_n(x,y,t) dxdy, \quad n=1, 2, 3, 4.
\end{equation*} \]

Рис. 3. Отношения концентраций ионов водорода и катиона соли $c_{3av}/c_{1av}$ (a) и концентраций ионов гидроксила и аниона соли $c_{4av}/c_{2av}$ (b). Рис. (c) — увеличение фрагмента рис. (b). Результаты расчета с учетом (сплошные линии) и без учета ЭК (пунктирные линии) и с учетом РДРВ для $c_0=10$ моль/м$^3$, $i=10$ А/м$^2$, $c_0=1$ моль/м$^3$, $i=1$ А/м$^2$, $c_0=0.1$ моль/м$^3$, $i=0.1$ А/м$^2$
[Figure 3. The relationship between the concentrations of hydrogen ions and the cation of the salt $c_{3av}/c_{1av}$ (a) and the concentrations of hydroxyl ions and the anion of the salt $c_{4av}/c_{2av}$ (b). Fig. (c) shows an enlarged fragment of Fig. (b). The calculation results are presented with (solid lines) and without (dashed lines) consideration of EC, and with consideration of RDRW DRRW for $c_0=10$ mol/m$^3$, $i=10$ A/m$^2$, $c_0=1$ mol/m$^3$, $i=1$ A/m$^2$, $c_0=0.1$ mol/m$^3$, $i=0.1$ A/m$^2$]

Из рис. 3 видно, что для всех рассматриваемых значений концентрации электролита $c_0$ с течением времени отношение концентраций ионов водорода и катиона соли $c_{3av}/c_{1av}$ убывает (и не превышает значения 0.001), а отношение концентраций ионов гидроксила и аниона соли $c_{4av}/c_{2av}$ возрастает до значений, способных оказывать значимое влияние на процесс массопереноса. При этом меньшей концентрации $c_0$ соответствует большее значение $c_{4av}/c_{2av}$, например, при $t=100$ с для $c_0=10$ моль/м$^3$ значение $c_{4av}/c_{2av} \approx 0.004$; для $c_0=1$ моль/м$^3$ — $c_{4av}/c_{2av} \approx 0.04$; для $c_0=0.1$ моль/м$^3$ — $c_{4av}/c_{2av} \approx 0.01$.

Таким образом, при протекании тока в слое раствора электролита у поверхности катионообменной мембраны влияние РДРВ проявляется возрастанием концентрации ионов гидроксила.

Скорость роста отношения $c_{4av}/c_{2av}$ в расчете с учетом ЭК и РДРВ для $c_0=10$ моль/м$^3$ увеличивается после достижения переходного времени, при этом максимальное значение $c_{4av}/c_{2av}$ в рассматриваемом временном диапазоне не превышает 0.004 (рис. 3, c). Поэтому при $t \leqslant \tau$ отличие значений скачка потенциала в системе не превышает 3 %, а при $t > \tau$ отличие мгновенных значений скачка потенциала может достигать 16 %, но диапазон колебаний скачка потенциала совпадает.

Для исследования влияния РДРВ на ЭК в расчетах с учетом и без учета РДРВ рассчитана усредненная по длине канала толщина электроконвективного вихревого слоя $d_{ec}$ (рис. 2, b, d, f). Граница электроконвективного слоя определялась по точкам, в которых разница среднеквадратического значения скорости и ее максимального значения $V_{\max}$ (на рассматриваемом временном отрезке $V_{\max} \approx 0.0044$ м/с) превышает 1 % от $V_{\max}$. Так как диапазон колебаний скачка потенциала при $t > \tau$ в расчетах с учетом и без учета РДРВ для $c_0=10$ моль/м$^3$ приблизительно совпадает, колебания $d_{ec}$ также находится приблизительно в одном диапазоне значений (рис. 2, b).

В расчете для концентрации $c_0=1$ моль/м$^3$ значение отношения концентраций $c_{4av}/c_{2av}$ на порядок выше, чем для $c_0=1$ моль/м$^3$ (рис. 2, b). Образование дополнительных носителей заряда в результате диссоциации молекул воды в таком количестве снижает скачок потенциала в слое раствора электролита (рис. 2, c), что уменьшает толщину электроконвективного вихревого слоя (рис. 2, d).

Электроконвективные вихри, в области у поверхности мембраны перемешивая обедненный раствор электролита с раствором из объема (рис. 4), в свою очередь, приводят к уменьшению скачка потенциала (рис. 2, a, c, e), в результате чего уменьшается концентрация ионов гидроксила (рис. 3, b, c). Таким образом, развитие ЭК в мембранной системе снижает интенсивность РДРВ.

Рис. 4. Распределение концентрации ионов $\mathrm{Na}^+$ (величина показана цветом) и течения раствора (белые линии) в момент времени $t = 140$ c, рассчитанные без учета (a) и с учетом РДРВ (b) для $c_0=1$ моль/м$^3$, $i=1$ А/м$^2$
[Figure 4. The distribution of $\mathrm{Na}^+$ ion concentrations (showed in color) and the flow of the solution (white lines) at the time $t = 140$ s, calculated without considering (a) and with consideration of (b) for $c_0=1$ mol/m$^3$, $i=1$ A/m$^2$]

4. Заключение

Разработана двумерная математическая модель переноса ионов в слое раствора бинарного электролита у поверхности ионообменной мембраны при сверхпредельных постоянных токах с учетом развития ЭК и РДРВ. Модель построена на основе системы связанных уравнений Нернста–Планка–Пуассона–Навье–Стокса и нового гальванодинамического граничного условия для потенциала. Переходное время, рассчитанное с использованием предлагаемой модели, совпадает с высокой точностью (отличие менее 0.03 %) с аналитической оценкой Санда для этой величины. С использованием разработанной модели впервые рассчитаны хронопотенциограммы мембранной системы с учетом влияния электроконвекции и реакции диссоциации/рекомбинации молекул воды. Структура рассчитанных хронопотенциограмм качественно совпадает со структурой экспериментальной кривой.

Показано, что отношение концентрации продуктов диссоциации воды к концентрации ионов соли определяет различное соотношение эффектов электроконвекции и диссоциации.

Рассмотрены следующие варианты соотношения эффектов электроконвекции и диссоциации молекул воды:

1. значимое влияние на массоперенос оказывает электроконвекция, влияние реакции диссоциации воды — незначительные;
2. электроконвекция и диссоциация оказывают значимое влияние на процессы переноса: образование дополнительных носителей заряда в результате диссоциации молекул воды снижает скачок потенциала в слое электролита, что уменьшает интенсивность электроконвекции; развитие электроконвекции, в свою очередь, замедляет процесс диссоциации молекул воды;
3. продукты интенсивной диссоциации молекул воды препятствует развитию электроконвекции.

Конкурирующие интересы. Конкурирующих интересов не имею.
Авторский вклад и ответственность. Я несу полную ответственность за предоставление окончательной версии рукописи в печать. Окончательная версия рукописи мною одобрена.
Финансирование. Исследование выполнено за счет гранта Российского научного фонда № 23–29–00534, https://rscf.ru/project/23-29-00534/.
Благодарность. Автор благодарен проф. М. Х. Уртенову за внимание к данной работе и полезные дискуссии.

Об авторах

Аминат Магометовна Узденова

Автор, ответственный за переписку.
Email: uzd_am@mail.ru
ORCID iD: 0000-0001-5951-9876
SPIN-код: 2810-2267
Scopus Author ID: 55821149800
http://www.mathnet.ru/person188273

кандидат физико-математических наук, доцент; доцент; каф. информатики и вычислительной математики

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

Дополнительные файлы

Доп. файлы

Действие

1. JATS XML

Скачать

2. Рис. 1. Двумерное сечение слоя раствора $\mathrm{NaCl}$ у поверхности катионообменной мембраны (КОМ). Схематически показаны потоки ионов соли $\mathrm{Na}^+$, $\mathrm{Cl}^-$, а также водорода $\mathrm{H}^+$ и гидроксила $\mathrm{OH}^-$, образующихся в расширенной области пространственного заряда (ОПЗ) в сверхпредельном режиме

Скачать (122KB)

Метаданные

3. Рис. 2. Хронопотенциограммы $d_{\varphi}'(t)$ и средние значения толщины электроконвективного вихревого слоя $d_{ec}(t)$, рассчитанные для концентрации $c_0=10$ моль/м$^3$, при плотности тока $i=10$ А/м$^2$ (a), (b); $c_0=1$ моль/м$^3$, $i=1$ А/м$^2$ (c), (d); $c_0=0.1$ моль/м$^3$, $i =0.1$ А/м$^2$ (e), (f). Расчеты с учетом (сплошные линии) и без учета (пунктирные линии) ЭК, с учетом (синие лини) и без учета (красные линии) РДРВ. На рис. (a) также приведена экспериментальная хронопотенциограмма для катионообменной мембраны Neosepta CMX в растворе $\mathrm{CuSO}_4$ 10 моль/м$^3$ (черная линия) при $i=10$ А/м$^2$, полученная в работе [34]

Скачать (524KB)

Метаданные

4. Рис. 3. Отношения концентраций ионов водорода и катиона соли $c_{3av}/c_{1av}$ (a) и концентраций ионов гидроксила и аниона соли $c_{4av}/c_{2av}$ (b). Рис. (c) — увеличение фрагмента рис. (b). Результаты расчета с учетом (сплошные линии) и без учета ЭК (пунктирные линии) и с учетом РДРВ для $c_0=10$ моль/м$^3$, $i=10$ А/м$^2$, $c_0=1$ моль/м$^3$, $i=1$ А/м$^2$, $c_0=0.1$ моль/м$^3$, $i=0.1$ А/м$^2$

Скачать (289KB)

Метаданные

5. Рис. 4. Распределение концентрации ионов $\mathrm{Na}^+$ (величина показана цветом) и течения раствора (белые линии) в момент времени $t = 140$ c, рассчитанные без учета (a) и с учетом РДРВ (b) для $c_0=1$ моль/м$^3$, $i=1$ А/м$^2$

Скачать (1MB)

Метаданные