Откуда берётся приближение: дискретизация и итерация

Урок разбирает два главных приёма, лежащих в основе почти любого численного алгоритма: дискретизацию и итерацию.

Дискретизация — замена непрерывного объекта (функции, области, времени) конечным набором значений в выбранных точках. Итерация — построение последовательности приближений, сходящейся к ответу.

Дискретизация: непрерывное → конечное

Компьютер не умеет работать с «настоящей» непрерывной функцией f(x) на отрезке — у отрезка бесконечно много точек, а памяти конечно. Поэтому первый шаг почти любого метода: выбрать конечную сетку точек и работать только со значениями функции в них. Хотим взять интеграл — разбиваем отрезок на полоски. Хотим найти производную — берём две близкие точки. Хотим решить дифференциальное уравнение — двигаемся по времени маленькими шагами.

Размер шага сетки h — главный параметр. Чем меньше h, тем точнее дискретное приближает непрерывное, но тем больше точек надо обработать. Это фундаментальный компромисс: точность против стоимости.

import math

# Дискретизируем sin(x) на [0, pi] и приблизим площадь суммой прямоугольников
def площадь(n):
    a, b = 0.0, math.pi
    h = (b - a) / n
    s = 0.0
    for i in range(n):
        x = a + (i + 0.5) * h   # середина i-го кусочка
        s += math.sin(x) * h
    return s

for n in [2, 4, 8, 16, 64, 256]:
    приближение = площадь(n)
    print(f"n={n:4}: интеграл = {приближение:.8f}  ошибка = {abs(приближение - 2.0):.2e}")

print("\nточное значение интеграла sin на [0,pi] = 2")

Вывод:

n=   2: интеграл = 2.22144147  ошибка = 2.21e-01
n=   4: интеграл = 2.05234431  ошибка = 5.23e-02
n=   8: интеграл = 2.01290909  ошибка = 1.29e-02
n=  16: интеграл = 2.00321638  ошибка = 3.22e-03
n=  64: интеграл = 2.00020081  ошибка = 2.01e-04
n= 256: интеграл = 2.00001255  ошибка = 1.25e-05

точное значение интеграла sin на [0,pi] = 2

Заметьте закономерность: при удвоении n ошибка падает примерно вчетверо. Это и есть «порядок точности» метода — позже мы научимся его измерять и предсказывать. А пока зафиксируем: дискретизация работает, и измельчая сетку, мы управляем точностью.

Итерация: приближение → лучшее приближение

Второй приём нужен, когда ответ нельзя получить прямой формулой даже после дискретизации. Тогда строят итерационный процесс: берут грубое начальное приближение x0 и по некоторому правилу получают x1, затем x2, и так далее. Если правило выбрано удачно, последовательность сходится к искомому значению.

Классический пример — вычисление квадратного корня методом, который знали ещё в Вавилоне. Чтобы найти √S, берём любое x > 0 и повторяем x ← (x + S/x) / 2 — среднее между x и S/x. Если x больше корня, то S/x меньше, и их среднее зажимается к истине.

def корень(S, x0=1.0):
    x = x0
    print(f"старт: x = {x}")
    for i in range(1, 7):
        x = (x + S / x) / 2        # вавилонское правило
        print(f"шаг {i}: x = {x:.15f}")
    return x

import math
S = 612.0
ответ = корень(S)
print(f"\nточное math.sqrt = {math.sqrt(S):.15f}")
print(f"совпало знаков: {ответ == math.sqrt(S)}")

Вывод:

старт: x = 1.0
шаг 1: x = 306.500000000000000
шаг 2: x = 154.248368678629703
шаг 3: x = 79.107997864354715
шаг 4: x = 43.422128682151481
шаг 5: x = 28.758162428779126
шаг 6: x = 25.019538536995718

точное math.sqrt = 24.738633753705962
совпало знаков: False

За шесть шагов из абсурдного старта x = 1 мы пришли к 25.02 при истинном 24.74. Ещё пара шагов — и совпадут все 15 знаков. Это поразительно быстро: число верных цифр здесь почти удваивается с каждым шагом (квадратичная сходимость — позже разберём, почему).

Как работает под капотом

Любой итерационный метод нуждается в двух вещах: правиле перехода x_{n+1} = g(x_n) и критерии остановки. Правило должно быть «сжимающим» — приближать к ответу. Критерий остановки чаще всего такой: остановиться, когда |x_{n+1} − x_n| стало меньше заданного порога eps, или когда выполнено максимальное число итераций (страховка от зацикливания, если метод вдруг расходится).

x = x0
for k in range(max_iter):
    x_new = g(x)
    if abs(x_new - x) < eps:   # достаточно близко — выходим
        return x_new
    x = x_new
raise RuntimeError("не сошлось за max_iter шагов")

Эта рамка повторяется в десятках методов: Ньютон, секущие, простая итерация, Якоби, Зейдель, степенной метод — все они отличаются только видом функции g.

Частые ошибки

• Бесконечно измельчать сетку. При очень малом h начинают доминировать ошибки округления (об этом — отдельный урок), и точность не растёт, а падает. Есть оптимальный h.
• Останавливать итерацию только по числу шагов. Лучше комбинировать: и порог по |x_{n+1}−x_n|, и лимит итераций как страховка.
• Считать, что итерация всегда сходится. Неудачное правило g или плохой старт могут давать расходимость или зацикливание — сходимость надо обосновывать.

Итоги

• Дискретизация заменяет непрерывный объект конечным набором значений; шаг h управляет точностью ценой объёма счёта.
• Итерация строит сходящуюся к ответу последовательность приближений по правилу x_{n+1} = g(x_n).
• Критерий остановки итерации: малое изменение |x_{n+1} − x_n| плюс лимит шагов как страховка.
• Эти два приёма — каркас, на который натянуты почти все методы курса.