realizar 0.3.2

//! IMP-600b: GPU vs SIMD Benchmark for Matrix-Vector Multiplication
//!
//! Popperian Falsification Test:
//! CLAIM: "trueno wgpu GPU CANNOT match llama.cpp cuBLAS performance for matvec"
//!
//! KEY INSIGHT: For LLM token generation, n=1, making matmul a MATVEC operation.
//! Research shows cuBLAS "hurts more than helps" for matvec due to latency overhead.
//!
//! Run: cargo run --example gpu_matvec_benchmark --features gpu

use std::hint::black_box;
use std::time::Instant;

// Scalar implementation for baseline
fn matvec_scalar(matrix: &[f32], vector: &[f32], rows: usize, cols: usize) -> Vec<f32> {
    let mut result = vec![0.0f32; rows];
    for i in 0..rows {
        let mut sum = 0.0f32;
        for j in 0..cols {
            sum += matrix[i * cols + j] * vector[j];
        }
        result[i] = sum;
    }
    result
}

// SIMD implementation using trueno
fn matvec_simd(matrix: &[f32], vector: &[f32], rows: usize, cols: usize) -> Vec<f32> {
    use trueno::Vector as TruenoVector;

    let v = TruenoVector::from_slice(vector);
    let mut result = vec![0.0f32; rows];

    for (i, res) in result.iter_mut().enumerate().take(rows) {
        let row_start = i * cols;
        let row_end = row_start + cols;
        let row = TruenoVector::from_slice(&matrix[row_start..row_end]);
        *res = row.dot(&v).unwrap_or(0.0);
    }
    result
}

fn main() {
    println!("=== IMP-600b: GPU vs SIMD Matvec Benchmark ===");
    println!("Popperian Falsification: Can trueno GPU achieve parity?\n");

    // Typical LLM dimensions
    let test_cases = [
        ("phi2 hidden", 2560, 2560),
        ("phi2 ffn", 10240, 2560),
        ("llama 7B hidden", 4096, 4096),
        ("llama 7B ffn", 11008, 4096),
        ("llama 70B hidden", 8192, 8192),
    ];

    println!(
        "{:<20} {:>10} {:>10} {:>10} {:>12} {:>12}",
        "Operation", "Rows", "Cols", "Scalar(µs)", "SIMD(µs)", "Speedup"
    );
    println!("{}", "-".repeat(80));

    for (name, rows, cols) in test_cases {
        // Create test data
        let matrix: Vec<f32> = (0..(rows * cols)).map(|i| (i as f32) * 0.0001).collect();
        let vector: Vec<f32> = (0..cols).map(|i| (i as f32) * 0.001).collect();

        let iterations = 100;

        // Warmup
        black_box(matvec_scalar(&matrix, &vector, rows, cols));
        black_box(matvec_simd(&matrix, &vector, rows, cols));

        // Benchmark scalar
        let start = Instant::now();
        for _ in 0..iterations {
            black_box(matvec_scalar(&matrix, &vector, rows, cols));
        }
        let scalar_us = start.elapsed().as_micros() / iterations as u128;

        // Benchmark SIMD
        let start = Instant::now();
        for _ in 0..iterations {
            black_box(matvec_simd(&matrix, &vector, rows, cols));
        }
        let simd_us = start.elapsed().as_micros() / iterations as u128;

        let speedup = scalar_us as f64 / simd_us as f64;

        println!(
            "{:<20} {:>10} {:>10} {:>10} {:>10} {:>12.2}x",
            name, rows, cols, scalar_us, simd_us, speedup
        );
    }

    println!("\n=== GPU Capability Check ===");

    #[cfg(feature = "gpu")]
    {
        use trueno::backends::gpu::GpuDevice;

        match GpuDevice::new() {
            Ok(device) => {
                println!("GPU device initialized successfully");

                // Benchmark GPU matvec for large dimensions
                println!("\n=== GPU Matvec Benchmark (Large Dimensions) ===");
                println!(
                    "{:<20} {:>10} {:>10} {:>10} {:>12}",
                    "Operation", "Rows", "Cols", "SIMD(µs)", "GPU(µs)"
                );
                println!("{}", "-".repeat(65));

                for (name, rows, cols) in &[
                    ("llama 7B hidden", 4096, 4096),
                    ("llama 70B hidden", 8192, 8192),
                ] {
                    let matrix: Vec<f32> =
                        (0..(rows * cols)).map(|i| (i as f32) * 0.0001).collect();
                    let vector: Vec<f32> = (0..*cols).map(|i| (i as f32) * 0.001).collect();

                    // SIMD timing
                    let iterations = 20;
                    let start = Instant::now();
                    for _ in 0..iterations {
                        black_box(matvec_simd(&matrix, &vector, *rows, *cols));
                    }
                    let simd_us = start.elapsed().as_micros() / iterations as u128;

                    // GPU timing via matmul (treating vector as 1-column matrix)
                    let mut result = vec![0.0f32; *rows];
                    let start = Instant::now();
                    for _ in 0..iterations {
                        device
                            .matmul(&matrix, &vector, &mut result, *rows, *cols, 1)
                            .expect("GPU matmul failed");
                        black_box(&result);
                    }
                    let gpu_us = start.elapsed().as_micros() / iterations as u128;

                    println!(
                        "{:<20} {:>10} {:>10} {:>10} {:>12}",
                        name, rows, cols, simd_us, gpu_us
                    );
                }

                println!("\n=== Falsification Analysis ===");
                println!("- Token generation uses MATVEC (n=1), not GEMM");
                println!("- Research shows cuBLAS 'hurts more than helps' for MATVEC");
                println!("- GPU launch overhead dominates for latency-bound operations");
                println!("- SIMD with 4-accumulator dot product is likely FASTER for matvec");
            },
            Err(e) => {
                println!("GPU not available: {}", e);
                println!("\nFallback analysis: SIMD-only performance for matvec");
            },
        }
    }

    #[cfg(not(feature = "gpu"))]
    {
        println!("GPU feature not enabled. Run with: cargo run --example gpu_matvec_benchmark --features gpu");
    }

    println!("\n=== Falsifiable Claims (IMP-600b) ===");
    println!("1. SIMD matvec should be 4-7x faster than scalar (PROVEN by 4-accum)");
    println!("2. GPU overhead makes it SLOWER than SIMD for small-medium matvec");
    println!("3. GPU only wins for very large batch operations (GEMM with n >> 1)");
    println!("4. trueno CAN achieve parity because matvec is latency-bound, not throughput-bound");
}