task1
This commit is contained in:
255
task1/uniform_search.py
Normal file
255
task1/uniform_search.py
Normal file
@@ -0,0 +1,255 @@
|
||||
"""Лекция 4 (стр. 15)
|
||||
|
||||
Задача. Найти пару {x*, f(x*)} для целевой функции вида
|
||||
f(x) = sqrt(x^2 + 9) / 4 + (5 - x) / 5
|
||||
при условии, что X = [-3, 8].
|
||||
Взять N = 10, eps_x = 0.05, eps_f = 0.001.
|
||||
Применить метод равномерного поиска.
|
||||
"""
|
||||
|
||||
import math
|
||||
from pathlib import Path
|
||||
|
||||
import matplotlib.pyplot as plt
|
||||
import numpy as np
|
||||
|
||||
# Глобальный счётчик обращений к оракулу
|
||||
oracle_calls = 0
|
||||
|
||||
|
||||
def f(x: float) -> float:
|
||||
global oracle_calls
|
||||
oracle_calls += 1
|
||||
return math.sqrt(x**2 + 9) / 4 + (5 - x) / 5
|
||||
|
||||
|
||||
a = -3
|
||||
b = 8
|
||||
N = 10
|
||||
eps_x = 0.05
|
||||
eps_f = 0.001
|
||||
|
||||
"""
|
||||
Количество итераций может быть найдено из условия:
|
||||
(2/N)^k * (b - a) / 2 <= eps_x
|
||||
"""
|
||||
print("Метод равномерного поиска")
|
||||
print("=" * 80)
|
||||
print(
|
||||
f"Расчётное количество итераций: {math.ceil(math.log(2 * eps_x / (b - a), 2 / N))}"
|
||||
)
|
||||
|
||||
# Создаём папку для графиков
|
||||
Path("uniform_search_plots").mkdir(exist_ok=True)
|
||||
|
||||
# Сохраняем начальные значения для графиков
|
||||
a_init = a
|
||||
b_init = b
|
||||
|
||||
# Метод равномерного поиска
|
||||
k = 0
|
||||
a_k, b_k = a, b
|
||||
prev_best_f = None
|
||||
best_x = None
|
||||
best_fx = None
|
||||
history = [] # для отладки/анализа
|
||||
|
||||
print("\nНачальный интервал: [{:.6f}, {:.6f}]".format(a_k, b_k))
|
||||
print("-" * 80)
|
||||
|
||||
while True:
|
||||
k += 1
|
||||
# Разбиваем [a_k, b_k] на N равных частей
|
||||
Delta = (b_k - a_k) / N
|
||||
xs = [a_k + i * Delta for i in range(N + 1)] # включая концы
|
||||
# Считаем значения в узлах x_1..x_{N-1}
|
||||
vals = [f(x) for x in xs]
|
||||
# Ищем минимум среди внутренних узлов (по лекции)
|
||||
i_min = min(range(1, N), key=lambda i: vals[i])
|
||||
x_min = xs[i_min]
|
||||
fx_min = vals[i_min]
|
||||
|
||||
# Сжимаем интервал к соседям минимума (x_{i-1}, x_{i+1})
|
||||
a_next = xs[i_min - 1]
|
||||
b_next = xs[i_min + 1]
|
||||
|
||||
# Обновляем лучшее известное
|
||||
best_x = x_min
|
||||
best_fx = fx_min
|
||||
|
||||
# Критерии остановки (из лекции): по точности по x и стабилизации по f
|
||||
Delta_next = b_next - a_next # длина нового интервала
|
||||
center_next = (a_next + b_next) / 2 # удобная оценка x*
|
||||
fx_center_next = f(center_next)
|
||||
|
||||
# По лекции можно использовать центр нового интервала как текущую оценку x*
|
||||
# и/или минимальный узел — обе оценки допустимы. Здесь — центр интервала,
|
||||
# чтобы гарантировать |x - x*| <= Delta_next/2.
|
||||
approx_x = center_next
|
||||
approx_fx = fx_center_next
|
||||
|
||||
if prev_best_f is None:
|
||||
prev_best_f = approx_fx
|
||||
|
||||
# Условия остановки
|
||||
stop_by_x = (Delta_next / 2) <= eps_x
|
||||
stop_by_f = abs(approx_fx - prev_best_f) <= eps_f
|
||||
|
||||
history.append(
|
||||
{
|
||||
"iter": k,
|
||||
"a_k": a_k,
|
||||
"b_k": b_k,
|
||||
"i_min": i_min,
|
||||
"x_min": x_min,
|
||||
"f(x_min)": fx_min,
|
||||
"a_next": a_next,
|
||||
"b_next": b_next,
|
||||
"Delta_next/2": Delta_next / 2,
|
||||
"approx_x": approx_x,
|
||||
"approx_fx": approx_fx,
|
||||
"xs": xs.copy(),
|
||||
"vals": vals.copy(),
|
||||
# "stop_by_x": stop_by_x,
|
||||
# "stop_by_f": stop_by_f,
|
||||
}
|
||||
)
|
||||
|
||||
# Построение графика для текущей итерации
|
||||
x_plot = np.linspace(a_init, b_init, 500)
|
||||
y_plot = [f(x) for x in x_plot]
|
||||
|
||||
plt.figure(figsize=(12, 7))
|
||||
plt.plot(x_plot, y_plot, "b-", linewidth=2, label="f(x)")
|
||||
|
||||
# Текущий интервал
|
||||
plt.axvline(a_k, color="red", linestyle="--", alpha=0.7, linewidth=1.5)
|
||||
plt.axvline(b_k, color="red", linestyle="--", alpha=0.7, linewidth=1.5)
|
||||
plt.axvspan(
|
||||
a_k,
|
||||
b_k,
|
||||
alpha=0.1,
|
||||
color="red",
|
||||
label=f"Текущий интервал [{a_k:.3f}, {b_k:.3f}]",
|
||||
)
|
||||
|
||||
# Точки разбиения
|
||||
plt.plot(
|
||||
xs, vals, "o", color="orange", markersize=6, alpha=0.6, label="Узлы разбиения"
|
||||
)
|
||||
|
||||
# Минимальная точка
|
||||
plt.plot(
|
||||
x_min,
|
||||
fx_min,
|
||||
"go",
|
||||
markersize=12,
|
||||
label=f"Минимум: x={x_min:.4f}, f(x)={fx_min:.4f}",
|
||||
zorder=5,
|
||||
)
|
||||
|
||||
# Следующий интервал
|
||||
plt.axvspan(
|
||||
a_next,
|
||||
b_next,
|
||||
alpha=0.15,
|
||||
color="green",
|
||||
label=f"Новый интервал [{a_next:.3f}, {b_next:.3f}]",
|
||||
)
|
||||
|
||||
plt.xlabel("x", fontsize=12)
|
||||
plt.ylabel("f(x)", fontsize=12)
|
||||
plt.title(
|
||||
f"Равномерный поиск - Итерация {k}\nN={N}, D={Delta:.4f}, D_next/2={Delta_next / 2:.4f}",
|
||||
fontsize=13,
|
||||
fontweight="bold",
|
||||
)
|
||||
plt.legend(fontsize=9, loc="best")
|
||||
plt.grid(True, alpha=0.3)
|
||||
plt.tight_layout()
|
||||
plt.savefig(f"uniform_search_plots/iteration_{k:02d}.png", dpi=150)
|
||||
plt.close()
|
||||
|
||||
# Вывод информации по итерации
|
||||
print(f"\nИтерация {k}:")
|
||||
print(f" Интервал: [{a_k:.6f}, {b_k:.6f}], D = {Delta:.6f}")
|
||||
print(f" Минимум найден в узле i={i_min}: x = {x_min:.6f}, f(x) = {fx_min:.6f}")
|
||||
print(
|
||||
f" Новый интервал: [{a_next:.6f}, {b_next:.6f}], D_next/2 = {Delta_next / 2:.6f}"
|
||||
)
|
||||
print(f" Аппроксимация: x ~= {approx_x:.6f}, f(x) ~= {approx_fx:.6f}")
|
||||
print(f" Критерии: D/2 <= eps_x? {stop_by_x}, |Df| <= eps_f? {stop_by_f}")
|
||||
|
||||
# Переход на следующую итерацию
|
||||
a_k, b_k = a_next, b_next
|
||||
prev_best_f = approx_fx
|
||||
|
||||
if stop_by_x and stop_by_f:
|
||||
break
|
||||
|
||||
# Дополнительная защита: если интервал перестал уменьшаться (численная стабильность)
|
||||
if Delta_next <= 1e-12:
|
||||
break
|
||||
|
||||
# Финальный ответ
|
||||
x_star = approx_x
|
||||
f_star = approx_fx
|
||||
|
||||
# Финальный график
|
||||
x_plot = np.linspace(a_init, b_init, 500)
|
||||
y_plot = [f(x) for x in x_plot]
|
||||
|
||||
plt.figure(figsize=(12, 7))
|
||||
plt.plot(x_plot, y_plot, "b-", linewidth=2, label="f(x)")
|
||||
plt.axvline(a_k, color="red", linestyle="--", alpha=0.7, linewidth=1.5)
|
||||
plt.axvline(b_k, color="red", linestyle="--", alpha=0.7, linewidth=1.5)
|
||||
plt.axvspan(
|
||||
a_k,
|
||||
b_k,
|
||||
alpha=0.15,
|
||||
color="red",
|
||||
label=f"Финальный интервал [{a_k:.6f}, {b_k:.6f}]",
|
||||
)
|
||||
plt.plot(
|
||||
x_star,
|
||||
f_star,
|
||||
"r*",
|
||||
markersize=20,
|
||||
label=f"x* = {x_star:.6f}, f(x*) = {f_star:.6f}",
|
||||
zorder=5,
|
||||
)
|
||||
|
||||
plt.xlabel("x", fontsize=12)
|
||||
plt.ylabel("f(x)", fontsize=12)
|
||||
plt.title(
|
||||
f"Равномерный поиск - Финальный результат\nВыполнено итераций: {k}",
|
||||
fontsize=14,
|
||||
fontweight="bold",
|
||||
)
|
||||
plt.legend(fontsize=10, loc="best")
|
||||
plt.grid(True, alpha=0.3)
|
||||
plt.tight_layout()
|
||||
plt.savefig("uniform_search_plots/final_result.png", dpi=150)
|
||||
plt.close()
|
||||
|
||||
print("\n" + "=" * 80)
|
||||
print(f"\nИтераций выполнено: {k}")
|
||||
print(f"x* ~= {x_star:.6f}")
|
||||
print(f"f(x*) ~= {f_star:.6f}")
|
||||
print(f"Длина конечного интервала / 2 <= eps_x? {(b_k - a_k) / 2 <= eps_x}")
|
||||
print(f"|Df| <= eps_f? {abs(approx_fx - prev_best_f) <= eps_f}")
|
||||
# print(f"\nОбращений к оракулу (вычислений функции): {oracle_calls}")
|
||||
print("\nГрафики сохранены в папке 'uniform_search_plots/'")
|
||||
|
||||
# При желании выведем краткую таблицу нескольких последних итераций
|
||||
print("\n" + "=" * 80)
|
||||
print("Сводная таблица:")
|
||||
print("-" * 80)
|
||||
for row in history:
|
||||
print(
|
||||
f"Итер {row['iter']:2d}: "
|
||||
f"[{row['a_k']:7.4f}, {row['b_k']:7.4f}] -> "
|
||||
f"x_min={row['x_min']:7.4f} (f={row['f(x_min)']:.4f}) -> "
|
||||
f"[{row['a_next']:7.4f}, {row['b_next']:7.4f}], "
|
||||
f"D/2={row['Delta_next/2']:.4f}"
|
||||
)
|
||||
Reference in New Issue
Block a user