2

I would like to have an implementation in PARI/GP for the calculation of

a_1 ^ a_2 ^ ... ^ a_n (mod m)

which manages all cases, especially the cases where high powers appear in the phi-chain.

Does anyone know such an implementation ?

max taldykin
  • 12,459
  • 5
  • 45
  • 64
Peter
  • 214
  • 1
  • 10

1 Answers1

0

Here's a possibility using Chinese remainders to make sure the modulus is a prime power. This simplifies the computation of x^n mod m in the painful case where gcd(x,m) is not 1. The code assumes the a_i are > 1; most of the code checks whether p^a_1^a_2^...^a_n is 0 mod (p^e) for a prime number p, while avoiding overflow.

\\ x[1]^x[2]^ ...^ x[#x] mod m, assuming x[i] > 1 for all i
tower(x, m) =
{ my(f = factor(m), P = f[,1], E = f[,2]);
  chinese(vector(#P, i, towerp(x, P[i], E[i])));
}

towerp(x, p, e) =
{ my(q = p^e, i, t, v);

  if (#x == 0, return (Mod(1, q)));
  if (#x == 1, return (Mod(x[1], q)));
  if (v = valuation(x[1], p),
    t = x[#x]; i = #x;
    while (i > 1,
      if (t >= e, return (Mod(0, q)));
      t = x[i]^t; i--);
    if (t * v >= e, return (Mod(0, q)));
    return (Mod(x[1], q)^t);
  );
  Mod(x[1], q)^lift(tower(x[^1], (p-1)*p^e));
}

For instance

? 5^(4^(3^2)) % 163 \\ direct computation, wouldn't scale
%1 = 158
? tower([5,4,3,2], 163)
%2 = Mod(158, 163)
K.B.
  • 861
  • 5
  • 14