As both a functional language and a language with support for arbitrary precision integers, Scheme is a beautiful language for describing these algorithms.
The Fermat primality algorithm and the Solovay-Strassen primality algorithm are, by themselves, straightforward.
The Fermat test
The Fermat primality test is based on Fermat's little theorem.
Fermat's little theorem implies that for any prime $p$ and any natural $a$ such that $1 \leq a < p$, it must be the case that
\[ a^{p-1} \mod p = 1\text. \]However, for any given composite number $n$ and some natural $a$ such that $1 \leq a < n$, it is usually not the case that:
\[ a^{n-1} \mod n = 1\text. \]Given a number $m$, we can repeatedly test whether $a^{m-1} \mod m = 1$ for many different values of $a$.
Each time we get back 1, we become increasingly confident that the natural $m$ is truly prime.
However, if we ever get back any value other than 1, we know immediately that $m$ is composite.
Warning: Charmichael numbers
It must be noted that there are composites, known as Charmichael numbers, that can pass this test for any value of $a$.
In general, Charmichael numbers are sufficiently rare that this test can be used in practice, with caution.
561 is the lowest Charmichael number.
The Solovay-Strassen test
The Solovay-Strassen test uses a slightly more complex test than the Fermat test, but it provides a tight, probabilistic guarantee on the result.
Euler's criterion implies that for any $a$ coprime with some odd prime $p$:
\[ a^{(p-1)/2} \mod p = \left( \frac{a}{p} \right) \][The notation $(\frac{a}{m})$ denotes the Jacobi symbol.]
Meanwhile, for a composite $m$, given an arbitrary natural $a$, there is less than a 50% chance that:
\[ a^{(m-1)/2} \mod m = \left( \frac{a}{m} \right) \]By successively testing the condition with different values of $a$, the probability that $m$ is composite falls by at least half with each successful test.
Thus, it becomes possible to choose a level of confidence that $m$ is prime and perform the necessary number of tests to achieve that.
The complex part of the algorithm is implementing the Jacobi symbol:
(define (jacobi a n)
(cond
((= n 1) 1)
((= a 1) 1)
((not (= (gcd a n) 1)) 0)
((and (= a 2)
(let ((n-mod-8 (modulo n 8)))
(cond
((or (= n-mod-8 1) (= n-mod-8 7)) 1)
((or (= n-mod-8 3) (= n-mod-8 5)) -1)))))
((> a n) (jacobi (modulo a n) n))
((even? a) (* (jacobi (/ a 2) n) (jacobi 2 n)))
((even? n) (* (jacobi a (/ n 2)) (jacobi a 2)))
(else (* (jacobi n a) (if (even? (/ (* (- a 1) (- n 1)) 4)) 1 -1)))))
Code
Warning: If you want to use this in practice, you'll also have to add a weak-key detection pass to make it reasonably secure.