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

УЧЕНЫЕ ЗАПИСКИ ЦАГИ

Том XXVI 1995 № 1-2

УДК 517.97

КОРРЕКЦИЯ КВАДРАТИЧЕСКОЙ АППРОКСИМАЦИИ НЕВЫПУКЛОЙ ФУНКЦИИ ДЛЯ ИТЕРАЦИОННОЙ ПРОЦЕДУРЫ ОПТИМИЗАЦИИ

*»

А. Э. Бунаков, Н. М. Гревцов, О. Е. Ефимов, И. О. Мельц, А. Б. Трубецкой

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

Методы оптимизации второго порядка могут быть эффективны в вычислительном отношении при достаточно небольших размерах вектора оптимизируемых параметров. Одним из важных достоинств методов второго порядка является их потенциально высокая точность. Это обстоятельство является особенно важным, когда требуется достаточно точно выявить качественные закономерности зависимости оптимальных решений от параметров задачи. Методы второго порядка применимы для минимизации выпуклых функций [1].

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

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

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

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

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

При использовании методов второго порядка для решения ряда задач оптимального управления в дискретной форме для каждого интервала дискретизации требуется нахождение оптимального вектора сравнительно небольшой размерности [2, 3]. Характерной в этом смысле является задача оптимизации траекторий самолета (применение методов второго порядка для решения этой задачи рассматривается в [3]). Для таких задач рассмотренные здесь подходы к построению шага итерационной процедуры являются полезными.

1. Квадратическая аппроксимация минимизируемой функции. При локальной квадратической аппроксимации дважды дифференцируемой скалярной функции / = /(х), х еИ" около точки х, исходной для построения шага Ах итерационной процедуры, приращение функции имеет следующий вид:

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

Если матрица С положительно определенная, необходимое условие минимума (1) есть

откуда следует, что улучшающее приращение вектора х определяется формулой

Это обычная формула метода Ньютона. Введение в эту формулу положительного множителя у, 0 < у < 1, регулирующего длину шага

(1)

где

Д х = - С 1 g.

(2)

Ах = -у С 1g,

обеспечивает убывание / если g ф О [4].

Если матрица С не является положительно определенной, то формула (2) не может быть использована для построения шага итерационной процедуры, даже если |С|*0. Так, например, для функции двух переменных, если

'*0 , с = '1 0 Л

1 о

имеем

\ 7 I 1

А/ = gx^xl + g2^x2+-^x{ -т-Ах|.

Очевидно, что из условий

дА/

ЗАХ!

= 0,

дА / д А х2

= 0

будут получены значения АХ1 И ДХ2, соответствующие минимуму ПО ДХ1 и максимуму по Дх2 функции А/ При формальном использовании (3) получим

Л/ = -У [ 1 - 1 (^1 - ^|)-

При g^ < даже при малых у > 0 получим А/> 0.

2. Построение улучшающего приращения Ах для задачи на безусловный экстремум. Рассмотрим построение улучшающего приращения Ах для задачи/(х) -> шш, х е Я”. В этом случае для А/имеем формулу (1). Осуществим преобразование системы координат у = 51' Ах [5, 6], где квадратная матрица 5 имеет своими столбцами линейно независимые единичные ортонормированные собственные векторы $1, / = 1 ,п матрицы С : 5 = (5Ь ..., 5Й).

Матрица 51 является ортогональной, т. е. 51-1 = 51', в силу чего

Ах = ЗУ Подставляя (4) в (1), получаем

. А/ = ё'у + ^у'Лу,

(4)

(5)

где = g'S, Л = Б'СБ.

Матрица Л является диагональной, а элементами ее диагонали являются собственные числа Х1г... ,Х„ матрицы С [5, 6], соответствующие собственным векторам ...,

Л =

Л,1 . 0

0

'•л у

Эти собственные числа являются действительными, так как матрица С симметричная.

Переход от (1) к (5) соответствует переходу к главным осям квадратической формы с ортами / = 1,я, образующими базис. При этом в (5) во втором слагаемом отсутствуют произведения разноименных компонентов вектора у.

Будем считать, что матрица

где Л] положительно определенная диагональная матрица, содержит все положительные собственные числа, а диагональная матрица Л2 — остальные собственные числа.

Представим теперь вектор у в виде двух блоков:

где у\ объединяет компоненты у, для которых собственные числа положительные, а У2 объединяет остальные компоненты вектора у. При этом (5) можно переписать в виде

где g[ =g'Sl, gj = &Б2, а матрицы и ^ являются соответствующими блоками матрицы £:

Отметим, что такому представлению матриц Л, Б и вектора у соответствуют определенным образом упорядоченные по номерам собственные числа >-1, А.2, . . . , собственные векторы І1, $2, • • • и компоненты у\, У2,... вектора у.

Произведем коррекцию квадратической аппроксимации (6) путем отбрасывания квадратичного слагаемого по уі-

Поскольку Лх положительно определенная матрица, то минимизация (8) в «выпуклом» направлении, т. е. по у\, даст формулу, аналогичную (3):

Очевидно, что теперь можно сконструировать подходящий шаг в «невыпуклом» направлении по У2 в соответствии с градиентной процедурой:

(6)

(7)

(8)

У\ =-УіЛ11Яі.

(9)

Уі - 2 І2-

(10)

В приведенных выше формулах У1 и У2 — положительные параметры, регулирующие длину шага в «выпуклом» и «невыпуклом» направлениях, которые подлежат выбору.

Используя (4) и (7), а также (9) и (10), получаем следующую конструкцию для улучшающего приращения Ах вектора х:

в которой сочетаются методы второго и первого порядков.

Рассмотрим другой возможный вариант коррекции (5), отличный от (8). Этот вариант основан на искусственной замене неположительных собственных чисел матрицы С положительными числами, которые

где А2 — диагональная матрица, заменяющая матрицу Лг-

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

(П)

нужно выбирать с учетом специфики конкретной задачи. При этом матрица А заменяется на матрицу

/ Кбадратичесхая аппвонситиия

х

Квадратическая аппроксимация функции в окрестности точки хо и построение Ах для скорректированной аппроксимации

Используя Л вместо Л, получаем у = - у Л 1 g и вместо (11)

ДХ = -уС-1£,

где

С = БАБ', С"1 =БА~1Б'.

3. Построение улучшающего приращения Ах для задачи с условиями связи. Рассмотрим теперь построение улучшающего приращения Ах вектора х для задачи вида

иш/, X = {х:0 = 0}, хеХ

где вектор-функция

П = П(х) =

задает условия связи, х е Я", р < п.

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

Пусть С2 = О(х) * 0 в точке х. Потребуем, чтобы на шаге Ах вектор П изменился на заданную величину АЛ:

ААх = АП,

(12)

где матрица А = дО./дх размерами р х п имеет своими строками компоненты градиентов а, функций со, (х) в точке х : а, = (Эсо/ (х)/йх)'. При этом А = (а\ ... ар).

Будем считать, что градиенты а, V/= 1 ,р линейно независимы в точке х. Поставим задачу определения Ах следующим образом: требуется определить Ах, обеспечивающий минимум

А/ = g'Ax + ■^■Ах'СДх

при условии (12).

Если матрица С является положительно определенной, то ее решение, обеспечивающее минимум, имеет вид

Ах = -С Хи + А'р),

где

р = -(АС -1А'У1 (АП + АС -1 £),

причем матрица АС ~1А’ размерами р х р при оговоренных допущениях является невырожденной.

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

Введем матрицу проекции градиента [1]

Р = / - А\АА)~1А

и осуществим построение вектора Ах следующим образом:

Ах = Д^ + Д2х, (13)

где

Ді* = Р%, А2х = РАП, (14)

вектор % подлежит определению и Р = А'(АА')~1.

Вектор Дрс является компонентом Ах, лежащим в пересечении линейных условий связи ААх = 0, которое представляет собой подпространство N1 размерности п - р пространства Лп. Вектор Д2х является компонентом Дх в подпространстве N2 размерности р пространства Яп. Базисными векторами пространства N2 могут быть ортонормированные при помощи процедуры Грама—Шмидта [5] векторы а\, ..., ар. Подпространства N1 и N2 взаимно ортогональны.

Можно показать, что

ЛД]Х = АР^ = 0, АА2х = АП. (15)

Введение конструкции (13), (14) с учетом (15) позволяет рассматриваемую задачу свести к задаче определения вектора % в (13), не привлекая условия (12), поскольку оно заведомо выполняется.

Подставляя (13), (14) в (1), получаем

А/ = т+^'рс^+й?, (16)

где

Ё = g + СА2х, (I = £'Д2х + Д2х'СД2х.

Представим теперь неизвестный вектор § через базисные векторы в N1 и N2.

Обозначим ортонормированные базисные векторы в N1 через щ, Ы2, —> ип-р- Эти базисные векторы являются собственными векторами матрицы проекции Р с собственными числами, равными единице. Действительно, так как Ри1 = щ V / = \,п - р, то Ри( = V / = 1,и - р, где

*.!= 1, 1, ... , К-р= 1- (17)

Обозначим ортонормированные базисные векторы в N2 через М„-р+ 1, и„.р+ 2, ■■■, и„. Эти базисные векторы являются собственными векторами матрицы проекции Р с собственными числами, равными нулю. Действительно, так как Ри( = О V / = п - р + \,п, то Л/, = V / = п - р + \,п, где .

*л-р+1 = °> ^п-р+2 = ••• > = °- (18)

Введем теперь матрицы Щ = (щ ... и„_р), и2 = (и„_р+1... ип) и матрицу и = (У\ I 1/2), которая является ортогональной: 1Г1 = II'.

Заметим, что

РЦ = иь ри2 = 0. (19)

Вектор % в новом базисе можно записать теперь в виде

5 = щ, (20)

где

Л =

Лі

Иг

и, в свою очередь, тц и т|2 ■— векторы координат в базисах N1 и N2 соответственно. В связи с этим (20) принимает вид

\ = Ялі + ^2Л2-Подставляя (20) в (16) с учетом (19) и (21), имеем

А/ = I'' Щ Лі + \ Л' и'РСРЩ + </.

(21)

(22)

Получим теперь выражение для матрицы 11'РСРи. Учитывая (17) и (18), имеем

'/

0

и ри -

Поскольку 1ЛГ = I,

и>рсри = (и'Ри)(и'Си)(и'РЦ) =

где матрица С\ размерами (и - р) х (п - р) есть блок матрицы

с = иси =

Сі 0

0 0

(А с2)

Используя этот результат в (22), получаем

АҐ = ї,и1 лі + -| лівілі + 4-

(23)

Остается теперь перейти к главным осям квадратичной формы в N1. Для этого полагаем

Л1 = Щ, (24)

где матрица IV = (м^ ... мп_р) имеет своими столбцами ортонормиро-ванные собственные векторы матрицы С1. Этим собственным векторам

соответствуют собственные числа щ, •••> Ч-я-р матрицы Сх и диаго-

нальная матрица

составленная из этих чисел. Подставив (24) в (23), получим квадратическую аппроксимацию Д/в N1 в главных осях:

Учитывая, что матрица М диагональная, нетрудно осуществить, если необходимо, коррекцию квадратической формы (25), аналогичную той, которая описана при переходе от (6) к (8), и получить конструкцию для улучшающего приращения Дрс. В этом случае при соответствующей нумерации ..., У^п-р и ць цл-р матрица М может быть представлена в виде

где М\ — диагональная матрица, содержащая положительные собственные числа, а М2 — диагональная матрица, содержащая остальные собственные числа.

Разбивая соответствующим образом вектор £ и матрицу IV на два блока:

Теперь на основании (14) с учетом (19), (20) и (24) получаем конструкцию

(25)

\у2;

из (25) после коррекции получаем

Д/ = АКі+!сіМ1Сі+АК2+4,

где И{ = £'1^, Аї=ї'и21У.

Форма (26) аналогична (8) и, аналогично (9), (11), имеем

Сі = - її Щ\, & = -чгЬ-

(26)

(27)

где, как и ранее, уі и у2 — константы, подлежащие выбору.

Если использовать другой вид коррекции, который рассмотрен в конце предыдущего раздела, то матрица М будет заменена на матрицу М. При этом формула (27) примет вид

Д]Х = -циЯУАГ1 IV' Щ g.

Это выражение можно записать в форме, соответствующей методу Ньютона:

где

Ахх=-уС lg,

С = ЩУУМЦГ' Щ, С1 = и^М-11Г' щ.

4. Минимизация функции двух переменных. Рассмотрим в качестве примера применение изложенного подхода для задачи минимизации функции двух переменных / = Дрс), х = (х1; *2)' сначала без дополнительных ограничений.

Компоненты вектора градиента g и элементы матрицы вторых производных С

8 =

, с = 'СЦ с\2 ^

<£2, ,с12 с22,

вычисляются по формулам

3/

д/

_ Э2/

£1 Яу ’ ^ ~ Рпг ’ ~ =„2 ’ -

ОХ\ 0X7 дх{

а2/

дV

дх]дх2 ’ °22 дх%

Корни характеристического уравнения, которые здесь обозначены через VI, ь2, определяются по формуле

.. с11+с22 ^ %2 =------5---- ±1

г С\1 ~с22 ")2 , „2 ч 2 ) +С»

Уравнение для определения собственных векторов, которые здесь обозначим через г\, г2 и которые соответствуют собственным числам матрицы С или корням характеристического уравнения оь у2. имеет вид

С11 - с\2 . с12 с22 - ь\

= 0, / = 1,2.

При с12 = 0 и сц = С22 корни характеристического уравнения кратные: VI — У2 = С и = С22. В качестве собственных векторов можно принять два любых единичных вектора, например,

'1> '0"

г\ = > Ъ = ,1,

При С\2 * 0 всегда у2ти собственные векторы могут быть вычислены по формулам

1\ =

СОвф

^БШф

(-вШф^ ^ = \> ^ СОвф )

где

вШф = -

си~°1

I 2 7 * СО&<р = -

+с12 1(С11-Ч) + с12

* с12

При расчете в соответствии с (11) требуется построить матрицы $1, 5г и Лх. Рассмотрим различные случаи их построения и приведем соответствующие формулы для расчета Ах.

1. При «1 > 0, иг > 0 полагаем =v1, к2 =и2, *1 =Г\, $2 =г2>

А-х 0 "1

- (*1> *2)> *^2 - А. -

О А.

2) 1-1 <

Из (11) имеем Ах = -у1 или в скалярной форме:

Л*1 = -У1 т^ + -г^ £1+ + &

( _2 _2 > (

511 , 21 Л + \

^2

/ \ Г .2 „2 > “

*11*12 , *21*22 а + *21 , *22 #2

1 я.1 Х2 , ^1

Дх2 = - у 1

2. При »1 > 0,' 02 £ 0 полагаем Ах = 1»х, Х2 =у2) =г15 52 =г2>

= $1, Б2 = 52, Л = А.х. ИЗ (11) ИМееМ Ах = - — $х *1 £ - У2*2 *^ £ 11111 в

Я.1

скалярной форме:

Л*1 = - Т1 (*11 «1 + 5П *12 82 ) - У 2 (*21Й + *21 *22 82),

К1

А^2 = - -р- (*хх *12 81 + *?2 £2) -У2 («2! *22 81 + 42 й)-Л-1

3. При <, 0, ь2 > 0 полагаем А^ = ь2, А.2 = »1, *1 = г2, $2 = г1; ^ = ^1, ^ = 52, Л = А^. Формулы для расчета Ах остаются теми же, что и в предыдущем случае.

4. При VI 5 0, у2 < 0 полагаем в (11) 5х = 0 или, что то же, Ух = 0, ^х = гх, $2 - г2, Б2 = (^х, д2)> после чего из (11) имеем

Дх = -у252^Я = -у2Я,

так как Б2Б2 = /.

При введении в эту же задачу ограничения типа равенства О = =П(хъ *2) = 0 требуется вычислить матрицы Щ, Щ, Щ и М\ в (27). Можно положить

/ 8ІПф (COS<p>

U\ = «1 = , U2 = «2 = •

l^-coscpj ^sm<p;

где

вІПф = -

а2

, coscp =

Н+а:

, ai = дП/дхх, а2 = дП/ дх2■

Матрица С\ для рассматриваемого примера является скаляром и вычисляется по формуле

Q = Сц sin2 ф - 2сП sin ф COS ф + С22 cos2 ф,

а ее собственный вектор и собственное число принимают значения

w= 1, ц= С\.

Матрица проекции Р вычисляется по формуле

Р =

sin2 ф - sin ф COS ф

—sin ф COS ф COS2 ф

При (х > 0 полагаем W\ = w, Щ = 0 и из (27) имеем

ДіХ = -Уі

вІПф

-совф) ц '

У1 ( • 2 Sin ф -вІПфСОвф

и ^-ЗІПфСОвф 2 COS ф

где

и

g = g + СРДф = g + СД2х

Р =

\Р2у

а1 а2 Р\ " 2 1 ’ Pi ~ 2 5"

При ДП = 0 полученная формула для Дрс является при уі = 1 аналогом обычной формулы метода Ньютона.

При ц ^ 0 полагаем Щ = 0, Щ = м> и из (27) имеем

ДіХ = -у2

' вІПф Л ,-С08ф,

1-1 • (sin ф|—COS ф)^ = — у 2Р & •

При АП = 0 эта формула соответствует стандартному методу проекции градиента.

Компонент А** для ц > 0 и ц 5 0 вычисляется по формуле

Л2х = РАф.

Работа выполнена при подцержке Российского фонда фундаментальных исследований.

ЛИТЕРАТУРА

1. Кюнци Г. П., Крелле В. Нелинейное щкираммирование.—М.: Советское радио.—1965.

2. Брайсон А., Хо Ю-Ш и. Прикладная теория оптимального управления.—М.: Мир.—1972.

3. Ефимов О. Е. Метод второго порядка оптимизации нелинейных динамических систем и его использование для оптимизации траекторий самолета // Ученые записки ЦАГИ.—1991. Т. 22, № 3.

4. Моисеев Н. Н., Иванилов Ю. П., Столяров Е. М. Методы

. оптимизации.—М.: Наука,—1978.

5. Ган'тмахер Ф. Р. Теория матриц.—М.: ГИТТЛ.—1953.

6. Стражева Т. В., Мелкумов В. С. Векторно-матричные методы в механике полета.—М.: Машиностроение,—1973.

7. Мельц И. О. Оптимизация управляющих функций и параметров динамических систем // Труды ЦАГИ.—1969. Вып. 1149.

8 Fletcher R., Powell М. J. D. Rapidly convergent descent method for minimization // Computer Journal.—1963. Vol. 6, N 2.

Рукопись поступила 6/ХИ 1993 г.