Оптимизация с линейными ограничениями в виде равенств и неравенств

В этой статье рассмотрен субпакет minbleic - оптимизатор, поддерживающий простые ограничения и линейные ограничения в виде равенств/неравенств. Этот субпакет пришел на смену устаревшему модулю minasa, который продолжает поддерживаться из соображений обратной совместимости. Алгоритм BLEIC (boundary, linear equality-inequality constraints) способен решать задачи вида

Содержание     1 Об алгоритме            MinBLEIC как метод активных множеств            Дополнительная информация     2 Приступая к работе            Выбор между аналитическим и численным градиентом            Масштаб переменных            Предобуславливатель            Условия остановки            Ограничения     3 Примеры

Об алгоритме

MinBLEIC как метод активных множеств

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

В итоге мы получаем метод, состоящий из двух чередующихся стадий:

• на первой стадии (внутренние итерации) мы решаем субпроблему с ограничениями-равенствами.
• на второй стадии (внешняя итерация) принимается решение об активации или деактивации ограничений (обычно в зависимости от знака множителей Лагранжа), после чего мы переходим к новому циклу внутренних итераций.

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

1. мы начинаем со стартовой точки, в которой все ограничения неактивны
2. в результате решения первой субпроблемы мы приходит в точку (0,1), в которой активируется ограничение (1).
3. решение второй субпроблемы приводит нас к точке (0,0), в которой мы активируем ограничение (3) и деактивируем ограничение (1)
4. наконец, третья субпроблема приводит нас к точке (1,0), в которой активируются сразу два ограничения - (2) и (3), после чего алгоритм прекращает работу

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

Дополнительная информация

В этом разделе мы рассмотрим ключевые свойства алгоритма, реализованного в субпакете minbleic. Ниже мы будем использовать следующие обозначения: N будет обозначать количество переменных в задаче, K - количество линейных ограничений общего вида (равенств и неравенств), Ki - количество линейных ограничений-неравенств.

Можно отметить следующие свойства алгоритма:

• в качестве базового алгоритма оптимизации используется метод сопряженных градиентов
• для обработки ограничений-равенств мы модифицируем целевую функцию, проецируя её на подпространство, задаваемое ограничениями
• если это необходимо, то ограничения-неравенства обрабатываются путем активации (т.е. интерпретации как равенств)
• особенностью алгоритма является то, что он не поддерживает смену предобуславливателя "на лету". Модификация предобуславливателя во время процесса оптимизации не окажет влияния на решение. Хотя при следующем запуске оптимизатора будет использоваться уже новый предобуславливатель.
• алгоритм прозрачным для пользователя образом добавляет к задаче служебные переменные так что все линейные ограничения-неравенства превращаются в линейные равенства плюс одно простое ограничение (неотрицательность) на служебную переменную. В результате к задаче с Ki ограничениями-неравенствами добавляется Ki служебных переменных.
• простые ограничения обрабатываются отдельно от линейных ограничений общего вида, как важный частный случай. Простые ограничения требуют меньше времени на обработку и всегда выполняются в точности - как в окончательном решении, так и на промежуточных стадиях. Например, при ограничении вида x≥0 мы никогда не будем вычислять функцию в точках с отрицательным x.
• линейные ограничения общего вида (равенства или неравенства) требуют большего времени на обработку, чем простые ограничения. Линейные ограничения - как равенства, так и неравенства - могут выполняться с небольшой погрешностью (порядка N·ε). Например, при ограничении вида x+y≤1 мы можем остановиться в точке, которая фактически лежит по другую сторону от линии x+y=1 (на расстоянии порядка N·ε от допустимого множества).
• наличие ограничений добавляет накладные расходы. Дополнительные вычисления требуются на двух стадиях: при активации/деактивации ограничений и при каждом вычислении функции. При активации/деактивации хотя бы одного ограничения требуется реортогонализировать всю матрицу ограничений, что требует O((N+Ki)·K ²) операций. При каждом вычислении функции требуется O(N) дополнительных операций на вектор простых ограничений, и O((N+Ki)·K) - на матрицу линейных ограничений.
• Вычислительные расходы не зависят от того, активно ограничение или нет. Единственная ситуация, когда расходы на обработку ограничения незначительны - это простое ограничение вида -∞<x или x<+∞.

Приступая к работе

Выбор между аналитическим и численным градиентом

Самый первый выбор, который вам предстоит сделать - Например, если вы ищете минимум функции f(x,y)=x ²+exp(x+y)+y ², то для работы алгоритма потребуются как значение функции в какой-то промежуточной точке, так и производные df/dx и df/dy. Откуда взять эти производные? Пользователям пакета ALGLIB доступно несколько вариантов:

1. пользователь вычисляет производные самостоятельно, обычно путем аналитического дифференцирования функции (т.н. аналитический градиент). Этот вариант является предпочтительным по двум причинам. Во-первых, значение градиента вычисляется с минимальной погрешностью. Во-вторых, очень часто для вычисления всех N компонент градиента требуется лишь в несколько раз больше времени, чем на вычисление одного значения функции. Это связано с тем, что знание внутренней структуры дифференцируемой функции часто позволяет вычислять вектор производных параллельно с вычислением самой функции.
2. градиент вычисляется пакетом ALGLIB путем численного дифференцирования, с использованием четырехточечной разностной формулы. При этом пользователь вычисляет только значение функции, оставляя вопросы, связанные с дифференцированием, на усмотрение пакета ALGLIB. Этот вариант удобнее предыдущего, потому что пользователю не надо тратить время на написание кода, вычисляющего производную функции. Это позволяет очень быстро создавать работающие прототипы. Однако, есть и недостатки. Во-первых, численный градиент менее точен, чем аналитический, что может замедлить скорость сходимости алгоритма. Во-вторых, численное дифференцирование требует 4*N вычислений функции для получения только одного значения градиента. В результате численное дифференцирование можно использовать только на задачах малой размерности (не более нескольких десятков). На задачах большей размерности алгоритм также будет работать, но очень, очень долго.
3. пользователь вычисляет производные с использованием автоматического дифференцирования. В результате оптимизатор получает быстрый аналитический градиент, а пользователь освобождается от необходимости дифференцировать функцию вручную. Пакет ALGLIB не поддерживает автоматическое дифференцирование, но мы можем рекомендовать несколько пакетов, которые можно использовать для того, чтобы вычислить градиент и передать его в ALGLIB:
   • AutoDiff для программистов на C#
   • FuncDesigner для программистов на Python
   • FADBAD++ для программистов на C++

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

• создаете оптимизатор с помощью функции minbleiccreate
• для запуска оптимизации используете вариант функции minbleicoptimize, принимающий в качестве параметра подпрограмму, вычисляющую значение функции и градиента (одновременно) в заданной точке. Если вы ошибочно выберете вариант, использующий только функцию, то оптимизатор сгенерирует исключение при первой же попытке вычислить градиент.

Если вы решили использовать численное дифференцирование, то вы:

• создаете оптимизатор с помощью функции minbleiccreatef. Эта функция отличается от своего брата-близнеца тем, что принимает дополнительный параметр - шаг дифференцирования. Дифференцирование осуществляется с фиксированным шагом, который, однако, может быть различным для различных переменных (в зависимости от их масштаба, установленного подпрограммой minbleicsetscale).
• для запуска оптимизации используете вариант функции minbleicoptimize), принимающий в качестве параметра подпрограмму, вычисляющую значение функции (только функции) в заданной точке. Если вы ошибочно выберете вариант, использующий градиент, то оптимизатор сгенерирует исключение.

Масштаб переменных

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

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

Предобуславливатель

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

Предобуславливатель нужен:

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

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

Пакет ALGLIB поддерживает несколько типов предобуславливателей:

• предобуславливатель "по умолчанию", который не делает ничего. Этот тип предобуславливателя может быть установлен при помощи функции minbleicsetprecdefault.
• диагональный предобуславливатель на основе приближенного Гессиана. Для такого предобуславливателя необходимо рассчитать диагональ приближенного Гессиана функции (не обязательно точного Гессиана), после чего вызвать функцию minbleicsetprecdiag. Диагональная матрица должна быть положительно определенной - если вы попытаетесь передать диагональ с нулем или отрицательным элементом, то алгоритм сгенерирует исключение. Этот тип предобуславливателя можно использовать для выпуклых функций, либо в ситуациях, когда можно гарантировать, что приближенный Гессиан будет положительно определенным .
• диагональный предобуславливатель на основе масштаба переменных. Этот предобуславливатель включается функцией minbleicsetprecscale. Его можно использовать, если ваши переменные сильно отличаются по своему масштабу (в сотни и более раз). В таких задачах основной сложностью является именно различный масштаб, создающий сложности для оптимизатора. Для того, чтобы такой предобуславливатель работал, необходимо сначала установить масштаб переменных (см. предыдущий раздел).

Условия остановки

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

Пользователю доступны три типа внутренних критериев остановки, которые устанавливаются вызовом функции minbleicsetinnercond:

• после снижения градиента F(x) до заданной величины
• после совершения достаточно малого шага
• после достаточно малого изменения функции на последнем шаге

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

Замечание #1
В общем случае нельзя гарантировать, что сработает именно тот критерий остановки, который вы установили. Например, алгоритм может сделать шаг, который приведет нас точно в минимум функции, и тогда сработает критерий, связанный с нулевым значением градиента - независимо от того, какие критерии были установлены. Возможны и другие ситуации, когда срабатывает не тот критерий, который вы установили (например, из-за погрешностей операций с плавточкой).

Замечание #2
Некоторые критерии остановки используют масштаб переменных, который должен быть установлен при помощи вызова отдельной функции (см. предыдущую секцию).

Для остановки внешних итераций, т.е. алгоритма в целом, доступны два типа критериев:

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

Ограничения

Алгоритм BLEIC поддерживает три типа ограничений:

• простые ограничения, т.е. ограничения вида l_i ≤x_i ≤u_i 
• линейные ограничения-неравенства, т.е. ограничения вида a₀ ·x₀ +...+a_N-1 ·x_N-1 ≥b или a₀ ·x₀ +...+a_N-1 ·x_N-1 ≤b
• линейные ограничения-равенства, т.е. ограничения вида a₀ ·x₀ +...+a_N-1 ·x_N-1 =b

Простые ограничения устанавливаются при помощи функции minbleicsetbc. Эти ограничения обрабатываются очень эффективно - N ограничений требует всего лишь O(N) дополнительных операций при каждом вычислении функции. Наконец, эти ограничения всегда выполняются в точности. Мы не вычисляем функцию в точках, координаты которых лежат за пределами интервала, задаваемого ограничениями. Результат оптимизации также всегда строго удовлетворяет ограничениям (находится либо во внутренней области, либо на границе).

Линейные ограничения общего вида могут быть равенствами или неравенствами. Эти ограничения устанавливаются при помощи функции minbleicsetlc. Линейные ограничения обрабатываются менее эффективно, чем простые: они требуют O((N+Ki)·K) дополнительных операций при каждом вычислении функции, где N - размерность задачи, K - число линейных ограничений, Ki - число ограничений-неравенств. При активации/деактивации хотя бы одного ограничения требуется реортогонализировать всю матрицу ограничений, что требует O((N+Ki)·K ²) операций. В отличие от простых ограничений, линейные ограничения могут выполняться с небольшой погрешностью (порядка N·ε). Например, при ограничении вида x+y≤1 мы можем остановиться в точке, которая фактически лежит по другую сторону от линии x+y=1 (на расстоянии порядка N·ε от допустимого множества).

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

Примеры

Для облегчения работы с алгоритмом BLEIC мы подготовили несколько примеров:

• minbleic_d_1 - этот пример демонстрирует оптимизацию с использованием простых ограничений
• minbleic_d_2 - этот пример демонстрирует оптимизацию с использованием линейных ограничений общего вида
• minbleic_ftrim - этот пример демонстрирует минимизацию функции с сингулярностями на границе области определения. Для того, чтобы алгоритм не застял на границе области, используется прием под названием 'подрезка функции'.
• minbleic_numdiff - этот пример демонстрирует минимизацию функции с использованием численного дифференцирования.

Также мы рекомендуем прочитать статью 'Оптимизация - советы и приемы', которая обсуждает ряд типичных проблем, возникающих при оптимизации, и способы их решения.

This article is intended for personal use only.

Скачать ALGLIB