Документация

Методы доверительной области для нелинейной минимизации

Многие методы, используемые в решателях Optimization Toolbox™, основаны на областях доверия, простой, но мощной концепции оптимизации.

Чтобы понять подход доверительной области к оптимизации, рассмотрим задачу неограниченной минимизации, минимизируйте f (x), где функция принимает векторные аргументы и возвращает скаляры. Предположим, что вы находитесь в точке x в n-пространстве и хотите улучшить, т.е. переместитесь в точку с меньшим значением функции. Основная идея состоит в том, чтобы аппроксимировать f более простой функцией q, которая разумно отражает поведение функции f в окрестности N вокруг точки x. Эта окрестность является областью доверия. Пробный этап s вычисляется путем минимизации (или приблизительно минимизации) над N. Это подпроблема области доверия,

Текущая точка обновляется до x + s, если f ( x + s ) < f (x); в противном случае текущая точка остается неизменной, и N, область доверия, уменьшается, и вычисление пробного этапа повторяется.

Ключевые вопросы при определении конкретного подхода области доверия к минимизации f (x) заключаются в том, как выбрать и вычислить аппроксимацию q (определенную в текущей точке x), как выбрать и изменить область доверия N и как точно решить подпроблему области доверия. В этом разделе рассматривается неограниченная проблема. В последующих разделах обсуждаются дополнительные осложнения из-за наличия ограничений на переменные.

В стандартном способе доверительной области ([48]) квадратичное приближение q определяется первыми двумя членами приближения Тейлора к F при x; окрестность N обычно имеет сферическую или эллипсоидальную форму. Математически обычно указывается подпроблема области доверия

где g - градиент f в текущей точке x, H - гессеновская матрица (симметричная матрица вторых производных), D - диагональная масштабная матрица, Δ - положительный скаляр, а ∥. ∥ - 2-норма. Существуют хорошие алгоритмы для решения уравнения 2 (см. [48]); такие алгоритмы обычно включают вычисление всех собственных значений H и процесс Ньютона, применяемый к светскому уравнению

Такие алгоритмы обеспечивают точное решение уравнения 2. Однако они требуют времени, пропорционального нескольким факторизациям Х. Поэтому для масштабных проблем необходим другой подход. Несколько аппроксимационных и эвристических стратегий, основанных на уравнении 2, были предложены в литературе ([42] и [50]). Подход аппроксимации, используемый в решателях Optimization Toolbox, заключается в ограничении подпроблемы области доверия двумерным подпространством S ([39] и [42]). После вычисления подпространства S работа по решению уравнения 2 является тривиальной, даже если необходима полная информация о собственном значении/собственном векторе (поскольку в подпространстве задача является только двумерной). Доминирующая работа теперь перешла к определению подпространства.

Двумерное подпространство S определяют с помощью процесса предварительно кондиционированного сопряженного градиента, описанного ниже. Решатель определяет S как линейное пространство, охватываемое s1 и s2, где s1 находится в направлении градиента g, а s2 является либо приближенным направлением Ньютона, т.е. решением

или направление отрицательной кривизны,

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

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

Эти четыре шага повторяются до сходимости. Измерение Δ области доверия корректируется в соответствии со стандартными правилами. В частности, он уменьшается, если стадия испытания не принимается, то есть f (x + s) ≥ f (x). Для обсуждения этого аспекта см. [46] и [49].

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

Метод предварительно кондиционированного сопряженного градиента

Популярным способом решения больших, симметричных, положительных определённых систем линейных уравнений Hp = -g является метод Предварительно кондиционированных сопряжённых градиентов (PCG). Этот итеративный подход требует возможности вычисления матрично-векторных произведений формы H· v, где v - произвольный вектор. Симметричная положительная определенная матрица М является предварительным условием для Н. То есть  М = ^C2, где ^C-1HC-1 является хорошо кондиционированной матрицей или матрицей с кластеризованными собственными значениями .

В контексте минимизации можно предположить, что матрица Гессена H симметрична. Однако H гарантированно является положительным определенным только в районе сильного минимизатора. Алгоритм PCG выходит, когда он сталкивается с направлением отрицательной (или нулевой) кривизны, то есть ^dTHd ≤ 0. Направление р выхода PCG является либо направлением отрицательной кривизны, либо приближенным решением для системы Ньютона Hp = -g. В любом случае p помогает определить двумерное подпространство, используемое в подходе «доверительная область», обсуждаемом в документе «Методы доверительной области для нелинейной минимизации».

Ограничения линейного равенства

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

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

где A - матрица m-на-n (m ≤ n). Некоторые решатели Optimization Toolbox предварительно обрабатывают А для удаления строгих линейных зависимостей методом, основанным на факторизации LU ^AT [46]. Здесь A принимается рангом m.

Метод, используемый для решения уравнения 5, отличается от неограниченного подхода двумя значимыми способами. Сначала вычисляется начальная осуществимая точка x0 с использованием разреженного шага наименьших квадратов, так что Ax0  = b. Во-вторых, алгоритм PCG заменяется на уменьшенные предварительно кондиционированные сопряженные градиенты (RPCG), см. [46], чтобы вычислить приблизительный уменьшенный шаг Ньютона (или направление отрицательной кривизны в нулевом пространстве A). Ключевой шаг линейной алгебры включает в себя решение систем вида

где $$\stackrel{A˜}{}$$ аппроксимирует A (малые ненулевые значения A устанавливаются в ноль при условии, что ранг не теряется), а C является разреженным симметричным положительно-определенным приближением к H, т.е. C = H. Подробнее см. [46].

Ограничения коробки

Проблема, связанная с рамкой, имеет вид

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

Масштабированный модифицированный шаг Ньютона возникает при изучении необходимых условий Куна-Такера для уравнения 7,

где

и вектор v (x) определяется ниже для каждого  1 ≤  i ≤ n:

Уравнение 8 нелинейной системы не дифференцируется везде. Недифференцируемость возникает, когда _vi = 0. Избежать таких моментов можно, сохранив строгую выполнимость, т.е. ограничив l < x < u.

Масштабированный модифицированный шаг Ньютона _sk для нелинейной системы уравнений, заданный уравнением 8, определяется как решение линейной системы.

на k-ой итерации, где

и

Здесь ^Сп играет роль якобианца | v |. Каждая диагональная составляющая диагональной матрицы Jv равна 0, -1 или 1. Если все компоненты l и u являются конечными,  Jv = diag (знак (g)). В точке, _где  gi = ₀, vi может быть не $${}_{\mathrm{дифференцируемым}.}^{}\mathrm{}$$В такой точке определяется Jiiv = 0. Недифференцируемость этого типа не вызывает беспокойства, поскольку для такого компонента не является существенным, какое _{значение} принимает vi. _{Кроме того}, | vi | все еще будет прерывистым в этой точке, но _{функция} | vi |· gi непрерывна.

Во-вторых, отражения используются для увеличения размера шага. (Один) этап отражения определяется следующим образом. Учитывая шаг p, который пересекает ограничивающее ограничение, рассмотрим первое ограничивающее ограничение, пересекаемое p; предположим, что это ограничение i-ой границы (i-я верхняя или i-я нижняя граница). Затем шаг отражения ^pR = p, за исключением i-ой составляющей, где ^pRi  = -_pi.

Введение

При ограниченной оптимизации общая цель состоит в преобразовании проблемы в более простую подпроблему, которая затем может быть решена и использована как основа итеративного процесса. Характеристикой большого класса ранних методов является преобразование ограниченной задачи в основную неограниченную задачу путем использования штрафной функции для ограничений, которые находятся вблизи или за пределами границы ограничения. Таким образом, ограниченная задача решается с помощью последовательности параметризованных неограниченных оптимизаций, которые в пределе (последовательности) сходятся к ограниченной задаче. Эти методы в настоящее время считаются относительно неэффективными и были заменены методами, которые были сосредоточены на решении уравнений Каруша-Куна-Такера (KKT). Уравнения KKT являются необходимыми условиями оптимальности для задачи ограниченной оптимизации. Если задача является так называемой выпуклой задачей программирования, то есть f (x) и _Gi (x), i = 1,..., m, являются выпуклыми функциями, то уравнения KKT являются как необходимыми, так и достаточными для глобальной точки решения.

Ссылаясь на GP (уравнение 1), уравнения Куна-Такера могут быть указаны как

в дополнение к исходным ограничениям в уравнении 1.

Первое уравнение описывает отмену градиентов между целевой функцией и активными ограничениями в точке решения. Для отмены градиентов необходимы множители Лагранжа (_{λ i}, i = 1,..., m), чтобы сбалансировать отклонения по величине целевой функции и градиентов ограничения. Поскольку в эту операцию отмены включены только активные ограничения, неактивные ограничения не должны включаться в эту операцию, и поэтому для них заданы множители Лагранжа, равные 0. Это косвенно указано в двух последних уравнениях Куна-Такера.

Решение уравнений ККТ составляет основу для многих алгоритмов нелинейного программирования. Эти алгоритмы пытаются вычислить множители Лагранжа напрямую. Ограниченные квазиньютоновские методы гарантируют сверхлинейную сходимость, накапливая информацию второго порядка относительно уравнений KKT, используя процедуру квазиньютоновского обновления. Эти методы обычно называют методами последовательного квадратичного программирования (SQP), поскольку подпроблема QP решается при каждой крупной итерации (также известной как методы итеративного квадратичного программирования, рекурсивного квадратичного программирования и ограниченной переменной метрики).

'active-set' алгоритм не является масштабным; см. Крупномасштабные и среднемасштабные алгоритмы.

Последовательное квадратичное программирование (SQP)

Методы SQP представляют собой современные методы нелинейного программирования. Шиттковский [36], например, реализовал и протестировал версию, которая превосходит любой другой проверенный метод по эффективности, точности и проценту успешных решений, по большому количеству проблем тестирования.

Основываясь на работах Биггса [1], Хана [22] и Пауэлла ([32] и [33]), метод позволяет близко имитировать метод Ньютона для ограниченной оптимизации так же, как это делается для неограниченной оптимизации. При каждой большой итерации делается аппроксимация гессена лагранжевой функции методом квази-ньютоновского обновления. Затем используется для генерации подпроблемы QP, решение которой используется для формирования направления поиска для процедуры поиска строки. Обзор SQP найден у Fletcher [13], Gill et al. [19], Пауэлл [35] и Шиттковский [23]. Общий способ, однако, изложен здесь.

Учитывая описание задачи в GP (уравнение 1), основная идея состоит в формулировании подпроблемы QP, основанной на квадратичной аппроксимации лагранжевой функции.

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

Подпроблема квадратного программирования (QP)

Эта подпроблема может быть решена с помощью любого алгоритма QP (см., например, Quadratic Programming Solution). Решение используется для формирования новой итерации

Параметр α k длины шага определяется соответствующей процедурой поиска строки, так что получается достаточное уменьшение функции качества (см. Обновление матрицы Гессена). Матрица Hk является положительной определённой аппроксимацией гессенской матрицы лагранжевой функции (уравнение 13). Hk может обновляться любым из квазиньютоновских методов, хотя метод BFGS (см. Обновление гессенской матрицы) представляется наиболее популярным.

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

Рассмотрим функцию Розенброка с дополнительным нелинейным ограничением неравенства, g (x),

Это было решено реализацией SQP в 96 итерациях по сравнению со 140 для неограниченного случая. На рис. 5-3, метод SQP для нелинейно ограниченной функции Розенброка показывает путь к точке решения x = [0,9072,0,8228], начиная с x = [-1,9,2,0].

Реализация SQP

Реализация SQP состоит из трех основных этапов, которые кратко рассматриваются в следующих подразделах:

• Обновление матрицы Гессена

• Решение для квадратного программирования

• Инициализация

• Функция поиска и оценки строк

Обновление Гессенской матрицы.  При каждой большой итерации положительное определённое квази-ньютоново приближение гессена лагранжевой функции, H, вычисляется с помощью метода BFGS, где _{λ i}, i = 1,..., m - оценка множителей Лагранжа .

где

Пауэлл [33] рекомендует сохранять гессенский позитив определенным, даже если он может быть положительным неопределенным в точке решения. Поддерживается положительный определенный гессен при условии, что $${}_{}^{}{}_{\mathrm{qkTsk}}$$ является положительным при каждом обновлении, и что H инициализируется положительной определенной матрицей. Если $${}_{}^{}{}_{\mathrm{qkTsk}}$$ не является положительным, _qk изменяется на поэлементной основе так, что $${}_{}^{}{}_{\mathrm{qkTsk}}\mathrm{>}$$0. Общая цель этой модификации состоит в том, чтобы как можно меньше искажать элементы qk, которые способствуют позитивному определенному обновлению. Поэтому в начальной фазе модификации самый отрицательный элемент qk * sk многократно уменьшается вдвое. Эта процедура продолжается $${}_{}^{\mathrm{до\; тех\; пор,\; пока}}{}_{}$$qkTsk не будет больше или равен малому отрицательному допуску. Если после этой $${}_{\mathrm{процедуры}}^{}{}_{}$$qkTsk все еще не является положительным, _{измените} qk, добавив вектор v, умноженный на постоянный скаляр w, то есть

где

и систематически увеличивать w до тех пор, пока $${}_{}^{}{}_{\mathrm{qkTsk}}$$ не станет положительным.

Функции fmincon, fminimax, fgoalattain, и fseminf все используют SQP. Если Display имеет значение 'iter' в options, то приводятся различные сведения, такие как значения функций и максимальное нарушение ограничения. Когда гессен должен быть изменен с использованием первой фазы предыдущей процедуры, чтобы сохранить его положительным определенным, то Hessian modified отображается. Если гессен должен быть снова модифицирован с использованием второй фазы описанного выше подхода, то Hessian modified twice отображается. Если подпроблема QP неосуществима, то infeasible отображается. Такие дисплеи обычно не вызывают беспокойства, но указывают на то, что проблема в высшей степени нелинейна и что сходимость может занять больше времени, чем обычно. Иногда сообщение no update отображается, указывая, что $${}_{}^{}{}_{\mathrm{qkTsk}}$$ почти равен нулю. Это может указывать на неправильную настройку проблемы или на то, что вы пытаетесь минимизировать неслагаемую функцию.

Решение для квадратного программирования.  При каждой основной итерации метода SQP решается задача QP следующей формы, где _Ai ссылается на iтретья строка матрицы А «m-на-n».

Метод, используемый в функциях Optimization Toolbox, представляет собой стратегию активного набора (также известную как метод проекции), аналогичную стратегии Gill et al., описанной в [18] и [17]. Он был изменен для проблем линейного программирования (LP) и квадратичного программирования (QP).

Процедура решения включает две фазы. Первая фаза включает вычисление возможной точки (если она существует). Вторая фаза включает в себя генерацию итеративной последовательности возможных точек, которые сходятся к решению. В этом способе $${\stackrel{}{}}_{}$$поддерживается активное множество A _ k, которое является оценкой активных ограничений (то есть тех, которые находятся на границах ограничения) в точке решения. Практически все алгоритмы QP являются методами активного набора. Этот момент подчеркивается, поскольку существует много различных методов, которые очень похожи по структуре, но которые описаны в широко различных терминах.

$${\stackrel{}{В}}_{}$$каждой итерации k обновляется k, и это используется для формирования основы для $${\stackrel{}{\mathrm{направления}}}_d$$^ k поиска. Ограничения равенства всегда остаются в активном $${\stackrel{}{}}_{\mathrm{наборе}}$$A _ k. Обозначение переменной $${\stackrel{}{}}_{}$$d ^ k используется здесь, чтобы отличить _ее от dk в основных итерациях метода SQP. $${\stackrel{}{}}_{}$$Направление поиска d ^ k вычисляется и минимизирует целевую функцию, оставаясь на любых активных границах ограничения. Возможное $${\stackrel{}{}}_{}$$подпространство для d ^ k формируется _из базиса Zk, столбцы которого ортогональны оценке $${\stackrel{}{}}_{}$$активного множества $${\stackrel{A}{}}_{}{\_}_{}k\mathrm{}$$(то есть A _ k = 0). Таким образом, направление поиска, которое формируется из линейного суммирования любой комбинации столбцов Zk, гарантированно остается на границах активных ограничений.

Матрица _Zk формируется из последних m-1  столбцов QR-разложения матрицы $${\stackrel{}{A}}_{\_}^{}$$kT, где l - число активных ограничений и l < m То есть Zk задается как

где

Как только _Zk найден, $${\stackrel{}{}}_{}$$ищут новое направление поиска d ^ k, которое минимизирует q (d), $${\stackrel{где}{}}_{}$$d ^ k находится в нулевом пространстве активных ограничений. То $${\stackrel{}{}}_{\mathrm{есть}}$$ d ^ k является линейной комбинацией _{столбцов} $${\stackrel{}{}}_{}{\mathrm{Zk}}_{:}$$d ^ k = Zkp для некоторого вектора p.

Тогда, если вы рассматриваете квадратичный как функцию p, заменяя d $${\stackrel{}{}}_{^}$$k, вы имеете

Дифференцирование этого относительно выходов p

∇q (p) называется спроецированным градиентом квадратичной функции, поскольку он является градиентом, спроецированным в подпространстве, определенном _Zk. Термин $${}_{}^{}{}_{}$$ZkTHZk называется проектируемым гессенским. Предполагая, что матрица Гессена H является положительной определённой (что и происходит в данной реализации SQP), то минимум функции q (p) в подпространстве, определяемом Zk, возникает, когда ∇q (p ) = 0, что является решением системы линейных уравнений

Затем выполняется шаг формы.

В каждой итерации из-за квадратичной природы целевой функции существует только два выбора длины шага α. Шаг единицы вдоль $${\stackrel{}{d}}_{^}$$ k является точным шагом к минимуму функции, ограниченной нулевым пространством $${\stackrel{}{}}_{\mathrm{A\; \_}}$$k. Если такой шаг может быть сделан без нарушения ограничений, то это решение QP ( уравнение 18). В противном случае шаг $${\stackrel{}{}}_{\mathrm{от}}$$ d ^ k до ближайшего ограничения меньше единицы, и новое ограничение включается в активный набор при следующей итерации. Расстояние до границ ограничения в любом $${\stackrel{}{}}_{}$$направлении d ^ k задается

которая определена для ограничений, не входящих в активный набор, и где направление $${\stackrel{}{d}}_{^}$$ k направлено к границе ограничения, т.е. $${}_{}{\stackrel{}{}}_{\mathrm{Aid}}\mathrm{^}\text{ k}>\mathrm{}0\mathrm{,\; i}=$$1,..., m.

Когда в активное множество включены n независимых ограничений, без расположения минимума, вычисляются множители Лагранжа λ _k, удовлетворяющие неингулярному множеству линейных уравнений

Если все элементы _{λ k} являются положительными, _xk является оптимальным решением QP (уравнение 18). Однако если какая-либо составляющая _{λ k} отрицательна, и составляющая не соответствует ограничению равенства, то соответствующий элемент удаляется из активного множества и ищется новый итерат.

Инициализация.  Для запуска алгоритма требуется выполнимая точка. Если текущая точка из метода SQP неосуществима, то можно найти точку, решив задачу линейного программирования

Обозначение _Ai указывает i-ю строку матрицы A. Можно найти осуществимую точку (если она существует) уравнения 26, установив значение x, которое удовлетворяет ограничениям равенства. Это значение можно определить, решив недо- или переопределенный набор линейных уравнений, сформированных из набора ограничений равенства. Если есть решение этой задачи, то переменная ослабления γ устанавливается на ограничение максимального неравенства в этой точке.

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

Если с помощью предшествующего метода LP найдена приемлемая точка, вводится основная фаза QP. Направление поиска $${\stackrel{}{d}}_{^}$$ k инициализируется с