1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56

pub fn extended_gcd(a: i64, b: i64, p: &mut i64, q: &mut i64) -> i64 {
    if b == 0 {
        *p = 1;
        *q = 0;
        a
    } else {
        let d = extended_gcd(b, a % b, q, p);
        *q -= a / b * *p;
        d
    }
}

pub fn chinese_remainder_theorem(b: &[i64], modulo: &[i64]) -> Option<(i64, i64)> {
    let (mut result, mut m) = (0, 1);
    for i in 0..b.len() {
        let (mut p, mut q) = (0, 0);
        let d = extended_gcd(m, modulo[i], &mut p, &mut q);
        if (b[i] - result) % d != 0 {
            return None;
        }
        let tmp = ((b[i] - result) / d * p) % (modulo[i] / d);
        result += m * tmp;
        m *= modulo[i] / d;
    }
    Some(((result % m + m) % m, m))
}

#[cfg(test)]
mod tests {
    use super::*;
    use rand;
    use rand::Rng;

    #[test]
    fn test_crt() {
        let mut rng = rand::thread_rng();
        let n = 10;
        let max_m = 100;
        for _ in 0..1000 {
            let ans = rng.gen::<u32>() as i64;
            let mut b = vec![0; n];
            let mut m = vec![0; n];
            for i in 0..n {
                m[i] = rng.gen::<u8>() as i64;
                m[i] %= max_m;
                m[i] += 1;
                b[i] = ans % m[i];
            }

            let (a, _) = chinese_remainder_theorem(&b, &m).unwrap();
            for i in 0..n {
                assert_eq!(a % m[i], b[i]);
            }
        }
    }
}