О быстром варианте алгоритма отжига Текст научной статьи по специальности «Математика»

УДК 519.626

О БЫСТРОМ ВАРИАНТЕ АЛГОРИТМА ОТЖИГА

А.С.Тихомиров

Институт электронных и информационных систем НовГУ, Alexey.Tikhomirov@novsu.ru

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

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

Fast versions of the simulated annealing algorithm are constructed. It is shown, that if the objective function is non-degenerate the number of its evaluations required to obtain the desired accuracy £ in the solution can be a slowly (namely, logarithmically) growing function as £ approaches zero.

Keywords: simulated annealing algorithm, random search, stochastic search, global optimization, stochastic optimization

1. Введение

Пусть целевая функция f: X a R (где, например, X = Rd ) принимает минимальное значение в единственной точке x*. Рассмотрим задачу поиска точки глобального минимума x* с заданной точностью є. Один из способов решения этой задачи состоит в применении алгоритмов случайного поиска экстремума функции (см. [1-15]). Такие методы давно и успешно используются при решении сложных задач оптимизации. Однако теоретических результатов о скорости сходимости этих алгоритмов мало (см. [3-6]). Данная работа посвящена исследованию скорости сходимости варианта, пожалуй, одного из самых знаменитых алгоритмов стохастической глобальной оптимизации — алгоритма отжига (simulated annealing).

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

функции, требуемого для достижения заданной точности е решения задачи, при стремлении е к нулю. Такой упрощенный подход объясняется тем, что теоретическое исследование скорости сходимости алгоритма отжига оказалось сложной задачей. Известные теоретические оценки скорости сходимости этого алгоритма или носят асимптотический характер или очень неточны (см. [4-7]). Кроме того, эти оценки показывают крайне медленную скорость сходимости алгоритма отжига. Поэтому основное внимание мы уделим изучению порядка скорости сходимости.

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

Метод оптимизации будем называть «быстрым», если число вычислений целевой функции, требуемое для достижения заданной точности е решения задачи, имеет логарифмический (медленный) порядок роста при стремлении е к нулю (т.е. имеет вид 0(|1пе |а) при а > 0 ). Метод оптимизации назовем «медленным», если он имеет гораздо более худшую

— степенную (т.е. 0(1/еа) при а > 0) зависимость требуемого числа вычислений целевой функции от £.

По известным теоретическим оценкам скорости сходимости (см., напр., [4-7]) алгоритм отжига является медленным методом оптимизации. В лучшем случае зависимость требуемого числа вычислений целевой функции от е имеет вид O(1/ £d), где d

— размерность пространства оптимизации. Такой же порядок зависимости от е имеет «чистый случайный поиск», т.е. простейший алгоритм случайного бросания точек в множество оптимизации, см. [17]. Процитируем: «The theoretical rate of convergence of the simulated annealing is very slow; this convergence is based on the convergence of the pure random search which is contained within the simulated annealing algorithms ...» [18].

Быстрые варианты алгоритма отжига получены в работах [14] и [15]. Данная работа продолжает и дополняет результаты статей [14] и [15]. В этой работе получена асимптотическая оценка вероятности «успеха» для одного из вариантов алгоритма отжига, и с ее помощью показано, что (при оптимизации невырожденных целевых функций) алгоритм отжига может быть быстрым, если его хорошо организовать. Результаты работы опровергают распространенное мнение о том, что с теоретической точки зрения алгоритм отжига является медленным методом оптимизации.

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

Назовем пространством оптимизации множество оптимизации X, снабженное метрикой р. Пусть ц

— мера на с-алгебре борелевских подмножеств X. Замкнутый шар радиуса r с центром в точке x обозначим как Br (x) = {y е X: р(x, y) < r} . Везде в дальнейшем будем считать, что выполнены следующие условия.

Условие а. 0 < д(Br (x)) < +да для всех x е X и любого r > 0 .

Условие б. д(Br (x)) = д(Br (y)) для всех

x, y е X и любого r > 0.

Основным практическим примером таких метрических пространств будет, конечно, d-мерное евклидово пространство Rd с какой-либо «обычной» метрикой р и мерой Лебега ц.

В дальнейшем всегда будет предполагаться, что целевая функция f: X a R ограничена снизу, измерима и удовлетворяет следующим условиям:

1) принимает минимальное значение в единственной точке x*;

2) непрерывна в точке x*;

3) неравенство inf{ f (x): x g Br (x*)} > f (x*) верно для любого r > 0 .

Отметим, что функции указанного класса могут быть многоэкстремальными в любой окрестности глобального минимума. Еще одно условие на целевую функцию f и определение невырожденности f будут введены ниже.

Опишем исследуемый алгоритм отжига {4 n }n>0 с помощью алгоритма моделирования (подробнее см. [4], [14], [15]).

Алгоритм 1

Шаг 1. 4 0 ^ x , n ^ 1.

Шаг 2. ?n ^ Pn (Sn-1, •).

ТТТ •, £ [ ^n , с вероятностью ^(Pn , An ),

Шаг 3. 4n

14n-1, с вероятностью 1 - p(en , Аn ).

Шаг 4. Если n < N, то n ^ n +1 и перейти к шагу 2, иначе — STOP.

Здесь x — начальная точка поиска, N — число шагов поиска. Обозначение « qn ^ Pn (4n-i, •) » читается как «получить реализацию случайной величины q n с распределением Pn (4n-1, •)». В соответствии со структурой алгоритма 1 распределения Pn (4n-1, •) будем называть пробными переходными функциями, а случайные величины qn — пробными точками. Вероятности p(P n, Д n), называемые вероятностями принятия пробных точек, задаются следующей формулой

[1, если Д „ < 0, p(Pn, Д n) = 1 ( „ . ) . 0

Iexp(-P nД n X если Д n > 0

где Дn = f (qn) - f (4n-1), а величины pn > 0 являются параметрами алгоритма.

Ограничимся рассмотрением алгоритма отжига, пробными переходными функциями которого служат равномерные распределения в шарах. Кроме того рассмотрим поведение алгоритма отжига только при достаточно больших значениях параметра Р n. Далее будет показано, что и такой просто устроенный вариант алгоритма отжига может быть быстрым.

Равномерное распределение в шаре Ba (x) радиуса a > 0 с центром в точке x е X обозначим через Ua (x, • ) , т.е. положим Ua (x, • ) = д( • П Ba (x))/Д^ (x)) . Пусть k — натуральное число, 0 < rk <... < r0,

aj = rj + r;-1 при j = 1,.,k . Пусть n^...,nk — натуральные числа, Nj = n1 +... + nj при j = 1,.,k и N0 = 0 . Рассмотрим класс поисков алгоритма 1, для которых общее число шагов N равно Nk , а пробные переходные функции Pn = Uaj при Nj-1 < n < Nj .

Числа k, rj, nj являются параметрами метода. Здесь k — число «этапов» поиска, а nj — число шагов на j-м этапе, причем радиусы aj шаров, соответствующих переходным функциям Ua j , постоянны для всех

шагов одного этапа.

При отыскании точки минимума x* с заданной точностью £ > 0 нас должно интересовать попадание поиска в шар B£ (x*). Может, однако, случиться так, что поиск, оказавшись в B£ (x*) на шаге n, выйдет из B£ (x*) на одном из последующих шагов. Чтобы избежать анализа таких эффектов, введем множества M£ = {xеB£(x*):f(x)<f(y) для любого ygB£(x*)}.

Рассмотрим величины 4„ = а^шіп{ f (40),..., f (4п)} . Случайная величина 4 „, попав в множество Mє, из него больше не выйдет. Поэтому мы будем изучать момент попадания поиска в множество МЕ, где є сохраняет смысл требуемой точности поиска. Соответственно мерой близости точки х к х* оказывается

не расстояние р( х, х*), а число

5( х) = ігі"{г > 0: х є Мг}.

Нам потребуется еще одно ограничение на поведение целевой функции £ а именно:

4) и г>о Мг = X .

Это условие гарантирует конечность величины 5(х) при любом х є X .

В задачах оптимизации важно учитывать информацию о поведении целевой функции. Ниже вся информация о целевой функции £ будет содержаться в виде коэффициента асимметрии ¥£ (г) = д(Мг) / д(Вг (х*)).

В силу условий, наложенных на целевую функцию, ¥£ (г) > 0 при всех г > 0. Функции, у которых

іішіпГ ¥£ (г) > 0 при г а 0, будут называться невырожденными. Подробнее наложенные на £ ограничения и свойства коэффициента асимметрии ¥ £ обсуждаются в [4, 12, 13].

3. Исследование скорости сходимости

Справедлива следующая теорема.

Теорема. Пусть целевая функция £: X а Я удовлетворяет условиям 1)-4). Пусть є > 0 и начальная точка поиска х є X фиксированы. Пусть {4 „ }„>0 — описанный алгоритм отжига, у которого гк < є , г0 > 5(х), а параметры в„ не зависят от номера шага „ (т.е. в„ = в ). Пусть 0 < т}- < д(М) и

Ф= д(Ва, (х*)). Тогда

1 -I 1 -

Ф1

lim P(4N єMє) >

р^+ад

n ^ k f f m Y

П

j=2

1-

1 --

j-1

1-

Ф

і У

. (*)

Покажем, что оценку (*) можно использовать для построения быстрого поиска. Рассмотрим пространство X = [0,1)й с ^-мерной мерой Лебега д и метрикой

рш (x,у) = тахтш{| хп -уп |Д-1 хп -уп |},

1<п<^

где х = (х1,..., хй) и у = (у1,...,уа). Для выбранной метрики пространство [0,1)а топологически является ^-мерным тором. Выбор метрики тора в кубе [0,1)й вместо обычной метрики позволяет избежать изучения «краевых эффектов», возникающих вблизи границы куба. В то же время для изучения порядков скорости сходимости случайного поиска выбор тора в качестве пространства оптимизации не является

принципиальным. Действительно, при х* е (0,1)а

малые окрестности точки x* в кубе и торе совпадают, и результаты скорости сходимости, полученные для тора, автоматически переносятся на куб. Выбор метрики тора сказывается здесь только на конструкции поиска.

Рассмотрим невырожденную целевую функцию f с коэффициентом асимметрии F f. Пусть F * — оценка снизу коэффициента асимметрии, и пусть inf F* > 0 . Пусть r0 = 1/ 2 . Выберем число 0 < q < 1 и положим rj = q} /(2 + 2q) при j > 1. Пусть 0 <е < 1/4 и k = min{ j: rj < е}. Зафиксируем число 0 < у < 1. В качестве значений параметров n j возьмем наименьшие натуральные числа, при которых значения соответствующих множителей оценки (*) превосходят

величину у17k (оценка (*) представляет собой произведение k множителей). Для описанного поиска верно неравенство lim P(4n е M£) > у, и кроме того

N = 0(| lne | ln | lne |). Таким образом, при всех достаточно больших значениях параметра в мы получаем быстрый вариант алгоритма отжига. Число шагов такого поиска, при котором достижение множества M £ гарантировано с вероятностью не меньшей чем у, логарифмически зависит от е и сравнимо по порядку зависимости от е с быстрыми методами локальной оптимизации.

В заключение приведем пример применения исследуемого варианта алгоритма отжига. Пример взят из известной книги [5] (см. с.219-220). Здесь

пространство оптимизации X = R2, x = (x1, x2), f (x) = f (x1,x2) = xj4 + xj2 + x1 x2 + x| . Функция f принимает минимальное значение в единственной точке x* = (0,0) и f (x*) = 0. Начальной точкой поиска выбрана x = (1,1) и f (x) = 4. Три столбца (T = 0,01, T = 0,1, T = 1) нижеследующей таблицы составляют данные таблицы 8.2 из [5], в последнем ее столбце (Алгоритм 1) представлены результаты поиска данной работы.

В таблице 8.2 даны средние значения функции f, полученные при 40 повторениях рассматриваемых алгоритмов отжига, причем значения параметра T интерпретируются как начальная температура. Отмечено, что значительные усилия были потрачены на эффективный выбор параметров алгоритма отжига.

При построении нашего поиска использовалась метрика рш (x, y) = max | xn - yn |. Число шагов на эта-

1<n<d

пе поиска равнялось 5. Радиусы шаров a}- образовывали геометрическую прогрессию. Начальный радиус a1 во всех случаях равнялся 1. Конечные радиусы ak

равнялись 0,01, 10-18 и 10-120 соответственно. Конечные радиусы выбирались близкими к ожидаемой точности поиска при аппроксимации по аргументу. Значения в j также образовывали геометрическую прогрессию и выбирались достаточно большими. Началь-

m

m

ное значение Р1 во всех случаях равнялось 50. Конечные значения вк равнялись 105, 1037 и 10240 соответственно. Конечные значения в к выбирались близкими к ожидаемой точности поиска при аппроксимации по значению. Здесь для оценки минимального значения / поиск повторялся 500 раз (число повторений было увеличено для повышения точности оценки).

Оценки минимального значения целевой функции / при N вычислений (шагов поиска)

N T = 0,01 T = 0,1 T = 1 Алгоритм 1

100 1,86 0,091 0,763 1,6 x 10-5

1000 0,0092 0,067 0,506 8,2 x 10-38

10000 0,00038 0,0024 0,018 7,5 x 10-237

Хорошо видно, что поиск данной работы оказался значительно точнее алгоритма отжига с параметрами работы [5]. Причем при больших значениях N преимущество огромно. Кроме того, никаких особых усилий по подбору параметров исследуемого поиска не прилагалось.

Отметим в заключение, что процедура моделирования предложенного варианта поиска очень проста. Эта процедура проще процедуры моделирования алгоритма отжига [5], основанного на моделировании многомерного нормального распределения, и лишь немного сложнее процедуры моделирования простейшего «чистого случайного поиска».

1. Ермаков С.М., Жиглявский А.А. // Теория вероятностей и ее применения. 1983. №1. С.129-136.

2. Ермаков С.М., Жиглявский А.А., Кондратович М.В. // Журн. вычисл. математики и мат. физики. 1989. Т.29. №2. С.163-170.

3. Жиглявский А.А., Жилинскас А.Г. Методы поиска глобального экстремума. М.: Наука, 1991. 248 с.

4. Zhigljavsky A., Zilinskas A. Stochastic Global Optimization. Berlin: Springer-Verlag, 2008. 262 s.

5. Spall J.C. Introduction to stochastic search and optimization: estimation, simulation, and control. New Jersey: Wiley, 2003. 618 p.

6. Spall J.C., Hill S.D., Stark D.R. Theoretical framework for comparing several stochastic optimization approaches // Probabilistic and randomized methods for design under uncertainty. London: Springer, 2006. P.99-117.

7. Yin G. // SIAM J. on Optimization. 1999. V.10. No.1. P.99-120.

8. Ingber L. // Math. Comput. Modelling. 1989. V.12. P.967-973.

9. Лопатин А.С. Метод отжига // Стохастическая оптимизация в информатике. 2005. Вып.1. С.133-149.

10. Граничин О.Н., Поляк Б.Т. Рандомизированные алгоритмы оценивания и оптимизации при почти произвольных помехах. М.: Наука, 2003. 291 с.

11. Абакаров А.Ш., Сушков Ю.А. Статистическое исследование случайного поиска // Математические модели. Теория и приложения. Вып.2. СПб.: Изд-во НИИХ СПбГУ, 2002. C.70-86.

12. Nekrutkin V.V., Tikhomirov A.S. // Acta Applicandae Mathematicae. 1993. V.33. P.89-108.

13. Тихомиров А.С., Некруткин В.В., Марковский монотонный поиск экстремума. Обзор некоторых теоретических результатов // Математические модели. Теория и приложения. Вып.4. СПб.: ВВМ, 2004. С.3-47.

14. Тихомиров А.С. О быстрых вариантах алгоритма отжига (simulated annealing) // Стохастическая оптимизация в информатике. 2009. Вып.5. С.65-90.

15. Тихомиров А.С. // Журн. вычисл. математики и мат. физики. 2010. Т.50. №1. С.24-37.

16. Жиглявский А.А., Жилинскас А.Г. Указ. соч. С.11.

17. Zhigljavsky A., Zilinskas A. Op.cit. Р.38.

18. Ibid. P.119.