Построение оптимальных регуляторов по нелинейным критериям качества на примере одной динамической системы Текст научной статьи по специальности «Математика»

Математическое моделирование. Оптимальное управление Вестник Нижегородского университета им. Н.И. Лобачевского, 2012, № 2 (1), с. 165-176

УДК 519.71+519.853.4

ПОСТРОЕНИЕ ОПТИМАЛЬНЫХ РЕГУЛЯТОРОВ ПО НЕЛИНЕЙНЫМ КРИТЕРИЯМ КАЧЕСТВА НА ПРИМЕРЕ ОДНОЙ ДИНАМИЧЕСКОЙ СИСТЕМЫ

© 2012 г. С.Ю. Городецкийх, А.С. Сорокин 2

1 Нижегородский госуниверситет им. Н.И. Лобачевского 2 ООО «МФИ Софт», Н. Новгород

gorosyu@gmail. com

Поступала в редакцаю 18.01.2012

Исследуется применение численных методов оптимизации для выбора параметров регулятора по нелинейным критериям качества с учетом нелинейных уравнений объекта. Предлагается двухуровневая процедура, на внутреннем уровне которой происходит построение модального или LQR-регулятора по линеаризованной модели при параметрах квадратичного или модального критерия, выбираемых численным методом, который на внешнем уровне решает задачу оптимизации по дополнительным нелинейным критериям. Работа выполнена на примере управления конкретным нелинейным объектом.

Ключевые слова: численная оптимизация, нелинейные критерии качества, LQR-регулятор, модальный регулятор, численное моделирование.

Введение

Классические подходы линейной теории автоматического регулирования, развитые в первой половине прошлого века, позволяли решать вопросы устойчивости вплоть до получения аналитических описаний областей устойчивости и обобщенной устойчивости [1, 2]. Это давало возможность конструктору выбирать конкретные параметры регулятора в пределах построенных областей, привлекая (в основном, на интуитивном уровне) дополнительную неформализованную информацию об объекте управления и дополнительных целях управления. Развитые за последние двадцать-тридцать лет новые подходы к синтезу оптимальных, модальных и робастных регуляторов [3-7] позволяют в значительной степени автоматизировать процесс их построения с обеспечением устойчивости как одного из минимальных требований. Однако оптимальность определяется только по отношению к критериям качества ограниченных видов и не позволяет непосредственно учитывать произвольные дополнительные требования. В ряде случаев используемые критерии качества сами обладают параметрами, выбор которых обычно не регламентируется и не всегда является очевидным, хотя существенно влияет на качество управления.

В то же время применение численных методов нелинейного математического программирования [8-11] предоставляет принципиальную возможность построения регуляторов, опти-

мальных по достаточно произвольным критериям качества, вычисляемым на траекториях нелинейной модели объекта. Однако использование этого подхода при оптимизации непосредственно по параметрам регулятора (в рамках выбранной его структуры) может привести к ряду проблем. Одна из них связана с невозможностью вычисления некоторых критериев качества при выходе параметров регулятора из области устойчивости. В то же время предварительное получение ограничений, описывающих область устойчивости, рассматривается как нежелательное действие, поскольку усложняет процесс подготовки к численному выбору параметров.

Компромиссом может являться сочетание построения в рамках линейной теории оптимальных регуляторов, автоматически обеспечивающих асимптотическую устойчивость, с нелинейной оптимизацией (по свободным параметрам использованного метода построения) дополнительных критериев качества, которые можно численно определять по исходной нелинейной математической модели объекта. Возможности такого подхода неочевидны, поскольку численно определяемые критерии могут обладать «плохими» свойствами с точки зрения гладкости или липшицевости, что ухудшит характеристики решаемых оптимизационных задач. В данной работе указанный подход исследуется экспериментально с использованием численных методов локальной оптимизации.

Рис. 1. Вид модели лабораторной установки «Портальный кран»

1. Изучаемый объект и неформализованная цель управления

Исследование выполнено на примере конкретной нелинейной динамической системы, в качестве которой использована математическая модель лабораторной установки с компьютерным управлением «Портальный кран», размещенной в учебно-исследовательской лаборатории кафедры теории управления и динамики машин факультета ВМК ННГУ.

В принимаемой идеализации установка представляет собой каретку (тележку) массой Мк с легкими (невесомыми) зубчатыми колесами радиуса г, движущуюся по горизонтально расположенным зубчатым рельсам, что исключает проскальзывание. На каретке установлен электромотор постоянного тока, вращающий колеса через редуктор с передаточным числом #р. Масса мотора Мм, момент инерции ротора J, омическое сопротивление активной катушки ротора R, а ее индуктивность L; питающее напряжение v(t). Далее при составлении математической модели будет учитываться совокупная эквивалентная масса каретки с электромотором и колесами, обозначаемая через М, а также совокупный коэффициент трения ^ при движении каретки вдоль рельсов.

К центру масс каретки на шарнире свободно подвешен невесомый абсолютно жесткий стержень длины I с прикрепленной к его концу точечной массой т, имитирующей перемещаемый «краном» груз. Угол отклонения стержня от вертикали - ф. В узле крепления стержня к каретке присутствует вязкое трение с коэффициентом hф. Общая схема установки показана на рис. 1.

(1)

Вначале каретка и груз покоятся в положении, соответствующем координате центра масс каретки x = 0 (ось x направлена вдоль рельсов). За счет регулирования величины напряжения v(t), подаваемого на электромотор, требуется обеспечить перемещение и стабилизацию груза в положении, соответствующем x = x (точка С).

Установка снабжена датчиками, измеряющими значения x(t), ф(у) и их первые производные. Используется линейный (в каждом из диапазонов ±п + 2nk отклонений по углу ф) регулятор по состоянию, формирующий входное напряжение

v(t) = kj ((ф + л) mod 2л - л) + ^ф +

+ k3( x - x*) + k4 x.

Значения физических параметров установки приведены в табл. 1.

Примем, что неформализованной целью является выбор параметров регулятора, обеспечивающий перемещение x-координаты груза в заданную окрестность требуемого значения x = = x и стабилизацию. При этом перемещение нужно выполнить за возможно меньшее время при амплитудном ограничении колебательности в зависимости A(t), где A(t) = x - x(t) -—&т(ф(ф. Должно учитываться ограничение на предельно возможное значение напряжения от источника питания |v(t)| < V+.

2. Математическая модель объекта регулирования

Система нелинейных уравнений движения разомкнутой системы рассматривается без учета индуктивности L активной катушки ротора, поскольку данная величина, являясь пренебрежимо малой, при своем включении в модель повысила бы дифференциальный порядок системы уравнений и привела к необходимости учета быстрых движений, что влечет увеличение погрешности при численном интегрировании системы уравнений. С учетом сказанного принятая модель объекта имеет вид:

mix cos ф + ml2 ф + mgl sin ф = -^ф,

(M + m)x + тіфcosф-тіф2 sinф = f (t)-hxx, f (t) = y((v(t) -(s/p)x)/ R),

(2)

Значения параметров модели лабораторной установки «Портальный кран»

Тaблaцa 1

Параметр м m l R g hф hx в Y є

Значение 1.073 0.019 0.641 2.6 9.81 0.0024 5.4 0.024 4.481 0.00767

Размерность, СИ кг кг м Ом кг-м/с2 Н-м-(рад/с)-1 Н-(м/с)-1 м/рад Н/А В-(рад/с)-1

где Д?) - эквивалентная сила линейной тяги, развиваемая электромотором с учетом редуктора и колес, Р = г#;1 - коэффициент пропорцио-

нальности в связи

х = Р®м (юМ - угловая

ско-

рость вращения ротора), у - коэффициент в зависимости X?) = у1(?) (I - ток в активной катушке ротора), є - коэффициент самоиндукции в активной катушке ротора, М = Мк + Мм +

г ) - совокупная эквивалентная масса

каретки и двигателя с учетом его эквивалентной инерционной массы.

Следует обратить внимание, что главный определитель Д системы (2), возникающий при ее разрешении относительно старших производных х, ф и имеющий вид Д = т12(т(со^'ф -

- 1) - М), всегда отрицателен и, следовательно, никогда не обращается в ноль. В конечном счете, отсюда (с учетом вида остальных элементов уравнений) вытекает выполнение достаточных условий дифференцируемости решений замкнутой системы (1), (2) по параметрам регулятора к7

(7 =1,4) на рассматриваемом далее (в разделах 5-7) конечном интервале времени и промежутке изменения ф е (-л,+л). Это свойство учитывается далее в разделе 6.

Линеаризованная в окрестности ф = 0 система (2), записанная относительно вектора фазовых переменных х = (хь х2, х3, х4)Т = (ф,ф,х-х*,х)Т, имеет вид

х = А • х + Ь • у(£),

где

' 0 1 0 0 > Г 0 ^

А = а 2і а22 0 а24 , ъ = Ъ2

0 0 0 1 0

5! а 42 0 а44 VЪ4 ,

(3)

(4)

а, = —

(М + т) g

(М + т)й

М£

Мт£

= _1_ 4 = М

а41 =—, а42 =-------------------------•

41 М 42 МТ

ув

--------+ hx

яр х

(5)

mg

Ъ2 = —■

у

Ъ 4 =

ув

--------+ hx

яр х

МК£’ 4 МК' (6)

Система (3) замыкается обратной связью у(?) = кТ х = к1 х1 + к2 х2 + к3 х3 + к4 х4. (7)

Значения параметров физической модели «Портальный кран» представлены в табл. 1.

3. Обсуждение некоторых подходов

линейной теории синтеза регуляторов

Современная теория управления предоставляет различные методы синтеза регуляторов. Остановимся кратко на тех, которые будут использованы в дальнейшем изложении.

LQR-регуляmор. При этом подходе коэффициенты к = (к\, к2, к3, к4)Т линейного регулятора (7) для объекта (3), (4) определяются из условия достижения минимума квадратичного функционала вида (8), вычисляемого в силу системы

да

3[у] = |(х^)Т Q х(?) + гу(?)2) Ж, (8)

0

где скалярный параметр г > 0, а матрица Q неотрицательно определена. Известно [3,4,6], что оптимальный с точки зрения функционала 3[у] закон регулирования в классе линейных регуляторов (7) характеризуется вектором параметров вида

к = —г-1ЪТКТ, (9)

где К - симметрическая матрица, являющаяся положительно определенным решением матричного алгебраического уравнения Риккати

Q — КЪг-1ЪТКТ + КА + АТКТ = 0. (10)

Существование такого решения гарантируется [6], если система (А; Ъ) является полностью управляемой.

Заметим, что поскольку в рассматриваемой задаче г - скаляр, его значение, без ограничения общности, может быть принято равным единице. При положительной определенности матрицы Q гарантируется [3] асимптотическая устойчивость замкнутой системы при использовании регулятора (7) с коэффициентами (9). Применяя разложение Холесского [8], симметрическую положительно определенную матрицу Q представим в виде

Q = LTDL, (11)

где D - диагональная матрица с положительной диагональю С77, а L - верхняя треугольная с единичной диагональю и элементами 1ц над ней. Таким образом, выбор функционала (8) при использовании представления (11) определяется набором из 10 параметров

Р = (с„ > 0 (7 = 14), 1и (7 = 1,3; 4 > ] > /)). (12)

Свободой их выбора можно распорядиться для оптимизации дополнительных нелинейных критериев качества переходных процессов по этим параметрам. Это использовано в разделе 6.

Заметим, что для значений параметров из табл. 1 условие полной управляемости системы (А; Ъ) из (3)—(6) выполнено, поскольку det(b АЪ А2Ъ А3Ъ) = 3794.14 ф 0.

а22

Модальный регулятор. Поскольку изучаемая динамическая система при рассматриваемых параметрах из табл. 1 является полностью управляемой, с использованием стандартных подходов модального управления [3] можно определить значения коэффициентов регулятора (7), обеспечивающие заданное размещение корней характеристического полинома замкнутой системы (3), (4), (7). А именно, если определить коэффициенты характеристического полинома матрицы А

а (р) =| рЕ— А| = р4 + а3р3 + а2р2 + а1 р + а0,

а также коэффициенты «желаемого» характеристического полинома системы (3), (4), (7)

g (р) = р4 + gзPІ + g2P1 + glP + go, то коэффициенты модального регулятора в каноническом базисе примут вид ki = аі — gi (і = 0,3), а в исходном базисе определятся как

k = Р~ = Р(а — g). (13)

Матрица преобразования Р, определяющая переход к каноническим переменным ~ = Рх, вычисляется по известной формуле

Р = ЯЯ ‘, (14)

где Я , Я - матрицы управляемости для системы (3), (4) в новом и старом базисах:

Я = Ъ АЪ А2Ъ А3~ ; Я = ||ЪАЪА2Ъ А3Ъ ||.

Матрица А и вектор Ъ являются результатом преобразования пары (А; Ъ) в нормальную форму и имеют вид

' 0 1 0 0 > Г 0 ^

0 0 1 0 0

А = , Ъ =

0 0 0 1 0

ч— а0 — а1 — а2 — а3 у V1 у

В качестве набора р свободных параметров метода можно использовать либо вектор коэффициентов g = ^0, gl, g2, g3)T «желаемого» характеристического полинома замкнутой системы (но здесь существует проблема определения адекватной области их варьирования), либо 4 параметра, определяющие требуемое размещение корней X 7 (7 = 1,4) этого полинома. Последний выбор предпочтительнее, поскольку позволяет ввести содержательные ограничения на область варьирования параметров р. Для обеспечения заданного качества протекания переходных процессов положения корней целесообразно варьировать в области комплексной плоскости, имеющей вид пересечения вертикальной полосы и конуса вида:

— л2 < ReX7 < -^1 < 0; |1тX7| < —уReX7,

(і = 1,4).

(15)

Как известно, это обеспечивает в линеаризованной системе экспоненциальное затухание переходных процессов с показателями, лежащими между -Ц\( и -^2^, а также затухание не

-2л/у г

хуже, чем в е ' раз, за период для колебательных составляющих (на каждой из собственных частот переходного процесса).

Имеющаяся свобода в выборе параметров р из (15) используется в разделах 5, 6 для оптимизации переходных процессов по дополнительным нелинейным критериям качества, зависящим от р.

Поскольку существует три варианта типов корней (если классифицировать по количеству комплексно сопряженных пар), то задача оптимизации по параметрам их размещения в области (15) распадается на три отдельные задачи с различными наборами параметров р. Первая задача соответствует наличию двух пар сопряженных корней, при этом принимается

р = (я.еX!,2, 11тX 1,2 |, ReX3,4,11тX3,4 |) (16) и соответствующим образом определяется область их варьирования

D=А ={—'Л2 <Р1,3 <—'Л: <0; \р2,\<—т^Р1,3}. (17) Вторая задача соответствует одной паре сопряженных корней Х3>4, а третья - случаю четырех действительных корней. При этом для второй задачи принимается

Р = (X!, X 2^е X 3,4,11т X 3,4 |), (18)

D = ^ ={ л 2 < Р1,2,3 < —Л1 < 0;

0 < Р 4 <—Р3}

а для третьей -

Р = (X1, X2 , Xз, X3) ,

D = А = {—Л 2 < Р1,2,3,4 <—Л1 < 0}. (21)

После получения решений каждой из трех указанных вспомогательных задач в качестве окончательного решения должно выбираться лучшее из них по целевому критерию.

Линейные матричные неравенства. Их использование дает альтернативный подход к реализации модального управления без возможности точного размещения каждого из корней характеристического уравнения замкнутой системы. Обеспечивается лишь их расположение в заданной области определенного вида. Требуемую область размещения (15) можно описать как LMI-область, т.е. множество точек комплексной плоскости z, порождаемое системой линейных матричных неравенств относительно переменных х = Re z и у = 1т z. Введем обозначение у = аг^(у).

Используя методику, описанную в [7], и опуская промежуточные построения, получим, что для размещения всех собственных чисел

(19)

(20)

критериев качества

матрицы С в соответствующей LMI-области необходимо и достаточно, чтобы нашлась симметрическая матрица X = ХТ, удовлетворяющая системе линейных матричных неравенств

М (С, X) < 0, X > 0. (22)

Знаки неравенств указывают на требуемую знакоопределенность матриц. Матрица М(С, X) для области (15) блочно-диагональна с блоками CX + XCT + 2л 1X; — CX — XCT — 2л 2 X и

(

(CX + XCT )sin у (CX - XCT )cos у

Л

. (23)

— (СХ — ХСТ )cos у (СХ + ХСТ )sin уу

В замкнутой линейной системе (3), (4) с регулятором и = ктх из (7) матрица С примет вид С = А + Ъкт, и первое неравенство в системе (22) сведется к линейному матричному неравенству по отношению к матрицам X > 0 и X = кТХ. Применение к полученной системе (22) вычислительных процедур пакета Ма^аЬ [12] позволяет найти одно из решений Х, X . Тогда искомое значение коэффициентов регулятора определяется как к* = (2*(Х*)- ')Т.

Таким образом, описанный выше подход ставит в соответствие области вида (15) единственное значение к параметров регулятора, что не позволяет непосредственно использовать аппарат линейных матричных неравенств в сочетании с оптимизацией дополнительных нелинейных критериев качества.

4. Дополнительные нелинейные критерии качества

Для реализации неформализованной цели построения регулятора, указанной в конце первого раздела статьи, необходимо ввести набор дополнительных критериев качества переходных процессов для замкнутой нелинейной модели объекта (1), (2), позволяющих количественно оценивать необходимые характеристики.

С точки зрения последующего применения локальной численной оптимизации желательна гладкость до второго порядка этих критериев по параметрам. Выполнить это требование, как правило, не удается, и необходимые критерии

могут оказаться даже разрывными. В этом случае важно, чтобы геометрические места точек разрывов и других нарушений гладкости, по крайней мере, не проходили через области, содержащие решения. Априори неизвестно, выполнено ли это требование, и тогда требуется предварительное численное исследование построенных критериев, поскольку возможность успешного применения методов численной оптимизации не является в этом случае очевидной и требует экспериментального подтверждения. В любом случае способ построения критериев должен преследовать цель максимального соблюдения указанных выше свойств.

Для рассматриваемой задачи было введено пять нелинейных критериев качества (T, H, hi, h2, ^max), описанных ниже, значения которых зависели от вектора параметров k регулятора (1). Параметры k в дальнейшем, как правило, определяются в виде функций k = k(p) от других параметров p, выбор которых рассматривался в разделе 2. Геометрический смысл первых четырех критериев поясняет рис. 2. Их значения находятся по виду зависимости Д(?) = x - x(t) --/’sin 9(t), вычисляемой в процессе численного интегрирования на заданном промежутке [0, T ] нелинейной системы (1), (2). Значение Д(0 определяет величину отклонения по оси x центра масс перемещаемого груза от заданной целевой точки. Критерий T оценивает время сходимости, а критерии H, h1, h2 введены для контроля немонотонности в зависимости Д(0.

Приведем формальное описание введенных критериев качества:

t*, если 3t* е [0,T*]: (3t < t*: |A(t)| > 8 > 0;

T = <! |A(t*)| = 8;Vt e(t*,T*]: |Д(0| <8) (24)

T в противном случае;

0, если Vt е[0, T*]: A(t) > 0,

H = <!- min{A(t): t е [0,T*] и A(t) < 0} (25)

в противном случае.

Обозначим через t точки локальных минимумов зависимости Д(0, формально определяемые условиями Д' (t) = 0, Д" (t) > 0. Тогда

max {Д(0 : t е [0,T*] & Д'(Г) = 0 &

& Д' (Г) > 0 & Д(Г) > 0}, если 3 Д(Г) > 0,

0, иначе.

h =

(2б)

Аналогично, через 7 обозначим точки локальных максимумов зависимости Д(7), формально определяемые условиями Д' (7) = 0, Д" (?) < 0 . Тогда

h2 =

max {A(f): t є [0,T*] & A'(t) = 0 &

& A'(t) < 0 & A(t) > О}, если З A(t) > О,

О, иначе.

(27)

Пятый критерий определим как наибольшее за время наблюдения значение модуля управляющего напряжения, формируемого регулятором (1), т.е.

^тах = таХ{| У(7)|: 7 е [0, Т *] }. (28)

Заметим, что критерий Н характеризует немонотонность в изменении Д(7) только при наличии эффекта перерегулирования, когда Д(7) принимает на части промежутка [0, Т*] отрицательные значения. При этом будет получено Н > > 0. Если же всегда Д(7) > 0, то Н = 0. В этом случае немонотонность контролируется только дополнительными критериями Ь и ^, которые принимают положительные значения при наличии у зависимости Д(7) локальных минимумов и максимумов в области значений Д(7) > 0 на [0, Т*].

Для применения численных методов важным является вопрос о свойствах введенных критериев как функций от определенных ранее параметров р. В разделе 2 было показано, что значения Д(7) непрерывно дифференцируемы по параметрам регулятора к, если угол отклонения ф(7) остается в пределах (-п, п). Однако сама структура критериев (24)-(28) может порождать нарушение их гладкости. Например, при гладком изменении по параметрам формы зависимости Д(7), приводящем к уменьшению значения максимума «Ь» (см. рис. 2) до значения, меньшего 5, произойдет скачкообразное уменьшение значения Т с «с» на «а», т.е. возникнет разрыв.

Критерий Н в силу своей структуры, по крайней мере, непрерывен. Критерии Ь и ^ могут изменять свои значения с нарушением непрерывности в случае исчезновения порождающих эти значения локальных экстремумов при гладком по параметрам изменении кривой Д(7). На приведенных в разделе 5 видах линий равного уровня критериев для ряда двумерных сечений наличие разрыва критерия ^ обнаруживается, например, на рис. 4. Непрерывная дифференцируемость ^ (или ^) может нарушаться при значениях параметров, когда локальный минимум (или максимум) с наибольшим положительным значением выравнивается по этому значению с другим локальным минимумом (максимумом).

Таким образом, следует ожидать кусочную непрерывную дифференцируемость предложенных критериев по параметрам к регулятора (1). Однако при решении рассматриваемых в разделах 5, 6 задач оптимизации в качестве па-

раметров критериев качества применялись не только коэффициенты к, но и значения р, вводимые при использовании LQR или модальных регуляторов. В последних случаях коэффициенты к определялись по параметрам р разного вида согласно (9)-(12) или (16), (17); (18), (19); (20), (21). Возникающая зависимость к = к(р) может оказывать дополнительное влияние на характер гладкости критериев по новым параметрам р.

Проанализируем возможность такого влияния. При использовании модального подхода связь к = к(р) с параметрамир из (16), (17); (18), (19); (20), (21), очевидно, будет гладкой до любого порядка, поскольку к линейно зависит, согласно (13), (14), от набора коэффициентов g полинома, выражаемых по теореме Виета через компоненты желаемых собственных чисел из (15), и не будет оказывать влияния на порядок гладкости дополнительных нелинейных критериев по параметрам (16), (18) и (20).

При использовании LQR-подхода коэффициенты к, получаемые из (9), гладко зависят от параметров р из (12), что вытекает из непрерывной дифференцируемости по элементам Q элементов матрицы К, являющейся единственным положительно определенным решением уравнения (10) в случае полной управляемости пары (А , Ъ). Последнее следует, например, из того, что матрицу К в (10) можно рассматривать как стационарное решение дифференциального матричного уравнения Риккати [4,6], непрерывная дифференцируемость решений которого по элементам Q следует из соответствующей теоремы для систем обыкновенных дифференциальных уравнений.

Таким образом, в рассматриваемой задаче связь к = к(р) не оказывает дополнительного влияния на гладкость критериев качества по параметрам р.

5. Постановка задач оптимизации и предварительный анализ структуры критериев качества

Оптимальные значения параметров р вычислялись из решения задач следующей структуры:

шп {т(р): Н(р) < Н+; \(р) < Ц;

peD

h2(p)<h2+; V(p) < V + }

(29)

где верхние ограничители Н, к+, к+ составляли 3% от уровня начального отклонения |А(0)|. В расчетах раздела 6 было выбрано А(0) = 1, что соответствовало абсолютным значениям ограничителей Н+, к+, к+ = 0.03. Исходя из параметров источника питания исследуемой лабораторной установки принято V+ = 12 В.

Были исследованы три вида задач (29) для определения оптимальных значений параметров регулятора к из (1). В первом случае к определялось зависимостью к = к(р), реализующей -регулятор (9)-(11), где вид вектора р соответствовал (12). Во втором - зависимостью к = = к(р), реализующей модальный регулятор (13), (14), где вектор р имел три вида (16), (18) и (20), что порождало три последовательно решаемые задачи с областями D видов (17), (19), (21) с последующим выбором лучшего из полученных результатов. В третьем случае полагалось к = р, т.е. оптимизация проводилась по прямым параметрам регулятора.

Для вычисления значений критериев Т(р), Н(р), М(р), h2(p), ^ах(р) выполнялись следующие действия:

a) для задач первого и второго видов по параметрам р с использованием линеаризованной модели (3)—(6) строился LQR или модальный регулятор с коэффициентами к = к(р), вычисляемыми по описанным выше правилам (внутренний уровень задачи оптимизации); для задач третьего типа полагалось к = р;

b) при найденных значениях к проводилось численное интегрирование на промежутке [0, Т] системы нелинейных уравнений (1), (2) и по рассчитанным значениям А(?) и v(t) вычислялись критерии (24)-(28).

Интегрирование выполнялось методом Рун-ге-Кутта 4-го порядка с использованием решателя ode45 системы МаїЬаЬ при заданной точности 106.

Для предварительного анализа характера зависимости критериев (24)-(28) от параметров р в задачах оптимизации трех видов были построены линии равного уровня основных критериев (24), (28) и (25), (27) в ряде сечений. Характерные виды показаны на рис. 3-8. Поведение критерия ^ не приводится, поскольку оно структурно подобно поведению ^, хотя с ним и не совпадает. На рисунках, соответствующих ограничиваемым критериям Н, ^, ^, V, темным тоном выделены подобласти в этих сечениях, где значения критериев удовлетворяют указанным в (29) ограничениям.

Рисунки 3, 4 соответствуют задаче (29) первого вида. Сечения построены по двум диагональным элементам матрицы Q, что соответствует параметрам р1 и р3 в (12), значения pj для] = = 5,...10 приняты равными нулю. В данном сечении задача выглядит одноэкстремальной. Следует обратить внимание, что с учетом используемых в (29) ограничений все критерии существенны, поскольку влияют на вид допустимого множества. У критерия ^ в окрестности положительной диагонали области наблюдается разрыв, причем можно ожидать, что оптималь-

ное значение р может оказаться именно на многообразии разрыва. Это может препятствовать работе локальных методов, и по результатам численного решения задача может выглядеть как многоэкстремальная.

Для задачи второго вида (использующей модальный регулятор) согласно данному выше описанию необходимо последовательно решить три задачи вида (29), в которых вектор параметров р принимает вид (16), (18) и (20) и варьируется в областях поиска (17), (19), (21). Решением является лучший результат из трех. На рис. 5, 6 показано характерное поведение в выбранном сечении основных критериев для первой из трех связанных подзадач, соответствующей параметрам (16) из области поиска (17). На рис. 5 (слева) присутствует значительная подобласть стационарности критерия Т, а также контур, при пересечении которого критерий Т быстро изменяется и, возможно, имеет разрывы. Также у Т наблюдается несколько локальных минимумов в сечении (один из них находится вне множества допустимых значений). Критерий ^ также имеет многообразие разрывов и подобласть стационарности. При численной оптимизации обнаруживается значительная зависимость результата от начальной точки поиска, что говорит о многоэкстремальности.

Для задач третьего вида оптимизация в (29) выполняется непосредственно по параметрам регулятора, т.е. к = р. Следует отметить, что поведение линий равного уровня критериев в различных сечениях для задачи этого типа существенно различается. Например, в сечениях по параметрам к], к2 в областях с большими значениями Н, ^ у этих критериев наблюдается значительное количество локальных минимумов. В сечении на рис. 7, 8 задача выглядит как одноэкстремальная. Можно обратить внимание на наличие потерь гладкости критерия Т на рис. 7 вдоль некоторых кривых, а также наличие у него зоны стационарности; на рис. 8 наблюдается разрывность у критерия ^.

6. Результаты вычислительных экспериментов

В данном разделе описаны условия и приведены характерные результаты вычислительных экспериментов по оптимальному выбору параметров регулятора посредством решения задач (29) методами локальной оптимизации системы МаЛаЬ. Процесс оптимального выбора параметров р с использованием общих методов нелинейного математического программирования образует внешний уровень решения задачи. Таким образом, в задачах первого и второго видов использованы двухэтапные процедуры.

200

180

160

140

120

100

80

60

40

20

-12-

-10-

-ІІ-

-ІІ-

Рис. 3. Линии равного уровня критериев Т и Vm

50 100 150 200

по параметрам Q11, Q33 при диагональной форме

матрицы <2 в -подходе (для сечения 222 = 0.5, 244 = 0.5)

Рис. 4. Линии равного уровня критериев Н и h2 по параметрам Qn, Q33 при диагональной форме матрицы 2 В Ь<2К-подходе (для сечения 022 = 0.5, б44 =0.5)

Рис. 5. Поведение критериев Т и Vmax в плоскости параметров Re Х12, |1т Х12| при использовании модального подхода (для сечения Яе Х3 4 = -5.6, |1т Х3 4| = 3.7)

Рис. 6. Поведение критериев Н и в плоскости параметров Re Х12, |1т Х12| при использовании модального подхода (для сечения Re Х34 = -5.6, |1т Х34| = 3.7)

30

25

20

15

10

-25-

-22-

-17-

-12-

ш

І2

-30 -25 -20 -15 -10 -5 0

Рис. 7. Зависимость критериев Т и Утяк непосредственно от параметров къ к3 регулятора (для сечения к2= 1 ,к4 = 5)

Рис. 8. Зависимость критериев Н и h2 непосредственно от параметров к1, к3 регулятора (для сечения к2 = 1, к4 = 5)

Таблица 2

Примеры численной оптимизации по элементам диагональной матрицы Q в двухэтапной процедуре с использованием LQR-подходa

№ Начальная диагональ Q Локально оптимальная диагональ Q Коэффициенты регулятора к Значения критериев

Т Н V

1 10; 0.01; 10; 0.01 110.54; 0.00; 144.00; 0.00 -8.0; 1.4; 12.0; 5.3 2.30 0.005 0.025 0.057 11.99

2 300; 10; 300; 10 288.50; 0.01; 143.90; 7.90 -13.9; 1.9; 11.9; 7.5 3.16 0.000 0.000 0.000 11.99

При численном решении задач (29) были выбраны следующие области поиска D.

Для задач первого вида (при использовании LQЯ-подхода) D задавалась в виде гиперинтервала в пространстве размерности N = 10:

D = [р—,р+] є Я10 (30)

при граничных значениях р~ = (10-1, 10-1, 10-1, 10-1, -10-3,., -10-3 ), р+ = (103,...,103). При использовании на внутреннем уровне модального регулятора (задачи второго вида) в качестве областей поиска D применялись множества из Я4 вида (17), (19), (21), где П1 = 0.1, П2 = 6, у = 1. В задачах третьего типа - множество в виде гиперинтервала в пространстве размерности N = 4: D = [Р ,Р+ ] є Я4 (31)

при граничных значениях Р+ = (+ 20,...,+20).

Указанные задачи (29) решались с использованием функции £тіпсоп системы MatLab [12]. Для учета ограничений в этой функции приме-

нялся либо метод активного набора в сочетании с квазиньютоновским методом решения вспомогательных задач, либо метод внутренней точки (барьерных функций) [8, 9]. В задачах первого типа лучшие результаты обеспечил метод внутренней точки. Параметры останова выбирались стандартными: максимальное нарушение ограничений То1Соп = 10-16, относительная норма текущего шага при останове TolX = 10-10. Как уже отмечалось выше, во всех расчетах целевое значение х принималось равным 1, что соответствовало Д(0) = 1.

Результаты численного решения задач представлены в табл. 2-5. Числа, выделенные в таблицах подчеркиванием, соответствуют найденным локально оптимальным значениям целевого критерия Т. Первые две таблицы соответствуют использованию двухэтапной процедуры с LQR'-регулятором. В табл. 2 оптимизация проводилась только по диагональным элементам

Таблица 3

Примеры численной оптимизации по элементам полной матрицы Q в двухэтапной процедуре с использованием LQR-подхода

№ Показатель Значение показателя

1 Начальное Q (10 0 0 0); (1 0 0); (10 0); (1)

Локально оптимальное Q (209.2535 -9.7812 -84.6526 -13.9990); ( 10.4357 4.5847 11.4506); (134.8462 6.4285); (13.2875)

Собственные числа Х12; Х3,4 (-5.5670 ± І4.0819; -0.7254 ± И.3688)

Коэффициенты регулятора к -14.6293 0.9688 11.6123 5.6851

Критерии: Т, Н, Ы1, Ы2, V 1.8228; 5.610-5; 0; 9.110-7; 11.6123

Число итераций 25

2 Начальное Q (10 1 0 0); (1 1 0); (10 1); (1)

Локально оптимальное Q (224.0763 5.2483 -84.2376 7.5215); (12.1744 2.1011 9.2310); (130.1759 -12.5088); (11.8363)

Собственные числа Х12; Х34 (-6.1133 ± і-3.2440; -0.7034 ± І 1.3606)

Коэффициенты регулятора к -15.2459 0.5699 11.4095 5.7159

Критерии: Т, Н, Ы1, Ы2, V 1.7628; 2.510-7; 1.910-5; 4.5-10-5; 11.4095

Число итераций 43

3 Начальное Q (10 -1 0 0); (1 -1 0); (10 -1); (1)

Локально оптимальное Q (192.7226 -20.6724 -79.0012 -28.5299); (7.2360 6.5206 8.6686); (131.0184 19.6216); (12.3023)

Собственные числа Х12; Х34 (-5.4573 ± І 3.9361; -0.7308 ± І 1.3985)

Коэффициенты регулятора к -13.7902 0.8942 11.4463 5.4388

Критерии: Т, Н, Ы1, Ы2, V 1.7444; 7.6-10-5; 0; 1.910-5; 11.4463

Число итераций 22

4 Начальное Q (10 -1 -1 -1); (1 -1 -1); (10 -1); (1)

Локально оптимальное Q (196.7297 -25.4903 -77.7772 -14.4188); (8.6338 7.7437 8.5840); (132.5863 11.4881); (10.8123)

Собственные числа Х12; Х34 (-5.6074 ± і-3.7874; -0.7235 ± И.3975)

Коэффициенты регулятора к -14.0431 0.8224 11.5146 5.5049

Критерии: Т, Н, Ы1, Ы2, V 1.7039; 5.210-5; 0; 2.510-5; 11.5146

Число итераций 26

Таблица 4

Примеры численной оптимизации в двухэтапной процедуре с использованием модального подхода

№ Показатель Значение показателя

1 Начальные значения собственных чисел Х12; Х34 (-1.0000 ± /'• 1.0000; -1.0000 ± И.0000)

Локально оптимальные значения собственных чисел (-4.8983 ± і-3.0993; -0.6350 ± і-1.2295)

Коэффициенты регулятора к -11.2374 -0.5035 6.5333 2.4440

Критерии: Т, Н, Ы1, Ы2, V 1.8913; 1.510-4; 3.110-9; 6.910-5; 6.5333

Число итераций 79

2 Начальные значения собственных чисел Х12; Х34 (-1.0000 ± і 0.0000; -0.1000 ± і 0.1000)

Локально оптимальные значения собственных чисел (-6.9597 ± і-0.9844; -0.5464 ± і 1.0097)

Коэффициенты регулятора к -18.1490 -1.3711 6.6121 3.5478

Критерии: Т, Н, Ы1, Ы2, V 2.4043; 4.810-5; 0; 2.5-10-4; 6.6121

Число итераций 72

3 Начальные значения собственных чисел Х12; Х34 (-2.0000 ± і 2.0000; -4.0000 ± і 4.0000)

Локально оптимальные значения собственных чисел (-5.5702 ± г4.0757; -0.7060 ± И.3767)

Коэффициенты регулятора к -14.5117 0.8578 11.5789 5.4920

Критерии: Т, Н, Ы1, Ы2, V 1.6132; 2.3-10-4; 1.310-8; 4.2-10-5; 11.5789

Число итераций 85

4 Начальные значения собственных чисел Х12; Х34 (-0.5000 ± /'• 1.0000; -2.0000 ± і 3.0000)

Локально оптимальные значения собственных чисел (-6.3174 ± і-3.7514; -0.6662 ± і-1.2504)

Коэффициенты регулятора к -17.6279 0.5298 11.0026 5.8613

Критерии: Т, Н, Ы1, Ы2, V 1.9098; 1.7-10-6; 4.910-9; 9.410-8; 11.0026

Число итераций 36

Таблица 5

Примеры прямой численной оптимизации по коэффициентам регулятора в одноэтапной процедуре

№ Начальные коэффициенты регулятора к Локально оптимальные коэффициенты регулятора к Значения критериев

Т Н *>1 Ы2 V

1 1; 1; 1; 1 -19.24; 0.70; 12.00; 6.63 -16.00; -1.15; 6.62; 3.29 -11.67; 1.90; 12.00; 5.60 1.89 3-10-7 2-10-7 2-10-7 12.00

2 -1; 1; 1; 1 2.28 2.03 3-10-4 2-10-8 3-10-4 6.61

3 -10; 5; 10; 5 3-10-4 3-10-4 3-10-4 12.00

1

0.8

0.6

0.4

0.2

0

-0.2

-0.4

т

\\xli

Л—-а

о

~г

Рис. 9. Вид переходных процессов по нелинейно преобразованной ошибке наведения Д(?) для нескольких локально оптимальных значений параметров регулятора (а - эксперимент 1 табл. 2, Ь - эксперимент 4 табл. 3, с - эксперимент 3 табл. 4, с! - эксперимент 1 табл. 5)

2000 4000 6000 8000 10000

2000 4000 6000 8000 10000

Рис. 10. Изменения основных критериев качества при переходе от одного локально оптимального размещения собственных чисел ^1,2,3,4 (эксперимент 1 табл. 4) к другому (эксперимент 3 табл. 4)

матрицы Q при нулевых значениях внедиаго-нальных элементов, а в табл. 3 - по всем ее свободным 10 элементам. Сопоставление результатов однозначно указывает на недостаточность варьирования только четырех диагональных элементов матрицы Q. Наблюдаемое различие результатов решения для различных начальных значений в табл. 3 свидетельствует о наличии многоэкстремальности в областях с малыми значениями целевого критерия.

Результаты, приведенные в табл. 4, 5, свидетельствуют о возможности решения задачи оптимизации параметров регулятора как с использованием двухэтапной процедуры, применяющей модальный подход, так и с помощью одноэтапной процедуры. Однако найденный в табл. 5 результат уступает предыдущим. Кроме того, прямая оптимизация непосредственно по параметрам регулятора значительно чаще приводит к попаданию в локальные минимумы с «плохими» значениями целевого критерия Т, чем при использовании двухэтапных процедур.

На рис. 9 приведен вид переходных процессов по отклонению Д(7) для лучших результатов из табл. 2-5. Графики построены в неравномерной шкале с использованием преобразования вида

~ Log((ст + |Д(7 )| Уст)

Д(7) = *\ ' ч \ ’ • ^п(Д(7)) , (32)

Log((ст + 1)/ст)

где было выбрано о = 0.005. Преобразование сохраняет значение единицы и растягивает область в окрестности нуля, позволяя наблюдать малые по модулю значения Д(7).

Для иллюстрации многоэкстремального поведения критериев в окрестности найденных оптимальных значений на рис. 10 на примере задачи (29), решаемой по параметрам (16) в области D = D! из (17), приведены графики изменения четырех основных критериев вдоль отрезка, соединяющего две оптимальные точки из экспериментов 1 и 4 табл. 4. Отрезок пройден с равномерным шагом, равным 10-4 его длины. По оси абсцисс указан номер текущей точки этого отрезка. Ограничиваемые критерии Н и ^, ^ изменяются вдоль отрезка немонотонно, причем ограничения на Гтах, Н и ^, ^ всюду выполнены. Целевой критерий Т вдоль построенного сечения имеет много небольших локальных минимумов, порожденных разрывами, относительная величина которых составляет около 1% от значения этого критерия. Поскольку в рассматриваемом случае наивысшая собственная частота в переходном процессе находится меж-

ду 3 и 4 рад-с-1, то возможные разрывы в значениях Т, порождаемые видом этого критерия (см. обсуждение в разделе 4), могут иметь величину, близкую к 2с (т.е. более 100% от оптимального значения). Поэтому наблюдаемые малые разрывы могут быть порождены только эффектами накопления погрешности с учетом высокой чувствительности значений критерия Т к малым ошибкам интегрирования. Заметим, что применяемые процедуры локальной оптимизации реагируют на эти малые локальные минимумы, что приводит к различным результатам решения задачи (29) локальными методами для разных начальных точек (см. табл. 2-5).

Обсуждение и выводы

В работе на примере конкретной динамической системы подтверждена возможность определения параметров регуляторов путем численной локальной оптимизации дополнительных нелинейных критериев качества и с учетом нелинейного характера модели объекта. Показаны преимущества предложенных двухэтапных процедур по сравнению с прямой оптимизацией по параметрам регулятора, определяемые тем, что при двухэтапном подходе выбор оптимального поведения переходных процессов происходит уже только на классах регуляторов, обеспечивающих асимптотическую устойчивость и приемлемое качество, что исключает ситуации численного интегрирования систем дифференциальных уравнений в условиях неустойчивости.

В вычислительных экспериментах выявлены признаки многоэкстремальности предложенных нелинейных критериев. Наличие основных локальных минимумов связано со свойствами самой задачи. Однако в подобластях со значениями целевого критерия, близкими к оптимальным, наблюдаются множественные локальные минимумы малой относительной глубины, порождаемые разрывами, предположительно вызываемыми влиянием вычислительных погрешностей численного интегрирования. Обнаружены подобласти чрезвычайно быстрого и значительного изменения значений в части критериев, что может свидетельствовать о наличии у

них разрывов, связанных со свойствами задачи. Разрывность препятствует использованию классических методов липшицевой глобальной оптимизации [9-11]. С учетом указанных особенностей к задаче могут быть применены методы с множественным значением константы Липшица (см., например, [11] и приведенную там библиографию) в модификациях, способных учитывать разрывные функциональные ограничения. Задача (29) может использоваться для тестирования новых методов указанного класса.

Работа выполнена в рамках бюджетной темы Н-12-Ю факультета ВМК.

Список литературы

1. Красовский А.А., Поспелов Г.С. Основы автоматики и технической кибернетики. М.: Госэнерго-издат, 1962. 600 с.

2. Неймарк Ю.И. Динамические системы и управляемые процессы. М.: Наука, 1978. 336 с.

3. Методы классической и современной теории автоматического управления: Учебник в 3-х т. Т. 2: Синтез регуляторов и теория оптимизации систем автоматического управления / Под ред. Н.Е. Егупова. М.: Изд-во МГТУ им. Н.Э. Баумана, 2000. 736 с.

4. Ройтенберг Я.Н. Автоматическое управление: Учебное пособие, изд. 2-е, перераб. и дополн. М.: Наука, 1978. 551 с.

5. Первозванский А.А. Курс теории автоматического управления. М.: Наука, 1986. 616 с.

6. Поляк Б.Т., Щербаков П.С. Робастная устойчивость и управление. М.: Наука, 2002. 303 с.

7. Баландин Д.В., Коган М.М. Синтез законов управления на основе линейных матричных неравенств. М.: Физматлит, 2007. 280 с.

8. Гилл Ф., Мюррей У., Райт М. Практическая оптимизация. М.: Мир, 1985. 507 с.

9. Городецкий С.Ю., Гришагин В.А. Нелинейное программирование и многоэкстремальная оптимизация. Н. Новгород: Изд-во ННГУ, 2007. 489 с.

10. Стронгин Р.Г. Численные методы в многоэкстремальных задачах (информационно-статистические алгоритмы). М.: Наука, 1978. 239 с.

11. Сергеев Я.Д., Квасов Д.Е. Диагональные методы глобальной оптимизации. М.: Физматлит, 2008. 352 с.

12. Ма1ЬаЪ R2006/2007/2008 + Slmulmk 5/6/7. Основы применения. 2-е изд., перераб. и доп. / В. П. Дьяконов. М.: Солон-Пресс, 2008. 800 с.

CONSTRUCTING OPTIMAL CONTROLLERS USING NONLINEAR PERFORMANCE CRITERIA ON THE EXAMPLE OF ONE DYNAMIC SYSTEM

S.Yu. Gorodetsky, A.S. Sorokin

The paper investigates the application of numerical optimization methods for selecting controller parameters by using nonlinear quality criteria and taking into account non-linear equations of the object. A two-level procedure is proposed. At the internal level, a modal or LQR controller is constructed using a linearized model with quadratic or modal criteria parameters selected by the numerical method. This method is used at the external level to solve the controller optimization problem using additional nonlinear objective and constraint functions. The study has been performed on a specific example of controlling a nonlinear object.

Keywords: numerical optimization, nonlinear performance criteria, LQR controller, modal controller, numerical simulation.