Name

number_powermod — Modular exponentiation.

Calling Sequence

   x = number_powermod ( a, n, m )
   
   

Parameters

a :

a 1x1 matrix of floating point integers, must be positive

n :

a 1x1 matrix of floating point integers, must be positive

m :

a 1x1 matrix of floating point integers, must be positive

x :

a 1x1 matrix of floating point integers

Description

Returns x = mod ( a^n, m ).

Uses a repeated squaring algorithm. The algorithm may fail if the number is so large that the square a^2 is out of range.

Some programming tricks are used to speed up the computation, and to allow computations in which A**N is much too large to store in a real word.

First, for efficiency, the power A**N is computed by determining the binary expansion of N, then computing A, A**2, A**4, and so on by repeated squaring, and multiplying only those factors that contribute to A**N.

Secondly, the intermediate products are immediately "mod'ed", which keeps them small.

For instance, to compute mod ( A**13, 11 ), we essentially compute 13 = 1 + 4 + 8 then A**13 = A * A**4 * A**8. This leads to mod ( A**13, 11 ) = mod ( A, 11 ) * mod ( A**4, 11 ) * mod ( A**8, 11 ).

Fermat's little theorem says that if P is prime, and A is not divisible by P, then ( A**(P-1) - 1 ) is divisible by P.

Examples

number_powermod( 2 , 3 , 7 ) // 1
number_powermod( 4 , 5 , 6 ) // 4

// This generates an error, as expected :
// intermediate terms are not in the range
a = 1028712535855
u = 34630287489
n = 4432676798593
number_powermod( a , u , n ) // 678384139479

   

Authors

Copyright (C) 2009 - 2010 - DIGITEO - Michael Baudin
Copyright (C) 2004 - John Burkardt