структуры и моделирование 2016. №2(38). С. 43-59
УДК 519.16
АЛГОРИТМЫ РЕШЕНИЯ ЗАДАЧИ ВЫБОРА ХАБА С ЗАДАННЫМ ЧИСЛОМ УЗЛОВ
А.В. Еремеев1
доцент, д.ф.-м.н., e-mail: [email protected] Э.А. Тарасенко2 магистрант, e-mail: [email protected]
1 Институт математики им. С.Л.Соболева Сибирского отделения РАН 2Омский государственный университет им Ф. М. Достоевского
Аннотация. Рассматривается задача выбора узлов хаба для рынка электроэнергии. Подмножество узлов электроэнергетической сети (хаб) используется для вычисления ценового индекса, который в каждый момент времени полагается равным среднему арифметическому цен электроэнергии по всем входящим в хаб узлам. В рассматриваемой задаче при известных значениях цен электроэнергии для каждого узла и каждого субъекта рынка за некоторый достаточно продолжительный период требуется выбрать хаб, состоящий из заданного числа узлов, с целью минимизации суммарного взвешенного среднеквадратического отклонения цен субъектов от ценового индекса хаба. Ввиду вычислительной сложности задачи для её решения предложены эвристические алгоритмы: метод локального поиска, генетический и гибридный генетический алгоритмы. Разработаны модификации предложенных алгоритмов для практически значимого варианта задачи с неопределёнными ценами в некоторые часы. Проведено тестирование алгоритмов на реальных данных и их сравнение с методами частично целочисленного программирования. Вычислительные эксперименты показали перспективность использования предложенных алгоритмов и выявили случаи, когда одни алгоритмы имеют преимущество перед другими.
Ключевые слова: поиск подмножества векторов, кластерный анализ, генетический алгоритм, локальный поиск, частично целочисленное программирование.
Введение
В современных рынках электроэнергии, основанных на узловом маржинальном ценообразовании, цена на электроэнергию не определяется единым значением, но зависит от узла электроэнергетической сети и периода времени [1]. При этом точное прогнозирование цены весьма проблематично, и участники рынка (оптовые продавцы и покупатели электроэнергии) заинтересованы в использовании ценовых индексов с целью снижения своих финансовых рисков
посредством заключения фьючерсных контрактов по цене этих индексов. Подмножество узлов электроэнергетической сети, называемое хабом, служит для вычисления ценового индекса, аппроксимирующего в некотором смысле цены электроэнергии субъектов рынка. Индекс хаба в каждый момент времени полагается равным среднему арифметическому цен электроэнергии по всем входящим в хаб узлам. Использование индекса хаба позволяет субъектам рынка заключать фьючерсные контракты по цене электроэнергии в хабе с целью снижения финансовых рисков, связанных с колебаниями цен на электроэнергию. При выделении хаба заданной мощности из найденного на предварительном этапе множества узлов (кластера) возникает задача выбора узлов хаба.
В работе предложены и реализованы алгоритм локального поиска, генетический и гибридный генетический алгоритм для задачи выбора узлов хаба. Предложена математическая модель для варианта задачи с неопределёнными ценами в некоторые часы. Проведён вычислительный эксперимент на задачах различной размерности, в ходе которого проведено сравнение результатов работы эвристических алгоритмов и пакета CPLEX 11.0.
1. Постановки задачи выбора узлов хаба и их свойства
В [4,5] сформулирована задача выделения заданного числа хабов из всего множества узлов энергетической сети как невыпуклая задача математического программирования. Данная модель является весьма сложной для стандартных методов невыпуклого программирования, поэтому вместо поиска ее точного решения может использоваться двухэтапный подход [4,5]: сначала с помощью методов кластерного анализа все множество узлов разбивается на кластеры, в которых узлы имеют «близкие» цены, а затем в каждом из кластеров выделяются хабы. Последняя задача является предметом исследования в настоящей работе.
1.1. Задача выбора узлов хаба
Рассмотрим математическую модель задачи выбора хаба с заданным числом узлов из множества узлов кластера. С учётом специфики рынка электроэнергии России будем полагать, что требуется минимизировать суммарное взвешенное квадратичное отклонение индекса хаба от цен субъектов рынка. Пусть цены электроэнергии для субъектов известны за весь анализируемый период (например, за прошедший год или несколько лет). Число узлов искомого хаба задано. Этот параметр выбирается достаточно большим, чтобы обеспечить предсказуемость динамики индекса хаба и ослабить влияние на него возможных удалений отдельных узлов из системы (на практике размер хабов, как правило, составляет несколько десятков или сотен узлов).
Параметрами модели являются следующие величины: п — число узлов заданного кластера;
т — число субъектов рынка, для которых выбирается хаб; N — число узлов искомого хаба;
Т — число часов в анализируемом периоде;
оц — цена ¿-го узла в час Ь (г = 1,..., п, Ь = 1,..., Т);
^ ^ 0 — суммарный объём покупки или продажи электроэнергии субъекта ] (] = 1,...,т), характеризующий значимость отклонения его цены от индекса хаба в момент времени Ь;
— цена электроэнергии, по которой её покупает или продаёт ]-й субъект рынка в час Ь = 1,..., т, Ь = 1,..., Т).
Переменные задачи: хг = 1, если г-й узел включается в хаб; иначе хг = 0 (г = 1,..., п); Oi — ценовой индекс хаба в час Ь (Ь = 1,..., Т).
Требуется минимизировать суммарное квадратичное отклонение индекса хаба от цен электроэнергии субъектов рынка с учётом весовых коэффициентов. Математическая модель выглядит следующим образом [4]:
T m
/ = - wt ^ min; (1)
t=1 j=1
1 П
ct = 1 = 1,...,T; (2)
i=1
n
ё Xi = N; (3)
i= 1
xi G {0,1}, i = 1,...N. (4)
Запишем (1)-(4) в виде модели частично целочисленного линейного программирования [4]. Введём обозначение: Ф^ = (ct — pjt)2. Тогда
$jt=(^N ёxicit—pj^j. (5)
Раскрывая скобки, получим:
1 n 2 n k-1
фjt = n^S (c2t— 2NcitPjt)xi + n^^EE cktcitxk xi + p2t. (6) i=1 k=1 i=1
Введём переменные zki, k = 1,..., n, l = 1,..., k — 1, подчинённые условиям
zk = xkxi, k = 1,..., n, l = 1,..., k — 1. (7)
Последнее равенство можно заменить на систему линейных ограничений:
Zki ^ xk, zk ^ xi, k = 1,..., n, l = 1,..., k — 1; (8)
zki ^ xk + xi — 1, k = 1,..., n, l = 1,..., k — 1. (9)
Запишем (1)-(4) с учётом (5)-(9).
n n k—1
Co + J] СгХ + У^ У^ Bk^Zki ^ min; (10)
i=1 k=l i=1 Zki ^ Xk + xi — 1, k = 1,..., n, l = 1,...,k - 1; (11)
n
^ Xi = N; (12)
¿=1
zki ^ Xi, k = 1,... ,n, l = 1,...,k — 1; (13)
X, G {0,1}, i = 1,...,n; (14)
zki G {0,1}, k = 1,..., n, l = 1,..., k — 1. (15)
Для целевой функции (10) ограничения (8) реализуются автоматически. Поэтому можно оставить только условия (9). Рассматривая zki как непрерывные переменные, получим модель частично целочисленного линейного программирования:
n n k—1
Co + ^ CiXi + У^ У^ BkiZki ^ min; (16)
¿=1 k=1 i=1 Zki ^ Xk + Xi — 1, k = 1,...,n, l = 1,..., k — 1; (17)
n
^ Xi = N; (18)
¿=1
X, G {0,1}, i = 1,...,n; (19)
zki ^ 0, k = 1,..., n, l = 1,..., k — 1. (20)
Из (6) следует, что коэффициенты в целевой функции имеют вид:
m T
co = J2Y1
j=11=1
1 m T
Ci = N2 ^ — ^С^К^ j = 1 t=1
2 m T Bki = CktCit^jt.
j=1 t=1
1.2. Вычислительная сложность задачи выбора узлов хаба
В [3] установлено, что задача выбора узлов хаба является NP-трудной в сильном смысле. Этот же вывод может быть получен исходя из связи рассматриваемой задачи с задачей выбора подмножества векторов с кратчайшим средним при ограничении на мощность (Subset with the Shortest Average with cardinality restriction, SSA) [2]. Последняя формулируется следующим образом.
Дано: множество Y = {ш,...,Уп} точек (векторов) из RT и натуральное число M.
Найти: подмножество C Ç Y такое, что |C| ^ M и
11I У^ У II —;► min.
|C|" U
Любая индивидуальная задача SSA сводится к последовательности из M индивидуальных задач выбора узлов хаба с теми же значениями n и T, что и в SSA, с узловыми ценами cit = yit, t = 1,...,T, с одним субъектом рынка, имеющим цены p1t = 0, t = 1,..., T и N = 1,..., M. Поскольку, как показано в [2], для задачи SSA не существует приближенных алгоритмов с какой-либо гарантированной оценкой относительной погрешности при P=NP, то такое же утверждение справедливо и для задачи выбора узлов хаба.
1.3. Задача выбора узлов хаба при неопределённых ценах в некоторые часы
На практике исходные данные не всегда содержат информацию о ценах в узлах и ценах субъектов рынка для каждого часа. Постановка, рассматриваемая в данном параграфе, отличается от предыдущей тем, что она учитывает возможность отсутствия в некоторые часы данных о ценах в узлах и ценах субъектов.
Заданы множества: N — множество узлов кластера, определённых в час ¿, |М| = Т — множество часов, когда определена цена электроэнергии субъекта ], Т- = 1Т1
Тз |Тзг
Требуется минимизировать суммарное квадратичное отклонение индекса ха-ба от всех цен покупки или продажи электроэнергии по субъектам рынка с учётом весовых коэффициентов. Построенная модель выглядит следующим образом:
- ^ min; (21)
j=i teTj
XiCit
ct = ^-, t =1,...,T; (22)
ieNt
Ë Xi = N ; (23)
i=1
X G {0,1}, i = 1,...,n; (24)
ct ^ 0, t = 1,...,T. (25)
В отличие от моделей, рассмотренных в п 1.1 и п. 1.2, данная модель содержит дробно-линейные ограничения, что существенно осложняет её решение, однако,
для практики эта постановка задачи имеет большое значение, поскольку данные о ценах узлов и субъектов по техническим причинам не всегда определены. Для решения поставленной задачи целесообразно применять эвристические алгоритмы.
2. Эвристические алгоритмы решения задачи выбора узлов хаба
Ввиду того, что задача выбора узлов хаба является NP-трудной в сильном смысле и неаппроксимируемой при Р=ЫР, для её решения целесообразно использовать эвристические алгоритмы, показавшие эффективность для решения задач с большой вычислительной сложностью. В настоящем разделе для решения задачи выбора узлов хаба (1)-(4) и варианта задачи с неопределёнными ценами в некоторые часы (21)-(25) предложены алгоритм локального поиска, генетический алгоритм и комбинация двух этих методов.
2.1. Алгоритм локального поиска
Пусть Б — множество допустимых решений. Окрестностью Ф(х0) С Б решения xo называется множество допустимых решений, в некотором смысле близких к данному. Идея алгоритма локального поиска состоит в следующем. Имеется допустимое решение задачи х. Строится окрестность Ф(х). Из окрестности выбирается решение х', дающее улучшение целевой функции. Различают несколько вариантов улучшения текущего решения: из окрестности может выбираться решение, дающее лучшее значение целевой функции, либо из окрестности выбирается первое улучшение целевой функции. После перехода в решение х' строится окрестность Ф(х'), и процесс продолжается. Если улучшить текущее решение не удаётся, то оно является локальным оптимумом.
Рассмотрим алгоритм локального поиска для задачи выбора узлов хаба. Расстоянием Хэмминга Л,(х,х') для булевых векторов х и х' называется число несовпадающих компонент. Окрестностью к-зи'ар вектора х называется множество таких точек х', что Л,(х,х') = 2к, где Л,(х,х') — расстояние Хэмминга для булевых векторов х и х'. Рассмотрим окрестность Ь^ар, в эту окрестность включаются те решения, которые получаются из текущего исключением из ха-ба одного узла и включением другого узла. Допустимым решением будет булев
п
вектор: х = (хь ... ,хп) такой, что ^ хг = N, где хг=1, если узел г включили в
г=1
хаб. В предложенном алгоритме в окрестности текущего решения выбирается решение, дающее наибольшее улучшение целевой функции.
При переходе от текущего решения х к новому допустимому решению х' используется быстрый пересчёт значения целевой функции для задачи (1)-(4):
Б(х') = Б(х) - Сг + С - X] Виг ^.
м:жи=1
В случае задачи (21)-(25) значение функции Б(х') вычисляется следующим образом:
Т т
Б(х') = ЕЕ(с(х') -
4=1 ¿=1
где
с(х0 — Сг4 + Сд
ct(x') =
2.2. Генетический алгоритм
Exi
¿eNí
Генетический алгоритм (ГА) представляет собой метод оптимизации, основанный на эволюционных принципах. Рассмотрим задачу оптимизации в общем виде:
max{f (x) : x Е X},
где X — пространство решений.
Популяцией П = ... , £г) численности l называется вектор простран-
ства B¿,B = {0,1}u, координаты которого называются генотипами особей. Целевая функция f исходной задачи заменяется в генетическом алгоритме на функцию приспособленности генотипа F. С помощью данной функции можно задать для каждой особи среднее число оставляемых ею потомков или вероятность того, что эта особь оставит потомка. Алгоритм начинает свою работу с некоторой начальной популяции П0. Шагом генетического алгоритма является переход от текущей популяции к следующей, то есть получение новой популяции ní+1 из П. В генетическом алгоритме новые решения строятся как комбинации элементов решений текущей популяции с использованием операторов селекции, скрещивания (кроссинговера) и мутации.
Оператором селекции Sel : B1 ^ {1,...,l} называется одноместный оператор, действие которого заключается в выборе номера родителя в популяции. Оператор селекции особей на пространстве популяций имеет то же значение, что и естественный отбор в природе.
Оператором скрещивания Cross : B х B ^ B х B называется двуместный оператор, действие которого заключается в построении генотипов потомков путём комбинации участков родительских генопипов.
Широко используемый одноточечный оператор кроссинговера [6] действует следующим образом. Пусть генотипы родителей имеют вид:
£ =(£ъ6,...,Ы,
П = (П1,П2,... ,Пи), тогда результат кросинговера (£',п') = Cross(£,n) формируется в виде:
£' = (£ъ&,..., £х,Пх+1,...,П«),
п' = (п1,п2,...,пх,£х+1 ,...,£u
где случайная позиция х выбрана с равномерным распределением от 1 до u — 1. Оператор кроссинговера применяется в алгоритме с вероятностью Pcross. В случае, когда оператор кроссинговера не применяется, родительские генотипы копируются в новую популяцию.
Оператором мутации называется одноместный оператор Mut : B ^ B, который вносит случайные изменения в генотип особи. Оператор мутации [6] в каждой позиции аргумента с заданной вероятностью Pmut заменяет её содержимое на другой элемент алфавита {0, 1}.
Предлагаемый генетический алгоритм для решения задачи выбора узлов хаба описывается следующим образом. Генотипом является вектор x G Bn,
n
причём X = N. Популяцией численности e является вектор П = ..., xe}. i=1
Функция приспособленности генотипа
T m
F (et) = —/(et) = — Y1 — Pjí)2wjt. t=1 j=1
Генетический алгоритм
1. Случайным образом строится начальная популяция П0.
2. Пока не превышено максимальное число итераций, выполняется:
2.1. Оператор селекции Se/'. Из текущей популяции выбираются k особей с равномерным распределением, из них выбираются две особи £ и n с наибольшим значением функции приспособленности.
2.2. Оператор скрещивания Cross'. Два потомка £' и n' строятся с помощью оператора кроссинговера, описываемого ниже.
2.3. Оператор мутации Mut'. К потомкам £' и n' применяется оператор мутации с вероятностью Pmut, который также описан ниже.
2.4. Из двух особей потомков £' и n' и двух особей родителей £ и n выбираются две особи с наибольшим значением функции приспособленности и помещаются в популяцию П'.
Рассмотрим оператор кроссинговера Cross', предложенный для рассматриваемой задачи. Особенностью данного оператора является то, что потомки
n
должны наследовать свойство родителей ^ x = N. Сохранение этого свойства
i=1
обеспечивается следующим правилом формирования потомков:
1) если в потомке, полученном в результате применения оператора кроссинговера, описанного в п. 2.2.1, имеется N' единиц, причём N' > N, то из множества позиций, содержащих единицы, с равномерным распределением выбираются N' — N позиций, которые заполняются нулями;
2) если в потомке, полученном в результате применения оператора кроссин-говера, имеется N' единиц, причём N' < N, то из позиций, в которых стоят
единицы в родительских особях, а в потомках стоят нули, случайным образом выбираются N - N позиций, которые заменяются единицами.
Рассмотрим оператор мутации Мм£', применяемый в предложенном алгоритме. Отличие данного оператора мутации от оператора Ми классического генетического алгоритма состоит в том, что сохраняется свойство родительских
п
решений = N. Для этого каждый элемент вектора с вероятностью Рти
г=1
заменяется на противоположный, при этом подсчитывается — количество замен 1 на 0 и — количество замен 0 на 1, затем — раз вы-
полняются обратные замены для — случайно выбранных позиций
из числа тех, которые подвергались замене.
Рассмотрим вычисление функции приспособленности генотипа для задач (1)-(4) и (21)-(25). Пусть х — генотип, для которого известно значение функции приспособленности Б(х), значение индекса хаба с4(х) и число узлов хаба ^(х),£ = 1,...,Т. Вычисляем Б(х') для задачи (1)-(4) следующим образом:
С(х0 = N ( С(х) + Е Сй - Е ^
У ге/х ге/-х
Тт
Б (х') = ЕЕ (С4(х') - ¿)2
4=1 ¿=1
Для задачи (21)-(25) вычисление Б(х') происходит следующим образом:
С(х) + Е Сц - Е Сг4 / 'ч _ ге/1 ге/-1
С4(х)= Ж4(х) + ЕП=1(Хг - Х) ,
Тт
Б(х') = ЕЕ(с4(х') - ¿)Ч,,
4=1 ¿=1
где /1 = (г : ж - ж» = 1}, /-1 = (г : ж - ж» = -1}.
2.3. Комбинация генетического алгоритма с алгоритмом локального поиска
Для решения задачи выбора узлов хаба (1)-(4) и варианта задачи с неопределёнными ценами в некоторые часы (21)-(25) был также предложен гибридный генетический алгоритм. Идея гибридных алгоритмов заключается в сочетании генетического алгоритма с некоторым другим алгоритмом решения, в котором учитывается специфика задачи. В данном случае в качестве второго метода использовался алгоритм локального поиска. На каждом поколении полученный потомок оптимизировался этим методом, после чего добавлялся в популяцию.
Рассмотрим схему гибридного генетического алгоритма для задачи выбора узлов хаба (1)-(4) и варианта задачи с неопределёнными ценами в некоторые
часы (21)-(25). Для каждой из двух полученных особей на итерации генетического алгоритма применяется алгоритм локального поиска, полученные локальные оптимумы помещаются в популяцию П'. Данная схема в англоязычной литературе носит название «memetic algorithm».
3. Вычислительный эксперимент
В настоящем разделе приведены результаты тестирования предложенных выше алгоритмов для задачи выбора узлов хаба. Для проведения экспериментов алгоритмы были реализованы на языке С++ и тестировались на ЭВМ Athlon II X4 с тактовой частотой процессора 2,9 ГГц и объёмом оперативной памяти 2 Гб. В качестве тестовых примеров использовались задачи, построенные с использованием исходных данных о ценах электроэнергии «PJM Interconnection» (США), доступных на сайте http://www.pjm.com, а также реальные данные рынка электроэнергии России.
Для решения задачи выбора узлов хаба (1)-(4) и варианта задачи с неопределёнными ценами в некоторые часы (21)-(25) реализованы алгоритм локального поиска, генетический алгоритм и гибридный генетический алгоритм.
В качестве тестовых использовались задачи, представленные в табл. 1. В данной таблице содержится описание параметров задач, таких как число узлов заданного кластера n, число субъектов рынка m, число узлов хаба N, число часов T, также в последнем столбце обозначена принадлежность соответствующего кластера США или России. Задачи 1-8 представляют данные для хабов США и содержат полную информацию о ценах в узлах и ценах в субъектах рынка, они соответствуют модели (1)-(4). Задачи 9-14 представляют данные из европейской и сибирской ценовых зон соответственно, в некоторые часы данные не содержат информацию о ценах в узлах или в субъектах, эти задачи соответствуют модели (21)-(25).
3.1. Выбор параметров генетического алгоритма
Процедура подбора параметров генетического алгоритма решения задачи выбора узлов хаба (1)-(4) проводилась по следующей схеме.
Этап 1. Подбор вероятности мутации Pmut с шагом 0.1 (таким образом рассматриваются все варианты Pmut = 0,0.1,..., 1), при этом вероятность запуска оператора кроссинговера Pcross = 0.
Этап 2. Подбор размера популяции. Рассматриваются варианты: 50, 100, 200, 300, 400, 500. При этом вероятность мутации фиксирована и равна лучшему значению параметра, полученного на этапе 1.
Этап 3. Подбор размера турнира k. Осуществляется путём запуска алгоритма с числом участников турнира k = 2,5,20.
На этапе 1 для задачи (1)-(4) относительное отклонение от лучшего известного значения функции оказалось наименьшим при вероятности мутации Pmut = 0.1. На этапе 2 лучшие значения целевой функции для всех рассматриваемых задач получены при размере популяции равном 100. Исключение со-
Таблица 1. Параметры тестовых примеров
№ п т N Т Географическое расположение
1 43 5 20 24 США, РЛМ
2 152 10 20 24 США, РЛМ
3 199 20 10 24 США, РЛМ
4 199 30 100 24 США, РЛМ
5 233 100 5 24 США, РЛМ
6 408 5 100 24 США, РЛМ
7 642 30 10 24 США, РЛМ
8 642 5 100 24 США, РЛМ
9 63 204 50 20472 Россия (европейская ценовая зона)
10 95 181 76 20472 Россия (европейская ценовая зона)
11 259 266 207 20472 Россия (европейская ценовая зона)
12 411 656 330 20472 Россия (европейская ценовая зона)
13 49 52 40 20472 Россия (сибирская ценовая зона)
14 80 52 64 20472 Россия (сибирская ценовая зона)
ставляла задача 6, для которой лучшие результаты были получены при размере популяции равном 50. На этапе 3 лучшее значение целевой функции для всех рассматриваемых задач получено при к = 20. Процедура подбора параметров генетического алгоритма для варианта задачи с неопределёнными ценами в некоторые часы дала аналогичные результаты.
3.2. Решение задач выбора узлов хаба
Рассмотрим решение задачи выбора узлов хаба (1)-(4). Данная задача была решена с помощью алгоритма локального поиска, генетического алгоритма, а также гибридного генетического алгоритма. Для каждого набора тестовых данных 1-8 алгоритм локального поиска применялся 20 раз, в качестве результата выбиралось лучшее из 20 решений.
Для гибридного генетического алгоритма применены те же настройки параметров, что и для генетического алгоритма. Для решения задач генетическим алгоритмом и гибридным генетическим алгоритмом задавалось максимальное время решения, равное времени решения задачи алгоритмом локального поиска. Результаты эксперимента для всех рассмотренных алгоритмов представлены для задач 1-8 в табл. 2. В таблице указано число узлов заданного множества, число субъектов рынка, число узлов хаба, а также полученные рассматриваемыми алгоритмами значения целевой функции и время счёта т в секундах. Лучшие известные значения целевой функции в этой и последующих таблица обозначены жирным шрифтом.
Таблица 2. Результаты эксперимента по сравнению эвристических алгоритмов
№ n m N Локальный Генетический Гибридный
поиск алгоритм алгоритм
f т, сек. f f
1 43 5 20 4.993 1 4.993 4.993
2 152 10 20 0.1857 1 0.1857 0.1857
3 199 20 10 0.0478 3 0.0478 0.0478
4 199 30 100 3.338 21 3.338 3.338
5 233 100 5 61.057 4 61.123 61.057
6 408 5 100 1315.440 309 1315.869 1315.468
7 642 30 10 4413.308 17 4413.279 4413.317
8 642 5 100 0.0049 582 0.0062 0.0052
Из представленной таблицы видно, что рассматриваемые алгоритмы получали одинаковое значение целевой функции для задач 1-4. На задачах 5, 6, 8 гибридный генетический алгоритм получал значения целевой функции меньше, чем генетический алгоритм за то же самое время, кроме того на этих задачах алгоритм локального поиска получал значение целевой функции меньше, чем гибридный генетический алгоритм. Значение функции, полученное генетическим алгоритмом для задачи 7 меньше, чем значение функции, полученное алгоритмом локального поиска и гибридным генетическим алгоритмом. Результаты проведённого вычислительного эксперимента позволяют сделать вывод, что алгоритм локального поиска наряду с гибридным генетическим алгоритмом успешно решают задачи малой размерности.
Для сравнения полученных результатов работы алгоритма локального поиска с коммерческим пакетом прикладных программ была составлена модель частично целочисленного линейного программирования (10)-(15) с помощью пакета GAMS 22.6, а также модель квадратичного программирования (1)-(4). Все тестовые задачи решены с помощью пакета CPLEX 11.0, входящего в состав GAMS. При решении задач квадратичного частично целочисленного программирования использовался параметр CPLEX miqcp=1 (решение задачи квадратичного программирования в каждом узле дерева ветвей и границ).
Решение задач с помощью пакета CPLEX 11.0 для каждой модели было проведено в два этапа:
1. Задавалось максимальное время решения, равное времени решения задачи алгоритмом локального поиска (rmax = tl).
2. Задавалось максимальное время решения, равное 30 минутам. В обоих случаях максимальное число итераций задавалось равным 108, но в ходе тестирования предел числа итераций не превышался.
Результаты экспериментов приведены в сводных таблицах (в табл. 3 — результаты сравнения алгоритма локального поиска и пакета СРЬЕХ в режиме частично целочисленного линейного программирования, в табл. 4 — того же алгоритма и пакета СРЬЕХ в режиме квадратичного программирования). В обеих таблицах для каждой задачи указано число узлов заданного множества, число субъектов рынка, число узлов в хабе, а также полученные алгоритмом локального поиска и пакетом СРЬЕХ 11.0 значения целевой функции и время счета т в секундах. Символ «*» обозначает оптимальные решения. В скобках указано время доказательства оптимальности пакетом СРЬЕХ.
Таблица 3. Сравнение алгоритма локального поиска и пакета СРЬЕХ в режиме частично целочисленного линейного программирования
№ п т N Локальный поиск СРЬБХ 11.0
! т, сек. ттах тЬ ттах = 1800С
! число итер. ! число итер.
1 43 5 20 4.993 1 20.337 3031 7.711 974854
2 152 10 20 0.1857 1 0.1874 412 0.1858 604440
3 199 20 10 0.0478 3 0.051 251 0.048 453946
4 199 30 100 3.338 21 3.478 1874 3.334 32189
5 233 100 5 61.057* 4 78.339 59 61.939 4396
6 408 5 100 1315.4* 309 3290.99 7789 3290.98 19555
7 642 30 10 4413.3 17 5245.243 96 5245.242 2011
8 642 5 100 0.0049 582 11.972 4396 0.0109 7234
Как видно из табл. 3, на первом этапе СРЬЕХ 11.0 во всех задачах получал решение со значением целевой функции, большим, чем алгоритм локального поиска. На втором этапе СРЬЕХ 11.0 за время т = 30 мин. не доказал оптимальность решения для задач 1-8 и получал значение целевой функции, большее, чем алгоритм локального поиска.
Как видно из табл. 4, для задач 1, 5, 6, 8 с помощью пакета СРЬЕХ 11.0 установлена оптимальность найденных решений. В задачах 5, 6 алгоритм локального поиска получил оптимальные решения. В задачах 2, 3, 5, 7 значения целевой функции, полученные на первом этапе пакетом СРЬЕХ, оказались больше значений, полученных алгоритмом локального поиска.
Проведённое сравнение алгоритма локального поиска и пакета СРЬЕХ показывает, что алгоритм локального поиска, наряду с пакетом СРЬЕХ в режиме квадратичного программирования, может эффективно применяться для решения задачи выбора узлов хаба. Однако использование коммерческого пакета
Таблица 4. Сравнение алгоритма локального поиска и пакета СРЬЕХ в режиме квадратичного
программирования
№ п т N Локальный поиск CPLEX 11.0
! т Ттах = 1800С
! !
1 43 5 20 4.993 1 4.992* (1 с) 4.992* (1 с)
2 152 10 20 0.1857 1 0.1858 0.1856
3 199 20 10 0.0478 3 0.0479 0.0478
4 199 30 100 3.338 21 3.338 3.338
5 233 100 5 61.057* 4 61.085 61.057* (13 с)
6 408 5 100 1315.4* 309 1315.4* (3 с) 1315.4* (3 с)
7 642 30 10 4413.3 17 5168.2 4413.3
8 642 5 100 0.0049 582 0.0039* (4 с) 0.0039* (4 с)
СРЬЕХ 11.0 в режиме частично целочисленного линейного программирования оказалось нецелесообразным.
3.3. Решение задач выбора узлов хаба при неопределённых ценах в некоторые часы
Для решения варианта задачи с неопределёнными ценами в некоторые часы (21)-(25) реализован алгоритм локального поиска, генетический алгоритм, а также гибридный генетический алгоритм, описанные в разделе 2. В качестве данных о ценах в узлах, ценах субъектов рынка и объёмах покупки или продажи электроэнергии используются реальные данные рынка электроэнергии России. В качестве тестовых задач используются четыре кластера из европейской ценовой зоны (задачи 11-14) и два кластера из сибирской ценовой зоны (задачи 9, 10) — см. http://www.mosenex.ru/. Алгоритм локального поиска запускается один раз, число итераций генетического алгоритма равно 2000. Максимальное время работы гибридного генетического алгоритма задаётся равным времени решения задачи генетическим алгоритмом. Локальный поиск в гибридном алгоритме запускается с вероятностью Р1оса = 0.2. Результаты экспериментов представлены в табл. 5, которая содержит информацию о числе узлов заданного множества, числе субъектов рынка, числе узлов хаба, а также полученные рассматриваемыми алгоритмами значения целевой функции и время счета т в секундах.
Во всех задачах генетический алгоритм и алгоритм локального поиска получали близкие значения целевой функции, однако, время работы алгоритма
Таблица 5. Решение задач выбора узлов хаба при неопределённых ценах в некоторые часы
№ п т N Локальный Генетический Гибридный
поиск алгоритм алгоритм
! т, сек. ! т, сек. !
9 63 204 50 2.5600е12 1917 2.5601е12 2213 2.5600е12
10 95 181 76 4.3826е12 6748 4.3839е12 3239 4.3839е12
11 259 266 207 2.3320е12 160650 2.3343е12 18677 2.3358е12
12 411 656 330 3.4374е12 518400 3.4386е12 71430 3.4386е12
13 49 52 40 6.1329е11 211 6.1329е11 454 6.1328е11
14 80 52 64 6.5781е11 970 6.5842е11 721 6.5781е11
локального поиска даже на одной итерации оказалось больше, чем время работы генетического алгоритма, исключение составляют задачи 9, 13. Для задач 9, 14 гибридный алгоритм получал такое же значение целевой функции, что и алгоритм локального поиска, для задачи 13 гибридный алгоритм получал меньшее значение функции, чем алгоритм локального поиска.
Проведённый вычислительный эксперимент показывает, что для варианта задачи с неопределёнными ценами в некоторые часы (21)-(25) может эффективно применяться алгоритм локального поиска и гибридный генетический алгоритм для задач небольшой размерности. Для задач большой размерности целесообразно применять генетический алгоритм, в то время как применение алгоритма локального поиска или гибридного генетического алгоритма нецелесообразно в связи с большой трудоёмкостью локального поиска (быстрый пересчёт значения целевой функции в данном случае не применим).
4. Заключение
В данной работе рассмотрены две постановки задачи выбора узлов хаба, предложена математическая модель варианта задачи с неопределёнными ценами в некоторые часы (21)-(25). Предложены и реализованы алгоритм локального поиска, генетический и гибридный генетический алгоритмы. Все алгоритмы были реализованы и протестированы на задачах различной размерности. Проведено сравнение алгоритмов между собой.
Полученные результаты позволяют сделать вывод о перспективности применения предложенных эвристических алгоритмов для приближенного решения задачи выбора узлов хаба большой размерности. В случае, когда число субъектов мало, точные решения могут быть получены пакетом СРЬЕХ в режиме квадратичного частично целочисленного программирования. Использование пакета СРЬЕХ в режиме частично целочисленного линейного программирования представляется нецелесообразным.
Для приближённого решения задачи с неопределёнными ценами в неко-
торые часы целесообразно применять генетический алгоритм, в то время как использование алгоритма локального поиска или гибридного генетического алгоритма оказалось менее эффективно в связи с большой трудоёмкостью локального поиска.
Благодарности
Работа выполнена при поддержке Программы фундаментальных научных исследований государственных академий наук на 2013-2020 годы, п. 1.5.1.5. «Исследование и решение задач комбинаторной оптимизации с использованием целочисленного программирования».
Литература
1. Давидсон М.Р., Догадушкина Ю.В., Крейнес Е.М., Новикова Н.М., Удальцов Ю.А., Ширяева Л.В. Математическая модель конкурентного оптового рынка электроэнергии в России // Изв. РАН. Теория и системы управления. 2004. № 3. С. 72-83.
2. Еремеев А.В., Кельманов А.В., Пяткин А.В. О сложности некоторых евклидовых задач оптимального суммирования. Доклады Академии наук. 2016. т. 468, N 4. C. 372-375.
3. Тарасенко Э.А. О сложности задачи выбора узлов хаба // Труды региональной научной студенческой конференции «Молодежь III тысячелетия». Омск : ОмГУ, 2010. С. 45-48.
4. Borisovsky P.A., Eremeev A.V., Grinkevich E.B., Klokov S.A., Vinnikov A.V. Trading Hubs Construction for Electricity Markets // In: Optimization in the Energy Industry / Kallrath J., Pardalos P.M., Rebennack S. Springer : Berlin, 2009. P. 27-51.
5. Hartshorn A., Hastings D. Midwest ISO Trading hubs development. Cambridge, MA : LECG, LLC, 2003. URL:~https://www.misoenergy.org/_layouts/ miso/ecm/redirect.aspx?id=30607 (дата обращения: 12.02.2016).
6. Holland J. H. Adaptation in Natural and Artifical systems. Ann Arbor, MI : University of Michigan Press, 1975.
7. Reeves C. Genetic algorithm for the operations researcher // INFORMS Journal on Computing. 1997. V. 9, № 3. P. 231-250.
ALGORITHMS FOR SOLVING THE HUB CONSTRUCTION PROBLEM WITH
GIVEN NUMBER OF HUB NODES
A.V. Eremeev1
Associate Professor, Dr. Sci. (Phys.-Math.), e-mail: [email protected]
E.A. Tarasenko2 MSci. Student, e-mail: [email protected]
1Omsk Branch of Sobolev Institute of Mathematics SB RAS 2Dostoevsky Omsk State University
Abstract. The trading hub construction problem for electricity markets is considered. A subset of nodes of an electricity grid (called a hub) is used for computing the price index, which at every moment is assumed to be equal to the arithmetic mean of electricity prices over all nodes of the hub. Given historical prices for all nodes and for all market participants over a sufficiently long period of time, the problem is to choose a hub consisting of a given number of nodes so as to minimize the sum of squared differences between the hub price and the prices of participants. In view of problem complexity, the heuristic algorithms are proposed for its solution: a local search method, a genetic algorithm and a memetic algorithm. Modifications of these algorithms were developed for a practically significant version of the problem where the prices for some hours are undefined. The algorithms are tested and compared to the integer programming methods on the real-life data. Computational experiments indicated prospectiveness of usage of the proposed algorithms and displayed the cases where some algorithms are in advantage to others.
Keywords: vectors subset search, cluster analysis, genetic algorithm, local search, mixed integer programming.