CC BY
13ВЕСТН. САМАР. ГОС. ТЕХН. УН-ТА. СЕР. ТЕХНИЧЕСКИЕ НАУКИ. 2016. № 2 (50)
Электротехника
УДК 621.365
КОМПЛЕКСНОЕ МОДЕЛИРОВАНИЕ И УПРАВЛЕНИЕ ПРОЦЕССОМ НЕПРЕРЫВНОГО ИНДУКЦИОННОГО НАГРЕВА ФЕРРОМАГНИТНЫХ ЗАГОТОВОК*
А.А. Базаров, А.И. Данилушкин, В.А. Данилушкин
Самарский государственный технический университет Россия, 443100, г. Самара, ул. Молодогвардейская, 244
Предложен подход к комплексному решению задачи моделирования и управления процессом непрерывного индукционного нагрева ферромагнитной загрузки до температур пластической деформации. Моделирование электромагнитных и тепловых процессов осуществляется в двумерной постановке методом конечных элементов с помощью программы COMSOL Multiphysics®. Для синтеза системы регулирования использована аналитическая модель объекта, содержащая инерционное звено и звено транспортного запаздывания. Анализ качества процессов управления осуществляется с помощью численных методов решения обыкновенных дифференциальных уравнений, встроенных в пакет БтиНпк. Для снижения амплитуды колебаний и времени переходного процесса предложено разделение индуктора на две независимые секции. Исследование динамических режимов системы с разными аппроксимирующими выражениями для звена запаздывания позволило найти параметры системы регулирования с компенсатором, наиболее корректно учитывающие особенности практической реализации. Приведены результаты моделирования системы с учетом распределенных параметров объекта.
Ключевые слова: индукционный нагрев, распределенная система, численное моделирование, управление, температурное поле, аппроксимация.
Индукционные нагревательные установки представляют собой сложные технологические объекты, характеризующиеся наличием распределенного управляющего воздействия и распределенной функции состояния, в качестве которой чаще всего рассматривают температурное распределение по объему нагреваемого изделия. В таких распределенных системах управляющие воздействия и функция состояния объекта зависят от пространственных координат и времени [1-6]. В этом отношении индукционные нагревательные установки как объекты управления существенно отличаются от объектов с сосредоточенными параметрами, так как состояние последних однозначно определяется заданием описыва-
* Работа поддержана грантом РФФИ № 15-08-03053/16 (904/15-17).
Александр Александрович Базаров (д.т.н., доц.), профессор кафедры «Электроснабжение промышленных предприятий».
Александр Иванович Данилушкин (д.т.н., проф.), профессор кафедры «Электроснабжение промышленных предприятий».
Василий Александрович Данилушкин (к.т.н.), доцент кафедры «Электроснабжение промышленных предприятий». 128
ющих процесс величин в какой-либо одной точке. Задачи управления объектами с распределенными параметрами оказываются качественно более сложными по сравнению с аналогичными задачами для объектов с сосредоточенными параметрами ввиду целого ряда принципиальных особенностей [2, 5, 6].
Решение практически важных задач эффективного управления индукционными нагревателями требует комплексного подхода, включающего разработку адекватных моделей взаимосвязанных электромагнитных и тепловых процессов и создание алгоритмов и систем управления, обеспечивающих качественное функционирование всего технологического комплекса [2, 5, 7-9]. Поставленные задачи должны решаться на основе проблемно-ориентированных математических моделей управляемых температурных полей, позволяющих создать основу для повышения эффективности, качества и экономичности технологического процесса. Полученные проблемно-ориентированные модели позволят разработать оптимальную конструкцию нагревателя, реализовать эффективные алгоритмы управления, обеспечить более высокую точность управления температурным распределением, высокое быстродействие и удобство регулирования.
Решение задач управления для объектов с распределенными параметрами наиболее эффективно при использовании аналитических методов описания пространственно-временных распределений управляемых величин. Непосредственное применение цифровых моделей как электромагнитного, так и температурного полей, безусловно, повышает точность расчетов и является незаменимым инструментом при решении задач проектирования конструктивных параметров установок, но при решении задач управления приводит к резкому возрастанию необходимого для расчетов компьютерного времени, обусловленного проведением достаточно сложных вычислительных процедур с многократным расчетом состояния объекта на различной стадии процесса управления. Учитывая сказанное, целесообразно применение в качестве основной модели аналитического решения, но уточнение параметров должно осуществляться с помощью численных моделей [5, 9].
В общем случае процесс индукционного нагрева рассматриваемого класса объектов описывается нелинейной взаимосвязанной системой уравнений Максвелла [11, 12,] и Фурье [13] соответственно для электромагнитного и теплового полей с соответствующими краевыми условиями:
rot H = +— ;div H = 0;
- Jè - 0)
rot E =--;div E = 0;
dt
ЯТ __'
c(T)g(T) — = div(l(T)grad T) - div EH . (2)
dt
Здесь H, E, B, D - векторы напряженностей магнитного и электрического полей, магнитной и электрической индукции, ô = E / р - плотность тока проводимости, c(T ), à(T), T ) - удельные значения теплоемкости, теплопроводности и плотности металла, T = T (r, x, t ) - температурное распределение в металле загрузки, r, x - радиальная м аксиальная координаты цилиндрической загрузки, t - время. Объемная плотность внутренних источников тепла, индуцируемых в загрузке, определяется дивергенцией вектора Пойнтинга П = -div EH [12].
В качестве основной расчетной модели удобна конечно-элементная формулировка. Моделирование электромагнитных и тепловых процессов осуществляется в двумерной постановке методом конечных элементов с помощью программы COMSOL Multiphysics®.
Решение задачи электромагнитного поля достигается использованием векторного магнитного потенциала {A} и скалярного электрического потенциала V, которые выражаются следующим образом:
М= rot\A\L (3)
M = (4)
Чтобы функция {a} была определена, нужно определить значение ее дивергенции. Для этого добавляется условие, которое называется калибровкой Кулона:
div{A}= 0. (5)
В результате получим следующую систему уравнений:
'¿J + }; (6)
rot\A}= B} (7)
div{A}= 0. (8)
Используя соотношение
rot (rot{A}j= grad(div{A }j- V2 {a} (9)
при p=const из (6) получим уравнение
V 2 {A}- 7Ws{i}-{7 }= 0. (10)
Уравнение Пуассона (10) дополняется граничными условиями Дирихле и Неймана на различных участках границы:
}= 0 на Sj; (11)
/
rot
= 0 на S2. (12)
дп
При построении конечно-элементной модели в результате ансамблирования получаем систему алгебраических уравнений:
М + /М" $+{?}= 0. (13)
Здесь [К], - матрицы жесткости и источников тепла. Решение данной задачи осуществляется итерационным методом.
Мощность внутренних источников тепла, характеризующих нагрев проводящих тел индукционной системы, вычисляется для каждого элемента по закону Джоуля - Ленца:
Í
М=1 (
рУе> = - I s
2
Vole ^ 0
E • E
dVol, (14)
*
где Е - величина, сопряженная к Е .
Рассмотрим постановку тепловой задачи для двумерного случая, на базе которой построим конечно-элементную формулировку. Переход к более простым случаям реализуется с помощью упрощений.
Первый закон термодинамики для объемных тел записывается в виде дифференциальных уравнений
С ^т1 „ Л
gc
T
--— + {v}— {l}—) + {l}— {q} = q, (15)
V-t 0
где {l}=|--—■ —- векторный оператор; {v}=|vx vy — - вектор, характеризующий скорость переноса тепла за счет перемещения сред (адвекция); {q}- вектор теплового потока; q - скорость образования тепла в конечном объеме. Полученные в результате преобразований матрицы жесткости [K ], теплоемкости
[С ] и источников {Q} с учетом замены временной производной -— конечно-
-t
разностным аналогом, объединяются в систему уравнений.
i^ + Y [к ]V }„+1 =VH-(1 -yMV }„ +{Q}„ (16)
VAt У VAt У
где At- временной шаг, n - номер шага, g- коэффициент, принимающий значения от 0 до 1. Последнее выражение переписывается в виде
[K^}n+1 =Q} •
Comsol позволяет осуществлять полноценное решение мультидисциплинар-ных задач, например электротепловых, но в данном случае такой подход очень затратен, так как тепловая задача, содержащая перемещающиеся среды, обладает низкой вычислительной устойчивостью. Решение проблемы заключается в уменьшении шага по времени до десятых долей секунды. С учетом достаточно большого числа элементов объем выходного файла измеряется в гигабайтах. В такой ситуации объединение моделей было бы неоправданным.
Важным является не просто решение тепловой задачи, но и моделирование, хотя бы в упрощенном виде, работы системы регулирования. Средства Comsol позволяют реализовать полноценное релейное управление с помощью функции Хевисайда (flc1hs(u1 - —d, scale)), имеющейся в библиотеке. Здесь переменная u1 соответствует температуре точки объекта, выбранной в качестве точки регулирования; —d - температура задания. Величина переменной scale выбирается исходя из требования к длительности перехода функции из одного состояния в другое.
Аналоговый регулятор ограничен в выборе функций, так как использование оператора Лапласа не предусмотрено. Можно осуществить построение модели на основе связки Simulink-Comsol, но этот путь для данной задачи также требует слишком больших временных затрат.
Решение задачи рассмотрим на примере индукционного нагревателя непрерывного действия для нагрева стальных колец под раскатку. Внешний диаметр кольца - 0,15 м, внутренний - 0,08 м. Скорость перемещения заготовок -0,005 м/c. Нагрев производится до температуры 1273 К. При достижении температуры 1023 К происходит резкое снижение магнитной проницаемости до единицы. Отклонение температуры по сечению заготовки на выходе из нагревателя не превышает 5 % от заданного значения. Для обеспечения возможности регулирования температуры предусмотрен 20%-ный запас по мощности.
В результате предварительного итерационного решения электромагнитной и тепловой задач получили следующие параметры индуктора и распределения температуры: активная мощность загрузки 270 кВт; активная мощность загрузки на участке до 1023 К - 170 кВт; активная мощность на участке с немагнитной загрузкой 80 кВт; длина ферромагнитного участка 0,5 м; длина немагнитного участка 1 м; отклонение температуры от заданного значения в заготовке на выходе ±50 К. Моделирование процесса нагрева движущейся загрузки при применении релейного регулятора с бесконтактным датчиком температуры показало значительные колебания на выходе из индуктора (рис. 1). Сравнивая полученный результат с нагревом загрузки в двухсекционном индукторе (рис. 2) [15] с раздельным управлением, нужно отметить отсутствие колебаний температуры на выходе из индуктора во втором случае. Это объясняется более слабой реакцией на снижение мощности в начале индуктора из-за процессов регулирования. На основе анализа сделан вывод о необходимости разделения индуктора на две независимые зоны с нагревом до точки магнитных превращений и с нагревом до требуемой по технологии температуры на выходе. На первом участке возможно применение разомкнутой системы при наличии стабилизации уровня напряжения. Частота напряжения может быть различной, от 50 Гц и выше. Проведенные расчеты для секций индуктора показали возможность применения для первой секции как частоты 1000 Гц, так и частоты 50 Гц. В первом случае проще обеспечить стабилизацию температуры на выходе первой секции за счет возможности плавного регулирования напряжения источника питания. Во втором случае, при использовании частоты 50 Гц, использование регуляторов, за исключением релейных, затруднительно. Отклонения сетевого напряжения приведут к изменению выходной температуры. Однако эти возмущения проще отрабатываются системой регулирования второй секции, так как имеют меньшую амплитуду и нерегулярный характер по сравнению с процессами в односекционном нагревателе (см. таблицу).
Параметры секций индуктора
Секция Р , Гц Ь, м Р2, кВт Ри , кВт и и, В w
1 50 0,5 170 233 210 15
1 1000 0,5 170 219 545 15
2 1000 1 80 119 696 27
Для второй секции имеется запас по напряжению, что обеспечивает добавку мощности 30 %.
Для синтеза системы управления рассмотрим упрощенную постановку задачи, когда можно пренебречь распределением температуры по толщине нагреваемого тела, а также переносом тепла по длине за счет теплопроводности. Теплообмен с окружающей средой рассмотрим в виде конвекции. Уравнение теплопроводности формулируется как для пластины или стержня:
су
гдТ (х, t) дТ (х, t)] а
V
+ V-
дt дх
+ -(Т (х, t)-Тс ) = б(х, t),
г
Г (17)
где а - коэффициент конвективного теплообмена; V - скорость; г - радиус загрузки.
Граничные условия:
T(x,0) = T0 (x), T(0,t) = g(t), X > 0, t > 0
(18)
1300 1200 1100 1000
^ 900
<o
£ 800 I
£ 700 s
P
600 500 400 300 200
0 100 200 300 400 500 600 700 800 900 1000
время, с
Рис. 1. Переходный процесс при нагреве движущейся загрузки при наличии участков
с разной магнитной проницаемостью: 1 - температура поверхности на выходе из индуктора; 2 - температура поверхности на выходе
из первого участка индуктора
1 Г
2
\
\
Приведем (17) к стандартному виду:
a дТМ + b ЭТМ + kT(x, t) = f (x, t).
dt
Здесь a = cg; b = cgv, k = a.
r
dx
(19)
Импульсная переходная функция для уравнения (19) с учетом граничных условий (18) будет иметь вид [4]
к
G(x, X, t) = l(x - X) exp
5[bt - a(x -X)].
Ь(х -X).
Выбрав любую удобную точку, чаще всего х=Ь, можно получить выражение для передаточной функции
W (x, X, p) = l(x - X)exp
ap + k
_ b(x -X).
(20)
Здесь p - оператор Лапласа. Для одноточечной системы
W (x, p) = 1 exp
ap + k b(x -X)_
= exp
ax
p
• exp
k b
x
(21) 133
b
1300 1200 1100 1000 900
(С
£ 800
700
£
ш
600 500 400 300 200
1 /
2,
3 /
................ ! I
...............
100 200 300 400 500 600
время, с
700
800
900 1000
Рис. 2. Переходный процесс при нагреве движущейся загрузки при наличии участков
с разной магнитной проницаемостью: 1 - температура поверхности на выходе из второй секции индуктора; 2 - температура поверхности на расстоянии 0,5 м от выхода из первой секции; 3 - температура поверхности на выходе из первой
секции индуктора
При использовании такого подхода нужно иметь в виду, что переход от распределенной системы к сосредоточенной сопровождается существенным упрощением модели процесса и может отразиться, за редкими исключениями, на качестве процесса регулирования.
Моделирование системы управления, содержащей звено запаздывания, часто удобнее проводить, используя аппроксимирующие выражения. Для этого существует несколько вариантов. Наименее трудоемкими являются аппроксимации звена чистого запаздывания разложением в ряд Паде второго и четвертого порядков. Аппроксимирующее выражение второго порядка имеет вид
2 6 12 Р --Р +
Ж (х, X, р) = ехр(- рт) =
т
2 6 12 Р +-Р + т т2
(22)
где т
ах
Т
На рис. 3 показана блок-схема аппроксимирующего выражения для звена запаздывания.
Параметры схемы определяются на основе величины времени запаздывания, являющейся функцией длины индуктора и скорости перемещения:
т = ах = = 200.
ь
еуу
Числовые значения коэффициентов к1 и k2 при т = 200 с будут равны
6
к1 = И = _12- = 0.0003; к 2 = -
.2
200
2
т 200
0.03.
Рис. 3. Блок-диаграмма аппроксимации транспортного запаздывания разложением в ряд Паде второго порядка
Сравнение переходных процессов для систем со звеном запаздывания из библиотеки 81шиНпк и аппроксимирующего выражения (см. рис. 3) показало, что имеются незначительные отличия в характере процессов [14, 15, 16]. В целом результаты близкие.
Рис. 4. Модель системы управления нагревом движущейся загрузки с компенсатором
запаздывания
Модель системы управления для объекта, содержащего звено запаздывания, представлена на рис. 4. Для самого объекта применено более точное представление блока транспортного запаздывания из библиотеки, а в схеме компенсатора
135
запаздывания использована схема разложения в ряд Паде. Такой подход точнее соответствует реальному положению, когда параметры процесса и его модели имеют отличия. Исследования, проведенные со схемами, когда сам объект и его модель в компенсаторе представлены одинаковыми блоками, показали более простую настройку регулятора для получения переходных процессов желаемого вида.
Для схемы, показанной на рис. 4, настройка параметров регулятора позволила получить приемлемые показатели. При коэффициенте усиления, равном 10, постоянной времени дифференцирования 25 с и постоянной времени интегрирования 0,05 с отклонение от заданной величины не превышает 5 % (рис. 5).
1200
о.
Sí
I 400 .................:...............i................:..................-
си :
ь ; : :
200 ...............i...............:...............¿..............-
0-¡-1-¡-
0 500 1000 1500 2000
время, с
Рис. 5. Переходный процесс в системе с запаздыванием при использовании компенсатора и ПИД-регулятора
Полученные результаты при моделировании процессов на численной модели (см. рис. 1, 2) и при построении переходных процессов в системе управления (см. рис. 5) иллюстрируют изложенный подход для определения передаточных функций и синтеза системы управления.
Выводы. На основе комплексного решения задач численного моделирования процессов нагрева ферромагнитноой загрузки до температур пластической деформации в индукционном нагревателе непрерывного действия и идентификации объекта с использованием аналитических моделей показана целесообразность разделения нагревателя на две секции. Параметры секций нагревателя определяются координатой точки магнитных превращений. Разработан алгоритм синтеза системы управления распределенным объектом с транспортным запаздыванием на базе аналитических моделей. Анализ переходных процессов для температурного распределения, выполненный с помощью численных методов решения обыкновенных дифференциальных уравнений, встроенных в пакет Sim-ulink, показал целесообразность управления процессом нагрева ферромагнитных заготовок в нагревателе непрерывного действия с разделением нагревателя на две секции в точке перехода к немагнитному состоянию. Использование компенсатора запаздывания позволило снизить амплитуду колебаний до приемле-
мых значений. Аппроксимирующее выражение в виде ряда Паде второго порядка обеспечивает допустимую погрешность.
Результаты исследования показали эффективность применения предложенного подхода к анализу и синтезу системы управления рассматриваемого класса объектов.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Бутковский А.Г. Структурная теория распределенных систем. - М.: Наука, 1977. - 320 с.
2. Рапопорт Э.Я. Оптимизация процессов индукционного нагрева металла. - М.: Металлургия, 1993. - 279 с.
3. Сиразетдинов Т.К. Оптимизация систем с распределенными параметрами. - М.: Наука, 1977. -4B0 с.
4. Рапопорт Э.Я. Структурное моделирование объектов и систем управления с распределенными параметрами. - М.: Высшая школа, 2003. - 299 с.
5. Базаров А.А., Данилушкин А.И., Зимин Л.С. Анализ задачи пространственно-распределенного управления индукционным нагревом колец на основе структурного метода для систем с распределенными параметрами // Вестник Самарского государственного технического университета. Сер. Технические науки. - 199B. - № 5. - С. 115-120.
6. Данилушкин А.И., Данилушкин И.А. Метод вторичных источников для моделирования электромагнитных процессов при индукционном нагреве // Вестник Самарского государственного технического университета. Сер. Физико-математические науки. - 199B. - № б. - С. 140-142.
7. Rapoport E.Y. Analytical Construction of Aggregated Controllers in Systems with Distributed Parameters // Journal of Computer and Systems Sciences International. - 2012. - Vol. 51. - No 3. -рp. 375-390. ИФ журнала в Web of Science 0,191.
B. Дилигенская А.Н., Рапопорт Э.Я. Идентификация пространственного распределения внутренних источников тепла в обратных задачах теплопроводности // Вестник Самарского государственного технического университета. Сер. Технические науки. - 2011. - № 4(32). - С. 157164. ИФ журнала в РИНЦ 0,097.
9. Базаров А.А. Моделирование процесса теплопроводности для задач синтеза систем управления в среде MATLAB // Вестник Самарского государственного технического университета. Сер. Технические науки. - 2005. - № 33. - С. 7-11.
10. Rudnev V. Handbook of Induction Heating. Marcel Dekker Inc. - New York, 2003, p. 277.
11. Немков В.С., Демидович В.Б. Теория и расчет устройств индукционного нагрева. - Л.: Энерго-атомиздат, 19BB. - 2B0 с.
12. Вайнберг А.М. Индукционные плавильные печи. - М.: Энергия, 1967. - 415 с.
13. Лыков А.В. Тепломассообмен: Справочник. - М.: Энергия, 197B. - 4B0 с.
14. Данилушкин А.И., Мостовой А.П. Анализ эффективности пусковых режимов двухсекционного индукционного нагревателя методического действия // Вестник Самарского государственного технического университета. Сер. Технические науки. - 2014. - № 4 (44). - С. 113-120.
15. Данилушкин А.И., Князев С.В., Мостовой А.П. Моделирование стационарного распределения температуры металла в проходном индукционном нагревателе // Вестник Иркутского государственного технического университета. - 2012. - № 9(бВ). - С. 41-4б.
16. Базаров А.А. Аналоговое моделирование процесса теплопроводности для задач синтеза управления в среде MATLAB // Математическое моделирование и краевые задачи: Труды Всероссийской научной конференции. - Самара: СамГТУ, 2004. - С. 29-31.
Статья поступила в редакцию 3 марта 2016 г.
INTEGRATED MODELING AND PROCESS CONTROL
OF CONTINUOUS INDUCTION HEATING OF FERROMAGNETIC
WORKPIECES
A.A. Bazarov, A.I. Danilushkin, V.A. Danilushkin
Samara State Technical University
244, Molodogvardeyskaya st., Samara, 443100, Russian Federation
The approach to the complex solution of the modeling problems and control process continuous induction heating ferromagnetic workpieces to plastic deformation temperature is suggested. Simulation of electromagnetic and thermal processes is carried out in a two-dimensional formulation by finite element method by using the COMSOL Multiphysics software. The analytical model of the plant comprising an inertial element and a transport delay element was used. Analysis of the control process quality is carried out by means of numerical methods of solving ordinary differential equations, built in Simulink® software. To reduce oscillation amplitude and transient time the separation of inductor into two independent sections is proposed. Investigation of dynamic modes of the system with various approximating expressions for the transport delay element allowed to find the control system parameters with a compensator taking into account the advantages of implementation more correctly. The simulation system results taking into account the distributed facility parameters are given.
Keywords: induction heating, distributed system, numerical simulation, control, temperature field, approximation.
Aleksandr A. Bazarov (Dr. Sci. (Techn.)), Professor.
Aleksandr I. Danilushkin (Dr. Sci. (Techn.)), Professor. Vasily A. Danilushkin (Ph.D. (Techn.)), Associate Professor.
