Journal of Siberian Federal University. Engineering & Technologies, 2017, 10(1), 59-73
УДК 669.713.7:004.9
Three-Dimensional Mathematical Modeling of Dynamics Interfaces Between Aluminum, Electrolytes and Reverse Zone
of Oxidized Metal Depending on the Potencial Distribution
Tatyana V. Piskazhovaa, Nadezhda P. Savenkovab, Sergey V. Anpilovb, Alexey V. Kalmykovb*, Fedor S. Zaitsevb and Fedor A. Anikeevb
aSiberian Federal University 79 Svobodny, Krasnoyarsk, 660041 Russia bLomonosov Moscow State University 1 Leninskie gory corp., Moscow, 119991, Russia
Received 03.11.2015, received in revised form 15.12.2015, accepted 23.05.2016
The mathematical model with a high detailed description of the studied processes is presented in current paper. Results of modeling Soderbergh's reduction cell for model task and for reduction cell with multiple anodes are also presented. System of the equations of Navier-Stokes is used for modeling of hydrodynamics of process of electrolysis. Distribution of electromagnetic fields is fitted to Maxwell's system of equations. Influence of distribution of electric potential over the anode on MHD-stability of process is considered and comparative analysis of numerical experiments is also given.
Keywords: electrolysis of aluminum, magnetohydrodynamics, MHD-stability.
Citation: Piskazhova T.V., Savenkova N.P., Anpilov S.V., Kalmykov A.V., Zaitsev F.S., Anikeev F.A. Three-dimensional mathematical modeling of dynamics interfaces between aluminum, electrolytes and reverse zone of oxidized metal depending on the potencial distribution, J. Sib. Fed. Univ. Eng. technol., 2017, 10(1), 59-73. DOI: 10.17516/1999-494X-2017-10-1-59-73.
© Siberian Federal University. All rights reserved Corresponding author E-mail address: [email protected]
*
Трёхмерное математическое моделирование динамики границы раздела сред алюминия, электролита и зоны обратного окисления металла в зависимости от распределения потенциала
Т.В. Пискажоваа, Н.П. Савенкова5, С.В. Анпилов5, А.В. Калмыков5, Ф.С. Зайцев5, Ф.А. Аникеев5
аСибирский федеральный университет Россия, 660041, Красноярск, пр. Свободный, 79 бМосковский государственный университет
им. М.В. Ломоносова Россия, 119991, Москва, Ленинские горы, корп. 1
В статье представлена новая математическая модель с высокой степенью детализации описания изучаемых процессов. Проведены расчёты электролизёра Содерберга для модельной задачи, а также расчёты для многоанодного электролизёра. Для описания гидродинамики процесса электролиза используется система уравнений Навье-Стокса. Распределение электромагнитных полей описывается системой уравнений Максвелла. Рассматривается влияние распределения потенциала по аноду на МГД-стабильность электролиза, а также даётся сравнительный анализ численных экспериментов.
Ключевые слова: электролиз алюминия, магнитогидродинамика, МГД-устойчивость.
Большинство исследований, посвященных электролизу алюминия, так или иначе сталкиваются с проблемой МГД-устойчивости процесса [1, 2]. Особый интерес представляет определение параметров процесса, при которых достигается устойчивая работа электролизёра. Существует множество критериев МГД-устойчивости, например критерии Селе, Бояревич-Ромерио. Однако каждый из критериев имеет свои недостатки, что существенно ограничивает возможность их практического применения. В настоящее время имеется большое количество математических моделей, описывающих динамику основных процессов в электролизной ванне. Однако обычно модели не имеют высокой степени адекватности и не позволяют с необходимой точностью самосогласованно описывать динамику физико-химических процессов, а также не способны работать в реальном времени.
В статье описана новая математическая модель с высокой степенью детализации описания изучаемых процессов. Проведены расчёты электролизёра Содерберга для модельной задачи, а также расчёты для многоанодного электролизёра. Для описания гидродинамики процесса электролиза используется система уравнений Навье-Стокса. Распределение электромагнитных полей описывается системой уравнений Максвелла. Для давления в электролизной ванне предложена отдельная модель, полученная на основе уравнений неразрывности и Навье-Стокса.
Основным преимуществом предлагаемого подхода является возможность изучения распределения поля скоростей, электромагнитных полей, температуры и линий тока как в горизонтальной, так и в вертикальной плоскостях. Данный подход позволяет детально иссле-
довать динамику границы раздела сред алюминий - электролит и границы зоны обратного окисления металла, что очень важно для определения эволюции во времени критического значения МПР.
Возникновение вертикальных вихревых структур в процессе промышленного электролиза может привести к существенным колебаниям зеркала металла, а в отдельных случаях и к перемешиванию фаз, снижению выхода первичного алюминия. В этой связи большой интерес представляет выявление зависимости возникающих вертикальных скоростей в металле от управляющих параметров процесса. В данной работе рассматривается влияние распределения потенциала по аноду на МГД-стабильность электролиза, а также даётся сравнительный анализ численных экспериментов.
Математическая модель
Рассмотрим подробнее трёхмерную трёхфазную модель алюминиевого электролизёра, основанную на системе уравнений Навье-Стокса и системе уравнений Максвелла. Как уже отмечалось, в основе предлагаемой модели лежит многофазный подход, который предполагает, что в каждом малом элементарном объёме присутствует и алюминий, и электролит, и газ. При этом их смесь занимает объём целиком, а каждая из фаз - некоторую часть этого объёма - АУ1 для алюминия, АУ2 для электролита и АУ3 для газа. В каждой точке объёма, занятого смесью, вводится макроскопическая скорость компоненты смеси V т, давление Р и объёмная доля фазы (объёмное содержание фазы) ат:
АУ, АУ2 АУъ а1 " АУ " АУ ,(Хъ ~ АУ '
при этом а! + а2 + а3 = 1.
Здесь и далее переменные с индексом т = 1 относятся к среде алюминия, с индексом т = 2 - к среде электролита, с индексом т = 3 - к среде газа.
Уравнение сохранения т-й компоненты смеси имеет вид
д(ХтРт- + ¿Ы(атРтО = Мт ,
dt
где рт - плотность, Мт - массовая скорость превращения соответствующей фазы. Уравнения импульсов каждой составляющей можно представить в виде
р v —- — —- —- — —- —
m m- + (v , V)(a p v ) = -a Vp + a Av +a p g + P +a F
v m' /v m> m m' m * ' m m m> m <-> m m n
dt
где - вязкости среды; Fm - объёмная плотность силы, обусловленная электромагнитным полем (сила Лоренца):
F~=— Г rot Н х Н
m
Pm - объёмная плотность силы трения (Стоксовой силы) между компонентами смеси за счёт вязких сил, которая определяется из соотношений (пример для 1-й фазы):
- 61 -
p = axa2K2 (Vj - V2) + ахаъКъ (v - V3),
где
К2 = + а2ц2, К3 = + аъцъ.
Для уравнения Навье-Стокса на границе области поставим условия прилипания: = 0.
Для поиска температуры расчётная область разбивается на две части: область жидких фаз и область твёрдых бортов ванны. Таким образом, получаем систему для поиска температуры в виде [3]
Z (amPm)^- + div T Z (v™)] = div(Лgrad(T)) + f - в области
i=1,2,3 Ot ^ m=1,2,3 J
жидких фаз,
dT
= {кgrad (Т)) - в области гарнисажа и настыли,
здесь с - теплоёмкость вещества; к - коэффициент температуропроводности; X - коэффициент теплопроводности; / - слагаемое, отвечающее за нагрев смеси из-за протекающего тока, определяется из закона Джоуля-Ленца:
Г = 1 • Е.
Граничные условия для уравнения на внешней границе теплоизоляции Т = 20 °С. Для моделирования давления разработана отдельная математическая модель
div
Р2
с граничными условиями
aj- Vp + ^ Vp + ^ Vp
Рз
= div
Z((v ,V)(a v ) + ц Av +a F ) + g
1 \ m * /\ mm/ r^m m mm I о
у m=1,2,3
dp dn
( _ _ _ _ __ _ : Z ((^ , V)K vm ) + ^Avm + Fm ) + g
Уравнение для потенциала электрического поля выводится из уравнения для дивергенции электрического поля и его связи с потенциалом. Это уравнение имеет вид
0,
где правая часть равна нулю в силу предположения об электрической нейтральности среды.
Потенциал на аноде считается равным некоей наперёд заданной функции Ф(х,у), на катоде - нулевым. В непроводящей части границы ванны примем нормальную производную потенциала равной нулю. В результате получим следующие граничные условия:
=°(x,y),H,am0a = 0
дф dn
= 0.
С учётом индуцированного электрического поля, вызванного движением среды, общее электрическое поле таково:
Вихревым электрическим полем, которое возникает при изменении магнитной индукции, пренебрегаем.
По напряжённости электрического поля можно вычислить плотность текущего в электролизной ванне тока
Магнитное поле является соленоидальным векторным полем и, следовательно, имеет векторный потенциал А. Для нахождения А воспользуемся калибровкой Кулона
Принимая во внимание теорему о циркуляции магнитного поля, пренебрегая током смещения и считая изменение электрической индукции по времени несущественной, можно получить соотношение для определения векторного потенциала магнитного поля
Численное решение трёхмерной системы уравнений в частных производных проводилось разностным методом и относительно новым для данной области методом сглаженных частиц (8РИ) [4, 5]. Метод 8РИ обладает рядом достоинств, которые делают его привлекательным для детального описания процесса электролитического получения алюминия. К ним относятся: отсутствие сетки по фазовым переменным, сведение задачи к решению системы обыкновенных дифференциальных уравнений, простое и естественное описание разделов многокомпонентных сред, источников и стоков, границ сложной геометрии, узколокализо-ванных движущихся объектов, например пузырей, эффективная реализация на параллельных компьютерах.
При численном решении системы уравнений разностным методом применяется метод разделения по физическим процессам. При таком подходе на очередном шаге по времени каждое уравнение решается отдельно, а их взаимосвязь учитывается на следующем шаге. В методе 8РИ разделения по физическим процессам с характерными временами одного масштаба не требуется и все уравнения для всех фаз решаются одновременно.
В представленном численном эксперименте приводится попытка установить зависимость скорости движения фаз от распределения потенциала по аноду. Пусть в начальный момент времени скорости обеих фаз равняются нулю, граница разделения фаз является плоской, без каких-либо возмущений. Изменение потенциала по аноду Ф(х,у) во всех расчётах происходит по оси ОУ, вдоль оси ОХ потенциал остаётся неизменным.
j = аЕ.
divA = 0.
И0 j = rot rotA = grad divA - ДА = -ДА .
Для векторного потенциала используется следующее граничное условие:
Оси ориентированы следующим образом:
• ОХ - горизонтально, вдоль короткой стороны ванны,
• ОУ - горизонтально, вдоль длинной стороны ванны,
• 02 - вертикально.
На приведённых ниже рисунках изображены результирующие поля скоростей в различных срезах. Граница раздела фаз представлена в срезе плоскостью У2, так как по оси ОХ профиль границы меняется несущественно.
Равномерное распределение потенциала по аноду
В данном расчёте потенциал электрического заряда берётся равномерно распределённым по всей поверхности анода, функция Ф(х,у) является постоянной.
На рис. 1 присутствуют ярко выраженные вихревые структуры как в области электролита, так и в области металла. Данные вихри расположены у краёв ванны, в тех местах, где заканчивается анод, именно наличие горизонтальных токов влияет на возникновение вертикальных скоростей. Они могут привести к колебанию границы раздела фаз у бортов электролизной ванны.
На рис. 2 ясно просматривается искривление первоначально горизонтальной границы раздела сред металл-электролит. Новая граница раздела окрашена зелёным цветом. Наблюдается поднятие границы раздела в центре ванны и ее оседание у бортов.
Линейное распределение потенциала по аноду
В данном расчёте потенциал электрического заряда линейно меняется от нулевого значения до максимума вдоль оси ОУ (т. е. на рис. 3 слева - нулевое значение потенциала, справа -максимальное).
На рис. 3 присутствует вертикальный вихрь, смещённый в правую сторону. Наиболее интенсивное движение в вертикальном направлении приходится на правую границу анода, где сосредоточен максимум потенциала. В данном эксперименте видно, что плавное изменение распределения потенциала по аноду ведёт к небольшому изменению вертикальной скорости, в то время как резкое изменение потенциала у правого борта ванны вызывает резкий рост градиента скорости по оси 02. Вертикальные срезы, смещённые относительно середины ванны,
Г" ■ ( 1 ' t ( : 1 t i ' t -1 V J f I * f t t t 1 I t t 1 t t 1 < 111 11, i 11 t 11 111 111 111 Г t 1 11 11 f t f t f t 11 11 1111 1111 1111 1111 1111 1111 11 11 11 11 11 11 11 11 11 t t t t t t ' ' ! I f'"T^ t t t i ■ , 1 ■t t t J ■ , 1 ' 1 1 ' i t t t ■ - , . , t t , . , ( V \ f 1 ' . , l i \ t I ' \ J Я
111 ::: : ■
i ; . : :;; 10 * -
2 4 6 8
У
Рис. 1. Срез поля скоростей плоскостью У2, плоскость проходит через середину электролизной ванны
Рис. 2. Колебание границы раздела сред: а - профиль; б - вид сверху и вид снизу
............„„„^^s/. .............-.----■.■.'"■'Г. 1
- ■ ' 1 - - - >'/ i ; : : : :..................' / / - ................ ...........— ........... * ' ' ' '•-/'•
|_J-1-1-1 i _|-i 1-1-1-1-1 —i— Г 1
2 4 6 8
Y
Рис. 3. Срез поля скоростей плоскостью У7, плоскость проходит через середину электролизной ванны
демонстрируют качественно схожие с рис. 3 картины, однако модули скоростей на них существенно меньше и практически затухают у бортов ванны.
На рис. 4 зелёным цветом обозначено отклонение границы раздела от горизонтального положения. В отличие от предыдущего численного эксперимента здесь мы можем наблюдать проседание границы раздела фаз только у одного борта, того, где потенциал на аноде имеет максимальное значение.
б
Рис. 4. Колебание границы раздела сред: а - профиль; б - вид сверху и вид снизу
«Параболическое» распределение потенциала
Рассмотрим численный эксперимент, в котором распределение потенциала электрического тока имеет параболическую зависимость: потенциал электрического заряда имеет максимальное значение в центре оси ОУ, к краям он уменьшается до нуля (т. е. на рис. 5 в центре максимальное значение, к бортам - нулевое).
На рис. 5 представлены две ярко выраженные вихревые структуры противоположного направления. Область наиболее интенсивного движения по оси 02 приходится на максимальное значение потенциала на аноде. Как и в предыдущих расчетах, вертикальное движение сохраняется по всему объёму и затухает к бортам ванны.
На рис. 6 зелёным и желтым цветами изображена граница раздела жидких фаз в смеси. В данном эксперименте граница наиболее сильно просела в центре ванны, в том месте, где потенциал на аноде максимален.
Распределение потенциала с вырезом центральной линии
В заключительном численном эксперименте потенциал электрического заряда распределён равномерно по всему аноду кроме полоски шириной в 26 см в центре, где потенциал равен нулю.
Интенсивность вертикальных скоростей на рис. 7 вызвана отсутствием потенциала в середине анода. Как и в предшествующих расчётах, вертикальные вихревые структуры могут
Рис. 5. Срез поля скоростей плоскостью У7, плоскость походит через середину электролизной ванны
Рис. 6. Колебание границы раздела сред: а - профиль; б - вид сверху и вид снизу
вызвать МГД-нестабильность. Ближе к бортам ванны интенсивность скоростей снижается в связи с трением о боковые стенки.
Представленные на рис. 7 распределения скоростей в горизонтальных срезах в средах электролита и металла похожи на поля скоростей рис. 5, однако имеют обратные направления и меньшую интенсивность.
Данное распределение потенциала по аноду влечёт наиболее сильные возмущения границы раздела сред алюминия и электролита. На рис. 8 зелёным и желтым цветами обозначено
Рис. 7. Срез поля скоростей плоскостью У2, плоскость проходит через середину электролизной ванны
отклонение границы раздела от горизонтального положения. В данном расчёте наблюдается ярко выраженное поднятие границы раздела в центре расчётной области, а также менее интенсивные колебания к бортам ванны.
Численное моделирование многоанодного электролизёра
При моделировании промышленного электролизёра с 22 анодами были взяты реальные начальные распределения поля скоростей и границы раздела сред алюминий-электролит.
Основной вынуждающей силой является сила Лоренца, которая зависит от плотности тока в электролизной ванне и распределения магнитного поля. На рис. 9 представлено распределение плотности тока в плоскости У2, взятое в середине ванны. С увеличением глубины можно наблюдать более равномерное распределение плотности тока вдоль оси ОУ, в то время как непосредственно под анодами плотность тока распределена неравномерно.
В связи с этим магнитное поле (рис. 10, левый) в алюминии имеет одну большую вихревую структуру в плоскости ХУ. Для среды электролита (рис. 10, правый) более характерны маленькие вихревые структуры, расположенные непосредственно под анодом.
На рис. 11 представлено поле скоростей в плоскости У2 в середине электролизной ванны. Данное распределение похоже на поле скоростей при равномерном распределении потенциала по аноду. Однако наблюдается неравномерность вертикальных компонент скоростей, которая обусловлена многоанодностью.
Рис. 9. Распределение плотности тока в плоскости У2 в середине ванны
В результате при численном моделировании наблюдается возникновение ряби на поверхности раздела алюминий-электролит, количество волн равно 10 и соответствует промежуткам между анодами (рис. 12).
Моделирование зоны обратного окисления позволяет точнее проследить динамику МГД устойчивости процесса электролиза. На рис. 13 представлена зона обратного окисления с объёмной долей газа в расплаве 0,2 и 0,05 соответственно.
Как видно из рисунков, наибольшую концентрацию газ имеет между анодами, однако в меньшем количестве он встречается и в нижних слоях электролита.
Выемка крайней пары анодов
Одним из технологических процессов, который может привести к возникновению МГД-нестабильности, является выемка крайней пары выгоревших анодов. При замене анодов происходит перераспределение тока, протекающего в ванне, и, как следствие, изменение электромагнитных полей. На рис. 14 представлена плотность тока после перераспределения (убрана крайняя правая пара анодов).
2 4 6 8
У
Рис. 11. Срез поля скоростей плоскостью У7, плоскость проходит через середину электролизной ванны
Рис. 13. Зона обратного окисления с объёмной долей газа в расплаве 0,2 и 0,05
Рис. 14. Распределение плотности тока в плоскости У7 в середине ванны
Магнитное поле существенно меняется после выемки крайних анодов (рис. 15), особенно в расплавленном металле. Данное изменение влечёт за собой распределение силы Лоренца и возникновение волны на поверхности раздела фаз алюминий-электролит.
В новом распределении поля скоростей (рис. 16) появляются компоненты горизонтальных составляющих под местом выемки анодов, это ведёт к возникновению волны и последующей МГД-нестабильности.
На рис. 17 представлены границы раздела фаз алюминий-электролит при замене выгоревшей пары анодов. Отсутствие ряби по всей поверхности границы в дальнейшем вызовет волну, направленную в положительном направлении оси ОУ. Данная волна вызовет нестабильную работу электролизёра.
Наибольшая концентрация пузырьков, как и в предыдущем численном эксперименте, располагается вокруг анодов, однако в месте изъятия анодов граница зоны обратного окисления проседает существенно меньше (рис. 18).
Разработанное программное обеспечение позволяет наглядно визуализировать динамические процессы внутри электролизной ванны.
Показано, что распределение потенциала по аноду сильно влияет на вертикальную компоненту поля скоростей. Неизотропность заряда вызывает отклонение границы раздела фаз
2 4 6 8
У
Рис. 16. Срез поля скоростей плоскостью У7, плоскость проходит через середину электролизной ванны
Рис. 17. Граница раздела фаз алюминий-электролит в момент времени 1=0 с и 1=2 с.
Рис. 18. Зона обратного окисления с объёмной долей газа в расплаве 0,2 и 0,05
от устойчивого положения, приводит к возникновению волн на поверхности раздела сред и в отдельных случаях может вызвать перемешивание сред. В данных численных экспериментах наблюдается следующая зависимость:
• при увеличении плотности заряда на аноде граница раздела фаз продавливается вниз;
• при уменьшении плотности заряда на аноде граница раздела фаз выгибается вверх;
• при равномерном распределении заряда по аноду граница менее всего отклоняется от положения равновесия.
В силу того что невозможно обеспечить равномерное распределение тока по аноду в электролизёре Содерберга, одноанодный электролизёр является менее устойчивым к изменению потенциала, чем многоанодный электролизёр, что подтверждается численными расчётами. Анализ проведённых численных экспериментов позволяет предположить, что наиболее оптимальным будет электролизёр с большим количеством анодов в виде штырей с одинаковыми потенциалами и с минимальным расстоянием между ними.
Список литературы
[1] Савенкова Н.П., Анпилов С.В., Кузьмин Р.Н. и др. Цветные металлы. Красноярск, 2011, с. 282-286 [Savenkova N.P., Anpilov S.V., Kuzmin R.N. et.al. Non-ferrous metals. Krasnoyarsk, 2011, с. 282-286 (in Russian)]
[2] Нигматулин Р.И. Основы механики гетерогенных сред. М.: Наука, 1978 [Nigmatulin R.I. Foundations of mechanics of heterogeneous media. Moskow, Nauka, 1978 (in Russian)]
[3] Белолипецкий В.М., Пискажова Т.В. Математическое моделирование процесса электролитического получения алюминия. Решение задач управления технологией. Красноярск: СФУ, 2013, 271 с. [Belolipetskii V.M., Piskazhova T.V. Mathematical modeling of the electrolytic obtaining of aluminum. The task of management technology. Krasnoyarsk, SibFU, 2013, 271 p. (in Russian)]
[4] Monaghan J. J. Rep. Prog. Phys, 2005, 68, 1703-1759.
[5] Zaitsev F.S. Mathematical modeling of toroidal plasma evolution. English edition. Moscow, MAKS Press, 2014, 688 p.