VAC 05.17.00 Химическая технология
VAC 05.23.00 Строительство и архитектура
GRNTI 20.53 Технические средства обеспечения информационных процессов
Is developed mathematical model of hysteresis loop H=f(B) kind. Method of formatting partly asymmetric cycle of magnetization reversal is proposed. There is realized model of three-phase transformer hysteresis loop in the Matlab Simulink system.
hysteresis loop, three-phase transformer, symmetrical and asymmetrical cycle of magnetization reversal
Моделирование петли гистерезиса, а так же частных петель для силовых трансформаторов и нелинейных индуктивностей в программном комплексе Matlab Simulink возможно с помощью блока Powergui. Однако получающиеся при этом частные петли по форме значительно расходятся с фактическими. На рис. 1, а представлено семейство частных симметричных петель, которые при изменении напряженности магнитного поля (тока) практически не изменяют своей ширины, в отличие от фактических петель, показанных на рис. 1, б. [3]
|
I, А
|
|
20 15 10 5 0 -5 -10 -15 -20
|
|
|
|
-1.5 |
|
|
|
-0.5 |
|
0 |
|
0.5 |
|
1 |
|
1.5 |
|
2 |
|
а б
Рис. 1. Семейство частных петель гистерезиса (а) при моделировании, (б) фактические
Следовательно, для исследований некоторых специальных режимов работы трансформаторов необходимо использовать более корректную модель гистерезиса. Причем эта модель должна формировать как частные симметричные, так и несимметричные циклы перемагничивания. На сегодняшний день таким требованиям отвечают модели Джилса – Аттертона, Джона Чана и др. [2]. Наиболее простой для реализации является модель Джона Чана, так как она использует только основные параметры петли гистерезиса, такие как коэрцитивная сила (Hc), остаточная индукция (Br) и индукция насыщения (Bs). Эти параметры для электротехнических сталей (ЭТС) легко найти, например, в справочнике [1].
Разработка гистерезисной модели сердечника для программного комплекса Matlab Simulink, правильно реализующей частные петли гистерезиса, является актуальной задачей.
В модели гистерезиса Джона Чана верхняя и нижняя кривые петли предельного цикла перемагничивания описываются соответственно уравнениями 1 и 2 [2].
, (1)
, (2)
где – магнитная проницаемость вакуума, Гн/м.
Вторым слагаемым в выражениях 1 и 2 можно пренебречь в силу малых значений.
Начальная кривая намагничивания (КН) при моделировании вычисляется как половина суммы верхней и нижней кривых по выражению 3:
. (3)
Таким образом, предельный цикл перемагничивания имеет вид, представленный на рис. 2.
|
Hc |
|
Bs |
|
Br |
|
B |
|
H |
|
Bup(H) |
|
Bdn(H) |
|
Bmed(H) |
Рис. 2. Петля гистерезиса и начальная кривая намагничивания
Однако в действительности предельный цикл – это один случай перемагничивания и чаще всего трансформатор работает по частным циклам, например, в режиме холостого хода получается симметричный частный цикл. В рассмотренной в статье [2] методике предлагается строить частные циклы путем переноса нижней и верхней ветвей предельного цикла соответственно вверх и вниз. При этом формируется траектория D-E-F-G-F, представленная на рис. 3, а. Участок DE получается путем переноса нижней предельной ветки на величину Bd. Для моделирования КН были выбраны параметры ЭТС марки 3414 [1].
Частный цикл перемагничивания будет симметричным при 0 ≤ Bd≤ Br. При Bd>Br участок D-E должен пройти выше точки A начала координат. В этом случае перенос нижней ветви предельного участка вверх приведет к ломаной несимметричной частной петле, так как точка E окажется на верхней ветви. Для получения ровного несимметричного цикла, когда, например, происходит изменение H от -40 до +120 А/м (рис. 3, б), была разработана методика, которая отличается от предложенной в статье [2] более удобной реализацией.
При увеличении напряженности до величины H = 0 (участок D-E) нижняя ветвь частного цикла формируется простым переносом нижней ветви предельного цикла как в предыдущем случае. Затем в точке H = 0 необходимо зафиксировать величину текущей индукции B. При дальнейшем увеличении H до +120 А/м необходимо сохранять пропорцию между расстояниями до предельной петли и до средней (начальной) кривой. Согласно рис. 3, б получается пропорция 4 из соответствующих отрезков, а кривые при изменении H от 0 до + Hmax и до - Hmax описываются соответственно уравнениями 5 и 6:
. (4)
, (5)
, (6)
где - значение индукции при H = 0.
Таким образом, образуется частный цикл по траектории D-E-K-Br-D.
|
|
|
|
A |
|
C |
|
D |
|
E |
|
F |
|
G |
|
K |
|
-120 |
|
H, А/м
|
|
1.5 1.2 0.9 0.6 0.3 0 -0.3 -0.6 -0.9 -1.2 -1.5
|
|
-100 |
|
-80 |
|
-60 |
|
-40 |
|
-20 |
|
0 |
|
20 |
|
40 |
|
60 |
|
80 |
|
100 |
|
120 |
|
|
Br |
|
E |
|
K |
|
A |
|
C |
|
D |
|
Br' |
|
E' |
|
A' |
|
H, А/м
|
|
-100 |
|
-80 |
|
-60 |
|
-40 |
|
-20 |
|
0 |
|
20 |
|
40 |
|
60 |
|
80 |
|
100 |
|
120 |
|
-120 |
|
1.5 1.2 0.9 0.6 0.3 0 -0.3 -0.6 -0.9 -1.2 -1.5
|
а б
Рис. 3. (а) симметричный, (б) несимметричный частный цикл
Перейдем к реализации рассмотренной модели петли гистерезиса в модели силового трехфазного трансформатора в Matlab Simulink. Из множества библиотечных моделей наиболее подходящей для внедрения, а так же для исследований, связанных с магнитной системой трансформатора, является «Three-phase Two-Winding Saturable Transformer – Core Type». В отличие от других моделей здесь обмотки трансформатора разных фаз магнитосвязаны. Особенностью структуры этой модели является реализация КН в виде функции H=f(B). То есть, для применения модели петли гистерезиса, необходимо получить из имеющихся уравнений B(H) уравнения вида H(B).
Воспользовавшись математическим вычислительным комплексом Maple, после простых преобразований уравнений 1 и 2 имеем:
, (7)
. (8)
Решением задачи по преобразованию начальной КН, описанной уравнением 3, являются восемь корней, три из которых 9, 10, 11 определяют ее характер:
, (9)
, (10)
. (11)
Причем каждое из полученных уравнений действует на своем интервале B, то есть:
(12)
Здесь Bmed(Hc) величина магнитной индукции, определяемая по формуле 3 при H = Hc:
. (13)
Общий вид обратной петли гистерезиса и начальной КН представлен на рис. 4, а. Отдельный вид начальной КН, образованной тремя графиками, представлен на рис. 4, б.
|
Hc |
|
Bs |
|
Br |
|
B |
|
H |
|
Hmed1(B) |
|
Hmed2(B) |
|
Hmed3(B) |
|
Hdn(B) |
|
Hup(B) |
|
Bmed(Hc) |
|
-Bmed(Hc) |
|
B |
|
H |
|
Hmed1(B) |
|
Hmed2(B) |
|
Hmed3(B) |
|
Bmed(Hc) |
|
-Bmed(Hc) |
а б
Рис. 4. Общий вид (а) обратной петли гистерезиса, (б) начальной кривой
Теперь, чтобы получить частный симметричный цикл перемагничивания необходимо при увеличении индукции B сместить влево кривую Hdn(B) на величину Bd, а при уменьшении B сместить вправо кривую Hup(B) на ту же величину. В результате формируется частный симметричный цикл перемагничивания по траектории D-E-F-G-F (рис. 5, а).
Для получения частных несимметричных петель применяется аналогичная методика, как и с прямой петлей. Однако теперь необходимо вычислять H. Поэтому когда индукция B увеличивается и H меньше 0 кривая Hup(B) просто смещается вправо (участок DE на рис. 5, б), а когда H начинает входить в интервал от 0 до +Hmax, то необходимо соблюдать соотношение 4. В итоге формируется цикл по траектории D-E-K-Br-D (рис. 5, б).
|
A |
|
D |
|
E |
|
F |
|
G |
|
K |
|
C |
|
Bd |
|
|
-1.5 |
|
B, Тл
|
|
120 100 80 60 40 20 0 -20 -40 -60 -80 -100 -120
|
|
-1.2 |
|
-0.9 |
|
-0.6 |
|
-0.3 |
|
0 |
|
0.3 |
|
0.6 |
|
0.9 |
|
1.2 |
|
1.5 |
|
Br |
|
E |
|
A |
|
C |
|
D |
|
Br' |
|
E' |
|
A' |
|
K |
|
-1.5 |
|
B, Тл
|
|
-1.2 |
|
-0.9 |
|
-0.6 |
|
-0.3 |
|
0 |
|
0.3 |
|
0.6 |
|
0.9 |
|
1.2 |
|
1.5 |
|
|
120 100 80 60 40 20 0 -20 -40 -60 -80 -100 -120
|
а б
Рис. 5. (а) Симметричный, (б) несимметричный частный цикл обратной петли
Получим уравнения для кривой при несимметричном цикле, когда величина H изменяется от 0 до +Hmax. Для этого выразим переменную H из уравнений 5 и 6. При преобразовании каждого из уравнений появляется по восемь корней. Однако только два определяют характер кривой. Причем уравнения 14 и 15 описывают поведение кривой в первом квадранте, а 16 и 17 – в третьем квадранте:
, (14)
, (15)
, (16)
, (17)
Общий вид части обратной петли гистерезиса и начальной КН в первом квадранте, а так же кривая частного несимметричного цикла, образованная двумя графиками по уравнениям 14 и 15, при изменении H от 0 до +Hmax, представлены на рис. 6.
|
Hc |
|
Bs |
|
Br |
|
B |
|
H |
|
Hmed(B) |
|
Bright01 |
|
Hright01 |
|
Hright1(B) |
|
Hdn(B) |
|
Hup(B) |
|
Hright2(B) |
|
A |
Рис. 6. Петля гистерезиса и частный несимметричный цикл в первом квадранте
Как видно из рис. 6 графики Hright1(B) и Hright2(B) пересекаются в точке А с координатами Bright01 и Hright01. Следовательно, общий график Hright(B) изменяется сначала по кривой Hright1(B), затем в точке А переходит на кривую Hright2(B). Координаты точки А в первом и третьем квадранте описываются соответственно уравнениями 18, 19, 20, 21:
, (18)
, (19)
, (20)
. (21)
Полученную модель обратной петли гистерезиса можно представить в виде S – функции, описанной в Matlab. Для ее применения выбранный трансформатор необходимо преобразовать, убрав блоки с КН, а вместо них поставить блоки с S – функцией.
Рассмотрим подробнее работу модели трехфазного трансформатора с разработанной петлей гистерезиса. Помимо B-H характеристики стали в данной модели учитывается геометрия магнитной цепи. Структура модели блока «Three-phase Two-Winding Saturable Transformer – Core Type» представлена на рис. 7. Выводы первичных обмоток обозначены индексом 1, а вторичных – 2. Магнитная связь между обмотками, которые располагаются на едином сердечнике, моделируется в подсистеме «Core». На вход подсистемы подаются значения магнитных потоков «Fluxes», а на выходе получаются магнитодвижущие силы МДС. Магнитные потоки в такой модели трансформатора вычисляются из напряжений «Voltages», которые снимаются через датчики напряжений «V», подключенных к сопротивлениям Rm. Для управления источниками тока «Imag», которые моделируют ток намагничивания в обмотках, полученные МДС преобразуются в токи «Currents».
Рис. 7. Структура модели трехфазного трансформатора
В подсистеме «Core», которая представлена на рис. 8, применяются подсистемы (блоки магнитного сопротивления) нелинейной магнитной цепи, состоящей из трех стержней («L1_1», «L1_2», «L1_3») и четырех частей ярем («L2_1», «L2_2», «L2_3», «L2_4»), а так же блоки воздушных элементов «Air», моделирующие потоки рассеяния и поток нулевой последовательности. Установленные датчики напряжения «FMM» измеряют величину МДС каждой из шести катушек, а результаты поступают на блок «State - Space», который решает дифференциальные уравнения обобщенной электрической машины.
Рис. 8. Структура подсистемы «Core»
В подсистемах нелинейной магнитной цепи стандартной модели «Three-phase Two-Winding Saturable Transformer – Core Type» используется нелинейный блок «Look – Up Table», в котором задается B – H характеристика КН трансформатора, рис. 9, а. Разработанная модель петли гистерезиса в виде S–функции представлена на рис.9,б блоком «Function1».
а б
Рис. 9. Структура подсистемы нелинейной магнитной цепи (а) с КН, (б) с петлей гистерезиса
Произведем количественную оценку работы моделей трансформаторов с гистерезисом и без него. Примем к рассмотрению силовой трансформатор типа ТМЭ–400/10 со следующими параметрами: Sн = 400 кВА, Uн1 = 6300 В, Iн1 = 38.49 А, Iхх = 1.52 %, ΔPхх = 820 Вт, ΔPкз = 6062 Вт. Схема и группа соединения обмоток D/Yн–11. Магнитопровод трансформатора выполнен из стали ЭТ-А-3408 ГОСТ 21427.1-83. Длины средних линий магнитопровода – междуосевое расстояние по соседним стержням L2 = 325 мм, по ярмам L1 = 635 мм. Площадь сечения магнитопровода S = 225.264 см2. Количество витков обмотки НН w2 = 27, ВН w1 = 425.
Рассмотрим пуск трансформатора на холостом ходу. Структурная схема модели представлена на рис. 10. Здесь модель трансформатора представлена подсистемой «D/Y», которая подключается к блоку «Трехфазный источник напряжения» через блок «Выключатель». Для пуска трансформатора подается сигнал с блока «Задание» на включение выключателя. Величина и форма фазных токов измеряется в блоке «Трехфазный измеритель U и I» и выводятся на блок «Осциллограф». Результаты измерений для модели трансформатора без петли гистерезиса и с петлей представлены на рис. 11 а и б.
Рис. 10. Структурная модель пуска трансформатора
|
|
t, c
|
|
150 100 50 0 -50 -100 -150
|
|
0 |
|
0.2 |
|
0.4 |
|
0.6 |
|
0.8 |
|
1 |
|
1.2 |
|
1.4 |
|
1.6 |
|
1.8 |
|
2 |
|
IA |
|
IB |
|
IC |
|
400 300 200 100 0 -100 -200 -300 -400
|
|
t, c
|
|
0 |
|
0.2 |
|
0.4 |
|
0.6 |
|
0.8 |
|
1 |
|
1.2 |
|
1.4 |
|
1.6 |
|
1.8 |
|
2 |
|
IA |
|
IB |
|
IC |
|
а б
Рис. 11. Осциллограммы токов в обмотках ВН в модели (а) без петли, (б) с петлей
Как видно из рис. 11, б бросок токов намагничивания при включении оказывается гораздо больше, чем в модели без петли гистерезиса, рис. 11, а. Для фазы А их величина составляет соответственно 375 А и 143 А. Помимо этого видно, что затухание токов в модели с гистерезисом происходит явно быстрее. За время 0.6 с ток на рис. 11, б изменяется от 375 А до 50 А (в 7,5 раз), в то время как ток на рис. 11, а за то же время изменяется от 143 до 50 А (в 3 раза). Это объясняется тем, что наличие петли гистерезиса увеличивает потери мощности на перемагничивание.
Об адекватности разработанной модели силового трансформатора свидетельствуют данные, полученные экспериментальным путем. В ходе эксперимента комплектное распределительное устройство КРУ подключало сразу два трансформатора ТМЭ-400/10 со схемой соединения обмоток D/Yн и Y/Yн. Осциллограммы суммарных фазных токов первичных обмоток снимались с трансформаторов тока (ТТ), установленных в КРУ. Результаты измерения для фазы А представлены на рис. 12.
|
400 320 240 160 80 0 -80 -160
|
|
|
t, c-2 |
|
0 |
|
0.5 |
|
1 |
|
1.5 |
|
2 |
|
2.5 |
|
3 |
Рис. 12. Осциллограмма тока в фазе А при эксперименте
Аналогично эксперименту структурная модель пуска (рис. 10) была дополнена вторым трансформатором. Для сравнения на рис. 13, б представлена осциллограмма такой же фазы А, как и в эксперименте. Как известно ТТ не пропускают постоянную составляющую измеряемого тока, поэтому на осциллограмме рис. 12 наблюдаются отрицательные полуволны, в отличие от смоделированных осциллограмм на рис. 13, б, где использовались датчики тока ДТ. Для оценки адекватности было проведено восстановление кривой первичного тока с помощью интегрирующего - форсирующего звена с равными постоянными времени, результат приведен на рис. 13, а.
|
|
t, c |
|
0 |
|
0.5 |
|
1 |
|
1.5 |
|
2 |
|
2.5 |
|
3 |
|
400 320 240 160 80 0 -80
|
|
t, c |
|
400 320 240 160 80 0 -80
|
|
|
0 |
|
0.5 |
|
1 |
|
1.5 |
|
2 |
|
2.5 |
|
3 |
а б
Рис. 13. Осциллограммы тока в фазе А (а) восстановленного, (б) смоделированного
Как видно на рис. 13 в обоих случаях за время 0,5 с ток затухает до величины 70 А. Отличие состоит только в первом броске тока. Так на рис. 13, а бросок тока соответствует 400 А, когда на рис. 13, б он достигает величины 320 А. Это объясняется тем, что применяемый ТТ в первый момент времени насыщается. Для количественной оценки всего переходного процесса были посчитаны амплитуды токов за 16 периодов и построены затухающие экспоненциальные зависимости без первого пульса. В результате оценки средних значений было выявлено, что экспериментальный ток затухает на 5.5 % быстрее, чем смоделированный, что говорит об адекватности имитационной модели.
Выводы:
- Разработана модель обратной петли гистерезиса по преобразованным уравнениям Джона Чана и др. для сердечников трансформаторов и нелинейных индуктивностей;
- Разработана и предложена методика формирования частных симметричным и несимметричных циклов перемагничивания стали сердечника;
- Получены рекомендации по применению модели петли гистерезиса в существующих моделях трансформаторов программного комплекса Matlab;
- Доказана адекватность математической модели.
1. Molotilov B.V., Mironov L.V., Petrenko A.G. i dr. Holodnokatannye elektrotehnicheskie stali: spravochnik. M.: Izd. Metallurgiya, 1989. 168 s.
2. Volodin V. Gisterezisnaya model' nelineynoy induktivnosti simulyatora LTspice // Silovaya elektronika. 2010. №1. S. 56-60.
3. Mihaylova M. Yu., Yakovlev A. O., Nesterov A. M. Teoreticheskie osnovy elektrotehniki: metodicheskie ukazaniya k vypolneniyu lab. rabot dlya studentov special'nosti 140211-Elektrosnabzhenie. Belgorod: iz-vo. BGTU im. V. G. Shuhova, 2012. 54 s.
