Оптимизация формы множества Парето в задачах многокритериального программирования

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

Конечномерной многокритериальной моделью мы будем называть математическую модель с N целевыми функциями fk(х,и) ^ max     k = [1,N] ,                            (1)

Ж подлежащими максимизации, на. обладающем внутренними точками множестве элементов х € Еп, удовлетворяющих условиям вида

Уі(х,и) 6 0      г = [1,т],                                  (2)

где и € Ө С Е^Т - вектор параметров модели. При этом предполагается, что функции f_k(х, и) и Уі(х, и) достаточно гладкие, то есть они имеют непрерывные частные производные требуемого порядка, по всем своим аргументам.

Некорректность в общем случае подобной постановки задачи очевидна, поскольку экстремальный для одной из целевых функций элемент ж, вообще говоря, не является таковым для других. При этом, однако, некоторую полезную информацию может дать последовательное решение следующих однокритериальных задач поиска экстремума на множестве (2) каждой из функций (1) в отдельности для к = [1,2V]:

f k (ж, и) ^ max

X при условиях                                       (3)

Уі(ж,и) 6 0      г = [1,т] .

Целевую функцию f k (ж, и) назовем улучшаемой в допустимой (то есть удовлетворяющей условиям (3)) точке жо, если существует другая допустимая точка жі, для которой fk (жі,и) >  fk (жо,и). Понятно, что решение задачи (3) при любом к = [1, V] является неулучшаемым, или «идеальным» с точки зрения целевой функции fk (ж, и).

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

В дальнейшем будем называть второе множество множеством паретовского типа или просто множеством Парето.

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

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

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

В данной статье мы воспользуемся вариантом подхода, в котором компромиссной является задача поиска на множестве (2) элемента, минимизирующего величину рассогласования целевых функций (то есть максимального по всей совокупности критериев, отклонения каждого из критериев от своего неулучшаемого «идеала»). Иными словами, это есть задача математического программирования следующего вида:

минимизировать р при условиях р >  0:

Уі(ж,и) 6 0 г = [1,т] ,

fk(ж,и) > Ғ*(и) -р к = [1,V] , решение которой мы будем обозначать как р**(и) и ж**(и). Здесь FJ*(u) = fk^k(и),и).

Задачу (4) естественно назвать двухуровневой параметрической задачей, поскольку формулировка ее условия содержит F J * ( u), к = [1, V], - решения однокритериальных задач (3), которые будем называть задачами первого уровня. При этом как в задачах первого, так и второго уровня вектор параметров и Е Ө предполагается фиксированным.

Понятно, что экстремальная величина рассогласования критериев в общем случае определяется свойствами множества Парето и является для этого множества некоторой количественной характеристикой, зависящей, вообще говоря, от вектора параметров и. Поэтому естественной представляется постановка для модели (1) - (2) оптимизационной задачи третьего уровня, например, такой:

оптимизировать по и Е Ө выражение р**(и),

решением которой будут являться вектор параметров и*** € Ө и чиело р*** = р** (и***).

Далее в настоящей статье мы рассмотрим возможный подход к решению задачи (5).

Метод решения

Рассмотрим теперь задачу отыскания в пространстве параметров стандартным (например, градиентным) методом поиска экстремума величины рассогласованности значений целевых функций многокритериальной модели (3) - (4) - (5).

Отметим, что специфика этой схемы заключается в том, что постановка задачи (5) верхнего (третьего) уровня включает р**(и) - зависимость, являющуюся решением задачи (4) второго уровня, условие которой, в свою очередь, содержит зависимости Ғ^и) = fk (жк(u),u) Vk € [1,2V], определяемые решениями задач (3) нижнего (первого) уровня.

При этом зависимости р**(и) и Fj(и) Vk € [1, V] в общем случае (для гладких функций fk (ж, и) и уг(ж, и) ) не являются дифференцируемыми функциями, то есть непосредственное использование каких-либо численных методов, основанных на применении тейлоровских аппроксимаций, оказывается невозможным.

Для преодоления этого затруднения предлагается воспользоваться методом гладких штрафных функций (см. [3]), чтобы получить достаточно гладкие аппроксимации зависимостей р**(и) и F^u) Vk € [1, V]. При этом будем предполагать, что штрафная функция Р ( t,s), штрафующая нарушение ограничения S 6 0, удовлетворяет следующим условиям:

• 1°. Функция Р ( t,s) им еет V t >  0 и Vs непрерывные производные по всем своим аргументам до второго порядка включительно.

• 2°. Для всех т >  0 и Vs выполнены неравенства

  д ² Р

  aS² >⁰.

  
  
  
  

ЭР

> 0; ds

3°. Р(t,s) > 0 Vs и любого т > 0, причем, кроме того, lim Р(т,s) = ^+”, S> 0 ’                         (7)

-^+0 v 7     [     0, s < 0.                                 v 7

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

Будем использовать для однокритериальных задач (3) вспомогательную функцию вида

Mw)=f_k (ж, и) - ^р Р (т,уМи) Vk € [1,V ],            (8)

г =1

в то время как достаточно гладкая штрафная функция Р (т, s) удовлетворяет условиям (6) и (7).

Как показано в [4], в качестве сглаженной аппроксимации жк(и) - точного решения каждой из задач нижнего уровня (3), можно принять Жк(и) - стационарную точку вспомогательной функции (8), определяемую системой уравнений дА

дж- = 0    V7 € [1’"]’ или, что то же самое,

df k ^дж з

-

РУ ВРв»

gдр! дЖТ, = ^{0    V 1 €} М.

Поскольку в условие задачи второго уровня (4) входят зависимости Fk(n) = fk (хк(п),п) Vk Е [1,-^], не являющиеся дифференцируемыми функциями своих аргументов, то для этих зависимостей также необходимо выбрать сглаженную аппроксим ацию.

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

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

Поскольку F k (п) = А к (т,Х к (п), п), то по правилу дифференцирования сложной функ-

ции имеем	дFk   дАк     дАк Эх а дп р = дп р + ^ дх , дпр   Vp Е ^^
что в силу (9) дает	дҒк  дАк -— =-— ( т,х к ( п ) ,п )     VP Е [1,т].                       (10) дпр   дпр

Отметим, что последнее упрощение было бы невозможным, если для F k (п) вместо сглаживающей аппроксимации А к (т,Х к (п),п) использовать, более естественную на первый взгляд, аппроксимацию f_k(х_к(п), п).

Рассмотрим теперь процедуру решения задачи второго уровня.

Запишем условие задачи (4) в следующем, более удобном для применения метода штрафных функций, виде:

максимизировать — р при условиях — р 6 0,

У і (х,п) 6 0 г = [1,ш] ,

Үк(р,х,п) 6 0    к = [1,^] ,(11)

где    Ү к (р,х,п) = F^u) — р — f k (х,п).

Решение этой задачи мы будем соответственно обозначать как р**(п) и х**(п).

Определим следующим образом для задачи (10) вспомогательную функцию

Nm

Е(т, р,х,п) = -р -Р (т, -р) — ^Р (т, үк(р,х,п)) — ^2Р (т Уі(х,п)), к=1

заменив предварительно в Үк (р,х,п) зависимость F k (п) на ее сглаженную аппроксимацию F k (п).

Условия стационарности вспомогательной функции (12) по совокупности переменных { р, х і , х 2 , . . . х_п } будут

дЕ	дР N дР
-- =	—1 + — + V — = 0,
др	др к =1 дү к
	(13)
дЕ	N дР дf k m дР ду і
. дх з	^           ^ а а = 0     V3 Е [1’п]‘ к =1 дҮ к дх з   і =1 д У і дх з

Пусть решения системы (13) суть р(п) и х(п), тогда в качестве сглаженной аппроксимации зависимости р**(п) можно использовать функцию Е(п) = —Е (т,р(п) ,х(п) ,и) . Найдем производные этой функции по компонентам вектора параметров п.

По правилу дифференцирования сложной функции многих переменных имеем

что, с учетом выражение

ЭЕ _ ЭЕ   ” ЭЕ Эж , ЭЕ Эр Эпр Эпр + Уж , Эпр + др Эпр    Р^ Е ,Т , ЭЕ       дЕ —— = 0 и —— = 0 VI Е [1,п] из (13), дает существенно более простое др         Эж , эё эе = ^—д,р(п),ж(п),п)      Vp Е [1 ,т ] .                      (14) эпр дпр

Наконец получим формулы для компонент градиента от Е(п) в терминах функций,

используемых в формулировке многокритериальной модели (4) - (5) и методе гладких штрафных функций.

Исходя из (12) получаем

УЕ _  * ЭР ЭҮк _^ ЭРЭр^ дпр        ЭҮк дпр      Эуі дпр , г      к=1             і=1 эҮк  эғ к dfk           эғ к где —— = —--—— , а значения —— находятся из равенств (10).

Эпр   Эпр   Эпр            Эпр

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

Пример решения задачи

Проиллюстрируем особенности такой постановки задачи следующим примером.

Рассмотрим достаточно наглядную многокритериальную математическую модель, в которой ж = ||жіж2жз||^т Е Е ³ - вектор независимых переменных, а вектор параметров п = || піП2^^т Е Е² .

Требуется максимизировать по ж при некотором п Е Ө функции

/1(ж,п) = Ж1 , /2(ж,п) = Ж2 , /з(ж,п) = жз при условиях: ж1 > 0 , ж2 > 0 и жз > 0 ,

О1(пі,П2)жі + О2(пі,П2)ж2 + Оз(пі,П2)жз 6 Ь(пі,П2) , где положительные функции аДп), а2(п), аз(п) и Ь(п) задаются условиями, приведенными далее.

Геометрическая интерпретация многокритериальной модели приведена на рис. 1. Отметим ее основные особенности. Допустимой областью модели (при допустимом фиксированном п) является прямоугольная пирамида ОАВС. Множество Парето совпадает с гранью АВС или является его частью. Предполагается, что точка А имеет координатное представление | п1 0 0 ||T. а точіса В - || 0п2 0 ||T.

При этом будем считать, что множество Ө в пространстве параметров задается условием: сумма длин отрезков ОА, ОВ и ОС постоянна и равна трем.

Применив стандартные методы аналитической геометрии, получим, что для совместности системы ограничений модели необходимо существование г > 0 такого, что ax(ui,U2) = п2Т ,

О2(пі,п2) = п і Т , аз(пі,п2) = ед , Ь(п 1 ,п₂) = п 1 п₂т .

Рис. 1. Геометрическая интерпретация модели (6) - (7)

Причем выполнение условия непустоты множества Ө может быть обеспечено выбором г = 3 — иі — и2

при ограничениях 0.1 6 и і 6 2.5 и 0.1 6 и 2 6 2.5.

Минимальная величина рассогласованности критериев в данном примере зависит от формы множества Парето, которое является треугольником АВС либо его частью. Графическое представление зависимости величины рассогласования целевых функций от параметров и і и и2 показано на рис. 2 и 3.

Рассмотрим подробнее свойства этой зависимости.

Решения задач первого уровня (3) при фиксированных и і и и 2 очевидны:

/1(ж*(и)) = U1 , /2(ж*(и)) = U2 , /з(ж* (и)) = 3 — U1 — U2 .

Следовательно, задача второго уровня (4) - минимизации рассогласованности критериев, будет иметь вид:

минимизировать по {$і, $2, $з, р} вели чину р при условиях р > 0, ж1 > 0, $2 > 0 и $з > 0, и2г $1 + и1г $2 + U1U2 $3 6 иіи2г ,

Ж1 >  U1 — р,

$ 2 > U 2 — р,

$3 >  г — р, г = 3 — и 1 — и₂ .

Ее решение будем по-прежнему обозначать р**(иі,и2).

Наконец, задача третьего уровня (5), для данного примера такова:

минимизировать по {иі, U2} р**(иі,и₂) при 0.1 6 иі 6 2.5 и 0.1 6 U2 6 2.5.

Из теории математического программирования известно, что свойства зависимости р**(иі,и₂) в первую очередь определяются тем, как множество ограничений модели типа «неравенство» разделяется на активные и неактивные, то есть первые из которых выполняются как равенства, а вторые — как строгие неравенства.

Это разделение зависит от значений параметров модели и его оптимальный вариант определяет решение задачи второго уровня.

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

Рис. 2. Величина |ОЛ| + |ОВ| + |ОС| постоянна (ЗВ-график)

Рис. 3. Величина |ОЛ| + |ОВ| + |ОС| постоянна (вид изолиний)

В этом случае должны быть активны последние пять ограничений из задачи второго уровня, и мы получаем следующую систему уравнений, позволяющую найти аналитиче-

ский вид зависимости р**(пі,п2) :

П2Г Х1 + П1Г Х2 + П1 П₂ Хз = П1П2Г ,

Х1 = п - р,

< Х2 = П2 - Р, Хз = Г - р ,

. Г = 3 - П1 - П2 , исключая из которой неизвестные хі, Х2, хз и г, получаем искомую зависимость р от щ и

П2 в виде

Р**(П1,П2) =

.--------1---------- п1 п2   3 - п1 - п2

Нетрудно проверить, 11-3 3 hT || 3 - 3 |Г и

что стационарными для р** (111,112) являются точки || 1 1 ||T, || 3 3 |T , в первой из которых, согласно критерию Сильвестра,

3,

эта функция имеет локальный максимум со значением а остальные являются недопу

стимыми, поскольку для них нарушены условия неотрицательности для х і , Х2, Хз и г, .

Полученная формула справедлива лишь в некоторой области, содержащейся в Ө. Анализ системы изолиний, изображенных на рис. 3, позволяет выделить пять областей с различ ными наборами активных ограничений. Границы между областями показаны светлыми линиями. Полученная выше формула справедлива только в области 4. В этой области множество Парето рассматриваемой модели состоит из внутренних точек треугольника АВС.

Вне области 4 формула для р**(пі,П2) другая. Например, для области 1 р**(пі,П2) на ходится из системы уравнений

Н2Г Х1 + И1Г Х₂ + 111 U2 Хз = П1П2Г ,

Х1 = 0,

< Х2 = П2 - р ,

Хз = Г - р,

Г = 3 - П1 - П2 , поскольку множество активных ограничений в ней другое: оно содержит условие хі = 0

вместо условия хі = пі - р. Проверьте самостоятельно, что в области 1

р**(П1,П2) =

-+;---1--- п 3 - П1 - П2

Стационарных точек у этой зависимости нет.

Для областей 2 и 3 рассуждения и результаты аналогичные. Множествами Парето в областях 1, 2 и 3 являются стороны треугольника АВС: ВС, АС и АВ соответственно. Наконец, отметим, что в области 5 система условий (2) противоречива.

При этом точное решение задачи верхнего (третьего) уровня имеет вид

«Г = 1, п2* = 1, р** = 3.

Опишем теперь метод решения задачи третьего уровня для рассматриваемого демонстрационного варианта многокритериальной модели. Мы имеем целевые функции: максимизировать по х = ||хі Х2 хз|т € Ез

/1(х,п) = Х1,

/2(х,п) = Х2,

/з(х,п) = хз, при условиях: х1 > 0, х2 > 0, хз > 0 , оі(и)хі + о2(п)х2 + аз(п)хз 6 Ь(п) , где

О1(п) = U2(3 — Щ — U2) , о2(п) = П1(3 — П1 — U2) , 02 (и) = П1П2 ,

Ь(п)   = п1п2(3 — п1 — п2) .

У1

У2

У з

У 4

-

-

-

Х1 ,

Х2 ,

Х1 ,

01x1 + о₂х₂ + 03x3 — b ,

а в качестве штрафной функции

возьмем Р (т, s) = т exp - . ^т

При этом

0.1 6 пі 6 2.5 п = «п1п2И 6 ө={п | 0.J 6 „2 6 2.5 }

Введем обозначения

Тогда вспомогательные функции однокритериальных задач (8) будут иметь вид

А_к (т, х, п) = х к — ^ т exp

3 =1

(-)

для к = 1, 2, 3 .

Условия стационарности этих вспомогательных функций по компонентам х задаются следующими системами уравнений:

ЭАк ж— = §kj + exp дхз

(У") — 03 exp (У4)

Vk = 1, 2, 3,      V; = 1,2,3 .

Наконец, производные от сглаженных аппроксимирующих функций F к (п) по компонентам вектора параметров п находятся из соотношений (10):

(^F к )^ = |^{Ак (т},^хк^(п),^п) = f⁴^ ^{+ х}к2^ ^{+ х}к31^{03 —} ^^ ^exp f ^^4^

р     и п р                 \     и п р        и п р        и п р и п р /       \ т /        (1о)

Vk = 1,2,3, Vp = 1,2 .

Теперь рассмотрим задачу (11) - оптимизации рассогласования критериев модели. В нашем случае эта задача имеет вид:

минимизировать р по совокупноети переменных {хі, х2, хз, р} при условиях х1 > 0, х2 > 0, хз > 0, р > 0

и

Оіхі + О2х2 + Озхз 6 b.

Кроме того, должны выполняться неравенства хк > F к — р Vk = 1, 2, 3.

Эту задачу будем также решать методом гладких штрафных функций с той же Р(т, s), для чего будет удобно ввести (помимо определенных ранее) обозначения Ук = Fк — р — хк Vk = 1, 2, 3. В этом случае вспомогательная функция (12) имеет вид

Е^{( т,х,} р ^{,п )} = ^—р

-

^т exp ( -^ ^— Е^{т e}xp (^Ук)^— Е^{т e}xp (7).

Условия стационарности вспомогательной функции (13) в случае (18) принимают вид системы уравнений

ЭЕ

Эр

ЭЕ дх к

exp

—1 + exp

— exp к=1

(Ү)

(Ү) + exp (^ -к exp (^

0,

0 ^к = 1, 2, 3.

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

Е (и) = —Е (т, р(и), х(и), и) .

Производные этой функции по компонентам вектора параметров и находятся по формулам (14) и для рассматриваемого примера имеют вид

^,Л'       У4 дУ4 3     Үк дҮк кЕ (и)) иР       7) • ди; + g7) • дй" =

= exp

( 7 )-(j5x,^ —

дЬ ди;

\   v³^     (Ү к    ^дҒ к

:) ⁺ к = ^exp(

^р = 1, 2 .

При этом значения производных от функций Ғ к (и) находятся по формулам (16).

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

Таблица. 1

t	И 1	П 2	Е	Р
0	0.700000000	1.600000000	0.580923855	0.545812501
1	0.700000000	1.200000000	0.633041421	0.596126653
2	0.900876900	1.159092100	0.654535635	0.621144210
3	0.941239800	1.000136600	0.660271356	0.626927551
4	0.981719000	1.018621200	0.661487390	0.628152833
5	1.000234400	1.000052500	0.661620557	0.628286990
6	1.000071000	0.999946090	0.661620583	0.628287016

Таблица. 2

t	N grad	W 1	W 2	s
0	0.231363725	-1.8973*10~-9	-1.000000000	0.400000000
1	0.131142557	0.979887531	-0.199550560	0.205000000
2	0.076755604	0.246115279	-0.969240563	0.164000000
3	0.050581430	0.909645878	0.415384614	0.044500000
4	0.010207109	0.706089900	-0.708122202	0.026222500
5	2.43444*10 -4	-0.837971176	-0.545714493	0.000195000
6	3.73575*10 -5	-0.922621240	0.385707203	0.000076500

Пусть t = 0, 1, 2, ... - номер итерации. Тогда эта схема будет описываться соотноше ниями иР+1 = иР + st ^р ,    гДе ^р = v---(Е(и^)   , р = 1, 2

Wgrad \       /Ир

Значение нормы градиента, вычисляется по обычной для ортонормированого базиса, фор муле

^grad = а величина st - шага по улучшающему направлению, для каждой итерации оценивается методом дихотомии.

Результаты расчетов приведены в табл. 1 и 2.

Замечания об использовании алгоритмов второго порядка

В заключение рассмотрим метод нахождения для функции Е (и) элементов матрицы Гессе, знание которых позволит использовать в процессе поиска, стационарных точек методы второго порядка.

Применяя к функции (20) правила, дифференцирования сложной функции, получаем

( е (и)) "   = X        / PpUq

д 2 Е    ” д2Е дх ,    д 2 Е др дирдиу + ~^ дирдх , диу + дирдр диу,    ^^   ^,Г.

Вторые частные производные вычисляются непосредственно в точке {х, р} , а первые про-дх, Эр изводные, т.е. —— и ——, находятся согласно известной теореме о неявных функциях из диу   диу системы линейных уравнений

< п д^Е дхі       д2Е + ^1 дрдхі диу   + дрдиу        ’ п д 2Е дхг       д^Е + Етт^+         = 0   ^3 = [1,п] , і=і Эх, Эхі диу     дх, диу которая, в свою очередв, получается при последовательном дифференцировании по переменным р и Xj j = [1,п] условий стационарности (13).

Наконец, отметим, что для вычисления производных (21) также необходимо знать значения вторых производных от функций F^(и). Эти значения могут быть найдены (аналогичным использованному выше) методом из формул (10) и условий (9) - стационарности функций Ак(т,х,и).