The standard trick for computing a^p modulo m
is to use successive square. The idea is to expand p
into binary, say
p = e0 * 2^0 + e1 * 2^1 + ... + en * 2^n
where (e0,e1,...,en)
are binary (0
or 1
) and en = 1
. Then use laws of exponents to get the following expansion for a^p
a^p = a^( e0 * 2^0 + e1 * 2^1 + ... + en * 2^n )
= a^(e0 * 2^0) * a^(e1 * 2^1) * ... * a^(en * 2^n)
= (a^(2^0))^e0 * (a^(2^1))^e1 * ... * (a^(2^n))^en
Remember that each ei
is either 0
or 1
, so these just tell you which numbers to take. So the only computations that you need are
a, a^2, a^4, a^8, ..., a^(2^n)
You can generate this sequence by squaring the previous term. Since you want to compute the answer mod m
, you should do the modular arithmetic first. This means you want to compute the following
A0 = a mod m
Ai = (Ai)^2 mod m for i>1
The answer is then
a^p mod m = A0^e0 + A1^e1 + ... + An^en
Therefore the computation takes log(p)
squares and calls to mod m
.
I'm not certain whether or not there is an analog for factorials, but a good place to start looking would be at Wilson's Theorem. Also, you should put in a test for m <= n
, in which case n! mod m = 0
.