Oregon Curriculum Network
Discovering Math with Python
Calling modulo arithmetic "clock arithmetic" is a way of calling attention to the finite cyclic nature of always taking remainders, factoring away the "modulus". What does that mean?
Clock Arithmetic is important to Number Theory, and Number Theory is important to crytography.
Cyptography, in the not so distant past, would not concern an average civilian, whereas today the HTTPS protocol in every mainstream web browser (Mozilla, Chrome, IE, Safari, Opera...) is capable of implementing what we called TLS, or Transport Layer Security.
The web browser client, and server, will handshake, meaning agree on what cryptographic strategy to use: first a public one, to exchange a symmetric key, and then a symmetric one, most likely AES at the time of this writing.
AES is the Advanced Encryption Standard, the winner of a global contest as adjudged by NIST. AES uses "clock arithemetic" internally in order to apply several layers of "mixing it up" by invertible methods. Enciphering is reverisble (the whole point) but only if the key is handy.
In Python, one of our primitive operators is %, for "modulo", which relates to the built-in function divmod( ). Lets see how % and divmod do their jobs:
In [1]:
(5 + 5) % 10 # no remainder
Out[1]:
In [2]:
divmod(10, 5) # 5 goes into 10 twice, no remainder
Out[2]:
In [3]:
19 % 10 # divide by 10, give the remainder
Out[3]:
Try doing a number of these yourself. Guess the answer before hitting the Shift-Enter key.
You may find your first years of formal schooling did not emphasize cyclic arithmetic, except in the two most common forms, of telling time and scheduling appointments. We understand both the 24 hour day, and the 365 day year, comprise a modulus such that we're always talking about our degree of offset into each one. Adding hours or days, per a timedelta object, keeps us within the scope of some calendar.
The unit circle of 360 degrees, or 2 pi radians, is also treated in a clock-like fashion. Going around thousands of degrees from some present position, never takes us beyond 360. We're confined to a finite domain. Cyclic phenomena more generally are accommodated by our clock arithmetic, of complex numbers, of e to imaginary powers.
There's a kind of "closure" in this picture, which might be resonating with you by now, as one of those properties of a group. Will we find groups in this Chapter? You bet. Groups also enjoy the symmetry of a circle in that each member is paired with another "180 degrees opposite" although sometimes the inverse might be the same as the self.
Again, we call this "clock arithmetic" because even if we say "20 hours from now", or "30 hours ago" we'll still be somewhere on the circle marked out into 12 intervals, each representing one hour. The yearly calendar is likewise a kind of modulus. No matter how many days we add, we're still somewhere between 0 and 365 days into some year.
We also call it "modulo arithmetic" meaning the computations are always vis-a-vis some fixed modulus. This word "modulus" is our inspiration for the M-numbers coded below.
Lets build a class, the instances of which will multiply and add per some fixed modulus, meaning we're always factoring out the modulus and keeping the remainder.
In [4]:
class M: # for "modulo"
modulus = 10 # class level
def __init__(self, val, modulus=None):
if modulus:
self.modulus = M.modulus = modulus # resettable
else:
self.modulus = M.modulus
self.val = val % M.modulus
def __add__(self, other):
if self.modulus != other.modulus:
raise ValueError
return M((self.val + other.val) % self.modulus)
def __mul__(self, other):
if self.modulus != other.modulus:
raise ValueError
return M((self.val * other.val) % self.modulus)
def __pow__(self, exp):
raise NotImplemented
def __repr__(self):
return "(" + str(self.val) + " mod " + str(self.modulus)+ ")"
a = M(8)
b = M(7)
print(a, b)
print("(8 mod 10) * (7 mod 10) = ", a * b)
print("(8 mod 10) + (7 mod 10) = ", a + b)
OK, everything seems to be working, even though we haven't implemented powering yet. Eventually we'd like to go pow(M(3), -1) to get the inverse of M(3), such that M(3) times its inverse equals the multiplicative identity M(1).
But wait, does every M-number, with modulus set to 10, have an inverse? We can check that easily. First, lets make a list of all 10 M-numbers, (0) through (9):
In [5]:
elems = [M(n) for n in range(10)]
elems
Out[5]:
Now we can do like a "times table" wherein we pick a single M-number, say M(5), and multiply it by every number in elems...
In [6]:
[M(5) * x for x in elems]
Out[6]:
Interesting. In ordinary arithematic, the times table for 5 goes 0, 5, 10, 15, 20, 25, 30... and so on. Factoring out the 10s, leaving only remainders, we get M(0) or M(5) as our two M-numbers. We have no way to reach M(1) and so M(5) has no inverse. We don't have a group yet.
Lets try M(2):
In [7]:
[M(2) * x for x in elems]
Out[7]:
Same thing. We cycle around and around, always stopping at the same stations, like a model train going in a circle. We never stop at station M(1). So M(2) has no inverse either. What about M(3)?
In [8]:
[M(3) * x for x in elems]
Out[8]:
Aha! Now we're getting somewhere. M(3) * M(7) returns M(1), so these two are inverses of one another. One fact to notice immediately is neither has any factors in common with 10. In fact, both are prime numbers in the ordinary integer sense. M(9) is not prime, but again, 9 has no factors in common with 10. So does M(9) have an inverse? Lets check:
In [9]:
[M(9) * x for x in elems]
Out[9]:
Indeed it does. M(9) * M(9) = M(1), meaning M(9) is its own inverse.
When positive integers have no factors in common, aside from 1, we say they're "relatively prime" or that they're "strangers" to one another. 3 and 10 are strangers, as are 9 and 10. Sometimes we write 7 | 10 = 1, meaning their greatest factor in common is 1. On the other hand, 6 | 10 = 2, as 2 divides into both without remainder.
What we're about to discover is the strangers to a modulus comprise a group, i.e. M(n) where n | N = 1 and M.modulus = N. A Cartesian product of such "minions" would demonstrate closure, inverse for everyone, naturally a neutral element, and associativity i.e. makes no difference if you go M(a) M(b) M(c) by starting with either M(a) M(b) or M(b) M(c). Commutativity is not a requirement for grouphood.
It'd sure be handy at this point, to have a function, gcd, that returns the greatest common divisor of two numbers. Then we could find all the strangers to 10, that are positive integers less than 10, pretty easily. They would have no common divisor with 10.
We call these strangers the "totatives" of 10, and the number of totatives, is called the "totient" of 10.
Note that we're using list comprehension syntax in this example.
The for loop inside the square brackets makes x be every integer from 0 to modulus - 1, but then most of these get filtered out, because of factors in common with the modulus.
In [10]:
import math
modulus = 10
[x for x in range(modulus) if math.gcd(x, modulus) == 1]
Out[10]:
In [11]:
modulus = 12
[x for x in range(modulus) if math.gcd(x, modulus) == 1]
Out[11]:
That was pretty easy. Right away we're getting totatives. We could be using set comprehension syntax instead.
Using Python's "list comprehension" syntax, which allows an if clause for filtering, we were able to get the totatives of 10 and 12 respectively. Both 10 and 12 have a totient of four, meaning each has four totatives.
What's true is that the M-numbers that are totatives of some modulus (say 10), collectively form a group under multiplication.
If we confine ourselves to strangers to the modulus, we'll have an inverse for every element, closure, associativity, and M(1) will be included. That's a quick way to get a group.
We don't have closure for addition though. If we want to include addition, along with its neutral element zero, then we'll need to make our modulus a prime number. Then we will be guaranteed a field, a Galois Field to be more precise, though in German they say it differently.
Before we jump to fields though, lets build a function for computing totatives and use that to build some groups.
A "set comprehension" is just like a list comprehension in terms of syntax, but the curly braces make it for a set, instead of a list.
Sets are unordered and their elements are always unique (no duplicates). Since we know totatives are unique, we might as well practice using a set object.
In [12]:
def totatives(n):
return {x for x in range(n) if math.gcd(x, n) == 1}
In [13]:
elems = {M(x, 12) for x in totatives(12)}
elems
Out[13]:
What we just did there was compose a group of M-number objects, of twelves four strangers. If you're using this Notebook live, here might be a good time to insert some code cells of your own, testing out whether the Cayley Table, i.e. the everything by everything multiplication table, really shows closure, and 1 in every row (proof of inverse element).
In [14]:
M.modulus = 100
elems = {M(x, 100) for x in totatives(100)} # set comprehension
elems # print M-number totatives in no special order -- it's a set
Out[14]:
Note that set objects cannot be ordered even if we might like them to be as by definition they're not sequences. Does Python have an OrderedSet corresponding to an OrderedDict?
If we confine ourselves to a set of strangers, we're safe in assuming there's always an inverse for any given element. That suggests a "brute force" way of finding any element's inverse: just multiply it by every element in the same set of strangers, until their product is 1 (the identity element). The function below accomplishes this:
In [15]:
def inverse(n : M, elems : set) -> M:
for m in elems:
if (n * m).val == 1:
return m
In [16]:
inverse(M(79), elems)
Out[16]:
In [17]:
M(19) * M(79)
Out[17]:
In [18]:
def totient(n):
"""how many totatives have we?"""
return len(totatives(n))
totient(1_000_000) # using _ helps us parse the number, one million
Out[18]:
The number one million has a totient of 400,000, meaning that many numbers from one to one million are co-prime to it.
We now have the technology (tool set) we need to add the __pow__ method to our M class. A negative power triggers finding an inverse of a, the -1 part, then we switch back to positively powering. We're saying $M(n)^{-e}$ equals $(M(n)^{-1})^{e}$.
This is how we treat exponents normally, e.g. 2 ** -2 equals (2 ** -1)**2 or 1/4.
In [19]:
2 ** -2 == (2 ** -1) ** 2
Out[19]:
For example, pow(M(79), -2) means find the inverse of M(79) e.g. M(79)**-1, and then raise the result to the 2nd power.
In [20]:
class M: # for "modulo"
"""
Version 2, with inverting and powering added
Version 3, not listed, is in groups.py and uses xgcd (next Chapter)
"""
modulus = 10 # class level
def __init__(self, val, modulus=None):
if modulus:
self.modulus = M.modulus = modulus # resettable
else:
self.modulus = M.modulus
self.val = val % M.modulus
def __add__(self, other):
if self.modulus != other.modulus:
raise ValueError
return M(self.val + other.val, self.modulus)
def __mul__(self, other):
if self.modulus != other.modulus:
raise ValueError
return M(self.val * other.val, self.modulus)
def __invert__(self):
elems = {M(x) for x in totatives(self.modulus)}
for m in elems:
if (self * m).val == 1:
return m
def __pow__(self, exp): # pow() and ** both trigger this method
output = self
if exp < 0:
exp = abs(exp) # make exp positive now
output = ~self # -1 means take the inverse
elif exp == 0:
output = M(1) # return identity if exp is 0
elif exp == 1:
output = self # return M if exp is 1
if exp > 1:
for _ in range(1, exp): # use __mul__ (already defined)
output *= self
return output
def __repr__(self):
return "(" + str(self.val) + " mod " + str(M.modulus)+ ")"
In [21]:
M.modulus = 100
a = M(79)
a_inv = pow(a, -1)
a_inv
Out[21]:
In [22]:
M(13) * M(13) * M(13)
Out[22]:
In [23]:
M(13)**3 # confirm this is another way of saying it (we're testing)
Out[23]:
In [24]:
M(13)**-1 # give me the inverse of M(13) please
Out[24]:
In [25]:
M(13) * M(77) # these must be inverses then
Out[25]:
Now that we have our inverse function, we can ask questions such as:
In [26]:
print("What is the inverse of M(5, 7)?", ~M(5,7))
print("What is the inverse of M(7, 5)?", ~M(7,5)) # M(7,5) is same as M(2,5)
print("What is the inverse of M(97, 100)?", ~M(97, 100))
Now that we have a powering method in place, lets take a look a Fermat's Little Theorem. The logic here is useful, in the sense that we sometimes get tripped up by IF / THEN statements.
IF a number p is prime, THEN M(A, p) ** p == A, where 0 < A < p. If some number n passes this test for a given A, where gcd(A, n) == 1, then that's evidence n might be prime, but not proof. If n in fact turns out to be composite, then a is called a Fermat Liar or Fermat Fooler.
Lets take a look. Is 511 prime? Lets find an A with no factors in common with 511.
In [27]:
from math import gcd
gcd(511, 25)
Out[27]:
In [28]:
n = M(25, 511) # n will power modulo 511
n ** 511
Out[28]:
Here we say 25 is a "witness" to the fact that 511 is not prime. If it were, it would pass this Fermat test. How about 513?
In [29]:
n = M(25, 513)
n ** 513
Out[29]:
Nope, that's not passing either. Remember we expect A back out again, of 25 in this case. 25 is a witness that 513 is not prime. How about 523 then?
In [30]:
n = M(19, 523)
n ** 523
Out[30]:
In [31]:
n = M(25, 523)
n ** 523
Out[31]:
In [32]:
n = M(31, 523)
n ** 523
Out[32]:
OK, at last it appears we've struck gold. All of our As are testifying that 523 is a prime. Could all of these be liars?
Lets introduce the set of Carmichael Numbers, precisely those composites which sneak through the Fermat Test for any suitable base 1 < A < n.
561 is the lowest Carmichael number.
In [33]:
liars = [M(n, 561) for n in (19, 25, 31) if gcd(n, 561) == 1]
liars
Out[33]:
In [34]:
[n**561 for n in liars] # all pass the Fermat primality test
Out[34]:
And yet 561 is a composite number, with prime factors (3, 11, 17).
Another test for primality, called the AKS Test, is more fool-proof. We'll need the coefficients of the polynomial (x - 1) ** p, where p is our candidate prime. Those coefficients may be read of Pascal's Triangle. If all coefficients but the first and last (which are 1) are divisible by p, then p is prime.
Lets try that out, first by writing a Python generator for successive rows of Pascal's Triangle:
In [35]:
def pascal():
row = [1]
while True:
yield row
row = [(i+j) for i,j in zip(row + [0], [0] + row)]
gen = pascal()
for i in range(10):
print(next(gen))
That seems to be working. So Pascal's Triangle could be at the heart of a new iterator for producing successive prime numbers.
In [36]:
def primes():
yield 2
p = 1
gen = pascal()
next(gen) # [1]
next(gen) # [1, 1]
while True:
p += 2 # 3, 5... check odd rows only, once 2 yielded
next(gen) # skip even row
r = set(next(gen)[1:-1]) # drop 1s on both ends and dedupe
if sum([coeff%p for coeff in r]) == 0:
yield p
ps = primes()
print(str([next(ps) for _ in range(100)]))
Homework:
Watch this Youtube on the Chinese Remainder Theorem (CRT). Be prepared to explain its utility. What Youtube do you find most elucidating, regarding the CRT?
So now it looks like we have another fun tool for exploring Group Theory, while learning Python at the same time!
Now lets take what we've learned in this chapter and apply it to an important topic: Public Key Cryptography.
Back to Chapter 3: A First Class
Chapter 5: Public Key Cryptography
Introduction