Problem 21 "Amicable numbers"

Let d(n) be defined as the sum of proper divisors of n (numbers less than n which divide evenly into n).
If d(a) = b and d(b) = a, where ab, then a and b are an amicable pair and each of a and b are called amicable numbers.

For example, the proper divisors of 220 are 1, 2, 4, 5, 10, 11, 20, 22, 44, 55 and 110; therefore d(220) = 284. The proper divisors of 284 are 1, 2, 4, 71 and 142; so d(284) = 220.

Evaluate the sum of all the amicable numbers under 10000.

問 21 「友愛数」

d(\(n\)) を \(n\) の真の約数の和と定義する. (真の約数とは \(n\) 以外の約数のことである. )

もし, d(\(a\)) = \(b\); かつ d(\(b\)) = \(a\); (\(a\) ≠ \(b\) のとき) を満たすとき, \(a\) と \(b\) は友愛数(親和数)であるという.

例えば, 220 の約数は 1, 2, 4, 5, 10, 11, 20, 22, 44, 55, 110 なので d(220) = 284 である.

また, 284 の約数は 1, 2, 4, 71, 142 なので d(284) = 220 である.

それでは10000未満の友愛数の和を求めよ.

struct AmicableNumberScanner {
    below: u32,
    _primes: Vec<u32>,
    _aliquot_sum: Vec<u32>,
}

impl AmicableNumberScanner {
    fn new(below: u32) -> Self {
        Self {
            below: below,
            _primes: vec![2, 3],
            _aliquot_sum: vec![0u32; below as usize],
        }
    }
    fn _divide_fully(&self, n: &mut u32, d: u32, side: &mut u32, sum: &mut u32) {
        if *n % d != 0 {
            return;
        }
        let mut exp = 0u32;
        while {
            *n /= d;
            exp += 1;
            *n % d == 0
        } {}
        *side = (*n as f32).sqrt() as u32;
        *sum *= (d.pow(exp + 1) - 1) / (d - 1);
    }
    fn _sum_of_divisors(&mut self, mut n: u32) -> u32 {
        let mut side = (n as f32).sqrt() as u32;
        let mut sum = 1u32;
        for &p in &self._primes {
            if p > side || n == 1 {
                break;
            }
            self._divide_fully(&mut n, p, &mut side, &mut sum);
        }
        if n != 1 {
            sum *= (n * n - 1) / (n - 1);
            self._primes.push(n);
        }
        sum
    }
    fn pair_sum(&mut self) -> u32 {
        for n in 4..self.below {
            self._aliquot_sum[n as usize] = self._sum_of_divisors(n) - n;
        }
        let mut sum = 0u32;
        for (i, v) in &self._aliquot_sum.enumerate() {
            let vsize = *v as usize;
            if vsize >= self._aliquot_sum.len() {
                continue;
            }
            if vsize == i {
                continue;
            }
            if self._aliquot_sum[vsize] == i as u32 {
                sum += *v;
            }
        }
        sum
    }
}

fn main() {
    let sum = AmicableNumberScanner::new(10_000).pair_sum();

    println!("{}", sum);
    assert_eq!(sum, 31626);
}

struct Index {
    i: usize,
    _ite: Box<dyn Iterator<Item = usize>>,
}

impl Index {
    fn increment(&mut self) {
        self.i += self._ite.next().unwrap();
    }
    fn new() -> Self {
        Self {
            i: 5,
            _ite: Box::new(vec![2usize, 4].into_iter().cycle()),
        }
    }
}

fn rule_out(sieve: &mut Vec<bool>, prime: usize) {
    for i in (prime * prime..sieve.len()).step_by(prime) {
        sieve[i] = false;
    }
}

fn primes(below: u32) -> Vec<u32> {
    let mut primes: Vec<u32> = vec![2u32, 3u32];
    let mut sieve = vec![true; below as usize];
    let sqrt = (sieve.len() as f32).sqrt() as usize;
    let mut index = Index::new();
    while index.i <= sqrt {
        if sieve[index.i] {
            primes.push(index.i as u32);
            rule_out(&mut sieve, index.i);
        }
        index.increment();
    }
    while index.i < sieve.len() {
        if sieve[index.i] {
            primes.push(index.i as u32);
        }
        index.increment();
    }
    primes
}

struct AmicableNumberScanner {
    below: u32,
    _primes: Vec<u32>,
    _checked: Vec<bool>,
}

impl AmicableNumberScanner {
    fn new(below: u32) -> Self {
        Self {
            below: below,
            _primes: primes(below),
            _checked: vec![false; below as usize],
        }
    }
    fn _divide_fully(&self, n: &mut u32, d: u32, side: &mut u32, sum: &mut u32) {
        if *n % d != 0 {
            return;
        }
        let mut exp = 0u32;
        while {
            *n /= d;
            exp += 1;
            *n % d == 0
        } {}
        *side = (*n as f32).sqrt() as u32;
        *sum *= (d.pow(exp + 1) - 1) / (d - 1);
    }
    fn _sum_of_divisors(&mut self, mut n: u32) -> u32 {
        let mut side = (n as f32).sqrt() as u32;
        let mut sum = 1u32;
        for &p in self._primes.iter() {
            if p > side || n == 1 {
                break;
            }
            self._divide_fully(&mut n, p, &mut side, &mut sum);
        }
        if n != 1 {
            sum *= (n * n - 1) / (n - 1);
        }
        sum
    }
    fn pair_sum(&mut self) -> u32 {
        let mut pair_sum = 0u32;
        for a in 4..self.below {
            if self._checked[a as usize] {
                continue;
            }
            let sum = self._sum_of_divisors(a) - a;
            if sum <= a {
                continue;
            }
            if sum >= self.below {
                continue;
            }
            let b = self._sum_of_divisors(sum) - sum;
            self._checked[sum as usize] = true;
            if a == b {
                pair_sum += a + sum;
            }
        }
        pair_sum
    }
}

fn main() {
    let sum = AmicableNumberScanner::new(10_000).pair_sum();

    println!("{}", sum);
    assert_eq!(sum, 31626);
}