multiadd.rs

#[bench]
fn multiadd(b: &mut test::Bencher) {
    use std::vec::Vec;

    fn batch_invert(v: &mut [Fp]) -> Fp {
        let mut tmp = Vec::with_capacity(v.len());

        let mut acc = Fp::one();
        for p in v.iter() {
            tmp.push(acc);
            acc = acc * p;
        }

        acc = acc.invert().unwrap();
        let allinv = acc;

        for (p, tmp) in v.iter_mut().rev().zip(tmp.into_iter().rev()) {
            let tmp = tmp * acc;
            acc = acc * (*p);
            *p = tmp;
        }

        allinv
    }

    const K: usize = 20;
    const N: usize = 1 << K;

    let mut vp = Vec::with_capacity(N);
    let mut cur = G1Projective::generator();
    for _ in 0..N {
        vp.push(cur);
        cur = cur.double();
    }
    let mut v = vec![G1Affine::generator(); N];
    G1Projective::batch_normalize(&vp[..], &mut v[..]);

    b.iter(|| {
        let mut v = v.clone();
        let mut n = N;
        while n > (1 << 7) { // switch to mixed addition when the inversions get too expensive
            assert_eq!(v.len(), n);
            let mut tmp = vec![Fp::zero(); n >> 1];
            for i in 0..(n >> 1) {
                tmp[i] = v[2 * i].x - v[2 * i + 1].x;
            }
            batch_invert(&mut tmp);
            for i in 0..(n >> 1) {
                let lambda = v[2 * i].y - v[2 * i + 1].y;
                let lambda = lambda * tmp[i];
                let lambda2 = lambda.square();
                let xr = lambda2 - v[2 * i].x - v[2 * i + 1].x;
                let yr = lambda * (v[2 * i + 1].x - xr) - v[2 * i + 1].y;
                v[i].x = xr;
                v[i].y = yr;
            }
            v.truncate(n >> 1);
            n = n >> 1;
        }
        G1Affine::from(v.into_iter().fold(G1Projective::identity(), |a, b| {
            b + a
        }))
    });
}