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

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

УДК 519.853.4

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

© 2012 г. К.А. Баркалов

Нижегородский госуниверситет им. Н.И. Лобачевского

Ъагка1оу@йдр. unn.ru

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

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

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

Введение

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

Алгоритмы из работ [1-4] основаны на предположении о липшицевости всех функционалов задачи, что типично и для многих других подходов (см., например, [5, 6]). Но если для функций, входящих в постановку задачи оптимизации, известны значения производных, их использование может существенно ускорить сходимость. Так, в [7-10] предложены эффективные алгоритмы, использующие значения производной целевой функции в задачах глобальной оптимизации без ограничений. В то же время остается открытым вопрос использования значений производных в алгоритмах условной глобальной оптимизации. Настоящая статья посвящена разработке и обоснованию схемы такого алгоритма.

Проводимое рассмотрение ограничено одномерным случаем, обсуждение алгоритма для

решения многомерных задач глобальной условной оптимизации планируется осуществить в последующих публикациях. Здесь лишь отметим, что обобщение нового алгоритма на многомерный случай может быть основано на сведении решения многомерной задачи к решению эквивалентной ей одномерной задачи (или же параллельному решению набора таких одномерных задач). Соответствующая редукция основана на использовании разверток единичного отрезка вещественной оси на гиперкуб. Роль таких разверток играют непрерывные однозначные отображения типа кривой Пеано (кривые, заполняющие пространство) и их обобщения, называемые «множественными развертками». Последовательные алгоритмы построения указанных отображений подробно изложены в [2, 11], вопросы параллельного использования множественных разверток обсуждаются в [12-14].

1. Постановка задачи

Рассмотрим одномерную задачу глобальной оптимизации

ф(х*) = тт{ф(х): х е [а,Ъ\ gj{x)<0, 1 <]<т}. (1.1)

Решение задачи (1.1) сводится к построению оценки х* е[а,й], отвечающей некоторому понятию близости к точному решению х (например, чтобы |х - х*| < е, где 8 > 0 есть заданная точность) на основе некоторого числа к значений функционалов задачи, вычисленных в точ-

ках области поиска.

Задача (1.1) может быть рассмотрена в постановке, когда целевая функция ф(х) (в дальнейшем обозначаемая также gm+1) и левые части ограничений gj(x), 1 <1 < т, определены и вычислимы лишь в соответствующих подобластях с [а, Ь], где

б: = [а, Ь], = {х е ф) < 0},

1 < 1 < т. (1.2)

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

С учетом условий (1.2) исходная задача (1.1) может быть представлена в виде

ф(х*) = шт^т+1(х): х е бт+\}. (1.3) Для целей дальнейшего изложения введем классификацию точек х из области поиска [а,Ь] с помощью индекса V = у(х). Указанный индекс

V определяется условиями

g](x) < 0, 1 < 1 < V - 1, gv(x) > 0, (1.4) где последнее неравенство несущественно, если

V = т + 1, и удовлетворяет неравенствам

1 < V = v(x) < т + 1. Данная классификация порождает функцию /(х) = gv(x), V = v(x), (1.5)

определенную и вычислимую всюду на [а, Ь]. Ее значение в точке х есть либо значение левой части ограничения, нарушенного в этой точке (случай, когда V < т), либо значение минимизируемой функции (случай, когда V = т + 1). Поэтому определение значения /(х) сводится к последовательному вычислению величин gj(x), 1 < < 1 < V = v(x), т.е. последующее значение g/+l(x) вычисляется лишь в том случае, когда gj(x) < 0. Процесс вычислений завершается либо в результате установления неравенства gj(x) > 0, либо в результате достижения значения v(x) = = т + 1. Описанная процедура, называемая испытанием в точке х, автоматически приводит к определению индекса V этой точки. Тройка значений

^ = gv(х) г' = gV(x), v = v( x), (1.6) порожденная испытанием в точке х е [а, Ь], на-зываетсярезультатом испытания.

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

х е Qv, х й ^+1, то вспомогательная задача, всегда имеющая решение, может быть записана в виде

gм = gм(хМ) = ш1пС?м(х): х е бмX (1.7) где М есть максимально возможное значение индекса, т.е.

1 <М = шах^(х): х е [а, Ь]} < т + 1. (1.8) Поскольку множество бм всегда непусто,

задача (1.7) всегда имеет решение. При М =

Л * *

= т + 1 решение х = хт+1 является также решением исходной задачи (1.1). При М < т + 1 необходимо выполняющееся неравенство

* г\

gм < 0 может использоваться как индикатор несовместимости ограничений.

Воспользуемся основной идеей индексного подхода, которая состоит в редукции условной задачи (1.7) к безусловной задаче

у(х ) = шт{у(х): х е [а, Ь]},

где

Гgv (х), V < М, х) = ] v * ^ (1.9)

I gм - gм , v= М. Как результат, функция у(х) будет дифференцируемой на каждом множестве 1 < V <

< М. Однако эта новая функция будет иметь разрывы первого рода на граничных точках множеств бу из (1.2), поэтому применять для ее минимизации известные алгоритмы [6-8] непосредственно нельзя. Требуется дополнить алгоритм новыми правилами учета разрывов функции у(х), к чему мы и переходим.

2. Построение миноранты

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

< 1 < т + 1, а именно

(2.1)

(х1) - g} (х2)1< К1 1х1 - х2 I,

1 < 1 < т +1, [х1,х2] с [а,Ь].

Другое предположение, которое является основой рассматриваемых в работе алгоритмов, - предположение о липшицевости первых про-

изводных функций задачи, т.е. выполнение ус-

ловия

Я (х1) - §', (х2 )1 < — |х1 - х2 I

1 < ] < т +1, [х1,х2] с [а,Ь]. (2.2) Основываясь на данном предположении, можно вывести эффективные решающие правила алгоритма. В самом деле, условие (2.2) с учетом разложения функции в ряд Тейлора до второго порядка приводит к неравенству

Я](х) >

> тах <

\ Я] (Х1) + Я' (Х1)( х - Х1) - 0.5-- (х - Х1)2 Я] (*2 ) - я', (*2 )( х2 - х) - 05— (Х2 - X)

где 1 <] < т + 1, а х е [х1, х2]. Данное неравенство позволяет получить нижнюю оценку значений минимизируемой функции на отрезке [ х1, х2 ] , которая достигается в точке

х = -

- (Я, (х2) - Я] (х1)) + (я) (х2)Х2 - я) (х1 )х1 ) -, (х2 - х1) + (Я(х2) - Я' , (х1 ))

,2 0.5-,(х22 - х2)

— (х2 - х1) + (Я' (х2) - Я' (х1)) и может быть записана в одной из двух форм: Я](х) > Я](х1) + я) (х1)( х - х1) -

- 0.5-, (х - х1)2

(2.3)

или

Я] (х) > Я] (х2) - Я'] (х2 )(х2 - х) -

- 0.5— (х2 - х)2. (2.4)

Любая из данных оценок (например, (2.3)) может быть использована, если ] = v(x1) = v(x2), т. е. когда на границах отрезка нарушено одно и то же ограничение, или же точки х1, х2 являются допустимыми. В случае если v(х1) ^ v(х2), то в качестве нижней оценки можно использовать либо неравенство

Я] (х) > Я] (х1 ) + Я] (х1 )(х2 - х1) - 0 5-] (х2 - Х^

] =v(х1) > м(х2 ) (т.е. когда у функции вычислено значение и производная на левой границе отрезка), либо неравенство

Я] (х) > Я] (х2) - Я'] (х2)(х2 - х1) - 0 5(х2 - X1)2,

] =v(х2) >V(х1) (т. е. когда у функции вычислены значение и производная на правой границе отрезка).

Принимая во внимание правило построения оптимизируемой «индексной» функции (1.9), в данных формулах из значения функции Я](х)

нужно вычесть яМ при ] = М. Заметим, ни мак*

симальный индекс М, ни значение ЯМ , ни значе-

ния констант Липшица —м, 1 < V < М, не являются известными заранее. Однако проблему можно преодолеть, используя вместо этих величин их адаптивные оценки, получаемые в процессе решения задачи на основании результатов испытаний. Эти нижние оценки могут быть использованы в алгоритме в качестве характеристик интервалов поиска: очередное испытание нужно проводить в интервале, которому соответствует наименьшее значение миноранты, т. е. в интервале с минимальной характеристикой.

3. Описание алгоритма

Первое испытание осуществляется в произвольной внутренней точке х1 е (а, Ь). Выбор точки х + , к > 1, любого последующего испытания определяется следующими правилами.

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

а = х0 < х1 < ... < хг < ... < хк< хк+1 = Ь, (3.1) и сопоставить им значения 2г = ям (х,) и А = ЯV (х), V = м(х,), 1 < г < к, из (1.6), вычисленные в этих точках; точки х0 = а и хк+1 = Ь введены дополнительно (значения 20 , 2к+1, 20, 2'к+1 не определены) для удобства последующих обозначений.

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

1м = {г: 1 < г < к, м = м(хг)}, 1 < м< т+1, (3.2) содержащих номера всех точек хг, 1 < г < к, имеющих индексы, равные одному и тому же значению V. Граничные точки х0 = а и хк+1 = Ь интерпретируются как имеющие нулевые индексы, и им сопоставляется дополнительное множество 10 = {0, к+1}.

Определить максимальное значение индекса М = тах{м = v(xг), 1 < г < к}. (3.3) Правило 3. Вычислить текущие нижние границы для неизвестных констант Липшица — первой производной функций ям, 1< V < т + 1,

цм=тах

тах

2' - аа

х - х,

2[-( Аг - 2,)+2' ( х,. - х1)]

(х - х}- )2 2[( 2- 2;)-А'(х - х])]

(х - х] )2

г, ] е м,

г > ]

(3.4)

2

+

Если множество Д, содержит менее двух элементов или если из (3.4) оказывается равным нулю, то принять = 1. Из (3.4) следует, что оценки являются неубывающими, начиная с момента, когда (3.4) порождает первое положительное значение

Правило 4. Для всех непустых множеств Д,, 1 < V < т + 1, вычислить оценки

. Г0, v< М,

(3.5)

(хг ): ' е Дv}, V= М.

Правило 5. Для каждого интервала (хг-1, хг), 1 < г < к + 1, вычислить характеристику Я(г), при этом если V = v(хг-1) = v(хг), то

я(г) = (ъ-1-4 ) + г\-1(хг-хг-1) - 0.5г^ (хг-хг-1)2,

где

Xi =

xt xt-1 — s,

g2(x) = x sin (2nx - 0.5).

В предположении частичной вычислимости дуги функций задачи, соответствующие областям Qj, 1 < j < 2, из (1.2), представлены на рис. 1. Для проверки правильности работы метода было оценено (перебором по равномерной сетке) решение задачи x = 2.0795, ф(х ) = 0.565.

Ф (x)

- (zi -zi-1) + (z'ixi -z'-1 xi-1) + О-^ГуЦу (xf -хг2_!)

rv^v (x - x -i) +(z'-i) '

если v = v(xi) > v(xi-1), то

R(i) _(z- z*) - A(x- xi-i) - 0.5fv^v(x- xi-i)2, если v = v(xi-1) > v(xi), то

R(i)_(z-- zv)+z'-1(x- xi-1) - °.5rv^v(x- x-)2.

Величины rv> 1, 1 < v < m + 1, являются параметрами алгоритма. Подходящий выбор значений rv позволяет использовать произведение

rv|iv как оценку констант Липшица Lv, 1 < v < < m + 1, для первых производных.

Правило 7. Определить интервал (xt-1, xt), которому соответствует минимальная характеристика

R(t) _ min{R(i): 1 < i < k + 1}. (3.6) Правило 8. Провести очередное испытание в серединной точке интервала (xt-1, xt), если индексы его концевых точек не совпадают, т. е.

xk+1 = Х+Х^, v( xt-1) *v( xt),

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

x _ xt _

_ - (zt - zt-1) + (ztxt - Z'-1xt-1) + °.5rvlv(xf - xf-1)

rvlv (xt - xt-1) + (z't - z't-1) '

v_v(xt-1) _v(xt). (3.7) Описанные правила можно дополнить условием остановки, прекращающим испытания, если

(3.8)

где t из (3.6) и s > 0 есть заданная точность.

В качестве иллюстрации рассмотрим задачу (1.1) при x £ [0.6, 2.2], m _ 2:

9(x) _ cos (18x - 3)sin (10x - 7) + 1.5, g1(x) _ exp (-x/2)sin (6x - 1.5),

v =1 v=2 v=3

Рис. 1

Данная задача была решена индексным методом без использования и с использованием значений производных; при этом параметры метода выбирались как гу = 2, 1 < V < 3, и е = = 10-5. Оба метода нашли глобально-оптимальное решение. Координаты 63 точек испытаний, осуществленных алгоритмом без учета производных, отмечены на рис. 1 тремя рядами вертикальных штрихов. Штрихи верхнего ряда соответствуют точкам с единичным индексом, второго - точкам, индексы которых равны 2; точки, отмеченные штрихами нижнего ряда, являются допустимыми. На рис. 2 аналогично отмечены координаты 35 точек испытаний, выполненных алгоритмом с использованием производных.

ф(х)

v=1 v=2 v=3

Рис. 2

b

4. Операционные характеристики и сравнение алгоритмов

Один из известных подходов к оценке эффективности методов безусловной глобальной оптимизации основан на численном решении этими методами всех задач из некоторой случайно генерируемой выборки большого объема. В работе [4] описан генератор, порождающий задачи с т ограничениями вида (1.1), на основе заданного случайного механизма Е, дающего функции/(х), х е [а, Ъ]. Предложено два таких генератора: GSH на основе функций Шекеля [16]

j=i

1 < Kj < 3, 0 < Aj, Cj < 10,

и GHL на основе функций Хилла [17]

14

Íhl (x) = X \Aj sin (2 jnx) + Bj cos (2jnx)), x e [0,1],

j=1

-1 < Aj, Bj < 1,

Й-.,

[ (x - Aj )2 + Cj ], x e[0,10],

(4.1)

(4.2)

параметры А], Б], С], К] являются независимыми случайными величинами, равномерно распределенными в указанных выше интервалах.

Примененная процедура оценки основана на операционных характеристиках из [15] и состоит в следующем. Пусть некоторая задача из рассматриваемой выборки решается с помощью алгоритма ^ При этом задаче сопоставляется порождаемая алгоритмом последовательность точек испытаний {х }. Указанная последовательность усекается (т.е. процесс решения прекращается) либо в связи с первым попаданием точки очередного испытания в заданную 8-окрестность решения х , либо в связи с тем, что в ходе выполнения заданного числа Ктах испытаний такое попадание не имело места. В проведенных численных экспериментах использовалось значение Ктах = 1000.

Результат решения всех задач выборки с помощью алгоритма S представляется функцией PS(к), характеризующей долю задач, для которых в ходе к шагов поиска имели место попадания точек испытаний в заданную 8-окрестность решения. Такую функцию будем называть операционной характеристикой алгоритма ^

Для решения серии задач применялись индексный метод без использования производных из [1] (обозначим его ИМ) и описанный выше индексный метод с использованием значений производных (обозначаемый как ИМП). При этом использовались параметры г,, 1 < V < 5, определяемые выражениями

10, ку< 20, [2, ку > 20,

где к^, есть количество вычисленных значений функционала gV (т.е. мощность множества Д,). Точность попадания в окрестность решения задавалась как 8 = 10-4(Ъ - а), где [а, Ъ] - область поиска.

Операционные характеристики для ИМ и ИМП, полученные на классах задач GSH и GHL, представлены соответственно на рис. 3 и 4. При этом нижняя разрывная кривая характеризует ИМ, а верхняя непрерывная - ИМП.

1.0

0.8

0.6 0.4 0.2 0

X

/ / /

1 / f

/ _ S

50 100

150 200 Рис. 3

250 300 350

✓ У

/ / / *

/ /

i J /

J / / S __/

1.0 0.8

0.6 0.4 0.2 0

50 100 150 200 250 300 350 Рис. 4

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

Дополнительную информацию о полученном эффекте ускорения могут дать средние числа кг вычислений значений г-й функций задачи для классов GSH и GHL, представленные в табл. 1 и 2.

Таблица 1

Таблица 2

ИМ ИМП

k1 179 114

k2 107 68

кз 68 43

к4 46 28

к5 29 17

кб 18 11

ИМ ИМП

кх 180 122

кг 106 81

кз 62 54

к4 41 37

к5 24 23

к6 15 14

Заключение

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

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

Работа выполнена при поддержке РФФИ (грант М 11-01-00682-а), совета по грантам Президента Российской Федерации (грант М НШ-1960.2012.9), а также ФЦП «Научные и научно-педагогические кадры инновационной России», соглашениеМ 14.В37.21.0393.

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

1. Стронгин Р.Г., Маркин Д.Л. Минимизация многоэкстремальных функций при невыпуклых ограничениях // Кибернетика. 1986. №4. С. 63-69.

2. Стронгин Р.Г. Поиск глобального оптимума. М.: Знание, 1990.

3. Стронгин Р.Г., Баркалов К.А. О сходимости индексного алгоритма в задачах условной оптимизации с е-резервированными решениями// Математические вопросы кибернетики. М.: Наука, 1999. С. 273-288.

4. Баркалов К. А., Стронгин Р. Г. Метод глобальной оптимизации с адаптивным порядком проверки ограничений // Ж. вычисл. матем. и матем. физ. 2002. Т. 42, №9. С. 1338-1350.

5. Евтушенко Ю.Г. Методы решения экстремальных задач и их применение в системах оптимизации. М.: Наука, 1982.

6. Pinter J. Global optimization in action (Continuous and Lipschitz Optimization: Algorithms, Implementations and Applications). Kluwer Acad. Publ., 1996.

7. Гергель В.П. Об одном способе учета значений производных при минимизации многоэкстремальных функций // Ж. вычисл. матем. и матем. физ. 1996. Т. 36, №6. С. 51-67.

8. Gergel V.P. A Global Optimization Algorithm for Multivariate Functions with Lipschitzian First Derivatives // Journal of Global Optimization. 1997. V. 10. P. 257-281.

9. Sergeyev Ya.D. A global optimization algorithm using derivatives and local tuning. ISI-CNR, Report 1, 1994.

10. Gergel V.P., Sergeyev Ya.D. Sequential and parallel algorithms for global minimizing functions with lipschitzian derivatives// Computers & Mathematics with Applications. 1999. V. 37 №4-5. P. 163-179.

11. Strongin R.G., Sergeyev Ya.D. Global optimization with non-convex constraints. Sequential and parallel algorithms. Dordrecht: Kluwer Academic Publishers, 2000.

12. Стронгин Р.Г. Параллельная многоэкстремальная оптимизация с использованием множества разверток // Ж. вычисл. матем. и матем. физ. 1991. Т. 31, №8. С. 1173-1185.

13. Gergel V.P., Strongin R.G. Parallel computing for globally optimal decision making on cluster systems // Future Generation Computer Systems. 2005. V. 21. №5. P. 673-678.

14. Стронгин Р.Г., Гергель В.П., Баркалов К.А. Параллельные методы решения задач глобальной оптимизации// Известия высших учебных заведений. Приборостроение. 2009. Т. 52. №10. С. 25-32.

15. Гришагин В. А. Операционные характеристики некоторых алгоритмов глобального поиска// Проблемы случайного поиска. Рига: Зинатне, 1978. №7. С. 198-206.

16. Shekel J. Test functions for multimodal search technique// Proc. 5-th Princeton Conf. Inform. Sci. Systems. Princeton Univ. Press, 1971. P. 354-359.

17. Hill J.D. A search technique for multimodal surfaces// IEEE Trans. Systems Science and Cybernetics. 1969. V. 5. №1. P. 2-8.

AN ALGORITHM FOR SOLVING CONSTRAINED GLOBAL OPTIMIZATION PROBLEMS

USING DERIVATIVES

K.A. Barkalov

An algorithm is proposed to solve multiextremal problems with nonconvex constraints. The algorithm uses the values of the first derivatives of the problem which are assumed to satisfy the Lipschitz condition. An index scheme is used to reduce the constrained problem to the unconstrained one. The scheme is based on a separate account of each constraint, the penalty function method not being used. The algorithm is substantiated and numerical experiments have confirmed its convergence and efficiency as compared with those ones which do not take into account the values of the first derivatives of the problem.

Keywords: multiextremal optimization, nonconvex constraints, characteristical algorithms.