SIMD e Vetorização em Zig: Processamento Paralelo de Dados

Escalar (1 elemento por vez):
  a[0] + b[0] = c[0]
  a[1] + b[1] = c[1]
  a[2] + b[2] = c[2]
  a[3] + b[3] = c[3]
  → 4 instrucoes

SIMD (4 elementos por vez):
  [a[0], a[1], a[2], a[3]] + [b[0], b[1], b[2], b[3]] = [c[0], c[1], c[2], c[3]]
  → 1 instrucao

const std = @import("std");

test "vetores SIMD basicos" {
    // Criar vetores de 4 floats
    const a: @Vector(4, f32) = .{ 1.0, 2.0, 3.0, 4.0 };
    const b: @Vector(4, f32) = .{ 5.0, 6.0, 7.0, 8.0 };

    // Operacoes aritmeticas — compilam para instrucoes SIMD
    const soma = a + b;     // [6, 8, 10, 12]
    const prod = a * b;     // [5, 12, 21, 32]
    const diff = b - a;     // [4, 4, 4, 4]

    try std.testing.expectEqual(@as(f32, 6.0), soma[0]);
    try std.testing.expectEqual(@as(f32, 12.0), prod[1]);
    try std.testing.expectEqual(@as(f32, 4.0), diff[2]);

    // Reducao: soma de todos os elementos
    const soma_total = @reduce(.Add, soma); // 6 + 8 + 10 + 12 = 36
    try std.testing.expectEqual(@as(f32, 36.0), soma_total);

    // Min e Max
    const minimo = @reduce(.Min, a); // 1.0
    const maximo = @reduce(.Max, b); // 8.0
    try std.testing.expectEqual(@as(f32, 1.0), minimo);
    try std.testing.expectEqual(@as(f32, 8.0), maximo);
}

/// Soma vetorizada de um array de floats
fn somaVetorizada(dados: []const f32) f32 {
    const vec_len = 8; // Processar 8 floats por vez (AVX2)
    var soma_vec: @Vector(vec_len, f32) = @splat(0.0);

    // Processar blocos de 8 elementos
    var i: usize = 0;
    while (i + vec_len <= dados.len) : (i += vec_len) {
        const bloco: @Vector(vec_len, f32) = dados[i..][0..vec_len].*;
        soma_vec += bloco;
    }

    // Reduzir o vetor para um escalar
    var resultado = @reduce(.Add, soma_vec);

    // Processar elementos restantes (tail)
    while (i < dados.len) : (i += 1) {
        resultado += dados[i];
    }

    return resultado;
}

test "soma vetorizada" {
    const dados = [_]f32{ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 };
    const resultado = somaVetorizada(&dados);
    try std.testing.expectApproxEqAbs(@as(f32, 55.0), resultado, 0.001);
}

fn dotProduct(a: []const f32, b: []const f32) f32 {
    std.debug.assert(a.len == b.len);

    const vec_len = 8;
    var soma_vec: @Vector(vec_len, f32) = @splat(0.0);

    var i: usize = 0;
    while (i + vec_len <= a.len) : (i += vec_len) {
        const va: @Vector(vec_len, f32) = a[i..][0..vec_len].*;
        const vb: @Vector(vec_len, f32) = b[i..][0..vec_len].*;
        soma_vec += va * vb; // Multiply-accumulate vetorizado
    }

    var resultado = @reduce(.Add, soma_vec);

    // Tail loop
    while (i < a.len) : (i += 1) {
        resultado += a[i] * b[i];
    }

    return resultado;
}

test "dot product SIMD" {
    const a = [_]f32{ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 };
    const b = [_]f32{ 10, 9, 8, 7, 6, 5, 4, 3, 2, 1 };
    // 1*10 + 2*9 + 3*8 + 4*7 + 5*6 + 6*5 + 7*4 + 8*3 + 9*2 + 10*1 = 220
    const resultado = dotProduct(&a, &b);
    try std.testing.expectApproxEqAbs(@as(f32, 220.0), resultado, 0.001);
}

/// Ajustar brilho de uma imagem RGBA (4 canais, 8-bit por canal)
fn ajustarBrilho(pixels: []u8, fator: f32) void {
    const vec_len = 16; // 16 bytes por vez (SSE)

    var i: usize = 0;
    while (i + vec_len <= pixels.len) : (i += vec_len) {
        // Carregar 16 bytes como vetor
        const dados: @Vector(vec_len, u8) = pixels[i..][0..vec_len].*;

        // Converter para float, multiplicar, converter de volta
        const float_dados: @Vector(vec_len, f32) = @floatFromInt(dados);
        const ajustado = float_dados * @as(@Vector(vec_len, f32), @splat(fator));

        // Clamp entre 0 e 255
        const zero: @Vector(vec_len, f32) = @splat(0.0);
        const max: @Vector(vec_len, f32) = @splat(255.0);
        const clamped = @min(@max(ajustado, zero), max);

        // Converter de volta para u8
        pixels[i..][0..vec_len].* = @intFromFloat(clamped);
    }

    // Tail loop para pixels restantes
    while (i < pixels.len) : (i += 1) {
        const val = @as(f32, @floatFromInt(pixels[i])) * fator;
        pixels[i] = @intFromFloat(@min(@max(val, 0.0), 255.0));
    }
}

test "ajustar brilho SIMD" {
    var pixels = [_]u8{ 100, 150, 200, 255 } ** 8;
    ajustarBrilho(&pixels, 1.5);

    try std.testing.expectEqual(@as(u8, 150), pixels[0]); // 100 * 1.5
    try std.testing.expectEqual(@as(u8, 225), pixels[1]); // 150 * 1.5
    try std.testing.expectEqual(@as(u8, 255), pixels[2]); // 200 * 1.5 = 300, clamped a 255
    try std.testing.expectEqual(@as(u8, 255), pixels[3]); // 255 * 1.5 = 382, clamped a 255
}

/// Calcular distancia entre N pares de pontos 3D
fn distanciasLote(
    ax: []const f32, ay: []const f32, az: []const f32,
    bx: []const f32, by: []const f32, bz: []const f32,
    resultado: []f32,
) void {
    const vec_len = 8;
    const n = ax.len;

    var i: usize = 0;
    while (i + vec_len <= n) : (i += vec_len) {
        const dx_vec: @Vector(vec_len, f32) = ax[i..][0..vec_len].* - @as(@Vector(vec_len, f32), bx[i..][0..vec_len].*);
        const dy_vec: @Vector(vec_len, f32) = ay[i..][0..vec_len].* - @as(@Vector(vec_len, f32), by[i..][0..vec_len].*);
        const dz_vec: @Vector(vec_len, f32) = az[i..][0..vec_len].* - @as(@Vector(vec_len, f32), bz[i..][0..vec_len].*);

        const dist_sq = dx_vec * dx_vec + dy_vec * dy_vec + dz_vec * dz_vec;

        // sqrt vetorizado
        resultado[i..][0..vec_len].* = @sqrt(dist_sq);
    }

    // Tail
    while (i < n) : (i += 1) {
        const dx = ax[i] - bx[i];
        const dy = ay[i] - by[i];
        const dz = az[i] - bz[i];
        resultado[i] = @sqrt(dx * dx + dy * dy + dz * dz);
    }
}