A voyage through Numbers

Exploring Number Theory for Programming

Ashish Ranjan
Programming and Algorithms, IITR
9 min readDec 18, 2017

--

The best way to describe Mathematics is probably “an intersection of art and science”. Through ancient times, artists and philosophers have found inspiration in mathematics and vice-versa. Mathematics is also commonplace in Computer Science research, as well as being heavily applied to graph algorithms and areas of computer vision.

Problems in competitive programming usually require a good knowledge of mathematics, especially Number Theory, Geometry, and Combinatorics. Here, we will focus on some of the basic concepts of Number Theory. Number Theory is sometimes called “The Queen of Mathematics”. It deals with operations on the set of Natural Numbers, sometimes involving Integers. It is pivotal for Information Security and Cryptography. In fact, a crucial aspect of cryptography is the inability of modern computers to factor large numbers quickly.

We’ll only be discussing concepts useful in competitive programming contests here. You don’t need any special knowledge of mathematics to understand these concepts, only a knack for solving problems is enough. But if you’re interested in learning more, you can read ‘Discrete Mathematics and its Applications’ and ‘Elementary Number Theory’ by Kenneth H. Rosen or ‘Number Theory for Computing’ by S.Y. Yan. The latter also discusses the application of Number Theory in cryptography in detail.

XKCD’s Number Line

Modulo Arithmetic

The modulus or modulo operator(denoted by ‘%’ in almost all programming languages) gives the remainder on dividing one number by another. E.g. 10%3 = 1, 74%7 = 4 etc.

Properties of the modulus operator :

Consider calculating (a^b)%m. It can be easily calculated using the following naïve approach.

Naive algo for exponentiation

This approach is very inefficient as it has an O(b) time complexity but we may even have to solve it for 10⁶<b≤10¹⁸.

Getting complex!

In such a scenario, we can use the highly efficient recursive function :

This approach is a major advancement from the naïve approach as it reduces the time complexity to O(log(b)).

Prime Numbers

Prime numbers are, in all likelihood, the most essential component of Number Theory. A prime number is a positive integer having exactly two factors : 1 and itself. Numbers which are not prime are called composite numbers.

Checking if a number is prime :

Given a number N we have to determine whether it is prime or not. This can be easily determined by checking if the number is divisible by any number greater than 1 and less than N.

The time complexity for this will be O(N) but we know that if d|N then (N/d) also divides N. So for any composite number there will be at least one divisor d such that 1≤d≤sqrt(N). So we only need to check divisibility of N with numbers less than or equal to sqrt(N).

So the time complexity of primality test is reduced from O(N) to O(√N).

Primality Testing Algorithm

You can try to solve this problem for practice.

Finding all prime numbers less than N: Sieve of Eratosthenes

To find all primes less than N, we can check the primality of every number but that would have a time complexity of O(N √N). Moreover, it would be redundant because we know that if a number is prime none of its multiples can be prime as it would be one of the divisors.

This fact is used in the famous prime generating algorithm called the Sieve of Eratosthenes :

Implementation of Sieve of Eratosthenes
C++ implementation

The time complexity of this sieve algorithm is O(N loglogN) which can be proved by using bounds for harmonic sums.

SMPF (Smallest Prime Factor)

When we require to factorize small numbers (≤10⁷) many number of times, then we can factorize a number in O(logN) time complexity after doing a sieve precomputation of O(N loglogN).

We can do so by doing a sieve such that sieve[n] stores the smallest prime factor of the number n for all n ≤ N.

The SMPF sieve

Proof of complexity:

The sieve part is O(N loglogN) since it is very similar to the prime sieve.

For factorization, at each step, we divide n by sieve[n]. Since sieve[n]≥2, n is reduced to at least n/2 at every step. Hence, n will become 1 in at most O(logN) divisions.

You can try to solve this problem for practice.

Finding primes in a range: The Segmented Sieve

The sieve of Eratosthenes is indeed a good approach to find the primes. However, if N is a very large number, its space complexity will be O(N) which may not fit in a cache memory.

Using a segmented sieve, we can do the same with an O(√N) space, the time complexity remaining the same.

In a segmented sieve, we use the simple sieve for calculating primes up to √N.

For the other remaining elements, the order in which we mark off the primes is changed. Instead of marking all multiples of a prime p in range (p, N], we proceed by marking the multiples of p for data that will fit in our memory and then continue for the next set of data.

In simple words, say if the limit of data that can fit in memory is B. Then we will apply simple sieve for all elements up to √N and then apply segmented sieve. We will find primes in range [sqrt(N), sqrt(N)+B) , then [sqrt(N)+B, sqrt(N)+2*B) and so on.

The segmented sieve also comes in handy when we need to find primes in the range [L,R] and R is large enough but we can create an array of size [R-L+1].

You can try to solve this problem for practice.

The problem PRIME1 can also be solved using a segmented sieve to get a better running time.

Finding the sum of divisors of all numbers up to N

A naive approach would be to find all the divisors of i (1≤ i ≤ N) which will take O(sqrt(i)) time for each i. The loose upper bound for this will be O(N).

Better approach :

Let Sum[i] denote the sum of all divisors of i. Sum[i] is initialised to 0 for all i.

Then for all j from 1 to n, we check for the multiples of j and add j to Sum[i] ,i.e for all multiples of j, increment their sum of divisors, Sum by j.

For each j , the inner loop iterates n/j times . The total time required is (n/1 + n/2 + n/3 + … + n/n) = n(1 + ½+ ⅓ + …+ 1/n ). The time complexity for above code is O(nloglogn).

You can try solving problems DIVSUM, NDIV and NFACTOR on SPOJ for practice.

Greatest Common Divisor (GCD)

Greatest Common Divisor or GCD of two numbers is defined as the largest positive number that divides both the numbers.

gcd(a,b) = {max(d) : d|a,d|b}

A brute-force approach to finding the GCD of two numbers would be iterating from the lower number to 1. The first number found to divide both numbers would be the GCD.

A brute force approach for GCD calculation

The time complexity of this function is O(min(a,b)). This can be improved by using what is known as Euclid’s Algorithm for GCD which states GCD (A,B) = GCD(B,A%B).

Say, for A= 2177 and B= 147

GCD (2177,147) = GCD (147,2177%147) = GCD(147,119)

= GCD (119,147%119) = GCD(119,28)

=GCD (28,119%28) = GCD(28,7)

=GCD( 7,28%7) = GCD(7,0) = 7

GCD of a positive integer and 0 is the positive integer itself.

So GCD is computed as :

Euclid’s GCD Algorithm

You can try to solve this problem for practice.

Extended Euclid Algorithm

Equations of the form Ax + By = C have an integer solution if and only if C is a multiple of GCD(A,B). This is easy to prove as the LHS is always divisible by GCD(A,B).

So, effectively if we find a solution for Ax + By = G [where G = GCD(A,B)] say (p,q) the solution for Ax+By = C will be (C.p/G,C.q/G).

Moreover, GCD(A,B) has a special property (Bezout’s Identity) that it can always be represented in the form of an equation, i.e., Ax + By = GCD(A, B). The Extended Euclidean Algorithm uses this fact to find integer solutions for the equation Ax+By = GCD (A,B). This algorithm gives us the coefficients (x and y) of this equation which will be later useful in finding the Modular Multiplicative Inverse. These coefficients can be zero or negative in their value.

Implementation:

Euclid’s Extended GCD Algorithm

If x and y (any solution) are known then general solution (X,Y) can be calculated as follows :

X = x + k*(b/GCD(a,b))

Y = y — k*(a/GCD(a,b))

The Time Complexity of Extended Euclid’s Algorithm is O(log(max(A,B))).

You can try to solve this problem for practice.

Euler’s Totient Function

Euler’s Totient Function (represented by phi(n) or 𝛷(n)) is a very important number-theoretic function with a lot of useful properties. It is the count of numbers less than or equal to n and coprime to n. More formally

𝛷(n) =Count{ i : 1<=i<=n and gcd(i,n) =1}

Plot of 𝛷(n) vs n

The function has the following properties :

Properties of the Euler Totient Function

The SPOJ Problem ETF requires us to calculate 𝛷(n) for n in the range [1,10⁶] for T test cases.

A naïve approach to this problem can be:

But under the given constraints it will give TLE as the complexity would be O(T*NlogN).

Instead, using the properties we can calculate 𝛷(n) with a factorization algorithm :

This solves the problem in O(T*sqrt(N)) time.

Can we do better? Yes, we can. We can compute the Euler totient function of all numbers from 1 to N in O(NlogN) complexity and then answer each test case using the precalculated value in O(1). This reduces the overall complexity of the problem to O(N logN + T).

Similar to the segmented sieve for primes, we can also design a segmented sieve for the Euler totient function. This is left as an exercise for the reader.

You can try solving problems ETFS and ETFD on SPOJ for practice.

Fermat’s Little Theorem

This theorem states that for a prime p and a natural number a coprime to p

a^(p-1) ≡ 1 (mod p)

Euler’s Theorem

This is the generalization of Fermat’s Little Theorem. It states that for coprime numbers a and n

a^𝛷(n) ≡ 1 (mod n)

Fibonacci Numbers

Mathematically the nth Fibonacci number is defined as

F[0]=0, F[1]=1, F[n] = F[n-1]+F[n-2] n≥2

The sequence looks like this:

{0, 1, 1, 2, 3, 5, 8, 13, 21…}

Fibonacci Grid

Some Properties

Properties of Fibonacci Numbers

To compute all Fibonacci numbers less than a given number N

This can be done with a time complexity O(N) using the recursive definition.

To compute a single Fibonacci Number, the Nth in series

We can find F[N] in O(log(N)) this way:

Matrix exponentiation for Fibonacci numbers

This formula can be proven using induction.

Hence, to solve our original problem we can find X^N, for which the pseudo code is same as exponentiation of any number(which is O(log(N))) except the values are going to be matrices for this problem and the operations addition and multiplication on matrices need to be defined as functions with input and output as 2x2 arrays.

You can try solving this for practice.

More problems for the reader’s practice: MAIN74, INVPHI, DCEPC200, GGD, POWFIB, Counting fractions, Totient maximum, Marbles, Rare Easy Problem

INS16G, INS16H, INS17M are problems set by members of Programming and Algorithms Group, IITR that you can practice. These problems were set for our annual flagship programming contest Insomnia. You can find all the problems set by the members here.

If you are new to competitive programming and need help setting up an editor here is a comprehensive guide by Ketan Gupta, one of our members, to set up Sublime Text for Competitive Programming.

Programming and Algorithms Group, commonly known as PAG, is a student-run group that fosters competitive programming under the Software Development Section at IIT, Roorkee. We have a forum open to all doubts pertaining to competitive programming for everyone.

--

--