Плавающая точка на уровне битов

Целые числа устроены честно: тридцать два бита, два миллиарда значений подряд, шаг всегда единица. Числа с плавающей точкой устроены иначе. Те же тридцать два бита, но значения расставлены неравномерно: густо у нуля, редко на краю, и в самом конце шкалы вместо чисел сидят бесконечность и NaN. Из этой одной кодировки вырастает всё, что про float обычно рассказывают как список странностей: почему сумма двух десятых не даёт три десятых, почему у нуля два представления, почему деление на ноль не падает. В этом уроке мы соберём формат снизу, с мини-версии на шесть бит, где всё множество значений помещается в одну таблицу, и дойдём до настоящих f16, f32, f64, f80 и f128 в Zig. А потом будем писать код, который видит float только как биты.

Цели урока

• Собрать формат с плавающей точкой из трёх полей: знак, порядок со смещением, мантисса.
• Понимать, почему порядок хранится со смещением, а не в дополнительном коде.
• Читать любой битовый шаблон как одно из пяти значений: ноль, денормаль, нормальное число, бесконечность, NaN.
• Знать правило округления к ближайшему чётному и уметь показать его на конкретной паре соседей.
• Видеть в денормалях не исключение, а способ закрыть дыру у нуля равномерной сеткой.
• Гонять @bitCast между u32 и f32, между u64 и f64, и понимать, что это переименование, а не преобразование.
• Писать операции над float целочисленным кодом: удвоение, деление пополам, перевод в целое и обратно.
• Понимать, откуда берётся накопленная ошибка представления и чем она обошлась батарее Patriot в Дахране.
• Различать @floatCast, @intFromFloat, @floatFromInt и знать, что Zig делает на краю диапазона в разных режимах сборки.

Идея: научная запись, только двоичная

В прошлом уроке целые числа вели себя предсказуемо: разряды идут подряд, шаг всегда единица, переполнение сворачивается по кругу. С дробями так не выйдет.

Дробь можно хранить с фиксированной точкой: договориться, что младшие восемь бит это доли, старшие двадцать четыре целая часть. Просто, быстро, и совершенно негодно, когда в одной программе живут и заряд электрона, и масса Солнца: фиксированная точка либо потеряет мелочь, либо не дотянется до крупного.

Плавающая точка решает это ровно так же, как научная запись решает ту же задачу на бумаге. Число раскладывается на три части.

• Знак: плюс или минус, один бит.
• Порядок: на сколько разрядов сдвинута двоичная точка, целое число.
• Мантисса: значащие цифры, дробь от единицы до двух.

Значение собирается по формуле «мантисса умножить на два в степени порядка», и точка буквально плавает: одна и та же мантисса при разном порядке даёт и микроскопическое, и гигантское число. Расплата за это тоже понятна сразу: значащих цифр фиксированное количество, поэтому чем больше число, тем крупнее шаг сетки. Двадцати четырёх значащих двоичных разрядов хватает, чтобы точно записать любое целое до шестнадцати миллионов, а дальше идут только чётные, потом только кратные четырём, и так далее.

Дальше весь урок про то, как эти три части упаковать в фиксированное число бит так, чтобы кодировка сама разбиралась и с краевыми случаями.

Мини-формат: шесть бит, и всё множество на ладони

Разбирать f64 сразу тяжело: там 2⁶⁴ значений, никакой таблицей их не охватишь. Поэтому начнём с игрушки. Правила возьмём настоящие, а ширины уменьшим: один бит знака, три бита порядка, два бита мантиссы. Шесть бит, шестьдесят четыре шаблона, и половина из них помещается в один экран.

Правило кодирования такое же, как в IEEE-754.

• Поле порядка хранит не сам показатель, а показатель плюс смещение. Для трёх бит смещение равно 2² минус 1, то есть 3. Показатель 0 записывается как 3, показатель минус 2 как 1.
• Мантисса хранит только дробную часть. Целая часть подразумевается: перед двоичной точкой стоит единица, которую никто не записывает. Это скрытая единица, бесплатный лишний бит точности.
• Два значения поля порядка зарезервированы. Всё нули означает ноль и денормали, все единицы означает бесконечность и NaN.

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

Напишем декодер мини-формата и напечатаем всю неотрицательную половину.

const std = @import("std");

/// Мини-формат: 1 бит знака, 3 бита порядка, 2 бита мантиссы.
const exp_bits = 3;
const mant_bits = 2;
const bias = (1 << (exp_bits - 1)) - 1; // 3
const exp_max = (1 << exp_bits) - 1; // 7
const mant_count = 1 << mant_bits; // 4

const Kind = enum { zero, denormal, normal, infinity, nan };

const Decoded = struct { kind: Kind, value: f64 };

fn decode(pattern: u6) Decoded {
    const sign: f64 = if (pattern >> (exp_bits + mant_bits) == 1) -1 else 1;
    const e: u32 = (pattern >> mant_bits) & exp_max;
    const m: u32 = pattern & (mant_count - 1);
    const frac = @as(f64, @floatFromInt(m)) / @as(f64, mant_count);

    if (e == exp_max) return if (m == 0)
        .{ .kind = .infinity, .value = sign * std.math.inf(f64) }
    else
        .{ .kind = .nan, .value = std.math.nan(f64) };

    if (e == 0) return .{
        .kind = if (m == 0) .zero else .denormal,
        .value = sign * frac * std.math.pow(f64, 2, 1 - bias),
    };

    const power = @as(f64, @floatFromInt(e)) - bias;
    return .{ .kind = .normal, .value = sign * (1 + frac) * std.math.pow(f64, 2, power) };
}

pub fn main(init: std.process.Init) !void {
    var buf: [4096]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    var counts = [_]u32{0} ** 5;
    var pattern: u6 = 0;
    while (pattern < 32) : (pattern += 1) {
        const d = decode(pattern);
        counts[@intFromEnum(d.kind)] += 1;
        try out.print("{b:0>6}  E={d} M={d}  {s:<12} {d}\n", .{
            pattern,
            (pattern >> mant_bits) & exp_max,
            pattern & (mant_count - 1),
            @tagName(d.kind),
            d.value,
        });
    }
    try out.print("\nсчёт: ноль {d}, денормали {d}, нормальные {d}, inf {d}, NaN {d}\n", .{
        counts[0], counts[1], counts[2], counts[3], counts[4],
    });
    try out.flush();
}

$ zig run mini.zig
000000  E=0 M=0  zero         0
000001  E=0 M=1  denormal     0.0625
000010  E=0 M=2  denormal     0.125
000011  E=0 M=3  denormal     0.1875
000100  E=1 M=0  normal       0.25
000101  E=1 M=1  normal       0.3125
000110  E=1 M=2  normal       0.375
000111  E=1 M=3  normal       0.4375
001000  E=2 M=0  normal       0.5
001001  E=2 M=1  normal       0.625
001010  E=2 M=2  normal       0.75
001011  E=2 M=3  normal       0.875
001100  E=3 M=0  normal       1
001101  E=3 M=1  normal       1.25
001110  E=3 M=2  normal       1.5
001111  E=3 M=3  normal       1.75
010000  E=4 M=0  normal       2
010001  E=4 M=1  normal       2.5
010010  E=4 M=2  normal       3
010011  E=4 M=3  normal       3.5
010100  E=5 M=0  normal       4
010101  E=5 M=1  normal       5
010110  E=5 M=2  normal       6
010111  E=5 M=3  normal       7
011000  E=6 M=0  normal       8
011001  E=6 M=1  normal       10
011010  E=6 M=2  normal       12
011011  E=6 M=3  normal       14
011100  E=7 M=0  infinity     inf
011101  E=7 M=1  nan          nan
011110  E=7 M=2  nan          nan
011111  E=7 M=3  nan          nan

счёт: ноль 1, денормали 3, нормальные 24, inf 1, NaN 3

Посмотри на эту таблицу минуту, в ней есть весь урок.

1. Значения идут строго по возрастанию, ровно в том же порядке, что и битовые шаблоны. Это работа смещения.
2. Шаг сетки не постоянный. От 1/4 до 1/2 шаг равен 1/16, от 1 до 2 шаг равен 1/4, от 8 до 14 шаг равен 2. На каждой границе порядка шаг удваивается, потому что мантисса делит очередной отрезок вида 2ᴱ до 2ᴱ⁺¹ на четыре равные части. Разрешение падает ровно там, где числа растут.
3. У нуля свой отдельный шаблон 000000, и у него есть близнец 100000 со знаковым битом. Минус ноль это законное значение формата.
4. Три шаблона внизу таблицы это денормали, и их шаг такой же, как у первых нормальных: 1/16. Дыры между нулём и наименьшим нормальным числом нет.
5. На верхнем краю четыре шаблона отданы не числам: один бесконечности и три под NaN. Это не потери, это плата за то, что все краевые случаи умещаются в ту же кодировку.

Виджет ниже это та же таблица, только живая. В режиме «биты и значение» щёлкай по любому биту и смотри, как меняется класс значения и куда переезжает точка на шкале. Ползунки меняют ширины полей, от двух до пяти бит порядка и от двух до четырёх бит мантиссы: добавь бит порядка, и края шкалы разъедутся, добавь бит мантиссы, и сетка станет вдвое гуще. Нижняя шкала это увеличенная окрестность нуля: там видно, что денормали продолжают ту же равномерную сетку, что и первые нормальные числа. Режим «округление до чётного» пригодится в следующем разделе.

Пять классов из одного правила

Обрати внимание, что в декодере нет ни одного «особого случая», написанного отдельно от формата. Все пять классов получаются из двух проверок поля порядка.

Поле порядка	Мантисса	Что это	Значение
все единицы	ноль	бесконечность	плюс или минус inf по знаковому биту
все единицы	не ноль	NaN	ответа нет
все нули	ноль	ноль	плюс или минус ноль
все нули	не ноль	денормаль	0.mmm умножить на 2 в степени 1 минус смещение
остальное	любая	нормальное	1.mmm умножить на 2 в степени поле минус смещение

Две тонкости, которые стоит проговорить вслух.

Первая: у денормалей показатель равен единице минус смещение, а не нулю минус смещение, хотя поле порядка нулевое. Это сделано специально, чтобы сетка денормалей продолжила сетку первых нормальных чисел без разрыва. Проверь по таблице: наибольшая денормаль 3/16, наименьшее нормальное 1/4 = 4/16, разница 1/16, ровно шаг денормалей.

Вторая: у денормалей нет скрытой единицы, перед точкой подразумевается ноль. Именно поэтому шаблон из одних нулей это ноль, а не единица умножить на два в какой-то степени.

Округление до ближайшего чётного

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

Почему при ничьей именно к чётному, а не всегда вверх, как учили в школе? Потому что «всегда вверх» это систематический сдвиг. Округли так миллион случайных значений, и сумма поедет вверх на половину шага, умноженную на число ничьих. Округление к чётному распределяет ничьи поровну между «вверх» и «вниз», и накопленного сдвига не возникает.

Возьми виджет выше, переключись в режим «округление до чётного» и набери там 9. В нашем мини-формате соседи это 8 и 10, число ровно посередине. Шаблон восьмёрки 011000, младший бит мантиссы ноль, шаблон десятки 011001, младший бит единица. Побеждает восьмёрка. Теперь набери 11: соседи 10 и 12, снова ничья, но теперь чётная мантисса у двенадцати (011010), и результат едет вверх. Одна и та же ничья решается в разные стороны, и это ровно то, что нужно.

Набери ещё 15. Наибольшее конечное значение формата это 14, следующий узел сетки был бы 16, середина между ними как раз 15, и правило чётности отправляет результат в бесконечность. Так переполнение оказывается не отдельной проверкой, а следствием того же правила округления.

В настоящем f32 всё то же самое, только сетка мельче. Мантисса там двадцать три бита плюс скрытая единица, то есть двадцать четыре значащих разряда, и все целые до 2²⁴ = 16777216 записываются точно. Сразу за этой границей шаг сетки становится равен двум, и ничьи начинаются на каждом нечётном числе.

const std = @import("std");

pub fn main(init: std.process.Init) !void {
    var buf: [4096]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    // Ничья: ровно посередине между двумя соседями по сетке f32.
    inline for (.{ 16777216, 16777217, 16777218, 16777219, 16777220 }) |n| {
        const x: f32 = @floatFromInt(@as(i32, n));
        const bits: u32 = @bitCast(x);
        try out.print("{d:>9} -> {d:>10}  0x{X:0>8}  младший бит мантиссы {d}\n", .{
            n, x, bits, bits & 1,
        });
    }

    // Знаменитая пара: сумма в f64 не попадает в узел сетки.
    const a: f64 = 0.1;
    const b: f64 = 0.2;
    try out.print("\n0.1 + 0.2 = {d}   {x}\n", .{ a + b, a + b });
    try out.print("0.3       = {d}          {x}\n", .{ @as(f64, 0.3), @as(f64, 0.3) });
    try out.print("биты суммы  0x{X:0>16}\n", .{@as(u64, @bitCast(a + b))});
    try out.print("биты 0.3    0x{X:0>16}\n", .{@as(u64, @bitCast(@as(f64, 0.3)))});
    try out.print("равны ли: {}\n", .{a + b == 0.3});
    try out.flush();
}

$ zig run rounding.zig
 16777216 ->   16777216  0x4B800000  младший бит мантиссы 0
 16777217 ->   16777216  0x4B800000  младший бит мантиссы 0
 16777218 ->   16777218  0x4B800001  младший бит мантиссы 1
 16777219 ->   16777220  0x4B800002  младший бит мантиссы 0
 16777220 ->   16777220  0x4B800002  младший бит мантиссы 0

0.1 + 0.2 = 0.30000000000000004   0x1.3333333333334p-2
0.3       = 0.3          0x1.3333333333333p-2
биты суммы  0x3FD3333333333334
биты 0.3    0x3FD3333333333333
равны ли: false

Верхняя половина это правило чётности на настоящих числах: 16777217 едет вниз к 16777216, а 16777219 вверх к 16777220, и оба раза у победителя младший бит мантиссы нулевой.

Нижняя половина это самый цитируемый пример про float, и теперь он перестаёт быть загадкой. Одна десятая в двоичной записи бесконечная периодическая дробь, как одна треть в десятичной. Ни 0.1, ни 0.2 не лежат в узлах сетки f64, оба литерала уже при компиляции округлились к ближайшим соседям. Сумма двух округлённых значений округлилась ещё раз и попала на один шаг выше, чем ближайший к 0.3 узел. Шестнадцатеричная запись показывает это буквально: мантиссы отличаются последней цифрой, …334 против …333. Один младший бит. Никакой ошибки в процессоре нет, есть сетка с конечным шагом и три честных округления подряд.

Настоящие форматы: от f16 до f128

Теперь то же самое, но в полный рост. Zig даёт пять типов с плавающей точкой, и все они устроены по разобранным правилам, отличаются только ширины полей. Обзорный взгляд на те же форматы со стороны прикладного кода у тебя уже был в уроке про float в Rust, а здесь мы смотрим внутрь. Соберём их параметры через comptime, чтобы числа брались из компилятора, а не с потолка.

const std = @import("std");
const math = std.math;

fn row(out: *std.Io.Writer, comptime T: type) !void {
    const exp = math.floatExponentBits(T);
    try out.print("{s:<5}{d:>5}{d:>9}{d:>10}{d:>9}{d:>8}{d:>9}   {d}\n", .{
        @typeName(T),
        @bitSizeOf(T),
        exp,
        math.floatMantissaBits(T),
        (@as(u32, 1) << (exp - 1)) - 1,
        math.floatExponentMin(T),
        math.floatExponentMax(T),
        @sizeOf(T),
    });
}

pub fn main(init: std.process.Init) !void {
    var buf: [4096]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;
    try out.print("{s}\n", .{"тип    бит  порядок  мантисса смещение    Emin     Emax   байт"});
    inline for (.{ f16, f32, f64, f80, f128 }) |T| try row(out, T);
    try out.flush();
}

$ zig run formats.zig
тип    бит  порядок  мантисса смещение    Emin     Emax   байт
f16     16        5        10       15     -14       15   2
f32     32        8        23      127    -126      127   4
f64     64       11        52     1023   -1022     1023   8
f80     80       15        64    16383  -16382    16383   16
f128   128       15       112    16383  -16382    16383   16

Читай строку f32 как описание нашей игрушки, только с другими числами: один бит знака, восемь бит порядка, двадцать три бита мантиссы, смещение 127. Показатель нормальных чисел от -126 до 127, потому что поля 0 и 255 зарезервированы под краевые классы.

Три наблюдения из таблицы.

f16 это те же правила на шестнадцати битах. Формат настолько маленький, что все его значения можно пересчитать перебором, и это отличный способ проверить, что понимание кодировки правильное.

const std = @import("std");

pub fn main(init: std.process.Init) !void {
    var buf: [1024]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    var zero: u32 = 0;
    var denormal: u32 = 0;
    var normal: u32 = 0;
    var infinity: u32 = 0;
    var nan: u32 = 0;

    // Шестнадцать бит это всего 65536 шаблонов: перебираем их все.
    var bits: u17 = 0;
    while (bits < 65536) : (bits += 1) {
        const pattern: u16 = @truncate(bits);
        const exp = (pattern >> 10) & 0x1F;
        const mant = pattern & 0x03FF;
        if (exp == 0x1F) {
            if (mant == 0) infinity += 1 else nan += 1;
        } else if (exp == 0) {
            if (mant == 0) zero += 1 else denormal += 1;
        } else normal += 1;
    }

    try out.print("все 65536 шаблонов f16:\n", .{});
    try out.print("  нули          {d:>6}\n", .{zero});
    try out.print("  денормали     {d:>6}\n", .{denormal});
    try out.print("  нормальные    {d:>6}\n", .{normal});
    try out.print("  бесконечности {d:>6}\n", .{infinity});
    try out.print("  NaN           {d:>6}\n", .{nan});
    try out.print("наименьшая денормаль {e}, наименьшее нормальное {e}, наибольшее {e}\n", .{
        std.math.floatTrueMin(f16), std.math.floatMin(f16), std.math.floatMax(f16),
    });
    try out.flush();
}

$ zig run f16census.zig
все 65536 шаблонов f16:
  нули               2
  денормали       2046
  нормальные     61440
  бесконечности      2
  NaN             2046
наименьшая денормаль 6e-8, наименьшее нормальное 6.104e-5, наибольшее 6.55e4

Больше двух тысяч шаблонов из шестидесяти пяти тысяч отданы под NaN. В f32 их 16777214, в f64 больше девяти квадриллионов. Мы ещё вернёмся к вопросу, зачем формату столько способов сказать «ответа нет».

f80 это расширенная точность x87, и в нём есть уникальная особенность: скрытой единицы нет, старший бит мантиссы записан явно. Поэтому мантисса шестьдесят четыре бита, а не шестьдесят три. Ещё занятно, что @sizeOf(f80) равен шестнадцати, а не десяти: компилятор выравнивает значение, и шесть байт остаются мусором.

f128 даёт сто двенадцать бит мантиссы, но на большинстве машин у него нет аппаратной поддержки, и операции идут через программные библиотеки. Это медленно, поэтому в системном коде он появляется редко.

@bitCast: float это просто биты

Всё, что мы разбирали, лежит в памяти как обычное машинное слово, и в Zig есть способ посмотреть на него так. @bitCast берёт значение одного типа и читает те же биты как другой тип той же ширины. Никакого преобразования не происходит, меняется только то, как их интерпретируют.

const std = @import("std");

/// Разобрать f32 на три поля и напечатать вместе с классом значения.
fn show(out: *std.Io.Writer, x: f32) !void {
    const bits: u32 = @bitCast(x);
    const sign = bits >> 31;
    const exp: u8 = @truncate(bits >> 23);
    const mant = bits & 0x007F_FFFF;
    const class = switch (exp) {
        0 => if (mant == 0) "ноль" else "денормаль",
        255 => if (mant == 0) "бесконечность" else "NaN",
        else => "нормальное",
    };
    try out.print("{e:>10}  0x{X:0>8}  {b} {b:0>8} {b:0>23}  {s}\n", .{
        x, bits, sign, exp, mant, class,
    });
}

pub fn main(init: std.process.Init) !void {
    var buf: [4096]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    const values = [_]f32{ 1.0, 0.5, 2.0, -2.0, 3.0, 0.1, 0.0, -0.0 };
    for (values) |v| try show(out, v);
    try show(out, std.math.floatTrueMin(f32));
    try show(out, std.math.inf(f32));
    try show(out, std.math.nan(f32));

    // Обратный ход: биты собраны руками, а прочитаны как число.
    const hand: u32 = 0b0_10000000_10000000000000000000000;
    try out.print("\nсобрано руками: 0x{X:0>8} даёт {d}\n", .{ hand, @as(f32, @bitCast(hand)) });
    try out.flush();
}

$ zig run bits.zig
       1e0  0x3F800000  0 01111111 00000000000000000000000  нормальное
      5e-1  0x3F000000  0 01111110 00000000000000000000000  нормальное
       2e0  0x40000000  0 10000000 00000000000000000000000  нормальное
      -2e0  0xC0000000  1 10000000 00000000000000000000000  нормальное
       3e0  0x40400000  0 10000000 10000000000000000000000  нормальное
      1e-1  0x3DCCCCCD  0 01111011 10011001100110011001101  нормальное
       0e0  0x00000000  0 00000000 00000000000000000000000  ноль
      -0e0  0x80000000  1 00000000 00000000000000000000000  ноль
     1e-45  0x00000001  0 00000000 00000000000000000000001  денормаль
       inf  0x7F800000  0 11111111 00000000000000000000000  бесконечность
       nan  0x7FC00000  0 11111111 10000000000000000000000  NaN

собрано руками: 0x40400000 даёт 3

Разбери несколько строк вручную, это лучшая проверка понимания.

Единица: поле порядка 01111111 это 127, минус смещение 127 даёт показатель 0. Мантисса нулевая, значит значащая часть равна 1.000. Итого 1.0 умножить на 2⁰.

Двойка: поле порядка на единицу больше, 10000000 это 128, показатель 1. Мантисса та же. Итого 1.0 умножить на 2¹. Отсюда сразу видно правило, которым мы скоро воспользуемся: умножение на два это прибавление единицы к полю порядка.

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

Тройка: показатель 1, мантисса 100…0, значащая часть 1.5. Полторы умножить на два это три. И обратный ход в конце листинга: те же биты, набранные руками через двоичный литерал с подчёркиваниями по полям, прочитались как 3.

Одна десятая: показатель -4, мантисса 10011001100110011001101. Видишь период 1001, оборванный на двадцать третьем бите, и последний разряд 1 вместо 0? Это округление к ближайшему в действии, литерал приехал на один шаг вверх от истинного значения.

Для f64 всё то же самое, парный тип это u64.

const std = @import("std");

pub fn main(init: std.process.Init) !void {
    var buf: [1024]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    // f64 раскладывается на u64 ровно так же, как f32 на u32.
    const x: f64 = 1.0;
    try out.print("f64 1.0 -> 0x{X:0>16}\n", .{@as(u64, @bitCast(x))});
    try out.print("f64 0.1 -> 0x{X:0>16}\n", .{@as(u64, @bitCast(@as(f64, 0.1)))});

    // У f80 подразумеваемой единицы нет: старший бит мантиссы записан явно.
    const y: f80 = 1.0;
    const u: u80 = @bitCast(y);
    try out.print("f80 1.0 -> 0x{X:0>20}, старший бит мантиссы {d}\n", .{ u, (u >> 63) & 1 });
    try out.print("@sizeOf(f80) = {d}, @bitSizeOf(f80) = {d}\n", .{ @sizeOf(f80), @bitSizeOf(f80) });
    try out.flush();
}

$ zig run f80.zig
f64 1.0 -> 0x3FF0000000000000
f64 0.1 -> 0x3FB999999999999A
f80 1.0 -> 0x3FFF8000000000000000, старший бит мантиссы 1
@sizeOf(f80) = 16, @bitSizeOf(f80) = 80

Сравни последнюю пару строк с предыдущим листингом: у f64 единица кодируется как 0x3FF0…, у f80 как 0x3FFF8000…, и восьмёрка в мантиссе это та самая явно записанная старшая единица, которой в остальных форматах нет.

@bitCast требует совпадения ширины: f32 и u32, f64 и u64, f80 и u80. Ошибиться нельзя, компилятор не пропустит.

Денормали закрывают дыру у нуля

Представь, что зарезервированного нулевого поля порядка нет и наименьшее нормальное число это край. Тогда между нулём и этим краем лежала бы пропасть шириной с само число, а сетка внутри отрезка от него до вдвое большего была бы в миллионы раз гуще. Это ломало бы код на ровном месте: два разных числа могли бы дать в разности ноль.

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

const std = @import("std");

fn line(out: *std.Io.Writer, label: []const u8, x: f32) !void {
    const bits: u32 = @bitCast(x);
    try out.print("{s} {e:>14}  0x{X:0>8}\n", .{ label, x, bits });
}

pub fn main(init: std.process.Init) !void {
    var buf: [4096]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    try line(out, "наименьшее нормальное ", std.math.floatMin(f32));
    try line(out, "наибольшая денормаль  ", std.math.floatMin(f32) - std.math.floatTrueMin(f32));
    try line(out, "наименьшая денормаль  ", std.math.floatTrueMin(f32));

    // Постепенная потеря значимости: число не падает в ноль сразу.
    var x: f32 = std.math.floatMin(f32);
    var steps: u32 = 0;
    while (x != 0) : (steps += 1) x /= 2;
    try out.print("\nот наименьшего нормального до нуля: {d} шагов деления пополам\n", .{steps});

    // Два соседних числа у самой границы нормальных.
    const a: f32 = std.math.floatMin(f32);
    const b: f32 = @bitCast(@as(u32, @bitCast(a)) + 1);
    try out.print("\na и b различны: {}\n", .{a != b});
    try line(out, "b - a                 ", b - a);
    try out.flush();
}

$ zig run denormals.zig
наименьшее нормальное   1.1754944e-38  0x00800000
наибольшая денормаль    1.1754942e-38  0x007FFFFF
наименьшая денормаль            1e-45  0x00000001

от наименьшего нормального до нуля: 24 шагов деления пополам

a и b различны: true
b - a                           1e-45  0x00000001

Смотри на битовые шаблоны в первых трёх строках. Наименьшее нормальное это 0x00800000, наибольшая денормаль ровно на единицу меньше, 0x007FFFFF. Никакого разрыва: шаблоны идут подряд, значения тоже. Граница между классами существует только в правиле декодирования, на шкале её не видно.

Второй блок это постепенная потеря значимости в действии: от наименьшего нормального до нуля двадцать четыре деления пополам, а без денормалей хватило бы одного.

Третий блок это то самое свойство, ради которого всё затевалось. Числа a и b соседи по сетке, они различны, и их разность тоже отлична от нуля, причём она равна наименьшей денормали. Без денормалей ближайшее к нулю доступное значение было бы 2⁻¹²⁶, разность туда не влезла бы и обнулилась, и получилось бы, что a != b, но a - b == 0. Такое поведение ломает любой численный алгоритм, который проверяет сходимость через разность.

Платить за денормали приходится скоростью: на многих процессорах операции над ними идут по медленному пути, в десятки раз дольше обычных. Поэтому в графике и в звуке иногда включают режим, где денормали принудительно обнуляются. Это осознанный размен точности у самого нуля на предсказуемое время.

Бесконечности и NaN

Верхний край поля порядка отдан не числам. Мантисса ноль это бесконечность, мантисса не ноль это NaN.

const std = @import("std");

test "бесконечность это край шкалы, а не ошибка" {
    const inf = std.math.inf(f32);
    try std.testing.expectEqual(@as(u32, 0x7F80_0000), @as(u32, @bitCast(inf)));
    try std.testing.expectEqual(inf, std.math.floatMax(f32) * 2);
    try std.testing.expectEqual(inf, @as(f32, 1.0) / 0.0);
    try std.testing.expectEqual(@as(f32, 0.0), @as(f32, 1.0) / inf);
    try std.testing.expect(inf > std.math.floatMax(f32));
}

test "NaN появляется там, где ответа нет" {
    try std.testing.expect(std.math.isNan(@as(f32, 0.0) / 0.0));
    try std.testing.expect(std.math.isNan(std.math.inf(f32) - std.math.inf(f32)));
    try std.testing.expect(std.math.isNan(std.math.inf(f32) * 0.0));
    try std.testing.expect(std.math.isNan(std.math.sqrt(@as(f32, -1.0))));
}

test "у NaN много шаблонов, и старший бит мантиссы отличает тихий от сигнального" {
    const quiet: u32 = @bitCast(std.math.nan(f32));
    try std.testing.expectEqual(@as(u32, 0x7FC0_0000), quiet);
    // Тихий NaN: старший бит мантиссы единица.
    try std.testing.expect(quiet & 0x0040_0000 != 0);
    // Сигнальный: тот же порядок, но старший бит мантиссы ноль.
    const signaling: f32 = @bitCast(@as(u32, 0x7F80_0001));
    try std.testing.expect(std.math.isNan(signaling));
}

test "NaN не равен ничему, и это единственное значение с таким свойством" {
    const nan = std.math.nan(f32);
    try std.testing.expect(nan != nan);
    try std.testing.expect(!(nan < 1.0));
    try std.testing.expect(!(nan > 1.0));
    try std.testing.expect(!(nan == nan));
}

$ zig test specials.zig
1/4 specials.test.бесконечность это край шкалы, а не ошибка...OK
2/4 specials.test.NaN появляется там, где ответа нет...OK
3/4 specials.test.у NaN много шаблонов, и старший бит мантиссы отличает тихий от сигнального...OK
4/4 specials.test.NaN не равен ничему, и это единственное значение с таким свойством...OK
All 4 tests passed.

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

NaN появляется там, где ответа нет даже бесконечного: ноль делить на ноль, бесконечность минус бесконечность, корень из отрицательного. Ключевое свойство одно: NaN не равен ничему, включая себя, поэтому x != x истинно ровно для NaN и ни для чего другого.

Зачем формату столько шаблонов NaN? Разряды мантиссы у NaN не используются кодировкой, и их отдали под полезную нагрузку: туда можно положить код ошибки или свою метку, а процессор в большинстве случаев перенесёт её из аргумента в результат. Отдельно стандарт различает тихий и сигнальный NaN по старшему биту мантиссы: тихий тихо распространяется по вычислению, сигнальный задуман как повод поднять исключение. Практика чаще всего использует тихий, и std.math.nan(f32) даёт именно его.

Ещё одна нагрузка на эту идею живёт в динамических языках: в мантиссу тихого NaN укладывают указатель или тег типа, и одно 64-битное слово хранит либо число, либо ссылку. Приём называют упаковкой в NaN, и он появится в разделе, где мы будем писать свою виртуальную машину.

Порядок битов совпадает с порядком чисел

Вернёмся к свойству, ради которого порядок хранится со смещением. Если взять неотрицательное число и прочитать его биты как беззнаковое целое, то сравнение целых даст тот же ответ, что сравнение чисел. Проверим не на словах.

const std = @import("std");

test "порядок битов совпадает с порядком значений для неотрицательных" {
    var previous: f32 = 0;
    var bits: u32 = 0;
    // Идём по битовым шаблонам подряд: 0, потом денормали, потом нормальные.
    while (bits < 0x7F80_0000) : (bits += 0x1_0000) {
        const x: f32 = @bitCast(bits);
        try std.testing.expect(x >= previous);
        previous = x;
    }
}

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

С отрицательными числами свойство переворачивается: там больший шаблон означает меньшее значение, потому что величина растёт, а знак делает её отрицательной. Плюс отдельно стоят два нуля, которые по значению равны, но по битам различны. Из этих трёх наблюдений собирается сравнение float чисто целочисленным кодом, и это одна из задач в конце урока.

Операции над битами: удвоение и половина

Если умножение на два это плюс единица к полю порядка, то целую арифметику над float можно писать, не имея типа float вообще. В задачах главы 2 CS:APP это отдельный жанр: функция принимает u32, возвращает u32, и внутри разрешены только целочисленные операции.

Начнём с ядра идеи.

const std = @import("std");

const SIGN: u32 = 0x8000_0000;
const EXP: u32 = 0x7F80_0000;
const FRAC: u32 = 0x007F_FFFF;

/// Удвоение для нормального числа, которое не переполнится:
/// прибавить единицу к полю порядка, мантиссу не трогать.
fn twiceNormal(uf: u32) u32 {
    return uf + (1 << 23);
}

test "для нормальных чисел удвоение это плюс один к порядку" {
    const values = [_]f32{ 1.0, 3.5, -2.25, 0.1, 1e10, -7.0 };
    for (values) |x| {
        const bits: u32 = @bitCast(x);
        const doubled: f32 = @bitCast(twiceNormal(bits));
        try std.testing.expectEqual(x * 2, doubled);
    }
}

test "поле порядка сдвигает всю запись, а мантисса остаётся на месте" {
    const one: u32 = @bitCast(@as(f32, 1.0));
    const two: u32 = @bitCast(@as(f32, 2.0));
    try std.testing.expectEqual(two, one + (1 << 23));
    try std.testing.expectEqual(@as(u32, 0x3F80_0000), one);
    try std.testing.expectEqual(@as(u32, 0x4000_0000), two);
}

test "у денормали сдвиг влево сам переносит бит в поле порядка" {
    // Наибольшая денормаль: порядок ноль, мантисса из одних единиц.
    const big_denormal: u32 = FRAC;
    const doubled = (big_denormal & ~SIGN) << 1;
    const x: f32 = @bitCast(big_denormal);
    try std.testing.expectEqual(x * 2, @as(f32, @bitCast(doubled)));
    // Результат уже нормальный: поле порядка стало единицей.
    try std.testing.expect(doubled & EXP == (1 << 23));
}

$ zig test twice.zig
1/3 twice.test.для нормальных чисел удвоение это плюс один к порядку...OK
2/3 twice.test.поле порядка сдвигает всю запись, а мантисса остаётся на месте...OK
3/3 twice.test.у денормали сдвиг влево сам переносит бит в поле порядка...OK
All 3 tests passed.

Третий тест стоит того, чтобы задержаться. Наибольшая денормаль это шаблон 0x007FFFFF: поле порядка нулевое, мантисса из одних единиц. Сдвинем весь шаблон без знака влево на один бит, и старший бит мантиссы переедет в младший бит поля порядка. Формально мы сделали сдвиг целого числа, а по смыслу получили ровно удвоение: денормаль превратилась в наименьшее нормальное число, и подразумеваемая единица встала на место. Кодировку так и проектировали, чтобы граница между классами не требовала отдельного кода.

Что осталось за кадром и достанется тебе в задаче: переполнение (когда прибавление единицы к полю порядка упирается в зарезервированное значение и результат обязан стать бесконечностью), NaN и бесконечность на входе, и деление пополам, где всё сложнее. При делении числа с полем порядка, равным единице, результат уезжает в денормали, один бит мантиссы уходит за край сетки, и здесь нужно то самое округление к ближайшему чётному, которое мы разбирали на девятке и одиннадцати.

Преобразования в целое и обратно

Вторая пара задач жанра это перевод между u32-шаблоном float и обычным целым. Идея в обе стороны одна: значащая часть числа это двадцать четыре бита, где двоичная точка стоит сразу после старшего разряда, а поле порядка говорит, на сколько её сдвинуть.

const std = @import("std");

const FRAC: u32 = 0x007F_FFFF;
const HIDDEN: u32 = 0x0080_0000; // подразумеваемая единица перед точкой

/// Отсечение дробной части для положительных чисел из безопасной середины
/// диапазона: порядок от 0 до 22, поэтому сдвиг всегда вправо и всегда влезает.
fn f2iCore(uf: u32) i32 {
    const e: i32 = @as(i32, @intCast((uf >> 23) & 0xFF)) - 127;
    if (e < 0) return 0; // модуль меньше единицы, целая часть ноль
    const significand = (uf & FRAC) | HIDDEN; // 24 бита: 1.mmm как целое
    return @intCast(significand >> @intCast(23 - e));
}

/// Перевод целого в биты числа одинарной точности для тех значений,
/// что помещаются в 24 бита мантиссы: округлять нечего.
fn i2fCore(x: u32) u32 {
    if (x == 0) return 0;
    const top: u32 = 31 - @clz(x); // позиция старшей единицы, это и есть порядок
    const shift: u5 = @intCast(23 - top);
    const frac = (x << shift) & FRAC; // единица уезжает за край, остаётся мантисса
    return ((top + 127) << 23) | frac;
}

test "отсечение дробной части в середине диапазона" {
    const values = [_]f32{ 0.0, 0.5, 0.99, 1.0, 1.5, 2.9, 42.75, 4194303.5 };
    for (values) |x| {
        const expected: i32 = @intFromFloat(x);
        try std.testing.expectEqual(expected, f2iCore(@bitCast(x)));
    }
}

test "целое в float, пока оно влезает в мантиссу" {
    var x: u32 = 1;
    while (x < (1 << 24)) : (x = x * 3 + 1) {
        const expected: f32 = @floatFromInt(x);
        try std.testing.expectEqual(@as(u32, @bitCast(expected)), i2fCore(x));
    }
    try std.testing.expectEqual(@as(u32, @bitCast(@as(f32, 16777215))), i2fCore(16777215));
}

test "за границей 24 бит точное представление кончается" {
    // 2^24 + 1 в мантиссу уже не влезает: нужен ещё один бит,
    // и настоящее преобразование округляет к ближайшему чётному.
    const n: u32 = (1 << 24) + 1;
    const rounded: f32 = @floatFromInt(n);
    try std.testing.expectEqual(@as(f32, 16777216), rounded);
    // i2fCore сюда не годится: у него нет ни округления, ни места под лишний бит.
    try std.testing.expect(31 - @clz(n) > 23);
}

$ zig test convert.zig
1/3 convert.test.отсечение дробной части в середине диапазона...OK
2/3 convert.test.целое в float, пока оно влезает в мантиссу...OK
3/3 convert.test.за границей 24 бит точное представление кончается...OK
All 3 tests passed.

В f2iCore работа делается одной строкой: восстановить двадцать четыре бита значащей части, поставив подразумеваемую единицу на место, и сдвинуть вправо на столько разрядов, сколько осталось от мантиссы после целой части. Сдвиг вправо и есть отсечение дробной части в сторону нуля, ничего дополнительно отбрасывать не надо. Обрати внимание, что округления здесь нет и быть не должно: перевод к целому обязан именно усекать, поэтому 2.9 даёт 2, а не 3.

В i2fCore идём обратно. @clz считает старшие нули, значит 31 - @clz(x) это позиция старшей единицы, и она же показатель степени двойки. Сдвигаем число влево так, чтобы эта единица встала на двадцать третий бит, маской срезаем её саму (она подразумевается) и собираем шаблон.

Обе функции живут в безопасной середине. За её краями начинается настоящая работа: у f2iCore это отрицательные числа, значения по модулю меньше единицы, значения, которые не влезают в тридцать два бита, а также NaN и бесконечности, для которых процессор возвращает специальный ответ. У i2fCore это целые, которым нужно больше двадцати четырёх бит: там появляется округление, и перенос при округлении может увеличить порядок на единицу. Третий тест показывает границу буквально: 2²⁴ + 1 приезжает в 2²⁴, ровно по правилу чётности.

Patriot: цена накопленной ошибки

Двадцать пятого февраля 1991 года в Дахране ракета Scud попала в казарму, где размещались американские военнослужащие. Погибли двадцать восемь человек, около сотни ранены. Батарея Patriot, прикрывавшая объект, цель не перехватила: система её обнаружила, но следящий радар искал её не там.

Причина оказалась в представлении времени. Внутренние часы системы считали десятые доли секунды целым числом тиков, а для расчётов время переводилось в секунды умножением на одну десятую. Коэффициент лежал в двадцатичетырёхбитном регистре с фиксированной точкой, а двоичная запись одной десятой бесконечна и периодична: 0.0001100110011001100110011… Регистр обрывал её, и каждый тик приносил маленькую ошибку. Ошибка была бы безобидна, если бы её не умножали на число тиков.

Повторим счёт.

const std = @import("std");

/// Регистр хранит 24 бита, но двоичная запись одной десятой начинается с трёх
/// нулей: 0.0001100110011... Значащих разрядов под саму дробь остаётся 21,
/// а бесконечный периодический хвост отбрасывается.
const KEPT: u6 = 21;
const SCALE: u64 = @as(u64, 1) << KEPT;
const STORED_TENTH: u64 = @intFromFloat(0.1 * @as(f64, @floatFromInt(SCALE)));

/// Часы умножали целое число десятых долей секунды на этот приближённый коэффициент.
fn secondsStored(ticks: u64) f64 {
    return @as(f64, @floatFromInt(ticks * STORED_TENTH)) / @as(f64, @floatFromInt(SCALE));
}

pub fn main(init: std.process.Init) !void {
    var buf: [4096]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    const stored = @as(f64, @floatFromInt(STORED_TENTH)) / @as(f64, @floatFromInt(SCALE));
    try out.print("двоичная запись 0.1:  0.0001100110011001100110011...\n", .{});
    try out.print("сохранено в регистре: 0.000110011001100110011\n", .{});
    try out.print("это в десятичной:     {d:.16}\n", .{stored});
    try out.print("ошибка одного тика:   {e}\n\n", .{0.1 - stored});

    try out.print("{s}\n", .{"часов     тиков   уход часов, с   промах по цели, м"});
    for ([_]u64{ 1, 8, 20, 50, 100 }) |h| {
        const ticks = h * 3600 * 10;
        const drift = @as(f64, @floatFromInt(ticks)) / 10.0 - secondsStored(ticks);
        try out.print("{d:>5}{d:>10}{d:>15.4}{d:>18.0}\n", .{ h, ticks, drift, drift * 1676.0 });
    }
    try out.flush();
}

$ zig run patriot.zig
двоичная запись 0.1:  0.0001100110011001100110011...
сохранено в регистре: 0.000110011001100110011
это в десятичной:     0.0999999046325684
ошибка одного тика:   9.536743164617612e-8

часов     тиков   уход часов, с   промах по цели, м
    1     36000         0.0034                 6
    8    288000         0.0275                46
   20    720000         0.0687               115
   50   1800000         0.1717               288
  100   3600000         0.3433               575

Батарея в Дахране работала без перезапуска около ста часов. Уход часов к этому моменту составил примерно треть секунды. Scud летит со скоростью около 1676 метров в секунду, значит следящий радар искал цель почти в шестистах метрах от того места, где она была. Окно, в котором система подтверждает захват, оказалось пустым, и цель списали как ложную.

Отдельно горько то, что ошибку уже знали. Часть кода к тому моменту исправили на более точный перевод времени, но исправили не везде. Пока обе ветки считали одинаково плохо, ошибка сокращалась при вычитании двух моментов времени, и система работала. Когда одна ветка стала точной, а другая осталась старой, разность перестала сокращаться, и ошибка вылезла наружу целиком. Обновлённая версия шла в Дахран и прибыла на следующий день после удара.

Три вывода, которые стоит унести.

Первое: ошибка представления не случайна и не мала, она систематична. Одно и то же отклонение, умноженное на миллионы повторов, растёт линейно.

Второе: время и деньги в дробных единицах лучше держать целыми. Считай в тиках, миллисекундах, копейках, а в дробь переводи только на выходе, когда показываешь человеку.

Третье: смешивать точную и приближённую версию одного вычисления опаснее, чем везде использовать приближённую. Согласованная ошибка часто сокращается, рассогласованная нет.

Явные преобразования в Zig и что происходит на краю

Zig не делает неявных числовых преобразований, поэтому каждый переход между мирами приходится называть по имени.

• @floatFromInt берёт целое и даёт float. Всегда определён: если целое не влезает в мантиссу, результат округляется.
• @intFromFloat берёт float и даёт целое, отсекая дробную часть в сторону нуля. Определён только если целая часть влезает в тип назначения.
• @floatCast меняет ширину float. Всегда определён: при сужении результат округляется, а слишком крупное значение становится бесконечностью, слишком мелкое нулём.

const std = @import("std");

pub fn main(init: std.process.Init) !void {
    var buf: [4096]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    // @floatCast сужает и округляет; за краем диапазона это бесконечность.
    const big: f64 = 1e300;
    const narrowed: f32 = @floatCast(big);
    try out.print("@floatCast(1e300) в f32: {e}\n", .{narrowed});
    const tiny: f64 = 1e-300;
    try out.print("@floatCast(1e-300) в f32: {e}\n", .{@as(f32, @floatCast(tiny))});
    const pi: f64 = 3.141592653589793;
    try out.print("@floatCast(pi) в f32: {d}, обратно в f64: {d}\n", .{
        @as(f32, @floatCast(pi)), @as(f64, @as(f32, @floatCast(pi))),
    });

    // @floatFromInt всегда определён: если не влезло, округляет.
    const n: i64 = 9007199254740993; // 2^53 + 1
    try out.print("\n@floatFromInt(2^53 + 1) в f64: {d}\n", .{@as(f64, @floatFromInt(n))});

    // @intFromFloat усекает в сторону нуля.
    for ([_]f64{ 2.9, -2.9, 0.5, -0.5 }) |x| {
        try out.print("@intFromFloat({d}) = {d}\n", .{ x, @as(i32, @intFromFloat(x)) });
    }

    // Безопасный путь для значения, которое может не влезть.
    try out.print("\nlossyCast(1e300 -> i32) = {d}\n", .{std.math.lossyCast(i32, 1e300)});
    try out.print("lossyCast(NaN -> i32)   = {d}\n", .{std.math.lossyCast(i32, std.math.nan(f64))});
    try out.flush();
}

$ zig run casts.zig
@floatCast(1e300) в f32: inf
@floatCast(1e-300) в f32: 0e0
@floatCast(pi) в f32: 3.1415927, обратно в f64: 3.1415927410125732

@floatFromInt(2^53 + 1) в f64: 9007199254740992
@intFromFloat(2.9) = 2
@intFromFloat(-2.9) = -2
@intFromFloat(0.5) = 0
@intFromFloat(-0.5) = 0

lossyCast(1e300 -> i32) = 2147483647
lossyCast(NaN -> i32)   = 0

Строка с числом пи заслуживает внимания: сужение до f32 и возврат обратно в f64 не даёт исходное значение. Обратно приезжает 3.1415927410125732, потому что после сужения младшие разряды мантиссы потеряны навсегда. Сужение это не переименование, а честная потеря данных, в отличие от @bitCast.

Теперь про край. @intFromFloat со значением, которое не влезает в тип назначения, это неопределённое поведение. Zig не оставляет это без присмотра: в отладочных режимах вставлена проверка. Посмотрим на все три режима сразу.

const std = @import("std");

/// Значение лежит в изменяемой переменной и читается через volatile,
/// поэтому компилятор не сворачивает вызов на этапе компиляции.
var input: f64 = 1e30;

fn toI32(x: f64) i32 {
    return @intFromFloat(x);
}

pub fn main(init: std.process.Init) !void {
    var buf: [512]u8 = undefined;
    var w = std.Io.File.stdout().writer(init.io, &buf);
    const out = &w.interface;

    const x = @as(*volatile f64, &input).*;
    try out.print("{d} -> {d}\n", .{ x, toI32(x) });
    try out.flush();
}

$ zig run -ODebug overflow.zig
thread 8295858 panic: integer part of floating point value out of bounds
/private/tmp/lesson08/overflow.zig:8:12: 0x10433e1f7 in toI32 (overflow)
    return @intFromFloat(x);

$ zig run -OReleaseSafe overflow.zig
thread 8296223 panic: integer part of floating point value out of bounds
/private/tmp/lesson08/overflow.zig:8:12: 0x102312983 in toI32 (overflow)
    return @intFromFloat(x);

$ zig run -OReleaseFast overflow.zig
1000000000000000000000000000000 -> 2147483647

$ zig run -OReleaseFast -target x86_64-macos overflow.zig
1000000000000000000000000000000 -> -2147483648

В Debug и ReleaseSafe программа падает с внятным сообщением: целая часть значения вне границ типа. В ReleaseFast проверки нет, и наружу вылезает то, что сделала машина. И вот тут самое поучительное: на процессоре aarch64 получилось 2147483647, а тот же исходник, собранный под x86-64, дал -2147483648. Инструкции преобразования на этих архитектурах ведут себя по-разному, одна насыщает результат по верхней границе, другая возвращает специальное значение. Неопределённое поведение здесь не абстракция из стандарта, а буквально разный ответ на одинаковый вопрос.

Если значение может не влезть, проверяй границы сам или бери std.math.lossyCast: он насыщает результат по краям типа и превращает NaN в ноль, зато определён всегда. В тонкости того, во что каждая из этих операций превращается на уровне инструкций, зачем там лишняя буква t в мнемонике и как процессор сравнивает числа с NaN, мы залезем в уроке про float и SIMD: там те же числа, но глазами процессора.

Практика

Две задачи в жанре домашних заданий главы 2 CS:APP. Правила в обеих одинаковые и жёсткие: аргументы и результаты это битовые шаблоны u32, тип f32 тебе недоступен, арифметики с плавающей точкой в решении быть не должно. Разрешены целочисленные операции, сдвиги, маски, сравнения, if и switch. Проверка на это текстовая: отдельный тест читает исходник и ищет имена типов с плавающей точкой и встроенные функции преобразования, поэтому не упоминай их даже в комментариях. Эталон внутри теста считается настоящими операциями над f32, так что расхождение хотя бы в одном бите видно сразу.

Первая задача про знак, модуль, сравнение и степень двойки. floatNegate возвращает биты числа с обратным знаком, floatAbsval биты его модуля, и обе обязаны вернуть аргумент нетронутым, если он NaN. floatLe отвечает на вопрос «первое меньше либо равно второму», и здесь пригодится всё из раздела про порядок битов: NaN выпадает из любого сравнения, два нуля равны друг другу, а у отрицательных чисел порядок шаблонов перевёрнут. fpwr2 собирает биты двойки в заданной целой степени и обязан правильно разложиться на четыре области: бесконечность сверху, нормальные числа, денормали и ноль снизу.

Вторая задача сложнее, потому что в ней появляется округление. floatTwice и floatHalf это удвоение и деление пополам до последнего бита, включая переход денормали в нормальное число, переполнение наибольшего конечного в бесконечность и округление к ближайшему чётному при уходе в денормали. floatF2i это отсечение дробной части в сторону нуля с ответом 0x8000_0000 для NaN, бесконечностей и всего, что не влезает в тридцать два бита со знаком: ровно так поступает и процессор. floatI2f идёт в обратную сторону и снова требует округления, потому что тридцать два бита целого в двадцать четыре бита мантиссы не помещаются.

Упражнения

Итоги

• Число с плавающей точкой это три поля в одном машинном слове: знак, порядок со смещением, мантисса без старшей единицы. Значение это значащая часть, умноженная на два в степени порядка.
• Смещение вместо дополнительного кода нужно для того, чтобы битовые шаблоны неотрицательных чисел шли в том же порядке, что и значения. На этом свойстве строится целочисленное сравнение float.
• Два значения поля порядка зарезервированы. Всё нули это ноль и денормали, все единицы это бесконечность и NaN. Пять классов значений получаются из одного правила, а не из списка исключений.
• Сетка представимых значений неравномерна: на каждой границе порядка шаг удваивается. Все целые до 2²⁴ точно представимы в f32, до 2⁵³ в f64, дальше начинаются пропуски.
• По умолчанию округление идёт к ближайшему, а при ничьей к соседу с чётной мантиссой. Это убирает систематический сдвиг при накоплении. Переполнение это не отдельная проверка, а тот же выбор ближайшего соседа, которым оказалась бесконечность.
• Денормали не имеют скрытой единицы и продолжают сетку первых нормальных чисел без разрыва. Благодаря им из различия двух чисел следует ненулевая разность.
• NaN не равен ничему, включая себя, поэтому x != x истинно ровно для него. Свободные биты мантиссы у NaN это полезная нагрузка, и на ней держится приём упаковки значений в NaN.
• @bitCast между f32 и u32 (или f64 и u64) это переименование битов, а не преобразование. На нём пишутся операции над float целочисленным кодом: умножение на два это плюс единица к полю порядка, деление пополам минус единица, а у денормалей то же самое делает обычный сдвиг.
• @floatFromInt и @floatCast определены всегда и при нехватке разрядов округляют. @intFromFloat усекает в сторону нуля и определён только в границах типа: за ними Debug и ReleaseSafe роняют программу, а ReleaseFast даёт разный ответ на разных архитектурах.
• Ошибка представления систематична, поэтому накапливается линейно. Часы Patriot ушли на треть секунды за сто часов работы, и этого хватило, чтобы радар искал цель в шестистах метрах от неё.

Дальше

Блок про биты и числа закрыт: ты знаешь, как в машинном слове лежат и целые, и дроби, и умеешь работать с обоими на уровне отдельных разрядов. Дальше начинается машинный уровень. Мы возьмём обычный Zig-код, снимем с него настоящий ассемблер x86-64 и научимся читать его как исходник: регистры и их части, форматы данных, что меняется между Debug и ReleaseFast в дизассемблере одной функции. Первый урок блока это путь от Zig до машинного кода, и инструменты оттуда останутся с тобой до конца раздела.