Математическая модель и численная схема расчёта электрических полей в гальванических ваннах с плоским токонепроводящим экраном Текст научной статьи по специальности «Математика»

ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ И

ПРОЦЕССЫ УПРАВЛЕНИЯ N. 3, 2021 Электронный .журнал, рег. Эл. N ФС77-394Ю от 15.04.2010 ISSN 1817-2172

http ://diffjo urn al. яр Ь и. 1 "и/ e-mail: jodiff@mail ru

Прикладные задачи

Математическая модель и численная схема расчёта электрических полей в гальванических ваннах с плоским токонепроводящим экраном

И. Ю. Пчелинцева, Ю. В. Литовка Тамбовский государственный технический университет,

e-mail: [email protected], [email protected]

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

Ключевые слова: уравнение Лапласа, метод Ньютона, критерий неравномерности, толщина покрытия.

1 Введение

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

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

Чтобы решить задачу вычисления толщины покрытия на детали-катоде, необходимо рассчитать распределение потенциалов в гальванической ванне из уравнения Лапласа.

Одним из наилучших является метод, согласно которому пространство гальванической ванны разбивается сеткой, и производные функции распределения потенциала электролита в объеме ванны заменяются их разностными аналогами. Полученная система нелинейных алгебраических уравнений решается методом простых итераций с использованием метода верхней релаксации с прогонкой по строке |1,2] по анодным и катодным плотностям тока. Однако у такой итеративной процедуры есть недостаток медленная сходимость к приближённому решению. В то же время, при решении задачи нанесения более равномерного покрытия из-за её высокой размерности уравнение Лапласа приходится решать тысячи раз. Если задать высокую точность, то такую задачу нельзя решить за приемлемое время. В связи с этим задача повышения скорости вычислений при решении уравнения Лапласа является актуальной.

Анализ итерационных методов показал, что вместо простых итераций целесообразно применить метод Ньютона, который, как известно 3], обладает квадратичной скоростью сходимости.

Цель данной работы разработка математической модели процесса получения покрытия с бесконечно тонкой изоляционной стенкой и численного метода решения её уравнений, обеспечивающего более высокую скорость решения поставленной задачи по сравнению с известными методами. Особенностью подхода, предлагаемого в работе, является то, что он позволяет построить динамически разностную схему по исходным данным для математического пакета Maxima в виде нелинейной системы алгебраических уравнений.

2 Математическая модель процесса

Построим математическую модель процесса нанесения покрытия на плоскую деталь в гальванической ванне с бесконечно тонкой плоской перегородкой-

фптт

Тс >Н1 еа< ч. ст е? 1Кс а с: оч ВС

им 1

А но Л

0.3 0.6 0.9 1.2 1.5 1.8 2.1 2.4 2.7

Рис. 1: Горизонтальное сечение гальванической ванны.

изолятором, в которой имеются поперечные отверстия прямоугольной формы. Предполагается, что ванна имеет форму параллелепипеда. Анод и катод являются плоскими и располагаются напротив друг друга вдоль соответствующих стен. Тогда мы имеем ванну, когда в горизонтальном сечении, перпендикулярном оси Ог7 не меняется конфигурация электрического поля внутри ванны, т.е. для любого такого сечения мы имеем картину, представленную на рис. 1. Размер ванны I х I.

Поперечный вид плоской перегородки изолятора представлен на рис. 2.

Рис, 2: Поперечный вид плоской перегородки-изолятора.

Конфигурация поперечных отверстий прямоугольной формы соответствует их расположению в продольном сечении ванны на рис. 1. Высота токонепро-водящего экрана равна высоте гальванической ванны, а длина равна/.

Предположим, что к аноду подведено напряжение Ц"а, к катоду О В. Стенки ванны являются изолирующими, т.е. градиент потенциала по нормали к их поверхности равен нулю.

Обозначим через <£>{х, у) потенциал электрического поля в точке с координатами (х,у). Запишем уравнение Лапласа, описывающее распределение потенциала в электролите гальванической ванны

д V + dV =0

дх2 ду2

а)

Далее опишем граничные условия для рис. 1. Для границы с изоляторами стенками они имеют вид:

д^р дх

(0;y)

= °'уG (0; дх

= 0, у G (0; l),

(i;y)

др ду

др ду

= 0, х G [0; xai), др

(x;0) ду

= 0, х G [0; xci), др

(x;l) ду

(x;0)

= 0, х G (Xar; l],

= 0, x G (xcr; l],

(x;l)

(2)

(3)

(4)

где xa/, xar, xc/ и xcr - координаты расположения краев анода и катода соответственно. Например, для рис. 1 xa/ = 0,6 дм, xar = 2,1 дм, хс/ = 0,3 дм, Xcr = 2, 4 дм ир и/= 2,7 дм.

Потенциалы поля на аноде и катоде связаны соотношениями

<р(х, 0) + Fa (ia(x)) = Ua, X G [Xa/; Xar], (5)

<p(x, /) + Fc (ic(x)) = 0, X G [Xci; Xcr], (6)

где Fa (ia) и Fc (ic) - так называемые функции анодной и катодной поляризации [5,6], ia, ic - соответственно, анодная и катодная плотности тока.

Для плоской бесконечно тонкой перегородки изолятора с поперечными отверстиями граничные условия аналогичны условиям для изоляторов стенок:

д^ ду

= 0, X G

(x;p)

x(k); x(k) p/ ' pr

, k = 1, n,

(7)

(к) (к)

где р - у-координата перегородки-изолятора, хр/ , хр/ - координаты расположения краёв к-ого отверстия, п - количество отверстий.

Поскольку ток течёт от большего потенциала к меньшему (положитель-

у

на аноде и катоде запишется как

«о = -Х

д^ ду

[ ; ] ' = д^ , x G [xal; xar], «c х T¡

(x;0) дУ

(x;l)

, x G [xc/; Xcr], (8)

гДе X удельная проводимость электролита.

Таким образом, мы построили математическую модель Ц) распре-

деления потенциала внутри гальванической ванны.

Толщина получаемого на катоде покрытия определяется из формулы

К

^(ж) = — гс(х)Д£,

где ж С [хс/; хсг], К - электрохимический эквивалент металла покрытия, р -плотность металла покрытия, Д£ - время нанесения покрытия.

3 Описание численного метода

Для численного решения нелинейной задачи Д) заменим производные

их разностными аналогами. Введём сетку по координатам х и у:

хг = (г - 1)Нх, Уз = (з - 1)Ну,

г = 1,Жх, з = 1,Ду,

где НХ и Ну - шаги сетки по х и у соответственно, ДХ и Ду - количество узлов. Производным, входящим в задачу ([]]) сопоставим соотношения:

( <г-1,з - 2(рг,з + <г+1,з + <г,з-1 - 2(рг,з + <г,з+1 = о

НХ

Н2

= 0, з = 2,Д - 1, N + 1,Ду - 1;

г = 2, ДХ - 1 & з = 2, - 1, + 1, N - 1, г = р(к),р[к) & з = Д & к = Т~п;

<2,з - <1,3

НХ

<г,2 - <г,1 НУ

<г,Жу - <г,Жу-1

= 0, з = 2,Д - 1, N + 1,Ду - 1;

= 0, г = 1, а/ - 1, аг + 1, ДХ;

Н

= 0, г = 1,с/ - 1, сг + 1,Дх;

у

Ну

= 0, г =1,р(к) - 1, рГк) + 1,Дх & к =1,п;

(к)

. „ / <г,2 - <г,1 > гг А • -

<г,1 + ^ -х-¡Г- ) - иа = 0, г = а/, аг;

НУ ,

<г,Жу - <г,Жу-Л . _

-^Н-^ ) =0, г = С/,сг,

где а/, аг, с/, сг - номера узлов сетки, соответствующих левому и правому

(к) (к)

краю анода и катода, р/, рГ номера узлов сетки, соответствующих краям отверстий перегородки-изолятора, п = 4 Например, для рис. 1 а/ = 7, аг = 22 с/ = 4 сг = 25, р((1) = 4 рГ1) = 8 Р((2) = 11 рГ2) = 14 Р((3) = 18,

р[3) = 20 Р((4) = 22, р4) = 26 N = 25.

Заметим, что количество уравнений и неизвестных в системе равно т = N • Жу.

Соберём неизвестные величины

<1,Ъ . . . , <1,Жу, <2,1, <2,2, . . . , <2Д„, . . . , <N^,1, ^^ . . . , <Мх,Му

в вектор и обозначим его как

Ф = [^х, ^2, ...,

Левая часть полученной системы является вектор-функцией Л от Ф. Таким образом, мы перепишем систему как

Л(Ф) = 0. (9)

Для численного решения нелинейной системы Ц) предлагается использовать метод Ньютона. Обозначим через Ф(0) вектор начального приближения. Выбор значений его компонентов, т.е. осуществляется следующим образом.

Поскольку в направлении от анода к катоду потенциал убывает, примем:

= иа, = 0, ^ = ТЖ

В промежуточных узлах сетки сделаем линейную интерполяцию, т.е.

(0) тт ^ - 3 (лп\

^ = и N,-1. (10)

Как известно (см., например, |3]), численная схема метода Ньютона имеет вид:

ф(г) = фСт-1) - (ф(г-1))] -1 Л (ф(г-1)) , (11)

где г = Т, 2, 3,... - номер итерации, У Л (Ф(г-1)) - матрица Якоби функции Л в точке Ф(г-1).

На сегодняшний день матрица Якоби в соотношении (11) при реализации метода Ньютона в математических пакетах не обращается, т.к. это требует существенных вычислительных затрат. Для нахождения значения приближения Ф(г) на текущей итерации выражение (11) можно привести к системе линейных уравнений вида

УЛ (ф(г-1)) Ф(г) = УЛ (Ф(г-1)) Ф(г-1) - Л (Ф(г-1)) ,

для решения которой используется наиболее эффективный метод Ы1-разложение (например, в |7] для решения задач теории упругости).

4 Пример расчёта

Проиллюстрируем разработанный алгоритм на примере расчёта процесса нанесения никелевого гальванического покрытия на плоскую пластину для электролита, содержащего №804 (2 моль/л) и Н3ВО3 (0,5 моль/л) при температуре 25 °С.

Для вычислительного эксперимента было выбрано положение и конфигурация токонепроводящего экрана, показанного на рис. 1 и 2.

Чтобы рассчитать распределение толщины получаемого на катоде покрытия в мкм, используем соотношение

я(х, Дг) = 0(х) • Дг = 100Кгс(х)Д£ = -100 —

Р Р ду

« 100 Дг,

Р Ну

Дг

(х;/)

где х € [хс/; хсг], г = с/, сг; К - электрохимический эквивалент никеля, г/(А^ч); р - плотность никеля, г/см3, Дг - время нанесения покрытия, ч. В нашем случае К = 1,09 г/(А^ч), р = 8,902 г/см3, х = 0, 515 (Ом^дм)-1 [5]. Далее положим, что время Дг = Дгса1с = 0, 5 ч.

Введём обозначение

¿(х) = q (х, Д^а1с) .

Заметим, что для раствора электролита поляризационные кривые приводятся в литературе в виде графиков, например, в работе 5, с. 275, 289].

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

Анализ графиков из |5] показал, что наилучшим является квадратичный вид аппроксимирующей зависимости.

Полученное нами аналитическое выражение для анодной поляризации:

^а(га) = -4, 267га + 5,867га

для га € [0; 1, 5], А/дм2; для катодной поляризации:

Fc(гс) = 0,883г2 - 2, 242гс (12)

для гс € [0;3], А/дм2.

Отрицательное значение изменения ^с(гс) потенциала катоде (в вольтах

уже учтено в функции (12) на приведённом отрезке изменения плотности

с

5 Результаты вычислений

Для отыскания приближённого решения системы мы использовали pea-лизацию метода Ньютона в математическом пакете Maxima £] с точностью £ = 10-2. Под точностью вычислений мы понимаем модуль разности векторов ф(г) и Ф(г-1) на соседних итерациях, т.е. когда выполнено условие

<

вычисления завершаются. Для автоматизации формирования скрипта для математического пакета с системой из большого числа уравнений разработана модификация программы |9] на языке С с учётом наличия в ванне изолятора стенки.

Для взаимодействия с пакетом авторы использовали перенаправление ввода/вывода с вызовом команды maxima под Linux. Каждая команда записывается в отдельной строке входного файла. Перед вызовом в пакете реализации метода Ньютона нужно выполнить команду display2d:f alse$, чтобы результаты вычислений выводить в виде строк. В самой команде mnewt on О определить левые части уравнений системы в символьном виде, вектор ф неизвестных, а также начальное приближение, сформированное правилу

(10). После того, как пакет произведёт вычисления, результат (учитывается из выходного текстового файла.

На рис. 3 приведены результаты вычислительного эксперимента для Ua = 3 В, l = 2, 7 да, Nx = = 28 и hx = = 0,1 дм. При этом для наглядности точки сетки соединены кубическим сплайном.

С учётом того, что мы работаем с плоским сечением гальванической ванны, для оценки равномерности полученного покрытия используем выражение

r =1 г ¿х, (i3)

L J Xcl ömin

где L - длина катода в продольном сечении (рис. 1),

ömin = min ö (x).

x£ [xcl ;Xcr]

Так как мы работаем с дискретным аналогом функции ö(x), заменим

Рис, 3: Распределение покрытия по поверхности катода, размерность по оси х - дм, по оси 6 - мкм.

интеграл в формуле (13) суммой вида

Я =

6тлд L

(¿(хг) 6тш) На

г=сг

6тт(сг с/ )Нх ^

6тт) На

(6(хг) - 6тт)

г=сг

г=сг

6тш (сг с/)

Для рис. 3 значение Я составило 0,024 или 2,4%, что считается удовлетворительным.

Заметим, что для проведения вычислений время Дг = Дгса1с было выбрано условно, чтобы проверить адекватность модели и вычислить значение показателя Я. На практике требуется подобрать такое время нанесения покрытия Дг = Дгргос, чтобы обеспечить заданное среднее значение толщины получаемого покрытия. Получим расчётную формулу для этого времени.

1

1

Среднее значение функции

i=Cl

— 1 fxcr 1 Сгл

6 — - 6(x)dx « (—-—

Wxcl (cr - C1 )hx i=C' Cr - Q

i=Cl

По полученным значениям 6(xi) определяем среднее значеиие 6 толщины рассчитываемого покрытия. Далее отметим, что

6 = в • Arcale,

6set — в • Atproc,

где в - среднее значение функции в(х) на отрезке [xc1; xcr], откуда

a t = 6set A t

Atproc — Atcalc. 6

6 Заключение

Полученные результаты показывают эффективность применяемого численного метода квадратичная скорость сходимости метода Ньютона на сетках с большим количеством узлов даёт выигрыш во времени примерно в 10 раз по сравнению с одним из наилучших численных методов для такого вида задач итерационного метода, описанного в работах 1,2]. При этом был использован численный метод Ньютона с реализацией в математическом пакете Maxima для решения системы нелинейных алгебраических уравнений, получаемых как разностный аналог задачи распределения потенциала в гальванической ванне. Данная задача включает в себя уравнение Лапласа и нелинейные краевые условия III-его рода на аноде и катоде.

Разработанные модель и эффективный метод решения её уравнений могут быть использованы для решения задачи оптимизации минимизации критерия (13) неравномерности покрытия.

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

1] Дутов А. В., Литовка Ю. В., Нестеров В. А., Соловьев Д. С., Соловьева И. А., Сыпало К. И. Поиск оптимального управления токовыми

режимами в гальванических процессах со многими анодами при разнообразии номенклатуры обрабатываемых изделий // Известия Российской академии наук. Теория и системы управления, 2019, №1. С. 78 88.

[2] Литовка Ю. В., Михеев В. В. Численный расчет электрического поля в гальванической ванне с биполярными электродами // Теоретические основы химической технологии, 2006, том 40, №3. С. 328 334.

[3] Демидович Б. 77., Марон И. А. Основы вычислительной математики. СПб.: Лань, 2006. 672 с.

[4] Pchelintseva I. Yu., Pchelintsev A. N., Litovka Yu. V. Modeling of metal distribution when coating flat metal plates in electroplating baths // International Journal of Numerical Modelling: Electronic Networks, Devices and Fields, 2021, vol. 34, iss. 2, e2830, 10 pp.

[5] Кудрявцев H. Т. Электролитические покрытия металлами. M.: Химия, 1979. 352 с.

[6] Жук Н. П. Курс теории коррозии и защиты металлов. М.: Альянс, 2014. 472 с.

[7] Толмачев А. В., Коновалов А. В., Партии А. С. Эффективность алгоритма LU-разложения с двухмерным циклическим распределением матрицы для параллельного решения упругопластической задачи // Программные продукты и системы, 2013, №3. С. 94 99.

[8] Maxima computer algebra system,

[9] Пчелинцева И. IO.. Пчелинцев А. II.. Литовка Ю. В. Свидетельство о государственной регистрации программы для ЭВМ №2020614929. Численное решение уравнения Лапласа для расчёта распределения электрического потенциала в гальванической ванне на базе математического пакета Maxima. 29.04.2020 г.

http://maxima.sourceforge.net/ги/

Mathematical model and numerical scheme for calculation of electric fields in galvanic baths with non-conductive screen

I. Yu. Pchelintseva, Yu. V. Litovka Tambov State Technical University, e-mail:

[email protected], [email protected]

Abstract. A mathematical model of the electric field in a electroplating bath with a flat anode and cathode, which has an infinitely thin flat insulator wall with transverse slots, is considered. Such a non-conductive screen is necessary for a more uniform coverage of the cathode detail. In this article, a transition is made to the difference analogue of the problem. A numerical method based on Newton's method for solving nonlinear algebraic equations is described, a computational experiment for 4 slits is carried out. The obtained results show the effectiveness of the applied numerical method.

Keywords: Laplace's equation, Newton's method, non-uniformity criterion, coating thickness.