float dot(const float* a, const float* b, size_t n) {
    __m256 acc0 = _mm256_setzero_ps();
    __m256 acc1 = _mm256_setzero_ps();
    size_t i = 0;
    for (; i + 16 <= n; i += 16) {
        acc0 = _mm256_fmadd_ps(_mm256_loadu_ps(a + i),
                               _mm256_loadu_ps(b + i), acc0);
        acc1 = _mm256_fmadd_ps(_mm256_loadu_ps(a + i + 8),
                               _mm256_loadu_ps(b + i + 8), acc1);
    }
    __m128 v = _mm_add_ps(_mm256_castps256_ps128(acc0),
                          _mm256_extractf128_ps(acc1, 1));
    v = _mm_add_ps(v, _mm_movehl_ps(v, v));
    v = _mm_add_ss(v, _mm_shuffle_ps(v, v, 0x1));
    float s = _mm_cvtss_f32(v);
    for (; i < n; ++i) s += a[i] * b[i];
    return s;
}

dot(float const*, float const*, unsigned long):
        vxorps        xmm0, xmm0, xmm0
        vxorps        xmm1, xmm1, xmm1
        xor           eax, eax
.Lloop:
        vmovups       ymm2, ymmword ptr [rdi + 4*rax]
        vmovups       ymm3, ymmword ptr [rdi + 4*rax + 32]
        vfmadd231ps   ymm0, ymm2, ymmword ptr [rsi + 4*rax]
        vfmadd231ps   ymm1, ymm3, ymmword ptr [rsi + 4*rax + 32]
        add           rax, 16
        lea           rcx, [rax + 16]
        cmp           rcx, rdx
        jbe           .Lloop
        vaddps        ymm0, ymm0, ymm1
        vextractf128  xmm1, ymm0, 1
        vaddps        xmm0, xmm0, xmm1
        vmovhlps      xmm1, xmm0, xmm0
        vaddss        xmm0, xmm0, xmm1
        vzeroupper
        ret