Square-and-Multiply Algorithm Proof
📂Number TheorySquare-and-Multiply Algorithm Proof
Algorithm
Given a natural number a,k,m, b≡ak(modm) can be computed as follows.
Step 1. Binary Expansion of k
Represent ui=0 or ui=1 as follows.
k=i=0∑rui2i=u0+2u1+⋯+2rur
Step 2.
Calculate a2r≡(a2r−1)2≡Ar−12≡Ar(modm) as follows:
aa2a4a8≡(a1)2≡(a2)2≡(a4)2≡A02≡A12≡A22≡A0(modm)≡A1(modm)≡A2(modm)≡A3(modm)⋮
Step 3.
Calculating A0u0⋅A1u1⋅⋯Arur yields:
ak≡A0u0⋅A1u1⋅⋯Arur(modm)
Explanation
Although the advancement of computers has made it possible to compute powers quickly, humans are never satisfied. The method of repeated squaring is a way to get the answer without directly computing immensely large numbers in modular operations. By dividing by (modm) after each square, the number falls below m, thus reducing the amount of computation.
Example
For instance, let’s calculate 7327(mod853). Taking the common logarithm of 7327 gives 327log7≈327⋅0.8450=276.315, which is a large number with a digit count of 277. It’s too big to multiply honestly and then divide, so let’s use the above algorithm. If you do the binary expansion of 327 as described in the algorithm, it looks like this:
327=256+64+4+2+1=28+26+22+21+20
If you calculate according to Step 2,
772747871673276471287256≡(71)2≡(72)2≡(74)2≡(78)2≡(716)2≡(732)2≡(764)2≡(7128)2≡7≡72≡492≡6952≡2272≡3492≡6752≡1232≡6282≡7(mod853)≡49(mod853)≡2401≡695(mod853)≡227(mod853)≡349(mod853)≡675(mod853)≡123(mod853)≡628(mod853)≡298(mod853)
Thus,
7327=7256⋅764⋅74⋅72⋅71=298⋅123⋅695⋅49⋅7≡286(mod853)
This method might still feel hard, but the amount of computation is quite reasonable compared to multiplying 327 times individually.
Proof
ak==≡au0+2u1+⋯+2rurau0a2u1⋯a2rurA0u0A1u1⋯Arur(modm)
■
Code
Below is an example code for calculating 7327(mod853) using the repeated squaring method, written in R language. Note that the power
option accommodates negative integers, and if mod
is a prime number, entering power=-1
allows you to find the multiplicative inverse of the given base
.
FPM<-function(base,power,mod)
{
i<-0
if (power<0) {
while((base*i)%%mod != 1) {i=i+1}
base<-i
power<-(-power)}
if (power==0) {return(1)}
if (power==1) {return(base%%mod)}
n<-0
while(power>=2^n) {n=n+1}
A<-rep(1,n)
A[1]=base
for(i in 1:(n-1)) {A[i+1]=(A[i]^2)%%mod}
for(i in n:1) {
if(power>=2^(i-1)) {power=power-2^(i-1)}
else {A[i]=1} }
for(i in 2:n) {A[1]=(A[1]*A[i])%%mod}
return(A[1])
}
FPM(7,327,853)
FPM(7,-1,11)