CC BY
46
УДК 51-71:541.13
МОДЕЛИРОВАНИЕ ПЕРЕНОСА БИНАРНОГО ЭЛЕКТРОЛИТА В КАНАЛЕ ОБЕССОЛИВАНИЯ ЭЛЕКТРОДИАЛИЗНОГО АППАРАТА В ПОТЕНЦИОСТАТИЧЕСКОМ РЕЖИМЕ
Коваленко Анна Владимировна
к.э.н., доцент кафедры прикладной математики
Уртенов Махамет Хусеевич д.ф.-м.н., профессор кафедры прикладной математики
Кубанский государственный университет, Краснодар, Россия
UDC 51-71:541.13
SIMULATION OF TRANSFER OF BINARY ELECTROLYTE IN THE CHANNEL OF DESALTING OF ELECTRO DIALYSIS APPARATUS IN POTENTIAL STATIC MODE
Kovalenko Anna Vladimirovna Cand.Econ.Sci., assistant professor
Urtenov Mahamet Khuseevich, Dr.Sci.Phys.-Math., professor
Kuban State University, Krasnodar, Russia
Хромых Анна Александровна
старший преподаватель кафедры прикладной
математики
Кубанский государственный технологический университет, Краснодар, Россия
Khromikh Anna Aleksandrovna senior teacher of chair of applied mathematics
Kuban State Technological University, Krasnodar, Russia
Статья посвящена численному анализу краевых задач для системы уравнений Нернста-Планка-Пуассона (НПП), применению этих краевых задач к моделированию и изучению массопереноса в канале обессоливания электродиализного аппарата. Предлагаются различные математические модели переноса ионов в потенциостатическом режиме в виде системы квазилинейных уравнений с частными производными. Выявлены основные закономерности возникновения и развития пространственного заряда в канале обессоливания электродиализного аппарата
Article is devoted to the numerical analysis of regional problems for system of the equations of Nernst-Plank-Puasson (NPP), to application of these regional problems to modeling and studying of mass transfer in the channel desalting of electro dialysis device. Various mathematical models of the transfer of ions in potential static mode in the form of system of the quasilinear equations with private derivatives are offered. The basic rules of occurrence and development of a spatial charge in the channel desalting of electro dialysis device are revealed
Ключевые слова: МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ, ПОТЕНЦИО-СТАТИЧЕСКИЙ РЕЖИМ, ЭЛЕКТРО-ДИАЛИЗНЫЙ АППАРАТ, ЭЛЕКТРО-МЕМБРАННЫЕ СИСТЕМЫ
Keywords: MATHEMATICAL MODELING, POTENTIALSTATIC MODE, ELECTRO DIALYSIS DEVICE, ELECTROMEMBRANE SYSTEMS
Введение
Электромебранные системы (электродиализные аппараты,
электромембранные ячейки и т.д.) используют, как правило, в двух разных режимах работы: потенциостатическом режиме, когда постоянным
поддерживается падение потенциала в цепи и гальваностатическом режиме, когда постоянным является ток в цепи. В данной статье
исследуется потенциостатический режим. Изучение массопереноса под действием внешнего электрического поля в электродных и мембранных системах является одной из основных классических фундаментальных задач современной электрохимии [1]. В частности, эта задача имеет важное значение для понимания явлений, протекающих в мембранных системах при реализации интенсивных токовых режимов. Такие режимы являются уже общепринятой практикой проведения электродиализа с ионообменными мембранами при получении высокочистой воды.
Краевые задачи для систем уравнений НПП ставились и изучались многими авторами. Однако в имеющихся работах обычно ограничивались постановкой и рассмотрением одномерных задач [3]. В данной статье рассматривается двумерный случай.
1. Постановка задачи Пусть Н и Ь - ширина и длина камеры обессоливания, соответственно, .х = 0- соответствует условной межфазной границе анионообменная мембрана/раствор, х = Н - соответствует условной межфазной границе катионообменная мембрана/раствор, у = 0 - входу, а у = Ь - выходу из
камеры обессоливания, V - заданная скорость прокачивания раствора. Скорость течения раствора имеет форму параболы Пуазейля или моделируется по другому.
Для математического моделирования явлений переноса для бинарного электролита используются следующие уравнения [1]:
(1)
дг г е0 = — ^ (^1С1 + ^2 С2 )
Шу]г, г = 1,2
(2)
(3)
(4)
Здесь V - градиент, А - оператор Лапласа, V - скорость течения раствора, р0 - характерная плотность раствора, Р - давление, С1, С2 -концентрации катионов и анионов в растворе, соответственно, г2 -зарядовые числа катионов и анионов, Д , 02 - коэффициенты диффузии катионов и анионов, соответственно, ф - потенциал электрического поля, е0 - диэлектрическая проницаемость электролита, — - постоянная Фарадея, Я - газовая постоянная, Т - абсолютная температура, г - время, п - коэффициенты кинематической вязкости
Эта система уравнений удобна для моделирования потенциостатического режима, так содержит уравнение для потенциала (3). После решения соответствующей краевой задачи для системы из 7 уравнений (1)-(3) для семи неизвестных функций ji, С1, г = 1,2, ф плотность тока 1 находится по формуле (4). Напряженность и
потенциал электрического поля связаны соотношением Е = -Vф.
Эти уравнения после некоторых преобразований можно записать в виде, удобном для численного решения:
2. Модель переноса симметричного бинарного электролита
Важным частным случаем является симметричный бинарный электролит. Когда г1 = -г2. Не нарушая общности рассуждения в этом случае можно положить г1 = 1 и г2 = -1, тогда уравнения (5), (6) упрощаются и принимают вид:
(5)
(6)
(7)
Т—'
—■ = 02АС2 - ёы(с^)- —О2сИуСVф) дг КТ
£„ Др = - Р (С, - С 2 )
(9)
3. Начальные и краевые условия
К системе декомпозиционных уравнений должны быть добавлены соответствующие краевые условия, которые ставятся в зависимости от целей конкретного исследования процессов переноса в электродиализных аппаратах, от режима их эксплуатации гальваностатического или потенциостатического, и т.д. В этой работе мы будем рассматривать потенциостатический режим, которому соответствует условие:
означающее, что величина падения потенциала в канале постоянна. Поскольку важно только величина падения потенциала, то, не нарушая общности можно положить, например ф(г,0,у) = Сф, тогда ф(г,Н, у) = 0.
Каждому значению падения потенциала сф соответствует своя плотность
тока
являющейся функцией от времени, но не зависящая от х.
Все граничные условия должны быть согласованы со свойствами мембран и с величиной йр («интенсивным» или «мягким» токовым
режимом).
Наряду с условием (7) будем использовать следующие граничные условия:
1). На поверхности анионообменной мембраны х = 0, у є[0, Ь], ї > 0 будем считать граничную концентрацию анионов равной фиксированному заряду внутри мембраны:
Кроме того, предположим анионообменную мембрану идеально селективной, т.е. непроницаемую для катионов:
р(ї, Н, у) - р(ї, 0, у) = йр = еоп8ї,
(10)
0
(11)
С2(ї,0,у) = С
(12)
ЭС +^-с, Эр
(1,0, у) = 0. (13)
Эх ЯТ 1 1 Эх ^
2). На поверхности катионообменной мембраны х = Н, у є [0, Ь], ї > 0 будем считать граничную концентрацию катионов равной фиксированному заряду внутри мембраны:
Сі(ї,Н,у) = Сы (14)
Кроме того, предположим катионообменную мембрану идеально селективной, то есть непроницаемой для анионов:
ЭС2 Р Эр
—-+—
(ї,Н,у) = 0. (15)
Эх ЯТ Эх _
3). На входе в рассматриваемую область у = 0, хє [0, Н], ї > 0 будем полагать:
С,(ї,х,0) = Сі0(ї,х), 1 = 1,2, (16)
= 0, (17)
Эу
4). На выходе из рассматриваемой области у = Ь, х є[0, Н], ї > 0 будем использовать «мягкие» условия на концентрации и на потенциал:
ЭС1 (ї, х, Ь)
эу
Эр(?, х, Ь)
: 0,1 = 1,2. (18)
0, (19)
эу
5). Начальные условия при ї = 0 примем согласованными, по возможности, с остальными граничными условиями:
С1 (0, х, у) = С, ,п (х, у), 1 = 1,2. (20)
Для функции р в качестве начального условия можно взять, например, линейную функцию, независящую от у :
( х Л
р(0,х,у) = 1 - - йр (21)
V Н)
Для решения краевой задачи применяется метод конечных элементов [2].
4. Характерные данные
Нами были проведены численные эксперименты для раствора КаС1 в широком спектре таких параметров, как начальная концентрация, межмембранное расстояние, длина канала, значение тока в цепи йр, и
определены основные закономерности распределения электрохимических (концентрация, напряженность и потенциал электрического поля, и т.д.) полей. Ниже представлены некоторые результаты численных
экспериментов при следующих входных параметрах: ширина канала обессоливания Н = 1 мм или меньше, длина канала Ь = 7 мм, средняя скорость вынужденного течения раствора варьировалась от -0 = 2 10-6 м/с до У0 = 10-2м/с, а падение потенциала от й^= 0.2В до йр = 0.6В, начальная концентрация раствора принимает значения от С0 = 10 моль/м до С0 = 100 моль/м3 , температура раствора Т = 293 К, начальная плотность раствора р 0 = 1002.5 кг/м, коэффициент кинематической вязкости п = 1006 • 10-9 м2/с, коэффициент диффузии катиона и аниона, соответственно, Д = 1.33 • 10-9м2/с и В2 = 2.05 • 10-9м2/с.
5. Переход к безразмерному виду и оценка безразмерных параметров
Для перехода к безразмерному виду используем «естественные» безразмерные величины (индекс и означает, что соответствующая величина безразмерная):
у(и) = Х (и) = у С (и) = С С (и) = Скт С (и) = Сат V (и) =
тт 5 У тт 5 1 5 кт ^ 5 ат ^ 5 тт 5
Н Н С0 С0 С0 -0
1(и) = Н 1, Ё(и) = Н—Ё, р(и) = р, Д(и) = Д-, где, например, г Д0С0 1 ЯТ ЯТ 1 Д0
Д0 = ц,или Д0 = 10-9 м2/с, I(и) = Н I, г(и) = ^г, Ь(и) = Ь.
0 1 0 Д0С0 Я Н Н
Кроме того, введем в рассмотрение безразмерные параметры:
dрP) = -—dр - безразмерное падение потенциала, Ре = —0Н - число Пекле,
£
ЯТ р
ЯТ£0 Н 2С0—;
2
V
Н
- Дебаевская длина [1].
Эти безразмерные параметры имеют следующие характерные
—
пределы изменения: у е [0,1], у е [0,7], d
(и)
^р = kdр меняется от 5 до
ят р р
25, число пекле линейно зависит от начальной скорости и при изменении меняется от 1 до 104, £ меняется от 10-4 до 10-10 в
зависимости от значений -0, С0, Н. Таким образом, £ можно считать малым параметром. В безразмерных переменных исходная система уравнений (1-4) имеет вид (индекс «и» для упрощения записи опущен):
У, = -2 Д гС г Vр - Д V С1 + РеС—, 1 = 1,2
дС,
Ре-
дг
^1у ]г, 1 = 1,2
£Др ( 21С1 + 2 2 С 2 )
I = ^ + 2 21г
Соответственно изменяются и системы уравнений (5-6) и (7-9), например система уравнений (5) примет следующий безразмерный вид:
Ре —^ = Д ДС, + 2 Д^С, Vр)-Pediv(C1V) дг
1 = 1,2
(21)
£ Др (21С1 + 22С2 ) (22)
Для системы уравнений (21), (22) безразмерные краевые условия примут вид:
1). При у = 0 :
р(г,0,у) = dр; С2(г,0,у) = Са
дС
+ 2,С, Р
ду ду
(г,0,у) = 0.
(23)
2
2). При х = 1:
<р(г,1,у) = 0; С/Г,!,у) = Ст; (+ z1C, ^\г,1,у) = 0. (24)
у Эх Эх у
3). При у = 0:
х,0) = 0, С (г, х, 0) = с 0 (г, х), , = 1,2, (25)
Эу ,
4). При у = Ь:
= 0, ОС,(г,х,Ь) = 0, , = 1,2. (26)
Эу Эу
5) При г = 0
(р(0,х,у) = (1 - х)dj , с, (0, х, у) = с,п(х, у), I = 1Д. (27)
Наличие малого параметра не позволяет решать краевую задачу (21)-(27) в исходных переменных при е меньше 10-2, поскольку в области погранслоя необходимо деление хотя бы с шагом порядка е/10.
Например, если е имеет порядок 10-4, то шаг должен быть порядка от 10-5 и меньше, количество элементов при этом будет порядка от 1010 и больше. Попытки вне погранслоя увеличить шаг приводит к «взрыву погрешностей» [4].
Для решения краевой задачи прие меньших порядков чем 10-2,
введем специальную растянутую систему координат и переменных:
-1 -1 -1 -1 -2
£ = хе 3, £ =е 3, д = уе 3, V = Ье 3, т = ге 3, (р = Ф,
1 11 —► —► — +— —► —► —
Е = Ое 3, С = 5.е 3, Ж = Уе 3.
’ II’
Тогда
ЭС, -3 Э5, ЭС, Э5, ЭС, Э5,
______е з____________________, ________,
Эг Эт Эх Э£ Эу Эд
Э2С, -3 Э25, Э2С, -3 Э25,
----~~е 3--------тт , ----= е 3---------
Эх2 ' Э£2 ’ Эу2 Эд2
Следовательно ДСг = е 3 Д^., йю (С р) = е 3 Ф),
_1 _2
Шу(СV) = е 3йу(БгЖ), Ар = е 3 ДФ .
В этих формулах дифференциальные операторы слева считаются в
переменных (х, у), а справа,- в (X, д).
Подставляя эти формулы в уравнения (21), (22) в растянутых
переменных получим уравнения:
Отметим, что внешний вид уравнений не изменился, за исключением того, что в новые уравнения малый параметр уже не входит. Область изменения переменных растянулась, т.к., например, при е = 10-3 прямоугольник [0,1]х[0,7] становится прямоугольником [0,10]х[0,70], и поэтому количество узлов увеличивается (в 100-210 раз), однако это значительно меньше начального варианта и уже вполне приемлемо. Растянутые переменные позволяют проводить расчеты вплоть до е = 10-6.
Краевые условия в растянутых переменных имеют вид:
*\п
Ре = ДАБ. + іД1сііу(Б1Ф)-РеШф^), і = 1,2
(28)
(29)
1) При Х = 0:
(30)
(31)
3) При д = 0:
°д
4) При д = ду.
ЭФ(Т,Х,^1) = 0, (Т,Х,^1) = 0, I = 1,2. (33)
^ дд ( )
5) При т = 0
1 _1 1 1 Ф(0,£,?) = еЦйр, 8,(0,£,£) = е~3С]„(еТ, е3Х), > = 1,2, (34)
6. Основные закономерности переноса ионов в канале обессоливания
6.1. Возникновение и развитие пространственного заряда. Структура канала обессоливания.
а) б)
Рисунок 1. Г рафик плотности распределения заряда р = 81 + 2 2 8 2 в
безразмерном виде при е = 10 _3 а) йр = 8, б) йр = 30
Постепенно увеличивая йр, получаем, пространственный заряд
возникает около мембран при запредельном режиме, причем при
дальнейшем увеличении падении потенциала область пространственного заряда увеличивается (рис.1). Кроме того, возле непосредственно вблизи мембран возникает погранслой, обусловленный влиянием граничного условия. Это влияние вне погранслоя практически отсутствует.
Из рисунка видно, что в ядре потока, т.е. в центральной части канала обессоливания электронейтральность раствора сохраняется с большой точностью (р » 0), а в области пространственного заряда процесс практически стационарен.
В зависимости от начального значения возможно возникновение и начального погранслоя.
Таким образом, структура канала обессоливания может быть достаточно сложна и канал может делиться на области электронейтральности, пространственного заряда, промежуточные слои и пограничные слои (диффузионный, диффузный, начальный).
6.2. Соленоидальность плотности тока
Численные расчеты показывают, что Шіу I » 0 в широких пределах
изменения Шр, например, Шіуі »10-3 при е = 10-3, Шіуі »10-6 при
ґ —►
е = 10- , поэтому можно ожидать Шіуі = О(е), е ® +0.
Таким образом, можно считать, что во всем канале Шіуі » 0, т.е. поле плотности тока является соленоидальной. Для объяснения этого интересного результата обратимся к уравнению (2):
дС ^
---- = -шу] , і = 1,2
ді Л
Умножим каждое уравнение на г. и сложим, тогда:
Шіуі = -
ді
В области электронейтральности с большой точностью выполняется условие электронейтральности г1С1 + г2 С2 = 0 , следовательно, с такой же
точностью Шіуі = 0. С другой стороны в области пространственного заряда из-за стационарности процессов переноса имеем
д(г\С\ + 2 2 С2) Л Л- 1 Л
-----11---2—— = 0, следовательно, опять ш уі = 0.
ді
Таким образом, поле плотности тока с достаточной точностью можно считать соленоидальным во всем канале обессоливания, хотя это выполняется в разных участках канала по различным причинам.
6.3. Основные закономерности распределения концентрации ионов
а)
б)
Рисунок 2 - Графики S1, S2 при t = 100, d = 16, e = 10"
а)
б)
в)
г)
3
Рисунок 3 - Графики 81, 82 а), б) - при ? = 464, = 16, £ = 10-4, в), г) -
при ? = 10000, 3,р= 25, £ = 10-6
6с 10с 20с
Рисунок 4 - Формирование начального погранслоя. Графики С1 в размерном виде в разные моменты времени в секундах при dj = 0.6В
Расчеты показывают, что концентрации в области
электронейтральности практически меняются по линейному закону, а в области пространственного заряда малы, причем одна из концентраций экспоненциально мала, а вторая почти постоянна. В области пространственного заряда устанавливается стационарный режим (рис.2-4 ).
С уменьшением £ изменения становятся более резкими, погранслои ярко выраженными. Процесс быстро стабилизуется, начальный погранслой
формируется за короткое время, зависящее от е, в ядре потока появляется «плато» с практически постоянным значением концентрации.
Заключение
1. Построена математическая модель нестационарного переноса бинарного электролита в разбавленных растворах в электромембранных системах с учетом пространственного заряда в потенциостатическом режиме, которая представляет собой краевую задачу для системы квазилинейных уравнений в частных производных.
2. Проведен численный анализ и установлены основные закономерности переноса в каналах обессоливания электромембранных систем в условиях совместного пространственного заряда и вынужденной конвекции. Изучены нестационарные и переходные процессы.
Работа выполнена при финансовой поддержке Федерального Агентства по образованию и науке РФ в рамках темы 1.4.08 Единого заказ/наряда.
Литература
1. Ньюмен Дж. Электрохимические системы. М.: Мир, 1977. 463 с.
2. Норри Д., де Фриз Ж. Введение в метод конечных элементов: Пер. с англ. М.: Мир, 1981. 304 с.
3. Уртенов, М.Х. Краевые задачи для систем уравнений Нернста-Планка-Пуассона. Краснодар, КубГУ, 1998. - 126 с.
4. Дулан Э., Миллер Дж., Шилдерс У. Равномерные численные методы решения задач с пограничным слоем. М.: Мир, 1983. 254 С.
