ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ
УДК 519.853.4+519.853.6
Городецкий Станислав Юрьевич
канд. физ.-мат. наук, доцент ФГАОУ ВО «ННГУ им. НИ. Лобачевского», г. Нижний Новгород, РФ Е-mail: gorosyu@gmail.com
О НОВОЙ МОДЕЛИ ПОВЕДЕНИЯ ЦЕЛЕВОЙ ФУНКЦИИ ДЛЯ ДИАГОНАЛЬНОЙ РЕАЛИЗАЦИИ
DIRECT-ПОДОБНЫХ МЕТОДОВ
Аннотация
Рассмотрен актуальный вопрос о построении нового диагонального обобщения метода DIRECT для определения глобальных минимумов липшицевых функций без оценивания констант Липшица. Получено непротиворечивое обоснование существующих диагональных реализаций, предложен альтернативный двухточечный метод, основанный на новых предположениях о поведении целевой функции. Новая модель функции может быть использована для построения прямых диагональных обобщений DIRECT на задачи с ограничениями-неравенствами. Проведенные в R1 численные эксперименты подтверждают эффективность предложенного подхода.
Ключевые слова
Глобальная оптимизация, липшицевы функции, метод DIRECT, множество констант Липшица, диагональная двухточечная схема, операционные характеристики.
Gorodetsky Stanislav Yurievich
Ph.D., associate professor Lobachevsky State University (UNN), Nizhny Novgorod, Russia
ABOUT THE NEW MODEL OF OBJECTIVE FUNCTION BEHAVIOR FOR DIAGONAL IMPLEMENTATION OF DIRECT-LIKE METHODS
Annotation
Relevant question about the construction of new diagonal generalization of DIRECT method for finding global minimums of Lipschitz functions without estimating of Lipschitz constants is reviewed. Noncontradictory justification of the existing diagonal implementations is suggested; alternative two points technique is proposed based on new assumptions about objective function behavior. New function model can be used for construction of immediate diagonal generalizations of DIRECT method to tasks with inequality constraints. Numerical experiments conducted in R1 confirm efficiency of the proposed approach.
Keywords
Global optimization, Lipschitz functions, DIRECT method, a set of Lipschitz constants, diagonal two points
scheme, operational characteristics.
Обзор состояния предметной области. В работе исследуются вопросы, связанные со специальной группой компонентных методов поиска глобальных минимумов х* целевой функции f(x) на гиперинтервале D:
f*=f(x*) = minf(x), (1)
XED
D = [a,b] ={ie Rn: Vi = 1, ...,N: at < xt < bt}. (2)
НАУЧНОЕ ПЕРИОДИЧЕСКОЕ ИЗДАНИЕ «CETERIS PARIBUS» №12/2016 ISSN 2411-717Х Целевая функция f(x) в (1) считается липшицевой по отношению к норме ||. || на D с неизвестным значением константы Липшица L, т.е. BL > 0, что
Vx',x" ED: \f(x')-f(x")\ <Щх' -х"Ц.
Выполнение условия Липшица является естественным для многих прикладных задач, а также позволяет строить конечные оценки значений целевой функции по конечному числу ее измерений, что является основой построения многих методов глобальной оптимизации. Один из первых методов липшицевой глобальной оптимизации был предложен Пиявским С.А. в работе 1967 года (более позднюю публикацию см. в [1]). Метод был основан на правиле размещения измерений в минимуме текущей поточечной миноранты липшицевой целевой функции. Вопросам разработки и реализации подобных методов посвящена обширная литература, весьма полно представленная (за период до 2008 года) в библиографии к обзорным разделам монографии [2]. Среди множества работ особо укажем на [3]. В ней предложен одноточечный центральный компонентный метод, основанный на последовательном разбиении исходного гиперинтервала D = [а, Ь] с RN на меньшие компоненты-гиперинтервалы Di = [щ, Ь{\ путем половинного деления по большему ребру того из них, который признавался наиболее приоритетным. При этом в каждой компоненте D i проводилось только одно измерение целевой функции f(c{) в ее центральной точке Ci = (щ + bi)/2, а в качестве ее приоритета R( = R(Di,L) использовалась нижняя оценка функции f на Di в целом, построенная по единственному измерению f(c{) в центральной точке:
R(Di,L)=f(ci)-Lllbi-aiH/2. (3)
Заметим, что в последующем авторами [3] вместо деления на два использовалась более экономичная схема деления на три [4], которая позднее стала применяться в работах других авторов. Эта схема деления иллюстрируется рис. 1, где точками выделены проведенные измерения целевой функции, серым фоном на левой части рисунка отмечена наиболее приоритетная компонента на текущем шаге.
Рисунок 1 - Пример компонентных разбиений по одноточечной центральной схеме деления на три на двух
последовательных итерациях
Использованный компонентный подход позволил (за счет упрощения вида нижней оценки функции по сравнению с [1]) построить методы с простой вычислительной реализацией. В работах Пинтера Я.Д. (Pinter J.D., см. [5] и др.) вместо указанной выше центральной одноточечной схемы использовались различные двухточечные диагональные компонентные схемы проведения испытаний и разбиения наиболее приоритетного гиперинтервала D( на 2N или на две части с проведением измерений функции в вершинах aj и bj, расположенных на главных диагоналях новых гиперинтервалов. При этом нижняя оценка значений функции вдоль главной диагонали приобретала вид, показанный на рис. 2:
R(Di,L) = (f(at) + f(b{))/2 -LHbt-atll/2, (4)
где значение константы L должно было удовлетворять очевидному условию
L>Li,
Li = \f(bi)-f(ai)\/Hbi-aiH
Рисунок 2 - Нижняя оценка функции вдоль главной диагонали при использовании двухточечной
диагональной схемы проведения измерений
В [5] было показано, что существует коэффициент завышения а > 1 для константы Липшица, при использовании которого в (4), т.е. при замене в (4) L на ab, выполняется оценка:
Ухе Di: f(x)>R(Dl,aL). (6)
При этом в [5] показано, что а > 2 , а в [2] этот результат уточнен, оказалось, что достаточно использовать а > V2. Таким образом [2, 5], при достаточном завышении L нижняя оценка значений функции на диагонали из (4) будет являться нижней оценкой значений функции по гиперинтервалу Di в целом.
В работе [6] было обнаружено (см. также изложение этого вопроса в [2]), что диагональные компонентные методы Пинтера Я.Д. обладают избыточностью в проведении испытаний в том смысле, что могут требовать повторного вычисления функции в точках, где ее значение ранее уже было определено. Кроме того, в [6, 2] описана специальная техника реализации схемы деления на три, при которой (за счет введения специального правила для выбора ориентаций используемых главных диагоналей в Di) удалось добиться размещения испытаний в таких вершинах Di, что стал возможным порядок обхода всего множества точек испытаний, порождающий непрерывную ломаную, называемую адаптивной диагональной кривой. Она составлена из главных диагоналей существующих компонент разбиения Di гиперинтервала D, специально ориентированных при проведении разбиений делимых гиперинтервалов (см. рис. 3). Вдоль этой кривой возникает одномерное сечение целевой функции многомерной задачи, что позволяет для выбора наиболее приоритетных компонент текущего разбиения применять правила вычисления приоритетов-характеристик подинтервалов между точками испытаний, разработанные для методов решения одномерных задач. Алгоритмы выполнения разбиений с согласованным выбором ориентаций диагоналей детально описаны в [2]. Для предотвращения повторных измерений значений функции в вершинах делимых гиперинтервалов, в том случае, когда точка такого измерения совпадает с одной из точек ранее выполненных измерений, в алгоритме, описанном в [2] предполагается использовать результат этого ранее выполненного вычисления целевой функции. Для этого в [2] предложена организация базы данных специальной структуры, позволяющая быстро находить результат нужного измерения.
Заметим, что основная часть методов липшицевой оптимизации для задачи (1)-(2) основана на оценивании неизвестного значения константы Липшица L по результатам проводимых измерений функции f, поскольку без знания этой константы невозможно оценить положение глобальных минимумов и минимально возможное значение f* из (1). Правила размещения точек новых измерений хк функции /, как правило, опираются на эти оценки.
Иной подход впервые был предложен в 1993 году в работе [7] авторов: Jones D.R., Perttunen C.D., Stuckman B.E. Полученный в [7] оригинальный метод, получивший название DIRECT (от DIviding RECtangles - разбиение прямоугольных областей, согласно [2, стр.62]), использует правило выбора новых точек измерений липшицевой функции без оценивания ее константы Липшица. При этом данный метод, так же как и методы с оцениванием значения константы, строит существенно неравномерное покрытие области поиска D точками испытаний, более плотное в подобластях с меньшими значениями целевой функции. За счет отсутствия оценок константы L в этом методе невозможно построить правило останова по точности.
\ А \ А
V А К7
/ А / А
Рисунок 3 - Пример построения компонентного разбиения исходного гиперинтервала по безизбыточной двухточечной диагональной схеме деления на три (выделены точки измерений целевой функции)
При решении задачи (1)-(2) DIRECT строит адаптивное компонентное разбиение исходного гиперинтервала D из (2) по указанной выше центральной одноточечной схеме деления на три (см. рис. 1). Измерения функции выполняются в центральных точках гиперинтервалов разбиения и в последующем запоминаются. Основу метода составляют следующие предположения о целевой функции (эти предположения определяют принятую методом DIRECT модель ее поведения). Считается, что константа Липшица L во всех подобластях построенного к -й итерации разбиения [Di}k множества D принимает одинаковые, но неизвестные значения, которые могут быть любыми от 0 до +те, т.е. L E [0, +те).
В DIRECT правило выбора гиперинтервалов для деления на очередном шаге основано на понятии доминирования одних гиперинтервалов разбиения над другими при гипотетическом значении константы Z, а также понятии потенциально оптимальных гиперинтервалов. Эти понятия будут необходимы для дальнейшего изложения.
Кратко поясним терминологию [7, 2]. Будем говорить, что при значении константы Липшица, равной L, гиперинтервал Dj доминирует гиперинтервал Di по отношению к L, если при этом значении Z нижняя оценка R(Dj, L) функции в Dj, вычисленная по измерениям целевой функции, проведенным в Dj, оказывается меньшей, чем нижняя оценка R(Di, 1), вычисленная при тех же условиях для Dj. Гиперинтервал Dj называют недоминируемым по отношению к L, если ни один из гиперинтервалов Di текущего разбиения [D(}k не доминирует Dj при значении Z. Для Dj при данном Z нижняя оценка целевой функции равна минимуму нижних оценок по всем гиперинтервалам текущего разбиения.
Гиперинтервал Dt называют недоминируемым, если существует значение константы L из диапазона от 0 до те, по отношению к которому он будет недоминируемым, т.е. если
3L E [0, те): R(Dt, L) < R(DS, L) V Ds E [Di]k. (7)
Гиперинтервалы для деления на очередной итерации выбираются, согласно принципам DIRECT [2, 7], только среди недоминируемых гиперинтервалов.
Правило выбора кандидатов на деление включает еще одно требование. Обозначим через fk наименьшее найденное значение функции f после проведения к итераций:
ft = rnin[f(x1).....f(xk)},
а также введем последовательность положительных чисел %. Дополнительное требование состоит в том, что нижняя оценка целевой функции делимого интервала должна улучшать найденное значение /£, по крайней мере, на величину Цк. Гиперинтервалы, удовлетворяющие (7) и этим дополнительным требованиям называют потенциально оптимальными.
Приведем формальное определение. Гиперинтервал Dt текущего разбиения [Di}k называют, потенциально оптимальным [2, 7], если для нижней оценки R(Dt,L) функции f на гиперинтервале Dt выполняется условие (7), а также дополнительное требование:
R(Dt,L)<fk:-Vk- (8)
На каждой итерации метода производится деление всех обнаруженных потенциально-оптимальных гиперинтервалов Dt.
НАУЧНОЕ ПЕРИОДИЧЕСКОЕ ИЗДАНИЕ «CETERIS PARIBUS» №12/2016 ISSN 2411-717Х В [7] показано, что можно получить простое правило выделения потенциально оптимальных гиперинтервалов, если представить Di точками с координатами (di/2;F{) на плоскости (d/2;F) при Fi = f(Ci). Если через такую точку провести прямую с тангенсом угла наклона, равным L, то значение R(Di, L) из (3) можно определить как ординату точки пересечения этой прямой с осью F, как показано на рис. 4.
Отсюда следует, что потенциально оптимальными будут только те гиперинтервалы Dt, соответствующие точки ( dt/2; Ft) которых лежат на нижнем правом участке границы выпуклой линейной оболочки следующего множества точек:
Рк = {(ds/2; Fs): Ds 6 [Di}k } U {(0; fk - Vk)}. (9)
На рис. 4 точки, соответствующие потенциально оптимальным гиперинтервалам, выделены черным. Для их нахождения может быть использован алгоритм Грехема, описанный в [8].
Рисунок 4 - Представление гиперинтервалов текущего разбиения точками плоскости (d/2;F) с выделением потенциально оптимальных гиперинтервалов в методе DIRECT
Простота и оригинальность идеи, положенной в основу DIRECT, привлекли внимание многих исследователей. За прошедшее с момента опубликования метода время было построено несколько его модификаций (см., например, библиографию в [2]), а также получены его обобщения для задач глобальной оптимизации с измененными постановками. В [9, 10] для одномерного и многомерного вариантов постановки задачи предложены методы для случая измерения не только функции, но и ее градиента, удовлетворяющего условию Липшица, значение константы в котором неизвестно и может изменяться на промежутке [0, ю). Построено несколько обобщений с использованием центральной схемы измерений для задач с липшицевыми ограничениями-неравенствами [11, 12]. Кроме того, предложены обобщения метода DIRECT, связанные с использованием не одного, а сразу двух измерений в каждом многомерном гиперинтервале Di. Эти измерения проводятся на концах одной из его главных диагоналей, входящих в адаптивную диагональную кривую, непрерывно связывающую, согласно правилам [2], точки проведенных испытаний. Приведенный обзор показывает, что исследования по DIRECT-подобным методам сохраняют свою актуальность.
Анализ существующего диагонального обобщения метода DIRECT. Кратко рассмотрим диагональное двухточечное обобщение метода DIRECT с безызбыточной стратегией, подробно описанное в [2]. Этот метод выполняет компонентное разбиение множества D на гиперинтервалы по описанной выше безызбыточной двухточечной диагональной схеме деления на три. В каждом новом гиперинтервале Di, возникающем за счет деления на три одного из гиперинтрвалов предыдущего разбиения, используется ровно два измерения целевой функции, размещаемых на концах одной из главных диагоналей, имеющей нужную ориентацию. Далее, вместо нижней оценки функции на Di, вида (3) используется соотношение (4), которое отличается от (3) только тем, что вместо значения функции в центральной точке Fi = f(ci) берется среднее арифметическое значений функции на концах этой диагонали:
Fi = (f(al)+f(bl))/2. (10)
НАУЧНОЕ ПЕРИОДИЧЕСКОЕ ИЗДАНИЕ «CETERIS PARIBUS» №12/2016 ISSN 2411-717Х
Это позволяет применить для двухточечной схемы прежнее геометрическое представление набора гиперинтервалов [D^ точками (d{/2; Fi) на плоскости (d/2; F) с применением для Fi новой формулы (10). При выделении на текущей итерации гиперинтервалов для их последующего деления диагональный метод в [2] формально использует геометрическое правило, полученное в обычном методе DIRECT, а именно, выделяет те гиперинтервалы, которые расположены на нижней правой границе выпуклой оболочки множества точек (9), в котором значения Fs вычисляются согласно (10). Построенный метод алгоритмически корректен и обладает хорошими поисковыми свойствами, однако, как показано ниже, его правило выбора делимых гиперинтервалов нельзя рассматривать как непосредственное обобщение принципов метода DIRECT на двухточечную диагональную схему, если не проведено специальное дополнительное обоснование данного метода.
Утверждение 1. Правило выбора делимых гиперинтервалов, применяемое в рассмотренном в [2] двухточечном диагональном геометрическом обобщении DIRECT (метод НМ-БС в [2]), формально не соответствует принципу выбора потенциально-оптимальных гиперинтервалов метода DIRECT [7], определяемому соотношениями (7), (8).
Доказательство. Действительно, принцип потенциальной оптимальности (7)-(8) требует сравнения при каждом значении L E [0, те) нижних оценок функции f для всех гиперинтервалов текущего разбиения [ D(}k. В то же время, используемые в существующем обобщении метода [2] выражения (4) являются нижними оценками функции на главной диагонали Di только в том случае, когда в (4) выполнено условие L>Li для Li из (5) (см. рис. 2).
Таким образом, если при выбранном конкретном значении L для части гиперинтервалов Di неравенство L > Li выполнено, а для другой части - нет, то на плоскости сравнения (d/2;F) при данном L не все гиперинтервалы могут быть представлены точками (di/2; F^), а только те из них, для которых при данном L верно условие L > L^ для Li из (5). Следовательно, в двухточечном диагональном методе множество точек на плоскости ( d/2; F), аналогичное представленному в (9) множеству Pk, будет зависеть от значения L, т.е. преобразуется в множество Pk(L), вид которого приведен ниже в (12). А именно, сравнивать между собой два гиперинтервала Di и Dj по значениям нижних оценок вида (4) при некотором значении L можно лишь в том случае, когда для обоих гиперинтервалов это значение не меньше пороговых значений Li и Lj. Поскольку это не учитывается в геометрическом подходе при построении метода НМ-БС в [2], мы имеем формальное несоответствие принципу построения DIRECT. Утверждение доказано.
Цели исследования. Требуется построить непротиворечивые обобщения на двухточечную диагональную схему основных понятий, определяющих принципы метода DIRECT, в частности, понятий недоминируемых и потенциально оптимальных гиперинтервалов, а также провести корректное обоснование метода НМ-БС из [2]. Также целью исследования является построение новой модели поведения целевой функции и получение на ее основе альтернативного варианта двухточечного диагонального обобщения метода DIRECT, проведение численной апробации на классах одномерных задач.
Непротиворечивое обобщение правил DIRECT на двухточечную диагональную схему. Построим корректные обобщения для понятий «недоминируемый» и «потенциально оптимальный» гиперинтервал при использовании ранее сделанных предположений о поведении целевой функции.
Определение 1. Гиперинтервал Dt будем называть недоминируемым для диагональной схемы, если существует значение константы L из диапазона от 0 до те, что для всех гиперинтервалов Ds текущего разбиения [D(}k, имеющих в (5) значение Ls < L, выполняются неравенства R(Dt, L) < R(DS, V), т.е.
3L E [0, те): Lt < L и R(Dt, L) = min [R(DS, L): V Ds E [Di}k, где Ls < L}. (11)
Определение 2. Гиперинтервал Dt текущего разбиения [Di}k назовем, потенциально оптимальным для диагональной схемы, если для гиперинтервала Dt выполняется условие недоминирования (11), а также требование (8).
Как уже отмечалось выше, при отображении множества гиперинтервалов Ds разбиения [D(}k на плоскости сравнения ( d/2; F) нужно учитывать, что теперь это множество, в отличие от (9), будет зависеть от значения L, при котором сопоставляются эти гиперинтервалы:
При этом значения Fs в (12) вычисляются согласно (10).
Утверждение 2. При использовании двухточечной диагональной схемы разбиения множество недоминируемых гиперинтервалов, выделенных на основе формальной геометрической трактовки определения (7), может не совпадать с множеством недоминируемых для диагональной схемы гиперинтервалов, построенных по корректному определению (11).
Доказательство. Достаточно построить контр пример. Рассмотрим одномерную задачу поиска глобального минимума на отрезке D = [0, 6]. После первой итерации к = 1 отрезок будет разбит на три интервала Dt = [0, 2], D2 = [2, 4], D3 = [4, 6]. Пусть функция f в точках измерений 0, 2, 4, 6 принимает значения равные 4, 2, 4, 4. Тогда, f = 2, и, согласно (5), минимальные значения константы Липшица для интервалов Dt, D2, D3 будут равны Lt = 1, L2 = 1 и L3 = 0. Поэтому диапазон возможных значений констант Липшица будет разбит на две части: [0, Lt) и [Lt, от).
Применим определение 1. Для первого диапазона [0, Lt) на плоскости сравнения (d/2;F) будет присутствовать только одна точка с координатами (1, 4), соответствующая интервалу D3, поскольку для остальных интервалов не существует нижних оценок функции при значениях L < Lt, т.к. такие значения L невозможны в Dt, D2. Таким образом, для диапазона L из [0, Lt) интервал D3 не доминируется другими (поскольку точек других интервалов на плоскости сравнения нет), поэтому он будет отнесен к недоминируемым. Для второго диапазона с L > Lt на плоскости сравнения будут присутствовать две точки с координатами (1, 4) и (1, 3). Первая соответствует интервалу D3, а вторая - двум интервалам: Dt и D2. Поскольку первые координаты точек совпадают, то недоминируемые интервалы будут соответствовать той из них, у которой вторая координата, отвечающая значению Fs, будет меньше. Таким образом, для диапазона L > Lt недоминируемыми будут интервалы Dt и D2. Следовательно, в приведенном примере определение 1 приведет к выделению всех трех гиперинтервалов, как недоминируемых.
Теперь применим геометрическую трактовку определения недоминируемости из (7). Поскольку при формальной трактовке этого определения в сравнении всегда участвуют все гиперинтервалы, в этом случае на плоскости сравнения сразу будет размещено две точки (1, 4) и (1, 3), соответствующие трем интервалам: D3 и Dt, D2. Поэтому в качестве недоминируемых будут выделены только два последних. Таким образом, результаты получились разными, утверждение доказано.
Заметим, что приведенные выше утверждения 1 и 2 не означают, что построенное в [2] двухточечное диагональное обобщение DIRECT (метод НМ-БС [2]) является некорректным. Констатируется лишь несоответствие геометрического подхода, примененного в [2] к обоснованию метода, принципам построения метода DIRECT. Приведенная ниже теорема 1 обеспечивает корректный способ обоснования метода.
Теорема 1. Использование корректного определения потенциальной оптимальности для диагональной схемы (11), (8) приведет к выделению тех же гиперинтервалов Dt для деления, что и использование геометрической трактовки формально некорректного для диагональной схемы определения потенциальной оптимальности из (7)-(8).
Доказательство. Рассмотрим применение определения 2 в случае, когда значения Fi определяются соотношением (5). Допустим, на -й итерации для некоторого гиперинтервала Dt условия (11) и (8) из определения 2 выполнены. Наименьшее значение константы Липшица, которое соответствует выполнению условия (8) для Dt обозначим через L(Dt) Введем обозначения:
(12)
Равенство (5) запишем в эквивалентном виде:
Вычитая из (13) равенство (14) и учитывая, что по определению для любого выполнено неравенство К < , получим:
0 < 2ГГп - 2% < (Ьфг) - Ь1)\\Ъ1 - м. Следовательно, для гиперинтервалов удовлетворяющих определению 2, всегда выполняются неравенства вида:
Дад > Ьг-
Это означает, что при выделении потенциально оптимальных для диагональной схемы гиперинтервалов Ог, согласно правилам (11), (8), дополнительные условия вида Ь3 < Ь, отличающие (11) от (7), не оказывают влияния на выбор этих гиперинтервалов, т.к. эти условия всегда будут выполнены. Теорема доказана.
Новая модель поведения целевой функции, альтернативный диагональный метод. Поскольку при использовании двухточечной диагональной схемы в каждом гиперинтервале текущего разбиения [0(}к имеются результаты двух измерений целевой функции в концевых точках щ и Ъ одной из главных диагоналей , то минимальное значение Ь константы Липшица в этом гиперинтервале, в общем случае, отлично от нуля и определяется соотношением (5). Таким образом, с каждым связывается свое минимальное значение константы ^ > 0.
Примем новую адаптивную модель поведения целевой функции х), предположив, что на текущей -й итерации f является липшицевой в каждом из гиперинтервалов 6 [0(}к со своим значением константы Липшица, которое представимо в виде:
Ь = Ьь+АЬ, (15)
где приращение АЬ одинаково для всех гиперинтервалов, является неизвестным и может принимать любые значения от 0 до т.е. АЬ 6 [0, от).
Такая модель названа адаптивной потому, что в процессе деления гиперинтервалов предположения о функции согласованно изменяются в связи с появлением новых значений Ь3 в новых гиперинтервалах .
Получим выражение для нижней оценки функции в с учетом (5) и (15), используя ее представление (4), а также новое обозначение Яф^ АЬ):
Иф0АЬ) = ЯфьЬ,+АЬ) = (Пад+ПЪ{))/2 - (Ь1 + АЬ) \\Ъ - а,\\/2 = = (/(ад + ПЪ{))/2 - тЪд-Пад\/\Ъ - щ\\ + АЬ) \\Ъ{ - а,\\/2. Выполняя элементарные преобразования, окончательно имеем:
АЬ): = тт ^(а) №)} - АЬ \\Ъ - а1 \\/2. (16)
Для полученных выше оценок нового вида (16) введем представление гиперинтервалов точками (с1{/2; Р{)) на плоскости сравнения (с1/2; F), определив координаты этих точек следующим образом:
^ = \\Ъ{-а\ Fi=тin{f(аi); f(Ъi)} (17)
Аналогично ранее изложенному (см. рис. 4), используя соответствующее (17) представление на новой плоскости сравнения, значение нижней оценки АЬ), отвечающее приращению константы Липшица АЬ, можно интерпретировать как координату точки пересечения с осью ординат луча, проведенного из точки (й[/2; Fi)) с тангенсом угла наклона, равным АЬ (рис. 5).
Введем понятие потенциальной оптимальности, соответствующее новым предположениям о целевой функции.
Определение 3. Гиперинтервал Ог текущего разбиения назовем, потенциально оптимальным
для диагональной схемы по отношению к адаптивной модели поведения целевой функции, если для гиперинтервала Ог выполняется аналог условия недоминирования (7), имеющий вид:
ЗАЬ 6 [0,от): Яф^АЬ) = тт{Яф3,АЬ): Ч036 ^}к]. (18)
а также аналог требования (8), когда для АЬ из (18) выполняется условие
Яфг,АЬ)<К-Пк. (19)
Правила (18)-(19) подобны правилам (7)-(8), если заменить в последних Ь на АЬ и вместо формулы (10) вычислять Fi, исходя из соотношения (17). В силу обнаруженной аналогии, для выделения на плоскости
НАУЧНОЕ ПЕРИОДИЧЕСКОЕ ИЗДАНИЕ «CETERIS PARIBUS» №12/2016 ISSN 2411-717Х сравнения точек (dt/2; Ft)), соответствующих гиперинтервалам Dt, удовлетворяющим (18)-(19), достаточно определить точки, лежащие на нижней правой границе выпуклой оболочки множества точек , где
Рк = {(ds/2; Fs): Ds Е {Di}k} U {(0; ft - Лк)}, (20)
а значения F в (20) вычисляются согласно (17). Для нахождения гиперинтервалов, удовлетворяющих определению 3, можно, как и раньше, использовать алгоритм Грэхема, описанный в [7] (см. рис. 5).
Рисунок 5 - Представление гиперинтервалов текущего разбиения точками в новой адаптивной модели
поведения функции, выделение гиперинтервалов по правилам определения 3 (отмечены черным)
Заметим, что вычисление координат Fi = Fi согласно (17) приведет к дополнительному нежелательному эффекту, связанному с тем, что значительно чаще (по сравнению с существующим методом из [2], использующим для вычисления координат Fi соотношение (10)) будут возникать ситуации, когда найденной, согласно (18)-( 19), потенциально оптимальной точке (dt/2 ; Ft)) на плоскости сравнения будет соответствовать сразу несколько гиперинтервалов Dt с одинаковыми длинами главных диагоналей и равными минимальными измеренными значениями функции на концах используемой диагонали: min {f(ai); f(bi)}. Значения функции во второй диагональной точке, равные max {f(ai); f(bi)} будут, как правило, отличаться.
Поэтому для окончательного выбора гиперинтервалов для деления на данной итерации (будем обозначать их как Dt*) необходимо ввести следующее дополнительное правило.
Дополнительное правило отбора. На множестве гиперинтервалов Dt, удовлетворяющих условиям (18)-(19), выделяются подгруппы с одинаковыми длинами диагоналей dt (на плоскости сравнения (d/2;F) каждой такой группе соответствует одна точка (dt/2 ; Ft)). В каждой такой группе для деления на данной итерации выделяются только такие гиперинтервалы Dt*, для которых значения max {f(at*); f(bt*)} минимальны в этой группе.
После введения описанных выше правил отбора гиперинтервалов Dt* для деления на текущей итерации альтернативный диагональный двухточечный метод можно кратко описать, используя ссылку на подробное изложение существующего двухточечного метода НМ-БС в [2, стр.260] с безызбыточной диагональной стратегией разбиения, использующего значения Fi, вычисляемые по формуле (10). Альтернативный метод отличается от существующего метода [2] двумя моментами. Во-первых, вместо определения множества делимых на текущей итерации гиперинтервалов с использованием геометрической трактовки правил (7)-(8), нужно применять правила (18)-(19) совместно с описанным выше дополнительным правилом отбора. Во-вторых, в методе из [2] значения цк традиционно определяются согласно правилу:
Лк = гкШ1 (21)
где в [2] выбрано Sk = 0,01. В предложенном альтернативном методе используется другое правило определения %:
Ik = ?kbfk, (22)
в котором значение A fk вычисляется по правилам, описанным в [12] (согласно этим правилам значение A fk, начиная с некоторой итерации, является постоянным), а значение £k может изменяться в процессе решения.
Изменения в правило определения ik внесены для того, чтобы последовательность точек измерений целевой функции, порождаемая методом, стала инвариантной к аддитивной постоянной добавке к значениям целевой функции. Существующий метод из [2], использующий (21), такой инвариантностью не обладает, поскольку сама добавка (21) не инвариантна к операции добавления константы. Альтернативный метод, использующий инвариантное по добавлению константы правило вычисления A fk из [12], полностью нечувствителен к добавлению константы к целевой функции.
Заметим, что из способа построения альтернативного метода непосредственно следует, что для него остается справедливой доказанная в [2, стр. 263] теорема 5.6, которая устанавливает сходимость и характер сходимости существующего метода НМ-БС [2].
Использование новой адаптивной модели поведения функции для построения диагонального варианта прямого обобщения метода DIRECT на задачи с ограничениями - неравенствами.
Укажем на специальную область применения предложенной в данной работе новой адаптивной модели функции. В [12] рассмотрено несколько способов обобщения метода DIRECT на задачи с ограничениями неравенствами при использовании одноточечной центральной схемы деления на три. Один из этих методов, названный в [12] ExDIR, основан на прямом обобщении принципов DIRECT на задачи с ограничениями. В этом методе используются нижние оценки не только целевой функции, но и функции обобщенного ограничения g(x) < 0, имеющие вид:
R(Di,L) = Gi-L\\bi-ai\\/2, (23)
где при использовании центрального измерения Gi = g(ci). При этом построение данного метода принципиально опирается на возможность определения по знаку Gi наличия в Di измерения, проведенного в допустимом множестве, где g(x) < 0.
Поэтому для построения диагонального аналога метода ExDIR нужно, чтобы в диагональной схеме значение Gi в (23) определялось формулой Gi = min{g(ai); g(bi)}, а не соотношением вида Gi = (d(.ai) + g(bi))/2. Это возможно, если при построении диагонального аналога метода ExDIR из [12] использовать новую адаптивную модель поведения для функции обобщенного ограничения.
Вычислительный эксперимент. Ниже представлены результаты экспериментального сопоставительного тестирования двух диагональных двухточечных обобщений метода DIRECT в одномерном случае, когда начальная область поиска D представляет собой отрезок в R1: D = [a, b]. В отличие от версии алгоритма НМ-БС, описанной в [2], при проведении экспериментов в обоих методах не использовался встроенный в описание [2] для НМ-БС механизм балансировки локальной и глобальной стратегий. Это сделано для того, чтобы в сравнении участвовали более простые варианты методов в том виде, как они рассматривались в данной статье. Кроме того, поправка ik, используемая в условиях (8) и (19) выбиралась в обоих методах одинаково, согласно (22), где базовое значение A fk определялось по алгоритмам, представленным в [12]. Значение £k = £ в (22) принималось постоянным, использовалось как параметр методов.
При такой постановке эксперимента два метода отличались только принятой моделью функции. Существующий метод из [2], использующий традиционную модель, выбирает интервалы Dt для деления на текущей итерации, используя в геометрической трактовке правил (7)-(8) нижние оценки функции вида (4). При этом на плоскости сравнения (рис. 4) по оси ординат откладываются значения Fi из (10) в виде среднего арифметического значений функции на концах интервала Di. Поэтому данный метод в результатах эксперимента помечен как метод типа Midi.
Альтернативный метод, использующий новую адаптивную модель функции, выбирает интервалы Dt* для деления на текущей итерации, используя аналогичные условия (18)-(19) с нижними оценками функции вида (16), а также введенное дополнительное правило отбора. При этом на плоскости сравнения (рис. 5) по
НАУЧНОЕ ПЕРИОДИЧЕСКОЕ ИЗДАНИЕ «CETERIS PARIBUS» №12/2016 ISSN 2411-717Х оси ординат откладываются значения Fj из (17) в виде наименьшего из значений функции на концах интервала Dj. Поэтому данный метод в результатах эксперимента будет помечаться как метод типа Min.
Эксперименты проводились на двух классах тестовых функций: Пинтера и Шекеля. Набор тестовых задач Пинтера включал 100 функций, Шекеля - 100 функций.
При тестировании на классе функций Пинтера выполнялось решение 100 задач вида:
min /5(х),
хЕ[-5, 5]
/¿(х) = 0,025(х — х^)2 + sm2((x — х^) + (х — xi!)2) + sjn2(x — х^),
где параметр функции выбирался случайно на отрезке [—5, 5] при использовании датчика с равномерным законом распределения.
При тестировании на классе функций Шекеля проводилось решение 100 задач вида:
min /5(х), *е[0,10]
/я(х) = E1=°1 l/(ci,s + (х — Qj,s)2),
где для каждого s = 1,... ,100 параметры Cj,s, ( i = 1,... ,10) выбирались случайно на промежутке (0, 2], а Q(,s, 0 = 1, .,10) - случайно на отрезке [0, 10] при использовании датчика с равномерным законом распределения.
Во всех случаях использовался стандартный датчик rand() из библиотеки С++, перед генерацией набора из 100 функций каждого класса выполнялась его инициализация вызовом srand(1).
На рис. 6 показан пример размещения точек измерений функции двухточечным методом типа Min для одной из функций Шекеля.
Рисунок 6 - Размещение к = 50 измерений альтернативным двухточечным методом типа Min с параметром
£ = 0,05 для функции Шекеля
На каждом наборе тестовых функций проводилось решение каждым из методов с использованием заданного числа к вычислений целевой функции. Очередная тестовая задача считалась решенной, если для оценки решения х£, полученной методом в этой задаче после к вычислений функции, выполнялось условие: |х£ — х| < (Ь — а)Д , где х* - точка глобального минимума решаемой тестовой задачи с функцией /¿(х), а Д - заданная относительная точность решения.
В результате расчетов были получены операционные характеристики в виде зависимостей р^ (Д) от к, где р^ - доля правильно решенных задач из тестовой выборки в результате выполнения к вычислений целевой функции при решении каждой из задач. Примеры вида операционных характеристик на классе тестовых функций Шекеля для диагональных методов двух типов (Midi и Min) приведены на рис. 7.
Рисунок 7 - Операционные характеристики методов Midi и Min для класса тестовых функций Шекеля при относительной точности А= 0,001 и параметре методов £ = 0,003
При построении операционных характеристик количество вычислений функции к варьировалось от 1 до значения ктах = 100. Основными показателями, которые фиксировались в экспериментах, являлись:
р*(А) = тах{рк(А): к = 1, ...,ктах}
- наибольшее достигнутое значение доли задач, решенных с заданной точностью не более чем за кта х вычислений функции;
к*(А) = min{k: рк(А) = р*(А); к < ктах}
- наименьшее количество вычислений функции, при котором было впервые достигнуто наибольшее значение доли правильно решенных задач.
Заметим, что значения р*(А) и к*(А) кроме заданной точности решения А существенно зависят от значения параметра £. Чтобы исключить это влияние, при каждом А величина £ выбиралась для методов Midi и Min наилучшим образом так, чтобы показатель к*(А) оказался минимально возможным для данного метода. Определяемые таким образом значения £ обозначены ниже в результатах эксперимента как £*(А).
Результаты проведенных расчетов представлены в табл. 1 и табл. 2.
Таблица 1
Сравнительные результаты численных экспериментов на классе из 100 тестовых функций
Шекеля при ктах = 150
А Тип метода - Midi Тип метода - Min
Р*(А) к*(А) £*(А) Р*(А) к*(А) £*(А)
0,01 1,00 38 0,100 1,00 42 0,050
0,001 1,00 92 0,003 1,00 53 0,003
0,0001 0,93 102 0 0,93 77 0
Таблица 2
Сравнительные результаты численных экспериментов на классе из 100 тестовых функций
Пинтера при ктах = 150
А Тип метода - Midi Тип метода - Min
Р*(А) к*(А) £*(А) Р*(А) к*(А) £*(А)
0,01 1,00 26 0,5000 1,00 30 0,02000
0,001 1,00 40 0,0001 1,00 40 0,00001
0,0001 1,00 60 0 1,00 64 0
Еще раз отметим, что оба метода применялись без использования механизма балансировки локальной и глобальной стратегий поиска, который был предложен в [2] для ускорения сходимости.
Сравнение результатов показывает, что ни один из сравниваемых методов не является абсолютно лучшим. В ряде ситуаций большую экономичность показывает метод типа Midi, а в другой части - метод
НАУЧНОЕ ПЕРИОДИЧЕСКОЕ ИЗДАНИЕ «CETERIS PARIBUS» №12/2016 ISSN 2411-717Х типа Min. Таким образом, с учетом отмеченной ранее области использования новой адаптивной модели функции, она может быть рекомендована в качестве основы двухточечных диагональных обобщений метода ExDir [12] в задачах с ограничениями.
Заключение. В данной работе проведен анализ существующего двухточечного диагонального геометрического обобщения метода DIRECT, показана противоречивость применяемого геометрического способа его обоснования, построен корректный способ вывода данного метода. Кроме того, предложена новая адаптивная модель поведения липшицевой целевой функции для DIRECT-подобных методов при наличии различных нижних ограничений на значение константы Липшица в гиперинтервалах разбиения и неограниченном сверху диапазоне ее значений. Построено новое альтернативное диагональное обобщение метода DIRECT, проведено сравнительное численное исследование нового метода с сушествующим, показано, что оба метода сопоставимы по эффективности. Список использованной литературы:
1. Пиявский, С.А. Один алгоритм отыскания абсолютного экстремума функций / С.А. Пиявский // Ж. вычисл. матем. и матем. физ. — 1972. — Т.12, № 4. — С. 888-896.
2. Сергеев, Я.Д. Диагональные методы глобальной оптимизации /Я.Д. Сергеев, Д.Е. Квасов // М.: Физматлит, 2008. — 352 с.
3. Евтушенко, Ю.Г. Метод половинных делений для глобальной оптимизации функций многих переменных / Ю.Г. Евтушенко, В.А. Ратькин // Изв. АН СССР. Техническая кибернетика. — 1987. — Т.1. — С. 119-127.
4. Evtushenko, Y.G. Deterministic global optimization / Y.G. Evtushenko, M.A. Potapov // Algorithms for Continuous Optimization: The State of the Art / Ed. by E.Spedicato. — Dordrecht: Kluwer Academic Publishers, — 1994. — NATO ASI Series. — P. 481-500.
5. Pinter, J.D. Global Optimization in Action. (Continuous and Lipschitz Optimization: Algorithms, Implementations and Applications) / J.D. Pinter. — Dordrecht: Kluwer Academic Publishers, 1996.
6. Sergeyev, Y.D. A new optimization algorithm using adaptive diagonal curves and its numerical testing / Y.D. Sergeyev, D.E. Kvasov // Tech. Rep. 1 — Institute of Systems and Informatics, Rende(CS), Italy: ISI-CNR, 2002.
7. Jones, D.R. Lipschitzian optimization without the Lipschitz constant / D.R. Jones, C.D. Perttunen, B.E. Stuckman // J. Optim. Theory and Appl. — 1993. — V. 79, № 1. — P. 157-181.
8. Препарата, Ф.Ф. Вычислительная геометрия. Введение / Ф.Ф. Препарата, М.И. Шеймос // Монография.
— М.: Мир, 1989. — 478 с.
9. Sergeyev, Ya.D. A univariate global search working with a set of Lipschitz constants for the first derivative / Ya.D. Sergeyev, D.E. Kvasov // Optimization Lettes. — 2009. — № 3. — P. 303-318.
10.Kvasov, D.E. Lipschitz gradients for global optimization in a one-point-based partitioning scheme / D.E. Kvasov, Ya.D. Sergeyev // Journal of Computational and Applied Mathematics. — 2012. — V. 236. — P. 4042-4054.
11.Городецкий, С.Ю. Обобщения метода Direct на задачи с функциональными ограничениями / С.Ю. Городецкий // Высокопроизводительные параллельные вычисления на кластерных системах: Материалы XII-ой Всероссийской конференции (26-28 ноября 2012, Н.Новгород). — Н.Новгород: Изд-во ННГУ, 2012. — С.101-105.
12.Городецкий, С.Ю. Несколько подходов к обобщению метода DIRECT на задачи с функциональными ограничениями / С.Ю. Городецкий // Вестник Нижегородского университета им. Н.И. Лобачевского. — 2013,
— № 6(1). — С. 189-215.
© Городецкий С.Ю., 2016