Алгоритмы квадратичного программирования

Определение квадратичного программирования

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

\min_{x} \frac{1}{2} x^{T} H x + c^{T} x

(1)

таким образом, что A·x ≤ b, Aeq·x = beq, l ≤ x ≤ u.

`interior-point-convex` `quadprog` Алгоритм

interior-point-convex алгоритм выполняет следующие шаги:

Предварительно решить/Пострешить
Сгенерируйте начальную точку
Корректор предиктора
Остановка условий
Обнаружение недопустимости

Примечание

Алгоритм имеет два пути выполнения кода. Это берет тот, когда матрица Гессиана, из которой H является обычной (полной) матрицей, удваивается, и это берет другой, когда H является разреженной матрицей. Для получения дополнительной информации разреженного типа данных, смотрите Разреженные матрицы (MATLAB). Обычно алгоритм быстрее для больших проблем, которые имеют относительно немного ненулевых условий, когда вы задаете H как sparse. Точно так же алгоритм быстрее для небольших или относительно плотных проблем, когда вы задаете H как full.

Предварительно решить/Пострешить

Алгоритм сначала пытается упростить проблему путем удаления сокращения и упрощения ограничений. Задачи, выполняемые во время предварительно решить шага, могут включать следующее:

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

Если алгоритм обнаруживает неосуществимую или неограниченную проблему, он останавливает и выпускает соответствующее выходное сообщение.

Алгоритм может прибыть в одну допустимую точку, которая представляет решение.

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

Для получения дополнительной информации смотрите Гулда и Тойнта [63].

Сгенерируйте начальную точку

Начальная точка x0 поскольку алгоритм:

Инициализируйте x0 к ones(n,1), где n количество строк в H.
Для компонентов, которые имеют обоих верхняя граница ub и нижняя граница lb, если компонент x0 не строго в границах, компонент установлен в (ub + lb)/2.
Для компонентов, которые имеют только один связанный, измените компонент при необходимости, чтобы лечь строго в связанном.
Сделайте шаг предиктора (см. Корректор Предиктора), с незначительными коррекциями для выполнимости, не полным шагом корректора предиктора. Это помещает начальную точку ближе в central path, не влекущий за собой издержки полного шага корректора предиктора. Для получения дополнительной информации центрального пути, смотрите Носедэла и Райта [7], страница 397.

Корректор предиктора

Разреженные и полные выпуклые внутренней точкой алгоритмы отличаются в основном по фазе корректора предиктора. Алгоритмы подобны, но отличаются по некоторым деталям. Для описания основного алгоритма смотрите Mehrotra [47].

Алгоритмы начинаются путем превращения линейных неравенств Ax <= b в неравенства формы Ax> = b путем умножения A и b-1. Это не имеет никакого влияния на решение, но делает проблему той же формы найденной в некоторой литературе.

Разреженный корректор предиктора
Полный корректор предиктора

Разреженный корректор предиктора. Подобно fmincon алгоритм внутренней точки, разреженный interior-point-convex алгоритм пытается найти точку, где условия Karush-Kuhn-Tucker (KKT) содержат. Для проблемы квадратичного программирования, описанной в Определении Квадратичного программирования, эти условия:

$\begin{matrix} H x + c - A_{e q}^{T} y - {\bar{A}}^{T} z = 0 \\ \bar{A} x - \bar{b} - s = 0 \\ A_{e q} x - b_{e q} = 0 \\ s_{i} z_{i} = 0, i = 1, 2, ..., m \\ s \geq 0 \\ z \geq 0. \end{matrix}$

Здесь

$\bar{A}$ расширенная линейная матрица неравенства, которая включает границы, записанные как линейные неравенства. $\bar{b}$ соответствующий линейный вектор неравенства, включая границы.
s является вектором слаксов, которые преобразуют ограничения неравенства в равенства. s имеет длину m, количество линейных неравенств и границ.
z является вектором множителей Лагранжа, соответствующих s.
y является вектором множителей Лагранжа, сопоставленных с ограничениями равенства.

Алгоритм сначала предсказывает шаг от формулы Ньютона-Raphson, затем вычисляет шаг корректора. Корректор пытается лучше осуществить нелинейное ограничение _sizi = 0.

Определения для шага предиктора:

_rd, двойная невязка:
$r_{d} = H x + c - A_{e q}^{T} y - {\bar{A}}^{T} z .$
_req, основная невязка ограничения равенства:
$r_{e q} = A_{e q} x - b_{e q} .$
_rineq, основная невязка ограничения неравенства, которая включает границы и слаксы:
$r_{i n e q} = \bar{A} x - \bar{b} - s .$
_rsz, невязка взаимозависимости:
_rsz = Sz.
S является диагональной матрицей слабых условий, z является матрицей столбца множителей Лагранжа.
_rc, средняя взаимозависимость:
$r_{c} = \frac{s^{T} z}{m} .$

На шаге Ньютона, изменениях в x, s, y и z, дают:

(\begin{matrix} H & 0 & - A_{e q}^{T} & - {\bar{A}}^{T} \\ A_{e q} & 0 & 0 & 0 \\ \bar{A} & - I & 0 & 0 \\ 0 & Z & 0 & S \end{matrix}) (\begin{matrix} Δ x \\ Δ s \\ Δ y \\ Δ z \end{matrix}) = - (\begin{matrix} r_{d} \\ r_{e q} \\ r_{i n e q} \\ r_{s z} \end{matrix}) .

(2)

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

Кроме того, чтобы поддержать положение “в центре” во внутренней части, вместо того, чтобы пытаться решить _sizi = 0, алгоритм берет положительный параметр σ и пытается решить

_sizi = _σrc.

quadprog _rsz замен в уравнении шага Ньютона с _rsz + ΔsΔz – _σrc 1, где 1 вектор из единиц. Кроме того, quadprog переупорядочивает уравнения Ньютона, чтобы получить симметричное, более численно устойчивую систему для вычисления шага предиктора.

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

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

$\bar{A}$ расширенная линейная матрица, которая включает и линейные неравенства и линейные равенства. $\bar{b}$ соответствующий линейный вектор равенства. $\bar{A}$ также включает условия для расширения векторного x со слабыми переменными s, которые поворачивают ограничения неравенства к ограничениям равенства:

$\bar{A} x = (\begin{matrix} A_{e q} & 0 \\ A & I \end{matrix}) (\begin{matrix} x_{0} \\ s \end{matrix}),$

где x ₀ средних значений исходный вектор x.

Условия KKT

\begin{matrix} H x + c - {\bar{A}}^{T} y - v + w = 0 \\ \bar{A} x = \bar{b} \\ x + t = u \\ v_{i} x_{i} = 0, i = 1, 2, ..., m \\ w_{i} t_{i} = 0, i = 1, 2, ..., n \\ x, v, w, t \geq 0. \end{matrix}

(3)

Чтобы найти решение x, слабыми переменными и двойными переменными к уравнению 3, алгоритм в основном рассматривает шаг Ньютона-Raphson:

(\begin{matrix} H & - {\bar{A}}^{T} & 0 & - I & I \\ \bar{A} & 0 & 0 & 0 & 0 \\ - I & 0 & - I & 0 & 0 \\ V & 0 & 0 & X & 0 \\ 0 & 0 & W & 0 & T \end{matrix}) (\begin{matrix} Δ x \\ Δ y \\ Δ t \\ \begin{array}{l} Δ v \\ Δ w \end{array} \end{matrix}) = - (\begin{matrix} H x + c - {\bar{A}}^{T} y - v + w \\ \bar{A} x - \bar{b} \\ u - x - t \\ \begin{array}{l} V X \\ W T \end{array} \end{matrix}) = - (\begin{matrix} r_{d} \\ r_{p} \\ r_{u b} \\ \begin{array}{l} r_{v x} \\ r_{w t} \end{array} \end{matrix}),

(4)

где X, V, W и T являются диагональными матрицами, соответствующими векторам x, v, w и t соответственно. Вектора невязок на ультраправой стороне уравнения:

_rd, двойная невязка
_rp, основная невязка
_rub, невязка верхней границы
_rvx, невязка взаимозависимости нижней границы
_rwt, невязка взаимозависимости верхней границы

Алгоритм решает уравнение 4 первым преобразованием его к форме симметрической матрицы

(\begin{matrix} - D & {\bar{A}}^{T} \\ \bar{A} & 0 \end{matrix}) (\begin{matrix} Δ x \\ Δ y \end{matrix}) = - (\begin{matrix} R \\ r_{p} \end{matrix}),

(5)

где

$\begin{array}{l} D = H + X^{- 1} V + T^{- 1} W \\ R = - r_{d} - X^{- 1} r_{v x} + T^{- 1} r_{w t} + T^{- 1} W r_{u b} . \end{array}$

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

Чтобы вывести уравнение 5 из уравнения 4, заметьте, что вторая строка уравнения 5 совпадает со второй строкой матрицы уравнения 4. Первая строка уравнения 5 прибывает из решения последних двух строк уравнения 4 для Δv и Δw и затем решения для Δt.

Чтобы решить уравнение 5, алгоритм следует за существенными элементами Олтмена и Гондсио [1]. Алгоритм решает симметричную систему разложением LDL. Как указано авторами, такими как Вэндербеи и Карпентер [2], это разложение численно устойчиво без любого поворота, быть быстрым - также.

Полный quadprog алгоритм корректора предиктора является в основном тем же самым как этим в linprog 'interior-point' алгоритм, но включает квадратичные условия также. Смотрите Корректор Предиктора.

Ссылки

[1] Олтмен, Анна и Х. Гондсио. Упорядоченные симметричные неопределенные системы в методах внутренней точки для линейной и квадратичной оптимизации. Методы оптимизации и программное обеспечение, 1999. Доступный для скачивания здесь.

[2] Vanderbei, R. J. и Т. Дж. Карпентер. Симметричные неопределенные системы для методов внутренней точки. Математическое программирование 58, 1993. стр 1–32. Доступный для скачивания здесь.

Остановка условий

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

$ρ = \max (1, ‖ H ‖, ‖ \bar{A} ‖, ‖ A_{e q} ‖, ‖ c ‖, ‖ \bar{b} ‖, ‖ b_{e q} ‖)$

Алгоритм останавливается, когда всем этим условиям удовлетворяют:

$\begin{matrix} {‖ r_{p} ‖}_{1} + {‖ r_{u b} ‖}_{1} \leq ρ TolCon \\ {‖ r_{d} ‖}_{\infty} \leq ρ TolFun \\ r_{c} \leq TolFun, \end{matrix}$

где

$r_{c} = \max_{i} (\min (| x_{i} v_{i} |, | x_{i} |, | v_{i} |), \min (| t_{i} w_{i} |, | t_{i} |, | w_{i} |)) .$

_rc по существу измеряет размер остаточных значений взаимозависимости xv и tw, которые являются каждым нулевые вектора в решении.

Обнаружение недопустимости

quadprog вычисляет оценочную функцию φ в каждой итерации. Оценочная функция является мерой выполнимости. quadprog остановки, если оценочная функция становится слишком большой. В этом случае, quadprog объявляет проблему быть неосуществимым.

Оценочная функция связана с условиями KKT для проблемы — смотрите Корректор Предиктора. Используйте следующие определения:

$\begin{matrix} ρ = \max (1, ‖ H ‖, ‖ \bar{A} ‖, ‖ A_{e q} ‖, ‖ c ‖, ‖ \bar{b} ‖, ‖ b_{e q} ‖) \\ r_{eq} = A_{eq} x - b_{eq} \\ r_{ineq} = \bar{A} x - \bar{b} - s \\ r_{d} = H x + c + A_{eq}^{T} λ_{eq} + {\bar{A}}^{T} {\bar{λ}}_{ineq} \\ g = \frac{1}{2} x^{T} H x + c^{T} x - {\bar{b}}^{T} {\bar{λ}}_{ineq} - b_{eq}^{T} λ_{eq} . \end{matrix}$

Обозначение $\bar{A}$ и $\bar{b}$ означает линейные коэффициенты неравенства, увеличенные с условиями представлять границы для разреженного алгоритма. Обозначение ${\bar{λ}}_{ineq}$ так же представляет множители Лагранжа для линейных ограничений неравенства, включая связанные ограничения. Это было названо z в Корректоре Предиктора, и $λ_{eq}$ был назван y.

Оценочная функция φ

$\frac{1}{ρ} (\max ({‖ r_{eq} ‖}_{\infty}, {‖ r_{ineq} ‖}_{\infty}, {‖ r_{d} ‖}_{\infty}) + g) .$

Если эта оценочная функция становится слишком большой, quadprog объявляет проблему быть неосуществимым и остановы с выходным флагом -2.

`trust-region-reflective` `quadprog` Алгоритм

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

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

\min_{s} {q (s), s \in N} .

(6)

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

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

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

\min {\frac{1}{2} s^{T} H s + s^{T} g таким образом , что ‖ D s ‖ \leq Δ},

(7)

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

$\frac{1}{Δ} - \frac{1}{‖ s ‖} = 0.$

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

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

H \cdot s_{2} = - g,

(8)

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

s_{2}^{T} \cdot H \cdot s_{2} < 0.

(9)

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

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

Сформулируйте двумерную подпроблему доверительной области.
Решите уравнение 7, чтобы определить испытательный шаг s.
Если f (x + s) <f (x), то x = x + s.
Настройте Δ.

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

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

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

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

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

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

Линейные ограничения равенства

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

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

\min {f (x) such that A x = b},

(10)

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

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

[\begin{matrix} C & {\tilde{A}}^{T} \\ \tilde{A} & 0 \end{matrix}] [\begin{matrix} s \\ t \end{matrix}] = [\begin{matrix} r \\ 0 \end{matrix}],

(11)

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

Ограничения поля

Ограниченная проблема поля имеет форму

\min {f (x) such that l \leq x \leq u},

(12)

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

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

{(D (x))}^{- 2} g = 0,

(13)

где

$D (x) = diag ({| v_{k} |}^{- 1 / 2}),$

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

Если _gi < 0 и _ui < ∞ затем _vi = _xi – _ui
Если _gi ≥ 0 и _li > – ∞ затем _vi = _xi – _li
Если _gi < 0 и _ui = ∞ затем _vi = –1
Если _gi ≥ 0 и _li = – ∞ затем _vi = 1

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

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

\hat{M} D s^{N} = - \hat{g}

(14)

в k th итерация, где

\hat{g} = D^{- 1} g = diag ({| v |}^{1 / 2}) g,

(15)

\hat{M} = D^{- 1} H D^{- 1} + diag (g) J^{v} .

(16)

Здесь ^Jv играет роль якобиана |v |. Каждый диагональный компонент диагонального матричного ^Jv равняется 0, –1, или 1. Если все компоненты l и u конечны, ^Jv = diag (знак (g)). В точке, где _gi = 0, _vi не может быть дифференцируемым. $J_{i i}^{v} = 0$ задан в такой точке. Недифференцируемость этого типа не является поводом для беспокойства, потому что для такого компонента это не значительно, какое значение _vi принимает. Далее, |_vi | все еще будет прерывисто в этой точке, но функции _|vi| · _gi непрерывен.

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

`active-set` `quadprog` Алгоритм

После завершения предварительно решить шага, active-set алгоритм продолжает в двух фазах.

Фаза 1 — Получает допустимую точку относительно всех ограничений.
Фаза 2 — Итеративно понижает целевую функцию при ведении списка активных ограничений и поддержании выполнимости в каждой итерации.

active-set стратегия (также известный как метод проекции) похожа на стратегию Джилла и др., описанная в [18] и [17].

Предварительно решите шаг

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

Алгоритм может прибыть в одну допустимую точку, которая представляет решение.

Для получения дополнительной информации смотрите Гулда и Тойнта [63].

Фаза 1 Алгоритм

В Фазе 1 алгоритм пытается найти точку x это удовлетворяет всем ограничениям без фактора целевой функции. quadprog запускает алгоритм Фазы 1 только если предоставленная начальная точка x0 неосуществимо.

Чтобы начаться, алгоритм пытается найти точку, которая выполнима относительно всех ограничений равенства, такова как X = Aeq\beq. Если нет никакого решения x к уравнениям Aeq*x = beq, затем остановы алгоритма. Если существует решение X, следующий шаг должен удовлетворить границам и линейным неравенствам. В случае никаких ограничений равенства устанавливает X = x0, начальная точка.

Запуск с X, алгоритм вычисляет M = max(A*X – b, X - ub, lb – X). Если M <= options.ConstraintTolerance, затем точка X выполнимо и остановы алгоритма Фазы 1.

Если M > options.ConstraintTolerance, алгоритм вводит неотрицательную слабую переменную γ для вспомогательной задачи линейного программирования

$\min_{x, γ} γ$

таким образом, что

$\begin{matrix} A x - γ \leq b \\ Aeq x = beq \\ lb - γ \leq x \leq ub + γ \\ γ \geq - ρ . \end{matrix}$

Здесь, ρ является ConstraintTolerance опция, умноженная на абсолютное значение самого большого элемента в A и Aeq. Если алгоритм находит γ = 0 и точка x, который удовлетворяет уравнениям и неравенствам, то x является выполнимой точкой Фазы 1. Если нет никакого решения вспомогательной задачи линейного программирования x с γ = 0, то проблема Фазы 1 неосуществима.

Чтобы решить вспомогательную задачу линейного программирования, алгоритм устанавливает γ ₀ = M + 1, устанавливает x ₀ = X, и инициализирует активный набор как фиксированные переменные (если таковые имеются) и все ограничения равенства. Алгоритм переформулирует линейные переменные p программирования, чтобы быть смещением x от текущей точки x ₀, а именно, x = x ₀ + p. Алгоритм решает задачу линейного программирования теми же итерациями, как это берет в Фазе 2, чтобы решить задачу квадратичного программирования с соответственно модифицированным Гессианом.

Фаза 2 Алгоритм

В терминах переменной d проблема

\begin{matrix} \min_{d \in ℜ^{n}} q (d) = \frac{1}{2} d^{T} H d + c^{T} d, \\ A_{i} d = b_{i}, i = 1, ..., m_{e} \\ A_{i} d \leq b_{i}, i = m_{e} + 1, ..., m . \end{matrix}

(17)

Здесь, _Ai относится к iстрока th m-by-n матричный A.

Во время Фазы 2, активного набора ${\bar{A}}_{k}$ , который является оценкой активных ограничений (те на границах ограничений) в точке решения.

Обновления алгоритма ${\bar{A}}_{k}$ в каждой итерации k, формируя основание для поискового направления _dk. Ограничения равенства всегда остаются в активном наборе ${\bar{A}}_{k}$ . Поисковое направление _dk вычисляется и минимизирует целевую функцию при оставлении на любых активных границах ограничений. Алгоритм формирует выполнимое подпространство для _dk от основания _Zk, столбцы которого являются ортогональными к оценке активного набора ${\bar{A}}_{k}$ (то есть, ${\bar{A}}_{k} Z_{k} = 0$ ). Поэтому поисковое направление, которое формируется из линейного суммирования любой комбинации столбцов _Zk, как гарантируют, останется на контурах активных ограничений.

Алгоритм формирует матричный _Zk из последнего m – столбцы l разложения QR матрицы ${\bar{A}}_{k}^{T}$ , где l является количеством активных ограничений и l < m. Таким образом, _Zk дают

Z_{k} = Q [:, l + 1 : m],

(18)

где

$Q^{T} {\bar{A}}_{k}^{T} = [\begin{matrix} R \\ 0 \end{matrix}] .$

После нахождения _Zk алгоритм ищет новое поисковое направление _dk, который минимизирует q (d), где _dk находится на пустом пробеле активных ограничений. Таким образом, _dk является линейной комбинацией столбцов _Zk: ${\hat{d}}_{k} = Z_{k} p$ для некоторого векторного p.

При просмотре квадратичного как функцию p путем заменения _dk, дает

q (p) = \frac{1}{2} p^{T} Z_{k}^{T} H Z_{k} p + c^{T} Z_{k} p .

(19)

Дифференциация этого уравнения относительно урожаев p

\nabla q (p) = Z_{k}^{T} H Z_{k} p + Z_{k}^{T} c .

(20)

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

Z_{k}^{T} H Z_{k} p = - Z_{k}^{T} c .

(21)

Алгоритм затем получает шаг формы

$x_{k + 1} = x_{k} + α d_{k},$

где

$d_{k} = Z_{k} p .$

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

$α = \min_{i \in {1, ..., m}} {\frac{- (A_{i} x_{k} - b_{i})}{A_{i} d_{k}}},$

который задан для ограничений не в активном наборе, и где направление _dk находится к границе ограничений, то есть, $A_{i} d_{k} > 0, i = 1, ..., m$ .

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

{\bar{A}}_{k}^{T} λ_{k} = c + H x_{k} .

(22)

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

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

Примечание

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

Иногда, active-set алгоритм может испытать затруднения при обнаружении, когда проблема неограниченна. Эта проблема может произойти, если направление неограниченного v является направлением где квадратичный термин v'Hv = 0. Для числовой устойчивости и включить факторизацию Холесского, active-set алгоритм добавляет маленький, строго выпуклый срок в квадратичную цель. Этот маленький срок заставляет целевую функцию быть ограниченной далеко от –Inf. В этом случае, active-set алгоритм достигает предела итерации вместо того, чтобы сообщить, что решение неограниченно. Другими словами, алгоритм останавливается с выходным флагом 0 вместо –3.

Ссылки

[1] Жабры, P. E. В. Мюррей, М. А. Сондерс и М. Х. Райт. Практическая процедура антициклического повторения для линейно ограниченной оптимизации. Математика. Программирование 45 (1), август 1989, стр 437–474.

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