Доклады Российской академии наук. Математика, информатика, процессы управления, 2020, T. 492, № 1, стр. 20-23

Проблемы глобальной оптимизации и многомерного интегрирования уже много лет привлекают большое внимание исследователей. Это связано с многочисленными приложениями таких задач в самых разных областях. На указанные темы опубликованы тысячи статей и книг, и это число постоянно растет. В качестве примера можно указать недавние публикации [1–3] в международном журнале Journal of Global Optimization и книгу [4], содержащую основательный обзор методов многомерного интегрирования.

Насколько известно автору, все работы по указанной проблематике так или иначе основаны на манипуляциях с множеством значений аргумента функции и/или ее различных аппроксимациях (дифференциальные методы, методы случайного поиска, аппроксимативные методы, методы снижения размерности аргумента, методы интегрирования Монте-Карло и др.).

Здесь предлагается иной подход, основанный на работе c множеством значений функции и сочетании идей лебеговской теории меры и интеграла с простым вариантом метода Монте-Карло. Приведенные ниже элементарные оценки и первые результаты вычислительных экспериментов (этому будет посвящена отдельная публикация) позволяют надеяться, что предлагаемый подход (без проблем реализуемый на современных компьютерах) окажется более эффективным, чем известные в литературе.

1. Предположения. Всюду далее используются следующие предположения:

а) исследуемая функция f(x) непрерывна и определена на компактной области $D \subset {{\mathbb{R}}^{d}}$;

б) априори известны числа (A, B) такие, что

(здесь m и M – вообще говоря, неизвестные значения глобального минимума и глобального максимума функции).

Если предположение а) не представляется обременительным, то предположение б) требует комментария. Компактность области D позволяет найти (пусть и весьма грубые) интервалы изменения каждой переменной аргумента. После этого, если функция задана аналитически, можно дать (опять-таки грубые) оценки чисел m и M. Если же речь идет о реальных задачах, то оценки для m и M можно, как правило, указать, исходя из содержательного смысла. Таким образом, нам представляется, что условие б) не является слишком ограничительным.

2. Глобальная максимизация. Критерием качества существующих методов поиска является степень близости значения функции в точке, найденной поиском, и значения глобального экстремума.

Предлагается новый критерий качества поиска глобального экстремума (для определенности, максимума). А именно, требуется найти такую точку $x* \in D$, для которой отношение

не будет превосходить заданного малого числа (и не важно, насколько отличается $f(x*)$ от значения глобального максимума). Здесь и далее $\mu ( \cdot )$ – мера Лебега.

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

Пусть ${\text{\{ }}{{\xi }_{n}}{\text{\} }}$ – последовательность независимых и равномерно распределенных на D случайных векторов. Рассмотрим последовательность независимых и одинаково распределенных случайных величин ${{y}_{n}}\mathop = \limits^{{\text{def}}} {\text{\{ }}f({{\xi }_{n}}){\text{\} }}$ с общей функцией распределения $F(t) = {\mathbf{P}}{\text{\{ }}f(\xi ) \leqslant t{\text{\} }}$. Введем следующий случайный процесс на [A, B]:

где $\mathbb{I}(A)$ – индикаторная функция множества A. Процесс ${{Z}_{n}}(t)$ – это эмпирическая функция распределения (ф.р.) случайной величины $f(\xi )$.

Предлагаемый метод поиска глобального максимума состоит в следующем. Разбиваем отрезок [A, B] на k равных отрезков точками $({{t}_{1}},.,{{t}_{k}})$, причем $A = {{t}_{1}} < {{t}_{2}} < \; \ldots \; < {{t}_{k}} = B$ и строим k-мерный вектор ${{\mathbb{Z}}_{n}} = \left( {({{Z}_{n}}({{t}_{1}}), \ldots ,\;{{Z}_{n}}({{t}_{k}})} \right)$. Введем k-мерный вектор $\mathbb{F} = \left( {(F({{t}_{1}}),\; \ldots ,\;F({{t}_{k}})} \right)$. В силу классических предельных теорем с вероятностью 1 $\left\| {{{\mathbb{Z}}_{n}} - \mathbb{F}} \right\|$ → 0 при $n \to \infty $. Далее в качестве нормы векторов будем использовать норму l₁, т.е. $\mathop {max}\limits_s \left| {{{Z}_{n}}({{t}_{s}}) - F({{t}_{s}})} \right|$.

Пусть $C \subset D$ – произвольное (измеримое) подмножество D. Пусть $\eta $ – случайный вектор, равномерно распределенный на D. Тогда

Поэтому, если мы укажем такое n*, что величина ${\mathbf{P}}\left\{ {f(\xi ) \geqslant {{y}_{{(n*)}}}} \right\}$, где ${{y}_{{(n*)}}}$ – максимальный член соответствующего вариационного ряда ${\text{\{ }}{{y}_{{(i)}}}{\text{\} }}_{{i = 1}}^{{n*}}$, будет достаточно мала, то, зная лебегову меру множества D, мы из (1) получим оценку для абсолютной величины критерия качества процесса поиска глобального максимума в предлагаемой постановке. При этом величина ${{y}_{{(n*)}}}$ будет оценкой значения глобального максимума, а соответствующий аргумент функции f будет оценкой положения этого максимума. Сама величина этой вероятности даст возможность оценить относительную (по отношению к мере всей области определения функции) долю “плохих” точек.

Пусть $\epsilon > 0$ – фиксированное (малое) число. Мы ставим задачу нахождения такой точки $x* \in D$, для которой

Процедура поиска состоит из двух этапов. На первом этапе ставится задача определения числа k отрезков разбиения множества [A, B]. Положим

и заметим, что ${{G}_{n}}({{t}_{k}}) = G({{t}_{k}}) = 0$ по определению при любых $n$ и k. Цель первого этапа поиска – найти такое k, что $G({{t}_{{k - 1}}}) < \epsilon $. Пусть ${{k}_{0}} = \left[ {\tfrac{1}{{4\epsilon }}} \right] + 1$ – начальное значение числа k ([a] – целая часть числа a). Пусть $n \geqslant n{\text{*}}({{k}_{0}},\epsilon ) = \frac{{k_{0}^{2}}}{{{{\epsilon }^{2}}}}$ и ${{G}_{n}}({{t}_{{k - 1}}}) < 2\epsilon $. Тогда при $n \geqslant n{\text{*}}({{k}_{0}},\epsilon )$ с вероятностью, большей, чем $(1 - \epsilon )$, имеем $G({{t}_{{k - 1}}}) < \epsilon $. В этом случае k₀ – то значение, которое нам требуется. В противном случае значение k требуется увеличить, и за конечное число шагов найдем такое k*, что с вероятностью, большей, чем $\left( {1 - \tfrac{1}{{4k{\text{*}}}}} \right)$, при $n \geqslant \tfrac{{{{{(k*)}}^{2}}}}{{{{\epsilon }^{2}}}}$ будет выполнено неравенство ${{G}_{n}}({{t}_{{k - 1}}}) < 2\epsilon $. Фиксируем это значение k*.

Предложение 1. Пусть для любого (достаточно малого) $\epsilon > 0$ число k* отрезков разбиения [A, B] выбрано в соответствии с описанной выше процедурой, количество $n$ сгенерированных независимых случайных векторов, равномерно распределенных на D, удовлетворяет условию $n \geqslant \mathop {\left( {\tfrac{{k{\text{*}}}}{\epsilon }} \right)}\nolimits^2 $, оценкой значения глобального максимума функции f(x) на D является величина максимального члена вариационного ряда выборки ${\text{\{ }}f(\xi ){\text{\} }}_{{i = 1}}^{n}$, а оценкой положения этого максимума – соответствующее значение аргумента. Тогда предлагаемый критерий качества поиска удовлетворяет оценке

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

3. Многомерное интегрирование. Ставится задача оценивания интеграла Римана J = $\int\limits_D^{} {f(x)dx} $, где функция f(x) и область D удовлетворяют условиям раздела 1. Как известно, если интеграл Римана существует, то он совпадает с интегралом Лебега. Положим

Пусть $\xi $ – случайный вектор с равномерным распределением на D и $F(t) = {\mathbf{P}}{\text{\{ }}f(\xi ) \leqslant t{\text{\} }}$. Тог-да, в соответствии с идеей интеграла Лебега (см. (1)),

Таким образом, оценка интеграла сводится к вычислению математического ожидания случайной величины $f(\xi )$. Если в результате генерации последовательности независимых и одинаково равномерно распределеных на D векторов ${\text{\{ }}{{\xi }_{i}}{\text{\} }}_{{i = 1}}^{n}$ построена (на разбиении $A = {{t}_{1}} < {{t}_{2}}\; \ldots \; < {{t}_{{k - 1}}} < {{t}_{k}}$ = B, см. раздел 1) эмпирическая функция распределения ${{Z}_{n}}(t)$, то оценка $\hat {J}$ интеграла J имеет вид

Обозначим $B - A = L$, $C = \left( {\left| A \right| + \tfrac{L}{2}} \right)$. Имеем

Предложение 2. Для любого (достаточно малого) $\epsilon > 0$, задавая число отрезков разбиения ближайшим целым к $\tfrac{{2L}}{{\epsilon \mu (D)}}$ и выбирая $n > \tfrac{{4{{L}^{2}}{{C}^{2}}}}{{\mu {{{(D)}}^{4}}{{\epsilon }^{5}}}}$, с вероятностью не менее $(1 - \epsilon )$ получаем оценку

Замечание 2. Если μ(D) неизвестно, то, погружая область D в некоторый параллелепипед (в силу компактности это можно сделать), можно найти оценку для $\mu (D)$ простым методом Монте-Карло, генерируя последовательность независимых и равномерно на параллепипеде распределенных случайных векторов и отбирая из этой последовательности те векторы, которые попадают в область D (ср. с Замечанием 1).