## Saturday, November 18, 2017

### Chinese Remainder Theorem Part 2 - Non Coprime Moduli

As promised on the last post, today we are going to discuss the "Strong Form" of Chinese Remainder Theorem, i.e, what do we do when the moduli in the congruence equations are not pairwise coprime. The solution is quite similar to the one we have already discussed in the previous post, so hopefully, it will be a lot easier to understand.

#### Prerequisite

If you haven't read my previous post, you can find it here: Chinese Remainder Theorem Part 1. I strongly suggest you read that one before proceeding onwards.

## Problem

Given two sequences of numbers $A = [a_1, a_2, \dots, a_n]$ and $M = [m_1, m_2, \dots, m_n]$, find the smallest solution to the following linear congruence equations if it exists:

$$x \equiv a_1 \text{(mod } m_1)\\ x \equiv a_2 \text{(mod } m_2)\\ \dots \\ x \equiv a_n \text{(mod } m_n)$$

Note: The elements of $M$ are not pairwise coprime

## Working With Two Equations

Before we go on to deal with $n$ equations, let us first see how to deal with two equations. Just like before, we will try to merge two equations into one.

So given the two equation $x \equiv a_1 \text{(mod }m_1)$ and $x \equiv a_2 \text{(mod }m_2)$, find the smallest solution to $x$ if it exists.

### Existance of Solution

Suppose, $gcd(m_1,m_2) = g$. In that case, the following constraint must be satisfied: $a_1 \equiv a_2 \text{(mod }g)$, otherwise there does not exist any solution to $x$. Why is that?

We know that, $x \equiv a_1 \text{(mod }m_1)$
or, $x - a_1 \equiv 0 \text{(mod }m_1)$
or, $m_1 | x - a_1$
Since $m_1$ divides $x - a_1$, any divisor of $m_1$ also divides $x - a_1$, including $g$.
$\therefore x - a_1 \equiv 0 \text{(mod }g)$.
Similarly, we can show that, $x - a_2 \equiv 0 \text{(mod }g)$ is also true.
$\therefore x - a_1 \equiv x - a_2 \text{(mod }g)$
or, $a_1 \equiv a_2 \text{(mod }g)$
or, $a_1 - a_2 \equiv \text{(mod }g)$
or $g | (a_1 - a_2)$

Therefore, $g$ must divide $a_1 - a_2$, otherwise, there is conflict and no solution exists. From this point, we assume that $g | a_1 - a_2$.

### Finding Solution

Just like last time, we use Bezout's identity to find a solution to $x$. From Bezout's identity, we know that there exists a solution to the following equations:

$m_1p + m_2q = g$
Since both side is divisible by $g$, we can rewrite it as:
$$\frac{m_1}{g}p + \frac{m_2}{g}q = 1$$

Using Extended Euclidean Algorithm, we can easily find values of $p$ and $q$. Simply call the function ext_gcd(m1/g, m2/g, p, q) and it should be calculated automatically.

Once we know the value of $p$ and $q$, we can claim that the solution to $x$ is the following:

$$x = a_1\frac{m_2}{g}q + a_2\frac{m_1}{g}p$$

Why is that? Look below.

$$x = a_1\frac{m_2}{g}q + a_2\frac{m_1}{g}p \\ \text{From Bezout's Identity, we know that } \frac{m_1}{g}p = 1 - \frac{m_2}{g}q \\ \text{so, } x = a_1\frac{m_2}{g}q + a_2(1 - \frac{m_2}{g}q) \\ \text{or, } x = a_2 + (a_1 - a_2)\frac{m_2}{g}q \\ \text{Since, } g | (a_1 - a_2) \\ x = a_2 + \frac{a_1 - a_2}{g}m_2q \\ \therefore x \equiv a_2 \text{(mod }m_2) \\ \text{Similarly, } x \equiv a_1 \text{(mod }m_1)$$

Hence, $x = a_1\frac{m_1}{g}q + a_2\frac{m_2}{g}p$ is a valid solution. It may not be the smallest solution though.

### Finding Smallest Solution

Let $x_1$ and $x_2$ be adjacent two solutions. Then the following are true:

$x_1 \equiv a_1 \text{(mod }m_1)$
$x_2 \equiv a_1 \text{(mod }m_1)$
$x_1 \equiv x_2 \text{(mod }m_1)$
$x_1 - x_2 \equiv 0 \text{(mod }m_1)$
$m_1 | x_1 - x_2$
Similarly, $m_2 | x_1 - x_2$
Since both $m_1$ and $m_2$ divides $x_1 - x_2$, we can say that $LCM(m_1, m_2)$ also divides $x_1-x_2$.

Since any two adjacent solution of $x$ is divisible by $L = LCM(m_1,m_2)$, then there cannot be more than one solution that lies on the range $0$ to $L-1$. Hence, the following is the smallest solution to $x$:

$$x \equiv a_1\frac{m_2}{g}q + a_2\frac{m_1}{g}p \text{(mod }LCM(m_1,m_2))$$

## Working with $n$ Equations

When there are $n$ equations, we simply merge two equations at a time until there is only one left. The last remaining equation is our desired answer. If at any point, if we are unable to merge two equations due to conflict, i.e, $a_i \not\equiv a_j \text{ (mod }GCD(m_i, m_j))$, then there is no solution.

## Chinese Remainder Theorem: Strong Form

Given two sequences of numbers $A = [a_1, a_2, \dots, a_n]$ and $M = [m_1, m_2, \dots, m_n]$, a solution to $x$ exists for the following $n$ congrunce equations:

$$x \equiv a_1 \text{(mod } m_1)\\ x \equiv a_2 \text{(mod } m_2)\\ \dots \\ x \equiv a_n \text{(mod } m_n)$$

if, $a_i \equiv a_j \text{ (mod }GCD(m_i,m_j))$ and the solution will be unique modulo $L = LCM(m_1, m_2,\dots, m_n)$.

## Code

The code is almost same as weak form, with slight changes in few lines.

/** Works for non-coprime moduli.
Returns {-1,-1} if solution does not exist or input is invalid.
Otherwise, returns {x,L}, where x is the solution unique to mod L
*/
pair<int, int> chinese_remainder_theorem( vector<int> A, vector<int> M ) {
if(A.size() != M.size()) return {-1,-1}; /** Invalid input*/

int n = A.size();

int a1 = A[0];
int m1 = M[0];
/** Initially x = a_0 (mod m_0)*/

/** Merge the solution with remaining equations */
for ( int i = 1; i < n; i++ ) {
int a2 = A[i];
int m2 = M[i];

int g = __gcd(m1, m2);
if ( a1 % g != a2 % g ) return {-1,-1}; /** No solution exists*/

/** Merge the two equations*/
int p, q;
ext_gcd(m1/g, m2/g, &p, &q);

int mod = m1 / g * m2; /** LCM of m1 and m2*/

/** We need to be careful about overflow, but I did not bother about overflow here to keep the code simple.*/
int x = (a1*(m2/g)*q + a2*(m1/g)*p) % mod;

/** Merged equation*/
a1 = x;
if (a1 < 0) a1 += mod; /** Result is not suppose to be negative*/
m1 = mod;
}
return {a1, m1};
}


Once again, please note that if you use int data type, then there is a risk of overflow. In order to keep the code simple I used int, but you obviously need to modify it to suit your needs. You can have a look at my CRT template on github from here: github/forthright48 - chinese_remainder_theorem.cpp

Complexity: $O(n\times log(L))$.

## Conclusion

And we are done. We can now solve congruence equations even when the moduli are not pairwise coprime.

The first time I learned Chinese Remainder Theorem, I had to read a paper. The paper was well written, but I never actually managed to understand how it worked. Therefore, I always feared Chinese Remainder Theorem and used it as a black box. Until last week, I didn't even know that CRT worked with moduli that are not coprime. I really enjoyed learning how it works and all the related proof. Hopefully, you enjoyed it too.

## Resources

1. forthright48 - Chinese Remainder Theorem Part 1 https://forthright48.blogspot.com/2017/11/chinese-remainder-theorem-part-1.html
2. Wiki - Chinese Remainder Theorem - https://en.wikipedia.org/wiki/Chinese_remainder_theorem
3. Brilliant.org - Chinese Remainder Theorem - https://brilliant.org/wiki/chinese-remainder-theorem/
4. Cut The Knot - Chinese Remainder Theorem - https://www.cut-the-knot.org/blue/chinese.shtml

## Wednesday, November 15, 2017

### Chinese Remainder Theorem Part 1 - Coprime Moduli

Second part of the series can be found on: Chinese Remainder Theorem Part 2 - Non Coprime Moduli

Wow. It has been two years since I published my last post. Time sure flies by quickly. I have been thinking about resuming writing again for a while now. Took me long enough to get back into the mood.

I am really excited about today's post. I have been meaning to write on Chinese Remainder Theorem (CRT for short) for over two years now. Hopefully, everyone will find it interesting and easy to understand.

Here are the pre-requisite that you need to complete for understanding the post properly:

## An Old Woman Goes to Market

I remember reading the following story the first time I tried learning about CRT. It's quite interesting in my opinion. Kind of brings out the essence of CRT in a simple to understand fashion.

An old woman goes to market and horse steps on her basket and crashes the eggs. The rider offers to pay the damages and asks her how many eggs she had brought. She does not remember the exact number, but when she had taken them out two at a time, there was one egg left. The same happened when she picked 2,3,4,5,6 at a time, but when she took seven at a time they came out even. What is the smallest number of eggs she could have had?

Hopefully, everyone understood the story above. Let me now rephrase it in mathematical notations.

Suppose, the old woman had $x$ eggs in her basket. Now she claims that when she took out eggs from the basket $2$ at a time, there was only $1$ egg remaining, meaning, $x \equiv 1 \text{(mod 2)}$. So basically, she is giving us various congruence equations. From her statement, we get the following linear congruence equations:

$$x \equiv 1 \text{(mod 2)} \\ x \equiv 1 \text{(mod 3)} \\ x \equiv 1 \text{(mod 4)} \\ x \equiv 1 \text{(mod 5)} \\ x \equiv 1 \text{(mod 6)} \\ x \equiv 2 \text{(mod 7)}$$

Now, we need to find the smallest value of $x$ which satisfies all of the above congruence equations. Chinese Remainder Theorem helps us in finding the value of $x$.

Before we move on, let us redefine the problem formally.

## Problem: Formal Definition

Given two sequences of numbers $A = [a_1, a_2, \dots, a_n]$ and $M = [m_1, m_2, \dots, m_n]$, find the smallest solution to the following linear congruence equations if it exists:

$$x \equiv a_1 \text{(mod } m_1) \\ x \equiv a_2 \text{(mod } m_2) \\ \dots \\ x \equiv a_n \text{(mod } m_n)$$

## Chinese Remainder Theorem: Weak Form

As mentioned before, we can use Chinese Remainder Theorem to solve the above problem described. On this part of the CRT series, we will look into a weaker form of Chinese Remainder Theorem, which is easier to understand and occurs more frequently.

The weak Chinese Remainder Theorem states the following:

Given two sequences of numbers $A = [a_1, a_2, \dots, a_n]$ and $M = [m_1, m_2, \dots, m_n]$, where all elements of $M$ are pairwise coprime, there always exists an unique solution to $x$ mod $L$, where $L = m_1 \times m_2 \times \dots \times m_n$, such that $x$ satisfies the following linear congruence equations:

$$x \equiv a_1 \text{(mod } m_1)\\ x \equiv a_2 \text{(mod } m_2)\\ \dots \\ x \equiv a_n \text{(mod } m_n)$$

So the weak form of Chinese Remainder Theorem has a constraint: members of the array $M$ must be pairwise coprime. What do we mean by that? This means that $\text{GCD}(m_i,m_j) = 1$ when $i \neq j$.

As long as this condition is satisfied, the weak form of CRT state that a solution to $x$ always exists which is "unique mod $L$". Let us see an example to understand what it means to be "unique mod $L$".

### Weak Form of CRT: Example

Suppose we are given the following three congruence equations:

$$x \equiv 3 (\text{mod } 5) \\ x \equiv 2 (\text{mod } 7) \\ x \equiv 2 (\text{mod } 8)$$

Here we have $A=[3,2,2]$, $M=[5,7,8]$ and $L = 5 \times 7 \times 8 = 280$. Using Weak form of Chinese Remainder Theorem, we can find that $x \equiv 58 \text{(mod } 280)$. How did we find this? We will see that later.

Before you continue, please verify yourself that this indeed satisfies the given congruence. Also, $x = 58$ is the smallest solution and there exists more than one solution. $x = 58 + 280 \times k$, where $k$ is any non-negative integer, also satisfies the equations. But as we said before, $x$ is unique when you mod all the solutions with $L$. This is what it means to have a solution unique to mod $L$.

### Weak Form of CRT: Finding a Solution

We are first going to see how to solve when there are just two equations. Once we know how to solve for two equations, we will then generalize the method for $n$ equations.

#### When there are just two equations

Let us first just consider two equation: $x \equiv a_1 \text{(mod }m_1)$ and $x \equiv a_2 \text{(mod }m_2)$. What we are going to do is merge these two equations into one equation. Here is how it works.

Since $gcd(m_1,m_2) = 1$, from Bezout's Identity, we know that there exists a solution to the following equation: $$m_1p + m_2q = 1 \tag{1}\label{1}$$

Using Extended Euclidean Algorithm, we can find the value of $p$ and $q$. If we know the value of $p$ and $q$, then we can say that:

$x = a_1 m_2 q + a_2 m_1 p \text{ (mod }m_1m_2) \tag{2}$

See below to understand why so.

$x = a_1 m_2 q + a_2 m_1 p$
$x = a_1 (1 - m_1 p) + a_2 m_1 p$ (Using $\eqref{1}$)
$x = a_1 - a_1 m_1 p + a_2 m_1 p$
$x = a_1 + (a_2 - a_1) m_1 p$
$\therefore x \equiv a_1 \text{(mod }m_1)$

Similarly, $x = a_1 m_2 q + a_2 m_1 p \equiv a_2 \text{(mod } m_2)$

Great! We now have a possible solution for $x$, but why did we mod the solution with $m_1m_2$?. The solution we found may or may not be the smallest. Like we have seen before, there could be infinite solutions, but there exists a solution which is unique to mod $m_1 \times m_2$. How so? I guess this is the perfect time to look into its proof.

##### Proof of Uniqueness

Since there could be infinite solutions, let $x_1$ and $x_2$ be two such solution. Hence we can say the following:

$x_1 \equiv a_1 \text{(mod } m_1)$
$x_2 \equiv a_1 \text{(mod } m_1)$
$\therefore x_1 \equiv x_2 \text{(mod } m_1)$
$x_1 - x_2 \equiv 0 \text{(mod } m_1)$
$\therefore m_1 | x_1 - x_2$

Similarly, we can show that $m_2 | x_1 - x_2$.

What can we do with all these information? The difference of any two solutions $x_1$ and $x_2$ is divisible by both $m_1$ and $m_2$. We also know that, $m_1$ and $m_2$ are coprime. Doesn't that mean, the difference is also divisible by $m_1 \times m_2$? Yes. That's exactly what we are going towards

$m_1m_2 | x_1 - x_2$
$x_1 - x_2 \equiv 0 \text{(mod } m_1m_2)$
$x_1 \equiv x_2 \text{(mod } m_1m_2)$

The next question I will ask is: if the difference between any two solutions, $x_1$ and $x_2$, is divisible by $m_1m_2$, then how many solution can there exist in the range $0$ to $m_1m_2 - 1$? Only one. And hence, the solution $x$ is unique to modulo $m_1m_2$.

Therefore we can say that the solution $x \equiv a_1 m_2 q + a_2 m_1 p \text{ (mod }m_1m_2)$ is unique and smallest. Hence, the two equation, $x \equiv a_1 \text{(mod }m_1)$ and $x \equiv a_2 \text{(mod }m_2)$, after merging becomes $x \equiv a_1 m_2 q + a_2 m_1 p \text{ (mod }m_1m_2)$

#### When there are $n$ equations

Since we know how to merge two equations into one, all we need to do is take the first two equations from $n$ equations and merge them. Then we are left with $n-1$ equations. Since all the elements of $M$ were pairwise coprime, the new modulus is also coprime. Hence, we can continue to merge the equations until there is only one equation left. The equation left is our answer.

### Weak form of CRT: Code

The following code implements the idea we discussed above. Note that I used int data type, but in most of the problems you will most likely require long long. I am sure you can adjust the code accordingly yourself.

/** Return {-1,-1} if invalid input.
Otherwise, returns {x,L}, where x is the solution unique to mod L
*/
pair<int, int> chinese_remainder_theorem( vector<int> A, vector<int> M ) {
if(A.size() != M.size()) return {-1,-1}; /** Invalid input*/

int n = A.size();

int a1 = A[0];
int m1 = M[0];
/** Initially x = a_0 (mod m_0)*/

/** Merge the solution with remaining equations */
for ( int i = 1; i < n; i++ ) {
int a2 = A[i];
int m2 = M[i];

/** Merge the two equations*/
int p, q;
ext_gcd(m1, m2, &p, &q);

/** We need to be careful about overflow, but I did not bother about overflow here to keep the code simple.*/
int x = (a1*m2*q + a2*m1*p) % (m1*m2);

/** Merged equation*/
a1 = x;
m1 = m1 * m2;
}
if (a1 < 0) a1 += m1; /** Result is not suppose to be negative*/
return {a1, m1};
}


Once again, please be careful about overflow when solving a problem. From my experience, I have seen that the value of $p$ and $q$ becomes large quickly and intermediate calculations no longer fit into long long variables. I use __int128 data type to avoid overflow issues. So if you get "Wrong Answer", it is most likely due to overflow.

Complexity: $O(n \times log(L))$

## Conclusion

With this, we are done with Weak Form of Chinese Remainder Theorem. We can now find a solution to congruence equations when the moduli are co-prime. On next post, we will see a stronger version of Chinese Remainder Theorem, where the moduli are not co-prime.

Edit: The next post can be found here: Chinese Remainder Theorem Part 2 - Non Coprime Moduli

## Resources

1. Wiki - Chinese Remainder Theorem - https://en.wikipedia.org/wiki/Chinese_remainder_theorem
2. Brilliant.org - Chinese Remainder Theorem - https://brilliant.org/wiki/chinese-remainder-theorem/
3. Cut The Knot - Chinese Remainder Theorem - https://www.cut-the-knot.org/blue/chinese.shtml

## Related Problems

That's all I could find on CRT. If you know about any other related problem, please feel free to add it as a comment.

## Tuesday, September 29, 2015

### Modular Inverse from 1 to N

We already learned how to find Modular Inverse for a particular number in a previous post, "Modular Multiplicative Inverse". Today we will look into finding Modular Inverse in a bulk.

## Problem

Given $N$ and $M$ ( $N < M$ and $M$ is prime ), find modular inverse of all numbers between $1$ to $N$ with respect to $M$.

Since $M$ is prime and $N$ is less than $M$, we can be sure that Modular Inverse exists for all numbers. Why? Cause prime numbers are coprime to all numbers less than them.

We will look into two methods. Later one is better than the first one.

## $O(NlogM)$ Solution

Using Fermat's little theorem, we can easily find Modular Inverse for a particular number.

$A^{-1} \:\%\: M = bigmod(A,M-2,M)$, where $bigmod()$ is a function from the post "Repeated Squaring Method for Modular Exponentiation". The function has complexity of $O(logM)$. Since we are trying to find inverse for all numbers from $1$ to $N$, we can find them in $O(NlogM)$ complexity by running a loop.
int inv[SIZE]; ///inv[x] contains value of (x^-1 % m)
for ( int i = 1; i <= n; i++ ) {
inv[i] = bigmod ( i, m - 2, m );
}
But it's possible to do better.

## $O(N)$ Solution

This solution is derived using some clever manipulation of Modular Arithmetic.

Suppose we are trying to find the modular inverse for a number $a$, $a < M$, with respect to $M$. Now divide $M$ by $a$. This will be the starting point.

$M = Q \times a + r$, (where $Q$ is the quotient and $r$ is the remainder)
$M = \lfloor \frac{M}{a} \rfloor \times a + (M \:\%\: a )$

Now take modulo $M$ on both sides.

$0 \equiv \lfloor \frac{M}{a} \rfloor \times a + (M \:\%\: a ) \:\:\:\text{(mod M )}$
$(M \:\%\: a ) \equiv -\lfloor \frac{M}{a} \rfloor \times a \:\:\:\text{(mod M )}$

Now divide both side by $a \times ( M \:\%\: a )$.

$$\frac{M \:\%\: a}{a \times ( M \:\%\: a )} \equiv \frac{- \lfloor \frac{M}{a} \rfloor \times a } { a \times ( M \:\%\: a ) } \:\:\:\text{(mod M)}$$
$$\therefore a^{-1} \equiv - \lfloor \frac{M}{a} \rfloor \times ( M \:\%\: a )^{-1} \:\:\:\text{(mod M)}$$

The formula establishes a recurrence relation. The formula says that, in order to find the modular inverse of $a$, we need to find the modular inverse of $b = M \:\%\: a$ first.

Since $b = M \:\%\: a$, we can say that its value lies between $0$ and $a-1$. But, $a$ and $M$ are coprime. So $a$ will never fully divide $M$. Hence we can ignore the possibility that $b$ will be $0$. So possible values of $b$ is between $1$ and $a-1$.

Therefore, if we have all modular inverse from $1$ to $a-1$ already calculated, then we can find the modular inverse of $a$ in $O(1)$.

## Code

We can now formulate our code.
int inv[SIZE];
inv[1] = 1;
for ( int i = 2; i <= n; i++ ) {
inv[i] = (-(m/i) * inv[m%i] ) % m;
inv[i] = inv[i] + m;
}

In line $2$, we set the base case. Modular inverse of $1$ is always $1$. Then we start calculating inverse from $2$ to $N$. When $i=2$, all modular inverse from $1$ to $i-1=1$ is already calculated in array $inv[]$. So we can calculate it in $O(1)$ using the formula above at line $4$.

At line $5$, we make sure the modular inverse is non-negative.

Next, when $i=3$, all modular inverse from $1$ to $i-1=2$ is already calculated. This is process is repeated till we reach $N$.

Since we calculated each inverse in $O(1)$, the complexity of this code is $O(N)$.

## Conclusion

I saw this code first time on CodeChef forum. I didn't know how it worked back then. I added it to my notebook and have been using it since then. Recently, while searching over the net for resources on Pollard Rho's algorithm, I stumbled on an article from Come On Code On which had the explanation. Thanks, fR0DDY, I have been looking for the proof.

## Reference

1. forthright48 - Modular Multiplicative Inverse
2. forthright48 - Repeated Squaring Method for Modular Exponentiation
3. Come On Code On - Modular Multiplicative Inverse