Number Theory

20pts

Python

Java

C++

Algorithms

Math

An brief introduction to modular arithmetic, prime numbers, and factoring in computer science

Modular Arithmetic

Let’s say that you tried to compute the 100,000th, or 10⁵, Fibonacci number. Soon you realize that the number is completely off, however.

This is because the number in question has 63 digits, a lot more than even the widest bit integer in C++ (128) can support. Therefore, many competitive programming problems will instead ask you for the result **mod** some number, usually either 10⁹ + 7 or 998353244. A number, such as x, **mod** another number, such as y, simply means the remainder of x divided by y. For example, 5 **mod** 2 would be 1. Most programming languages use the % operator to denote mod.

We use **modular arithmetic** to deal with such numbers and find the result in mod efficiently. Instead of doing calculations with extremely large integers directly, we do calculations with their remainder, or mod some number, according to the following rules:

(a + b) mod c = ((a mod c) + (b mod c)) mod c

(a - b) mod c = ((a mod c) - (b mod c)) mod c

(a * b) mod c = ((a mod c) * (b mod c)) mod c

(a - b) mod c = ((a mod c) - (b mod c)) mod c

(a * b) mod c = ((a mod c) * (b mod c)) mod c

(Note that the pattern does __not__ apply to division)

With negative numbers, the modulo operation may give a negative output, resulting in inconsistent behavior in the code. To combat this, we can do ((a mod c) + c) mod c to ensure the result is positive.

Modular Inverse

A modular inverse of a number can be thought of as its **reciprocal**. More formally, define a modular inverse of **x** as **x**^{⁻¹} such that x * x^{⁻¹} **≡** 1 mod c for some divisor c. If c is prime (which it usually is), then it can be proven that x^{⁻¹} always exists for any x between 1 and c - 1, inclusive. Finding the modular inverse of a number can be done by finding x^{c - 2} and using binary exponentiation in the process. This can only be done if c is prime, but the vast majority of the time it will be. This process is also encapsulated in the templates following with the / operator.

As the operator suggests, modular inverses are useful for division calculations. For example, if we wanted to compute 3 divided by 5 modulo 11, we’d calculate 3 * (modular inverse of 5), which would be 3 * 9 = 27 ≡ 5 mod 11.

Modnum Templates

Modnum templates, such as Ecnerwala’s, can make your code much more readable and less prone to bugs compared to coding add, multiply functions etc. yourself. If that’s not an option, however, for contests such as USACO which don’t allow templates, we recommend trying to memorize or understand part of Benq’s short template.

Prime Numbers

Prime numbers are the building blocks of all integers. The fundamental theorem of arithmetic tells us that every number has a unique prime factorization.

Check if one number is prime

To check that n is prime, we iterate from 2 to √(n), checking if they divide n. If none of them divide n, n is prime, otherwise, it’s composite. This takes O(√(n)) time.

Code:

```
def is_prime(x):
if x <= 1:
return False
i = 2
while i * i <= x:
if x % i == 0:
return False
i+=1
return True
```

Code:

```
static boolean isPrime(int x) {
if (x <= 1) {
return false;
}
for (int i = 2; i * i <= x; i++) {
if (x % i == 0) {
return false;
}
}
return true;
}
```

Code:

```
bool is_prime(int x) {
if (x <= 1) {
return false;
}
for (int i = 2; i * i <= x; i++) {
if (x % i == 0) {
return false;
}
}
return true;
}
```

Finding primes from 1 through n

Say we are trying to find all prime numbers from 1 to n. We keep an boolean array prime[n]. Let prime[i]=1 indicate that i is a prime and prime[i]=0 indicate that i is composite. We iterate from 2 to n. If our current number i is marked as prime, we go through the array and mark all the other multiples of i as composite. This way, after iterating through all the numbers, all the primes will have value 1 and all composites will have value 0 (we ignore 1, since it is neither prime nor composite). The algorithm takes O(n log n) time, and is known as the Sieve of Eratosthenes.

Code:

```
prime = [True for i in range(N + 1)]
prime[0] = False
prime[1] = False
for i in range(2, N + 1):
if not prime[i]:
continue
for j in range(i * 2, N + 1, i):
prime[j] = False
```

Code:

```
boolean[] prime = new boolean[N + 1];
prime[0] = false;
prime[1] = false;
for (int i = 2; i <= N; i++) {
if (!prime[i]) {
continue;
}
for (int j = i * 2; j <= N; j++) {
prime[j] = false;
}
}
```

Code:

```
vector<bool> prime(N + 1);
prime[0] = false;
prime[1] = false;
for (int i = 2; i <= N; i++) {
if (!prime[i]) {
continue;
}
for (int j = i * 2; j <= N; j++) {
prime[j] = false;
}
}
```

Factors

Quick Prime Factorization

Suppose we want to compute the prime factorization for every number between 2…N. Usually this would take around O(√(N)) per query, but we can speed this up immensely by precomputing the highest factor of every number to a runtime of O(log N) per query with O(N log N) preprocessing.

Code:

```
hfac = [1 for i in range(N + 1)]
for i in range(2, N + 1):
if hfac[i] == 1:
for j in range(i, N + 1, i):
hfac[j] = i
# 60 -> [2, 2, 3, 5]
def factors(x):
ret = []
while x != 1:
ret.append(hfac[x])
x = int(x / hfac[x])
ret.reverse()
return ret
```

Code:

```
int[] hfac = new int[N + 1];
Arrays.fill(hfac, 1);
for (int i = 2; i <= N; i++) {
if (hfac[i] == 1) {
for (int j = i; j <= N; j += i) {
hfac[j] = i;
}
}
}
// 60 -> [2, 2, 3, 5]
static List
``` factors(int x) {
List ret = new ArrayList<>();
while (x != 1) {
ret.add(hfac[x]);
x /= hfac[x];
}
Collections.reverse(ret);
return ret;
}

Code:

`vector` hfac(N + 1, 1);
for (int i = 2; i <= N; i++) {
if (hfac[i] == 1) {
for (int j = i; j <= N; j += i) {
hfac[j] = i;
}
}
}
// 60 -> [2, 2, 3, 5]
vector factors(int x) {
vector ret;
while (x != 1) {
ret.push_back(hfac[x]);
x /= hfac[x];
}
reverse(ret.begin(), ret.end());
return ret;
}

You can play around with all the code we've used in this article on Replit:

Join our discord server for updates, announcements, future event planning, to ask questions, or just for fun!

Copyright ©2021 Howard County Hour of Code

Copyright ©2021 Howard County Hour of Code