CC BY
19Научный вестник НГТУ. - 2012. - № 3(48)
УДК 621.313.333
Численное моделирование процесса пуска асинхронного двигателя*
А.С. ЗАМЧАЛКИН, В.А. ТЮКОВ
Продемонстрирована возможность расчета динамических режимов асинхронных двигателей на базе численного решения нелинейного и нестационарного уравнения для вектор-потенциала магнитного поля в двумерной постановке с учетом движения ротора и наличием внешних электрических цепей в программном пакете COMSOL Multiphysics. Представлены результаты численного моделирования процесса прямого пуска асинхронного двигателя с нагрузкой на валу вентиляторного типа.
Ключевые слова: асинхронный двигатель, численное моделирование, электромагнитное поле, пусковой режим.
ВВЕДЕНИЕ
На сегодняшний день, несмотря на более чем вековую историю развития теории асинхронных двигателей (АД), до конца не изучены многие вопросы, касающиеся неустановившихся режимов работы АД. Практика электромашиностроительных предприятий ясно показывает, что существующие методы расчетов АД зачастую дают существенные расхождения с экспериментом даже при расчете простейших случаев установившегося режима. Наиболее ярко это проявляется в крупных АД с нетрадиционной геометрией активного объема и значительными электромагнитными нагрузками. Методики предприятий для расчета неустановившихся режимов, как правило, представляют собой математические модели, основанные на сложившейся в электромеханике системе допущений, которая предполагает переход к эквивалентной двухфазной системе с сосредоточенными параметрами, а не редко и к однофазной эквивалентной схеме замещения [1]. Параметры этих моделей определяются при помощи полуэмпирических зависимостей, причем, чем более необычным или нагруженным является активный объем, тем более неопределенными становятся параметры этих моделей.
Возросшее число публикаций по переходным процессам в АД свидетельствует об интенсивной работе над совершенствованием старых и созданием новых методов [2], позволяющих уже на этапе проектирования точно определить характеристики двигателя при различных режимах работы, максимально возможные токи и электромагнитные моменты, их длительность и зависимость от конкретных особенностей конструкции. Существенное развитие вычислительной техники и программного обеспечения для конечно-элементных расчётов сложных научно-технических задач сегодня позволяет проводить подобные расчеты с минимальными допущениями, а главное - минимизировать потребность в использовании экспериментальных данных. Наиболее популярными на сегодняшний день программными продуктами для расчета электрических машин методом конечных элементов, являются ANSYS и Ansoft Maxwell.
В данной работе при моделировании процесса пуска АД серии 4АРМАк-400/6000 УХЛ4 использовался программный продукт COMSOL Multiphysics 3.5a разработки шведской компании COMSOL. Этот комплекс принципиально отличается от других программных продуктов, предназначенных для решения задач математической физики, тем, что типовые формы уравнений в частных производных обладают наибольшей общностью и имеются неограниченные
* Получена 14 декабря 2011 г.
возможности их комбинирования при построении систем уравнений. Уникальной чертой пакета является способность использовать возможности МАТЬЛВ.
1. УРАВНЕНИЯ ЭЛЕКТРОМАГНИТНОГО ПОЛЯ В АД
Уравнения электромагнитного поля Максвелла, применительно к низкочастотным электромагнитным устройствам, где можно пренебречь током смещения и свободными зарядами, выраженные через векторный и скалярный потенциалы, принимают вид [3]
с— + V х (ц-Уг"V х А) - сУ х (V х А) + аУУ = Je , (1)
8 А
—+Ух(ц0-цг "УхА)-сух(УхА) + сУК =„е
где А - векторный магнитный потенциал, У - скалярный электрический потенциал, Je -вектор плотности тока сторонних источников, Цо - магнитная проницаемость вакуума, -относительная магнитная проницаемость (зависящая от модуля вектора магнитной индукции, если речь идет о нелинейных средах), с - удельная электрическая проводимость, V - вектор скорости.
Расчет трехмерного нестационарного электромагнитного поля в активном объеме АД с учетом вращения ротора представляет собой чрезвычайно сложную задачу даже для современных персональных компьютеров. Поэтому, если предположить, что магнитное поле остается неизменным по длине активного объема, и пренебречь его искажением вблизи торцевых зон, полный ток, а следовательно, и векторный магнитный потенциал будут иметь только одну ненулевую г -компоненту. В этом случае задача сводится к расчету плоскопараллельного магнитного поля в поперечном сечении АД О, для которого уравнение (1.1) принимает вид
с(В)-1 VAZ) + сУ ^ = сДУ + Je , (2)
где В - модуль вектора магнитной индукции равный В =
8А ? ( 8Аг . дтл
+1--- | , ДУ - раз-
8х
г
8у
ность потенциалов на длине I.
Независимо от используемого пакета конечно-элементных расчетов процесс создания двухмерной модели АД примерно одинаков. Разделим расчетную область О на подобласти
О^, Ой , Ог^ , О*"", О5.
Подобласти О^, ОК и О5 принадлежат магнитопроводу статора, ротора и воздушному пространству соответственно. Поскольку пренебрежение потерями в стали шихтованного маг-нитопровода не оказывает существенного влияния на процесс пуска, удельную электрическую проводимость в этих областях можно принять равной нулю. В результате дифференциальное
уравнение (2) для подобластей О^, ОК и О5, в которых отсутствуют токи и поле является безвихревым, преобразуется к виду
V•(ц-Уг (В)-1 VAz ) = 0, (3)
где для подобластей О5, принимаем равной единице.
Магнитопровод как статора, так и ротора обладает шихтованной конструкцией и зачастую радиальными вентиляционными каналами. В двумерной постановке задачи неравенство
относительной магнитной проницаемости в аксиальном направлении задать не представляется возможным. Поэтому учет шихтовки и наличие радиальных вентиляционных каналов необходимо проводить путем корректировки кривой (В). Иными словами, магнитной проницаемости ставить в соответствие не индукцию В, которая определяется свойствами стали при данном (как это было бы в массивном стальном теле), а индукцию на ЛВ меньшую, где ЛВ вычисляется как
лв = _ В (_ 1)(Цг - 1) ,
^тМг
где kст - коэффициент заполнения магнитопровода сталью.
Таким образом, для подобластей QS и ОК вводится скорректированная кривая (В1) .
Подобласти П*± принадлежат катушкам обмотки статора (I - номер катушки, «+» -
направление намотки совпадает с направлением оси г, «-» - не совпадает). Поскольку размеры поперечного сечения элементарного проводника катушки обычно сопоставимы с толщиной скин-слоя при промышленной частоте, влиянием вихревых токов можно пренебречь, принимая удельную электрическую проводимость для этих областей равной нулю. В результате дифференциальное уравнение (2) для подобластей О™ примет вид
V24 =-Цо(X, >>), (X, >>) = , (4)
где Nk - число витков в катушке, Sk - суммарная площадь поперечного сечения эффективных проводников в катушке, ±1^ - мгновенное значение тока параллельной ветви фазы статора (¡^ обычно является неизвестным, так как его значение зависит от Аг и внешних элементов, связанных с обмоткой; « ± » - соответствует направлению намотки для подобласти ). Области, занятые обмоткой ротора типа «беличьей клетки»
¡ - номер стержня),
должны рассматриваться с учетом влияния вихревых токов, следовательно, для этих областей уравнение (2) принимает следующий вид:
(сА ЛК ^
^ 4 =
г ст
с 1ст у
(5)
где ЛКст - разность потенциалов стержня обмотки ротора на длине /ст, равной длине магнитопровода ротора (обычно также является неизвестной величиной), сст - удельная электрическая проводимость материала стержня ротора.
2. СВЯЗЬ УРАВНЕНИЙ ЭЛЕКТРОМАГНИТНОГО ПОЛЯ С УРАВНЕНИЯМИ ЭЛЕКТРИЧЕСКИХ ЦЕПЕЙ В АД
Активные сопротивления и индуктивности рассеяния торцевых частей обмотки статора и ротора оказывают существенное влияние на технические характеристики АД, однако по понятным причинам не могут быть напрямую учтены в двухмерной полевой модели, что зачастую значительно занижает точность расчета. Например, для рассматриваемого в статье АД серии 4АРМАк-400/6000 УХЛ4, омическое сопротивление вылетов стержня ротора (рис. 2.2, а) составляет 60 % от общего омического сопротивления стержня, что значительно влияет на
пусковые токи и моменты АД. В этом случае, учесть конечную длину активного объема, а также замыкание токов через торцевые части обмоток с конечным полным сопротивлением позволяет совместное решение нелинейного нестационарного дифференциального уравнения (2) для вектор-потенциала магнитного поля в двумерной постановке и систем дифференциальных уравнений электрической цепи для статора и ротора. Эквивалентные схемы замещения приведены на рис. 2.1.
Рис. 2.1. Эквивалентная схема замещения для цепи статора (а) и ротора АД (б)
Параметры АД серии 4АРМАк-400/6000 УХЛ4 представлены в таблице. Далее, при расчете, пренебрегалось индуктивностью рассеяния вылетов стержней ротора.
Параметры АД серии 4АРМАк-400/6000 УХЛ4
Рн , кВт Р cosф , о. е. Л , % /, Гц , о. е. Я , Ом. La, Гн
400 1 0.908 95.06 50 0.008 1.04 0.018
Nk а Яст , мкО Lk , мкГн ^, мкО Z^ Z2 /ст, м
15 1 40 0.64 9.4 42 34 0.43
2.1. Связь уравнений для статора
В общем случае ток параллельной ветви фазы статора ¡^ зависит от приложенного к обмотке напряжения и , падения напряжения в активном сопротивлении R, индуктивности рассеяния торцевой части обмотки статора и наведенной в обмотке электродвижущей силы е^ (ЭДС), которую можно выразить через вектор-потенциал магнитного поля как
ег = 4*
X] —- ] —-dQ i=1 Sk ^ i=1 Sk ^
(6)
где / - длина проводника в сердечнике статора, щ - количество катушек в параллельной
пр ветви.
Для рассматриваемого АД, в соответствии с эквивалентной схемой замещения (рис. 2.1, а), можно записать систему уравнений:
иА - ив + ^+ Л • 1А -^ • —в - Л • - + ев = 0
Шв Шс
ив - ис + А, • — + Л • ¡в - А, • — - Л • к - ев + ес = 0 =
(7)
/А + ¡в + с = 0
где Л - активное сопротивление фазы обмотки статора с учетом лобовых частей и индуктивность рассеяния лобовых частей обмотки статора.
Индуктивность можно определить из трехмерной задачи магнитостатики или рассчитать по известной приближенной формуле
г2 ,
Ьа = 2Мю-15Хл ,
РЯ
где г - число витков в фазе, р - число пар полюсов, я - число пазов на полюс и фазу, /д - условная длина поля рассеяния, Xл - коэффициент магнитной проводимости лобового рассеяния.
Однако при расчете по приближенной формуле не будут учитываться взаимоиндуктивные связи лобовых частей всех катушек обмотки.
2.2. Связь уравнений для ротора
В соответствии с (5) полный ток, протекающий в стержне ротора, может быть определен как:
= I
(
А
дt
- + ст„
ЛУ„
\
d □ = - |
Сст ^О + ^
ст /
дt
/
ст
=-|
дt г„„
откуда разность потенциалов стержня на длине /с
дА
Л^т = Гст/ст + Гст I °ст _Д d□ = Гст/ст + ес
□Г
(8)
где 5ст - площадь поперечного сечения стержня ротора, гст - омическое сопротивление части стержня находящегося в пазу, ест - ЭДС наведенная в /-м стержне.
Для короткозамкнутой обмотки ротора, рассматриваемого АД, электрическая схема замещения которой представлена на рис. 2.2, б, составим систему уравнений по первому и второму закону Кирхгофа:
1,2
¡1 ' Я1 + ¡1,2 ' 1,2 + А,2 ' —Т--¡2 ' Я2 - ¡1,2 ' Г1,2 - Ll.
1,2
=е2-е1;
¡2 ' Я2 + ¡2,3 ' Г2,3 + —2,3
'2,3
¡3 ' Я3 ¡2,3 ' Г2,3 —2,3
2,3
= е3-е2;
¡¡-1 ' Я,-1 + ¡¡- 1,г ' Г-1,г + --1,г
V-1,г'
-¡г ' Я,-
2 2 2 M(¡^^
- ^-1,, • Г-1; - L¡-1, —— = e¡-1- e¡;
С4/1 2 ^^/¡-1 ¡ ^^/¡ 1
¡1,2 ' Г1,2 + А,2 ' —Г" + ... + ¡г-1л • Г-1,/ + А- 1л--Т^ + ¡¡,1 ' Г,1 + А,1 ~Т" = 0;
2 2 2 ¡1,2 ' Г1,2 + М,2
'1,2 .2 2 т2 ЫЧ-1,¡ .2 2 т2 "'¡,1 п
— +... + ¡¡-у • rг—1,¡ + 4-у • —— + ¡¡,1' Г,1 - Lг,l • — = 0;
¡1,2 = ¡2,3 ¡2; ¡2,3 = /3,¡ ¡3; /¡-1,¡ = /1,¡ /¡;
¡1,2 = -¡1 + ¡и ; ¡2,3 = -¡2 + ¡1,2; ¡3,¡ = "¡3 + ¡2,3; ¡¡ = ¡и - ¡¡-М; ^ = 1...Z2,
(9)
где Я = гст + Яу - сумма гст и активного сопротивления вылетов стержня Яу ; — -1 ¡, г-1 ¡ -
индуктивность рассеяния и активное сопротивление короткозамыкающего кольца на один стержень ротора.
В общем случае активные сопротивления Яу , г-1 ¡ и индуктивность рассеяния -1 ¡ зависят от величины и частоты протекающего тока. Поэтому, сначала необходимо провести расчет трехмерной модели торцевой зоны ротора (рис. 2.2, а) для определения функциональных зависимостей Яу = / (¡ст, 5) , г-ц = / (¡к, 5) и L¡-l¡ (¡к, 5), где ¡к - ток короткозамыкаю-щего кольца, 5 - скольжение.
2
б
а
Рис. 2.2. Торцевая зона обмотки ротора (а) и картина магнитного поля в поперечном сечении АД (б) Рассчитанные зависимости входят в систему уравнений (9) в виде интерполяционных полиномов. Полученные из решения системы уравнений (9) токи // соответствуют полным токам стержней ротора, что позволяет, используя уравнение (8), определить разность потенциалов для каждой из областей ОГГ.
3. ОПРЕДЕЛЕНИЕ ЭЛЕКТРОМАГНИТНОГО МОМЕНТА И МГНОВЕННОЙ ЧАСТОТЫ ВРАЩЕНИЯ
На каждом временном слое по найденному распределению вектор-потенциала Аг ( х, у, t) в поперечном сечении активного объема АД определяем вх -, ву -компоненты вектора магнитной индукции и Нх -, Ну -компоненты напряженности магнитного поля во всем объеме
электрической машины. Мгновенный электромагнитный момент, действующий на ротор, в соответствии с [4], легко рассчитать по найденному распределению вх, ву, Нх, Ну, используя тензор натяжения Максвелла в виде
Mz = /сТX• Tny -y• Tnx)dS,
s у '
Tnx = - 2 (BXHX + ByHy) nx + (nxHx + nyHy) Bx Tny = - 2 (BxHx + ByHy) ny +(nxHx + nyHy) By
(10)
ny 2 у x x
где x и y - координаты точки интегрирования, Tn , Tn - x - и y -компоненты вектора на-
x y
тяжения соответственно, nx, ny - x - и y -компоненты внешней нормали к поверхности S
охватывающей объем ротора.
Мгновенную угловую скорость ротора юг находим из решения уравнения механического равновесия
d^ =(Mz (t)-Мст) (11)
dt JuP '
где Jпр - приведенный момент инерции, Мст - момент нагрузки.
Вращение моделируется с использованием технологии деформируемой конечноэлемент-ной сетки «Arbitrary Lagrangian-Eulerian (ALE)» доступной в COMSOL Multiphysics. При этом статор и ротор представляют собой два раздельных геометрических объекта, которые входят в сборку и имеют независимые конечноэлементные сетки, что позволяет на каждом временном шаге перестраивать лишь специальный слой сетки на границе между статором и ротором. Дифференциальные уравнения для подобластей ротора решаются во вращающейся системе координат, неподвижной относительно ротора, а уравнения для подобластей статора - в неподвижной глобальной системе координат, поэтому член cv • VAz, учитывающий движение
стержней обмотки ротора в магнитном поле в уравнении (5), отсутствует. Интерполяционные уравнения связи на границе статора и ротора обеспечивают непрерывность решения в глобальной системе координат, связанной со статором. Вращение конечноэлементной сетки ротора задается уравнением трансформации
(12)
где хк , ук, х^, у5 - координаты узлов сетки ротора во вращающейся и глобальной системе координат.
хя соз(ю,1) - )
_ Уя _ _ _
4. РЕЗУЛЬТАТЫ МОДЕЛИРОВАНИЯ ПРОЦЕССА ПУСКА АД
Результаты математического моделирования процесса прямого пуска АД серии 4АРМАк-400/6000 УХЛ4 (поперечное сечение представлено на рис. 2.2, б представлены на рис. 4.1 и рис. 4.2. Расчет выполнялся при условии пуска электродвигателя с номинальным напряжением на выводах фаз равным 3464 В и нагрузкой вентиляторного типа.
Рис. 4.1. Электромагнитный момент АД 4АРМАк-400/6000 УХЛ4 при прямом пуске и момент сопротивления ротора насоса ПЭ 150-63
Рис. 4.2. Ток фазы А обмотки статора АД 4АРМАк-400/6000 УХЛ4 при прямом пуске и угловая
скорость ротора
На приведенных осциллограммах можно заметить, что кратность ударного момента превышает кратность критического момента, а ударного тока - кратность пускового тока. Также отчетливо видно наличие высокочастотных составляющих в кривой электромагнитного момента, которые достигают значительных величин на этапе пуска вплоть до номинальной угловой скорости вращения. Эти составляющие обусловлены главным образом наличием пазов на статоре и роторе. Мгновенные частоты и амплитуды отдельных составляющих могут быть легко вычислены любым методом для анализа нестационарных процессов, например, применяя вейвлет-преобразование или преобразование Гильберта - Хуанга [5]. Анализ приведенной кривой изменения электромагнитного момента в процессе пуска АД показал, что изменение частоты наиболее значительной высокочастотной составляющей происходит по известному закону
( Z Л
/,= ^(1-5)±р ; р = 2; k = 1,
I Р )
при этом амплитуда пульсаций изменяется от 700 до 35 Нм.
ВЫВОДЫ
Использование метода конечных элементов и современного программного обеспечения, в частности COMSOL Multiphysics, позволяет выстраивать математические модели АД, обладающие высокой степенью адекватности протекающим в электрической машине физическим процессам.
Приведенный выше метод позволяет с достаточной степенью точности описать работу АД в неустановившихся режимах, рассчитать необходимые характеристики АД (пусковые, рабочие токи и моменты; потери, полезную и потребляемую мощность) при любых видах нагрузки и напряжения питания.
СПИСОК ЛИТЕРАТУРЫ
[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 / editors, N.E. Huang, S.S.P. Shen. - World Scientific Publishing Co. Pte. Ltd. 5 Toh Tuck. Link, Singapore 596224.
Замчалкин Александр Сергеевич, магистр, аспирант кафедры электромеханика Новосибирского государственного технического университета. Основное направление научных исследований - переходные процессы в крупных асинхронных машинах. Имеет 10 публикаций. E-mail: zamchalkin-alex@mail.ru.
Тюков Валентин Андреевич, доктор технических наук, профессор кафедры электромеханика Новосибирского государственного технического университета. Основное направление научных исследований -вопросы теории электромеханического преобразования энергии. Имеет более 1 00 публикаций. E-mail: tyukov@em.power.nstu.ru.
A.S. Zamchalkin, V.A. Tyukov
Numerical modeling ofprocess 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.
Key words: induction motor, numerical modeling, electromagnetic field, starting mode.
