Также на странице курса.
где $S \subseteq \mathbb{R}^n$, $f_j: S \rightarrow \mathbb{R}, \; j = 0,\ldots,m$, $g_k: S \rightarrow \mathbb{R}, \; k=1,\ldots,p$
Все функции как минимум непрерывны.
Важный факт</span>: задачи нелинейной оптимизации
в их самой общей форме являются численно неразрешимыми!
если $x^*$ точка локального минимума дифференцируемой функции $f(x)$, тогда $$ f'(x^*) = 0 $$
если $x^*$ точка локального минимума дважды дифференцируемой функции $f(x)$, тогда
$$ f'(x^*) = 0 \quad \text{и} \quad f''(x^*) \succeq 0 $$пусть $f(x)$ дважды дифференцируемая функция, и пусть точка $x^*$ удовлетворяет условиям
$$ f'(x^*) = 0 \quad f''(x^*) \succ 0, $$тогда $x^*$ является точкой строго локального минимума функции $f(x)$.
Замечание: убедитесь, что Вы понимаете, как доказывать эти
результаты!
Но ведь $x^*$ неизвестна!
Тогда
\begin{align*} & \|x_{k+1} - x_k \| = \|x_{k+1} - x_k + x^* - x^* \| \leq \\ & \|x_{k+1} - x^* \| + \| x_k - x^* \| \leq 2\varepsilon \end{align*}Аналогично для сходимости по функции,
однако иногда можно оценить $f^*$!
Замечание: лучше использовать относительные изменения
этих величин!
Например $\dfrac{\|x_{k+1} - x_k \|_2}{\| x_k \|_2}$
Определение: оракулом называют некоторое абстрактное
устройство, которое отвечает на последовательные вопросы
метода
Аналогия из ООП:
Концепция чёрного ящика
Метод доверительных областей: фиксируется допустимый размер области по некоторой норме $\| \cdot \| \leq \alpha$ и модель целевой функции, которая хорошо её аппроксимирует в выбранной области.
Далее производится поиск направления $h_k$, минимизирующего модель целевой функции и не выводящего точку $x_k + h_k$ за пределы доверительной области
Для заданного класса задач сравнивают следующие величины:
Линейная (геометрическая прогрессия) $$ \| x_{k+1} - x^* \|_2 \leq Cq^k, $$ где $q \in (0, 1)$ и $ 0 < C < \infty$
Сверхлинейная $$ \| x_{k+1} - x^* \|_2 \leq Cq^{k^p}, $$ где $q \in (0, 1)$, $ 0 < C < \infty$ и $p > 1$
Оптимальным методам и нижним оценкам будет, возможно,
посвящён отдельный семинар или часть домашнего задания.
In [3]:
%matplotlib inline
import matplotlib.pyplot as plt
USE_COLAB = False
if not USE_COLAB:
plt.rc("text", usetex=True)
import numpy as np
C = 10
alpha = -0.5
q = 0.9
num_iter = 7
sublinear = np.array([C * k**alpha for k in range(1, num_iter + 1)])
linear = np.array([C * q**k for k in range(1, num_iter + 1)])
superlinear = np.array([C * q**(k**2) for k in range(1, num_iter + 1)])
quadratic = np.array([C * q**(2**k) for k in range(1, num_iter + 1)])
plt.figure(figsize=(12,8))
plt.semilogy(np.arange(1, num_iter+1), sublinear,
label=r"Sublinear, $\alpha = -0.5$")
plt.semilogy(np.arange(1, num_iter+1), superlinear,
label=r"Superlinear, $q = 0.5, p=2$")
plt.semilogy(np.arange(1, num_iter+1), linear,
label=r"Linear, $q = 0.5$")
plt.semilogy(np.arange(1, num_iter+1), quadratic,
label=r"Quadratic, $q = 0.5$")
plt.xlabel("Number of iterations, $k$", fontsize=28)
plt.ylabel("Error rate upper bound", fontsize=28)
plt.legend(loc="best", fontsize=26)
plt.xticks(fontsize = 28)
_ = plt.yticks(fontsize = 28)
Что дают теоремы сходимости
Что НЕ дают теоремы сходимости
Мораль: нужно проявлять разумную осторожность
и здравый смысл!
Методы нулевого порядка: оракул возвращает только значение функции $f(x)$
Методы первого порядка: оракул возвращает значение функции $f(x)$ и её градиент $f'(x)$
Методы второго порядка: оракул возвращает значение функции $f(x)$, её градиент $f'(x)$ и гессиан $f''(x)$.
Вопрос: существуют ли методы более высокого порядка?
Идея из информатики первого семестра:
делим отрезок $[a,b]$ на две равные части
пока не найдём минимум унимодальной функции.
Тогда $$ |x_{K+1} - x^*| \leq \frac{b_{K+1} - a_{K+1}}{2} = \left( \frac{1}{2} \right)^{\frac{N-1}{2}} (b - a) \approx 0.5^{K} (b - a) $$
In [5]:
def binary_search(f, a, b, epsilon, callback=None):
c = (a + b) / 2.0
while abs(b - a) > epsilon:
# Check left subsegment
y = (a + c) / 2.0
if f(y) <= f(c):
b = c
c = y
else:
# Check right subsegment
z = (b + c) / 2.0
if f(c) <= f(z):
a = y
b = z
else:
a = c
c = z
if callback is not None:
callback(a, b)
return c
In [6]:
def my_callback(a, b, left_bound, right_bound, approximation):
left_bound.append(a)
right_bound.append(b)
approximation.append((a + b) / 2.0)
In [8]:
import numpy as np
left_boud_bs = []
right_bound_bs = []
approximation_bs = []
callback_bs = lambda a, b: my_callback(a, b,
left_boud_bs, right_bound_bs, approximation_bs)
# Target unimodal function on given segment
f = lambda x: (x - 2) * x * (x + 2)**2 # np.power(x+2, 2)
# f = lambda x: -np.sin(x)
x_true = -2
# x_true = np.pi / 2.0
a = -3
b = -1.5
epsilon = 1e-8
x_opt = binary_search(f, a, b, epsilon, callback_bs)
print(np.abs(x_opt - x_true))
plt.figure(figsize=(10,6))
plt.plot(np.linspace(a,b), f(np.linspace(a,b)))
plt.title("Objective function", fontsize=28)
plt.xticks(fontsize = 28)
_ = plt.yticks(fontsize = 28)
Идея:
делить отрезок $[a,b]$ не на две равные насти,
а в пропорции "золотого сечения".
Оценим скорость сходимости аналогично методу дихотомии:
$$ |x_{K+1} - x^*| \leq b_{K+1} - a_{K+1} = \left( \frac{1}{\tau} \right)^{N-1} (b - a) \approx 0.618^K(b-a), $$где $\tau = \frac{\sqrt{5} + 1}{2}$.
In [9]:
def golden_search(f, a, b, tol=1e-5, callback=None):
tau = (np.sqrt(5) + 1) / 2.0
y = a + (b - a) / tau**2
z = a + (b - a) / tau
while b - a > tol:
if f(y) <= f(z):
b = z
z = y
y = a + (b - a) / tau**2
else:
a = y
y = z
z = a + (b - a) / tau
if callback is not None:
callback(a, b)
return (a + b) / 2.0
In [10]:
left_boud_gs = []
right_bound_gs = []
approximation_gs = []
cb_gs = lambda a, b: my_callback(a, b, left_boud_gs, right_bound_gs, approximation_gs)
x_gs = golden_search(f, a, b, epsilon, cb_gs)
print(f(x_opt))
print(f(x_gs))
print(np.abs(x_opt - x_true))
In [14]:
plt.figure(figsize=(10,6))
plt.semilogy(np.arange(1, len(approximation_bs) + 1), np.abs(x_true - np.array(approximation_bs, dtype=np.float64)), label="Binary search")
plt.semilogy(np.arange(1, len(approximation_gs) + 1), np.abs(x_true - np.array(approximation_gs, dtype=np.float64)), label="Golden search")
plt.xlabel(r"Number of iterations, $k$", fontsize=26)
plt.ylabel("Error rate upper bound", fontsize=26)
plt.legend(loc="best", fontsize=26)
plt.xticks(fontsize = 26)
_ = plt.yticks(fontsize = 26)
In [8]:
%timeit binary_search(f, a, b, epsilon)
%timeit golden_search(f, a, b, epsilon)
In [16]:
f = lambda x: np.sin(np.sin(np.sin(np.sqrt(x))))
x_true = (3 * np.pi / 2)**2
a = 2
b = 60
epsilon = 1e-8
plt.plot(np.linspace(a,b), f(np.linspace(a,b)))
plt.xticks(fontsize = 28)
_ = plt.yticks(fontsize = 28)
In [17]:
left_boud_bs = []
right_bound_bs = []
approximation_bs = []
callback_bs = lambda a, b: my_callback(a, b,
left_boud_bs, right_bound_bs, approximation_bs)
x_opt = binary_search(f, a, b, epsilon, callback_bs)
print(np.abs(x_opt - x_true))
In [23]:
left_boud_gs = []
right_bound_gs = []
approximation_gs = []
cb_gs = lambda a, b: my_callback(a, b, left_boud_gs, right_bound_gs, approximation_gs)
x_gs = golden_search(f, a, b, epsilon, cb_gs)
print(np.abs(x_opt - x_true))
In [24]:
plt.figure(figsize=(8,6))
plt.semilogy(np.abs(x_true - np.array(approximation_bs, dtype=np.float64)), label="Binary")
plt.semilogy(np.abs(x_true - np.array(approximation_gs, dtype=np.float64)), label="Golden")
plt.legend(fontsize=28)
plt.xticks(fontsize=28)
_ = plt.yticks(fontsize=28)
plt.xlabel(r"Number of iterations, $k$", fontsize=26)
plt.ylabel("Error rate upper bound", fontsize=26)
Out[24]:
In [13]:
%timeit binary_search(f, a, b, epsilon)
%timeit golden_search(f, a, b, epsilon)