CC BY
47
Джанунц Гарик Апетович E-mail: janunts@inbox.ru Тел: 8- 918-506-90-24
Romm Yakov Evseevich
Taganrog State Pedagogical Institute E-mail: romm@List.ru
48, Initsiativnaia, Taganrog, 347936. Phone: 88634 60-18-99
Dzhanunts Garik Apetovich
E-mail: janunts@inbox.ru Phone: 8- 918-506-90-24
УДК 517.91: 518.1
КОМПЬЮТЕРНЫЙ АНАЛИЗ УСТОЙЧИВОСТИ СИСТЕМ ЛИНЕЙНЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ С ПРИЛОЖЕНИЕМ К ОЦЕНКЕУСТОЙЧИВОСТИ СИНХРОННОГО ГЕНЕРАТОРА
Изложен метод компьютерного анализа устойчивости систем линейных ОДУ, определяющий необходимые и достаточные условия устойчивости при ограничениях общего вида на основе разностных схем. Представлены результаты компьютерного моделирования устойчивости систем и численный эксперимент применительно к анализу устойчивости синхронного генератора, работающего на сеть большой мощности. Анализ реализуется на персональном компьютере в режиме реального времени.
Метод; анализа; время.
THE COMPUTER ANALYSIS OF A STABILITY OF SYSTEMS OF LINEAR DIFFERENTIAL EQUATIONS WITH APPLICATION TO AN ESTIMATION STABILITIES OF THE SYNCHRONOUS GENERATOR
The method of the computer analysis of a stability of an ODE linear systems, defining necessary and sufficient conditions for stability is stated at restrictions of a general view on the basis of difference schemes. Outcomes of computer modelling of a stability of systems and numerical experiment with reference to the analysis of a stability of the synchronous generator working on a web of the big potency are presented. The analysis is realised on the personal computer in a condition of real time.
Method; analysis; time.
Рассматривается задача Коши для линейной однородной системы дифференциальных уравнений
. . , . .
Ya.E. Romm, S.G. Bulanov
Y (t0) = Y0
где Y - вектор-функция, координаты которой - искомые функции yj, y2,..., yn от независимой переменной t; A (t) - матрица коэффициентов, n х n . Предполагается, что все функции aij(t) i,j = 1,...,n определены, непрерывны и непрерывно дифференцируемы на отрезке [t0, t ], при любом выборе t = const, t е 110, ^) ; Y(t0) = Y0 - начальный вектор. Всюду ниже используется канониче-
n
ская норма матрицы IIAII = max у \aik \ и согласованная с ней норма вектора
11 11 1< i < n1 1
к=1
|| Y || = max |yk |. Предполагается, что для системы (1) выполнены все условия
1< к < n
существования и единственности решения на [ t0, ^ ). Эти же условия предполагаются выполненными для каждого элемента множества возмущённых решений Y = Y (t), соответствующих возмущённому начальному вектору Y (t0) = Y0, по крайней мере, для некоторого 8 > 0 и любых Y0, таких, что
|~0 - Yc||<8. (2)
Пусть R обозначает область, включающую множество всех решений за-(1) , (2),
R :{ t0 < t <-; Y(t), Y(t): - Y„ || <8 }. (3)
(3) (1)
смысле Ляпунова. Традиционное определение устойчивости заимствовано из [1] с упрощениями [2] на случай рассматриваемого множества R.
В данных условиях для VT е [t0, ^) на отрезке [t0, T ] выполнено:
о | lA (t )| < c , где c = const зависит от T ;
2) если fk(ti, Y ) = ak 1( ti) У1(1) +... + akn(ti) Уп (l), то \fk (t, Y (t ):>| < cJ,
| fk- (t, Y(t)) | < c1, где c1 = const зависит от T и Y0, k = 1,...,n.
(1)
Yi+1 = Yi + hA (ti)Yi, i = 0,1,..., или
Yi +1 = (E + hA (ti)) Yi, (4)
h - .
В дальнейшем при любом выборе t = const, tе[t0, ^) , h и i всегда
предполагаются связанными следующими соотношениями:
t ____ t
t = ti+1, h = , i = 0,1,.... (5)
i + 1
(5)
~+1 = (E + hA (ti)) Yi , (6)
где Y(t0) = Y0, Y0 из (2). При этом всюду ниже предполагается, что Y (t) не
R , , Y(t) (3).
Для возмущения с учетом остаточного члена на каждом шаге получим [2]:
~+1 _ Y+1 = (E + hA (ti))(~ _ Yi) + Qei , (7)
где \ Qe ; ^ cjh2, С = const [2]. Рекуррентное преобразование (7) влечёт
i+1
П <E + hA(tt_f))(~o _ Yo) + Lt
=0
где
i i_k
Li = ZП(E + hA(ti_))QEk_1 + Qeі .
(В)
(9)
£ = 1 Ы 0
Лемма 1. В рассматриваемых условиях имеет место соотношение
1г = О (к). (10)
Доказательство леммы приводится в [2, 3]. Отсюда вытекает
1. 1
Нш Ь, = 0 . (11)
к ^0
Предельный переход в равенстве (8) при любом г из (5) влечет равенство
г
~(г) - у (г) = Иш П ( е + кА (г,_,))(~ - у0) + liш ьг.
к^0 к^-0
? = 0
Стремление к к нулю равносильно стремлению г к бесконечности. Отсюда, с учётом (11) для любого г е [г0, ^) выполнено
~(t)_Y(t) = ИіпП(E + hA(ti_))(~o _Yo).
(12)
f=o
h (12) i : h ( i ) =
i = 0,1,..., поэтому бесконечное
(12) , каждое из которых отличается от предыдущего не только добавлением нового , -лей параметра к в каждом из сомножителей но вого частичного произведения.
Практическая программная реализация условий устойчивости будет выполняться с постоянным шагом на фиксированном промежутке. Это делается по двум причинам: во-первых, экспериментально не обнаруживается разница результатов моделирования при изменении шага в наборе сомножителей, во, .
(12) :
Теорема 1. Чтобы в рассматриваемых условиях решение задачи (1) было
,
limП (E + hA (ti_f))
=0
< c1 = const
(13)
для У г е [ г0, го). Решение асимптотически устойчиво тогда и только тогда, когда
(13) , ,
г
Нш П (Е + кА (г, _/))
(14)
Раздел I. Математические методы синтеза систем
Мультипликативная форма выражений под знаком предела предоставляет возможность запрограммировать их вычисление в виде цикла по числу матричных сомножителей. Это влечет возможность компьютерного анализа устойчивости по значению нормы текущего произведения матриц из левой части.
Если матрица А в (1) не зависит от времени то условия (13), (14), соответственно, примут вид:
1ІШ В
І ^го
І +1
< С1 = сош!
(15)
ДЛЯ
V/ є [і0, го), где В = Е + НА , і +1 = 2к
1ІШ В
І +1
В этом случае условия устойчивости отличаются тем, что не требуют информации о характеристическом многочлене матрицы и о его корнях.
Условия устойчивости, аналогичные (13), (14), строятся и на основе разностных методов более высокого порядка [2]. В частности, для метода Эйлера-Коши имеет место аналог условий (13), (14) [2, 3].
(1) ,
1іш П (Е +1( А ( /і ) + А ( /і _, + Н)(Е + НА ( /і ))))
Ы 0
< 'о2 = сош!
(16)
для V/ е [/0, го). Решение асимптотически устойчиво тогда и только тогда, когда
выполнено (16) и, кроме того, имеет место соотношение ■ \ и
Нш П (Е + - (А () + А (/г _,+ И)(Е + И А (/ _,))))
( = 0
При использовании метода Рунге-Кутта условия устойчивости примут
:
1іш П (Е + ~(Ри _, + 2Р2 і _, + 2Рзі _, + Р4 І _,))
=0
< ~3 = COnst , ДЛЯ V/ є [ і0, го ),
1іш П (Е + - (Р_, + 2Р2І_, + 2РзІ_, + Р4І_, ))
=0
Рц = А(/і ),
Н Н
Р2І = А(/і + -)(Е + —А(/і )),
Н
Н
Н,
Н
Рзі = А(/і + -)(Е + -А(/і + -)(Е + -А(/і )))
2
2
2'
2
И И И И
Р4 г = А( /г + И)(Е + И А( + —)( Е + — А( + -)( Е + - А( )))).
Использование данных условий обеспечивает более высокую достоверность анализа устойчивости при более существенных ограничениях на функции (1) [2, 3].
І ^го
Моделирование устойчивости в условиях теоремы 1 проводится на Delphi 7. При накоплении частичного произведения на каждом шаге цикла вычисляется и через некоторое количество шагов к выводится на печать норма текущего .
program Stability_Ejler;
{$APPTYPE CONSOLE} uses SysUtils;
const h=0.000001; n=4; t0 = 0; TT=100; type matr=array[1..n,1..n] of extended;
var A,B,C:matr; i,j,l:integer; t,s,s0,norma:extended; k:longint; s01:array[1..n] of extended; procedure matrinput (t:extended;var A:matr); begin
a[1,1]:=0; a[1,2]:=1; a[1,3]:=cos(t); a[1,4]:=sin(t);
a[2,1]: = -1; a[2,2]:=0; a[2,3]:=sin(t); a[2,4]: = -
cos(t);
a[3,1]: = -cos(t); a[3,2]: = -sin(t); a[3,3]:=0; a[3,4]:=1;
a[4,1]: = -sin(t); a[4,2]:=cos(t); a[4,3]: = -1 a[4,4]:=0;
end;
procedure metod_E (var A:matr); begin
matrinput (t,A);
for i:=1 to n do for j:=1 to n do
begin a[i,j]:=h*a[i,j]; if i=j then a[i,j]:=a[i,j]+1; end; end;
procedure norma_s; begin
for i:=1 to n do
begin s0:=0; for j:=1 to n do s0:=s0+abs(a[i,j]);
s01[i]:=s0; end;
norma:=s01[1]; for i:=2 to n do
if s01[i]>norma then norma:=s01[i];
end;
begin
t:=t0; k:=0; metod_E (A); repeat t:=t+h; metod_E (B); for i:=1 to n do for j:=1 to n do begin
s:=0; for l:=1 to n do s:=s+a[i,l]*b[l,j]; c[i,j]:=s; end; for i:=1 to n do for j:=1 to n do a[i,j]:=c[i,j]; norma_s; k:=k+1; if k>=1250000 then begin write(norma:20); k:=0; end
until t>=TT; writeln; for i:=1 to n do begin for j:=1 to n do write(' ',a[i,j]:20); writeln; end; readln; end.
Приближение к П (E + hA (tj)) реализуется так, что шаг h не стре-
i = 0
мится к нулю, а берется достаточно малым и фиксированным. Как показал экс, . Пример 1. Дана система уравнений [4]:
f ,
Xj = + x2 + x3 cos t + x4 sin t,
x2 = -xj + x3 sint- x4cost,
x3 =-xj cos t - x2 sin t + + x4,
x4 = -xj sin t + x2 cos t - x3.
Нулевое решение этой системы устойчиво [4], но не асимптотически как слева, так и справа - в силу ограниченности решения на [ t0, ^ ) и на
( , to ].
Результаты работы программы Stability_Ejler следующим образом отражают этот ожидаемый результат:
1.8060 1.6695 1.8499 1.9344 1.6324 1.7364 1.4973 1.7628 1.4162 1.5230
1.6788
1.4056 1.8337 1.4703 1.4526 1.6665 1.7839 1.5841 1.7526 1.6986 1.8427
1.6939
Колебание значений нормы на рассматриваемом промежутке ограничено константой. В соответствии с критерием (13) система устойчиво справа.
С другой стороны, компьютерный анализ устойчивости слева влечет:
1.8060 1.6695 1.8499 1.9344 1.6324 1.7364 1.4973 1.7628 1.4162 1.5230
1.6788
1.4056 1.8337 1.4703 1.4526 1.6665 1.7839 1.5841 1.7526 1.6986 1.8427
1.6939
Значения нормы ограничены, что в соответствии с критерием (13) означает устойчивость решения слева.
В случае постоянства матрицы А процесс нахождения частичного произведения ускорится, если вместо перемножения матриц В возводить в квадрат произведение на текущем шаге. Если в (15) положить /' +1 = 2к, то степень
2к
В находится путём к умножений матрицы на себя. Таким образом,
к = 1о§2 0 +1).
Ниже приводится система линейных ОДУ, моделирующая работу синхронного генератора работающего на сеть большой мощности [2, 3].
Ed:
= 10-
f- 0.5610D 1 0.6793 0.6099 0 0.4948 0.5463 0 - 0.9520
0 -13.7658 1.4409 0 3.6163 1.1781 0 8.5472
0 -15.5076 -150.1554 0 -12.6793 38.9205 0 42.4023
0 - 6.5352 -1.1714 - 2.0723D2 0.9552 2.2156 0 5.4592
0 5.6334 0.4076 0 -16.5675 1.4111 0 - 4.2309
0 -3.8073 52.б270 0 -13.1829 -15б.9117 0 -38.8349
0 2.9781 3.9766 0 -10.6238 - 4.7247 - 4.4063D3 - 5.2010
10000 0 0 -10000 0 0 0 0 10000 0 0 0 0 0 -10000 0
- 0.7494
- 3.3161
- 21.4333
- 2.3385 10.1170 б8.5987 10.711б
0
0
«1
q2 Eq
d2 Ed
«2
x q3 Eq
Ed3
«3
^12
) ^13 13
2
E
12
Результаты анализа устойчивости исследуемой системы для случая В1 = в2 = В3 = 1.0 о.е. на промежутке [0, 109] с шагом |И| = 10-14 имеют вид:
3.8Е+0000 6.6Е+0000 1.2Е+0001 2.3Е+0001 4.4Е+0001 7.7Е+0001
8.5Е+0001
72.4Е-0127 2.0Е-0254 4.3Е-0509 5.8Е-1018 3.8Е-2036 2.1Е-4072
0.0Е+0000
В соответствии с условиями устойчивости результат трактуется как асим-.
Для случая В1 = В 2 = В3 = 0 на промежутке [0, 109 ] с шагом |И| = 10-14:
1.1Е+0000 1.2Е+0000 1.4Е+0000 1.7Е+0000 2.4Е+0000 3.8Е+0000 6.6Е+0000
7.7Е+0001 4.6Е+0001 4.7Е+0001 4.2Е+0000 3.0Е+0000 1.3Е+0000
1.4Е+0000
Ограниченное колебание значений нормы в соответствии с условием (15) соответствует неасимптотической устойчивости.
Анализ устойчивости данной системы выполняется за 0,1 с на Бе1рЫ 7, персональный компьютер Рейшт 4 [2, 3]. Согласно приводимой таблице анализ может осуществляться в реальном времени для матриц большой размерности ( . 1).
1
Время компьютерного анализа устойчивости систем большой размерности
Размерность матрицы 9 х 9 18 х 18 36 х 36 72 х 72 144х144 288х288
- полнения программного анализа 0,1 с 0,3 с 1 с 3 с 7 с 60 с
Необходимо принять во внимание возможность адаптации программы к параллельной вычислительной системе. Ядро программы включает лишь умножение текущей матрицы самой на себя. Умножение матриц обладает естествен,
пары матриц имеет вид Т (п3) = О (1о§2 п). Число умножений в данной про-
2 к
, где т таково, что в результате покрывается весь промежуток [ /0, Т ] приближенного решения системы по методу Эйлера. Иными словами т • И = Т — (0. Отсюда следует, что с точностью до целого
к = 1о§2 ^ Т ^ ^ ^. В итоге, временная сложность рассматриваемого анализа в
максимально параллельной форме без учета обмена составит
Т — t
Т(п3) = О(1оя2п) х 1оя2! -----0
,
.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Чезари Л. Асимптотическое поведение и устойчивость решений обыкновенных дифференциальных уравнений. - М.: Мир, 1964. - 478 с.
2. Ромм Я.Е., Буланов СТ. Метод компьютерного анализа устойчивости систем линейных дифференциальных уравнений / ТГПИ. - Таганрог, 2009. - 119 с. ДЕП в ВИНИТИ 30.04.09. № 268 - В2009.
3. . . -
ния устойчивости систем линейных дифференциальных уравнений на основе матричных мультипликативных преобразований разностных схем / Диссертация на соискание ученой степени кандидата технических наук. - Таганрог: ТРТУ, 2006. - 232 с.
4. . . -
М.: Наука, 1971. - 534 с.
Ромм Яков Евсеевич
Таганрогский государственный педагогический институт E-mail: romm@List.ru
347936, г. Таганрог, ул. Инициативная, д. 48. Тел: 88634 60-18-99
Буланов Сергей Георгиевич
E-mail: romm@List.ru Тел: 8-909-43-69-543
Romm Yakov Evseevich
Taganrog State Pedagogical Institute E-mail: romm@List.ru
48, Initsiativnaia, Taganrog, 347936. Phone: 88634 60-18-99
Bulanov Sergei Georgievich
E-mail: romm@List.ru Phone: 8-909-43-69-543
УДК 681.51
В.В. Соловьев, В.И. Финаев
ПОСТАНОВКА ЗАДАЧИ СИНТЕЗА УПРАВЛЕНИЯ СЛОЖНОЙ СИСТЕМОЙ В УСЛОВИЯХ АПРИОРНОЙ НЕОПРЕДЕЛЕННОСТИ
Выполнена классификация неопределенностей в условиях функционирования сложной системы. Показано применение принципа адаптации для управле-, . Неопределенность; адаптация.
