# W4 Notes: Number Theory 2

[ToC]

## Generating Primes with the Sieve of Eratosthenes

The [*Sieve of Eratosthenes*](https://en.wikipedia.org/wiki/Sieve_of_Eratosthenes) generates all primes up to $n$.

The algorithm first marks all the integers in $[1, n]$ prime, then iteratively takes each prime (ignoring non-primes) and marks all its multiples non-prime.

### Implementation

vector<bool> isprime; // Declare globally so it's accessible anywhere.

void sieve(int n) {
    // Mark all numbers as prime.
    isprime.assign(n + 1, true);
    // Handle trivial cases: 0 and 1 are not prime.
    isprime[0] = isprime[1] = false;
    // Loop over all candidates.
    for (int i = 2; i <= n; ++i) {
        // If `i` is already marked non-prime, continue with next candidate.
        if (!isprime[i]) continue;
        // If `i` is still marked prime, it is indeed a prime.
        // Loop over all multiples `j` of `i` s.t. `j > i` and mark them non-prime.
        for (int j = (i + i); j <= n; j += i) {
            isprime[j] = false;
        }
    }
}

Time complexity: $O(n \log{\log{n}})$

:::info
**Note:** On the surface, it may seem that the complexity should be $O(n \log{n})$. But notice that we are entering the inner loop only when we find a prime at the outer loop. Hence, the actual complexity comes from the distribution of primes. Notice that $O(n \log{\log{n}})$ is only **better** than $O(n \log{n})$, not worse. In fact, it's almost as good as $O(n)$.
:::

### Better implementation

The above implementation can be further optimized with a few changes:

1. Mark all even integers greater than $2$ (i.e., the multiples of $2$), non-prime, and never iterate over them again in the main loop.
1. **[important]** In the main loop, iterate with $i$ only up to $\sqrt n$ (i.e., run the loop only while $i^2 \leq n$ is true).
:::spoiler Think: Why does it work? *(Click for spoiler)*
If $m$ is a composite number, then its smallest prime factor is less than or equal to $\sqrt m$. Remember how we could generate all factors of an integer by iterating up to its square-root?
:::
1. **[important]** In the inner loop, start iterating with $j$ from $i^2$ instead of $2i$.
:::spoiler Think: Why does it work? (*Click for spoiler)*
If $p$ is a prime, then $p^2$ is the smallest integer such that **only** $p$ can factorize it. So, by the time we have found a prime $i$ in the main loop, all composite numbers less than $i^2$ will already have been marked non-prime by at least one prime smaller than $i$.
:::

All these changes don't improve the time complexity, but they reduce the actual number of instructions noticeably.

vector<bool> isprime;

void sieve(int n) {
    isprime.assign(n + 1, true);
    isprime[0] = isprime[1] = false;
    for (int j = 4; j <= n; j += 2) {
        isprime[j] = false;
    }
    for (int i = 3; (i * i) <= n; i += 2) {
        if (!isprime[i]) continue;
        for (int j = (i * i); j <= n; j += (i + i)) {
            isprime[j] = false;
        }
    }
}

:::info
**Do it yourself:** Instead of boolean primality, sieve to generate the *smallest* prime factor of each natural number up to a certain number.
:::

### Storing primes

If required, we can even store the sieved primes in any container, e.g., a vector.

vector<bool> isprime;
vector<int> primes;

void sieve(int n) {
    // . . .
    // insert code for sieve here
    // . . .
    primes.push_back(2);
    for (int i = 3; i <= n; i += 2) {
        if (isprime[i]) {
            primes.push_back(i);
        }
    }
}

## Modified Sieve for Number and Sum of Divisors

We can apply the technique used in Sieve of Eratotosthenes to find the number and sum of divisors of all natural numbers up to $n$.

### Sieve for NoD

vector<int> nod;

void sieve_nod(int n) {
    nod.assign(n + 1, 0);
    for (int i = 1; i <= n; ++i) {
        for (int j = i; j <= n; j += i) {
            ++nod[j];
        }
    }
}

The complexity: $O(n\log{n})$

### Sieve for SoD

vector<int> sod;

void sieve_sod(int n) {
    sod.assign(n + 1, 0);
    for (int i = 1; i <= n; ++i) {
        for (int j = i; j <= n; j += i) {
            sod[j] += i;
        }
    }
} }
}

The complexity: $O(n\log{n})$

## Euclidean Algorithm: GCD and LCM

The [Euclidean algorithm](https://en.wikipedia.org/wiki/Euclidean_algorithm) computes the $\gcd$ ([greatest common divisor](https://en.wikipedia.org/wiki/Greatest_common_divisor)) of two integers.

Once we have computed the $\gcd$ of two integers, their $\lcm$ ([least common multiple](https://en.wikipedia.org/wiki/Least_common_multiple)) can then be computed using this relation:

$$
\DeclareMathOperator*{\lcm}{lcm}
a \cdot b = \gcd{(a,b)} \cdot \lcm{(a,b)}\\
\Rightarrow \lcm{(a,b)} = \frac{a \cdot b}{\gcd{(a,b)}}
$$

### GCD Implementation

int gcd(int a, int b) {
    // base case
    if (b == 0) return a;
    // general case
    return gcd(b, a % b);
}

Time complexity: $O(\log \min{(a, b)})$

:::info
**Note:** Mathematically speaking, it is not a necessity for $a$ and $b$ to be non-negative. But we know that the `%` operator does not give a positive value for negative numbers, even though the mathematical definition of remainder requires it to be non-negative. If we need to find $\gcd(a, b)$ where either of $a$ and $b$ is negative, we can simply pass its absolute value as the argument according to this relation:
$\begin{align*}
&\gcd(a, b) \\
= &\gcd(-a, b) \\
= &\gcd(a, -b) \\
= &\gcd(-a, -b)
\end{align*}$
:::

### LCM Implementation

int lcm(int a, int b) {
    return (a / gcd(a, b)) * b; // divide first to avoid overflow 
}