ЭВОЛЮЦИОННЫЕ АЛГОРИТМЫ ДЛЯ РЕШЕНИЯ ЗАДАЧИ ОПТИМАЛЬНОГО УПРАВЛЕНИЯ

Полный текст

Численные методы решения задачи оптимального управления в основном заключаются в редукции задачи к задаче нелинейного программирования и использовании различных градиентных методов для поиска оптимального решения [1]. Следует отметить, что градиентные методы находят точки с нулевым значением градиента. Этому требованию удовлетворяют не только точки локального минимума, но и седловые точки, которых у произвольной функции большой размерности может быть больше, чем точек локального минимума. Например, если рассматривать произвольную функцию в пространстве размерности p, то для точки с нулевым значением градиента может быть 2 p вариантов знаков для вторых производных, среди которых только один вариант, когда все вторые производные больше нуля, т.е. данная точка есть локальный минимум, один вариант, когда все производные меньше нуля, т.е. точка - локальный максимум и 2 p - 2 вариантов, когда данная точка является седлом. Градиентные методы как правило эффективно работают для выпуклых функций, у которых отсутствуют седловые точки и локальный минимум совпадает с глобальным. В задаче оптимального управления, приведенной к задаче нелинейного программирования, достаточно сложно определить свойства функционала и гарантировать его выпуклость даже в ограниченной области пространства поиска. Современные поисковые эволюционные алгоритмы [2] не требуют от целевой функции специальных свойств, и находят точки локального минимума, причем в процессе поиска могут находить несколько точек локального минимума и в итоге определяют среди них наилучшую. До недавнего времени единственным препятствием применения эволюционных алгоритмов для решения многомерных задач оптимизации была слабость вычислительной техники. Современные компьютеры позволяют эффективно использовать эволюционные алгоритмы для решения задачи оптимального управления. Цель данной работы заключается в подтверждении высокой эффективности эволюционных алгоритмов в задачах поиска оптимального управления. В работе рассматривается решение задачи оптимального управления популярными эволюционными алгоритмами: генетическим алгоритмом, алгоритмом дифференциальной эволюции, методом роя частиц и алгоритмом летучих мышей. Для сравнения решаем также эту же задачу алгоритмами случайного поиска и градиентного наискорейшего спуска. Результаты сравниваем по средним значениям функционалов. Генетический алгоритм Он был разработан в 80-х годах прошлого столетия [3]. В генетическом алгоритме каждую компоненту вектора параметров q = [q 1 …q 8 ] Т кодируем двоичным кодом Грея. Для кодирования используем c бит под целую часть числа и d бит под дробную часть числа. Всего код одного вектора из восьми параметров содержит 8(c + d) бит. В результате код каждой компоненты вектора параметров определяет положительное число g i , i = 1, …, 8, в диапазоне от 0 до 2 c . В алгоритме первоначально генерируем случайно H двоичных кодов из 8(c + d) бит. Предполагаем, что сгенерированные коды представляют собой коды Грея Далее вычисляем значения функционалов для каждого кода возможного решения. Для этой цели первоначально переводим код Грея каждого возможного решения S i = (s i 1 , …, s i 8(c+d) ) в двоичный код B i = (b i 1 , …, b i 8(c+d) ) по формуле Переводим двоичный код в десятичный для каждых (c + d) бит Вычисляем значения вектора параметров q i = [q i 1 …q i 8 ] Т возможного решения i по формуле Формируем множество оценок целевых функций В алгоритме выполняем одноточечное скрещивание. Отбираем случайно два кода возможных решений Вычисляем вероятность скрещивания по формуле Находим случайно точку скрещивания σ ∈ {1, …, 8(c + d)} и выполняем скрещивание. Получаем два новых кода возможных решений Выполняем операцию мутации с заданной величиной вероятности p m . Для кода нового возможного решения S H+1 находим случайно точку мутации μ ∈ {1, …, 8(c + d)} и меняем случайно компоненту S μ H+1 ∈ {0, 1}. Тоже повторяем для второго нового возможного решения S H+2 . Оцениваем новые возможные решения по значению минимизируемого функционала, получаем f H+1 , f H+2 . Если оценка нового возможного решения лучше наихудшей оценки во множестве возможных решений, то заменяем наихудшее решение на новое возможное решение. S w ← S H+i , f w ← f H+i если f w > f H+i , где f w = max{f 1 , …, f H }, i = 1, 2. После выполнения K операций скрещивания алгоритм останавливаем и находим наилучшее возможное решение во множестве. Алгоритм дифференциальной эволюции Алгоритм дифференциальной эволюции впервые был представлен в 1995 году авторами К. Прайсом (K. Price) и Р. Сторном (R. Storn) [4]. В алгоритме исполь- зуется адаптивная мутация особей, которая зависит от распределенности особей популяции в пространстве поиска. При сильной распределености особей, алгоритм производит существенные вариации, а при малой распределенности производятся лишь небольшие изменения особей популяции. Классическая схема алгоритма дифференциальной эволюции выглядит следующим образом: если доступная информация о системе не позволяет выбрать начальное приближенное решение, то начальная популяция s 0 заполняется случайным образом. Далее на каждой итерации в текущей популяции s k для каждой особи s создается особь-потомок s′ на основе трех родительских особей s 1 , s 2 , s 3 , выбранных из этой популяции случайным образом (s 1  s 2  s 3 ). Для этого хромосому q′ формируют путем линейной комбинации хромосом q 1 , q 2 , q 3 , которые являются вещественными векторами: где F - скалярный коэффициент влияния. После формирования особей s′ проверяется их приспособленность и, если она выше чем у родительской особи s, то новая особь заменяет родительскую в популяции: При формировании новой хромосомы q′ может возникнуть ситуация, когда один или несколько ее компонентов (генов) окажутся за пределами допустимых значений (q i - ; q i + ). В этом случае вычисляется новое случайное значение q, которым заменяется компонент в хромосоме. Также существует модификация алгоритма дифференциальной эволюции, в котором в качестве базовой особи для мутации используется не случайно выбранная из текущей популяции, а наиболее приспособленная особь s best . Тогда формирование новой хромосомы q′ производится по правилу: где λ - дополнительный коэффициент для повышения эффективности использования наиболее приспособленной особи. Метод роя частиц Метод роя частиц основан на имитации социально-поведенческих моделей организованных групп [5]. Метод строится на тех же принципах, на основе которых птицы в стае и рыбы в косяке ведут себя удивительно синхронно, двигаясь в том или ином направлении как единое целое. Авторами метода являются Дж. Кеннеди (J. Kennedy) и Р. Эберхарт (R. Eberhart), которые предложили идею канонического метода роя частиц в 1995 году. В методе роя частиц распределенные в пространстве параметров задачи оптимизации частицы играют роль агентов. Частицы перемещаются в пространстве параметров и изменяют направление и скорость своего движения на основе определенных правил. На каждой итерации вычисляется соответствующее положению частицы значение целевой функции. Также частице доступна информация о наилучшей частице из числа соседних с ней частиц и информация о лучшем собственном значении на предыдущих итерациях. На основе этой информации определяется правило изменения положения и скорости частицы в пространстве параметров. Новое положение частицы s i в момент времени t определяется вектором ее координат q i , а ее скорость - вектором v i : где P n [a; b] - n-мерный вектор псевдослучайных значений, равномерно распределенных на интервале [a; b]; q b i,t - вектор координат частицы с лучшим собственным значением на предыдущих итерациях; q g,t - вектор координат лучшей соседней частицы; α, β, γ - свободные параметры алгоритма со следующими рекомендуемыми значениями: α = 0,7298; β = γ = 1,4962. Алгоритм летучих мышей Алгоритм разработан в 2010 году [6; 7]. Летучие мыши обладают уникальными средствами эхолокации, которая используется для обеспечения полетов в темноте и обнаружения добычи. В процессе поиска добычи летучие мыши генерируют сигналы, имеющие частоту ω i и громкость a i . Частота сигналов и их громкость может изменяться в заданных пределах - [ω min ; ω max ] и [0; 1] соответственно. Также может изменяться частота повторения сигналов r i ∈ [0; 1]. На первом шаге метода случайным распределением задаются соответственно начальные значения q i,0 , скорость v i,0 , частота сигнала ω i,0 , громкость сигнала a i,0 и частота повторения сигнала r i,0 . Определяется вектор координат X b агента, доставляющего глобально лучшее решение. Выполняется миграционная процедура по следующей схеме: где P 1 [-1; 1] - случайное число, равномерно распределенное на интервале [-1; 1]. Далее реализуется случайный локальный поиск по следующей схеме: где P n [-1; 1] - n-мерный вектор псевдослучайных значений, равномерно распределенных на интервале [-1; 1]; a - текущее среднее значение громкостей всех агентов популяции: Схема локального поиска повторяется либо до получения вектора q′ i,t , при котором значение фитнес-функции φ(q′ i,t ) лучше, чем при q i,t , либо до достижения терминального числа повторов схемы. На последнем этапе реализуется эволюция параметров a i и r i по правилу: где b a ∈ (0; 1), b r > 0 - свободные параметры алгоритма, рекомендуемые значения которых равны 0,9. Алгоритм случайного поиска Генерируем случайный вектор параметров q˜ = [q˜ 1 …q˜ 8 ] Т по формуле где ξ - случайная величина, ξ ∈ [0; 1]. Вычисляем оценку вектора f ˜ по значению минимизируемого функционала. Далее генерируем новый вектор q = [q 1 …q 8 ] Т по формуле (1) и вычисляем его оценку f. Если оценка нового вектора лучше f < f ˜ , то заменяем вектор и его оценку, q˜ ← q, f ˜ ← f. Повторяем генерацию вектора и вычисление его оценки R раз. В результате получаем наилучший вектор параметров q˜ = [q˜ 1 …q˜ 8 ] Т . Алгоритм наискорейшего градиентного спуска В соответствии со стратегией алгоритма наискорейшего градиентного спуска [8] генерируем случайный вектор параметров q˜ = [q˜ 1 …q˜ 8 ] Т по формуле (1). Вычисляем значение оценки для полученного вектора параметров Вычисляем вектор градиента целевой функции в точк где δq - заданная малая положительная величина Вычисляем значение второго вектора параметров q = [q 1 …q 8 ] Т в направлении антиградиента где первоначально M = M 0 , M 0 - заданное положительное целое число. Проверяем ограничения для второго вектора параметров q = [q 1 …q 8 ] Т . Если где δm - заданное положительное целое число, δm << M 0 . Если M > 0, то пересчитываем значение второго вектора параметров по формуле (2), иначе завершаем вычисления и считаем решением вектор параметров q˜ = [q˜ 1 …q˜ 8 ] Т . Повторяем вычисления по формулам (2), (4), пока выполняются условия (3). В результате получаем второй вектор параметров q = [q 1 …q 8 ] Т . Вычисляем значение минимизируемого функционала для второго вектора параметров Вычисляем значения векторов в точках золотого сечения При нарушении условия (10) определяем первый вектор параметров. Если f α < f β , то q˜ = q α , f ˜ = f α , иначе q˜ = q β , f ˜ = f β . Повторяем вычисления по формулам (5)-(10). Всего выполняем все вычисления заданное количество K 1 раз. Результатом вычислений будет последний вектор параметров q˜ = [q˜ 1 …q˜ 8 ] Т . Задача параметрического оптимального управления мобильным роботом Задана математическая модель мобильного робота: Заданы ограничения на управление: Заданы начальные условия: Заданы терминальные условия: Задан функционал Необходимо найти управление в классе кубических полиномов где значения коэффициентов полиномов ограничены Найденные оптимальные коэффициенты кубических полиномов должны обеспечивать минимальное значение заданному функционалу Сравнительный эксперимент Для проведения сравнительного эксперимента для каждого возможного решения формировались полиномы для управления в виде (16), (17), интегрировалась система дифференциальных уравнений (11)-(13) из заданных начальных условий (14) и вычислялась оценка по значению функционала (15). В вычислительном эксперименте использовались следующие параметры модели (11)-(13): x 1 0 = 10, x 2 0 = 10, x 3 0 = 0, x 1 f = 0, x 2 f = 0, x 3 f = 0, u 1 - = -10, u 2 - = -10, u 1 + = 10, u 2 + = 10, t + = 4, ε = 0,02, q i - = -8, q i + = 8, i = 1, …, 8. В генетическом алгоритме использовали следующие значения параметров: H = 512, K = 2 14 = 16384, c = 4, d = 12, p m = 0,7. В алгоритме случайного поиска количество сгенерированных векторов параметров было R = 40000. В алгоритме наискорейшего градиентного спуска использовали следующие параметры: q = 0,0001, M 0 = 4096, δm = 8, K 1 = 128. В методе роя частиц размер популяции составлял H = 16, число итераций - W = 2048. В методе, инспирированным летучими мышами, размер популяции составлял H = 32, число итераций - W = 1024. В алгоритме дифференциальной эволюции размер популяции составлял H = 64, число итераций - W = 1024, число скрещиваний на каждой итерации - S = 32. Для каждого алгоритма проводилось по 10 испытаний. Вычисления выполняли на компьютере с процессором Intel ® Core™, i7-2640M, CPU@2.80 GHz. Результаты вычислений, представляющие собой продолжительность перехода из начального в конечное положение для каждого алгоритма (таблица) показали, что все эволюционные алгоритмы решают рассмотренную задачу оптимального управления приблизительно одинаково и лучше, чем алгоритм случайного поиска. Алгоритм наискорейшего градиентного спуска при решении данной задачи работает неустойчиво и не находит в большинстве случаев приемлемого решения. Среди эволюционных алгоритмов наилучшие результаты показа генетический алгоритм. Таблица Результаты вычислительного эксперимента [The results of computational experiments] Испытание GA DE PSO BIA RS FGD 1 4,9349 6,5507 4,8288 4,3877 5,9813 9,7561 2 4,1093 4,6874 4,8866 4,7422 6,1006 5,4075 3 4,6475 4,7433 4,4649 5,3366 5,2015 18,2190 4 4,4244 5,8050 4,6277 4,6821 5,7782 11,3907 5 4,2987 5,8152 4,4758 4,4747 4,6759 22,5127 6 4,0434 6,4671 4,7797 4,1769 5,2597 28,3650 7 4,1218 4,6943 4,8279 4,5676 5,0867 15,1068 8 4,6303 4,9991 4,9260 5,4031 5,0530 17,1744 9 4,2148 4,6895 4,7935 4,3658 5,9370 30,7676 10 4,4533 5,0282 4,6172 5,5366 5,1588 13,9237 Сред. 4,3878 5,3480 4,7228 4,7673 5,4233 17,2624

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