## Coverage Report

#### Created: 2026-01-25 15:05

next uncovered line (L), next uncovered region (R), next uncovered branch (B)

    /home/noah/src/realizar/src/inference/simd.rs

    Line
    Count
    Source
    1
    //! SIMD-accelerated operations for inference
    2
    //!
    3
    //! Provides high-performance primitive operations using trueno's SIMD backend.
    4
    //! All operations are designed for cache efficiency with tiled implementations.
    5
    //!
    6
    //! ## Operations
    7
    //!
    8
    //! - [`simd_matmul`] - Matrix-vector multiplication with SIMD dot products
    9
    //! - [`simd_dot`] - SIMD-accelerated dot product
    10
    //! - [`simd_add`] - Vector addition
    11
    //! - [`simd_mul`] - Element-wise multiplication
    12
    //! - [`simd_silu`] - SiLU activation (x * sigmoid(x))
    13
    //! - [`simd_gelu`] - GELU activation (approximate)
    14
    //! - [`simd_softmax`] - Numerically stable softmax
    15
    //!
    16
    //! ## Performance
    17
    //!
    18
    //! Uses trueno's Vector::dot for all dot products, enabling:
    19
    //! - AVX2/SSE on x86
    20
    //! - NEON on ARM
    21
    //! - WASM SIMD in browsers
    22
    //! - Scalar fallback everywhere else
    23
    24
    use trueno::Vector;
    25
    26
    /// Tile size for cache-efficient tiled matmul
    27
    const TILE_SIZE: usize = 64;
    28
    29
    /// SIMD-accelerated matrix-vector multiplication
    30
    ///
    31
    /// Uses trueno's optimized SIMD backend for maximum performance.
    32
    /// Falls back to scalar for non-SIMD architectures.
    33
    ///
    34
    /// # Arguments
    35
    ///
    36
    /// * `input` - Input vector of length `in_dim`
    37
    /// * `weight` - Weight matrix stored row-major [out_dim × in_dim]
    38
    /// * `in_dim` - Input dimension
    39
    /// * `out_dim` - Output dimension
    40
    ///
    41
    /// # Returns
    42
    ///
    43
    /// Output vector of length `out_dim`
    44
    ///
    45
    /// # Example
    46
    ///
    47
    /// ```
    48
    /// use realizar::inference::simd_matmul;
    49
    ///
    50
    /// // 2x3 matrix times 3-vector = 2-vector
    51
    /// let input = vec![1.0, 2.0, 3.0];
    52
    /// let weight = vec![
    53
    ///     1.0, 0.0, 0.0,  // row 0: extracts x
    54
    ///     0.0, 1.0, 0.0,  // row 1: extracts y
    55
    /// ];
    56
    /// let output = simd_matmul(&input, &weight, 3, 2);
    57
    /// assert_eq!(output.len(), 2);
    58
    /// ```
    59
    #[must_use]
    60
    9
    pub fn simd_matmul(input: &[f32], weight: &[f32], in_dim: usize, out_dim: usize) -> Vec<f32> {
    61
        // Convert to trueno types for SIMD acceleration
    62
    9
        let input_vec = Vector::from_slice(input);
    63
    64
        // Compute each output element using SIMD dot product
    65
    9
        let mut output = vec![0.0; out_dim];
    66
    67
        // Use tiled approach for better cache utilization
    68
    11
        for tile_start in (0..out_dim)9.step_by9(TILE_SIZE) {
    69
    11
            let tile_end = (tile_start + TILE_SIZE).min(out_dim);
    70
    71
    275
            for row in tile_start11..tile_end11 {
    72
    275
                let row_start = row * in_dim;
    73
    275
                let row_end = row_start + in_dim;
    74
    275
                let row_vec = Vector::from_slice(&weight[row_start..row_end]);
    75
    275
                output[row] = input_vec.dot(&row_vec).expect("dot product failed");
    76
    275
            }
    77
        }
    78
    79
    9
        output
    80
    9
    }
    81
    82
    /// SIMD-accelerated dot product
    83
    ///
    84
    /// Uses trueno's SIMD backend for vectorized computation.
    85
    ///
    86
    /// # Example
    87
    ///
    88
    /// ```
    89
    /// use realizar::inference::simd_dot;
    90
    ///
    91
    /// let a = vec![1.0, 2.0, 3.0];
    92
    /// let b = vec![4.0, 5.0, 6.0];
    93
    /// let result = simd_dot(&a, &b);
    94
    /// assert!((result - 32.0).abs() < 1e-5);
    95
    /// ```
    96
    #[inline]
    97
    #[must_use]
    98
    239
    pub fn simd_dot(a: &[f32], b: &[f32]) -> f32 {
    99
    239
        Vector::from_slice(a)
    100
    239
            .dot(&Vector::from_slice(b))
    101
    239
            .expect("dot product failed")
    102
    239
    }
    103
    104
    /// SIMD-accelerated vector addition (a += b)
    105
    ///
    106
    /// # Example
    107
    ///
    108
    /// ```
    109
    /// use realizar::inference::simd_add;
    110
    ///
    111
    /// let mut a = vec![1.0, 2.0, 3.0];
    112
    /// let b = vec![4.0, 5.0, 6.0];
    113
    /// simd_add(&mut a, &b);
    114
    /// assert_eq!(a, vec![5.0, 7.0, 9.0]);
    115
    /// ```
    116
    #[inline]
    117
    5
    pub fn simd_add(a: &mut [f32], b: &[f32]) {
    118
    15
        for (x, y) in a5.iter_mut5().zip5(b5.iter5()) {
    119
    15
            *x += y;
    120
    15
        }
    121
    5
    }
    122
    123
    /// SIMD-accelerated element-wise multiplication (a *= b)
    124
    ///
    125
    /// # Example
    126
    ///
    127
    /// ```
    128
    /// use realizar::inference::simd_mul;
    129
    ///
    130
    /// let mut a = vec![1.0, 2.0, 3.0];
    131
    /// let b = vec![4.0, 5.0, 6.0];
    132
    /// simd_mul(&mut a, &b);
    133
    /// assert_eq!(a, vec![4.0, 10.0, 18.0]);
    134
    /// ```
    135
    #[inline]
    136
    6
    pub fn simd_mul(a: &mut [f32], b: &[f32]) {
    137
    17
        for (x, y) in a6.iter_mut6().zip6(b6.iter6()) {
    138
    17
            *x *= y;
    139
    17
        }
    140
    6
    }
    141
    142
    /// SIMD-accelerated SiLU activation (x * sigmoid(x))
    143
    ///
    144
    /// Also known as Swish activation: f(x) = x / (1 + exp(-x))
    145
    ///
    146
    /// # Example
    147
    ///
    148
    /// ```
    149
    /// use realizar::inference::simd_silu;
    150
    ///
    151
    /// let mut data = vec![0.0, 1.0, -1.0];
    152
    /// simd_silu(&mut data);
    153
    /// assert!((data[0] - 0.0).abs() < 1e-5);  // silu(0) = 0
    154
    /// assert!((data[1] - 0.7311).abs() < 0.01);  // silu(1) ≈ 0.731
    155
    /// ```
    156
    #[inline]
    157
    8
    pub fn simd_silu(data: &mut [f32]) {
    158
    16
        for x in data8.iter_mut8() {
    159
    16
            *x = *x / (1.0 + (-*x).exp());
    160
    16
        }
    161
    8
    }
    162
    163
    /// SIMD-accelerated GELU activation (approximate)
    164
    ///
    165
    /// Uses the tanh approximation:
    166
    /// GELU(x) ≈ 0.5 * x * (1 + tanh(sqrt(2/π) * (x + 0.044715 * x³)))
    167
    ///
    168
    /// # Example
    169
    ///
    170
    /// ```
    171
    /// use realizar::inference::simd_gelu;
    172
    ///
    173
    /// let mut data = vec![0.0, 1.0, -1.0];
    174
    /// simd_gelu(&mut data);
    175
    /// assert!((data[0] - 0.0).abs() < 1e-5);  // gelu(0) = 0
    176
    /// assert!((data[1] - 0.8413).abs() < 0.01);  // gelu(1) ≈ 0.841
    177
    /// ```
    178
    #[inline]
    179
    9
    pub fn simd_gelu(data: &mut [f32]) {
    180
        // Approximate GELU: 0.5 * x * (1 + tanh(sqrt(2/pi) * (x + 0.044715 * x^3)))
    181
        const SQRT_2_OVER_PI: f32 = 0.797_884_6; // sqrt(2/π)
    182
        const COEF: f32 = 0.044715;
    183
    184
    12
        for x in data9.iter_mut9() {
    185
    12
            let x3 = *x * *x * *x;
    186
    12
            let inner = SQRT_2_OVER_PI * (*x + COEF * x3);
    187
    12
            *x = 0.5 * *x * (1.0 + inner.tanh());
    188
    12
        }
    189
    9
    }
    190
    191
    /// SIMD-accelerated softmax with numerical stability
    192
    ///
    193
    /// Uses the max-subtraction trick to prevent overflow:
    194
    /// softmax(x)_i = exp(x_i - max(x)) / sum(exp(x_j - max(x)))
    195
    ///
    196
    /// # Example
    197
    ///
    198
    /// ```
    199
    /// use realizar::inference::simd_softmax;
    200
    ///
    201
    /// let mut data = vec![1.0, 2.0, 3.0];
    202
    /// simd_softmax(&mut data);
    203
    ///
    204
    /// // Probabilities should sum to 1
    205
    /// let sum: f32 = data.iter().sum();
    206
    /// assert!((sum - 1.0).abs() < 1e-5);
    207
    ///
    208
    /// // Largest input should have largest probability
    209
    /// assert!(data[2] > data[1]);
    210
    /// assert!(data[1] > data[0]);
    211
    /// ```
    212
    105
    pub fn simd_softmax(data: &mut [f32]) {
    213
    105
        if data.is_empty() {
    214
    2
            return;
    215
    103
        }
    216
    217
        // Find max for numerical stability
    218
    103
        let max_val = data.iter().copied().fold(f32::NEG_INFINITY, f32::max);
    219
    220
        // Compute exp(x - max) and sum
    221
    103
        let mut sum = 0.0;
    222
    253
        for x in data103.iter_mut103() {
    223
    253
            *x = (*x - max_val).exp();
    224
    253
            sum += *x;
    225
    253
        }
    226
    227
        // Normalize
    228
    103
        if sum > 0.0 {
    229
    103
            let inv_sum = 1.0 / sum;
    230
    253
            for x in data103.iter_mut103() {
    231
    253
                *x *= inv_sum;
    232
    253
            }
    233
    0
        }
    234
    105
    }
    235
    236
    // ============================================================================
    237
    // BF16/F16 SIMD Conversion (T-QA-021 Optimization)
    238
    // ============================================================================
    239
    240
    /// Fast BF16→F32 conversion using bit manipulation
    241
    ///
    242
    /// BF16 is a truncated F32 (same exponent, fewer mantissa bits).
    243
    /// Conversion is just a 16-bit left shift.
    244
    ///
    245
    /// # Arguments
    246
    ///
    247
    /// * `input` - Raw BF16 bytes (2 bytes per value)
    248
    ///
    249
    /// # Returns
    250
    ///
    251
    /// F32 vector with converted values
    252
    ///
    253
    /// # Performance
    254
    ///
    255
    /// This implementation uses SIMD on x86_64 with AVX2 support,
    256
    /// processing 8 BF16 values in parallel.
    257
    ///
    258
    /// # Example
    259
    ///
    260
    /// ```
    261
    /// use realizar::inference::simd_bf16_to_f32;
    262
    ///
    263
    /// let bf16_bytes = half::bf16::from_f32(1.5).to_le_bytes();
    264
    /// let f32_vals = simd_bf16_to_f32(&bf16_bytes);
    265
    /// assert!((f32_vals[0] - 1.5).abs() < 0.01);
    266
    /// ```
    267
    #[must_use]
    268
    18
    pub fn simd_bf16_to_f32(input: &[u8]) -> Vec<f32> {
    269
    18
        let count = input.len() / 2;
    270
    18
        if count == 0 {
    271
    1
            return Vec::new();
    272
    17
        }
    273
    274
        #[cfg(all(target_arch = "x86_64", target_feature = "avx2"))]
    275
        {
    276
            simd_bf16_to_f32_avx2(input, count)
    277
        }
    278
    279
        #[cfg(not(all(target_arch = "x86_64", target_feature = "avx2")))]
    280
        {
    281
    17
            bf16_to_f32_fast(input, count)
    282
        }
    283
    18
    }
    284
    285
    /// AVX2-accelerated BF16→F32 conversion
    286
    #[cfg(all(target_arch = "x86_64", target_feature = "avx2"))]
    287
    fn simd_bf16_to_f32_avx2(input: &[u8], count: usize) -> Vec<f32> {
    288
        use std::arch::x86_64::*;
    289
    290
        let mut output = vec![0.0f32; count];
    291
        let chunks = count / 8;
    292
        let remainder = count % 8;
    293
    294
        // SAFETY: AVX2 target_feature is required by cfg, input bounds checked by chunks calculation,
    295
        // output vector pre-allocated to count elements
    296
        unsafe {
    297
            for i in 0..chunks {
    298
                let in_offset = i * 16;
    299
                let out_offset = i * 8;
    300
    301
                // Load 8 BF16 values (16 bytes)
    302
                let bf16_bytes = _mm_loadu_si128(input.as_ptr().add(in_offset) as *const __m128i);
    303
    304
                // Unpack lower 4 BF16 to F32 (zero-extend and shift left by 16)
    305
                let lo = _mm_unpacklo_epi16(bf16_bytes, _mm_setzero_si128());
    306
                let lo_shifted = _mm_slli_epi32(lo, 16);
    307
    308
                // Unpack upper 4 BF16 to F32
    309
                let hi = _mm_unpackhi_epi16(bf16_bytes, _mm_setzero_si128());
    310
                let hi_shifted = _mm_slli_epi32(hi, 16);
    311
    312
                // Store results
    313
                _mm_storeu_ps(
    314
                    output.as_mut_ptr().add(out_offset),
    315
                    _mm_castsi128_ps(lo_shifted),
    316
                );
    317
                _mm_storeu_ps(
    318
                    output.as_mut_ptr().add(out_offset + 4),
    319
                    _mm_castsi128_ps(hi_shifted),
    320
                );
    321
            }
    322
        }
    323
    324
        // Handle remainder with scalar
    325
        let remainder_start = chunks * 8;
    326
        for i in 0..remainder {
    327
            let offset = (remainder_start + i) * 2;
    328
            let bits = u16::from_le_bytes([input[offset], input[offset + 1]]) as u32;
    329
            output[remainder_start + i] = f32::from_bits(bits << 16);
    330
        }
    331
    332
        output
    333
    }
    334
    335
    /// Fast scalar BF16→F32 conversion using bit manipulation
    336
    #[cfg(not(all(target_arch = "x86_64", target_feature = "avx2")))]
    337
    17
    fn bf16_to_f32_fast(input: &[u8], count: usize) -> Vec<f32> {
    338
    17
        let mut output = Vec::with_capacity(count);
    339
    567
        for chunk in input17.chunks_exact17(2) {
    340
    567
            let bits = u16::from_le_bytes([chunk[0], chunk[1]]) as u32;
    341
    567
            output.push(f32::from_bits(bits << 16));
    342
    567
        }
    343
    17
        output
    344
    17
    }
    345
    346
    /// Scalar fallback for non-AVX2 platforms
    347
    #[cfg(all(target_arch = "x86_64", target_feature = "avx2"))]
    348
    fn bf16_to_f32_fast(input: &[u8], count: usize) -> Vec<f32> {
    349
        let mut output = Vec::with_capacity(count);
    350
        for chunk in input.chunks_exact(2) {
    351
            let bits = u16::from_le_bytes([chunk[0], chunk[1]]) as u32;
    352
            output.push(f32::from_bits(bits << 16));
    353
        }
    354
        output
    355
    }
    356
    357
    /// Fast F16→F32 conversion using the half crate
    358
    ///
    359
    /// Unlike BF16, F16 has a different exponent bias and requires
    360
    /// proper conversion (not just bit shifting).
    361
    ///
    362
    /// # Arguments
    363
    ///
    364
    /// * `input` - Raw F16 bytes (2 bytes per value)
    365
    ///
    366
    /// # Returns
    367
    ///
    368
    /// F32 vector with converted values
    369
    #[must_use]
    370
    2
    pub fn simd_f16_to_f32(input: &[u8]) -> Vec<f32> {
    371
    2
        input
    372
    2
            .chunks_exact(2)
    373
    7
            .map2(|chunk| {
    374
    7
                let bits = u16::from_le_bytes([chunk[0], chunk[1]]);
    375
    7
                half::f16::from_bits(bits).to_f32()
    376
    7
            })
    377
    2
            .collect()
    378
    2
    }
    379
    380
    /// SIMD-accelerated BF16 dot product
    381
    ///
    382
    /// Computes dot product of two BF16 vectors without full conversion.
    383
    /// Converts small chunks at a time to keep F32 data in L1 cache.
    384
    ///
    385
    /// # Arguments
    386
    ///
    387
    /// * `a` - First BF16 vector (raw bytes)
    388
    /// * `b` - Second BF16 vector (raw bytes)
    389
    ///
    390
    /// # Returns
    391
    ///
    392
    /// Dot product as F32
    393
    #[must_use]
    394
    2
    pub fn simd_bf16_dot(a: &[u8], b: &[u8]) -> f32 {
    395
        const CHUNK_SIZE: usize = 64; // 64 BF16 values = 128 bytes, fits in L1
    396
    397
    2
        let count = a.len().min(b.len()) / 2;
    398
    2
        let mut sum = 0.0f32;
    399
    400
    5
        for chunk_start in (0..count)2.step_by2(CHUNK_SIZE) {
    401
    5
            let chunk_end = (chunk_start + CHUNK_SIZE).min(count);
    402
    5
            let byte_start = chunk_start * 2;
    403
    5
            let byte_end = chunk_end * 2;
    404
    5
    405
    5
            // Convert chunk to F32
    406
    5
            let a_f32 = simd_bf16_to_f32(&a[byte_start..byte_end]);
    407
    5
            let b_f32 = simd_bf16_to_f32(&b[byte_start..byte_end]);
    408
    5
    409
    5
            // Compute dot product of chunk using SIMD
    410
    5
            sum += simd_dot(&a_f32, &b_f32);
    411
    5
        }
    412
    413
    2
        sum
    414
    2
    }
    415
    416
    /// SIMD-accelerated BF16 matmul
    417
    ///
    418
    /// Computes matrix-vector product with BF16 weights.
    419
    /// Uses batch conversion to minimize conversion overhead.
    420
    ///
    421
    /// # Arguments
    422
    ///
    423
    /// * `input` - F32 input vector
    424
    /// * `weight_bf16` - BF16 weight matrix (raw bytes, row-major)
    425
    /// * `in_dim` - Input dimension
    426
    /// * `out_dim` - Output dimension
    427
    ///
    428
    /// # Returns
    429
    ///
    430
    /// F32 output vector
    431
    ///
    432
    /// # Performance
    433
    ///
    434
    /// This function batch-converts BF16 rows to F32 in tiles to amortize
    435
    /// conversion overhead and improve cache utilization. The tile size is
    436
    /// chosen to fit in L2 cache (~256KB per tile).
    437
    #[must_use]
    438
    2
    pub fn simd_bf16_matmul(
    439
    2
        input: &[f32],
    440
    2
        weight_bf16: &[u8],
    441
    2
        in_dim: usize,
    442
    2
        out_dim: usize,
    443
    2
    ) -> Vec<f32> {
    444
        // Batch conversion: convert all BF16 weights to F32 once
    445
        // This is more efficient than row-by-row conversion because:
    446
        // 1. Amortizes function call overhead
    447
        // 2. Better SIMD utilization (longer vectors)
    448
        // 3. Memory prefetching works better
    449
        //
    450
        // Trade-off: Uses 2x memory (BF16 + F32), but much faster
    451
    2
        let weight_f32 = simd_bf16_to_f32(weight_bf16);
    452
    453
        // Now use optimized F32 matmul
    454
    2
        simd_matmul(input, &weight_f32, in_dim, out_dim)
    455
    2
    }
    456
    457
    /// SIMD-accelerated BF16 matmul with streaming conversion
    458
    ///
    459
    /// This variant uses row-by-row conversion for lower memory usage
    460
    /// at the cost of performance. Use for very large matrices that
    461
    /// don't fit in memory when fully converted.
    462
    ///
    463
    /// # Arguments
    464
    ///
    465
    /// * `input` - F32 input vector
    466
    /// * `weight_bf16` - BF16 weight matrix (raw bytes, row-major)
    467
    /// * `in_dim` - Input dimension
    468
    /// * `out_dim` - Output dimension
    469
    ///
    470
    /// # Returns
    471
    ///
    472
    /// F32 output vector
    473
    #[must_use]
    474
    0
    pub fn simd_bf16_matmul_streaming(
    475
    0
        input: &[f32],
    476
    0
        weight_bf16: &[u8],
    477
    0
        in_dim: usize,
    478
    0
        out_dim: usize,
    479
    0
    ) -> Vec<f32> {
    480
        const TILE_SIZE: usize = 64;
    481
    482
    0
        let mut output = vec![0.0f32; out_dim];
    483
    0
        let input_vec = Vector::from_slice(input);
    484
    485
    0
        for tile_start in (0..out_dim).step_by(TILE_SIZE) {
    486
    0
            let tile_end = (tile_start + TILE_SIZE).min(out_dim);
    487
    488
    0
            for row in tile_start..tile_end {
    489
    0
                let row_byte_start = row * in_dim * 2;
    490
    0
                let row_byte_end = row_byte_start + in_dim * 2;
    491
    0
                let row_bf16 = &weight_bf16[row_byte_start..row_byte_end];
    492
    0
    493
    0
                // Convert row to F32
    494
    0
                let row_f32 = simd_bf16_to_f32(row_bf16);
    495
    0
                let row_vec = Vector::from_slice(&row_f32);
    496
    0
    497
    0
                output[row] = input_vec.dot(&row_vec).expect("dot product failed");
    498
    0
            }
    499
        }
    500
    501
    0
        output
    502
    0
    }
    503
    504
    // ============================================================================
    505
    // EXTREME TDD: Comprehensive Tests
    506
    // ============================================================================
    507
    508
    #[cfg(test)]
    509
    mod tests {
    510
        use super::*;
    511
    512
        // ------------------------------------------------------------------------
    513
        // simd_matmul Tests
    514
        // ------------------------------------------------------------------------
    515
    516
        #[test]
    517
    1
        fn test_simd_matmul_identity() {
    518
            // 3x3 identity matrix
    519
    1
            let input = vec![1.0, 2.0, 3.0];
    520
    1
            let identity = vec![
    521
                1.0, 0.0, 0.0, // row 0
    522
                0.0, 1.0, 0.0, // row 1
    523
                0.0, 0.0, 1.0, // row 2
    524
            ];
    525
    1
            let output = simd_matmul(&input, &identity, 3, 3);
    526
    1
            assert_eq!(output, vec![1.0, 2.0, 3.0]);
    527
    1
        }
    528
    529
        #[test]
    530
    1
        fn test_simd_matmul_projection() {
    531
            // 2x3 projection matrix
    532
    1
            let input = vec![1.0, 2.0, 3.0];
    533
    1
            let weight = vec![
    534
                1.0, 1.0, 1.0, // row 0: sum
    535
                1.0, 0.0, -1.0, // row 1: x - z
    536
            ];
    537
    1
            let output = simd_matmul(&input, &weight, 3, 2);
    538
    1
            assert_eq!(output.len(), 2);
    539
    1
            assert!((output[0] - 6.0).abs() < 1e-5); // 1+2+3 = 6
    540
    1
            assert!((output[1] - (-2.0)).abs() < 1e-5); // 1-3 = -2
    541
    1
        }
    542
    543
        #[test]
    544
    1
        fn test_simd_matmul_expansion() {
    545
            // 4x2 expansion matrix
    546
    1
            let input = vec![1.0, 2.0];
    547
    1
            let weight = vec![
    548
                1.0, 0.0, // row 0: x
    549
                0.0, 1.0, // row 1: y
    550
                1.0, 1.0, // row 2: x+y
    551
                1.0, -1.0, // row 3: x-y
    552
            ];
    553
    1
            let output = simd_matmul(&input, &weight, 2, 4);
    554
    1
            assert_eq!(output.len(), 4);
    555
    1
            assert!((output[0] - 1.0).abs() < 1e-5);
    556
    1
            assert!((output[1] - 2.0).abs() < 1e-5);
    557
    1
            assert!((output[2] - 3.0).abs() < 1e-5);
    558
    1
            assert!((output[3] - (-1.0)).abs() < 1e-5);
    559
    1
        }
    560
    561
        #[test]
    562
    1
        fn test_simd_matmul_large_tiled() {
    563
            // Test that tiling works for large matrices
    564
    1
            let in_dim = 128;
    565
    1
            let out_dim = 256;
    566
    128
            let input1: Vec<f32>1 = (0..in_dim)1.map1(|i| i as f32).collect1();
    567
    568
            // Create a simple weight matrix (diagonal-ish)
    569
    1
            let mut weight = vec![0.0; out_dim * in_dim];
    570
    128
            for i in 0..out_dim1.min1(in_dim1) {
    571
    128
                weight[i * in_dim + i] = 1.0;
    572
    128
            }
    573
    574
    1
            let output = simd_matmul(&input, &weight, in_dim, out_dim);
    575
    1
            assert_eq!(output.len(), out_dim);
    576
    577
            // First `in_dim` outputs should equal inputs
    578
    128
            for i in 0..in_dim1 {
    579
    128
                assert!((output[i] - i as f32).abs() < 1e-5);
    580
            }
    581
            // Remaining outputs should be zero
    582
    128
            for i in in_dim1..out_dim1 {
    583
    128
                assert!((output[i]).abs() < 1e-5);
    584
            }
    585
    1
        }
    586
    587
        #[test]
    588
    1
        fn test_simd_matmul_empty() {
    589
    1
            let input: Vec<f32> = vec![];
    590
    1
            let weight: Vec<f32> = vec![];
    591
    1
            let output = simd_matmul(&input, &weight, 0, 0);
    592
    1
            assert!(output.is_empty());
    593
    1
        }
    594
    595
        // ------------------------------------------------------------------------
    596
        // simd_dot Tests
    597
        // ------------------------------------------------------------------------
    598
    599
        #[test]
    600
    1
        fn test_simd_dot_basic() {
    601
    1
            let a = vec![1.0, 2.0, 3.0];
    602
    1
            let b = vec![4.0, 5.0, 6.0];
    603
    1
            let result = simd_dot(&a, &b);
    604
    1
            assert!((result - 32.0).abs() < 1e-5); // 1*4 + 2*5 + 3*6 = 32
    605
    1
        }
    606
    607
        #[test]
    608
    1
        fn test_simd_dot_orthogonal() {
    609
    1
            let a = vec![1.0, 0.0];
    610
    1
            let b = vec![0.0, 1.0];
    611
    1
            let result = simd_dot(&a, &b);
    612
    1
            assert!((result).abs() < 1e-5);
    613
    1
        }
    614
    615
        #[test]
    616
    1
        fn test_simd_dot_self() {
    617
    1
            let a = vec![3.0, 4.0];
    618
    1
            let result = simd_dot(&a, &a);
    619
    1
            assert!((result - 25.0).abs() < 1e-5); // 3^2 + 4^2 = 25
    620
    1
        }
    621
    622
        #[test]
    623
    1
        fn test_simd_dot_negative() {
    624
    1
            let a = vec![1.0, -1.0];
    625
    1
            let b = vec![-1.0, 1.0];
    626
    1
            let result = simd_dot(&a, &b);
    627
    1
            assert!((result - (-2.0)).abs() < 1e-5);
    628
    1
        }
    629
    630
        #[test]
    631
    1
        fn test_simd_dot_large() {
    632
    1
            let n = 1024;
    633
    1
            let a: Vec<f32> = vec![1.0; n];
    634
    1
            let b: Vec<f32> = vec![1.0; n];
    635
    1
            let result = simd_dot(&a, &b);
    636
    1
            assert!((result - n as f32).abs() < 1e-3);
    637
    1
        }
    638
    639
        // ------------------------------------------------------------------------
    640
        // simd_add Tests
    641
        // ------------------------------------------------------------------------
    642
    643
        #[test]
    644
    1
        fn test_simd_add_basic() {
    645
    1
            let mut a = vec![1.0, 2.0, 3.0];
    646
    1
            let b = vec![4.0, 5.0, 6.0];
    647
    1
            simd_add(&mut a, &b);
    648
    1
            assert_eq!(a, vec![5.0, 7.0, 9.0]);
    649
    1
        }
    650
    651
        #[test]
    652
    1
        fn test_simd_add_zeros() {
    653
    1
            let mut a = vec![1.0, 2.0, 3.0];
    654
    1
            let b = vec![0.0, 0.0, 0.0];
    655
    1
            simd_add(&mut a, &b);
    656
    1
            assert_eq!(a, vec![1.0, 2.0, 3.0]);
    657
    1
        }
    658
    659
        #[test]
    660
    1
        fn test_simd_add_negative() {
    661
    1
            let mut a = vec![1.0, 2.0, 3.0];
    662
    1
            let b = vec![-1.0, -2.0, -3.0];
    663
    1
            simd_add(&mut a, &b);
    664
    1
            assert_eq!(a, vec![0.0, 0.0, 0.0]);
    665
    1
        }
    666
    667
        // ------------------------------------------------------------------------
    668
        // simd_mul Tests
    669
        // ------------------------------------------------------------------------
    670
    671
        #[test]
    672
    1
        fn test_simd_mul_basic() {
    673
    1
            let mut a = vec![1.0, 2.0, 3.0];
    674
    1
            let b = vec![4.0, 5.0, 6.0];
    675
    1
            simd_mul(&mut a, &b);
    676
    1
            assert_eq!(a, vec![4.0, 10.0, 18.0]);
    677
    1
        }
    678
    679
        #[test]
    680
    1
        fn test_simd_mul_ones() {
    681
    1
            let mut a = vec![1.0, 2.0, 3.0];
    682
    1
            let b = vec![1.0, 1.0, 1.0];
    683
    1
            simd_mul(&mut a, &b);
    684
    1
            assert_eq!(a, vec![1.0, 2.0, 3.0]);
    685
    1
        }
    686
    687
        #[test]
    688
    1
        fn test_simd_mul_zeros() {
    689
    1
            let mut a = vec![1.0, 2.0, 3.0];
    690
    1
            let b = vec![0.0, 0.0, 0.0];
    691
    1
            simd_mul(&mut a, &b);
    692
    1
            assert_eq!(a, vec![0.0, 0.0, 0.0]);
    693
    1
        }
    694
    695
        #[test]
    696
    1
        fn test_simd_mul_negative() {
    697
    1
            let mut a = vec![2.0, 3.0];
    698
    1
            let b = vec![-1.0, -2.0];
    699
    1
            simd_mul(&mut a, &b);
    700
    1
            assert_eq!(a, vec![-2.0, -6.0]);
    701
    1
        }
    702
    703
        // ------------------------------------------------------------------------
    704
        // simd_silu Tests
    705
        // ------------------------------------------------------------------------
    706
    707
        #[test]
    708
    1
        fn test_simd_silu_zero() {
    709
    1
            let mut data = vec![0.0];
    710
    1
            simd_silu(&mut data);
    711
    1
            assert!((data[0]).abs() < 1e-5); // silu(0) = 0
    712
    1
        }
    713
    714
        #[test]
    715
    1
        fn test_simd_silu_positive() {
    716
    1
            let mut data = vec![1.0];
    717
    1
            simd_silu(&mut data);
    718
            // silu(1) = 1 / (1 + exp(-1)) ≈ 0.7311
    719
    1
            assert!((data[0] - 0.7311).abs() < 0.01);
    720
    1
        }
    721
    722
        #[test]
    723
    1
        fn test_simd_silu_negative() {
    724
    1
            let mut data = vec![-1.0];
    725
    1
            simd_silu(&mut data);
    726
            // silu(-1) = -1 / (1 + exp(1)) ≈ -0.2689
    727
    1
            assert!((data[0] - (-0.2689)).abs() < 0.01);
    728
    1
        }
    729
    730
        #[test]
    731
    1
        fn test_simd_silu_large_positive() {
    732
    1
            let mut data = vec![10.0];
    733
    1
            simd_silu(&mut data);
    734
            // silu(10) ≈ 10 (sigmoid(10) ≈ 1)
    735
    1
            assert!((data[0] - 10.0).abs() < 0.01);
    736
    1
        }
    737
    738
        #[test]
    739
    1
        fn test_simd_silu_large_negative() {
    740
    1
            let mut data = vec![-10.0];
    741
    1
            simd_silu(&mut data);
    742
            // silu(-10) ≈ 0 (sigmoid(-10) ≈ 0)
    743
    1
            assert!((data[0]).abs() < 0.01);
    744
    1
        }
    745
    746
        #[test]
    747
    1
        fn test_simd_silu_batch() {
    748
    1
            let mut data = vec![0.0, 1.0, -1.0, 2.0, -2.0];
    749
    1
            simd_silu(&mut data);
    750
    1
            assert!((data[0]).abs() < 1e-5);
    751
    1
            assert!(data[1] > 0.0);
    752
    1
            assert!(data[2] < 0.0);
    753
    1
            assert!(data[3] > data[1]); // monotonic for x > 0
    754
    1
        }
    755
    756
        // ------------------------------------------------------------------------
    757
        // simd_gelu Tests
    758
        // ------------------------------------------------------------------------
    759
    760
        #[test]
    761
    1
        fn test_simd_gelu_zero() {
    762
    1
            let mut data = vec![0.0];
    763
    1
            simd_gelu(&mut data);
    764
    1
            assert!((data[0]).abs() < 1e-5); // gelu(0) = 0
    765
    1
        }
    766
    767
        #[test]
    768
    1
        fn test_simd_gelu_positive() {
    769
    1
            let mut data = vec![1.0];
    770
    1
            simd_gelu(&mut data);
    771
            // gelu(1) ≈ 0.841
    772
    1
            assert!((data[0] - 0.841).abs() < 0.01);
    773
    1
        }
    774
    775
        #[test]
    776
    1
        fn test_simd_gelu_negative() {
    777
    1
            let mut data = vec![-1.0];
    778
    1
            simd_gelu(&mut data);
    779
            // gelu(-1) ≈ -0.159
    780
    1
            assert!((data[0] - (-0.159)).abs() < 0.01);
    781
    1
        }
    782
    783
        #[test]
    784
    1
        fn test_simd_gelu_large_positive() {
    785
    1
            let mut data = vec![3.0];
    786
    1
            simd_gelu(&mut data);
    787
            // gelu(3) ≈ 3 (tanh approaches 1)
    788
    1
            assert!((data[0] - 3.0).abs() < 0.01);
    789
    1
        }
    790
    791
        #[test]
    792
    1
        fn test_simd_gelu_large_negative() {
    793
    1
            let mut data = vec![-3.0];
    794
    1
            simd_gelu(&mut data);
    795
            // gelu(-3) ≈ 0 (tanh approaches -1)
    796
    1
            assert!((data[0]).abs() < 0.01);
    797
    1
        }
    798
    799
        #[test]
    800
    1
        fn test_simd_gelu_symmetry_breaking() {
    801
            // GELU is NOT symmetric: gelu(-x) != -gelu(x)
    802
    1
            let mut pos = vec![1.0];
    803
    1
            let mut neg = vec![-1.0];
    804
    1
            simd_gelu(&mut pos);
    805
    1
            simd_gelu(&mut neg);
    806
    1
            assert!((pos[0] + neg[0]).abs() > 0.1); // sum should not be 0
    807
    1
        }
    808
    809
        // ------------------------------------------------------------------------
    810
        // simd_softmax Tests
    811
        // ------------------------------------------------------------------------
    812
    813
        #[test]
    814
    1
        fn test_simd_softmax_sums_to_one() {
    815
    1
            let mut data = vec![1.0, 2.0, 3.0];
    816
    1
            simd_softmax(&mut data);
    817
    1
            let sum: f32 = data.iter().sum();
    818
    1
            assert!((sum - 1.0).abs() < 1e-5);
    819
    1
        }
    820
    821
        #[test]
    822
    1
        fn test_simd_softmax_preserves_order() {
    823
    1
            let mut data = vec![1.0, 2.0, 3.0];
    824
    1
            simd_softmax(&mut data);
    825
    1
            assert!(data[2] > data[1]);
    826
    1
            assert!(data[1] > data[0]);
    827
    1
        }
    828
    829
        #[test]
    830
    1
        fn test_simd_softmax_uniform() {
    831
    1
            let mut data = vec![1.0, 1.0, 1.0];
    832
    1
            simd_softmax(&mut data);
    833
            // Should be uniform distribution
    834
    4
            for &x3 in &data {
    835
    3
                assert!((x - 1.0 / 3.0).abs() < 1e-5);
    836
            }
    837
    1
        }
    838
    839
        #[test]
    840
    1
        fn test_simd_softmax_empty() {
    841
    1
            let mut data: Vec<f32> = vec![];
    842
    1
            simd_softmax(&mut data);
    843
    1
            assert!(data.is_empty());
    844
    1
        }
    845
    846
        #[test]
    847
    1
        fn test_simd_softmax_single() {
    848
    1
            let mut data = vec![5.0];
    849
    1
            simd_softmax(&mut data);
    850
    1
            assert!((data[0] - 1.0).abs() < 1e-5);
    851
    1
        }
    852
    853
        #[test]
    854
    1
        fn test_simd_softmax_numerical_stability() {
    855
            // Large values that would overflow without max subtraction
    856
    1
            let mut data = vec![1000.0, 1001.0, 1002.0];
    857
    1
            simd_softmax(&mut data);
    858
    1
            let sum: f32 = data.iter().sum();
    859
    1
            assert!((sum - 1.0).abs() < 1e-5);
    860
    3
            assert!1(data.iter()1.all1(|&x| x.is_finite()));
    861
    1
        }
    862
    863
        #[test]
    864
    1
        fn test_simd_softmax_negative() {
    865
    1
            let mut data = vec![-1.0, -2.0, -3.0];
    866
    1
            simd_softmax(&mut data);
    867
    1
            let sum: f32 = data.iter().sum();
    868
    1
            assert!((sum - 1.0).abs() < 1e-5);
    869
            // Order reversed: -1 > -2 > -3
    870
    1
            assert!(data[0] > data[1]);
    871
    1
            assert!(data[1] > data[2]);
    872
    1
        }
    873
    874
        #[test]
    875
    1
        fn test_simd_softmax_temperature_effect() {
    876
            // Larger differences should give more peaked distribution
    877
    1
            let mut narrow = vec![1.0, 2.0, 3.0];
    878
    1
            let mut wide = vec![1.0, 10.0, 100.0];
    879
    880
    1
            simd_softmax(&mut narrow);
    881
    1
            simd_softmax(&mut wide);
    882
    883
            // Wide should be more peaked (largest value dominates)
    884
    1
            assert!(wide[2] > narrow[2]);
    885
    1
        }
    886
    887
        // ------------------------------------------------------------------------
    888
        // Integration Tests
    889
        // ------------------------------------------------------------------------
    890
    891
        #[test]
    892
    1
        fn test_matmul_then_activation() {
    893
    1
            let input = vec![1.0, 2.0];
    894
    1
            let weight = vec![
    895
                1.0, 1.0, // sum: 3
    896
                -1.0, 1.0, // diff: 1
    897
            ];
    898
    1
            let mut output = simd_matmul(&input, &weight, 2, 2);
    899
    1
            assert!((output[0] - 3.0).abs() < 1e-5);
    900
    1
            assert!((output[1] - 1.0).abs() < 1e-5);
    901
    902
    1
            simd_gelu(&mut output);
    903
            // gelu(3) ≈ 3, gelu(1) ≈ 0.841
    904
    1
            assert!((output[0] - 3.0).abs() < 0.01);
    905
    1
            assert!((output[1] - 0.841).abs() < 0.01);
    906
    1
        }
    907
    908
        #[test]
    909
    1
        fn test_residual_connection() {
    910
    1
            let input = vec![1.0, 2.0, 3.0];
    911
    1
            let weight = vec![
    912
                0.1, 0.0, 0.0, 0.0, 0.1, 0.0, 0.0, 0.0, 0.1, // 0.1 * I
    913
            ];
    914
    1
            let proj = simd_matmul(&input, &weight, 3, 3);
    915
    916
    1
            let mut residual = input.clone();
    917
    1
            simd_add(&mut residual, &proj);
    918
    919
            // residual = input + 0.1 * input = 1.1 * input
    920
    1
            assert!((residual[0] - 1.1).abs() < 1e-5);
    921
    1
            assert!((residual[1] - 2.2).abs() < 1e-5);
    922
    1
            assert!((residual[2] - 3.3).abs() < 1e-5);
    923
    1
        }
    924
    925
        #[test]
    926
    1
        fn test_gated_activation() {
    927
            // SwiGLU style: gate * up
    928
    1
            let mut gate = vec![0.0, 1.0, 2.0];
    929
    1
            let up = vec![1.0, 2.0, 3.0];
    930
    931
    1
            simd_silu(&mut gate);
    932
    1
            simd_mul(&mut gate, &up);
    933
    934
            // gate[0] = silu(0) * 1 = 0
    935
    1
            assert!((gate[0]).abs() < 1e-5);
    936
            // gate[1] = silu(1) * 2 ≈ 0.7311 * 2 ≈ 1.46
    937
    1
            assert!((gate[1] - 1.46).abs() < 0.05);
    938
    1
        }
    939
    940
        // ------------------------------------------------------------------------
    941
        // BF16/F16 Conversion Tests (T-QA-021)
    942
        // ------------------------------------------------------------------------
    943
    944
        #[test]
    945
    1
        fn test_simd_bf16_to_f32_empty() {
    946
    1
            let result = simd_bf16_to_f32(&[]);
    947
    1
            assert!(result.is_empty());
    948
    1
        }
    949
    950
        #[test]
    951
    1
        fn test_simd_bf16_to_f32_single() {
    952
            // BF16 representation of 1.0: 0x3F80
    953
    1
            let bf16_bytes = half::bf16::from_f32(1.0).to_le_bytes();
    954
    1
            let result = simd_bf16_to_f32(&bf16_bytes);
    955
    1
            assert_eq!(result.len(), 1);
    956
    1
            assert!((result[0] - 1.0).abs() < 1e-6);
    957
    1
        }
    958
    959
        #[test]
    960
    1
        fn test_simd_bf16_to_f32_various_values() {
    961
    1
            let values = [0.0f32, 1.0, -1.0, 0.5, 2.0, -0.5, 100.0, -100.0];
    962
    1
            let mut bf16_bytes = Vec::new();
    963
    9
            for &v8 in &values {
    964
    8
                bf16_bytes.extend_from_slice(&half::bf16::from_f32(v).to_le_bytes());
    965
    8
            }
    966
    967
    1
            let result = simd_bf16_to_f32(&bf16_bytes);
    968
    1
            assert_eq!(result.len(), values.len());
    969
    970
    8
            for (i, (&expected, &actual)) in values1.iter1().zip1(result.iter()1).enumerate1() {
    971
                // BF16 has limited precision, allow some tolerance
    972
    8
                let tol = expected.abs().max(1.0) * 0.01;
    973
    8
                assert!(
    974
    8
                    (actual - expected).abs() < tol,
    975
    0
                    "Value {} mismatch: expected {}, got {}",
    976
                    i,
    977
                    expected,
    978
                    actual
    979
                );
    980
            }
    981
    1
        }
    982
    983
        #[test]
    984
    1
        fn test_simd_bf16_to_f32_large_batch() {
    985
            // Test with more than 8 values to verify SIMD remainder handling
    986
    1
            let count = 17; // 2 SIMD chunks (8+8) + 1 remainder
    987
    1
            let mut bf16_bytes = Vec::with_capacity(count * 2);
    988
    17
            for i in 0..count1 {
    989
    17
                let v = i as f32 * 0.1;
    990
    17
                bf16_bytes.extend_from_slice(&half::bf16::from_f32(v).to_le_bytes());
    991
    17
            }
    992
    993
    1
            let result = simd_bf16_to_f32(&bf16_bytes);
    994
    1
            assert_eq!(result.len(), count);
    995
    996
    17
            for i in 0..count1 {
    997
    17
                let expected = i as f32 * 0.1;
    998
    17
                let tol = 0.01;
    999
    17
                assert!(
    1000
    17
                    (result[i] - expected).abs() < tol,
    1001
    0
                    "Index {} mismatch: expected {}, got {}",
    1002
                    i,
    1003
                    expected,
    1004
    0
                    result[i]
    1005
                );
    1006
            }
    1007
    1
        }
    1008
    1009
        #[test]
    1010
    1
        fn test_simd_f16_to_f32_single() {
    1011
    1
            let f16_bytes = half::f16::from_f32(1.0).to_le_bytes();
    1012
    1
            let result = simd_f16_to_f32(&f16_bytes);
    1013
    1
            assert_eq!(result.len(), 1);
    1014
    1
            assert!((result[0] - 1.0).abs() < 1e-3);
    1015
    1
        }
    1016
    1017
        #[test]
    1018
    1
        fn test_simd_f16_to_f32_various_values() {
    1019
    1
            let values = [0.0f32, 1.0, -1.0, 0.5, 2.0, -0.5];
    1020
    1
            let mut f16_bytes = Vec::new();
    1021
    7
            for &v6 in &values {
    1022
    6
                f16_bytes.extend_from_slice(&half::f16::from_f32(v).to_le_bytes());
    1023
    6
            }
    1024
    1025
    1
            let result = simd_f16_to_f32(&f16_bytes);
    1026
    1
            assert_eq!(result.len(), values.len());
    1027
    1028
    6
            for (expected, actual) in values1.iter1().zip1(result.iter()1) {
    1029
                // F16 has limited precision
    1030
    6
                assert!(
    1031
    6
                    (actual - expected).abs() < 0.01,
    1032
    0
                    "Mismatch: expected {}, got {}",
    1033
                    expected,
    1034
                    actual
    1035
                );
    1036
            }
    1037
    1
        }
    1038
    1039
        #[test]
    1040
    1
        fn test_simd_bf16_dot_basic() {
    1041
            // Create two simple BF16 vectors
    1042
    1
            let a_vals = [1.0f32, 2.0, 3.0, 4.0];
    1043
    1
            let b_vals = [1.0f32, 1.0, 1.0, 1.0];
    1044
    1045
    1
            let mut a_bytes = Vec::new();
    1046
    1
            let mut b_bytes = Vec::new();
    1047
    4
            for (&a, &b) in a_vals1.iter1().zip1(b_vals1.iter1()) {
    1048
    4
                a_bytes.extend_from_slice(&half::bf16::from_f32(a).to_le_bytes());
    1049
    4
                b_bytes.extend_from_slice(&half::bf16::from_f32(b).to_le_bytes());
    1050
    4
            }
    1051
    1052
    1
            let result = simd_bf16_dot(&a_bytes, &b_bytes);
    1053
            // 1+2+3+4 = 10
    1054
    1
            assert!((result - 10.0).abs() < 0.1);
    1055
    1
        }
    1056
    1057
        #[test]
    1058
    1
        fn test_simd_bf16_dot_large() {
    1059
            // Test with many chunks to exercise chunked processing
    1060
    1
            let n = 256; // 4 chunks of 64
    1061
    1
            let mut a_bytes = Vec::with_capacity(n * 2);
    1062
    1
            let mut b_bytes = Vec::with_capacity(n * 2);
    1063
    1064
    256
            for i in 0..n1 {
    1065
    256
                let v = ((i % 10) as f32) * 0.1;
    1066
    256
                a_bytes.extend_from_slice(&half::bf16::from_f32(v).to_le_bytes());
    1067
    256
                b_bytes.extend_from_slice(&half::bf16::from_f32(1.0).to_le_bytes());
    1068
    256
            }
    1069
    1070
    1
            let result = simd_bf16_dot(&a_bytes, &b_bytes);
    1071
            // Sum of 0.0, 0.1, 0.2, ..., 0.9, 0.0, 0.1, ... (26 complete cycles)
    1072
            // Each cycle sums to 0+0.1+0.2+...+0.9 = 4.5
    1073
            // 256/10 = 25 full cycles + 6 remainder (0.0+0.1+0.2+0.3+0.4+0.5 = 1.5)
    1074
    1
            let expected = 25.0 * 4.5 + 1.5;
    1075
    1
            assert!(
    1076
    1
                (result - expected).abs() < 1.0,
    1077
    0
                "Large dot: expected ~{}, got {}",
    1078
                expected,
    1079
                result
    1080
            );
    1081
    1
        }
    1082
    1083
        #[test]
    1084
    1
        fn test_simd_bf16_matmul_identity() {
    1085
            // 3x3 identity in BF16
    1086
    1
            let input = vec![1.0f32, 2.0, 3.0];
    1087
    1
            let mut weight_bytes = Vec::new();
    1088
    1
            let identity = [[1.0f32, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]];
    1089
    4
            for row3 in &identity {
    1090
    12
                for &v9 in row {
    1091
    9
                    weight_bytes.extend_from_slice(&half::bf16::from_f32(v).to_le_bytes());
    1092
    9
                }
    1093
            }
    1094
    1095
    1
            let output = simd_bf16_matmul(&input, &weight_bytes, 3, 3);
    1096
    1
            assert_eq!(output.len(), 3);
    1097
    1
            assert!((output[0] - 1.0).abs() < 0.01);
    1098
    1
            assert!((output[1] - 2.0).abs() < 0.01);
    1099
    1
            assert!((output[2] - 3.0).abs() < 0.01);
    1100
    1
        }
    1101
    1102
        #[test]
    1103
    1
        fn test_simd_bf16_matmul_projection() {
    1104
            // 2x4 projection matrix
    1105
    1
            let input = vec![1.0f32, 2.0, 3.0, 4.0];
    1106
    1
            let weight = [
    1107
    1
                [1.0f32, 1.0, 1.0, 1.0], // row 0: sum = 10
    1108
    1
                [1.0, 0.0, 0.0, -1.0],   // row 1: 1-4 = -3
    1109
    1
            ];
    1110
    1
            let mut weight_bytes = Vec::new();
    1111
    3
            for row2 in &weight {
    1112
    10
                for &v8 in row {
    1113
    8
                    weight_bytes.extend_from_slice(&half::bf16::from_f32(v).to_le_bytes());
    1114
    8
                }
    1115
            }
    1116
    1117
    1
            let output = simd_bf16_matmul(&input, &weight_bytes, 4, 2);
    1118
    1
            assert_eq!(output.len(), 2);
    1119
    1
            assert!(
    1120
    1
                (output[0] - 10.0).abs() < 0.1,
    1121
    0
                "Sum: expected 10, got {}",
    1122
    0
                output[0]
    1123
            );
    1124
    1
            assert!(
    1125
    1
                (output[1] - (-3.0)).abs() < 0.1,
    1126
    0
                "Diff: expected -3, got {}",
    1127
    0
                output[1]
    1128
            );
    1129
    1
        }
    1130
    }