УДК 519.872: 519.622.2
С. М. Тиховод, И. О. Афанасьева, Т. М. Корнус
Метод компьютерного моделирования установившихся периодических электромагнитных процессов
Разработан численный метод расчета установившихся периодических несинусоидальных процессов в электрических цепях. В основу методики расчета положена аппроксимация производной решения дифференциальных уравнений при помощи полиномов. Метод проверен на конкретных примерах.
В настоящее время до 40% вырабатываемой электроэнергии перед потреблением подвергается пре-образованию[1]. Тиристорные и транзисторные преобразователи, как правило, вырабатывают питающие напряжения, отличающиеся по форме от синусоидальных. Уточненный расчет установившихся периодических электромагнитных процессов для различных нагрузок при питании от преобразователей с несинусоидальной формой выходного напряжения востребован практикой и актуален.
В настоящее время основным и наиболее распространенным методом расчета установившихся несинусоидальных процессов в электрических цепях является метод, основанный на разложении несинусоидальных электродвижущих сил (ЭДС) в ряд Фурье и расчете всех гармоник тока, вызванных по отдельности соответствующими гармониками ЭДС, а затем - на последующих расчетах комплексным (символическим) методом [2]. В этом случае несинусоидальные токи и напряжения выражаются суммами бесконечных тригонометрических рядов и не всегда могут быть вычислены без компьютерной техники. При другом известном способе [2] установившийся процесс находят численным методом [3] путем расчета переходного процесса до его полного затухания. Этот способ при медленно затухающих переходных процессах требует значительного времени расчета, что дополнительно приводит к значительной накапливаемой ошибке. В работе [4] предложен метод расчета установившегося режима путем вычислений переходного процесса со специальным подбором начальных условий, где начальные условия приближаются к оптимальным условиям путем многократного итерационного процесса. Это также требует значительных затрат компьютерного времени.
Целью данной статьи является разработка более быстродействующего (по сравнению с известными методами) численного метода расчета установившихся периодических несинусоидальных токов и напряжений в электрических цепях.
Суть метода изложим на простом примере. Рассмотрим во временной области ^ < t < tN линейное дифференциальное уравнение с постоянными коэффициентами, которое получено на основании закона напряжений Кирхгофа для одноконтурной цепи, содержащей последовательно соединенные активное сопротивление И, индуктивность L и источник ЭДС в^):
т di
L — + Ri = e(t) dt
(1)
где i(t) - функция изменения тока, протекающего через индуктивный элемент (т. е. i(t) - переменная состояния электрической цепи).
Допустим, что производную решения можно аппроксимировать полиномом p(t) (/У-1)-ой степени:
di 2 N_1
i (t) =— = p(t) = a1 + a2t + a3t +... + aNt (2) dt
В работе [5] такое допущение выполнено для частного случая при N = 4. Рассмотрим общий случай, когда N - произвольное целое число. Интервал изменения аргумента разобьем на N - 1 отрезков точками t1, t2 ,... tN. Для аппроксимирующего полинома (2) зададим дополнительное условие, чтобы в точках tk деления интервала изменения аргумента выполнялось равенство:
i (tk) = P(tk),
(3)
где к = 1, 2,...М - номер опорной точки.
Если условие (3) записать для каждой опорной точки tк , то получим систему линейных алгебраических уравнений порядка N:
a1 + a2t1 + a3t12 +... + aNt1N 1 = i' (t1); a + a2t2 + a3t22 +... + a^2N _1 = i' fo);
aN + a2tN + a3tN + + aNtN = i (tN )
В матричной форме система (4) имеет вид:
= V • A = I'
(4)
1 tl t12 • • tN_1 1 '1 a1 " i'(t1) "
1 t t22 • • t n_1 a2 = i'(t2)
1 tn tn' • ■ t n_1 'n . .an. i'' (tn ).
(5)
где
© С. М. Тиховод, И. О. Афанасьева, Т. М. Корнус 2009 р.
V - матрица Вандермонда,
A =
aL " i'(ti) "
. I =
_aN _ i ('n )_
Пусть номер к-го отрезка, на которые разделен рассматриваемый интервал изменения аргумента, совпадает с номером точки деления tk, расположенной слева отрезка. Проинтегрируем выражение (3) на к-ом отрезке.
ik+i - ik =
'k+1 = | p(t)dt
(6)
tk
Li[ + Ri1 = e(t1); Li'2 + Ri2 = e(t2);
(10)
Представим систему (10) в матричной форме:
L
L
" ii "
i2
+
_iN-1 _
R
R
R
Подстановка в (6) выражения (2) после интегрирования дает:
'к+1 _'к - у ('к+1 _ Ч ) + у(' к+1 _1 к ) +
i1 " e(ti) "
i2 = e('2)
_iN-1 _ e(tN-i)_
(11)
a
, "N itN ,N \
+--(t k+1 - t k )
N
(7)
где к = 1,..., N - 1.
Если уравнение (7) записать для всех N точек, то получим систему алгебраических уравнений:
'2 - '1 =^2 - t1) + - t12) + ••• + ol- (t2N - tn );
i
El 1
al(* 2 t 2) 2
a2 (t 2 - t 2) 2
n
a1 a2 2 2 an n n
l3 - l2 = ~r (t3 - t2) + (t3 -12 ) + L (t3 - t2 );
ln ln-1 = j (tN tN-1) + 2(tn tn-1) + + nn (tN tn-1)
В матричной форме система (8) имеет вид:
.(8)
-1 1 -1 1
-1 1
t2 -11
t2 -tl
t2 -t2 '3 '2
tN 'N-1 'n 'N-1
+N tN '2 - 'l
N
NN t2 - t1
N
tN - tN
tN tN -1
N
(9)
Если для N-1 точек также записать исходное дифференциальное уравнение (1), то получим следующую систему линейных уравнений:
Введем в рассмотрение вектор размерности, равной 3!:
У -О1, а2,---аы,'2 ,'1, 7, (12) который можно записать в виде:
Yi =
(13)
где вектор I имеет вид:
(14)
Используя вектор неизвестных У1, объединим уравнения (5), (9) и (11) в одну систему уравнений:
Mi • Yi = Fi
(15)
Коэффициенты указанной системы уравнений приведем в табл. 1. Для рассматриваемого установившегося процесса последние две строки в табл. 1 соответствуют условию периодичности получаемого решения:
i(ti) = i(tN); i'(ti) = i' ('n )
(16)
В результате этого получим систему из 3N уравнений с количеством неизвестных, равным 3М В качестве единственного решения для этой системы уравнений определим значения коэффициентов а1, в2,...,в!
X
L
X
2
I
N
2
N
1
2
a
'3 - '2
a
2
a
N
1
2
Таблица 1. Структура матриц уравнения (15)
A I' I Правая часть
а1 а2 aN i'1 i N i'1 iN
1 t1 t N -1 t\ -1
1 t2 t N -1 2
1 tN t N -1 N -1
t2 - tx 12 - 12 »2 ^ tN tN t2 t1 1 -1
1 2 N
t3 — t2 12 -12 l3 2 tN tN t3 2 1
1 2 N
tN — tN-1 12 -12 lN 'N-1 tN tN N N-1 -1
1 2 N
L R e(ti)
e(t2)
L R e(tN-i)
-1 1 0
-1 1
гов интегрирования, что привело бы, в свою очередь, к формированию полиномов высокого порядка. Это могло бы вызвать потерю точности аппроксимации решения полиномами. Поэтому ограничимся только рассмотрением решения для установившегося периодического процесса.
Рассмотрим частный случай, когда длины всех отрезков одинаковы и равны h, т. е. при: tk+1 - tk = h. В этом случае структура матрицы, приведенной в табл. 1, примет вид, показанный в табл. 2. Если коэффициенты (сопротивление R и индуктивность L) дифференциального уравнения (1) являются функциями от тока и времени L(t, i), R(t, i), то уравнение (1) является нелинейным. В этом случае матричное уравнение (15) следует решать, используя итерационные вычислительные процессы.
Обратим внимание на то, что решение систем дифференциальных уравнений производится аналогичным образом, как и решение одного уравнения. Для каждой искомой функции переменной состояния
аппроксимирующего полинома (2), значения искомой функции ¡(() и ее производной ¡'(() в опорных точках.
Если в уравнении (7) верхний предел интегрирования принять переменным параметром t, то получим следующее общее решение для искомого процесса:
'(I) ='к + - %) + у С2 -12 к) +
+ ^2 -'2к)• + ^(ГМ -ГМк) (17)
3 N К '
Таким образом, определив коэффициенты полинома и значения тока в опорных точках, получим решение во всех произвольных промежуточных точках любого из N отрезков в области изменения независимой переменной t.
Заметим, что предложенным способом можно рассчитывать не только установившиеся процессы, но и переходные. Однако, в таком случае, как правило, пришлось бы использовать большое число ша-
Таблица 2. Структура матриц уравнения (15) при постоянном шаге интегрирования
A I' I Правая часть
а1 а 2 aN i'1 i N h i2 iN
1 0 0 0 -1
1 h hN
1 (2h)2 (2h)N
1 ((N-1)h)1 ((N-1)h)N -1
h h2 2 hN N 1 -1
h 3 £ 2 (2h)N - hN N 1 -1
h (2hf - hN N [(N - 1)h] -[(N - 2)h] N 1 — 1
L R e(h)
R e(2h)
L R e((N-1)h)
-1 1 0
-1 1
моники выше третьей, то приемлемой погрешности нельзя добиться ни при каких значениях N.
Для устранения этого недостатка предложено использовать аппроксимирующий полином не для всего периода, а только для его части. Для этого период разбиваем на К временных участков и на каждом участке записываем матричное уравнение, аналогичное (15). Для узлов, в которых граничат смежные участки, запишем уравнения связи:
¿(Кк+1) = , к); 1
¿•Си+1)=¿'(^к)], (18)
где первый индекс обозначает номер узловой точки на участке, а второй индекс - номер участка. Таким образом, уравнения (18) задают условия, что в узлах смежности участков приравниваются соответственно токи и их производные, определенные в конце и в начале для смежных участков. При этом условие периодичности задает, что при установившемся периоди-
задается своя подматрица М.. и подвектор правой части Р.
Для исследования разработанного метода разработана компьютерная программа и выполнен расчет установившегося режима в одноконтурной цепи, содержащей последовательно соединенные элементы И, ^ и источник ЭДС вЦ). Расчет выполнен предложенным численным методом и для сравнения - точным комплексным методом. Установлено, что для случая, когда ЭДС источника содержит первую и третью гармоники одинаковой амплитуды, результаты расчета предложенным численным и комплексным методами совпадают при выбранном количестве опорных точек N е"16 на период первой гармоники. При увеличении N до 45 погрешность расчета уменьшается, однако при дальнейшем увеличении количества опорных точек N > 45) эта погрешность начинает возрастать. Это связано с потерей точности аппроксимации решения полиномами высокого порядка. Установлено, что если ЭДС источника содержит гар-
ческом процессе ток и его производная в конце периода повторения совпадают с соответствующими значениями в начале этого периода:
/(/ц) = i(tN, к ); 1 /'(/1Д) = Г (tN к)| ■
(19)
Объединив в одну систему уравнения (15), записываемые для каждого временного участка, а также граничные условия (18) и (19) и, получим матричное алгебраическое линейное уравнение вида:
М • У = Е, (20)
где Р - вектор-столбец, определяющий ЭДС всех узловых точек на протяжении всего периода; У - вектор решений для всего периода; М1 - подматрица, структура которой определена в табл. 2, но для которой последние две строки изменены:
M =
M1
: ; 0001 -1000
; 0001 ; -1000
M1
0001 -1000
0001 -1000
M1
; ; -1000 0001
; -1000 ; 0001
(21)
Согласно данному алгоритму в системе МаАаЬ разработана компьютерная программа и выполнен пример расчета тока в цепи, описываемое уравнением (1). Результаты расчета показаны на рис. 1. В примере расчета задана ЭДС, содержащая первую (с периодом, равным 0,02 с) и девятую гармоники одинаковой амплитуды:
e(t ) = 10[sin(rot) + sin(9rat)],
а также следующие параметры: R = 1,2 Ом; L = 0,05 Гн ; количество опорных точек на участке: N = 4; количество участков на период: K = 12.
При увеличении количества участков K погрешность вычислений уменьшается, и можно добиться относительной погрешности, имеющей значение меньше 0,01 %■ Повышение до N > 5 степени интерполяционного полинома к существенному уменьшению погрешности вычислений не приводит. Поэтому на практике достаточно принять N = (3-5). Отслеживание процессорного времени расчета с помощью команд tic/toc показало, что предложенный метод оказывается более быстродействующим (примерно на 15 %) по сравнению с методом разложения в ряд Фурье, если число гармоник в функции ЭДС равно двадцати. При этом относительное быстродействие возрастает с увеличением числа гармоник.
Покажем возможности данного метода при моделировании электрических процессов в двигателе постоянного тока, питаемом напряжением сложной формы (рис. 2), полученным на выходе трехфазного тиристорного выпрямителя с мостовой схемой выпрямления (при угле управления а = 25 эл. град.). Рассмотрим один период повторения напряжения меж-
Рис. 1. Результаты расчета установившегося режима в одноконтурной цепи предложенным методом (показано сплошной линией) и точным комплексным методом (показано звездочками)
ду моментами времени t1 и {2. При этом усложним форму напряжения наложением широтно-импульс-ной модуляции (рис. 3).
Расчет выполнен для двигателя типа 2ПБ180М мощностью 12 кВт. В расчете приняты: суммарное сопротивление (сглаживающего реактора и якоря двигателя И = 0,282 Ом); суммарная индуктивность (реактора и якоря двигателя L = 0,012 Гн); противо-ЭДС двигателя Ея = 228 В. Зависимость от времени установившегося тока показана на рис. 4.
и, в
■— s К r-s
\ \ ;
\ : \ S
;............ Л ! .....I-
|............
; .....I
1
; .... : I !
.....;...... .....j- ■ i
; ...........;............ ........... ............ ......... 1...........1.....
:
don? rroiu unrc ппт а 01 | т | П015 rV
Рис. 2. Зависимость от времени напряжения, полученного на выходе тиристорного выпрямителя
и, В
тщ
•ffitr 350 300 250 200 150; 100 50 0
% С kl и3
Рис. 3. Форма напряжения на одном периоде повторения после наложения широтно-импульсной модуляции
Г. А
Рис. 4. Временная зависимость установившегося тока двигателя на одном периоде повторения напряжения
Выводы
1. Разработанный численный метод позволяет выполнять расчеты установившихся периодических процессов в электрической цепи при воздействии ЭДС произвольной формы (с точностью, определяемой принятым числом разбиений периода ЭДС).
2. Время расчета при использовании предложенного метода оказывается меньше на 15% по сравнению с методом разложения в ряд Фурье, если число гармоник в функции ЭДС равно двадцати. При этом
скорость расчета при использовании предложенного метода возрастает с увеличением числа гармоник.
Перечень ссылок
1. Зиновьев Г С. Основы силовой электроники. В 2 ч. Ч. 1. Учебник / Зиновьев Г. С. - Новосибирск : НГТУ, 1999. - 199 с.
2. Теоретичн основи електротехыки : пщручник. У 3 т. Т. 1 : Устален режими лУйних електричних кт iз зосередженими параметрами / Бойко В. С., Бойко В. В., Видолюб Ю. Ф. та н ; за заг ред. I. М. Чижен-ка, В. С. Бойка. - К. : 1ВЦ «Полтехшка», 2004. -272 с.
3. Ортега Дж. Введение в численные методы решения дифференциальных уравнений : [пер. с англ.] / Дж. Ортега, У. Пулл. - М. : Наука, 1986. - 288 с.
4. Самотий В. В. Аналiз усталених ферорезонанс-них режимiв трифазних трансформаторiв з вра-хуванням пстерезизу / В. В. Самотий, Ю. I. Грицюв // Техшчна електродинамка. - 1996. - №5. -С. 59-62.
5. Тиховод С. М. Совершенствование численных методов расчета электромагнитных процессов в сложных нелинейных электрических и магнитных цепях / С. М. Тиховод // Електротехтка та елек-троенергетика. - 2007. - № 1. - С. 56-60.
Поступила в редакцию 09.12.08 г.
После доработки 22. 12.08 г.
Розроблений чисельний метод розрахунку сталих пер'юдичних несинусодних процесв в електричних колах. У основу методики розрахунку покладена апроксимаця похiдноí рiшення диферен-цальних рiвнянь за допомогою полiномiв. Метод перевiрений на конкретних прикладах.
The numerical calculation of steady-state periodical non-sinusoidal processes in electric circuits is developed. Approximation of derivative of differential equations solution is the basis of calculation method by means of polynomials. Present calculation is tested by the help of the concrete instances.
УДК 621.313
И. А. Орловский, Ю. В. Голянчук
Расчет на рекуррентных нейронных сетях математических моделей асинхронного двигателя при питании от автономного инвертора напряжения
Предложена методика расчета и примеры расчета моделей асинхронного двигателя (АД) на полиномиальных рекуррентных нейронных сетях (ПРНС), исходя из экспериментальных данных о режимах его работы при питании от автономного инвертора напряжения с широтно-импуль-сной модуляцией. Методом имитационного моделирования исследована точность предложенной методики и рассчитаны внутренние параметры АД из весовых коэффициентов ПРНС.
Введение
Качественное управление в электроприводе (ЭП) с нелинейными и изменяющимися в процессе работы параметрами является сложной и актуальной задачей [1]. Осуществление желаемых настроек системы управления ЭП может быть достигнуто, если известна математическая модель ЭП и выполняется идентификация его внутренних параметров в процессе © И. А. Орловский, Ю. В. Голянчук 2009 р.
работы. Идентификация параметров сложных объектов может выполняться для отдельных его узлов. Одним из основных узлов ЭП переменного тока является асинхронный двигатель (АД). Как правило, в частотно-регулируемых современных ЭП переменного тока питание АД осуществляется от автономного инвертора напряжения (АИН) с широтно-импульсной модуляцией (ШИМ) выходного напряжения.