УДК 621.313.333
А.С. Замчалкин, В.А. Тюков
Численное моделирование процесса пуска асинхронного двигателя
Продемонстрирована возможность расчета динамических режимов асинхронных двигателей на базе численного решения нелинейного и нестационарного уравнений для вектор-потенциала магнитного поля в двумерной постановке с учетом движения ротора и наличием внешних электрических цепей в программном пакете COMSOL Multiphysics. Представлены результаты численного моделирования процесса прямого пуска асинхронного двигателя с нагрузкой на валу вентиляторного типа.
Ключевые слова: асинхронный двигатель, численное моделирование, электромагнитное поле, пусковой режим.
На сегодняшний момент, несмотря на более чем вековую историю развития теории асинхронных двигателей (АД), до конца не изучены многие вопросы, касающиеся неустановившихся режимов работы АД. Практика электромашиностроительных предприятий ясно показывает, что существующие методы расчетов АД зачастую дают существенные расхождения с экспериментом даже при расчете простейших случаев установившегося режима. Наиболее ярко это проявляется в крупных АД с нетрадиционной геометрией активного объема и значительными электромагнитными нагрузками. Методики предприятий для расчета неустановившихся режимов, как правило, представляют собой математические модели, основанные на сложившейся в электромеханике системе допущений, которая предполагает переход к эквивалентной двухфазной системе с сосредоточенными параметрами, а нередко и к однофазной эквивалентной схеме замещения [1]. Параметры этих моделей определяются при помощи полуэмпирических зависимостей, причем, чем более необычным или нагруженным является активный объем, тем более неопределенными становятся параметры этих моделей.
Возросшее число публикаций по переходным процессам в АД свидетельствует об интенсивной работе над совершенствованием старых и созданием новых методов [2], позволяющих уже на этапе проектирования точно определить характеристики двигателя при различных режимах работы, максимально возможные токи и электромагнитные моменты, их длительность и зависимость от конкретных особенностей конструкции. Существенное развитие вычислительной техники и программного обеспечения для конечно -элементных расчётов сложных научно-технических задач сегодня позволяет проводить подобные расчеты с минимальными допущениями, а главное - минимизировать потребность в использовании экспериментальных данных. Наиболее популярными на сегодняшний момент программными продуктами для расчета электрических машин методом конечных элементов, являются ANSYS и Ansoft Maxwell.
В данной работе при моделировании процесса пуска АД 4АРМАк-400/6000 УХЛ4 использовался программный продукт COMSOL Multiphysics разработки шведской компании COMSOL. Этот комплекс принципиально отличается от других программных продуктов, предназначенных для решения задач математической физики, тем, что типовые формы уравнений в частных производных обладают наибольшей общностью, и имеются неограниченные возможности их комбинирования при построении систем уравнений. Уникальной чертой пакета является способность использовать возможности MATLAB.
Уравнения электромагнитного поля в АД
Уравнения электромагнитного поля Максвелла, применительно к низкочастотным электромагнитным устройствам, где можно пренебречь током смещения и свободными зарядами, выраженные через вектор-потенциал электромагнитного поля и скалярный потенциал электрического поля, принимают вид [3]
(1)
где А - векторный магнитный потенциал, V - скалярный электрический потенциал, Зе - вектор плотности тока сторонних источников, цо - магнитная проницаемость вакуума, цг - относительная магнитная проницаемость (зависящая от модуля вектора магнитной индукции, если речь идет о нелинейных средах), ст - удельная электрическая проводимость, V - вектор скорости.
Расчет трехмерного нестационарного электромагнитного поля в активном объеме АД с учетом вращения ротора представляет собой чрезвычайно сложную задачу даже для современных персональных компьютеров. Поэтому, если предположить, что магнитное поле остается неизменным по длине активного объема и пренебречь его искажением вблизи торцевых зон, полный ток, а, следовательно, и векторный магнитный потенциал будет иметь только одну ненулевую 7 - компоненту. В этом случае задача сводится к расчету плоскопараллельного магнитного поля в поперечном сечении АД О, для которого уравнение (1) принимает вид:
где В - модуль вектора магнитной индукции, Д V - разность потенциалов на длине I.
Независимо от используемого пакета конечно -элементных расчетов процесс создания двухмер-
пространству соответственно. Поскольку пренебрежение потерями в стали шихтованного магнито-провода не оказывает существенного влияния на процесс пуска, удельную электрическую проводимость в этих областях можно принять равной нулю. В результате дифференциальное уравнение (2)
Магнитопровод как статора, так и ротора обладает шихтованной конструкцией и зачастую радиальными вентиляционными каналами. В двумерной постановке задачи неравенство относительной магнитной проницаемости в аксиальном направлении задать не представляется возможным. Поэтому учет шихтовки и наличие радиальных вентиляционных каналов необходимо проводить путем корректировки кривой цг (В). Иными словами, магнитной проницаемости цг ставить в соответствие не индукцию В , которая определяется свойствами стали при данном цг (как это было бы в массивном стальном теле), а индукцию, на ДВ меньшую, где ДВ вычисляется как: ДВ = -В(кст-1)(цг-1)/кстцг , где кст - коэффициент заполнения магнитопровода сталью. Таким
Б К
образом, для подобластей О и О вводится скорректированная кривая цг (В1) .
Подобласти ОУ~ принадлежат катушкам обмотки статора (г - номер катушки, «+» - направление намотки совпадает с направлением оси 7, «-» - не совпадает). Поскольку обычно размеры поперечного сечения элементарного проводника катушки сопоставимы с толщиной скин-слоя при промышленной частоте, влиянием вихревых токов можно пренебречь, принимая удельную электрическую проводимость для этих областей равной нулю. В результате дифференциальное уравнение
(2) для подобластей О™ примет вид
где N - число витков в катушке, Бк - суммарная площадь поперечного сечения эффективных проводников в катушке, ±іу- - мгновенное значение тока параллельной ветви фазы статора (і^ обычно является неизвестным, так как его значение зависит от Л2 и внешних элементов, связанных с об-
(2)
ной модели АД примерно одинаков. Разделим расчетную область О на подобласти ОБ , О^, О™+ , О^г, О^. Подобласти ОБ , О^ и О^ принадлежат магнитопроводу статора, ротора и воздушному
Б К 5
для подобластей О , О и О , в которых отсутствуют токи и поле является безвихревым, преобразуется к виду
где для подобластей О^ цг принимаем равной единице.
<іг^га<іЛ2 =-щ^е2 (х,у); • (X у) = ±і/^к 1 Бк,
моткой; « ± » - соответствует направлению намотки для подобласти О™±).
Области, занятые обмоткой ротора типа «беличьей клетки» П]1'Г ( і - номер стержня), должны рассматриваться с учетом влияния вихревых токов, следовательно, для этих областей уравнение (2) принимает следующий вид:
<ііу§га<і Az = цостСт
дA* ДГсл
д,
1с
(3)
где ДVСT - разность потенциалов стержня обмотки ротора на длине /ст, равной длине магнитопро-вода ротора (обычно также является неизвестной величиной), стст - удельная электрическая проводимость материала стержня ротора.
Связь уравнений электромагнитного поля с уравнениями электрических цепей в АД
Активные сопротивления и индуктивности рассеяния торцевых частей обмотки статора и ротора оказывают существенное влияние на технические характеристики АД, однако по понятным причинам не могут быть напрямую учтены в двухмерной полевой модели, что зачастую значительно занижает точность расчета. Например, для рассматриваемого в статье АД 4АРМАк-400/6000 УХЛ4 омическое сопротивление вылетов стержня ротора (рис. 2, а) составляет 60% от общего омического сопротивления стержня, что значительно влияет на пусковые токи и моменты АД. В этом случае учесть конечную длину активного объема, а также замыкание токов через торцевые части обмоток с конечным полным сопротивлением позволяет совместное решение нелинейного нестационарного дифференциального уравнения (2) для вектор-потенциала электромагнитного поля в двумерной постановке и систем дифференциальных уравнений электрической цепи для статора и ротора. Эквивалентные схемы замещения приведены на рис. 1.
Г©24-
и
к и
я -АА\— и
■е*
■&
■е5
а б
Рис. 1. Эквивалентная схема замещения для цепи статора (а) и ротора АД (б)
Параметры АД 4АРМАк-400/6000 УХЛ4 представлены в таблице Далее, при расчете, пренебре-галось индуктивностью рассеяния вылетов стержней ротора.
Параметры АД
4АРМАк-400/6000 УХЛ4
Рн, кВт Р СОБф, о. е. П, % Ї, Гц 8н , о. е. м О Ь<7 , Гн
400 1 0,908 95,06 50 0,008 1,04 0,018
а Ест, мкО Ьк, мкГн Гк , мкО ¿1 ¿2 /ст, м
15 1 40 0,64 9,4 42 34 0,43
Связь уравнений для статора
В общем случае ток параллельной ветви фазы статора гу зависит от приложенного к обмотке напряжения и, падения напряжения в активном сопротивлении Я, индуктивности рассеяния торцевой части обмотки статора и наведенной в обмотке электродвижущей силы еу (ЭДС), которую можно выразить через вектор-потенциал магнитного поля как
е/ = 1пр
где /пр - длина проводника в сердечнике статора, щ - количество катушек в параллельной ветви.
Для рассматриваемого АД, в соответствии с эквивалентной схемой замещения (рис. 1, а), можно записать систему уравнений:
uA -uB + La —А+R A -La —~— R'yB -eA + eB - 0;*A + уБ + C - 0 dt dt
T diB . die n . n
ub -ue + Lct + R'iB -; R'-eB + ee - 0-
dt dt
где R - активное сопротивление фазы обмотки статора с учетом лобовых частей и La - индуктивность рассеяния лобовых частей обмотки статора.
Индуктивность La можно определить из трехмерной задачи магнитостатики или рассчитать по известной приближенной формуле
La- 2R)/g X л w2/ pq ,
где w - число витков в фазе, p - число пар полюсов, q - число пазов на полюс и фазу, /g - условная длина поля рассеяния, X л - коэффициент магнитной проводимости лобового рассеяния. Однако при расчете по приближенной формуле не будут учитываться взаимоиндуктивные связи лобовых частей всех катушек обмотки.
Связь уравнений для ротора
В соответствии с (3) полный ток, протекающий в стержне ротора, может быть определен как
*ст - J i-aCT A+ ст-АЧd^ = - J aCTdAZdQ+^s^^ =- J aCTdAZdQ+AV^,
qwr V dt 1ст J Qwr dt 4т Qwr dt Гст
откуда разность потенциалов на длине /ст:
дА
Агст - гстіст + гст J аст -д dQ- гстіст + єст , (4)
Qwr t
где 5ст - площадь поперечного сечения стержня ротора; гст - омическое сопротивление части стержня, находящегося в пазу; ест - ЭДС, наведенная в i-м стержне.
Для короткозамкнутой обмотки ротора, рассматриваемого АД, электрическая схема замещения которой представлена на рис. 1, б, составим систему уравнений по первому и второму закону Кирхгофа:
*
di'i 2 * * * di'i 2
i1'R1 + i1,2 'r1,2 + L1,2 ~t i2 'R2 -i1,2 'r1,2 -L1,2 ^~t~- e2 -e1';i1,2 -i2,3 -i2;i2,3 -i3,i -i3;
dt dt
*
di2 3 * * * di2 3
i2 • r2 + i2,3 'r2,3 + L2,3 —T i3 'R3 -i2,3 '*2,3 -L2,3 —"T--e3 -e2';
dt dt
*
dii 1 i * * * dii 1 i
ii-1' Ri-1 + ii-1,i * ri-1,i + Li-1,i dt ii' Ri - ii-1,i * ri-1,i - L-1,i dt = ei-1 - ei; ii-1,i - i1,i - *; (5)
di1,2 dii- 1,i r dii,1 a * *
*1,2 'r1,2 + L1,2 —"7“ +... + ii-1,i'ri-1,i + Li-1,i — + ?i,1'ri,1 + Li,1 '~T~-0;*1,2 --i1 + 71,i;
dt dt dt
* * *
* * * d?12 * * * dii-1 i * * * dii 1
i1,2 'r1,2 + L1,2 —7“ +... + ii-1,i'ri-1,i + Li-1,i Г^ + *і,1'ri,1 -Li,1 '~t~-0;
dt dt dt
*2,3- -*2 + *1,2; *з,і - _i3 + *2,3; ч - *1,i- **-1,i; *'-1... Z2,
где Ri - гст + Ry - сумма гст и активного сопротивления вылетов стержня Ry; Ly- у, r- у - индуктивность рассеяния и активное сопротивление короткозамыкающего кольца на один стержень ротора.
В общем случае активные сопротивления Ry, r- у и индуктивность рассеяния Ly- у зависят
от величины и частоты протекающего тока. Поэтому сначала необходимо провести расчет трехмерной модели торцевой зоны ротора (рис. 2, а) для определения функциональных зависимостей Ry - f (іст,5), Г-1,г - f (ік ,5) и Li-1,у (ік ,s), где ік - ток короткозамыкающего кольца, 5 - скольжение.
Рассчитанные зависимости входят в систему уравнений (5) в виде интерполяционных полиномов. Полученные из решения системы уравнений (5) токи іу соответствуют полным токам стержней ротора, что позволяет, используя уравнение (4), определить разность потенциалов для каждой из областей Qwr.
а б
Рис. 2. Торцевая зона обмотки ротора (а) и картина магнитного поля в поперечном сечении АД (б)
Определение электромагнитного момента и мгновенной частоты вращения
На каждом временном слое по найденному распределению вектор-потенциала Л2 (х,у,і) в поперечном сечении активного объема АД определяем Вх, Ву компоненты вектора магнитной индукции и Нх, Ну компоненты напряженности магнитного поля во всем объеме электрической машины. Мгновенный электромагнитный момент, действующий на ротор, в соответствии с [4], легко рассчитать по найденному распределению Вх, Ву, Нх, Ну , используя тензор натяжения Максвелла
в виде
Мг = /стс))(х-ТПу - у •ТПх;
ТПх = —"2(ВХНх + ВуНу у^ПХ + {^ХНх + ПуНу )ВХ ; Тпу = —"2(ВхНх + ВуНу )у + {^хНх + пуНу )Ву,
где х и у - координаты точки интегрирования, Тп , Тп - х и у компоненты вектора натяжения
соответственно, n
ЧУ
- x и y компоненты внешней нормали к поверхности S, охватывающей
объем ротора. Мгновенную угловую скорость ротора юг находим из решения уравнения механического равновесия:
d ar dt
= ( (t)- Mс (av ))//
где Jnp - приведенный момент инерции, Мс - момент нагрузки.
Вращение моделируется с использованием технологии деформируемой конечно-элементной сетки «Arbitrary Lagrangian-Eulerian (ALE)», доступной в COMSOL Multiphysics. При этом статор и ротор представляют собой два раздельных геометрических объекта, входящих в сборку, и имеют независимые конечно-элементные сетки, что позволяет на каждом временном шаге перестраивать лишь специальный слой сетки на границе между статором и ротором. Дифференциальные уравнения для подобластей ротора решаются во вращающейся системе координат, неподвижной относительно ротора, а уравнения для подобластей статора - в неподвижной глобальной системе координат, поэтому член <rv • gradAz, учитывающий движение стержней обмотки ротора в магнитном поле в уравнении (3) отсутствует. Интерполяционные уравнения связи на границе статора и ротора обеспечивают непрерывность решения в глобальной системе координат, связанной со статором. Вращение конечно-элементной сетки ротора задается уравнением трансформации
(6)
xR cos(art) -sin(art) xS
_ yR _ sin(a rt) cos(a rt) _ ys _
где х^ , ук , х$, у$ - координаты узлов сетки ротора во вращающейся и глобальной системе координат.
Результаты моделирования процесса пуска АД
Результаты математического моделирования процесса прямого пуска АД 4АРМАк-400/6000 УХЛ4 (поперечное сечение представлено на рис. 2, б) приведены на рис. 3 и 4. Расчет выполнялся при условии пуска электродвигателя с номинальным напряжением на выводах фаз, равным 3464 В, и нагрузкой вентиляторного типа.
3500 3000 2500 2000 1500 й 1000 ^ 500 0
-500 -1000 -1500
0 0,5 1 , ч 1,5 2 2,5
/(сек)
Рис. 3. Электромагнитный момент АД и момент сопротивления ротора насоса ПЭ 150-63
Рис. 4. Ток фазы А обмотки статора АД и угловая скорость ротора
На приведенных осциллограммах можно заметить, что кратность ударного момента превышает кратность критического момента, а ударного тока - кратность пускового тока. Также отчетливо видно наличие высокочастотных составляющих в кривой электромагнитного момента, которые достигают значительных величин на этапе пуска вплоть до номинальной угловой скорости вращения. Эти составляющие обусловлены главным образом наличием пазов на статоре и роторе. Мгновенные частоты и амплитуды отдельных составляющих могут быть легко вычислены любым методом для анализа нестационарных процессов, например применяя вейвлет-преобразование или преобразование Гильберта-Хуанга [5]. Анализ приведенной кривой изменения электромагнитного момента в процессе пуска АД показал, что изменение частоты наиболее значительной высокочастотной составляющей происходит по известному закону / = кр~^12 (1 - 5) ±в; в = 2; к = 1, при этом амплитуда пульсаций изменяется от 700 до 35 Нм.
Заключение
Использование метода конечных элементов и современного программного обеспечения, в частности СОМБОЬ МиШрЬу8Ю8, позволяет выстраивать математические модели АД, обладающие высокой степенью адекватности протекающим в электрической машине физическим процессам.
Приведенный выше метод, позволяет с достаточной степенью точности описать работу АД в неустановившихся режимах, рассчитать необходимые характеристики АД (пусковые, рабочие токи и моменты; потери, полезную и потребляемую мощность) при любых видах нагрузки и напряжения питания.
Литература
1. Копылов И.П. Математическое моделирование электрических машин / И.П. Копылов. - М.: Высш. шк., 2001. - 327 с.
2. Персова М.Г. О новом подходе к проектированию электрических машин на основе численного моделирования / М.Г. Персова, Ю.Г. Соловейчик, З.С. Темлякова и др. // Электротехника. -2007. - № 9. - С. 15-21.
3. Туровский Я. Техническая электродинамика / Я. Туровский. - М.: Энергия, 1974. - 488 с.
4. Иванов-Смоленский А.В. Электромагнитные силы и преобразование энергии в электрических машинах / А.В. Иванов-Смоленский. - М.: Высш. шк., 1989. - 312 с.
5. Huang N.E. The Hilbert-Huang transform and its applications / N.E. Huang, S.S. Shen. - Singapore: World Scientific Publishing Co. Pte. Ltd. - 2005. - 325 p.
Замчалкин Александр Сергеевич
Магистр, аспирант каф. «Электромеханика»
Новосибирского государственного технического университета (НГТУ)
Тел.: (383-3) 15-39-42
Эл. почта: [email protected]
Тюков Валентин Андреевич
Д-р техн. наук, профессор каф. «Электромеханика» НГТУ
Тел.: (383-3) 15-39-42
Эл. почта: [email protected]
Zamchalkin A.S., Tyukov VA.
Numerical modeling of process of start-up of the induction motor
The possibility of calculation of dynamic modes of induction motor on the basis of the numerical solution of the nonlinear and non-stationary equation for magnetic field vector-potential in two-dimensional statement, taking into account movement of a rotor and presence of external electric circuits in software package COMSOL Multiphysics is shown. Results of numerical modeling of process of start-up of the induction motor with loading on a shaft fan-driven type are presented.
Keywords: induction motor, numerical modeling, electromagnetic field, starting mode.