为何Rust编译器无法自动向量化该浮点点积实现？

先看这个简单的点积实现：

pub fn dot(a: &[f32], b: &[f32]) -> f32 {
    a.iter().zip(b.iter()).fold(0.0, |acc, (x, y)| acc + x * y)
}

用rustc 1.68加参数-C opt-level=3 -C target-feature=+avx2,+fma编译后，得到的是循环展开的标量汇编：

.LBB0_5:
        vmovss  xmm1, dword ptr [rdi + 4*rsi]
        vmulss  xmm1, xmm1, dword ptr [rdx + 4*rsi]
        vmovss  xmm2, dword ptr [rdi + 4*rsi + 4]
        vaddss  xmm0, xmm0, xmm1
        vmulss  xmm1, xmm2, dword ptr [rdx + 4*rsi + 4]
        vaddss  xmm0, xmm0, xmm1
        vmovss  xmm1, dword ptr [rdi + 4*rsi + 8]
        vmulss  xmm1, xmm1, dword ptr [rdx + 4*rsi + 8]
        vaddss  xmm0, xmm0, xmm1
        vmovss  xmm1, dword ptr [rdi + 4*rsi + 12]
        vmulss  xmm1, xmm1, dword ptr [rdx + 4*rsi + 12]
        lea     rax, [rsi + 4]
        vaddss  xmm0, xmm0, xmm1
        mov     rsi, rax
        cmp     rcx, rax
        jne     .LBB0_5

这里不仅没用到AVX2的向量寄存器，连FMA指令都没合并，完全是标量运算。但把f32换成i32后，就能得到预期的向量化汇编：

.LBB0_5:
        vmovdqu ymm4, ymmword ptr [rdx + 4*rax]
        vmovdqu ymm5, ymmword ptr [rdx + 4*rax + 32]
        vmovdqu ymm6, ymmword ptr [rdx + 4*rax + 64]
        vmovdqu ymm7, ymmword ptr [rdx + 4*rax + 96]
        vpmulld ymm4, ymm4, ymmword ptr [rdi + 4*rax]
        vpaddd  ymm0, ymm4, ymm0
        vpmulld ymm4, ymm5, ymmword ptr [rdi + 4*rax + 32]
        vpaddd  ymm1, ymm4, ymm1
        vpmulld ymm4, ymm6, ymmword ptr [rdi + 4*rax + 64]
        vpmulld ymm5, ymm7, ymmword ptr [rdi + 4*rax + 96]
        vpaddd  ymm2, ymm4, ymm2
        vpaddd  ymm3, ymm5, ymm3
        add     rax, 32
        cmp     r8, rax
        jne     .LBB0_5

核心原因：浮点数运算的结合律限制 #

整数加法严格满足结合律，不管运算顺序怎么调整，最终结果都完全一致。但浮点数加法受限于IEEE 754标准的精度要求，运算顺序会影响最终结果——比如(a + b) + c和a + (b + c)可能因为舍入误差产生细微差异。

LLVM默认会严格遵循IEEE 754的语义，不会通过改变运算顺序来做向量化优化。而点积的fold操作是按顺序累积相加：acc = acc + x1*y1; acc = acc + x2*y2; ...，如果改成SIMD向量化，会先计算多个x*y的向量，再把向量内的元素相加，最后再和累积值相加，这相当于改变了加法的顺序，LLVM不敢默认这么做。

而整数运算没有这个精度顾虑，所以LLVM可以安全地将循环向量化，用AVX2的向量寄存器一次处理8个i32元素，大幅提升效率。

怎么让f32点积也能自动向量化？ #

如果你可以接受浮点数运算顺序改变带来的微小精度差异，可以通过以下方式开启优化：

1. 启用fast-math优化
   编译时添加参数-C fast-math，或者更精细地控制：-C llvm-args=-enable-fp-contract=fast。fast-math会允许LLVM忽略IEEE 754的部分语义限制，包括重新排序浮点运算，这样LLVM就会生成向量化和FMA指令。

2. 手动使用SIMD intrinsics
   直接用core::arch::x86_64里的AVX2/FMA指令手动实现，比如：

use core::arch::x86_64::*;

#[target_feature(enable = "avx2,fma")]
pub unsafe fn dot_simd(a: &[f32], b: &[f32]) -> f32 {
    assert_eq!(a.len(), b.len());
    let mut sum = _mm256_setzero_ps();
    let mut i = 0;
    let n = a.len() & !7; // 对齐到8个元素的倍数
    while i < n {
        let va = _mm256_loadu_ps(a.as_ptr().add(i));
        let vb = _mm256_loadu_ps(b.as_ptr().add(i));
        sum = _mm256_fmadd_ps(va, vb, sum);
        i += 8;
    }
    // 归约向量中的8个元素
    let mut tmp = _mm256_hadd_ps(sum, sum);
    tmp = _mm256_hadd_ps(tmp, tmp);
    let res = _mm_cvtss_f32(_mm256_castps256_ps128(tmp));
    // 处理剩余元素
    let mut res_scalar = res;
    while i < a.len() {
        res_scalar += a[i] * b[i];
        i += 1;
    }
    res_scalar
}

3. 使用并行迭代库
   比如rayon的par_iter，并行归约会自动处理向量化，同时允许重新排序运算顺序：

use rayon::prelude::*;

pub fn dot_par(a: &[f32], b: &[f32]) -> f32 {
    a.par_iter().zip(b.par_iter())
        .map(|(x, y)| x * y)
        .sum()
}

这种方式既简洁又能获得向量化的性能提升，同时sum()的并行实现会处理好归约的精度问题（虽然还是和标量fold有细微差异，但大多数场景可以接受）。

内容的提问来源于stack exchange，提问作者benjamin-lieser