Запрещённая математика в твоём autograd: бесконечно малые, дуальные числа и нестандартный анализ

TL;DR

Когда вы пишете loss.backward(), ваш autograd делает то, что 200 лет считалось математической ересью: оперирует бесконечно малыми как настоящими числами. В 1960 году Абрахам Робинсон формализовал эту «ересь» в виде нестандартного анализа. Forward-mode автодифференцирование, на котором держатся JAX, PyTorch и пол-индустрии — это его обрезанная версия. В этой статье разберём гиперреалы и монады, реализуем дуальные числа в коде.

Проблема, о которой не говорят

Откройте любой учебник термодинамики. Найдите там первое начало:

dU=δQ−δA

Один значок прямой, другой — кривой. Спросите автора учебника, чем δQ отличается от dU — вам начнут объяснять про «полный» и «неполный» дифференциал. Спросите математика-аналитика, что вообще такое dx — он начнёт рассказывать про предельные переходы.

Проблема в том, что с точки зрения стандартного матанализа dU и δQ — это один и тот же объект. Оба определяются через предел. И тот, и другой стремится к нулю. Если у нас в руках только аппарат пределов, мы не можем сказать, чем интегрирование первого по замкнутому контуру отличается от интегрирования второго.

Физик это интуитивно понимает: δQ — маленький кусочек тепла, зависящий от пути. dU — приращение функции состояния, не зависящий. Но математически в стандартном анализе разница затаскивается в форму записи и определение интеграла, а не сидит в самих объектах. Дифференциалы как объекты в обычном анализе попросту отсутствуют — только пределы.

Краткая история запрещённой идеи

Лейбниц в 1684 году ввёл в анализ бесконечно малые — числа, которые больше нуля, но меньше любого положительного действительного. Производную он определял буквально: dy/dx, отношение двух бесконечно малых. Ньютон делал что-то похожее, кстати говоря.

Это работало. С помощью этого аппарата Эйлер посчитал такое, что до сих пор переоткрывают. Но был один нюанс: бесконечно малые непонятно как определить. Если ε > 0, но ε < r для любого положительного действительного r — то что это вообще такое за объект?

Беркли накинулся с философской критикой («призраки исчезнувших величин»). В XIX веке Коши, а потом более строго Вейерштрасс заменили бесконечно малые на пределы и (ε, δ)-определения. Получилось строго. Получилось тяжело. Получилось неудобно для физиков. Бесконечно малые формально из математики изгнали — точно так же, как когда-то изгоняли ноль и отрицательные числа.

Так оно и существовало до 1960 года.

В 1960-м Абрахам Робинсон опубликовал работу, в которой строго, в рамках теории моделей, построил расширение действительных чисел — гиперреалы *R — содержащее бесконечно большие и бесконечно малые числа как полноправные элементы поля. Бесконечно малые перестали быть «недочислами» и превратились в обычные числа, с которыми можно делать всё то же, что и с обычными.

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

Что такое гиперреалы

Достаточно знать следующее.

*R — это поле, которое содержит обычные действительные числа R и при этом содержит:

• бесконечно большие числа H, такие что H > n для любого натурального n;

• бесконечно малые числа ε = 1/H, такие что 0 < |ε| < 1/n для любого натурального n.

Внутри действуют все обычные алгебраические правила. Можно складывать, умножать, делить. 2ε — тоже бесконечно малое. ε² — ещё меньше, чем ε. H + ε — бесконечно большое. И так далее.

Главное свойство: для любого конечного гиперреала x (то есть не бесконечно большого) существует единственное действительное число st(x), бесконечно близкое к нему. Это число называется стандартной частью x. Например, st(2 + ε) = 2, а st(3 + 5ε - 7ε²) = 3.

Вокруг каждой действительной точки клубится облако бесконечно близких к ней гиперреалов — это и есть монада этой точки в смысле Лейбница. Монада — это окрестность, состоящая из самой точки и всех бесконечно малых отклонений от неё. Топологическая интуиция тут принципиально другая: вместо «точки сколь угодно близко» появляется «целое облако точек, каждая из которых уже неотличима от центра».

Аксиома Архимеда не работает — и это фича

Аксиома Архимеда: для любых двух положительных чисел a и b найдётся натуральное n, такое что n·a > b. В обычной математике это очевидно: складывай a достаточно много раз — рано или поздно перешагнёшь любое b.

В *R это не работает. Если ε бесконечно мало, а b = 1, то n·ε остаётся бесконечно малым для любого натурального n. Никогда не дотянет до единицы. Это не баг — это базовое свойство нестандартного континуума.

И, кстати, это не такая уж экзотика. Когда в физике одновременно требуют dt → 0 и t → ∞, и волнуются, как именно эти пределы связаны — фактически вопрос про не-архимедову структуру.

Производная без пределов

В стандартном анализе:

В нестандартном анализе:

где ε — конкретное бесконечно малое, а st(…) берёт стандартную часть.

Никаких пределов. Никаких дельта-эпсилон

Для f(x) = x²:

(f(x + ε) - f(x)) / ε
  = ((x + ε)² - x²) / ε
  = (x² + 2xε + ε² - x²) / ε
  = (2xε + ε²) / ε
  = 2x + ε

st(2x + ε) = 2x. Готово. Никакого предела не брали — только алгебра в гиперреалах плюс взятие стандартной части. Производная играет роль объекта.

Дуальные числа: бесконечно малые, в которые поверил инженер

Полный нестандартный анализ — штука с серьёзной теоретико-модельной начинкой, и реализовать его на компьютере целиком невозможно. Но есть его обрезанная версия, которая прекрасно реализуется и которой пользуется буквально вся индустрия.

Это дуальные числа.

Дуальное число — это пара (a, b), которую обычно записывают как a + bε, где ε² = 0, но ε ≠ 0.

Это как комплексные числа, только вместо i² = -1 мы постулируем ε² = 0. Получается коммутативное кольцо (не поле, но нам и не надо).

Сложение и умножение работают по правилам:

(a + bε) + (c + dε) = (a + c) + (b + d)ε
(a + bε) · (c + dε) = ac + (ad + bc)ε + bd·ε²
                    = ac + (ad + bc)ε        ← ε² = 0

И теперь фокус. Подставим x + ε в любую гладкую функцию и формально разложим в ряд Тейлора:

f(x + ε) = f(x) + f'(x)·ε + ½ f''(x)·ε² + ...
         = f(x) + f'(x)·ε                   ← всё, что с ε², обнуляется

Вычисляем функцию в дуальной точке — на выходе получаем пару (значение, производная) бесплатно. Никаких разностных схем, ни потери точности из-за деления маленького на маленькое. Точная производная, посчитанная одним прогоном алгоритма.

Реализация в NumPy

Делаем класс дуальных чисел:

import numpy as np

class Dual:
    """Дуальное число: real + dual·ε, где ε² = 0."""
    __slots__ = ("real", "dual")

    def __init__(self, real, dual=0.0):
        self.real = np.asarray(real, dtype=float)
        self.dual = np.asarray(dual, dtype=float)

    def __repr__(self):
        return f"Dual({self.real}, {self.dual}ε)"

    def __add__(self, other):
        other = other if isinstance(other, Dual) else Dual(other)
        return Dual(self.real + other.real, self.dual + other.dual)

    def __radd__(self, other):
        return self + other

    def __sub__(self, other):
        other = other if isinstance(other, Dual) else Dual(other)
        return Dual(self.real - other.real, self.dual - other.dual)

    def __rsub__(self, other):
        return Dual(other) - self

    def __mul__(self, other):
        other = other if isinstance(other, Dual) else Dual(other)
        return Dual(
            self.real * other.real,
            self.real * other.dual + self.dual * other.real,
        )

    def __rmul__(self, other):
        return self * other

    def __truediv__(self, other):
        other = other if isinstance(other, Dual) else Dual(other)
        return Dual(
            self.real / other.real,
            (self.dual * other.real - self.real * other.dual) / (other.real ** 2),
        )

    def __pow__(self, n):
        # n — обычное число
        return Dual(self.real ** n, n * self.real ** (n - 1) * self.dual)

    def __neg__(self):
        return Dual(-self.real, -self.dual)

Дальше — элементарные функции. Для каждой мы вручную записываем её разложение в дуальной арифметике, используя её производную:

def sin(x):
    if isinstance(x, Dual):
        return Dual(np.sin(x.real), np.cos(x.real) * x.dual)
    return np.sin(x)

def cos(x):
    if isinstance(x, Dual):
        return Dual(np.cos(x.real), -np.sin(x.real) * x.dual)
    return np.cos(x)

def exp(x):
    if isinstance(x, Dual):
        e = np.exp(x.real)
        return Dual(e, e * x.dual)
    return np.exp(x)

def log(x):
    if isinstance(x, Dual):
        return Dual(np.log(x.real), x.dual / x.real)
    return np.log(x)

И обёртка для вычисления производной:

def derivative(f, x):
    """Производная f в точке x — одним прогоном."""
    return f(Dual(x, 1.0)).dual

Тестируем на чём-нибудь нетривиальном:

def f(x):
    return sin(x ** 2) * exp(-x) + log(1 + x ** 2)

x0 = 1.3
print("f'(1.3) через дуальные числа:    ", derivative(f, x0))

Сверяем с центральной разностью:

h = 1e-6
fd = (f(x0 + h) - f(x0 - h)) / (2 * h)
print("f'(1.3) через центральную разность:", fd)

Числа совпадут до плавающей точности. Только дуальная версия посчитала точно (до ошибок округления самой арифметики), а центральная разность — приближённо, с ошибкой O(h²) плюс чувствительностью к выбору h.

Где это реально полезно за пределами «вау-факта»

Несколько мест, где гиперреальная (или дуальная) интуиция работает на практике.

Численные методы. Когда вы сравниваете центральную разность с дуальными числами, они дают одинаковый ответ при сопоставимой стоимости — но разностная схема ломается на функциях, чувствительных к ошибкам округления (вычитание близких чисел в числителе). Дуальные не ломаются. Это используют в CFD-симуляциях для построения якобиана и в оптимизации, где нужны очень точные градиенты.

Чувствительностный анализ. В моделях с большим числом параметров иногда хочется получить одновременно значение функции и её производную по одному выбранному направлению — без полного градиента. Forward-mode AD — это ровно оно, и он дешевле обратного режима, когда выходов много, а входов мало.

Гессианы и высшие производные. Если вложить дуальные числа сами в себя — Dual(Dual(...)) — получается машинерия для второй производной. Расширение на «многоуровневые ε» даёт высшие порядки. По сути это вторая аппроксимация нестандартного анализа.

Дискретные модели пространства-времени. В подходах с минимальной длиной (петлевая квантовая гравитация, причинные множества, разные дискретные модели) не-архимедова структура возникает естественно. Гиперреалы дают язык, на котором это удобно описывать без перехода к «всё-таки в пределе непрерывно».

Парадоксы стандартной теории множеств. Парадокс Банаха–Тарского — как разрезать шар на пять кусков и собрать из них два таких же шара — это сигнал, что с обычной аксиоматикой не всё благополучно на бесконечности. Альтернативные конструкции вещественной прямой — повод задуматься, какая из них вообще описывает физическую реальность, а не идеализированный континуум.

К чему это все?

Бесконечно малые два с половиной века считались либо метафорой, либо ошибкой. Физики ими тихо пользовались, математики делали вид, что их не существует, и переписывали всё через пределы. В 1960 году Робинсон показал, что они не метафора. Это нормальные числа в нормальном поле, просто в большем, чем R.

А в 2026-м ваш autograd считает на их обрезанной версии градиенты для трансформеров.

Иногда математика возвращается с того света неожиданным образом.