python - Numpy matrix power/exponent with modulo?

Question

Is it possible to use numpy's linalg.matrix_power with a modulo so the elements don't grow larger than a certain value?

Answer 1

In order to prevent overflow, you can use the fact that you get the same result if you first take the modulo of each of your input numbers; in fact:

(M**k) mod p = ([M mod p]**k) mod p,

for a matrix M. This comes from the following two fundamental identities, which are valid for integers x and y:

(x+y) mod p = ([x mod p]+[y mod p]) mod p  # All additions can be done on numbers *modulo p*
(x*y) mod p = ([x mod p]*[y mod p]) mod p  # All multiplications can be done on numbers *modulo p*

The same identities hold for matrices as well, since matrix addition and multiplication can be expressed through scalar addition and multiplication. With this, you only exponentiate small numbers (n mod p is generally much smaller than n) and are much less likely to get overflows. In NumPy, you would therefore simply do

((arr % p)**k) % p

in order to get (arr**k) mod p.

If this is still not enough (i.e., if there is a risk that [n mod p]**k causes overflow despite n mod p being small), you can break up the exponentiation into multiple exponentiations. The fundamental identities above yield

(n**[a+b]) mod p = ([{n mod p}**a mod p] * [{n mod p}**b mod p]) mod p

and

(n**[a*b]) mod p = ([n mod p]**a mod p)**b mod p.

Thus, you can break up power k as a+b+… or a*b*… or any combination thereof. The identities above allow you to perform only exponentiations of small numbers by small numbers, which greatly lowers the risk of integer overflows.