Двойственная функция
\begin{equation*} \begin{split} g(\mu) & = -b^{\top}\mu + \inf_x(f(x) + \mu^{\top}Ax) \\ & = -b^{\top}\mu - \sup_x((-A^{\top}\mu)^{\top}x -f(x)) \\ & = -b^{\top}\mu - f^*(-A^{\top}\mu) \end{split} \end{equation*}Двойственная задача
$$ \max_\mu -b^{\top}\mu - f^*(-A^{\top}\mu) $$Подход 1: найти сопряжённую функцию и решить безусловную задачу оптимизации
Трудности
или
$$ \begin{bmatrix} f' & A^{\top} \\ A & 0 \end{bmatrix} \begin{bmatrix} x^{\\*} \\ \mu^{\\*} \end{bmatrix} = \begin{bmatrix} 0 \\ b \end{bmatrix} $$Подход 2: решить нелинейную в общем случае систему методом Ньютона.
Вопрос: в каком случае система окажется линейной?
Из условий оптимальности имеем
$$ \begin{bmatrix} f''(x) & A^{\top} \\ A & 0 \end{bmatrix} \begin{bmatrix} v \\ w \end{bmatrix} = \begin{bmatrix} -f'(x) \\ 0 \end{bmatrix} $$Шаг метода Ньютона определён только для невырожденной матрицы!
Упражнение. Посчитайте за сколько итераций метод Ньютона сойдётся для квадратичной функции с ограничениями типа равенств.
Важно: начальная точка должна лежать в допустимом множестве!
def NewtonEqualityFeasible(f, gradf, hessf, A, b, stop_crit, line_search, x0, tol):
x = x0
n = x.shape[0]
while True:
newton_matrix = [[hessf(x), A.T], [A, 0]]
rhs = [-gradf(x), 0]
w = solve_lin_sys(newton_matrix, rhs)
h = w[:n]
if stop_crit(x, h, gradf(x), **kwargs) < tol:
break
alpha = line_search(x, h, f, gradf(x), **kwargs)
x = x + alpha * h
return x
Получим выражение для значения
$$ f(x) - \inf_v(\hat{f}(x + v) \; | \; A(x+v) = b), $$где $\hat{f}$ - квадратичная аппроксимация функции $f$.
Для этого
$$ \langle h^{\top} \rvert \cdot \quad f''(x)h + A^{\top}w = -f'(x) $$с учётом $Ah = 0$ получаем
$$ h^{\top}f''(x)h = -f'(x)^{\top}h $$Тогда
$$ \inf_v(\hat{f}(x + v) \; | \; A(x+v) = b) = f(x) - \frac{1}{2}h^{\top}f''(x)h $$Вывод: величина $h^{\top}f''(x)h$ является наиболее адекватным критерием остановки метода Ньютона.
Сходимость метода аналогична сходимости метода Ньютона для задачи безусловной оптимизации.
Теорема Пусть выполнены следующие условия
Тогда, метод Ньютона сходится к паре $(x^*, \mu^*)$ линейно, а при достижении достаточной близости к решению - квадратично.
где $r_p(x, \mu) = Ax - b$ и $r_d(x, \mu) = f'(x) + A^{\top}\mu$
или более подробно
$$ \begin{bmatrix} f''(x) & A^{\top}\\ A & 0 \end{bmatrix} \begin{bmatrix} z_p\\ z_d \end{bmatrix} = - \begin{bmatrix} r_p(x, \mu)\\ r_d(x, \mu) \end{bmatrix} = - \begin{bmatrix} f'(x) + A^{\top}\mu\\ Ax - b \end{bmatrix} $$def NewtonEqualityInfeasible(f, gradf, hessf, A, b, stop_crit, line_search, x0, mu0, tol):
x = x0
mu = mu0
n = x.shape[0]
while True:
z_p, z_d = ComputeNewtonStep(hessf(x), A, b)
if stop_crit(x, z_p, z_d, gradf(x), **kwargs) < tol:
break
alpha = line_search(x, z_p, z_d, f, gradf(x), **kwargs)
x = x + alpha * z_p
mu = mu + alpha * z_d
return x
def linesearch(r, x, mu, z_p, z_d, c, beta):
alpha = 1
while np.linalg.norm(r(x + alpha * z_p, mu + alpha * z_d)) >= (1 - c * alpha) * np.linalg.norm(r(x, mu)):
alpha *= beta
return alpha
Результат аналогичен результаты для допустимой начальной точки
Теорема. Пусть
Тогда сходимость метода линейна при удалении от решении и квадратичная при достаточном приближении к решению.
где $f_i$ - выпуклые и дважды непрерывно дифференцируемы, $A \in \mathbb{R}^{p \times n}$ и $\mathrm{rank} \; A = p < n$.
Предполагаем, что задача строго разрешима, то есть выполняется условие Слейтера.
где $I_-$ - индикаторная функция
$$ I_-(u) = \begin{cases} 0, & u \leq 0\\ \infty, & u > 0 \end{cases} $$Проблема. Теперь целевая функция - недифференцируема.
In [1]:
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
x = np.linspace(-2, 0, 100000, endpoint=False)
plt.figure(figsize=(10, 6))
for t in [0.1, 0.5, 1, 1.5, 2]:
plt.plot(x, -t * np.log(-x), label=r"$t = " + str(t) + "$")
plt.legend(fontsize=20)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlabel("u", fontsize=20)
Out[1]:
называется логарифмическим барьером. Её область определения - множество точек, для котороых ограничения типа неравенств выполняются строго.
Упражнение. Найдите градиент и гессиан $\phi(x)$
для $\lambda = \lambda^*(t)$ и $\mu = \mu^*$.
def BarrierMethod(f, x0, t0, tol, alpha, **kwargs):
x = x0
t = t0
while True:
x = SolveBarrierProblem(f, t, x, **kwargs)
if m * t < tol:
break
t *= alpha
return x
Простой поиск допустимой точки
\begin{equation*} \begin{split} & \min s\\ \text{s.t. } & f_i(x) \leq s\\ & Ax = b \end{split} \end{equation*}Похож на барьерный метод, но
При некоторых предположениях о начальной точке и начальном значении $\mu$, можно показать, что для достижения $\mu_k \leq \varepsilon$ потребуется
$$ \mathcal{O}\left(\sqrt{n}\log \left( \frac{1}{\varepsilon}\right)\right) $$итераций
Доказательство и все детали можно найти тут или тут