Многопараметрическая идентификация конструктивных параметров методом
объединенного принципа максимума А.А. Костоглотов, А.И. Костоглотов, С.В. Лазаренко, Д.С. Андрашитов
Введение. Синтез оптимального управления механическими системами в настоящее время базируется на методах классического вариационного исчисления, динамического программирования Р. Беллмана и принципа максимума Л.С. Понтрягина [4,5]. Однако эти методы дают условия оптимальности для рассматриваемого данного момента времени, а проблема синтеза решалась в ряде работ с использованием гипотез об управлениях, линиях переключения, функциях Беллмана и др. [4,6]. Для решения проблемы синтеза управления применялись также подходы, основанные на методах декомпозиции и использовании заданной программной траектории [7], на использовании управлений по обратной связи [8] и введении функции А.М. Ляпунова.
Если принцип максимума применить к признаку истинного движения Гамильтона -Остроградского [3], то условия оптимальности получаются для конечного промежутка времени, то есть рассматривается задача синтеза оптимального управления. Решение такой экстремальной задачи получается в форме объединенного принципа максимума из условия максимума функции обобщенной мощности [9,10]. Здесь синтез оптимального управления строится на основе анализа структуры фазового пространства. При этом устанавливается аналитическая зависимость управлений от фазовых координат и связь с целевым функционалом.
1. Теорема объединенного принципа максимума [9,10].
Пусть задан целевой функционал
tk
J = | F(q, q)dt ^ min . (1)
*0
Строится расширенный функционал [2,9,10]
tk
J„, =ЦЛ(Г + A) + F\dt, (2)
*0
— n
где X - неопределенный множитель Лагранжа, T = — ^ askqsqk - кинетическая энергия,
2 s,k=—
üsk - коэффициент инерции, q = \ql,...,q], q = \ql,...,qn] - вектора обобщенных координат и
n tk
скоростей, A = ^JQsdqs - работа обобщенных сил, Q = \Q,...,Qn] - вектор обобщенных
s=1 ,
*0
сил, зависящий от вектора идентифицируемых параметров z = \zy,..., zn ].
Теорема. Для того чтобы обобщенная сила Q(q, q, z, t) ^ GQ и соответствующая ей
траектория (q, q) e R2n доставляли минимум расширенному функционалу (2), необходимо выполнить условия максимума для обобщенной мощности
n
ф(^ q, Q X) = max £ \XQs(^ q) + Vs q, (3)
q^gq s=—
где X = const > 0, а на концах траектории t = t0, t = tk выполняются условия трансверсальности
X( A - T) + F = 0, (4)
5F
здесь V =-------фиктивная сила, зависящая от формы целевого функционала.
Sqs
Доказательство. Пусть к расширенному функционалу (2) последовательно применено асинхронное и игольчатое варьирование [3,4]. Для произвольной обобщенной силы Q ^ ^ асинхронная вариация функционала будет иметь выражение
д/„ =[мг + a) + f ]-д(|;; + zj
І ST S], + ^~ Sq, + Q,Sj,1 + V,Sj,
s dq , J
(5)
где ¿¿¡,, - синхронные вариации обобщенных координат и скоростей;
Мех( = ¿1 ^ + ¿ехіАї - асинхронная вариация функционала.
Интегрирование по частям первого слагаемого под знаком интеграла и замена в граничных условиях синхронной вариации так, чтобы асинхронная вариация равнялась нулю Ац =&],; + цАї = 0, откуда = —цАї, преобразует выражение (5) к виду
І
Zii^Sq.dt SqJtk-Zj
Sq, Sq, lt0
ST
,=1 t ti
=l t
'dA\ dT _ ,ST _
— I-----Sq, + A-----------Sq,
v dt J dq , dq ,
dt =
n k
= -A2T•Дt|tk - Zj
'dX\ ST _ ,ST _
— I----Sq, + A-----------Sq,
v dt J Sq , Sq ,
(6)
dt.
Из выражений (5) и (6) следует A = const [2], а при преобразовании краевых условий применена теорема Эйлера об однородных функциях [3]. С учетом преобразований первая асинхронная вариация (5) приводит к условиям трансверсальности (4) и выражению
n k
J=Z j
s = 1
A -
d ST ST
dt Sq Sq,
+ ^- + Q,
+V
5qsdt > 0.
(7)
Пусть из допустимой области Се выбрана другая обобщенная сила, но полученная из первой игольчатым варьированием ^ = Q + <Х), ф 0 при / е [г, т + 8т] [4].
Асинхронная вариация функционала для этой обобщенной силы запишется аналогично (7)
n lk
J»*=Zj
,=i t,.
A
d ST ST
\
■ + ■
dt Sq , Sq ,
+ Qs
+ V„
Sq, Л.
(8)
В силу произвольности синхронные вариации можно получить одинаковыми ^ при ? = т + 8 .
Из сравнения (7) и (8) получается вторая асинхронно-игольчатая вариация функционала
А2 ^ = /„ — / =
n lk
Z jjA
d S(T - T) + SJ^ + Q-Qs )
+ (V,- V, )\Sqsdt.
(9)
ж дд. дд.
Отрезок [¿0, ^ ] можно разделить на три части. На полуоткрытом интервале I [¿0 ,т) произвольная и варьированная обобщенные силы совпадают, поэтому А2= 0. На ограниченном замкнутом интервале II [г, т + М & ф Q, но в силу малости интервала 8 = О(<), Т< — Т = О(<) [1]. Вторая вариация функционала определяется соотношением
п т+,81 п
Д2= X |— Q.)+К. — К=ТШ.— Q.)+К — К)81.з. (10)
5=1 т .=1
На полуоткрытом интервале (т + Зt, гк ] <2, = Q и выражение под знаком интеграла (9) будет равно нулю
0
=l
0
d_ 6(Te-T) d(Te-T) =0 dt dqs dqs '
Это уравнение Лагранжа второго рода для возмущенного движения с начальными условиями
t = r + s5t, Sq(T + sSt) = -qsSt, Sq(T + sSt) = -qsSt.
Вторая вариация функционала будет
n ^
¿J*,,,, =Z jV - v Sqd
s=1 t+eSt
или
d ,? =t V - V У (t) • t = T + sSt, HJ^T + St) = J,. (11)
dt s=1
При предельном переходе s ^ 0, Vs- Vs ^ 0, qs^ q, qs ^ q, Sq ^ 0, Sq ^ 0, St - произвольно. Из (10) с учетом (1) получается
Д2 J n
Hm =-£[Л(0, -Q.) + V - Vqa2 > 0. (12)
s^0 s
s s=1
Если обобщенная сила Q доставляет минимум функционалу (2), то из (12) вытекает теорема объединенного принципа максимума (3)
n n
ф(ч, q, Q,x) = Z (Щ + V.)q. = + Vs. q. (13)
s=1 Q^GQ s=1
Это условие не нарушается вдоль траектории, так как в соответствии с (11) и условиями предельного перехода &JextIII = const вдоль траектории t е(т + sSt, t*.].
Из теоремы объединенного принципа максимума легко выяснить, что множество, на котором функция Ф(ч, q, Q, Л) достигает максимума, определяется совпадением знаков сомножителей sign(ЛQs + V ) = signqs или их пропорциональностью (Л(25 + V ) = ßs (q, q)q, где ßs (q, q) - синтезирующая знакопостоянная функция. Теорема позволяет с точностью до функции ßs (q, q) определить искомые обобщенные силы
Qs =Л-1\Мз(^q)qs -Vs]s = 1,n . (14)
Обратная подстановка обобщенной силы (14) в условие (13) устанавливает знакоотрицательность синтезирующей функции
n
^ q ал)=^ßsq2s <0.
s=1
2. Построение синтезирующей функции. Согласно (14) равенства
Qs =Л-1\Мз (^ q)qs - Vs ] = 0, s =1, n (15)
определяют в фазовом пространстве гиперповерхности переключения управления. Так как на этой гиперповерхности обобщенная сила равна нулю, условия трансверсальности (4)
преобразуют в условие постоянства обобщенного кинетического потенциала в данный
момент времени
L(q, q, t) = ЛТ(q, q) - F(q, q) = l = const. (16)
Это уравнение представляет собой поверхность гиперболического параболоида в фазовом
пространстве переменных Лагранжа qs, qs (s = 1, n ) .
Преобразование Лежандра [9] функции Ь(д, д, г) по переменным д (5 = 1, п ) есть функция Гамильтона, представляющая поверхность эллипсоида в переменных Гамильтона
У*, Р5(5 = 1 п )
А
Н(д,р,г) = ЛТ + Р = -Л ^ АкР*Рк + Р = к = овтг
(17)
в которой величины д выражены через д, р, г при помощи уравнений
дЬ -
Р5 =— , 5 = 1, П
дд*
для обобщенных импульсов, при этом при проведении преобразования величины д, г играют роль параметров. Здесь А1к - алгебраическое дополнение элемента а&к гессиана кинетического потенциала
д2Ь
дд5ддк
= Ф 0 .
5,к=1
Качественный анализ поверхности гиперболического параболоида (16) устанавливает, что ее главные сечения суть параболы Р = —21,ЛТ = 21, направленные в разные стороны. Сечения гиперболического параболоида плоскостью I = 0 суть семейство прямых ЛТ — Р = 0, которых на поверхности гиперболического параболоида бесчисленное множество. Плоскости I Ф 0 пересекают гиперболический параболоид по гиперболам.
Качественный анализ поверхности эллипсоида (17) устанавливает, что ее сечения являются эллипсами.
Из анализа так же следует, что если центр симметрии эллипса является терминальной точкой в задаче управления, то траектории, ведущие в нее, должны быть сопряжены к линии эллипса и являются прямыми, гиперболами, эллипсами. Они также называются линиями переключения.
Построение синтезирующей функции л в связи с этим проводится в два этапа. В начале выражение (15) подставляется в уравнение Лагранжа
Л
ё дТ дТ дА
- дд5 дд,5 дд_5
= &5д5 — К, 5 = - п
что позволяет получить ¡л,, - угловой коэффициент касательной к эллипсу.
Исключение с помощью условий трансверсальности (4) фиктивной силы
дР / дТ дА ^
V =----= Л
^5
^5 дд.
5 У
преобразует уравнение Лагранжа на поверхности переключения к виду
Л
-г
дА
дд.
= ля*.
(18)
5 У
Но так как на поверхности переключения -------= 0, а по (18) ЛЛд = V то из (18) получаются
дЧ5
два соотношения
Л5 = —Л
ёР5
и Л-^ = — V, 5 = 1, п
-г
Первое устанавливает, что функция ¡Лц является модулем углового коэффициента касательной к фазовым траекториям д (г) на поверхности переключения Н(д, Р,г) или
п
Ь(д, р, г) с коэффициентом деформации Л ; второе является уравнением этих траекторий на этой же поверхности.
Равенства (16) и (17) показывают, что поверхности переключения при фиксированном времени являются изоэнергетическими. Откуда по уравнениям Уиттекера определяется угловой коэффициент Д на этой поверхности
-Р5
= —Л
дН
дЧ5
дН
дР 5
V
п л
У АкРк
И
Согласно С - свойству измерительных функций Н.Н. Лузина, синтезирующая
функция Д (д, д) измерима, так как может быть сделана непрерывной Д =
_ — VI
на
а
множестве сколь угодно малой меры а .
В зависимости от начальных условий изоэнергетические поверхности (16), (17) образуют семейство, вырождающееся в точку (вершина или фокус) при I ^ 0, к ^ 0. Поэтому всякая траектория (линия переключения), сопряженная указанным изоэнергетическим поверхностям, должна иметь свое направление в сторону фокуса. Но это возможно, если модуль углового коэффициента касательной к линии переключения Д (д, д) и модуль углового коэффициента касательной к изоэнергетической поверхности
Д (д, д) находятся в соотношении
Д -д =
1
, м > 0 .
(19)
Вследствие этого для семейства линий переключения можно записать синтезирующую функцию общей для всех перечисленных линий формулой
Л5 (^ д) =—Л
-Р5
-д,
\д5
му,\+
, 5 = 1, п
а
А закон изменения обобщенной силы на фазовой траектории истинного движения
получит вид
а=л—
— \д5\д5
ММ+'
5 = 1, п
(20)
где Л,М8 ,а - константы, определяемые из решения краевой задачи управления.
Согласно (19) соответствующее семейство касательных к линиям переключения
^=—, 5 = т,п,
ёд. му 5 му,У ,
где величины Рх, дв ,V и производные вычисляются вдоль линии переключения.
Справедливость разработанной теории подтверждается так же совпадением предлагаемых решений с решениями задач известными методами в частных случаях. В простейших случаях уравнение (20) интегрируется. Пусть р = д, V = д , М,, = М . Тогда уравнение линии переключения
д + С = ±кдм,
где С , к - константы, определяемые из краевых условий.
При М = 2, С = 0 - это уравнение параболы, проходящей через начало координат, соответствует задаче о переводе фазовой точки в начало координат [11]. Возможны два решения для уравнений. Пусть а = и . Тогда:
.V
1. Если траектория точки совпадает с линией переключения, то
|q|q
u = X 1
£lM\q\ + ss q_
u G G„
что методом принципа максимума Л.С. Понтрягина получить нельзя;
2. Если траектория точки пересекает линию переключения, то
к\д\м 1 д
u = X 1
u G G„
что в точности совпадает с решениями Л.С. Понтрягина и А.Т. Фуллера [4, 11, 12].
3. Задача идентификации. Рассматривается динамическая система, движение которой подчиняется принципу Гамильтона- Остроградского
*к
S'R = J (ST -S'A)dt,
*0
n
S'A = Z QsSqs - элементарное приращение работы.
s=1
Уравнение наблюдения имеет вид
У = H q *) .
Требуется определить такие постоянные параметры z. e[z,...,zn], j = 1,m , что бы
достигался минимум целевого функционала невязки
*к 1 *к
Ji = J F (q, y, z)dt = - J (y-H )T • (y - H )dt ^ min . (21)
*0 *0
4. Метод решения. Идентификация параметров z производится в два этапа. На первом этапе в соответствии с заданным функционалом (21) по основной теореме объединенного принципа максимума [2, 9, 10]
n
Ф(^ q N, X) = max Z [xNs +(ys-Hs)],
q^gq s=1
вычисляется значение обобщенной силы N = [N,...,N], реализующей наблюдаемое движение y = [y-,..., Уп ]
qs|qs|
N =х
+ (ys-Hs )
s = 1, n
+ E,
На втором этапе идентификация параметров проводится на основе оптимального закона движения qs = -H~*(y) путем сопоставления обобщенной силы Ns с истинной силой Q = Q (q, q, t, z), явная зависимость которой от конструктивных параметров известна. Для этого конструктивные параметры выбирают такими, чтобы среднеквадратическая ошибка I за характерный (произвольный) отрезок времени ^ -10 между обобщенной силой и заданной обобщенной силой была бы минимальной
к 2 1 = ¡\Qs(q,q,Zt) -Ns(q,q,Zt)] dt ^ min • (22)
-s '
*0
Так как обобщенные координаты и наблюдения теперь известны как функции времени: д = д(1), д = д(1), у = у{1), то вместо минимизации функционала I
отыскивается минимум функции I(zj). При небольшом количестве неизвестных zj и их линейном вхождении в обобщенную силу а можно решать систему
& п •
----= 0, j = 1, т ,
дzi
а в более сложных случаях применить программу поиска экстремума.
Пример 1. Рассматривается идентификация параметров жесткости с = 5,064 и сопротивления Ь = 1 динамической системы, математическая модель которой при п = 1 имеет вид
д + Ьд + сд = 4sin(3t),
г0 = 0 у(^0) = 1, КО = °.
При уравнении наблюдения у(г) = д(г) целевой функционал (21) записывается в следующей форме
1 гк
А = -1 (д-у)2 ж.
Уравнение оптимальной траектории, реализующей наблюдаемое движение, синтезированное на основе объединенного принципа максимума на первом этапе идентификации [7] записывается так
допт| дот.
= л
М|Чопш- у| + ‘
' - (допт - у)
г 0, допт0 у0, Чопт^ 0; ^ 4, доптk доптk 0,
где Л-1 = 1200, М = 1.3, а = 80.
Оценка конструктивных параметров определялась по уравнениям (23), которые для рассматриваемого случая в развернутой форме имеют вид
д1 иле
— = аис + аиЬ - А = 0, дс
д1 и А л
— = а2!с + а22Ь - А2 = 0, дЬ
11=| Ч 2&, а12 = а21 =| д • д&, а22 = | д 2ёг
где а11 = J д
г0
гк
А =-]Л1
д|д|
м|д - у| +
£
■-(д - у)
дЖг, А2 — —^ Л
д|д|
м|д - у| +
£
■-(д - у)
джг.
Результаты моделирования наблюдаемой траектории у , оптимальной траектории д , реализующей наблюдаемое движение, а также траектории объекта с идентифицируемыми параметрами д показаны на рис. 1.
у
Чопт
Р I лс* Л
0
0
0
Ошибки оценки показаны на рис. 2, где обозначено: Зопт = у - допт, 5д = допт - дгЛ,
5 = ^опт - Чш.
0.04
0.024
5д0
0.008
5д
¿Ч -------- - 0.008
- 0.024
. ♦ ,
1 ■, , , _■ ч ”|
0 : 5 ' '
•
- 0.04
t.
г
Рис. 2
Ошибки оценки траектории
В результате получены следующие оценки конструктивных параметров: с = 5.065, Ь = 1.007. Относительная погрешность идентификации параметров составляет соответственно 5с = 0,1% и 5Ь = 0,7%.
Пример 2. Рассматривается идентификация параметров динамической системы
д + Ьд + сд + ае~ыд3 = 4 бш t,
где Ь = 3 , с = 2, а = 2, к = 0,1.
Уравнение наблюдения имеет вид у{1) = д{1) .
Для поиска экстремума функционала применим процедуру половинного деления. Для этого интеграл (39) представим в виде
1к
I ==1
- Ьд - Сд - ае д + 4бш t -Л
д|д|
М\д - у + ,
где Л-1 = 1000, М = 1000, £ = 800.
Полученные оценки параметров приведены в Табл.1.
--(д - у)
л2
Таблица 1 - Оценки конструктивных параметров
Ь с а к
Истинные параметры 3 2 2 0.1
Оценки 3.002 1.996 2 0.1
Результаты сравнения фазовых траекторий действительного д^) и оцениваемого на первом этапе алгоритма движения допт^) представлены на рис.3. Результаты сравнения наблюдаемой траектории у^) и траектории с идентифицируемыми параметрами д^а ^) представлены фазовыми портретами на рис. 4.
У
допт'
2
1
- 2 - 1 0 1
- 1
- 2
- 3 - 4
у _
дтс1-
2
1
2 - 1 0 1
- 1
- 2
- 3 - 4
УЛопт
Рис. 3
Фазовые портреты оптимальной оценки траектории по наблюдениям и наблюдаемой траектории
У^ '
Рис. 4
Фазовые портреты наблюдаемой траектории и траектории с идентифицированными параметрами
2
2
Относительная погрешность оценки параметров составляет соответственно 5Ь = 0.2%, 5с = 0.4%, 5а = 5к = 0%.
Заключение. Новый метод идентификации конструктивных параметров объектов на основе предложенного метода объединенного принципа максимума обладает универсальностью, а синтезируемые на его основе алгоритмы отличаются минимумом вычислительных затрат и простотой. Его применение обеспечивает высокую точность расчетов, что подтверждается результатами численного моделирования: относительная погрешность идентификации параметров не превышает 0,7%. Построение метода не требует введения дополнительных гипотез в отличии от получивших в настоящее время распространение. При этом решить задачу идентификации можно даже в случае нелинейности входящих параметров в уравнение движения и функционал.
Литература
1. Сейдж Э.П., Мелса Дж.Л. Идентификация систем управления. - М.: Наука, 1974, 248 с.
2. Костоглотов А.А., Костоглотов А.И. Лазаренко С.В., Шевцова Л.А. Синтез оптимального управления на основе объединенного принципа максимума.// Известия ВУЗ. Сев. Кав. Регион, Технические науки, 2010, №2 (154), С. 31-38.
3. Лурье А. И. Аналитическая механика. М.: Государственное издательство физико -математической литературы, 1961. - 453 с.
4. Понтрягин Л.С., Болтянский В.Г., Гамкрелидзе Р.В., Мищенко Е.Ф. Математическая теория оптимальных процессов - М.: Наука, 1971, 384 с.
5. Беллман Р. Динамическое программирование.- М.: ИЛ, 1960, 400 с.
6. Наумов Г.В. Построение кривой переключения для задач оптимального управления с учащающимися переключениями. // Изв. РАН. ТиСУ. 2003. №3. С. 46-51.
7. Пятницкий Е.С. Принцип декомпозиции в управлении механическими системами. // Докл. АН СССР. 1988.т.300, №2. С. 300-303.
8. Ананьевский И.М. Непрерывное управление по обратной связи возмущенными механическими системами. // ПММ. 2003. т.67.вып. С 163-178.
9. Костоглотов А.А., Костоглотов А.И., Лазаренко С.В. Объединенный принцип максимума в задаче синтеза оптимального управления нелинейными системами. // Автоматика и вычислительная техника. 2007. № 5, С.52-61.
10. Костоглотов А.А., Костоглотов А.И., Лазаренко С.В. Объединенный принцип максимума в задачах оценки параметров движения маневрирующего летательного аппарата. // Радиотехника и электроника. 2009. т.54. № 4, С. 450-457.
11. Fuller A.T. Study of an optimal non-linear control system. // Journal of Electronics and Control.1963. №1(15). pp. 63-71.
12. Kostoglotov A.A. Solution of Fuller's problem on the basis of the joint Pontryagin -Hamilton - Ostrogradskii principle. // Automatic Control and Computer Sciences, 2007, Vol. 41, No. 4 pp. 179 - 187.