Метод Стронгина

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

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

Требуется найти точку $x^{*}\in [a;b]$ такую, что $f(x^{*})=\min \{f(x):x\in [a;b],g_j(x)⩽0,1⩽j⩽m\}$. Предполагается, что функции $f(x)$ и $g_j(x),j=\overline{1,m}$ удовлетворяют условию Липшица на отрезке $[a;b]$.

Обозначим $g_{m+1}(x)=f(x)$, тогда для $j=\overline{1,m+1}$ выполняются следующие неравенства:

$|g_j(x+\mathrm{\Delta }x)-g_j(x)|⩽L_j\mathrm{\Delta }x,a⩽x+\mathrm{\Delta }x⩽b,$

где $L_j⩾0$ — константы Липшица.

Описание схемы учета ограничений

Пусть $Q_0=[a;b]$. Ограничение, имеющее номер $j$, выполняется во всех точках области $Q_j=\{x\in [a;b]:g_j(x)⩽0\}$, которая называется допустимой для этого ограничения. При этом допустимая область $Q$ исходной задачи определяется равенством:

$Q=\bigcap _{j=0}^m{Q}_j.$

Испытание в точке $x\in [a;b]$ состоит в последовательном вычислении значений величин $g_1(x),\dots ,g_{\nu }(x)$, где значение индекса $\nu$ определяется условиями:

$x\in Q_j,0⩽j<\nu ,x\notin Q_{\nu }.$

Выявление первого нарушенного ограничения прерывает испытание в точке $x$. В случае, когда точка $x$ допустима, то есть $x\in Q$ испытание включает в себя вычисление всех функций задачи. При этом значение индекса принимается равным величине $\nu =m+1$.

Пара $\nu =\nu (x),z=g_{\nu }(x)$, где индекс $\nu$ лежит в границах $1⩽\nu ⩽m+1$, называется результатом испытания в точке $x$.

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

$\psi (x^{*})=\underset{x\in [a;b]}{\min }\psi (x),$

$\psi (x)=\{\begin{array}{cc}{g}_{\nu }(x)/L_{\nu }& \nu <M,\\ (g_M(x)-g_M^{*})/L_M& \nu =M.\end{array}$

Здесь $M=\max \{\nu (x):x\in [a;b]\}$, а $g_M^{*}=\min \{g_M(x):x\in \bigcap _{i=0}^{M-1}{Q}_i\}$.

В силу определения числа $M$, задача отыскания $g_M^{*}$ всегда имеет решение, а если $M=m+1$, то $g_M^{*}=f(x^{*})$.

Дуги функции $\psi (x)$ липшицевы на множествах $\bigcap _{i=0}^j{Q}_i,0⩽j⩽M-1$ с константой 1, а сама $\psi (x)$ может иметь разрывы первого рода на границах этих множеств.

Несмотря на то, что значения констант Липшица и величина $g_M^{*}$ заранее неизвестны, они могут быть оценены в процессе решения задачи.

Описание метода

Пусть $x^0=a,x^1=b$. Индексы концевых точек считаются нулевыми, а значения $z$ в них не определены. Первое испытание осуществляется в точке $x^3=(a+b)/2$. Выбор точки $x^{k+1},k⩾3$ любого последующего испытания определяется следующими правилами:

1. Перенумеровать точки $x^0,\dots ,x^k$ $k$ предшествующих испытаний нижними индексами в порядке увеличения значений координаты: $a=x_0<\dots <x_i<\dots <x_k=b$ и сопоставить им значения $z_i=g_{\nu }(x_i),\nu =\nu (x_i),i=\overline{1,k}$.
2. Для каждого целого числа $\nu ,1⩽\nu ⩽m+1$ определить соответствующее ему множество $I_{\nu }$ нижних индексов точек, в которых вычислялись значения функций $g_{\nu }(x)$:

   $I_{\nu }=\{i:\nu (x_i)=\nu ,1⩽i⩽k\},1⩽\nu ⩽m+1.$

   Также определить максимальное значение индекса $M=\max \{\nu (x_i),1⩽i⩽k\}.$

3. Вычислить текущие оценки для неизвестных констант Липшица:

   ${\mu }_{\nu }=\max \{|g_{\nu }(x_i)-g_{\nu }(x_j)|/(x_i-x_j):i,j\in I_{\nu },i>j\}.$

   Если множество $I_{\nu }$ содержит менее двух элементов или если значение ${\mu }_{\nu }$ оказывается равным нулю, то принять ${\mu }_{\nu }=1$.

4. Для всех непустых множеств $I_{\nu },\nu =\overline{1,M}$ вычислить оценки

   $z_{\nu }^{*}=\{\begin{array}{cc}\min \{g_{\nu }(x_i):x_i\in I_{\nu }\}& \nu =M,\\ -{\epsilon }_{\nu }& \nu <M,\end{array}$

   где вектор с неотрицательными координатами ${\epsilon }_R=({\epsilon }_1,\dots ,{\epsilon }_m)$ называется вектором резервов.

5. Для каждого интервала $(x_{i-1};x_i),1⩽i⩽k$ вычислить характеристику

   $R(i)=\{\begin{array}{cc}{\mathrm{\Delta }}_i+\frac{(z_i-z_{i-1}{)}^2}{(r_{\nu }{\mu }_{\nu }{)}^2{\mathrm{\Delta }}_i}-2\frac{z_i+z_{i-1}-2z_{\nu }^{*}}{r_{\nu }{\mu }_{\nu }}& \nu =\nu (x_i)=\nu (x_{i-1}),\\ 2{\mathrm{\Delta }}_i-4\frac{z_{i-1}-z_{\nu }^{*}}{r_{\nu }{\mu }_{\nu }}& \nu =\nu (x_{i-1})>\nu (x_i),\\ 2{\mathrm{\Delta }}_i-4\frac{z_i-z_{\nu }^{*}}{r_{\nu }{\mu }_{\nu }}& \nu =\nu (x_i)>\nu (x_{i-1}),\end{array}$

   где ${\mathrm{\Delta }}_i=x_i-x_{i-1}$.

   Величины $r_{\nu }>1,\nu =\overline{1,m}$ являются параметрами алгоритма. От них зависят произведения $r_{\nu }{\mu }_{\nu }$, используемые при вычислении характеристик в качестве оценок неизвестных констант Липшица.

6. Определить интервал $(x_{t-1};x_t)$, которому соответствует максимальная характеристика $R(t)=\max \{R(i),1⩽i⩽k\}$.
7. Провести очередное испытание в середине интервала $(x_{t-1};x_t)$, если индексы его концевых точек не совпадают:

   $x^{k+1}=\frac{1}{2}(x_t+x_{t-1}).$

   В противном случае провести испытание в точке

   $x^{k+1}=\frac{1}{2}(x_t+x_{t-1})-\frac{z_t-z_{t-1}}{2r_{\nu }{\mu }_{\nu }},\nu =\nu (x_t)=\nu (x_{t-1}),$

   увеличить $k$ на 1.

8. Если $x_t-x_{t-1}<\epsilon$ ($\epsilon >0$ — заданная точность метода), то прекратить выполнение алгоритма, иначе перейти на шаг 1.

Достаточные условия сходимости

Пусть исходная задача оптимизации имеет решение $x^{*}$ и выполняются следующие условия:

1. каждая область $Q_j,j=\overline{1,m}$ представляет собой объединение конечного числа отрезков, имеющих положительную длину;
2. каждая функция $g_j(x),j=\overline{1,m+1}$ удовлетворяет условию Липшица с соответствующей константой $L_j$;
3. компоненты вектора резервов удовлетворяют неравенствам $0⩽2{\epsilon }_{\nu }<L_{\nu }(\beta -\alpha )$, где $\beta -\alpha$ — длина отрезка $[\alpha ;\beta ]$, лежащего в допустимой области $Q$ и содержащего точку $x^{*}$;
4. начиная с некоторого значения $k$ величины ${\mu }_{\nu }$, соответствующие непустым множествам $I_{\nu }$, удовлетворяют неравенствам $r_{\nu }{\mu }_{\nu }>2L_{\nu }$.

Тогда верно следующее:

1. точка $x^{*}$ является предельной точкой последовательности $\{x^k\}$, порождаемой методом при $\epsilon =0$ в условии остановки;
2. любая предельная точка $x^0$ последовательности $\{x^k\}$ является решением исходной задачи оптимизации;
3. сходимость к предельной точке $x^0$ является двухсторонней, если $x^0\ne a,x^0\ne b$.

Модификации метода

Параллельная модификация

Общая схема последовательного метода выглядит следующим образом:

1. Упорядочить точки предшествующих испытаний в порядке возрастания их координат: $a=x_0<\dots <x_i<\dots <x_k=b$.
2. Вычислить для каждого интервала $(x_{i-1};x_i),1⩽i⩽k$ характеристику $R(i)$.
3. Определить интервал $(x_{t-1};x_t)$, которому соответствует максимальная характеристика $R(t)=\max \{R(i),1⩽i⩽k\}$.
4. Провести следующее испытание в точке $x^{k+1}=d(t)\in (x_{t-1};x_t)$, где $d(t)$ — правило размещения точки следующего испытания в интервале с номером $t$.
5. Проверить выполнение критерия остановки $x_t-x_{t-1}<\epsilon$.

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

Схема параллельного алгоритма:

1. Упорядочить точки предшествующих испытаний в порядке возрастания их координат: $a=x_0<\dots <x_i<\dots <x_k=b$.
2. Вычислить для каждого интервала $(x_{i-1};x_i),1⩽i⩽k$ характеристику $R(i)$.
3. Характеристики интервалов упорядочить по убыванию: $R(i_1)>\dots >R(i_k)$.
4. Для всех интервалов с номерами $i_1,\dots ,i_p$ провести испытания в точках $x^{k+j}=d(i_j)\in (x_{i_j-1};x_{i_j}),j=\overline{1,p}$.
5. Проверить выполнение критерия остановки: $\mathrm{\exists }j,1⩽j⩽p:x_{i_j}-x_{i_j-1}<\epsilon$.

Такая схема распараллеливания целесообразна, если проведение испытания (то есть вычисление функций задачи) — трудоемкий процесс.

Модификация для решения задач c гёльдеровыми функциями

Метод достаточно просто обобщается на случай, когда функции $g_j(x),j=\overline{1,m+1}$ удовлетворяют условию Гёльдера с показателем $1/n$, где $n\in \mathbb{N}$, то есть

$|g_j(x+\mathrm{\Delta }x)-g_j(x)|⩽H_j(\mathrm{\Delta }x{)}^{1/n},a⩽x+\mathrm{\Delta }x⩽b$.

На шаге 3 значения ${\mu }_{\nu }$ вычисляются по формуле:

${\mu }_{\nu }=\max \{|g_{\nu }(x_i)-g_{\nu }(x_j)|/(x_i-x_j{)}^{1/n}:i,j\in I_{\nu },i>j\}.$

На шаге 5 ${\mathrm{\Delta }}_i=(x_i-x_{i-1}{)}^{1/n}$.

На шаге 7 в случае совпадения индексов концевых точек

$x^{k+1}=\frac{1}{2}(x_t+x_{t-1})-\mathrm{sgn}(z_t-z_{t-1})\frac{|z_t-z_{t-1}{|}^n}{2r_{\nu }{\mu }_{\nu }^n},\nu =\nu (x_t)=\nu (x_{t-1}).$

На шаге 8 критерий остановки принимает вид $(x_t-x_{t-1}{)}^{1/n}<\epsilon$.

Замечания

• Параметры $r_{\nu }$ отвечают за надежность метода. Чем больше их значения, тем больше итераций метода требуется для достижения заданной точности и тем вероятнее выполнение условия сходимости 4. Если устремить все $r_{\nu }$ к бесконечности, то $R(i)={\mathrm{\Delta }}_i$, то есть метод превращается в перебор по равномерной сетке.
• Использование ненулевого вектора резервов позволяет ускорить сходимость метода, однако при этом необходимо оценить возможность выполнения условия сходимости 3.
• Одномерный метод может быть применен для решения многомерных задач без ограничений. Многомерная задача на множестве $S=\{(x_1,\dots ,x_n)\in {ℝ}^n:a_i⩽x_i⩽b_i,i=\overline{1,n}\}$ представляется в виде

$\underset{(x_1,\dots ,x_n)\in S}{\min }f(x_1,\dots ,x_n)=\underset{a_1⩽x_1⩽b_1}{\min }\underset{a_2⩽x_2⩽b_2}{\min }\dots \underset{a_n⩽x_n⩽b_n}{\min }f(x_1,\dots ,x_n).$

Для решения задачи $\underset{a_1⩽x_1⩽b_1}{\min }\varphi (x_1)$, где $\varphi (x_1)=\underset{a_2⩽x_2⩽b_2}{\min }\dots \underset{a_n⩽x_n⩽a_n}{\min }f(x_1,\dots ,x_n)$ можно использовать одномерный алгоритм, но, чтобы вычислить значение функции $\varphi (x_1)$, необходимо решить задачу оптимизации размерности $n-1$.

Если $n=2$, то задача $\underset{a_2⩽x_2⩽b_2}{\min }f(x_1,x_2)$ решается одномерным методом (значение переменной $x_1$ при этом зафиксировано), иначе к ней также применяется процедура снижения размерности. Такой способ решения многомерных задач довольно трудоемкий, поэтому на практике применим при $n⩽5$. Наличие нелинейных функциональных ограничений может привести к потере липшицевости во вспомогательных одномерных задачах.

Литература

1. Баркалов К. А., СтронгинР. Г. Метод глобальной оптимизации с адаптивным порядком проверки ограничений // Ж. вычисл. матем. и матем. физ. — 2002. — Т. 42. — № 9. — стр. 1338—1350.
2. Городецкий С. Ю., Гришагин В. А. Нелинейное программирование и многоэкстремальная оптимизация. — Нижний Новгород: Издательство Нижегородского Университета, 2007.
3. Стронгин Р. Г. Численные методы в многоэкстремальных задачах (информационно-статистические алгоритмы). — Шаблон:М.: Наука, 1978. — 240 с.
4. Sergeyev Ya. D., Grishagin V. A. Sequential and parallel algorithms for global optimization // Optimization Methods and Software, 3:1-3, 1994, pp. 111—124.
5. Маркин Д. Л., Стронгин Р. Г. Метод решения многоэкстремальных задач с невыпуклыми ограничениями, использующий априорную информацию об оценках оптимума // Ж. вычисл. матем. и матем. физ., 27:1 (1987), стр. 56—62.

Ссылки

• [1] - реализация метода на языке C++.
• [2] - реализация на языке C++ модификации метода метода для решения многокритериальных многомерных задач.