MatrixPower with Modulus

Question

I have a matrix that will be multiplied with itself by an extremely large amount of times under $Z/pZ$. The matrix itself contains small numbers, so MatrixPower[Mod[mat, p], pow, vec] won't work; and taking $\bmod$ after calculation is infeasible since the result is extremely large. Can I have similar functionality like some other function's Modulus -> p with MatrixPower?

score 13 · Accepted Answer · edited Aug 07 '16 at 06:10

13

Using an undocumented function:

Algebra`MatrixPowerMod[{{2, 3, 1}, {5, 2, 4}, {0, 3, 2}}, 4, 6]

{{1, 0, 5}, {1, 1, 5}, {3, 3, 4}}

Mod[MatrixPower[{{2, 3, 1}, {5, 2, 4}, {0, 3, 2}}, 4], 6]

{{1, 0, 5}, {1, 1, 5}, {3, 3, 4}}

Unfortunately, the undocumented function does not support the action form of MatrixPower[].

edited Aug 07 '16 at 06:10

Alexey Popkov

61,809
7
149
368

answered Aug 05 '16 at 13:04

J. M.'s missing motivation

124,525
11
401
574

This is exactly what I want, PE502 part 2 solved:) – vapor Aug 05 '16 at 13:10
Another one for (805). @happyfish There are over 500 PE problems now? That's a lot to work through! I stopped when there were only 330 or so, and I was getting repeatedly stumped by ones in the 130's. How many have you solved? – Mr.Wizard Aug 05 '16 at 13:26
2

@Mr.Wizard I have done nearly 200, most of them are simple ones. But 502 is on the top of their difficulty list and I am nearly done :) – vapor Aug 05 '16 at 16:29
1

This function also support negative exponent (modular inverse). – user202729 Jan 02 '20 at 01:53

score 5 · Answer 2 · edited Aug 07 '16 at 01:03

5

Bressoud & Wagon define MatrixPowerMod on page 34 of their book, A Course in Computational Number Theory.

MatrixPowerMod[a_, n_, m_] :=
   Block[{b = a, d = IntegerDigits[n, 2]},
         Do[
            b = Mod[b.b, m];
            If[d[[i]] == 1, b = Mod[b.a, m]],
            {i, 2, Length[d]}];
         b]

An example involving Fibonacci numbers is the following:

AbsoluteTiming[MatrixPowerMod[{{0, 1}, {1, 1}}, 10^8, 10^10]]
(* {0.005576, {{1300390626, 7760546875}, {7760546875, 9060937501}}} *)

edited Aug 07 '16 at 01:03

J. M.'s missing motivation

124,525
11
401
574

answered Aug 06 '16 at 15:17

KennyColnago

15,209
26
62

2

This is, of course, the "Russian peasant" algorithm for exponentiation. Here is a compact implementation: Fold[Mod[If[#2 == 1, a.#, #].#, m] &, a, Rest[IntegerDigits[n, 2]]]. – J. M.'s missing motivation Aug 06 '16 at 15:21

MatrixPower with Modulus

2 Answers2

Linked