CC BY
49
УДК 519.854.33
И. С. Масич
ПРИБЛИЖЕННЫЕ АЛГОРИТМЫ ПОИСКА ГРАНИЧНЫХ ТОЧЕК ДЛЯ ЗАДАЧИ УСЛОВНОЙ ПСЕВДОБУЛЕВОЙ ОПТИМИЗАЦИИ
Исследуется применение схемы поиска граничных точек для решения задачи оптимизации псевдобулевых функций с ограничениями на переменные. Рассматривается несколько приближенных алгоритмов, основанных на этой схеме. Сравнение эффективности алгоритмов между собой и с универсальными поисковыми алгоритмами проводится на соответствующей практической задаче большой размерности.
Методы оптимизации можно разделить на алгебраические методы, требующие задания целевой функции и ограничений на переменные в явном виде, и поисковые методы, для применения которых целевая функция и ограничения могут быть заданы алгоритмически, а также для случаев, когда функции заданы аналитическими выражениями сложного вида.
Поисковые методы не используют функций явного вида, но при этом они могут учитывать различные свойства функций. Эти свойства могут определяться исходя из специфики решаемой задачи или с помощью некоторой процедуры идентификации свойств, которая на основе пробной выборки проверяет гипотезу о наличии определенного свойства.
Универсальные поисковые алгоритмы, такие как генетический алгоритм (ГА), эволюционный алгоритм, алгоритм имитации отжига, могут применяться для решения любых задач оптимизации. Но для некоторых практических задач они не находят приемлемого решения. Это происходит при наличии сложной системы ограничений на переменные и, как следствие, малой допустимой области. Для таких случаев необходима разработка алгоритмов, учитывающих специфику решаемой задачи.
В данной статье рассматриваются задачи условной псевдобулевой оптимизации, т. е. задачи оптимизации вещественных функций булевых переменных с ограничениями на переменные.
Приведем определения, необходимые для описания работы алгоритмов оптимизации [1]:
- псевдобулевой функцией будем называть вещественную функцию на множестве булевых переменных: / : В2П ^ Я1, где В2 = {0,1}; В2п = В2 хВ2 х...хВ2;
- точки X1, X2 е В2И назовем к-соседними, если они отличаются значением к координат, где к = 1, п ; /-соседние точки будем называть просто соседними;
- множество Ок (X), к = 0, п, всех точек ВП, к-сосед-них к точке X, назовем к-м уровнем точки X;
- множество точек Ш(X0,X1) = {X0, X\ ...,X1} сВ2п
- это путь между точками X0 и X1, если V ' = 1, к, I точка X' является соседней к X'-1;
- множество А с В2П назовем связным множеством, если VX0, X1 е А существует путь Ш (X0, X1) с А ;
точка X* е ВП, для которой
/(X ) < /(X), VX е О1^ ), является локальным минимумом псевдобулевой функции/;
- псевдобулевую функцию, имеющую только один локальный минимум, будем называть унимодальной на множестве В2п функцией;
- унимодальную функцию/назовем монотонной на Bn2, если VXk е Ok (X*), k = 1, n :
f (Xk-1) < f (Xk), VXk-1 е Ok-i (X*) n Oi (Xk).
Пример. Полином от булевых переменных
m _____
фСх^..^хп) = ХajПX-,
j=1 ,е|3 j
где вj с {1,..., п}, является унимодальной и монотонной псевдобулевой функцией при aj > 0, j = 1, m .
Постановка задачи. Рассмотрим задачу следующего вида:
C(X) ^ max , (1)
X eS с в2
где C(X) - монотонно возрастающая от X0 псевдобу-левая функция; S с ВП - некоторая подобласть пространства булевых переменных, определяемая заданной системой ограничений, например:
Aj (X) < Hj, j = im . (2)
В общем случае множество допустимых решений S является несвязным множеством.
Множество допустимых решений задачи (1) обладает следующими свойствами:
- точка Y е A является граничной точкой множества A, если существует X е O1 (Y), для которой X й A ;
- точку Y е Oj (X0) n A будем называть крайней точкой множества A, с базовой точкой X0 е A, если
VX е O1 (Y) I Oi+1 (X0) выполняется X й A ;
- ограничение, определяющее подобласть пространства булевых переменных, будем называть активным, если оптимальное решение задачи условной оптимизации не совпадает с оптимальным решением соответствующей задачи оптимизации без учета ограничения.
Рассмотрим другие свойства множества допустимых решений [2]:
- если целевая функция является монотонной унимодальной функцией, а ограничение активно, то оптимальным решением задачи (1) будет точка, принадлежащая подмножеству крайних точек множества допустимых решений S с базовой точкой X0, в которой целевая функция принимает наименьшее значение:
C(X0) = min C(X);
X еВ2
- если функция ограничения (2) является унимодальной функцией, то множество допустимых решений S задачи (1) представляет собой связное множество;
- количество крайних точек связного множества допустимых решений задачи (1) s < smax = С[п2], где [п/ 2]
- целая часть числа п / 2. В случае s = smax все крайние точки принадлежат O^2 (X ) при четном n и O(n-1)/,2 (X ) (или O(n+1y2(X0)) при нечетном n.
Переход к безусловной оптимизации. Одним из способов учета ограничений в условных задачах является составление обобщенной штрафной функции. Рассмотрим обобщенную функцию
m
F (X) = C (X) - r Y max{0, Aj (X) - Hj}
(3)
j=1
для задачи (1) с ограничениями (2).
В этом случае имеет место следующее утверждение [3]:
Утверждение. Крайние точки множества допустимых решений S исходной задачи (1) с ограничениями (2) с базовой точкой X являются точками локальных максимумов для обобщенной штрафной функции (3) при параметре г, удовлетворяющем условию
> C (Z) - C (X')
r >---------------------
m
Y max{0, Aj (Z) - Hj}
j=1
для любой крайней точки X' е Ok (X0) и любой точки
Z е Ok+i(X0) n Oi(X').
Таким образом, задача (1) идентична задаче
F(X) ^ max , (4)
XеВП
в которой количество локальных максимумов (в случае связности множества допустимых решений) s < smax = 2]. В общем случае, когда допустимое мно-
жество не является связным, число локальных максимумов теоретически может достигать 2й-1.
Задача (4) решается известными поисковыми методами: локального поиска, генетическим алгоритмом, алгоритмом имитации отжига и т. д.
Существенным недостатком такого подхода является потеря свойства монотонности целевой функции C (X). При введении даже простых (например, линейных) ограничений обобщенная функция будет полимодальной немонотонной функцией с экспоненциальным числом точек локальных максимумов.
Приближенные алгоритмы поиска граничных точек. Для любой эвристики поиска граничных точек будем рассматривать два алгоритма: прямой и двойственный. Прямой алгоритм начинает поиск из допустимой области и движется по пути возрастания целевой функции, пока не найдет крайнюю точку допустимой области. Двойственный же алгоритм ведет поиск в недопустимой области по пути убывания целевой функции, пока не найдет некоторого допустимого решения.
Общая схема прямого алгоритма поиска следующая:
- шаг 1. Полагаем X1 = X0, i = 1;
- шаг 2. В соответствии с правилом выбираем X+1 е Ot (X0) n O1 (X) n S. Если таких точек нет, то переходим на шаг 3; иначе i = i +1 и повторяем шаг 2;
- шаг 3. Xopt = Xt+1.
Общая схема двойственного алгоритма поиска такова:
- шаг 1. Полагаем X1 е On (X0), i = 1;
- шаг 2. В соответствии с правилом выбираем X+1 е On-j (X0) n O1(Xi). Если Xt+1 е S , то переходим на шаг 3; иначе i = i +1 и повторяем шаг 2;
- шаг 3. X opt = Xi+1.
Таким образом, прямой алгоритм находит крайнюю точку допустимой области, в то время как двойственный
алгоритм находит граничную точку, которая может и не являться крайней. Поэтому для задач со связной областью допустимых решений прямой алгоритм в среднем находит лучшее решение, чем двойственный алгоритм. Но если мы применим прямой алгоритм для задачи с несвязной допустимой областью, то полученное в результате решение может быть далеким от оптимального, так как на пути возрастания целевой функции допустимые и недопустимые решения будут чередоваться. Для таких случаев более пригоден двойственный алгоритм, для которого это чередование не играет никакой роли. Для улучшения решения, полученного двойственным алгоритмом, рекомендуется применять соответствующий прямой алгоритм. Такое улучшение на практике оказывается весьма значительным.
Рассмотренные далее алгоритмы поиска граничных точек отличаются друг от друга лишь правилом выбора следующей точки на шаге 2 общих схем, прямого и двойственного алгоритмов поиска.
Будем использовать следующие правила поиска:
- правило 1 - случайный поиск граничных точек (СПГ). Точка Xi+1 выбирается случайным образом. Каждая точка на следующем шаге может быть выбрана с равной вероятностью. При решении практических задач эти вероятности могут быть не равными, а вычисляться на основании специфики задачи до начала поиска;
- правило 2 - Гриди-алгоритм. Точка Xt+1 выбирается по условию
М( X+1) = max ^( X1),
j
где Xj е Of (X0) n O1 (X) n S для прямого алгоритма и XJ е On-i (X0) n O1(Xi) для двойственного алгоритма.
Функция М(X) выбирается исходя из специфики задачи, например: целевая функция X(X) = C(X), удельная ценность X(X) = C(X)/A(X) (для одного ограничения) и т. д.;
- правило 3 - адаптивный случайный поиск граничных точек (АСПГ). Точка Xi+1 выбирается случайным образом в соответствии с вектором вероятностей
Pi = (p, p2pJ ),
где J - количество точек, из которых производится выбор;
i UXj) —
Pj = J , j=1, J >
Y M( x; )
l=1
здесь Xj е Oi (X0) n O1 (X) n S для прямого алгоритма и Xj е Onч (X0) n O1 (X) для двойственного алгоритма. Алгоритм АСПГ является дополнением к Гриди-алгоритму;
- правило 4 - модифицированный случайный поиск граничных точек (МСПГ). Точка Xi+1 выбирается по условию
М( X,+1) = max М( Xr),
r ___
где Xr - точки, выбранные по правилу 1, r = 1,R , здесь R - задаваемый параметр алгоритма.
Гриди-алгоритм является регулярным алгоритмом, т. е. при повторном запуске из одной и той же точки он выдает одинаковое решение. Остальные алгоритмы можно запускать несколько раз и из найденных решений выб-
рать наилучшее. Трудоемкость каждого запуска алгоритма ограничена:
. п(п +1)
T <■
2
но средняя трудоемкость Гриди-алгоритма и АСПГ значительно выше остальных, так как на каждом шаге они просматривают все соседние точки следующего уровня, в отличие от СПГ и МСПГ, которые в двойственной схеме просматривают на каждом шаге лишь одну и R точек соответственно.
Рассмотрим применение описанных алгоритмов на одной из практических задач большой размерности с несвязным множеством допустимых решений.
Задача оптимизации загрузки производственных мощностей литейных отделений. В литейных отделениях (ЛО) производится продукция различного вида. В каждом ЛО имеется специализация по видам продукции, которую могут производить на его литейных машинах (ЛМ). Имеется некоторое количество заказов на производство продукции. Каждому заказу соответствует объем, вид продукции и срок выполнения. Вид продукции характеризуется производительностью за смену (всего три смены в сутки). Замена производимой на ЛМ продукции требует ее перенастройки, занимающей одну смену. Необходимо загрузить производственные мощности ЛО таким образом, чтобы выполнялись все заказы, продукция производилась более равномерно во времени и число перенастроек оборудования ЛМ было минимальным.
Входные данные:
- I - количество дней планирования;
- 3 - количество ЛО;
- Кр - количество ЛМ в j-м ЛО, / = 1,3;
- Ь - количество заказов на продукцию, производимую на ЛМ (и соответствующее число видов продукции);
- У1 - производительность 1-го вида продукции за смену на ЛМ, I = 1, Ь ;
- Тг - срок выполнения 1-го заказа (на производство 1-го вида продукции) на ЛМ, I = 1, Ь ;
- Ж1 - объем 1-го заказа (на производство 1-го вида продукции) на ЛМ, I = 1, Ь ;
г - характеризует специализацию ЛО:
1, если ЛМ _/-го ЛО способны
2 р = - производить продукцию I -го вида,
0 в противном случае;
- а - коэффициент жесткости ограничения по требованию равномерности производства продукции по дням, 0 <а<1.
Для построения модели введем следующие булевы переменные:
у = {Уцы}е В; где В ={о, 1}, вп = в х в хЛх в2
- множество булевых переменных;
X = {Х/т } е вп ;
1, если в /-й день на к -й ЛМ / -го ЛО
- производится продукция I -го вида,
0 в противном случае;
1, если в /-й день на к-й ЛМ / -го ЛО начинается производство продукции I -го вида, 0 в противном случае;
У ijkl
ijkl
Общая размерность булева вектора У (и X)
п = I • Ь^К/.
/=1
Дополнительные замечания:
1) / ^ У/И V/, /, к, I;
2) / = У •(1 - У-1,/И) V/, /, к, I
(У0/м = 0 V/, к, I). ^ (5)
Оптимизационная модель будет иметь следующий вид:
1) целевая функция и основные ограничения:
С(Х) ^ шт, (6)
4(7) > Щ, I = , (7)
А/(У) >а-Г1, /= Ц, ае (0,1), ^ = Х^1, (8)
где С (X) = Х X X Xj = X X X Ху#« I1 - У/-1, jki); i=1 j=1 k=1 l=1 i=1 j=1 k=1 l=1
T J K
4(Y) = XXXV • yjkI(2 + y-1,m);
i=1 j=1 k=1
J Kj L
Ai (Y) = XxX Vl ' yijkl (2 + yi-1, jkl );
j=1 k=1 l=1 Уоjkl = 0 V j, k, l;
2) дополнительные ограничения:
LL
X yjki <1 (X j < 1), i = Ц, j = , k = iTKj, (9)
. j
l=1 l=1
Ут < ^ (xk < ^), i=й,
j = 1, J, k = 1, Kj , l = 1, L .
(10)
Модель обладает следующими свойствами:
- имеется два пространства булевых переменных (обозначим их Вх и В7), соответствующих векторам X и У. Каждой точке У е В7 соответствует единственная точка X е Вх, компоненты которой определяются соотношением (5). При этом точке X е Вх может соответствовать несколько точек У е ВУ (с разным значением функций ограничений);
- целевая функция (6) является линейной и унимодальной монотонной в пространстве Вх с точкой минимума X0 = (0,..., 0). В пространстве ВУ целевая функция является квадратичной и унимодальной немонотонной с точкой минимума У0 = (0,..., 0);
- функции ограничений (7) и (8) в пространстве ВУ являются квадратичными и унимодальными монотонными псевдобулевыми функциями с точкой минимума У0 = (0,..., 0). В пространстве Вх функции ограничений однозначно не определены;
- множество допустимых решений в пространствах
Вх и ВУ ограничено сверху IX Кр -м уровнем точки
/=1
минимума (X0 и У0) согласно ограничению (9). В пространстве ВУ этот уровень соответствует тому, что в любой день на каждой литейной машине производится продукция;
- множество допустимых решений в общем случае является несвязным множеством (в пространстве ВУ).
Таким образом, решение задачи полностью определяется переменными У, чего нельзя сказать про переменные X. Но целевая функция от X обладает хорошими конструктивными свойствами, в связи с чем поиск оптиму-
l=1
ма по X является более эффективным, чем по У. Но так как функции ограничений (7), (8) от X не определены, то следует находить значения этих функций от соответствующей точки У. Таких точек может быть несколько:
X ^ У„ У/,..., Ун , причем в некоторых из них решение может быть допустимым, а в некоторых нет. А так как функции ограничений здесь монотонны, то для конкретной X следует выбирать У, принадлежащую наиболее возможному уровню (с большими значениями функций):
/ \
"V н У1/к1
К и/’кл
где УЬ = (у1111, ..., уиК3Ь )
Ниже приведен один из алгоритмов такого преобразования (алгоритм 1): _____
- шаг 1. Полагаем, что N]к = 0 , / = 1,3, к = 1, Кр,
/ = 1; _________________ _________ _______
- шаг 2. Для / = 1,3, к = 1, Кр, I = 1, Ь выполняем условие: если хш = 1, то N и = I ■_____
У = тах_
У, ,И=1, н
- шаг 3. Для і = 1,3, к = 1, К] , I = 1, Ь выполняем условие: если N^ = І, то у цЫ = 1, иначе уіт = 0;
- шаг 4. Если і < I, то і = і +1 и переходим на шаг 2.
В то же время решение У, полученное из найденного
наилучшего вектора Хор4 таким способом, может соответствовать ситуации, когда количество выпускаемой продукции намного превосходит требуемое значение (что, в общем-то, не противоречит построенной модели, но может влиять на равномерность загрузки производственных мощностей, которая оптимизируется на следующем этапе). Поэтому после завершения первого этапа поиска Уор4 следует определять по правилу
УоР4 = агЕ
тіп
і, і,к,І
к
УіікІ
¥к :4(ГИ )>Щ ,1=1,Ь
В этом случае алгоритм преобразования несколько отличается от предыдущего (алгоритм 2): ____
- шаг 1. Полагаем, что N]к = 0, / = 1,3, к = 1, Кр;
- шаг 1а. Полагаем, что у/]Ы = 0, / = 1,I, / = 1,3,
к = 1, К , / = 1; ’
- шаг 2. Для і = 1,3, к = 1, К] , І = 1, Ь выполняем условие: если х ік1 = 1, то Nіk = І;
- шаг 3. Для і = 1,3, к = 1, Ку , , = 1, Ь выполняем
условие: если Nіk = І и 4 (У) < , то у ік1 = 1;
- шаг 4. Если і < I, то і = і +1 и переходим на шаг 2.
Здесь в шаг 3 добавлено условие Л] (У) < , а также
введен шаг 1а для возможности такого вычисления. На-
добности в таком преобразовании во время поиска нет, оно необходимо только для определения конечного Уор4.
Для решения задачи оптимизации загрузки производственных мощностей использовались двойственные алгоритмы: СПГ, Гриди-алгоритм и МСПГ. Алгоритм АСПГ не рассматривался из-за его чрезмерной трудоемкости при многократном запуске. Единичный же запуск АСПГ очень редко может выдать решение лучше, чем решение Гриди-алгоритмом. Стартовой точкой поиска являлась точка безусловного минимума целевой функции X0 = (0,..., 0). Найденное решение улучшалось соответствующим прямым алгоритмом.
Кроме того, задача решалась генетическим алгоритмом. При реализации ГА была выбрана схема, которая эффективно работала при неоднократном решении других задач комбинаторной оптимизации [4; 5].
Результаты экспериментов показали, что для данной задачи наиболее эффективными алгоритмами по точности и трудоемкости являются Гриди-алгоритм и МСПГ. Остальные алгоритмы при жестких ограничениях на переменные вообще не находят допустимого решения. Это связано со спецификой задачи, в которой имеется большое количество различных ограничений и, как следствие, относительно малая допустимая область.
Осредненные результаты решения 10 задач месячного планирования загрузки производственных мощностей. приведены в таблице. Средние значения входных данных: I = 31; 3 = 3 ; К1 = 12 ; К2 = 9 ; К3 = 7 ; Ь = 36 ; а = 0,5 ; V1 е [40, 50]; ^ е [20, 25 000]. При этом общая размерность булева вектора п = 31 248 .
Количество запусков L алгоритмов выбиралось таким образом, чтобы трудоемкость приблизительно равнялась трудоемкости единичного запуска Гриди-алгоритма. В данном случае трудоемкость Т ~ 8 -10 .
Модифицированный случайный поиск граничных точек по сравнению с Гриди-алгоритмом является более гибкой процедурой, так как позволяет выбирать параметры L и R, влияющие на трудоемкость алгоритма и точность решения. Гриди алгоритм такой возможности не предоставляет, и при больших размерностях время поиска может быть чрезмерно велико. Что касается эффективности, то при примерно одинаковой трудоемкости точность найденных решений различается несущественно.
Таким образом, алгоритмы поиска граничных точек показали высокую эффективность при решении задачи условной псевдобулевой оптимизации с несвязной допустимой областью. Наиболее эффективными для рас-
Алгоритм Количество запусков L Найденное решение СО^
спг 2 000 Не найдено
Гриди 1 49
МСПГ, Л=1000 12 47
МСПГ, Л=100 60 52
МСПГ, Л=10 200 Не найдено
ГА* - Не найдено
* Параметры ГА: турнирная селекция с величиной турнира 5, размер популяции 100, наибольшее число поколений 8 000, вероятность мутации 0,000 1.
смотренной задачи являются двойственные алгоритмы МСПГ и Гриди.
Библиографический список
1. Антамошкин, А. Н. Регулярная оптимизация псев-добулевых функций / А. Н. Антамошкин. Красноярск : Изд-во Краснояр. ун-та, 1989. 284 с.
2. Антамошкин, А. Н. Не улучшаемый алгоритм условной оптимизации монотонных псевдобулевых функций / А. Н. Антамошкин, И. С. Масич // Исследовано в России [Электронный ресурс]. Электрон. журн. 2004.
№ 64. С. 703-708. Режим доступа : http://zhurnal. ape.relarn.ru/articles/2004/064.pdf. Загл. с экрана.
3. Масич, И. С. Решение задач условной псевдобу-левой оптимизации посредством обобщенных штрафных функций / И. С. Масич // Вестн. молодых ученых. Сер. «Прикладная математика и механика». 2004. N° 1 (4). С. 63-69.
4. Goldberg, D. E. Genetic algorithms in search, optimization, and machine learning / D. E. Goldberg // Reading, MA : Addison-Wesley, 1989.
5. Schwefel, H.-P. Evolution and Optimum Seeking /
H.-P. Schwefel. New York : Whiley Publ., 1995.
I. S. Masich
THE HEURISTIC ALGORITHMS OF BOUNDARY POINTS SEARCH FOR AN CONSTRAINT PSEUDO-BOOLEAN OPTIMIZATION PROBLEM
Application of the boundary point search scheme for solving an pseudo-Boolean optimization problem with variable constraints is investigated. We consider some heuristic algorithms based on that scheme. Comparison of algorithm efficiency with each other and an universal search algorithm is made for an appropriate real-world problem with large dimension.
