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

УДК 517.977.58

БИЛИНЕЙНЫЕ ДИНАМИЧЕСКИЕ СИСТЕМЫ: ИССЛЕДОВАНИЕ ИТЕРАТИВНЫХ МЕТВДОВ ОПТИМИЗАЦИИ

О.В. Батурина

Рассмотрена задача оптимального управления динамическими системами, описываемыми билинейными дифференциальными уравнениями и квадратичным критерием оптимальности. Для такой задачи приведено сравнение эффективности метода Кротова и градиентного метода последовательного улучшения управления на основе их численной реализации.

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

ВВЕДЕНИЕ

Основополагающий результат теории оптимального управления — принцип максимума Л.С. Пон-трягина [1], который сводит решение задачи оптимального управления к решению краевой задачи для системы обыкновенных дифференциальных уравнений. Алгоритмические трудности на этом пути возрастают экспоненциально с увеличением размерности задачи. Поэтому наряду с уравнениями принципа максимума широкое применение нашли методы непосредственной оптимизации управления, в частности методы его последовательного улучшения [2—6].

Данная статья посвящена экспериментальному исследованию и сравнительному анализу итерационных методов решения задачи оптимального управления динамическими системами, которые описываются билинейными уравнениями. Этот класс задач играет существенную роль в теории оптимального управления. Так, к нему относится класс актуальных задач управления квантовыми системами [7]. Теоретическое исследование свойств рассматриваемой в статье задачи и методов ее решения выполнено в работе [8].

1. ПОСТАНОВКА ЗАДАЧИ

Рассматривается задача оптимального управления:

т

/(V) = {/0[?, х(0, и(?)]<* + Д(х(Т)) ^ шт, (1)

о

х = / [?, х(?), и(?)], х(0) = х0, х е Е”,

и е и с Я, (2)

т~'”

где Е — евклидово пространство размерности п с

”

произведением ху = ^ х 1у1, число Т> 0, вектор х0,

I = 1

вектор-функция /(?, х, и), функции Д(х), /0(?, х, и) и компактное множество и заданы. Множество пар вектор-функций V = (х(-), и(-)), определенных на отрезке [0, Т] и удовлетворяющих условиям (2), будем называть множеством допустимых процессов и обозначать через Б. При этом функция и(?) — кусочно-непрерырывна, х(?) — непрерывна и кусочно-дифференцируема. Предполагается, что множество Б непустое. Выделим из данной задачи подзадачу улучшения.

Пусть имеется допустимый неоптимальный процесс у0 = (х0(?), и0(?)) е Б. Требуется найти процесс V = (х(?), и(?)) е Б, такой что /(у) < 1(у0). Повторяя операцию улучшения, построим улучшающую последовательность допустимых процессов {^}с Б, для которой /(^ + 1) < Д^).

Рассмотрим частный случай задачи (1), (2): т

/(V) = |ри2(?)Л + (х(Т), Хх(Т)) ^ шт, (3)

0

х = (А + иВ)х, х(0) = х0, и е [ а, а], (4)

где X, А и В — матрицы размерности п х п, в 1 0.

т

Слагаемое | ри2^)^ в данной задаче имеет регуля-

0

ризующий смысл. Главная цель — минимизация функции Д(х(Т)) = (х(Т), Хх(Т)), при этом часть ресурса управления расходуется на минимизацию

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

2. МЕТОДЫ РЕШЕНИЯ ЗАДАЧИ УЛУЧШЕНИЯ

В данной статье рассматриваются градиентный метод [5, 6] и метод Кротова (глобальный метод по терминологии работы [4]) решения поставленной задачи.

Введем в рассмотрение конструкции из работы [4]:

Н(?, у, х, и) = у (А + Ви)х — ри ,

Я(?, х, и) = Н(?, у(?), х, и) + ^х,

С(х) = (х, Хх) + у(Т )х — у(0)х0, и (?, х) = аг§ шахЯ(?, х, и), и е [—а, а],

где у — вектор пространства Е”, функция у(?) — решение сопряженной системы

ДД х0(?), и0(?)) = ау/^ + (Ат + Вти0)у(?) = 0, (5)

СХ(х,(Т)) = у(Т) + Д^Т)) =

= у(Т) + 2Хх0(Т) = 0. (6)

Здесь и далее Ях, (?х, Дх и Яи — частные производные по соответствующим переменным.

2.1. Градиентный метод

Для унификации записи опишем градиентный метод в терминах метода Кротова.

Известны управление и0(?) и, в силу уравнений процесса и начальных условий (4), траектория х0(?). Пусть частная производная Яи(?, х0(?), и0(?)) не равна нулю тождественно на отрезке [0, Т]. Новое управление и(?) будем искать достаточно близким к управлению и0(?) так, чтобы знак разности А/ = /(V) — /(у0) совпадал со знаком ее главной линейной части

+

б/ = Gx(x0(T))6x(T) - J(Rx(t, X0(t), u0(t))6x(t)

0

+ Ru(t, x0(t), u0(t))6u(t))dt,

где 8х = х(?) — х0(?), 8и = и(?) — и0(?).

Определим семейство управлений и(?, е) = = и0 + Еу(?), т. е. 8и = Еу(?), с параметром е > 0, соответствующее ему семейство траекторий х(?, е) и

семейство процессов у(е) = (х(?, е), и(?, е)). Существует е > 0 такое, что /(у) < /(у0).

Алгоритм улучшения.

Шаг 1. Решая задачу Коши (5), (6), находим сопряженную функцию у(?) и соответственно имеем

б/ = -є JRu(t, X0(t), u0(t))y(t)dt.

(У)

Шаг 2. Функцию у(?) выбираем так, чтобы правая часть выражения (7) была положительной.

Шаг 3. Выбираем произвольное s > 0 и новое управление u(t, s), из системы (4) определяем траекторию x(t, s).

Шаг 4. Вычисляем новое значение функционала I(v(s)), соответствующее найденным управлению и траектории. Если улучшения не произошло, то уменьшаем s и повторяем шаг 3.

Основные недостатки градиентного метода:

— необходимость подбирать настроечный параметр s на каждой итерации, если выбрать s постоянным и достаточно малым для получения улучшения на каждой итерации, то скорость сходимости существенно снижается.

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

2.2. Метод Кротова

Известны начальное управление u0(t) и начальная траектория x0(t), полученная из системы (4).

Алгоритм улучшения, основанный на достаточных условиях оптимальности, выглядит следующим образом [4].

Шаг 1. Решая задачу Коши (5), (6) с известными управлением u0 и траекторией x0, строим сопряженную функцию y(t).

Шаг 2. Замыкая систему (4) управлением U (t, x), находим траекторию x(t) и соответствующее управление u(t) = U (t,x(t)).

Улучшение функционала гарантировано в том случае, если функция F(x) вогнутая. В этом случае функция R(t, x, u) выпуклая по x, а функция G(x) — вогнутая, т. е. выполнение соотношений

Rx(t, *0(<), M0(t)) = 0, Gx(x0(T)) = 0

обеспечивает выполнение условий теоремы Кротова

R(t, x0(t), u0(t)) = minR(t, x(t), u0(t)), G(x0(T)) = max G(x(T)).

0

Кроме того, существует актуальный класс задач [7], в котором невогнутая функция Р(х) = (х, Хх) может быть модифицирована таким образом, что глобальный алгоритм улучшения становится применимым.

Рассмотрим частный случай системы (4), когда матрицы А и В имеют следующую блочную структуру:

А =

0 Ра

-Ра 0

В =

0 Рв -Рв 0

где РА и Рв — симметричные матрицы размерности п/2* п/2.

В этом случае, независимо от выбора управления, значение произведения (х, х) остается постоянным, поэтому если функция Д(х) не является вогнутой, то ее следует преобразовать путем прибавления к матрице X слагаемого (—аЕ) с достаточно большим значением а > 0, где Е — единичная матрица. В работе [8] показана возможность такого преобразования. Будем считать, что матрица X в данном случае такова, что Д(х) — вогнутая функция.

Алгоритмы решения задачи улучшения методом Кротова при р = 0 и р > 0 существенно различаются. Использование регуляризующего слага-т

емого р |и2(?)Л при р > 0 упрощает итерационные

0

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

Рассмотрим эти два варианта задания целевого функционала.

Пусть р > 0. Функция Н является квадратичной функцией по и с отрицательным коэффициентом при и2. Обозначим К(?, х) = у(?)Вх(?). Новое управление, построенное на втором шаге алгоритма, будет иметь вид:

ї (ґ, х) =

К

-а, если — < -а,

2 в ,

к к

если -а < — < а,

2в

2в

к

а, если — > а. 2 в

В этом случае правая часть системы (4) при

и(?) = и (?, х) является непрерывной и решение задачи Коши не вызывает сложностей.

Пусть р = 0. Функция Н линейна по и, коэффициент при и равен К.

Исходя из этого, управление, отвечающее шагу 2 описанного алгоритма, будет иметь вид:

ї (ґ, х) =

-а, если К < 0,

любое из [-а, а], если К = 0, а, если К > 0.

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

В этом случае найденное управление будет состоять из нескольких участков, на каждом из которых оно будет принимать либо одно из своих граничных значений, либо такие значения, чтобы обеспечить продолжение движения вдоль поверхности К(?, х) = 0. Введем понятие особого режима. Управление называется особым на отрезке [0, Т ], если на этом отрезке К(?, х(?)) = 0 для всех ?. Решая уравнение аК/Л = 0, получаем

и(0 = и0(0 + ( (АВ -ВА )х( t) , у ( 0 ) .

(Вх( О, у( 0)

Таким образом, на участках особого режима следует решать систему

х = Га+Г и0( о +( (АВ -ВА)х, у ( t) л В)*

^ ^ (В2х,у(0) У *

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

Возможно другое решение задачи улучшения при р = 0, когда решение не включает в себя особое управление. В этом случае

ї (ґ, х) =

-а, если К < 0, а, если К > 0.

При решении задачи данным методом найденное управление принимает только граничные значения. Участки особого режима аппроксимируются частыми переключениями управления с —а на а и обратно (так называемый скользящий режим). При этом значение К(?, х) будет близко к нулю. При уменьшении шага интегрирования число этих переключений возрастает. Полученное таким образом управление достаточно сложно использовать на практике, поэтому на таких отрезках более целесообразно особое управление, полученное при помощи дифференцирования равенства К(?, х) = 0.

3. ПРИМЕРЫ

Далее рассматривается конкретная задача оптимального управления, для которой реализуются описанные итеративные методы решения. Исследовано изменение скорости сходимости методов при варьировании исходных данных задачи: отрезка времени Т, числа a, задающего границы управления, параметра регуляризации р.

Рассмотрим задачу (3), (4) со следующими исходными данными:

L =

-1 0 0 0

0 -2 0 0

0 0 -3 0

1 0 0 0 -1

Л

B =

0 0 -11

0 0 12

1 -10 0

-1 -2 0 0 )

A =

Л

0 0 1 -2'

0 0 -2 -1

-1 2 0 0

2 1 0 0 )

ґ 1 1

1 -1 1 )

В данной задаче матрица L отрицательно определенная, поэтому функция F(x(T)) — вогнутая, и к задаче может быть применен глобальный метод улучшения. Кроме того, матрицы A и B имеют блочную струк-

Таблица 1

Варьирование времени Т

Номер итерации Т = 0,5 Т = 3 Т = 6

0 -6,5945 -6,6382 -7,4850

1 -7,2889 -8,8029 -8,5175

2 -7,7260 -10,8684 -8,9996

3 -8,0689 -11,0557 -9,7743

4 -8,3702 -11,1628 -10,6686

5 -8,6404 -11,5634 -11,2690

6 -8,8816 -11,6864 -11,3286

7 -9,0869 -11,7846 -11,9899

8 -9,2653 -11,826 -11,9928

9 -9,4140 -11,8179 -11,9945

Таблица 2

Варьирование параметра р

Номер итерации р = 0 в = 0,1 в = 0,5 в = 1

0 -6,5945 -6,5945 -6,5945 -6,5945

1 -7,2889 -7,2578 -7,1564 -6,9970

2 -7,7260 -7,7260 -7,5375 -7,2881

3 -8,0689 -8,0535 -7,8494 -7,5315

4 -8,3702 -8,3592 -8,1508 -7,7623

5 -8,6404 -8,6208 -8,4272 -7,9935

6 -8,8816 -8,8097 -8,6833 -8,2229

7 -9,0869 -9,0710 -8,9322 -8,4407

8 -9,2653 -9,2533 -9,1448 -8,6361

9 -9,4140 -9,4152 -9,3589 -8,8017

Таблица 3

Варьирование допустимых границ управления

Номер итерации а = 2 а = 3 а = 6 а = 10

0 -6,5645 -6,5645 -6,5645 -6,5645

1 -7,2103 -7,2889 -7,3406 -8,2363

2 -7,5404 -7,7260 -7,8619 -8,9838

3 -7,7803 -8,0689 -8,3341 -9,2730

4 -7,9787 -8,3702 -8,7236 -9,5488

5 -8,1422 -8,6404 -9,0462 -9,6718

6 -8,2789 -8,8816 -9,2672 -9,8144

7 -8,3916 -9,0869 -9,4532 -9,9429

8 -8,4823 -9,2653 -9,6175 -10,0456

9 -8,5551 -9,4140 -9,9166 -10,1637

- -8,8241 -10,1606 -11,7656 -11,9317

туру и выражение (х, x) остается постоянным для различных способов задания управления. Поэтому метод Кротова применим при решении задачи минимизации функционала Ду) — а(х, Ex), а > 0.

В ходе расчетов интервал времени [0, Т ] аппроксимировался равномерной сеткой с шагом 0,0005, при этом значения управления, траектории и сопряженной траектории в узлах хранились в массивах. Их значения в промежуточных точках аппроксимировались линейными сплайнами. Начальное управление было задано постоянным и0 = 1 на всем отрезке [0, Т ]. При решении задачи Коши применялся метод Рунге — Кутта четвертого порядка. Задача решалась методом Кротова. Варьировалась правая граница промежутка времени Т (табл. 1), при этом и = [—3, 3], р = 0. Здесь и далее в клетках таблиц указаны значения функционала I.

Можно видеть, что при увеличении отрезка времени значение целевого функционала, которого удается достичь, становится меньше. При этом с увеличением Т оптимального значения можно достичь на более ранних итерациях. При больших значениях Т достигается абсолютный минимум квадратичной формы (х, Lx) на сфере (х, х) = (х0, х0). В данной задаче эта величина принимает значение —12.

Рассмотрим изменение процедуры улучшения функционала при варьировании параметра регуляризации р (табл. 2). Конечный момент времени Т = 0,5, границы управления задаются числом а = 3. При уменьшении р процесс улучшения становится более эффективным и позволяет точнее приблизиться к оптимальному значению функционала для данной задачи. Этот способ позволяет снизить трудоемкость вычислений, но ухудшает эффективность метода.

Рассмотрим как изменяется итерационный процесс при варьировании допустимых границ управления, при этом Т = 0,5, р = 0. В последней строке табл. 3 приведено наименьшее значение функционала, которого удается достичь. Отметим, что при увеличении границ управления оптимальное значение целевого функционала уменьшается.

В табл. 4 проведено сравнение решений задачи методом Кротова и градиентным методом. Конечный момент времени Т = 0,5, границы управления задаются числом а = 3, параметр регуляризации р = 0. Можно

заметить, что значения функционала I = —8,8 методом Кротова удается достичь на шестой итерации, а градиентный метод приводит к нему только на одиннадцатой.

Увеличим допустимые границы управления, пусть а = 10, при этом Т = 0,5, р = 0 (табл. 5). При расширении отрезка допустимых значений управления до [—10, 10] метод Кротова остается также более эффективным, значение функционала —9,6 при его примене-

Таблица 4

Сравнение метода Кротова и градиентного метода

Номер итерации Метод Кротова Градиентный метод

0 —б,5б45 -б,5б45

1 —7,2889 -б,8024

2 —7,72б0 -7,023б

3 —8,0б89 -7,2537

4 —8,3702 -7,4873

5 —8,б404 -7,7188

б -8,881б -7,9433

7 —9,08б9 —8,15б9

8 —9,2б53 -8,35б8

9 —9,4140 -8,541б

Рис. 1. График K(t, x) при переходе на особый режим

Таблица 5

Сравнение метода Кротова и градиентного метода при увеличении допустимых границ управления

Номер итерации Метод Кротова Градиентный метод

0 -б,5б45 -б,5б45

1 -8,23б3 -7,0158

2 -8,9838 -7,4758

3 -9,2730 -7,9341

4 -9,5488 -8,3537

5 -9,б718 -8,7139

6 -9,8144 —9,011б

7 -9,9429 -9,2534

8 -10,045б -9,448б

9 -10,1бЗ7 -9,б0бЗ

Таблица 6

Модификация функции F(x)

Номер итерации I = (x, Lx) I = (x, (L — 2E)x)

0 —7,7б01 -7,7б01

1 —8,4033 -8,1487

2 —9,9071 -9,7277

3 -10,8199 -10,8279

4 -11,2210 -11,4753

5 -11,4234 — 11,51б5

6 -11,5449 — 11,727б

7 -11,7511 —11,8034

8 — 11,82б7 -11,8498

9 — 11,8б77 -11,8422

10 -11,9100 -11,9421

Рис. 2. График управления

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

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

Рассмотрим применение метода Кротова к задаче с измененным функционалом, полученным с помощью прибавления к матрице L матрицы (—аЕ). В рассматриваемой задаче матрица L отрицательно определенная, поэтому при изменении критерия оптимальности вогнутость функции F(x) сохраняется. Новое управление, полученное в результате улучшения, на каждом шаге не единственно, но при этом оптимальное значение функционала не должно измениться.

Пусть Т = 4, а = 3, р = 0. В качестве начального управления было выбрано постоянное управление и0 = 0,3 на всем отрезке [0, Т] (табл. 6). Видно, что при изменении функционала путем прибавления к матрице L матрицы (—аЕ) процедура улучшения остается эффективной, но при различных а может изменяться скорость улучшения функционала.

Опишем более подробно одну итерацию метода Кротова при наличии особого режима. Начальное управление выберем постоянным и0 = 0,3 на всем отрезке [0, Т ], р = 0, а = 3, Т = 0,5.

Начиная с ? = 0,035, реализуется особый режим, где К(1, х) принимает значения, близкие к нулю, но с разными знаками. На рис. 1 более подробно выделен отрезок участка особого режима. При уменьшении шага

CONTROL SCIENCES № б • 2G1G

Таблица 7

Модификация матрицы L для невогнутой функции F(x)

Номер итерации Градиентный метод с матрицей L Метод Кротова с матрицей L Метод Кротова с матрицей Lj

0 б,5945 б,5945 б,5945

1 5,7352 8,5095 5,5292

2 5,3442 9.3013 4,8725

3 5,087б 7,7084 4,4702

4 4,8б75 9,1927 4,2824

5 4,б843 8,0539 4,1819

б 4,5385 9,1579 4,1319

7 4,4255 8,081б 4,1003

8 4,3392 9,1503 4,0807

9 4,2735 8,1003 4,0б8б

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

График управления, полученный с помощью глобального метода улучшения, представлен на рис. 2.

При Ї < 0,035 значение Щ, х) > 0, и* = аі§шахН(ґ, х) = = 3, и є [-3, 3]. При Ї > 0,035 реализуется особый режим.

Рассмотрим задачу, в которой функция F(x) невогнутая. Положим

L =

1 0 0 0 А 0 2 0 0 0 0 3 0 V 0 0 0 1

Т = 0,5, а = 3, р = 0, начальное управление и0 = 1 на всем отрезке времени [0, 0,5]. Решая задачу с модифицированной матрицей

L, = L - 3E =

-2

0

0

0

л

0 0 0 -10 0 0 0 0 0 0 -2 J

получим улучшение исходного функционала с матрицей L (табл. 7).

В исходной задаче с невогнутой функцией F(x) глобальный метод не дает улучшения, но в модифицированной задаче он показывает лучший результат, чем градиентный метод. Оптимальное значение для данной задачи I = 4,03. Градиентным методом оно достигается на 43-й итерации, а глобальный метод с модифицированной матрицей L1 позволяет получить его на 15-й.

ЗАКЛЮЧЕНИЕ

Экспериментально исследована задача оптимального управления динамическими системами, которые описываются билинейными дифференциальными уравнениями и квадратичным критерием оптимальности. Рассмотрен пример, когда

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

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

ЛИТЕРАТУРА

1. Математическая теория оптимальных процессов / Л.С. Понтрягин, В.Г. Болтянский, Р.В. Гамкрелидзе, Е.Ф. Мищенко. — М.: Наука, 1976. — 384 с.

2. Черноусько Ф.Л., Колмановский В.Б. Вычислительные и приближенные методы оптимального управления // Итоги науки и техники / Сер. Мат. анал., 14, ВИНИТИ, М., 1977. — С. 101—166.

3. Федоренко Р.П. Приближенное решение задач оптимального управления. — М.: Наука, 1978. — 488 с.

4. Krotov V.F. Global Methods in Optimal Control Theory. — N.-Y.: Marcel Бескег Inc., 1996. — 382 p.

5. Kelley H. J. Gradient theory of optimal flight paths // ARS Journal. — 1960. — N 30. — P. 947—954.

6. Энеев Т.М. О применении градиентного метода в задачах теории оптимального управления // Космические исследования. — 1966. — Т. 4, № 5. — С. 651—669.

7. Кротов В.Ф. Управление квантовыми системами и некоторые идеи оптимального управления // Автоматика и телемеханика. — 2009. — № 3. — С. 15—23.

8. Батурина О.В., Булатов А.В., Кротов В.Ф. К оптимизации линейных систем с управляемыми коэффициентами // ArXiv 1004.3237v1 [math.OC] 19 Apr 2010. — 18 с.

Статья представлена к публикации членом редколлегии

В.Ю. Рутковским.

Батурина Ольга Владимировна — вед. инженер,

Институт проблем управления им. В.А. Трапезникова РАН,

г. Москва, ®(495) 334-91-59, И [email protected].