Числа с плавающей точкой и SIMD

До сих пор мы читали ассемблер целочисленного кода: аргументы в rdi и rsi, результат в rax, арифметика через add и imul. Числа с плавающей точкой живут по другому адресу. У них свой регистровый файл, свои инструкции и своя логика сравнения, где NaN ломает привычный порядок. В этом уроке снимем с Zig-кода настоящий ассемблер float-операций, разберём, как знак и модуль делаются битовой маской, а сравнение с NaN уходит в отдельную ветку. Потом шагнём дальше книги: те же регистры шире целого слова, и один YMM держит восемь float сразу. @Vector в Zig это способ занять всю ширину и посчитать восемь результатов одной инструкцией.

Цели урока

• Знать, что float живут в отдельном файле регистров XMM и YMM, а не в rax и его родне.
• Читать скалярные инструкции addsd, mulsd, movsd и понимать суффиксы sd и ss.
• Снимать преобразования cvtsi2sd, cvttsd2si, cvtss2sd с Zig-операторов @floatFromInt, @intFromFloat, @floatCast и по ассемблеру восстанавливать сигнатуру.
• Видеть, что float-константа лежит в .rodata и грузится со смещением от rip.
• Понимать смену знака и модуль как xorps и andps с битовой маской, а не как арифметику.
• Читать сравнение ucomisd и знать, почему NaN проверяют через флаг чётности PF.
• Собирать @Vector(4, f64) и @Vector(8, f32) и видеть, во что они компилируются под AVX2.
• Понимать @reduce(.Add, v) как дерево свёрток и знать, почему обычный цикл суммирования не векторизуется.

Идея: у float свой регистровый файл

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

Целочисленный код весь урок про машинный уровень крутился вокруг регистров общего назначения: rax, rdi, rsi и так далее. Float туда не попадают. У процессора есть второй набор регистров, XMM: шестнадцать регистров xmm0 до xmm15, каждый 128 бит. С расширением AVX те же регистры видны как YMM, по 256 бит. System V ABI кладёт float-аргументы в xmm0 и дальше, а результат возвращает в xmm0, ровно как целые возвращаются в rax.

Одна инструкция работает либо над всем регистром (пакетно, packed), либо только над младшим числом (скалярно, scalar). Урок начнём со скалярных операций, где занят только младший слот, а в конце дойдём до пакетных, где заняты все.

Скалярные операции: addsd и его семья

Возьмём три функции: сложение двух f64, полином x * a + b и сложение двух f32. Каждая помечена export, чтобы оптимизатор не встроил её и оставил как отдельную функцию с C-совместимой сигнатурой.

// Скалярные операции с f64 и f32.
export fn add_d(a: f64, b: f64) f64 {
    return a + b;
}
export fn poly_d(x: f64, a: f64, b: f64) f64 {
    return x * a + b;
}
export fn add_s(a: f32, b: f32) f32 {
    return a + b;
}

Снимаем ассемблер. Напоминаю рецепт из урока про машинный код: -femit-asm печатает Intel-синтаксис, поэтому AT&T берём из objdump.

$ zig build-obj scalar.zig -O ReleaseFast -target x86_64-linux \
      -fomit-frame-pointer -femit-bin=scalar.o
$ objdump -d --no-show-raw-insn -M att scalar.o

add_d:
    addsd   %xmm1, %xmm0        # xmm0 = xmm0 + xmm1
    ret

poly_d:
    mulsd   %xmm1, %xmm0        # xmm0 = x * a
    addsd   %xmm2, %xmm0        # xmm0 = (x * a) + b
    ret

add_s:
    addss   %xmm1, %xmm0        # xmm0 = xmm0 + xmm1, но 32 бита
    ret

Разбери по пунктам.

1. Аргументы пришли не в rdi и rsi, а в xmm0, xmm1, xmm2. У poly_d это x, a, b по порядку. Результат уходит в xmm0.
2. addsd это add scalar double: сложить два f64 в младших 64 битах. mulsd это умножение. Суффикс sd значит f64.
3. У add_s мнемоника addss, add scalar single: то же самое, но f32 в младших 32 битах. Суффикс ss значит f32.
4. poly_d это две инструкции без единого перехода: сначала mulsd, потом addsd. Ровно та же форма, что целочисленный x * a + b, только в другом файле регистров.

Полная семья скалярных арифметических инструкций: addsd, subsd, mulsd, divsd для f64 и те же с суффиксом ss для f32. Пересылка это movsd и movss. Здесь пересылки нет, потому что всё уже в нужных регистрах, но она появится, как только понадобится константа из памяти.

Преобразования: cvt между целыми и float

Целое и float это разные представления одного числа, и перевод между ними это отдельная инструкция, а не переливание битов. В Zig перевод всегда явный: @floatFromInt и @intFromFloat. Добавим @floatCast для смены ширины float.

// Преобразования между целыми и числами с плавающей точкой.
export fn i_to_d(n: i64) f64 {
    return @floatFromInt(n);
}
export fn d_to_i(x: f64) i64 {
    return @intFromFloat(x);
}
export fn u_to_d(n: u32) f64 {
    return @floatFromInt(n);
}
export fn s_to_d(x: f32) f64 {
    return @floatCast(x);
}
export fn d_to_s(x: f64) f32 {
    return @floatCast(x);
}

i_to_d:                         # i64 -> f64
    cvtsi2sd    %rdi, %xmm0
    ret

d_to_i:                         # f64 -> i64
    cvttsd2si   %xmm0, %rax
    ret

u_to_d:                         # u32 -> f64
    movl        %edi, %eax      # обнулить старшие 32 бита rax
    cvtsi2sd    %rax, %xmm0     # и перевести как знаковое 64-битное
    ret

s_to_d:                         # f32 -> f64
    cvtss2sd    %xmm0, %xmm0
    ret

d_to_s:                         # f64 -> f32
    cvtsd2ss    %xmm0, %xmm0
    ret

Мнемоники читаются как формула. cvtsi2sd это convert scalar integer to scalar double: целое из регистра общего назначения в f64. Обратно cvttsd2si это convert with truncation scalar double to scalar integer. Лишняя буква t это усечение: дробная часть отбрасывается в сторону нуля, а не округляется. Это ловушка: @intFromFloat(2.9) даёт 2, @intFromFloat(-2.9) даёт -2. cvtss2sd и cvtsd2ss переводят между f32 и f64.

Отдельно посмотри на u_to_d. Аргумент это u32, беззнаковое 32-битное. У инструкции cvtsi2sd нет беззнакового варианта, она читает знаковое целое. Компилятор обходит это так: movl %edi, %eax копирует младшие 32 бита и заодно обнуляет старшие 32 (любая запись в 32-битный регистр обнуляет верх), а потом cvtsi2sd %rax переводит уже 64-битное знаковое, где старший бит гарантированно ноль. Так беззнаковое u32 честно становится положительным f64, даже если старший бит был единицей.

Реверс: восстанови сигнатуру по cvt

Здесь работает тот же приём, что в упражнениях книги 3.51 и 3.52: по ассемблеру восстановить типы. Мнемоника cvt несёт в себе обе стороны преобразования, а обходные движения вокруг неё выдают знаковость и ширину. Смотри на u_to_d выше: если бы аргумент был i32, компилятор написал бы movslq %edi, %rax (расширение со знаком) или прямо cvtsi2sd %edi. Обнуление старших через movl вместо расширения со знаком это подпись беззнакового 32-битного.

Разбери ещё один листинг. Что за функция?

mystery:
    cvttsd2si   %xmm0, %eax
    ret

Пойдём с конца. Результат в eax, это 32-битный регистр, значит возвращается i32. Вход в xmm0, cvttsd2si с усечением, значит на входе f64. Инструкция усекающая, значит это @intFromFloat, а не округление. Сигнатура: fn(x: f64) i32 с телом return @intFromFloat(x);. Разница с d_to_i только в ширине результата: там rax и i64, здесь eax и i32.

Константа в памяти грузится от rip

Float-литерал нельзя закодировать внутри инструкции так же, как целое $42. У float-операций нет формы с непосредственным операндом. Поэтому константа лежит в секции .rodata, а инструкция ссылается на неё адресом.

export fn scale(x: f64) f64 {
    return x * 3.14159;
}

scale:
    mulsd   .LC0(%rip), %xmm0     # умножить на константу из .rodata
    ret

Второй операнд mulsd это не регистр, а память: .LC0(%rip), адрес относительно счётчика команд. Такая rip-относительная адресация это стандартный способ x86-64 дотянуться до глобальных данных. objdump с флагом релокаций показывает, что там лежит: символ .LCPI2_0 из .rodata.cst8. Посмотрим байты этой константы.

$ objdump -s -j .rodata.cst8 masks.o
Contents of section .rodata.cst8:
 0000 6e861bf0 f9210940                    n....!.@

Восемь байт 6e 86 1b f0 f9 21 09 40 это и есть 3.14159 в формате f64, little-endian. Проверить можно битовым разбором из урока про float: старший байт 0x40 даёт экспоненту около единицы, дальше мантисса. Никакой магии: literal превратился в восемь байт в секции только для чтения, а mulsd читает их прямо из памяти.

Знак и модуль это битовая маска

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

export fn negate(x: f64) f64 {
    return -x;
}
export fn absval(x: f64) f64 {
    return @abs(x);
}

negate:
    xorps   .LC0(%rip), %xmm0     # перевернуть знаковый бит
    ret

absval:
    andps   .LC1(%rip), %xmm0     # сбросить знаковый бит
    ret

Ни subsd, ни mulsd, а битовые xorps и andps над регистром и маской из памяти. Маски лежат в .rodata.cst16.

$ objdump -s -j .rodata.cst16 masks.o
Contents of section .rodata.cst16:
 0000 ffffffff ffffff7f ffffffff ffffff7f  ................
 0010 00000000 00000080 00000000 00000080  ................

Первая маска по смещению 0 это 0x7fffffffffffffff: все биты единицы, кроме старшего. andps с ней оставляет число как есть, но гасит знаковый бит, и результат всегда положительный. Это @abs. Вторая маска по смещению 0x10 это 0x8000000000000000: только старший бит. xorps с ней переворачивает знак и не трогает остальное. Это унарный минус.

Две детали. Первая: мнемоники кончаются на ps, а не pd, хотя число это f64. xorps и andps это упакованные битовые операции: они не смотрят на числа как на float, а просто гоняют биты. Вариант ps на один байт короче, чем pd, поэтому компилятор берёт его. Для чистого перебора битов разницы нет. Вторая: маска шириной 16 байт (128 бит), хотя число 8 байт. Так удобнее выравнивать константу, а лишние биты всё равно не используются, ведь регистр читается как один f64.

Проверим, что маска и правда равна арифметике. Переворот старшего бита должен совпасть с -x, а сброс с @abs.

const std = @import("std");

// Смена знака это переворот старшего бита, ровно то, что делает xorps с маской.
fn negBits(x: f64) f64 {
    const sign_mask: u64 = 0x8000_0000_0000_0000;
    return @bitCast(@as(u64, @bitCast(x)) ^ sign_mask);
}

// Модуль это сброс старшего бита, ровно то, что делает andps с маской.
fn absBits(x: f64) f64 {
    const abs_mask: u64 = 0x7fff_ffff_ffff_ffff;
    return @bitCast(@as(u64, @bitCast(x)) & abs_mask);
}

test "переворот знакового бита совпадает с унарным минусом" {
    const values = [_]f64{ 3.5, -2.0, 0.0, 1e300 };
    for (values) |v| {
        try std.testing.expectEqual(@as(u64, @bitCast(-v)), @as(u64, @bitCast(negBits(v))));
    }
}

test "сброс знакового бита совпадает с @abs" {
    const values = [_]f64{ 3.5, -2.0, -1e300 };
    for (values) |v| {
        try std.testing.expectEqual(@abs(v), absBits(v));
    }
}

$ zig test masks_check.zig
1/2 masks_check.test.переворот знакового бита совпадает с унарным минусом...OK
2/2 masks_check.test.сброс знакового бита совпадает с @abs...OK
All 2 tests passed.

Ручной переворот бита через @bitCast и ^ даёт те же биты, что -x, а сброс через & те же, что @abs. Компилятор не выдумывает трюк, он просто знает про представление знака.

Сравнение: ucomisd и флаг чётности для NaN

Сравнение float опаснее целочисленного, потому что среди значений есть NaN, который не равен ничему, даже себе. Инструкция сравнения ucomisd кладёт результат не в отдельный регистр, а во флаги, и специально для NaN задействует флаг чётности PF. Возьмём несколько предикатов.

export fn lt(a: f64, b: f64) bool {
    return a < b;
}
export fn ge(a: f64, b: f64) bool {
    return a >= b;
}
export fn is_nan(x: f64) bool {
    return x != x;
}

lt:                             # a < b
    ucomisd %xmm0, %xmm1        # сравнить b с a
    seta    %al                 # a < b это b > a
    ret

ge:                             # a >= b
    ucomisd %xmm1, %xmm0
    setae   %al
    ret

is_nan:                         # x != x
    ucomisd %xmm0, %xmm0        # сравнить x сам с собой
    setp    %al                 # взять флаг чётности
    ret

ucomisd сравнивает два f64 и ставит флаги ZF, PF, CF, а set* превращает флаг в булев байт. Для a < b компилятор сравнивает b с a и берёт условие above (seta), потому что регистр флагов после ucomisd устроен как после беззнакового целочисленного сравнения. Но самое интересное это is_nan. Выражение x != x истинно ровно тогда, когда x это NaN. Компилятор сравнивает x сам с собой и берёт setp: флаг чётности. ucomisd ставит PF в единицу, когда любой операнд неупорядочен, то есть NaN. Так проверка на NaN это одна инструкция сравнения плюс чтение одного флага.

Теперь посмотрим, зачем нужен переход по чётности. Равенство float в ветке требует проверить два флага сразу.

extern fn on_equal() void;
extern fn on_diff() void;
// Равенство в ветке с вызовом: NaN идёт в on_diff, для этого нужен jp.
export fn route(a: f64, b: f64) void {
    if (a == b) on_equal() else on_diff();
}

route:
    ucomisd %xmm1, %xmm0
    jne     .diff              # значения разные по ZF: в ветку различия
    jp      .diff              # неупорядочены (NaN): туда же
    jmp     on_equal           # равны и не NaN: хвостовой переход
.diff:
    jmp     on_diff

Вот зачем нужен jp, переход по флагу чётности. Равенство a == b для float это два условия: значения совпадают (ZF установлен) и оба не NaN. Компилятор проверяет их по очереди: jne уводит в ветку различия, если значения разные, а jp уводит туда же, если операнды неупорядочены. Только когда оба перехода не сработали, управление доходит до on_equal. Это классический след из книги: сравнение float с NaN всегда идёт через jp, потому что NaN обязан выпасть из равенства, каким бы ни был второй операнд.

Соберём поведение NaN в один тест, чтобы увидеть его на значениях, а не только в ассемблере.

const std = @import("std");

test "NaN не равен ничему, включая себя" {
    const nan = std.math.nan(f64);
    try std.testing.expect(nan != nan);
    try std.testing.expect(!(nan < 0));
    try std.testing.expect(!(nan > 0));
    try std.testing.expect(!(nan == nan));
}

test "минус ноль равен плюс нулю по значению, но не по битам" {
    const pos: f64 = 0.0;
    const neg: f64 = -0.0;
    try std.testing.expect(pos == neg);
    try std.testing.expect(@as(u64, @bitCast(pos)) != @as(u64, @bitCast(neg)));
}

$ zig test nan_check.zig
1/2 nan_check.test.NaN не равен ничему, включая себя...OK
2/2 nan_check.test.минус ноль равен плюс нулю по значению, но не по битам...OK
All 2 tests passed.

Все сравнения с NaN ложны, а != истинно. И отдельная тонкость: минус ноль и плюс ноль равны по значению, но их биты различаются знаковым битом. Поэтому сравнение и битовое равенство это не одно и то же.

Дальше книги: один регистр, восемь чисел

До сих пор мы занимали только младший слот XMM-регистра. Но регистр широкий: 128 бит XMM или 256 бит YMM. Если положить туда не одно число, а несколько, одна инструкция посчитает их все разом. Это SIMD: одна инструкция, много данных.

Книга SIMD откладывает, а мы дойдём. В Zig векторный тип это @Vector(N, T): N значений типа T, лежащих в регистре как дорожки. Обычная арифметика над вектором работает поэлементно: дорожка результата зависит только от своей дорожки входов.

const V4 = @Vector(4, f64);
const V8 = @Vector(8, f32);
// Поэлементное сложение четырёх f64.
export fn vadd4(a: V4, b: V4) V4 {
    return a + b;
}
// Поэлементное умножение восьми f32.
export fn vmul8(a: V8, b: V8) V8 {
    return a * b;
}

Чтобы компилятор выпустил AVX2 и занял весь YMM, укажем уровень микроархитектуры x86_64_v3. Это baseline, где AVX2 гарантированно есть.

$ zig build-obj simd.zig -O ReleaseFast -target x86_64-linux \
      -mcpu=x86_64_v3 -fomit-frame-pointer -femit-bin=simd.o
$ objdump -d --no-show-raw-insn -M att simd.o

vadd4:                          # 4 × f64
    vaddpd  %ymm1, %ymm0, %ymm0    # одна инструкция на все четыре
    ret

vmul8:                          # 8 × f32
    vmulps  %ymm1, %ymm0, %ymm0    # одна инструкция на все восемь
    ret

Смотри, что произошло. vadd4 это одна инструкция vaddpd: add packed double, сложить четыре f64 попарно за раз. Регистр теперь ymm0, полные 256 бит. vmul8 это vmulps: mul packed single, восемь умножений f32 одной инструкцией. Префикс v это AVX-форма с тремя операндами: два источника и приёмник по отдельности. Суффикс ps это packed single (f32), pd это packed double (f64).

Что было бы без AVX. Четыре f64 это 256 бит, а XMM всего 128, поэтому вектор пришлось бы держать в двух регистрах и складывать двумя addpd, а сам вектор ABI передал бы через стек. AVX даёт регистр нужной ширины, и всё умещается в одну инструкцию. Разница в цене прямая: одна операция вместо нескольких плюс возня со стеком.

Покрути дорожки в виджете. Выбери ширину (четыре f64 или восемь f32), операцию (сложение, умножение, максимум) и жми «шаг»: результат заполняется дорожка за дорожкой, и видно, что дорожка i зависит только от i, а соседи в неё не смотрят. Внизу переключатель «скаляр против вектора» на разных N: он считает, сколько инструкций нужно скалярному коду и сколько векторному. На тысяче f64 вектор шириной четыре делает вчетверо меньше.

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

const std = @import("std");
const V4 = @Vector(4, f64);

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

    const a: V4 = .{ 1.0, 2.0, 3.0, 4.0 };
    const b: V4 = .{ 10.0, 20.0, 30.0, 40.0 };

    const sum = a + b;
    const prod = a * b;
    const hi = @max(a, b);

    try out.print("a + b   = {d}\n", .{sum});
    try out.print("a * b   = {d}\n", .{prod});
    try out.print("max     = {d}\n", .{hi});
    try out.print("reduce +  = {d}\n", .{@reduce(.Add, a)});
    try out.print("reduce max= {d}\n", .{@reduce(.Max, b)});
    try out.flush();
}

$ zig run vec_demo.zig
a + b   = { 11, 22, 33, 44 }
a * b   = { 10, 40, 90, 160 }
max     = { 10, 20, 30, 40 }
reduce +  = 10
reduce max= 40

a + b, a * b и @max(a, b) возвращают вектор: по дорожкам. А вот @reduce возвращает одно число, схлопывая весь вектор. Это уже не поэлементная операция, и ассемблер у неё другой.

Горизонтальная сумма: @reduce как дерево

Сложить дорожки между собой сложнее, чем сложить два вектора. Поэлементное сложение независимо, а горизонтальная сумма перемешивает дорожки внутри регистра. Процессор не умеет сложить все дорожки одной инструкцией, поэтому @reduce(.Add, v) разворачивается в дерево.

const V4 = @Vector(4, f64);
export fn hsum4(a: V4) f64 {
    return @reduce(.Add, a);
}

hsum4:                                  # сумма четырёх дорожек f64
    vshufpd     $1, %xmm0, %xmm0, %xmm1   # переставить дорожки 0 и 1
    vaddsd      %xmm1, %xmm0, %xmm1       # сложить пару в нижней половине
    vextractf128 $1, %ymm0, %xmm0         # достать верхние 128 бит
    vaddsd      %xmm0, %xmm1, %xmm1       # прибавить
    vshufpd     $1, %xmm0, %xmm0, %xmm0   # снова переставить
    vaddsd      %xmm0, %xmm1, %xmm0       # последнее сложение
    ret

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

Ту же структуру видно, если суммировать длинный срез вектором. Копим в аккумуляторе из четырёх дорожек, а в конце сворачиваем @reduce.

const V4 = @Vector(4, f64);
// Аккумулятор из четырёх независимых частичных сумм: одна vaddpd на итерацию.
export fn sum_vec(ptr: [*]const f64, len: usize) f64 {
    var acc: V4 = @splat(0);
    var i: usize = 0;
    while (i + 4 <= len) : (i += 4) {
        const chunk: V4 = ptr[i..][0..4].*;
        acc += chunk;
    }
    var tail: f64 = @reduce(.Add, acc);
    while (i < len) : (i += 1) tail += ptr[i];
    return tail;
}

sum_vec:
    ...
.loop:
    vaddpd  (%rdi,%rcx,8), %ymm0, %ymm0    # одна vaddpd на четыре элемента
    ...
    jbe     .loop
    vshufpd ...                            # затем дерево @reduce
    vaddsd  ...
    vextractf128 ...

Тело цикла это одна vaddpd: четыре элемента за инструкцию, четыре независимые частичные суммы в четырёх дорожках. Когда данные кончились, дерево свёртки складывает четыре частичные суммы в одну. Это в четыре раза меньше сложений, чем скалярный цикл. Но заметь: мы сами написали вектор. Компилятор не превратил обычный цикл s += x в это. Почему, в следующем разделе.

Когда векторизация не срабатывает

Логично ждать, что компилятор сам увидит цикл суммирования и разложит его по дорожкам. Иногда так и есть, но не с плавающей точкой по умолчанию, и причина глубже оптимизатора.

// Сумма среза: сложение f64 не ассоциативно, поэтому по умолчанию цикл не векторизуется.
export fn sum_strict(ptr: [*]const f64, len: usize) f64 {
    var s: f64 = 0;
    for (ptr[0..len]) |x| s += x;
    return s;
}

sum_strict:
    ...
.loop:
    vaddsd  (%rdi,%rcx,8),  %xmm0, %xmm0    # scalar, дорожка одна
    vaddsd  0x8(%rdi,%rcx,8), %xmm0, %xmm0  # и каждое сложение
    vaddsd  0x10(%rdi,%rcx,8), %xmm0, %xmm0 # зависит от предыдущего
    ...
    jne     .loop

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

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

const std = @import("std");

fn sumSerial(xs: []const f64) f64 {
    var s: f64 = 0;
    for (xs) |x| s += x;
    return s;
}

fn sumLanes(xs: []const f64) f64 {
    const V4 = @Vector(4, f64);
    var acc: V4 = @splat(0);
    var i: usize = 0;
    while (i + 4 <= xs.len) : (i += 4) {
        const chunk: V4 = xs[i..][0..4].*;
        acc += chunk;
    }
    var tail: f64 = @reduce(.Add, acc);
    while (i < xs.len) : (i += 1) tail += xs[i];
    return tail;
}

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

    var data: [1024]f64 = undefined;
    for (&data, 0..) |*slot, i| slot.* = 1.0 / @as(f64, @floatFromInt(i + 1));

    const serial = sumSerial(&data);
    const lanes = sumLanes(&data);
    try out.print("serial = {x}\n", .{serial});
    try out.print("lanes  = {x}\n", .{lanes});
    try out.print("равны ли биты: {}\n", .{serial == lanes});
    try out.print("разница: {e}\n", .{serial - lanes});
    try out.flush();
}

$ zig run assoc.zig
serial = 0x1.e096558f169dfp2
lanes  = 0x1.e096558f169dcp2
равны ли биты: false
разница: 2.6645352591003757e-15

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

Как разрешить переупорядочивание. Первый способ ты уже видел: написать вектор руками, как в sum_vec. Тогда порядок сложений выбираешь ты, а не компилятор. Второй способ: сказать компилятору, что можно жертвовать точным порядком ради скорости, встроенной функцией @setFloatMode(.optimized) в теле функции. Это Zig-эквивалент флага быстрой математики: компилятор получает право переупорядочивать и, где выгодно, векторизовать. Пользуйся осознанно: результат в младших битах поплывёт, а NaN и бесконечности начнут вести себя иначе, поэтому в коде, где важна точная воспроизводимость, этот режим не включают.

Шаг проекта: zt учится числам с плавающей точкой

Всё, что было в этом уроке, для дизассемблера из прошлых шагов сплошная темнота. Вот ztPolyD, то есть poly_d из начала урока, x * a + b, в выводе zt после урока 14:

$ nm fixtures/step_17.o > step_17.nm
$ ./zig-out/bin/zt disasm --syms=step_17.nm fixtures/step_17.o | sed -n '/<ztPolyD>/,/ 1ae:/p'
00000000000001a0 <ztPolyD>:
     1a0: 55                           	pushq	%rbp
     1a1: 48 89 e5                     	movq	%rsp, %rbp
     1a4: f2                           	(bad)
     1a5: 0f                           	(bad)
     1a6: 59                           	popq	%rcx
     1a7: c1 f2 0f                     	shll	$0xf, %edx
     1aa: 58                           	popq	%rax
     1ab: c2                           	(bad)
     1ac: 5d                           	popq	%rbp
     1ad: c3                           	retq
     1ae: 66 90                        	nop

Две инструкции, mulsd и addsd, превратились в семь строк, и среди них два popq, которых в программе нет. Хуже всего строка 1a7: байты c1 f2 0f честно декодируются как сдвиг, и в листинге появляется shll на месте умножения. Если бы ты искал в этом коде ошибку, ты искал бы её не там.

Префикс, который стал частью опкода

Посмотри на байты двух инструкций урока рядом: f2 0f 58 c1 это addsd %xmm1, %xmm0, а f3 0f 58 c1 это addss %xmm1, %xmm0. Опкод один и тот же, 0f 58. Разница в первом байте, и этот байт тебе знаком: f2 и f3 это префиксы-повторители repne и rep из строковых инструкций, а 66 из шага 10 меняет ширину операнда на 16 бит. Когда Intel добавляла SSE, свободных опкодов почти не осталось, и старые префиксы получили перед 0f новый смысл: они выбирают тип данных.

Буквы в хвосте мнемоники это тот же выбор словами: p значит упакованные (packed), s значит скаляр (scalar), s и d в конце значат одинарную и двойную точность. Скалярные инструкции урока трогают только младший слот регистра, упакованные все слоты сразу, а кодируются они одним опкодом.

Остальное устройство инструкции старое. За опкодом идёт ModRM, поле reg это номер регистра %xmm, поле r/m это второй %xmm или адрес в памяти в любой из одиннадцати форм, бит REX.R и бит REX.B добавляют к номерам четвёртый бит, и получаются %xmm8 до %xmm15. Курсору не нужно ничего нового: он и так запоминает префиксы в Prefixes, а readModRm принимает ширину и для регистровой формы подставляет её. Ширина .xmm заведена в registers.zig ещё в уроке 9, с пометкой, что имена появятся здесь. Это и есть первая правка шага:

const names_xmm = [16][]const u8{
    "xmm0", "xmm1", "xmm2",  "xmm3",  "xmm4",  "xmm5",  "xmm6",  "xmm7",
    "xmm8", "xmm9", "xmm10", "xmm11", "xmm12", "xmm13", "xmm14", "xmm15",
};

и в Register.name вместо unreachable:

            .xmm => names_xmm[register.index],

Таблица на четыре колонки

Всё остальное в новом файле src/x86/ops_sse.zig. У каждого опкода здесь не имя, а четыре имени, по одному на префикс, и форма, которая говорит, откуда брать операнды:

//! Подмножество SSE: пересылки, арифметика, логика, сравнения и
//! преобразования над регистрами `%xmm0` до `%xmm15`.
//!
//! Главная новость этого файла в том, что байты 0x66, 0xF3 и 0xF2 здесь
//! перестают быть префиксами и становятся частью опкода: они выбирают тип
//! данных. Один и тот же `0F 58` без префикса это `addps` (четыре f32),
//! с 0x66 `addpd` (два f64), с 0xF3 `addss` (один f32), с 0xF2 `addsd`
//! (один f64). Поэтому у каждого опкода не одно имя, а четыре.

const std = @import("std");

const cursor_mod = @import("cursor.zig");
const instruction = @import("instruction.zig");

const Cursor = cursor_mod.Cursor;
const Instruction = instruction.Instruction;
const Size = instruction.Size;

/// Четыре варианта опкода в порядке «без префикса, 0x66, 0xF3, 0xF2».
/// `null` значит, что такого сочетания нет.
const Names = [4]?[]const u8;

/// Откуда берутся операнды.
const Form = enum {
    /// `op xmm/m, %xmm`: результат в регистре из поля reg.
    load,
    /// `op %xmm, xmm/m`: запись в память или в регистр из поля r/m.
    store,
    /// `cvtsi2sd`: источник это целый регистр или память.
    from_integer,
    /// `cvttsd2si`: приёмник это целый регистр.
    to_integer,
    /// `movq %rdi, %xmm0`: биты целого регистра переезжают в xmm как есть.
    gpr_to_xmm,
    /// `movq %xmm0, %rax`: и обратно.
    xmm_to_gpr,
};

const Entry = struct { names: Names, form: Form };

/// Четыре имени по одному корню: `add` даёт `addps`, `addpd`, `addss`, `addsd`.
fn allFour(comptime root: []const u8) Names {
    return .{ root ++ "ps", root ++ "pd", root ++ "ss", root ++ "sd" };
}

/// Только упакованные варианты: логика и выровненные пересылки над
/// скалярами не бывают, для них хватает упакованной формы.
fn packedOnly(comptime root: []const u8) Names {
    return .{ root ++ "ps", root ++ "pd", null, null };
}

fn entry(opcode: u8) ?Entry {
    return switch (opcode) {
        0x10 => .{ .names = .{ "movups", "movupd", "movss", "movsd" }, .form = .load },
        0x11 => .{ .names = .{ "movups", "movupd", "movss", "movsd" }, .form = .store },
        0x14 => .{ .names = packedOnly("unpckl"), .form = .load },
        0x15 => .{ .names = packedOnly("unpckh"), .form = .load },
        0x28 => .{ .names = packedOnly("mova"), .form = .load },
        0x29 => .{ .names = packedOnly("mova"), .form = .store },
        0x2a => .{ .names = .{ null, null, "cvtsi2ss", "cvtsi2sd" }, .form = .from_integer },
        0x2c => .{ .names = .{ null, null, "cvttss2si", "cvttsd2si" }, .form = .to_integer },
        0x2d => .{ .names = .{ null, null, "cvtss2si", "cvtsd2si" }, .form = .to_integer },
        // Сравнение без префикса смотрит на f32, с 0x66 на f64. Флаги
        // выставляются как у беззнакового cmp, а NaN поднимает флаг чётности.
        0x2e => .{ .names = .{ "ucomiss", "ucomisd", null, null }, .form = .load },
        0x2f => .{ .names = .{ "comiss", "comisd", null, null }, .form = .load },
        0x51 => .{ .names = allFour("sqrt"), .form = .load },
        0x54 => .{ .names = packedOnly("and"), .form = .load },
        0x55 => .{ .names = packedOnly("andn"), .form = .load },
        0x56 => .{ .names = packedOnly("or"), .form = .load },
        0x57 => .{ .names = packedOnly("xor"), .form = .load },
        0x58 => .{ .names = allFour("add"), .form = .load },
        0x59 => .{ .names = allFour("mul"), .form = .load },
        0x5a => .{ .names = .{ "cvtps2pd", "cvtpd2ps", "cvtss2sd", "cvtsd2ss" }, .form = .load },
        0x5c => .{ .names = allFour("sub"), .form = .load },
        0x5d => .{ .names = allFour("min"), .form = .load },
        0x5e => .{ .names = allFour("div"), .form = .load },
        0x5f => .{ .names = allFour("max"), .form = .load },
        0x6e => .{ .names = .{ null, "movd", null, null }, .form = .gpr_to_xmm },
        0x7e => .{ .names = .{ null, "movd", null, null }, .form = .xmm_to_gpr },
        else => null,
    };
}

/// Номер варианта по обязательному префиксу. 0xF2 и 0xF3 сильнее 0x66:
/// если встретились оба, тип данных задаёт повторитель.
fn variant(cursor: Cursor) usize {
    if (cursor.prefixes.repne) return 3;
    if (cursor.prefixes.repe) return 2;
    if (cursor.prefixes.operand_size) return 1;
    return 0;
}

/// Опкоды SSE с ведущим 0x0F. `null`, если байт или сочетание с префиксом
/// не из нашего подмножества.
pub fn decode(cursor: *Cursor, opcode: u8) ?Instruction {
    const found = entry(opcode) orelse return null;
    const name = found.names[variant(cursor.*)] orelse return null;
    // Целый операнд у преобразований и у movd: REX.W делает его 64-битным.
    const integer: Size = if (cursor.prefixes.rexW()) .qword else .dword;

    return switch (found.form) {
        .load => {
            const modrm = cursor.readModRm(.xmm) orelse return null;
            return cursor.two(name, null, modrm.rm, xmm(cursor.*, modrm.reg));
        },
        .store => {
            const modrm = cursor.readModRm(.xmm) orelse return null;
            return cursor.two(name, null, xmm(cursor.*, modrm.reg), modrm.rm);
        },
        .from_integer => {
            const modrm = cursor.readModRm(integer) orelse return null;
            // По регистру ширина видна из имени, по памяти нет: тогда
            // LLVM дописывает суффикс, `cvtsi2sdq (%rdi), %xmm1`.
            const suffix: ?u8 = switch (modrm.rm) {
                .memory => integer.suffix(),
                else => null,
            };
            return cursor.two(name, suffix, modrm.rm, xmm(cursor.*, modrm.reg));
        },
        .to_integer => {
            const modrm = cursor.readModRm(.xmm) orelse return null;
            return cursor.two(name, null, modrm.rm, .{ .register = cursor.register(modrm.reg, integer) });
        },
        .gpr_to_xmm => {
            const modrm = cursor.readModRm(integer) orelse return null;
            return cursor.two(movName(integer), null, modrm.rm, xmm(cursor.*, modrm.reg));
        },
        .xmm_to_gpr => {
            const modrm = cursor.readModRm(integer) orelse return null;
            return cursor.two(movName(integer), null, xmm(cursor.*, modrm.reg), modrm.rm);
        },
    };
}

fn xmm(cursor: Cursor, index: u4) instruction.Operand {
    return .{ .register = cursor.register(index, .xmm) };
}

/// Четыре байта едут командой `movd`, восемь командой `movq`.
fn movName(size: Size) []const u8 {
    return if (size == .qword) "movq" else "movd";
}

test "префикс выбирает тип данных у одного и того же опкода" {
    const names = entry(0x58).?.names;
    try std.testing.expectEqualStrings("addps", names[0].?);
    try std.testing.expectEqualStrings("addpd", names[1].?);
    try std.testing.expectEqualStrings("addss", names[2].?);
    try std.testing.expectEqualStrings("addsd", names[3].?);
}

test "у логики нет скалярного варианта" {
    try std.testing.expectEqual(@as(?[]const u8, null), entry(0x57).?.names[3]);
}

Разберём, что здесь сделано и почему.

1. Имена склеиваются на этапе компиляции. allFour("add") превращается в четыре строки, как withPrefix из шага 12. Для логики есть packedOnly: xorps и andpd существуют, а xorsd нет. Урок показал почему: знак и модуль считаются маской по всему регистру, отдельная скалярная версия логике не нужна. Байты f2 0f 57 поэтому дают null и в листинге станут (bad), так же как у objdump.
2. У 0f 2e колонки сдвинуты. ucomiss кодируется без префикса, ucomisd с 66, а f2 и f3 ему не положены. Таблица на четыре колонки это выдерживает без особых случаев: в колонках скаляров просто null.
3. Порядок префиксов в variant. Встретив и 66, и f2, процессор берёт тип из повторителя. На практике компилятор такого не выдаёт, но правило в одну строку дешевле, чем разбираться, откуда в выводе взялся addpd на месте addsd.
4. Шесть форм, а не две. Большинство опкодов это load: результат в регистре из поля reg. Пересылки movsd и movaps бывают в обе стороны, и направление, как у целого mov, выбирает младший бит опкода: 10 читает, 11 пишет. Преобразования особенные, у них один операнд целый. У cvtsi2sd источник это %rdi или %eax, и его ширину задаёт REX.W, а приёмник %xmm. У cvttsd2si наоборот. Для этого в from_integer и to_integer одна сторона читается с шириной integer, другая с шириной .xmm.
5. Суффикс у cvtsi2sd только по памяти. По регистру ширина видна из имени: %eax это четыре байта, %rax восемь. По адресу (%rdi) её не видно, и LLVM дописывает суффикс: cvtsi2sdq (%rdi), %xmm1. Форматтер уже умеет печатать суффикс после мнемоники, нужно только его выставить.
6. movq между целым и %xmm. Опкоды 66 0f 6e и 66 0f 7e перекладывают биты как есть, без преобразования. Это тот самый @bitCast между f64 и u64 на уровне машины, и компилятор вставляет его всякий раз, когда ты смотришь на биты числа. Четыре байта это movd, восемь с REX.W это movq.

В src/x86/decoder.zig импорт const sse = @import("ops_sse.zig"); и две строки в двухбайтной таблице, перед 0xaf:

        // Числа с плавающей точкой: у этих опкодов префиксы 0x66, 0xF3 и 0xF2
        // выбирают тип данных, поэтому таблица своя, в `ops_sse.zig`.
        0x10, 0x11, 0x14, 0x15, 0x28, 0x29, 0x2a, 0x2c...0x2f => sse.decode(cursor, opcode),
        0x51, 0x54...0x5a, 0x5c...0x5f, 0x6e, 0x7e => sse.decode(cursor, opcode),

В src/root.zig в структуре x86 строка pub const ops_sse = @import("x86/ops_sse.zig");, чтобы встроенные тесты модуля попали в прогон.

Опкоды перечислены явно, диапазоном 0x10...0x7e их не заменить. В этих промежутках живут другие инструкции: 0x1f это многобайтный nop из шага 11, 0x40...0x4f это cmovcc из шага 12. Если отдать весь диапазон в ops_sse, он вернёт null на чужие опкоды, и cmovel превратится в (bad).

Последняя правка шага не в декодере, а в форматтере. Константа 3.14159 в ztScale грузится инструкцией f2 0f 59 05 00 00 00 00: адресация от %rip со смещением ноль, потому что адрес константы в .rodata впишет компоновщик. LLVM печатает такой адрес как mulsd (%rip), %xmm0, без нуля. Форматтер из урока 10 напечатал бы 0x0(%rip): нулевое смещение он опускает, только когда в адресе есть база или индекс. В прошлых фикстурах от %rip адресовался один movq 0x10(%rip), %rbx с ненулевым смещением, поэтому разница всплыла только сейчас:

--- a/src/x86/formatter.zig
+++ b/src/x86/formatter.zig
@@ -54,9 +54,11 @@
     // Замена сегмента идёт впереди всего адреса.
     if (memory.segment) |segment| try out.print("%{s}:", .{segment});
 
-    // Нулевое смещение не печатается, если адрес и без него полный.
+    // Нулевое смещение не печатается, если адрес и без него полный. До
+    // компоновки так выглядит любое обращение к данным: `(%rip)`, а на месте
+    // смещения нули, которые ждут своего перемещения.
     const show_displacement = memory.displacement != 0 or
-        (memory.base == null and memory.index == null);
+        (memory.base == null and memory.index == null and memory.rip_target == null);
     if (show_displacement) {
         if (memory.displacement < 0) {
             try out.print("-0x{x}", .{@abs(memory.displacement)});

Смещение без базы и индекса по-прежнему печатается всегда: 0x0 в movq %fs:0x0, %rax из урока 12 это настоящий абсолютный адрес, и выбросить его значит потерять операнд целиком.

Чего в шаге нет. Векторные инструкции AVX из второй половины урока, vaddpd и vmulps, кодируются совсем иначе: вместо префиксов и 0f у них двух- или трёхбайтный префикс VEX, который упаковывает в себя и REX, и тип данных, и третий регистр-операнд. Это отдельный декодер размером с половину нашего, и в курсе он не понадобится. Без флага -mcpu компилятор под x86_64 выдаёт только SSE2, поэтому весь код урока, собранный по умолчанию, zt теперь читает.

Фикстура и тесты шага

Фикстура это функции урока плюс одна функция на ассемблере с формами, которых в них не нашлось: адреса в памяти, регистры %xmm8 и выше, преобразования по памяти, movq и movd:

//! Фикстура шага 17: скалярное подмножество SSE.
//!
//! Собирается так:
//!   zig build-obj fixtures/step_17.zig -O ReleaseFast \
//!       -target x86_64-linux-musl -femit-bin=fixtures/step_17.o

/// Функции из урока: сложение, умножение, преобразования, константа
/// от `%rip`, знак и модуль маской, сравнения с флагом чётности.
export fn ztAddD(a: f64, b: f64) f64 {
    return a + b;
}

export fn ztPolyD(x: f64, a: f64, b: f64) f64 {
    return x * a + b;
}

export fn ztAddS(a: f32, b: f32) f32 {
    return a + b;
}

export fn ztIToD(n: i64) f64 {
    return @floatFromInt(n);
}

export fn ztDToI(x: f64) i64 {
    return @intFromFloat(x);
}

export fn ztUToD(n: u32) f64 {
    return @floatFromInt(n);
}

export fn ztSToD(x: f32) f64 {
    return @floatCast(x);
}

export fn ztDToS(x: f64) f32 {
    return @floatCast(x);
}

export fn ztScale(x: f64) f64 {
    return x * 3.14159;
}

export fn ztNegate(x: f64) f64 {
    return -x;
}

export fn ztAbsval(x: f64) f64 {
    return @abs(x);
}

export fn ztLt(a: f64, b: f64) bool {
    return a < b;
}

export fn ztIsNan(x: f64) bool {
    return x != x;
}

/// Сумма четырёх чисел одного вектора: без AVX вектор из четырёх f64
/// занимает два регистра xmm, и половинки складываются скалярно.
export fn ztHsum4(a: @Vector(4, f64)) f64 {
    return @reduce(.Add, a);
}

/// Формы, которых нет в функциях выше: память вместо регистра, старшие
/// регистры xmm8 до xmm15, пересылка между целым регистром и xmm.
export fn ztSseForms() callconv(.naked) void {
    asm volatile (
        \\movsd (%rdi), %xmm0
        \\movsd %xmm1, 0x8(%rsp)
        \\movsd %xmm1, %xmm0
        \\movss (%rdi), %xmm2
        \\movss %xmm2, (%rsi)
        \\movapd %xmm0, %xmm9
        \\movaps (%rax), %xmm1
        \\movaps %xmm15, 0x10(%rsp)
        \\movupd (%rdi,%rcx,8), %xmm3
        \\addsd 0x8(%rdi), %xmm0
        \\subsd %xmm1, %xmm0
        \\mulsd %xmm8, %xmm10
        \\divsd %xmm1, %xmm0
        \\sqrtsd %xmm0, %xmm0
        \\minsd %xmm1, %xmm0
        \\maxsd %xmm1, %xmm0
        \\subss %xmm1, %xmm0
        \\mulss %xmm1, %xmm0
        \\divss %xmm1, %xmm0
        \\addpd %xmm1, %xmm0
        \\mulps %xmm1, %xmm0
        \\xorps %xmm0, %xmm0
        \\xorpd %xmm1, %xmm1
        \\andpd %xmm1, %xmm0
        \\andnps %xmm1, %xmm0
        \\orps %xmm1, %xmm0
        \\ucomisd %xmm1, %xmm0
        \\ucomiss (%rdi), %xmm0
        \\comisd %xmm1, %xmm0
        \\cvtsi2sd %eax, %xmm0
        \\cvtsi2sdq (%rdi), %xmm1
        \\cvtsi2ss %rax, %xmm2
        \\cvttsd2si %xmm0, %eax
        \\cvttss2si %xmm0, %rax
        \\cvtsd2si %xmm1, %rcx
        \\cvtss2sd (%rdi), %xmm0
        \\unpcklpd %xmm1, %xmm0
        \\movq %rdi, %xmm0
        \\movq %xmm0, %rax
        \\movd %edi, %xmm1
        \\movd %xmm1, %eax
        \\ret
    );
}

zig build-obj fixtures/step_17.zig -O ReleaseFast \
    -target x86_64-linux-musl -femit-bin=fixtures/step_17.o
objdump -d fixtures/step_17.o > fixtures/step_17.objdump.txt

Запись в fixtures/root.zig:

pub const step_17 = Fixture{
    .object = @embedFile("step_17.o"),
    .objdump = @embedFile("step_17.objdump.txt"),
};

Тесты в tests/step_17.zig:

//! Шаг 17: подмножество SSE, скалярная арифметика с плавающей точкой.
//!
//! У этих инструкций байты 0x66, 0xF3 и 0xF2 это не префиксы, а часть
//! опкода: они выбирают тип данных. `0F 58` без префикса это `addps`,
//! с 0x66 `addpd`, с 0xF3 `addss`, с 0xF2 `addsd`. Регистры `%xmm`
//! кодируются в ModRM теми же четырьмя битами, что и целые.

const std = @import("std");

const fixtures = @import("fixtures");
const zt = @import("zt");

const support = @import("support.zig");

const decode = zt.x86.decoder.decode;

test "фикстура шага разбирается целиком" {
    try support.expectNoBadInstructions(fixtures.step_17);
}

test "вывод совпадает с objdump построчно" {
    try support.expectMatchesObjdump(fixtures.step_17);
}

test "предыдущие шаги остались зелёными" {
    try support.expectMatchesObjdump(fixtures.step_09);
    try support.expectMatchesObjdump(fixtures.step_10);
    try support.expectMatchesObjdump(fixtures.step_11);
    try support.expectMatchesObjdump(fixtures.step_12);
    try support.expectMatchesObjdump(fixtures.step_13);
    try support.expectMatchesObjdumpNamed(fixtures.step_14, fixtures.step_14_nm);
}

test "один опкод, четыре типа данных" {
    try expectOne(&.{ 0x0f, 0x58, 0xc1 }, "addps\t%xmm1, %xmm0");
    try expectOne(&.{ 0x66, 0x0f, 0x58, 0xc1 }, "addpd\t%xmm1, %xmm0");
    try expectOne(&.{ 0xf3, 0x0f, 0x58, 0xc1 }, "addss\t%xmm1, %xmm0");
    try expectOne(&.{ 0xf2, 0x0f, 0x58, 0xc1 }, "addsd\t%xmm1, %xmm0");
}

test "пересылка туда и обратно отличается младшим битом опкода" {
    try expectOne(&.{ 0xf2, 0x0f, 0x10, 0x07 }, "movsd\t(%rdi), %xmm0");
    try expectOne(&.{ 0xf2, 0x0f, 0x11, 0x4c, 0x24, 0x08 }, "movsd\t%xmm1, 0x8(%rsp)");
}

test "REX расширяет номера xmm так же, как номера целых регистров" {
    try expectOne(&.{ 0xf2, 0x45, 0x0f, 0x59, 0xd0 }, "mulsd\t%xmm8, %xmm10");
    try expectOne(&.{ 0x44, 0x0f, 0x29, 0x7c, 0x24, 0x10 }, "movaps\t%xmm15, 0x10(%rsp)");
}

test "константа грузится от rip" {
    // f2 0f 59 05 00 00 00 00: mulsd по адресу, который впишет компоновщик.
    try expectOne(&.{ 0xf2, 0x0f, 0x59, 0x05, 0x00, 0x00, 0x00, 0x00 }, "mulsd\t(%rip), %xmm0");
}

test "преобразования между целым регистром и xmm" {
    try expectOne(&.{ 0xf2, 0x48, 0x0f, 0x2a, 0xc7 }, "cvtsi2sd\t%rdi, %xmm0");
    try expectOne(&.{ 0xf2, 0x0f, 0x2a, 0xc0 }, "cvtsi2sd\t%eax, %xmm0");
    // По памяти ширину целого из имени не понять, и появляется суффикс.
    try expectOne(&.{ 0xf2, 0x48, 0x0f, 0x2a, 0x0f }, "cvtsi2sdq\t(%rdi), %xmm1");
    try expectOne(&.{ 0xf2, 0x48, 0x0f, 0x2c, 0xc0 }, "cvttsd2si\t%xmm0, %rax");
    try expectOne(&.{ 0xf2, 0x0f, 0x5a, 0xc0 }, "cvtsd2ss\t%xmm0, %xmm0");
    try expectOne(&.{ 0x66, 0x48, 0x0f, 0x6e, 0xc7 }, "movq\t%rdi, %xmm0");
    try expectOne(&.{ 0x66, 0x0f, 0x7e, 0xc8 }, "movd\t%xmm1, %eax");
}

test "сравнение: ucomisd ставит флаги, дальше обычные jcc и setcc" {
    try expectOne(&.{ 0x66, 0x0f, 0x2e, 0xc1 }, "ucomisd\t%xmm1, %xmm0");
    try expectOne(&.{ 0x0f, 0x9a, 0xc0 }, "setp\t%al");
}

test "у логики нет скалярного варианта, такие байты не наши" {
    try std.testing.expect(decode(&.{ 0xf2, 0x0f, 0x57, 0xc0 }, 0).bad);
}

/// Прогоняет байты через декодер и форматтер и сверяет мнемонику с операндами.
fn expectOne(code: []const u8, expected: []const u8) !void {
    const decoded = decode(code, 0);
    try std.testing.expectEqual(code.len, decoded.length());

    var buffer: [256]u8 = undefined;
    var out: std.Io.Writer = .fixed(&buffer);
    try zt.x86.formatter.format(&out, decoded);
    try std.testing.expectEqualStrings(expected, out.buffered());
}

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

И 17 в project_steps в build.zig.

Прогон

$ zig build test -Dstep=17 --summary all
Build Summary: 3/3 steps succeeded; 10/10 tests passed
test success
+- run test 10 pass (10 total) 16ms MaxRSS:3M

Вывод с подписями против objdump:

$ ./zig-out/bin/zt disasm --syms=step_17.nm fixtures/step_17.o > zt17.txt
$ objdump -d fixtures/step_17.o | tail -n +5 > od17.txt
$ diff od17.txt zt17.txt
5c5
<        a: f2 0f 10 c1                  	movsd	%xmm1, %xmm0            # xmm0 = xmm1[0],xmm0[1]
---
>        a: f2 0f 10 c1                  	movsd	%xmm1, %xmm0
39c39
<       96: 66 0f 14 c1                  	unpcklpd	%xmm1, %xmm0            # xmm0 = xmm0[0],xmm1[0]
---
>       96: 66 0f 14 c1                  	unpcklpd	%xmm1, %xmm0
55c55
<       ca: 66 0f 15 d1                  	unpckhpd	%xmm1, %xmm2            # xmm2 = xmm2[1],xmm1[1]
---
>       ca: 66 0f 15 d1                  	unpckhpd	%xmm1, %xmm2
58c58
<       d6: 66 0f 15 c0                  	unpckhpd	%xmm0, %xmm0            # xmm0 = xmm0[1,1]
---
>       d6: 66 0f 15 c0                  	unpckhpd	%xmm0, %xmm0
87c87
<      114: 0f 54 05 00 00 00 00         	andps	(%rip), %xmm0           # 0x11b <ztAbsval+0xb>
---
>      114: 0f 54 05 00 00 00 00         	andps	(%rip), %xmm0  # 0x11b
95c95
<      124: 0f 57 05 00 00 00 00         	xorps	(%rip), %xmm0           # 0x12b <ztNegate+0xb>
---
>      124: 0f 57 05 00 00 00 00         	xorps	(%rip), %xmm0  # 0x12b
103c103
<      134: f2 0f 59 05 00 00 00 00      	mulsd	(%rip), %xmm0           # 0x13c <ztScale+0xc>
---
>      134: f2 0f 59 05 00 00 00 00      	mulsd	(%rip), %xmm0  # 0x13c

Расходятся только комментарии, и обе их разновидности предсказуемы. У адресации от %rip objdump подписывает посчитанный адрес именем, а zt печатает голый адрес: подписи в комментариях появятся вместе с читателем ELF. А у movsd между регистрами и у unpck objdump рисует, какие слоты куда переехали: xmm0 = xmm1[0],xmm0[1] значит, что младший слот %xmm0 взят из %xmm1, а старший остался своим. Это не часть инструкции, а подсказка читателю, и из неё видно важное: скалярный movsd между регистрами пишет только младшие восемь байтов, а старшие оставляет приёмнику.

Теперь @reduce из урока можно прочитать своим инструментом:

$ ./zig-out/bin/zt disasm --syms=step_17.nm fixtures/step_17.o | sed -n '/<ztHsum4>/,/ e2:/p'
00000000000000b0 <ztHsum4>:
      b0: 55                           	pushq	%rbp
      b1: 48 89 e5                     	movq	%rsp, %rbp
      b4: 48 83 e4 e0                  	andq	$-0x20, %rsp
      b8: 48 83 ec 20                  	subq	$0x20, %rsp
      bc: 66 0f 28 4d 10               	movapd	0x10(%rbp), %xmm1
      c1: 66 0f 28 45 20               	movapd	0x20(%rbp), %xmm0
      c6: 66 0f 28 d1                  	movapd	%xmm1, %xmm2
      ca: 66 0f 15 d1                  	unpckhpd	%xmm1, %xmm2
      ce: f2 0f 58 d1                  	addsd	%xmm1, %xmm2
      d2: f2 0f 58 d0                  	addsd	%xmm0, %xmm2
      d6: 66 0f 15 c0                  	unpckhpd	%xmm0, %xmm0
      da: f2 0f 58 c2                  	addsd	%xmm2, %xmm0
      de: 48 89 ec                     	movq	%rbp, %rsp
      e1: 5d                           	popq	%rbp
      e2: c3                           	retq

Без AVX вектор из четырёх f64 не помещается в один регистр, поэтому он приезжает через стек и грузится двумя movapd в %xmm1 и %xmm0, по два числа в каждый. unpckhpd %xmm1, %xmm2 кладёт старшее число %xmm1 в младший слот %xmm2, addsd складывает с младшим. Ещё один unpckhpd и два addsd добавляют вторую пару. Порядок сложения виден по регистрам: ((a0 + a1) + a2) + a3, строго слева направо, как ты и ожидал бы от @reduce над f64 без разрешения переставлять слагаемые.

Практика

Соберём понятия сравнения и специальных значений в одну задачу. Функция findRange(x: f32) i32 определяет класс числа по IEEE-754: NaN, отрицательное, ноль, положительная субнормаль, нормальное, плюс бесконечность. Порядок проверок важен: NaN надо ловить первым, потому что любое сравнение с ним ложно, а минус бесконечность обязана попасть в отрицательные. Субнормаль это положительное число меньше наименьшего нормального, про денормали ты уже читал. Пиши минимумом сравнений, каждый класс проверь на границе.

Упражнения

mystery:
    movl        %edi, %eax
    cvtsi2ss    %rax, %xmm0
    ret

Итоги

• Float живут в отдельном файле регистров XMM (128 бит) и YMM (256 бит). Аргументы приходят в xmm0 и дальше, результат уходит в xmm0.
• Скалярные операции занимают младший слот регистра: addsd, subsd, mulsd, divsd для f64, те же с суффиксом ss для f32. Пересылка это movsd и movss.
• Преобразования это отдельные инструкции: cvtsi2sd из целого в f64, cvttsd2si из f64 в целое с усечением в сторону нуля, cvtss2sd и cvtsd2ss между ширинами. Беззнаковое u32 расширяется через movl перед переводом.
• Float-константа лежит в .rodata и грузится со смещением от rip, потому что у float-операций нет непосредственного операнда.
• Смена знака это xorps с маской 0x8000000000000000, модуль это andps с маской 0x7fffffffffffffff. Это перебор битов, а не арифметика.
• Сравнение ucomisd ставит флаги. Проверка на NaN идёт через флаг чётности PF: x != x компилируется в ucomisd плюс setp, а равенство в ветке защищают от NaN переходом jp.
• @Vector(N, T) кладёт N чисел в регистр как дорожки. Поэлементные операции под AVX2 это одна инструкция: vaddpd на четыре f64, vmulps на восемь f32. Дорожка i зависит только от дорожки i.
• @reduce(.Add, v) схлопывает вектор в число деревом из шаффлов и сложений, потому что процессор не складывает все дорожки одной инструкцией.
• Обычный цикл s += x над f64 не векторизуется: сложение float не ассоциативно, порядок меняет результат в младших битах, и компилятор держит строгий порядок. Чтобы получить дорожки, суммируй в @Vector-аккумулятор сам или разреши переупорядочивание через @setFloatMode(.optimized).
• В шаге проекта zt disasm научился подмножеству SSE: префиксы 66, f3 и f2 перед 0f выбирают тип данных, поэтому у каждого опкода четыре имени, а регистры %xmm кодируются в ModRM так же, как целые.

Дальше

Ты прочитал ассемблер float и увидел, что векторные регистры это не экзотика, а способ занять всю ширину машины. Дальше в разделе мы напишем ассемблер сами: не читать чужой, а класть свои байты в память и исполнять их. Начнём с inline-ассемблера в Zig, системного вызова без стандартной библиотеки и функции, целиком написанной на ассемблере в отдельном файле. А инструменты, которыми ты снимал листинги весь этот блок, скоро станут частью собственного дизассемблера: он уже умеет целочисленные инструкции, и SSE-подмножество из этого урока ложится в него следующим шагом.