CC BY
85
ISSN 1992-6502 (P ri nt)_
2014. Т. 18, № 2 (63). С. 237-242
Ъыьмт QjrAQnQj
ISSN 2225-2789 (Online) http://journal.ugatu.ac.ru
УДК 519.62
Идентификация динамических систем
на основе нелинейного матричного преобразования Ли
1 2 А. Н. Иванов , П. М. Кузнецов
1 a.ivanov@fz-juelich.de, 2 peekay.spbu@gmail.com
1 Научно-исследовательский центр Юлих (Р. Германия 1,2 ФГБОУ ВПО «Санкт-Петербургский государственный университет» (СПбГУ)
Поступила в редакцию 16.01.2014
Аннотация. Предлагается метод идентификации нелинейных динамических объектов, основанный на построении матричного отображения. Особое внимание уделяется производительности подхода. Конечный вычислительный алгоритм записывается в виде операций сложения и умножения над числовыми матрицами. Приводится общая классификация подходов решения систем обыкновенных дифференциальных уравнений, на основе которых строятся алгоритмы идентификации. Приведен алгоритм матричного интегрирования систем обыкновенных дифференциальных уравнений. Предлагаемый подход заметно выигрывает по производительности у традиционных методов при реализации его на высокопроизводительных вычислительных системах.
Ключевые слова: параллельные методы решения ОДУ; идентификация динамических систем.
ВВЕДЕНИЕ
Под динамической системой будем понимать объект, изменяющийся по времени, который можно описать системой обыкновенных дифференциальных уравнений (ОДУ)
— = F (t, X), dt
(1)
где X - вектор состояния системы (набор параметров, описывающих объект), Г - векторная функция, описывающая динамику изменения вектора состояния, удовлетворяющая условиям существования и единственности решения задачи Коши. При заданном начальном состоянии системы Х0 = X (¿о), уравнение (1) полностью определяет динамику изменения вектора X. Причем в новый момент времени Х1 состояние системы определяется однозначно.
Традиционным способом решения уравнения (1) является использование пошаговых методов интегрирования [1]. Например, простейший метод Эйлера выглядит как
¿о Хо
¿1 Х^Щ, Хо)((-(о),
ХЫ ХЫ=Г(х, ХЫ-1)(ХЫ-ХЫ-1).
Идея пошаговых методов интегрирования заключается в том, что в каждый новый момент времени вектор состояния системы вычисляется по строго определенной формуле, которая является некоторой аппроксимацией истинного решения. Алгоритм носит итеративный характер: решение задачи Коши происходит последовательным вычислением промежуточных состояний системы. При изменении начального состояния процесс вычисления повторяется заново.
Под идентификацией в этом случае понимается алгоритм поиска такой функции Г, что задаваемая ею система дифференциальных уравнений будет переводить заданное множество начальных состояний в заданное множество конечных состояний через заданный интервал времени (см. рис. 1).
'X1 XN
X02 > < X2 > < XN
XK, Ук
Рис. 1. Пошаговое интегрирование
<
>
ПОСТРОЕНИЕ ОТОБРАЖЕНИЯ РЕШЕНИЯ ДЛЯ СИСТЕМЫ ОДУ
К преимуществам подходов решения систем ОДУ, построенных на основе методов пошагового интегрирования, следует отнести возможность построения высокоточных схем высокого порядка [2]. К недостаткам можно отнести рост глобальной ошибки на длительных интервалах моделирования и низкую вычислительную производительность. Существует иная идеология решения систем ОДУ, основанная на построении отображения. При удовлетворении системой (1) условий теоремы Пикара решение определяется однозначно, и, следовательно, существует некоторое отображение, позволяющее вычислить новое состояние системы в заданный момент времени в зависимости от начального состояния и времени
X=я (I, X0, д.
(2)
Представим, что функция Я нам известна, или, по крайней мере, известна некоторая ее аппроксимация. Рассмотрим теперь задачу перевода множества начальных состояний системы в начальный момент времени во множество конечных состояний за заданный временной интервал (рис. 2). Формулу (2) в этом случае можно записать в общем виде
Хы = я о X0.
(3)
Отметим, что оператор Я переводит теперь любой элемент из множества начальных состояний в элемент из множества конечных состояний за заданный интервал времени. Оператор Я однозначно определяется самой системой. Чтобы вычислить новое состояние через заданный интервал времени, нужно применить оператор к начальному состоянию системы. При использовании пошаговых методов пришлось бы вычислять промежуточные состояния в моменты времени Хь Х2,..., Хм-1. При использовании отображения новое состояние вычисляется за одну итерацию.
X1
X
X?
+ Ы (
■+М (X02)-
+Ы (X? )-
XI
X
> X?
Рис. 2. Применение отображения
Под идентификацией динамической системы в этом случае понимается алгоритм поиска такого отображения Я, которое будет перево-
дить заданное множество начальных состояний в заданное множество конечных состояний через заданный интервал времени (см. рис. 2).
Известной сложностью является поиск самого отображения. Софус Ли, основатель теории непрерывных групп преобразований [3], показал, что в общем случае алгоритм построения отображения по сложности соответствует пошаговому интегрированию. Однако существует ряд подходов, применяемых для оценки этого отображения. К ним следует отнести построение отображения на основе непрерывных преобразований групп Ли [4], использование дифференциальной алгебры [5] для вычисления производных и применение матричного формализма [6].
В данной работе приведен алгоритм численной реализации матричного интегрирования, позволяющего оценивать эволюцию (динамику) отображения Я в заданном интервале времени. Преимуществом матричного подхода является линейность всех операций в конечном виде алгоритма. Нелинейная динамика системы описывается при помощи линейных операций сложения и перемножения числовых матриц.
Понятие идентификации системы, данное во введении, можно соотнести с параметрической идентификацией, используемой в классической теории управления. Понятие идентификации на основе отображения возникает, например, в случаях оценки локального преобразования, которому подвергается некоторый объект, и которое переводит его из одного состояния в другое.
В случае численного матричного интегрирования систем ОДУ под идентификацией следует понимать комбинацию двух методов. При этом целью подбора параметров дифференциального уравнения (1) является соответствие пошаговой динамике всего отображения (3).
ПОСТРОЕНИЕ МАТРИЧНОГО ОТОБРАЖЕНИЯ
Предположим, что систему дифференциальных уравнений (1) и ее общее решение (2) можно разложить в ряд Тейлора
^ = ±Рк (') X[к
Ш к=0
п
X =£ Як (Г)X0к
(4)
(5)
где п - заданный порядок нелинейности; Рк(?), Як(0 - матрицы коэффициентов разложения правых частей уравнения и решения в ряд Тейлора. Под оператором Xм понимается кронеке-ровская степень вектора Х с редуцированием
>
>
к=0
размерности [6]. Уравнение (4) описывает динамическую систему. Матрицы Рк задаются коэффициентами, варьирование которых изменяет динамику системы. Матрицы Як задают отображение (3) в матричном виде и зависят от времени
Пример. Пусть имеется динамическая система, заданная в соответствии с формулой (4)
Ш I х ЛV у
0 - 2
у
'х ^
ху
V У2,
тогда ее решением в виде матричного отображения будет являться соотношение
х'
(
са%{42 ■ г) 8т(л12 ■ г)/42 -41ът(л[2 ■ г) соъ(л[2 ■ г)
х0 у0 2
У0
где (х0; у0) - вектор состояния в начальный момент времени ¿0 = 0.
Здесь матричное отображение приведено в аналитической форме. Коэффициенты второго порядка нелинейности опущены в виду громоздкости. Данное решение может быть проверено непосредственным дифференцированием и сравнением результатов с изначальной системой ОДУ. Для произвольной системы (4) процесс построения такого отображения может быть затруднительным, поэтому будем искать отображения (5) численно, используя идеологию пошагового интегрирования, но примененную к отображению. Дифференцируя уравнение (5) по времени
§ = (г) X 0
лг к=0 лг
и сравнивая полученное соотношение с (4)-(5)
лX
— = р (г) + р (г) X + Р2 (г) X[2] + ... + Рр (г) X[ р] =
Ш
= Р0(г) +
+ Р1(г)(я0 (г) + Я^г)X0 + ... + Як (г)Xo[k])+ (г)Я>(г) + Я1(г^0 + ... + Як (г^ )[2] +
+ Р,
+ ... +
+Р
(г)Я (г) + я1 (г)X0 +... + Як (г)Xo[k])р],
' Р V ) Я0() + Я1( )X 0 + .•• + Як (
можно получить систему уравнений
ля0(г) л
п
= £ Рк (г )(Як (г))
к=0
ШЯ] (г) _ р ,
Рк (у № ] )т'
(6)
где j = 1, ..., п; вектор Х(0 задается равенством (2); под ()т понимается операция транспонирования. Система (6) является системой обыкновенных дифференциальных уравнений относительно коэффициентов отображения Я
— = ~(г, Р). лг
(7)
Заметим, что правые части уравнения (7), во-первых, вычисляются в матричном виде, а во-вторых, не зависят от вектора начального состояния Х0. Уравнение (7) полностью задается системой (4) и может быть разрешено некоторым пошаговым методом интегрирования:
¿0 ^ Я(10),
ь ^ яг),
¿1
N
^ Я( г N).
Алгоритм построения решения матричной динамической системы через заданные интервалы времени запишется в виде следующей последовательности шагов:
1. Задать параметры динамической системы - коэффициенты матриц Рк.
2. Пошагово вычислить эволюцию матричного отображения Я:
a. В начальный момент времени отображение тождественное (переводит начальное множество состояний само в себя).
b. В каждый момент времени отображение Я(^) вычисляется по заданной формуле на основе значений матриц Рк и Я(^г-1).
В случае матричного отображения следует особенно отметить тот факт, что можно применять оператор отображения (3) сразу ко множеству начальных состояний системы
{XI; X2;...; XК != Я о {XX2;...; Xк 1
I"' N ' ' N ^ I 0 ' 0 ' ' 0 Р
причем в случае матричного отображения это соотношение запишется в виде Хы = R ° Х0, где каждый столбец матриц Х1 и Х0 соответствует вектору состояния.
Алгоритм идентификации может быть представлен в виде блок-схемы, изображенной на рис. 3.
х
+
х
+
0
+
Рис. 3. Схема идентификации динамической системы методом матричного интегрирования
В данном случае отслеживается динамика изменения состояний. На основе поступающих данных (новых состояний) корректируется динамическая система (коэффициенты матриц Рк), которая порождает новую последовательность отображений. Уравнение (4) в силу своей нелинейности может описывать нелинейную динамику (например, включать в себя движение системы координат). Производительность метода при этом зависит целиком от порядка нелинейности и размерности задачи.
ВЕРИФИКАЦИЯ АЛГОРИТМА НА МОДЕЛЬНЫХ ЗАДАЧАХ
В качестве модельных выбраны некоторые наиболее хорошо изученные задачи нелинейной динамики.
Модель Лотки-Вольтерра описывает класс уравнений типа «хищник-жертва», соответствующих динамике популяции двух видов. В линейном случае решение представляет собой гармонический осциллятор.
Рис. 4. Модель Лотки-Вольтерра
В нелинейном случае уравнения имеют вид
х = - у - ху, У = х + ху.
На рис. 4 представлен фазовый потрет движения, соответствующий начальной точке хо=0.5, уо=0.5, уже лежащей в области влияния нелинейных эффектов, что видно по искажению эллипса.
Несколько более сложный пример представляет собой система
х = х(4 - х - у), у = У( х - 2)
которая имеет три неподвижные точки {(х, у): (0, 0), (4, 0), (2, 2)}. Наиболее интересно поведение решения в окрестности точки (2, 2). Проводя линейный анализ, можно заключить, что эта точка является устойчивым фокусом, спиральные кривые решения при этом закручены против часовой стрелки. На рис. 5 изображены траектории движения, соответствующие различным начальным точкам.
Осциллятор Ван-дер-Поля является общеизвестным примером системы, обладающей предельным циклом:
х = у,
У = -А,(1 - х2) у + а 2 х.
Линейный анализ устойчивости указывает, что неподвижная точка (0;0) является в случае X2 < 4ю2 неустойчивым фокусом. Однако с ростом х и у нелинейный член Хх2у начинает доминировать, и можно предположить затухание колебаний обратно по направлению к началу координат.
Рис. 5. Устойчивый фокус
Так как в окрестности нуля точки движутся от начала координат к центру, а вдали - в обратном направлении, то непрерывность решений требует существования предельного цикла, представляющего собой замкнутую траекторию. Решения, начинающиеся как внутри нее, так и снаружи, притягиваются ею и никогда ее не пересекают. На рис. 6 изображены результаты численного моделирования системы для значения X = ю = 1.
Рис. 6. Осциллятор Ван-дер-Поля
Для численного анализа и сравнения результатов применялись матричное отображение 4-го порядка нелинейности (см. рис. 4-6) и метод Рунге-Кутты 4-го порядка. Заметим, что при использовании как пошагового, так и матричного интегрирования получаются качественно одинаковые изображения фазовых плоскостей.
ЗАКЛЮЧЕНИЕ
В работе представлен метод построения систем идентификации динамических объектов на основе построения матричного нелинейного
отображения. Рассмотренный подход обладает рядом преимуществ по сравнению с пошаговыми алгоритмами интегрирования, а именно -позволяет строить решение системы ОДУ сразу в виде отображения. Матричный подход успешно применяется при моделировании сложных систем управления пучками частиц [6] на длительных временных интервалах. Дальнейшее развитие метода авторами видится в расширении сферы его применения. Так, например, особый интерес представляет применение подхода в области компьютерного зрения, где массовый параллелизм задач может удачным образом отразиться на матричной записи алгоритма идентификации. Отметим, что теории групп Ли уже используются для алгоритмов оценки движения на изображениях (см., напр., [7]). Однако в настоящее время разработанные подходы не оценивают динамику системы и основываются на способах оценки локального отображения.
СПИСОК ЛИТЕРАТУРЫ
1. Haier E., Lubich C. Geometric numerical integration. Nitherlands: Springer-Verlag, 2006. 644 p. [ E. Haier and C. Lubich, Geometric numerical integration. Nitherlands: Springer-Verlag, 2006. ]
2. Oevel W., Sofroniou M. Symplectic Runge—Kutta schemes II: classification of symmetric methods [Электронный ресурс]. URL: http://citeseerx.ist.psu.edu/viewdoc/ summary?doi=10.1.1.46.5060.pdf (дата обращ. 01.09.2013). [ W. Oevel and M. Sofroniou (2013, Sep. 01). Symplectic Runge — Kutta schemes II: classification of symmetric methods [Online]. Available: http://citeseerx.ist.psu.edu/viewdoc/sum-mary?doi=10.1.1.46.5060. ]
3. Полищук Е. М. Софус Ли. Л.: Наука, 1983. 212 с. [ E. M. Polischyuk, Sophus Lie, (in Russian). Leningrad: Science,
1983. ]
4. Dragt A. J., Douglas D. R. Particle tracking Lie algebraic methods // Computing in Accelerator Design and Operations.
1984. V. 215. P. 122-127. [ A. J. Dragt and D. R. Douglas, "Particle tracking Lie algebraic methods," Computing in Accelerator Design and Operations, vol. 215, pp. 122-127, 1984. ]
5. Hubbard J., Lundell M. A first look at differential algebra [Электронный ресурс]. URL: www.math.washington.edu/ ~bel05/Research/Hubbard_Lundell-A_First_Look_at_Different ial_Algebra.pdf. (дата обращения 01.09.2013). [ J. Hubbard and M. Lundell, (2013, Sep. 01). A first look at differential algebra [Online]. Available: www.math.washington.edu/~bel05/ Research/Hubbard_Lundell-A_First_Look_at_Differential_ Algebra.pdf. ]
6. Андрианов С. Н. Динамическое моделирование систем управления пучками частиц. СПб.: Изд-во СПб. унта, 2004. 376 с. [ S. N. Andrianov, Dynamical modeling of charged particles beam control systems, (in Russian). Saint-Petersburg: Publishing of SPbU, 2004. ]
7. Xu Q., Ma. D. Applications of Lie algebra to computer vision: a brief survey // Int. Conf. Systems and Informatics. 2012. P. 2024-2029. [ Q. Xu and D. Ma, "Applications of Lie algebra to computer vision: a brief survey," in Proc. of Int. Conf. Systems and Informatics, pp. 2024-2029, 2012. ]
ОБ АВТОРАХ
ИВАНОВ Андрей Николаевич, докторант науч. центра, асп. каф. комп. моделир. и многопроц. систем. М-р инф. тех-нол. (СПбГУ, 2012). Готовит дис. о нелинейном матричном интегрировании при моделир. ускорителей заряженных частиц.
КУЗНЕЦОВ Павел Михайлович, асп. каф. комп. моделир. и многопроц. систем. М-р инф. технол. (СПбГУ, 2012). Иссл. в обл. парал. методов для задач компьютерного зрения.
METADATA
Title: Dynamical systems identification based on nonlinear matrix Lie presentation.
Authors: A. Ivanov1, P. Kuznetcov2
Affiliation:
1 Institute for Nuclear Physics, Forschungszentrum Juelich GmbH (FZJ), Germany.
2 Saint-Petersburg State University (SPbU), Russia.
Email: 1 a.ivanov@fz-juelich.de.
Language: Russian.
Source: Vestnik UGATU (scientific journal of Ufa State Aviation Technical University), vol. 18, no. 2 (63), pp. 237-242, 2014. ISSN 2225-2789 (Online), ISSN 1992-6502 (Print).
Abstract: A map method for identification of nonlinear dynamical systems is proposed. The approach is based on matrix representation of map. This allows to implement all operation via multiplication and addition of numerical matrices. In the article a set of the numerical methods for ordinary differential equations is compared with matrix integration. In compare with step-by-step integration schemes the given approach has great advantage in performance and can be easy implemented in high-performance computational systems.
Key words: Parallel methods in ODE solving; dynamicak system identification.
About authors:
IVANOV, Andrei Nicolaevich, Postgrad. (PhD) Student, Dept. of Computer Modelling and multiprocessing systems. PhD student, Institute for Nuclear Physics, Forschungszentrum Juelich GmbH, Germany. Master of Information Technology (SPbU, 2012). KUZNETCOV, Pavel Michailovich, Postgrad. (PhD) Student, Dept. of Computer Modelling and multiprocessing systems. Master of Information Technology (SPbU, 2012).
