Многопараметрическая оптимизация многомодальных функций

УДК 621.438

С. В. Черных

МНОГОПАРАМЕТРИЧЕСКАЯ ОПТИМИЗАЦИЯ МНОГОМОДАЛЬНЫХ ФУНКЦИЙ

Описан ход разработки программного средства для многопараметрической оптимизации на базе генетических алгоритмов. Дается краткое описание предметной области и объясняется особенность поставленной задачи. Указаны основные этапы разработки и достигнутый на сегодня уровень исполнения и структура программы. Приводятся примеры расчетов.

The course of development of a software for multiparametrical optimization on the basis of genetic algorithms is described. A brief description of a subject domain is given and a feature of the problem is explained. The basic development cycles and achieved on a today level of execution and structure of the program are specified. Examples of accounts are resulted.

Ключевые слова: генетические алгоритмы, многомодальные функции, модель эволюции, функция полезности, фронт Парето, программная реализация.

Key words: genetic algorithms, multimodal functions, model of evolution, function of fitness, front Pareto, program realization.

Эволюционное моделирование (ЭМ) зародилось в 60—70-х гг. прошлого века в результате работ трех групп исследователей. Поскольку применяемые методы заимствовали некоторые природные процессы, ЭМ было отнесено к бионическому направлению искусственного интеллекта. Достаточно полное и вместе с тем компактное описание четырех основных составляющих ЭМ дано в работе [1].

Из значительного объема материала, степени детализации и большого количества вариантов, посвященных генетическим алгоритмам (ГА), видно, что именно эта составляющая ЭМ является на сегодня наиболее развитой и широко применяемой на практике (по сравнению с эволюционным программированием, эволюционными стратегиями и генетическим программированием), что подтверждается общим обзором литературы за последние 10 лет. В РГУ им. И. Канта возникло объяснение такой популярности ГА, как следствие наиболее близкого к человеку соотношения применения операторов кроссинговера и мутации при образовании потомков.

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

рисунке 1, бросают горсть песчинок, которые случайным образом распределяются по ее поверхности и приклеиваются, образуя начальную популяцию. На основе значений их координат Х и У, переведенных в двоичный код, для каждой особи формируется математическая хромосома. (В реальности именно из разброса данных хромосом формируется начальная популяция.)_____________________________________________

1. Применение эволюционного моделирования для решения задач многопараметрической оптимизации

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

Рис. 1. Пример многомодальной функции

Далее представлен способ описания сути работы ГА для определения глобального максимума такой функции. При этом вводятся понятия основных генетических операторов, необходимых для понимания дальнейшего материала.

Сначала на ландшафт функции, изображенной на

Вестник Российского государственного университета им. И. Канта. 2010. Вып. 10. С. 94—103.

По их координатам Ъ для каждой особи вычисляется значение' целевой функции (ЦФ): чем больше Ъ, тем больше ЦФ. Составляется таблица значений ЦФ особей в порядке ее возрастания— ранжирование. Особи, чьи значения ЦФ стоят в верхней половине таблицы, допускаются к размножению — оператор селекции (ОС).

Размножение происходит виртуальным дистанционно-информационным образом.

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

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

Повторяется ранжирование, и нижняя половина популяции

(с худшими значениями ЦФ) уничтожается — оператор редукции (ОР). После него размер популяции возвращается к исходному.

Далее опять следуют ОС, ОК, ОМ и ОР. Каждое новое поколение песчинок, состоящее из некоторого количества предков и некоторого числа потомков, оказывается в среднем размещено на ландшафте функции выше предыдущего. То есть песчинки от поколения к поколению оказываются (не перемещаясь физически) все ближе к максимумам функции — локальным и глобальному. Причем, поскольку мутация перемещает особи еще и по Х и У, то в конце концов все особи оказываются около глобального максимума. Когда одна первая особь достигает его (то есть максимальное значение ЦФ по популяции перестает расти), считается, что задача решена. Описанный процесс реализует только дарвиновскую модель эволюции; существует много разновидностей моделей ОС, ОК, ОМ, ОР.

Таким образом, ГА работает с популяцией, применяет операторы ОС, ОК, ОМ и ОР [1], что можно выразить соотношением

ГА = ГА <ОК, ОМ, ОР, ОС>

Наиболее часто ГА применяется для МР-комбинаторных задач: задачи коммивояжера, оптимального размещения изделий на европоддоне, раскладки заготовок на прессе, распиловки бревен и раскроя заготовок, оптимизации портфеля заказов и управления ресурсами предприятия, управления многопродуктовыми складами и торговыми сетями, например [2], подбору характеристик конструкционных деталей или элементов изделий. Самая известная научная школа в России по ГА — профессора В. М. Курейчика из Таганрога, основное направление которой — применение ГА в радиоэлектронике для трассировки электронных плат и проектирования СБИС.

Возможно применение ГА и в мультидисциплинарных предметных областях, которое и будет рассматриваться дальше. В частности, в [3] описана постановка задачи оптимального проектирования двигательного блока системы коррекции, ориентации и стабилизации искусственного спутника Земли. Она сводится к нахождению максимума многомерной целевой функции Я:

К = я(Рк, рь, Ты, 5^, ер),

где ¥К — критерии (требования) по функциональности; Рь — критерии по прочности; Ты — критерии по тепловому режиму; — требования по радиационной стойкости; ЕР — критерии по

конфигурации магнитного поля (если в двигательном блоке применяются стационарные плазменные двигатели СПД). В этом случае мы имеем дело и с мультидисциплинарной, и с многопараметрической оптимизацией.

2. Математическая постановка задачи оптимизации

Пусть Х — и-мерный вектор (варьируемых) параметров задачи. Будем полагать, что X с Яп, где Яп, п > 1, — и-мерное арифметическое пространство (пространство параметров).

Параллелепипедом допустимых значений вектора параметров назовем непустой параллелепипед

П = {X | хі є [х+, х+ ], і є [1, и]}.

--------------------------------------------------------------------W-------------------

На вектор X могут быть дополнительно наложеної к ограничении, составляющих множество

D = {XI gi(x)є [0, ^ г є [1, к]},

где gi (x) — ограничивающие функции. Известно, что если все функции gi (x) непрерывны, то множество D является замкнутым. Замкнутое множество DX = ПоD называется множеством допустимых значений вектора параметров.

Пусть Ф(X) = (Фj(X), Ф2(X),..., Фп(X)) — векторный критерии оптимальности, онределенный

на множестве П. Будем полагать, что 0(X) є Rm, где Rm, m > 1, — m-мерное арифметическое пространство (пространство критериев).

На частные критерии оптимальности Ф1 (X), Ф2 (X),..., Фп (X) могут быть наложены критериальные ограничения вида

ф: Ф,

где <1*, Ф( — некоторые вещественньге константы.

В общем случае лицу, принимающему решение (ЛПР), может быть, желательно уменьшить значения одних частных критериев оптимальности и увеличить значения других. Однако, поскольку задача максимизации критерия легко сводится к задаче минимизации обратного критерия, будем полагать, что желательна минимизация всех частных критериев оптимальности Ф^), Ф2(X),..., Ф„(X).

Таким образом, задача многокритериальной оптимизации формулируется следующим образом: при условии (2) найти min Ф^) = Ф^*), где вектор X* — решение задачи многокритериальной

X^x

оптимизации.

Для множества, в которое векторный критерий оптимальности Ф^) отображает множество DX, введем обозначение D0 и назовем его критериальным множеством.

Обозначим Y = Ф(Х), и пусть \л < Y2, если Y1 Ф Y2, у] $ у], і є [1, т\. Будем говорить, что вектор Y1 из критериального множества доминирует (по Парето) над вектором Y2 из того же множества,

если Y1 Y2.

Неформально множество Парето DJ, поставленной задачи многокритериальной оптимизации можно определить как совокупность векторов Y є D0, среди которых нет доминируемых. Формально множество Парето определяется так: DJ, = {Y* є Оф | {У є Оф | Y' Y*} = 0}.

Если Ф(X) є Dф, то говорят, что вектор X — эффективный по Парето. Множество векторов, которое порождает множество Парето, будем обозначать D*X и называть эффективным по Парето. Формально можно записать Dф =Ф(DX). Отметим, что в системе «Парето» параллелепипед П имеет в основном «технологический» смысл: определяет диапазоны изменения случайных чисел, формирующих расчетную сетку во множестве DX, а также обеспечивает ограниченность множества DX.

Множество D может принадлежать одному из следующих иерархически вложенных «классов множеств»: многосвязные множества D1; односвязные множества D2; невыпуклые множества D3; выпуклые множества D4; множества D5 заданные нелинейными ограничениями; множества D6, заданные линейными ограничениями (произвольный выпуклый многогранник); регулярное множество D7.

К классу множеств D71 относятся параллелепипеды с ребрами, параллельными координатным осям пространства параметров. Для определения параллелепипеда, принадлежащего этому множеству, должно быть задано к = 2п констант a, bi таких, что D = {X | хг є [аг, Ьг ], і є [1, 2n]}.

К классу множеств D7 2 принадлежат симплексы, ортогональные ребра которых параллельны

координатным осям пространства параметров. Для определения симплекса, принадлежащего этому множеству, должна быть задана к = п + 1 константа аг, і є [1, п + 1], и h такие, что

D = {X | xt -at є [0, да), (x1 -а1) + (х2 -а2) +... + (хп -ап) є (-да, h], і є [1, n]}.

--------------------------------------------4------------------Ш-----------------------

Введенная иерархия классов множеств обладает тем свойством, что если О с Оі, і є [1, 6], то, вообще говоря, Ох = П п О с О. Класс множеств О7, легко видеть, этим свойством не обладает: если множество О регулярно, то пересечение его с параллелепипедом П может быть нерегулярным множеством.

Классы множеств имеют так же, как параллелепипед П, «технологический» смысл. Например, если некоторая программа ориентирована на использование в качестве множества Ох множества, принадлежащего классу множеств Оі, то она не может быть использована с классом множеств ОІ+х, ОІ+2 и так далее. Аналогично, если программа ориентирована на множество Ох класса О7і, то она не может быть использована с множеством Ох класса О7), ] ф і. Заметим, что класс множества, которому принадлежит множество Ох, заведомо может быть неизвестен. В этом случае этот класс может быть уточнен в процессе решения задачи многокритериальной оптимизации.

Итак, при оптимизации существуют два подхода: подход ЦФ (функции полезности, функции пригодности, функции фитнеса), уже упоминавшейся выше, и подход фронта Парето.

Построение фронта Парето в результате работы ГА объясняется, например, в [4]. Приводится пример наглядно-осмысляемого двумерного фронта. Этот подход хорош тем, что ЛПР видит расширенную информацию — форму фронта в конкретном числовом диапазоне изменения параметров. Анализ изображений всех комбинаций двумерных фронтов Парето для многокритериальной задачи позволяет ЛПР понять полную картину взаимозависимости параметров и визуально определить главный максимум. Либо использовать эту информацию как предварительную для построения функции полезности и дальнейших исследований уже с ее помощью.

Подход ЦФ хорош тем, что дает однозначное определение оптимального решения как глобального максимума этой функции. С другой стороны, структура такой функции может быть достаточно сложна, и не очень понятно, из каких соображений ее конструировать. Пример такой сложной ЦФ представлен в работе [5]. В ней рассматривается проектирование многоэлементных антенн Яги-Уда, известных также под названием волновой канал (тех, что стоят на крышах наших домов). Использовалось 15-мерное проектное пространство, а хромосома восьмиэлементной антенны имела вид

Ь1О1І2О2ЬзОзЬ4О4Ь5О5ЬбОбЬ7О7Ь8, где Ьі — длина элемента антенны с номером І; Оі — расстояние между элементами антенны с номерами і и і + 1.

На кодирование каждого элемента хромосомы отводилось 6 бит информации, общая длина хромосомы составила 90 бит. При этом ЦФ этой антенны и нескольких иных приведенных в статье представляли собой сложные формулы, включающие такие переменные, как коэффициент усиления антенны, уровень боковых лепестков, рабочую полосу частот, коэффициент помехозащищенности, коэффициент стоячей волны в линии передачи, активную и реактивную составляющие входного сопротивления антенны, волновое сопротивление линии питания и сумму членов функции диаграммы направленности. Эти переменные составляли различные комбинации в виде сумм, разностей, дробей, знаков модуля с числовыми коэффициентами и штрафными коэффициентами, высчитываемыми по сложным формулам.

В ходе разработки изначально поставленная задача создания надстройки над САЕ-системой №БА была скорректирована в сторону создания программного средства (ПС), универсального для практического применения как шодеРгопНег, с одной стороны (то есть стыкуемого с основными САЕ-пакетами), причем коммерческой версии доступной всем. С другой стороны, ПС должно давать возможность реализовывать не только различные варианты возможных значений генетических операторов, но и разные модели ГА, а также последние изыскания в данной научной области типа нечетких ГА и недарвиновских моделей эволюции.

3. Тестирование созданного программного средства на модельной оптимизационной задаче

Обычно эффективность работы систем с ГА определяют на функциях типа Расстригина и подобных показанной на рисунке 1.

В этом случае в качестве эталонного средства использовалась система шодеБгопЙег [6], а в качестве реальной физической задачи — оптимизация геометрических характеристик

коромысла механизма газораспределения двигателя. Модель сшдавалась в системе АЫБУБ МиШрЬу81с8 на основе чертежа реальной детали конструкции (рис. 2).

Рис. 2. Чертеж коромысла газораспределительного механизма

Материал детали — сталь 45ХН со следующими характеристиками: модуль Юнга Е равен 200000 МПа, коэффициент Пуассона равен 0,26. Рассматривается напряженно-деформированное состояние детали; в качестве параметра оптимизации принимается максимальное напряжение по конструкции.

После создания конечно-элементной модели было произведено ее закрепление жестко по внутреннему краю втулки. Давления P1 и Р2 действующие на плечи коромысла, рассчитаны таким образом, чтобы уравновешивать друг друга (рис. 2). Давление, приложенное к носку коромысла (со

стороны клапана), Р2 = 48,0769 МПа. Давление, приложенное к хвостовику коромысла (со стороны штанги), Р1 = 25 МПа.

Далее из меню возможных вариантов расчетов ANSYS выбиралось вычисление эквивалентных напряжений по теории Von Mises^. Результаты расчета до проведения оптимизации приведены на рисунке 3. Область с максимальным эквивалентным напряжением находится в районе хвостовика коромысла и указана длинной белой стрелкой. Значение давления в ней достигает 359 МПа.

В ходе оптимизации ищется вариант конструкции с минимальным давлением. При этом параметры коромысла (рис. 2) меняются в пределах, указанных в таблице 1.

При этом на каждом шаге по одному из параметров программа оптимизации вызывает входной файл ANSYS, его модифицирует и запускает новый расчет.

Таблица 1

Задание пределов геометрических параметров, мм

Переменная L1 L2 T1 T2 T3

Номинал 20 10 20 16 8

Верхний предел 25 15 20 18 8

Нижний предел 15 8 10 8 6

В результате оптимизационного поиска приблизительно через 5 минут и 100 поколений находится особь, дающая уменьшение давления в проблемном месте до 138 МПа (рис. 4). Эта особь имеет параметры, приведенные в таблице 2.

Рис. 4. Распределение напряжений

В шодеБгопНег использовался МОСА-11 с популяцией в 10 особей, 10 поколений. Вероятность направленного кроссовера 0,5, вероятность селекции 0,05 и вероятность мутации 0,1. В расчетах с помощью нашего программного средства тоже было 10 особей, 10 поколений, обычный ГА с двухточечным кроссинговером. Вероятности те же, что и при расчете с помощью шодеБгопИега.

Таблица 2

Задание пределов геометрических параметров, мм

Переменная 11 12 71 72 73

Номинал 20 10 20 16 8

Верхний предел 25 15 20 18 8

Нижний предел 15 8 10 8 6

Оптимизированное шodeFrontiera 15 11 10 13 8

Оптимизированное нашей программой 16,6 9 10,2 16,2 7,9

Результаты расчетов (уменьшение напряжения от количества поколений) представлены на рисунках 5 и 6. На графике, выдаваемом тодеБгопЙег, показано среднее значение по популяции, а на нашем — наилучшее. Поэтому нижний график более монотонный. Видно, что для такого простого (с точки зрения поверхности отклика) случая однодисциплинарной, однокритериальной и пятипараметрической оптимизации скорость расчета с помощью алгоритма МОСА-ІІ, имеющего специализированный направленный кроссовер, фактически не отличается от расчета с помощью простого двухточечного кроссовера.

Рис. 5. Уменьшение максимального напряжения по конструкции в течение первых 100 поколений при оптимизации шodeFrontiera

Рис. 6. Уменьшение максимального напряжения по конструкции в течение первых 100 поколений при оптимизации нашей программой

Заключение

Итак, разрабатываемое ПС, включающее важные свойства инженерных пакетов (в том числе стыковка с интерфейсными файлами САЕ-систем №БА, АМБУБ) и имеющее научно-исследовательские свойства (возможность варьировать разновидности ГА, различные модели эволюции и создавать гибридизацию с нечеткими системами и нейронными сетями) создано в первом приближении и прошло апробацию на ряде тестовых задач и в сравнении с ПС мирового уровня.

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

1. Методы генетического поиска. Таганрог, 2002.

2. Черных С. В. Решение задачи управления товарными запасами с применением генетических алгоритмов в среде Matlab / / Сборник статей III Всерос. науч.-техн. конф. Пенза, 2005.

3. Толетель О. В. Эволюционное моделирование двигательных блоков космических аппаратов // Междунар. науч.-техн. конф. «Интеллектуальные системы» (AIS'04) и «Интеллектуальные САПР» (CAD-04). М., 2004. Т. 1. С. 86-90.

4. Колесников А. В. Гибридные интеллектуальные системы. Теория и технология разработки. СПб., 2001.

5. Зинченко Л. А. Олейник М. П., Сорокин С. Н. Эволюционное проектирование антенн Яги-Уда с улучшенными характеристиками // Междунар. науч.-техн. конф. «Интеллектуальные системы» (AIS'03) и «Интеллектуальные САПР» (CAD-03). М., 2003.

6. Esteeo. URL: www. esteco. com.

Об авторе

Сергей Владимирович Черных — ассист., РГУ им. И. Канта, chemych@inbox.ru.

Author

Sergey Chernikh — assistant, IKSUR, e-mail: chernych@inbox.ru.