scipy/sprs random sparse matrices multiplication

main.py

import time
import numpy as np
from scipy.sparse import csr_matrix


def make_matrix(shape, nnz):
    i = np.random.randint(0, shape[0] - 1, nnz)
    j = np.random.randint(0, shape[1] - 1, nnz)
    data = np.random.randn(nnz)

    return csr_matrix((data, (i, j)), shape=shape)


def main():
    # shape = (1000000, 1000000)
    # shape = (1500, 2500)
    shape = (15000, 25000)
    # shape = (150000, 250000)
    # shape = (1500000, 2500000)
    # shape = (15000000, 25000000)

    nnzs = [
        100, 500, 1000, 2000, 3000, 4000, 5000, 7500, 10000, 20000, 30000, 40000, 50000,
        60000, 70000, 80000, 90000, 100000, 125000, 150000, 175000, 200000, 300000, 400000,
        500000, 750000, 1000000, 2000000, 3000000,
    ]

    elapseds = []
    # nnzs = []

    for nnz in nnzs:
        m1 = make_matrix(shape, nnz)
        m2 = make_matrix(shape[::-1], nnz)

        print('\n`m1` shape={}, nnz={}'.format(m1.shape, m1.nnz))
        print('`m2` shape={}, nnz={}'.format(m2.shape, m2.nnz))

        t0 = time.monotonic()
        mm = m1 @ m2
        t1 = round((time.monotonic() - t0) * 1000)
        print(" -> m1 * m2 (mm_nnz={}): {:.2f} msec".format(mm.nnz, t1))
        elapseds.append(t1)
        # nnzs.append(m1.nnz)

    print("\n`m1 * m2` nnz: {}".format(nnzs))
    print("\n`m1 * m2` elapsed msec: {}".format(elapseds))


if __name__ == '__main__':
    np.random.seed(2846293)
    main()

main.rs

use std::time::Instant;

use ndarray::{Array1};
use ndarray_rand::RandomExt;
use ndarray_rand::rand_distr::{Uniform, Normal};

use sprs::{TriMat, CsMat, Shape};


fn make_matrix(shape: Shape, nnz: usize) -> CsMat<f64> {
    let rows = Array1::<usize>::random((nnz,), Uniform::new(0, shape.0 - 1));
    let cols = Array1::<usize>::random((nnz,), Uniform::new(0, shape.1 - 1));
    let data = Array1::<f64>::random((nnz,), Normal::new(0.0, 1.5).unwrap());

    let mut m = TriMat::new(shape);

    for ((&i, &j), &v) in rows.iter().zip(cols.iter()).zip(data.iter()) {
        m.add_triplet(i, j, v);
    }

    m.to_csr()
}


fn main() {
    // let shape: Shape = (1500, 2500);
    let shape: Shape = (15000, 25000);
    // let shape: Shape = (150000, 250000);
    // let shape: Shape = (1500000, 2500000);

    let nnzs: &[usize] = &[
        100, 500, 1000, 2000, 3000, 4000, 5000, 7500, 10000, 20000, 30000, 40000, 50000,
        60000, 70000, 80000, 90000, 100000, 125000, 150000, 175000, 200000, 300000, 400000,
        500000, 750000, 1000000, 2000000, 3000000,
    ];

    let mut elapseds = Vec::<u128>::new();

    for &nnz in nnzs {
        let m1 = make_matrix(shape, nnz);
        let m2 = make_matrix((shape.1, shape.0), nnz);

        println!("\n`m1` shape={:?}, nnz={}", m1.shape(), m1.nnz());
        println!("`m2` shape={:?}, nnz={}", m2.shape(), m2.nnz());
        let now = Instant::now();
        let mm = &m1 * &m2;
        let elapsed = now.elapsed().as_millis();
        println!(" -> m1 * m2 (nnz={}): {} msec", mm.nnz(), elapsed);
        elapseds.push(elapsed);
    }

    println!("\n`m1 * m2` nnz: {:?}", nnzs);
    println!("`m1 * m2` elapsed msec: {:?}", elapseds);
}