Important note: this article is outdated! Go to http://zsiciarz.github.io/24daysofrust/ for a recent version of all of 24 days of Rust articles. The blogpost here is kept as it is for historical reasons.

When I start learning a new programming language, I like to code at least several solutions to Project Euler problems. These are very math-oriented and may not be the best introduction to general purpose programming, but it's a start. Anyway, it's just fun to solve them! (...and way more fun to solve them in a fast way and not by brute force.)

A lot of Project Euler problems involve prime numbers in some way. These include finding nth prime, efficient factorization or checking whether some curious number is prime or not. You could of course write these mathematical procedures yourself, which is also an educational activity. But I'm lazy. I set out to find some ready-made code and stumbled upon the slow_primes library by Huon Wilson. Incidentally this was the first external dependency I ever used in a Rust program, long before crates.io.

By the way don't let the name fool you, it's not that slow. As Huon says:

Despite the name, it can sieve the primes up to 10⁹ in about 5 seconds.

So let's see what's in there, shall we?

Prime sieve

The first thing to do is to create a sieve (see Wikipedia on Sieve of Eratosthenes for a detailed explanation of the algorithm). We need to set an upper bound on the sieve. There's a clever way to estimate that bound (see the docs for estimate_nth_prime) but for simplicity I'll hardcode it now to 10000.

Let's actually check some numbers for primality:

extern crate slow_primes;
use slow_primes::Primes;

fn main() {
    let sieve = Primes::sieve(10000);
    let suspect = 5273u;
    println!("{} is prime: {}", suspect, sieve.is_prime(suspect)); // true
    let not_a_prime = 1024u;
    println!("{} is prime: {}", not_a_prime, sieve.is_prime(not_a_prime)); // guess
}

How about finding 1000th prime number?

let n = 1000u;
match sieve.primes().nth(n - 1) {
    Some(number) => println!("{}th prime is {}", n, number),
    None => println!("I don't know anything about {}th prime.", n),
}

The primes() method returns an iterator over all prime numbers generated by this sieve (2, 3, 5, 7...). Iterators in Rust have a lot of useful methods; the nth() method skips over n initial iterations, returning the nth element (or None if we exhausted the iterator). The argument is zero-based, so to find 1000th prime we need to pass 999 to nth().

Factorization

Factorization is a way to decompose a number into its divisors. For example, 2610 = 2 * 3 * 3 * 5 * 29. Here's how we can find it out with slow_primes API:

println!("{}", sieve.factor(2610));

When we run this, we'll get:

$ cargo run
Ok([(2, 1), (3, 2), (5, 1), (29, 1)])

What is this? Let's have a look at the result type of factor():

type Factors = Vec<(uint, uint)>;
fn factor(&self, n: uint) -> Result<Factors, (uint, Factors)>

Looks a bit complicated, but remember the Result type. The Ok variant wraps a vector of pairs of numbers. Each pair contains a prime factor and its exponent (how many times it appears in the factorization). In case of an error we'll get a pair (leftover value, partial factorization).

We can use factorization to find the total number of divisors (including compound ones). This is very important in number theory (although for reasons that are outside the scope of this blog).

Consider the following function:

fn num_divisors(n: uint, primes: &Primes) -> Option<uint> {
    use std::iter::MultiplicativeIterator;
    match primes.factor(n) {
        Ok(factors) => Some(factors.into_iter().map(|(_, x)| x + 1).product()),
        Err(_) => None,
    }
}

The trick is to multiply all prime factor exponents, incremented before multiplication. See the explanation at Maths Challenge for the curious. So when we call the function on our 2610 example, we'll get Some(24) as a result.

println!("{}", num_divisors(2610, &sieve));

Further reading

• Divisor function
• Miller-Rabin primality test
• GMP algorithms page
• consecutive prime sum - an interesting Project Euler problem

Code examples in this article were built with rustc 0.13.0-nightly and slow_primes 0.1.4.

The header photo (taken by me) shows Tamka street in Warsaw, Poland.