Competitive Programming in Haskell: modular arithmetic, part 2

In my last post I wrote about modular exponentiation and egcd. In this post, I consider the problem of solving modular equivalences, building on code from the previous post.

Solving linear congruences

A linear congruence is a modular equivalence of the form

.

Let’s write a function to solve such equivalences for . We want a pair of integers and such that is a solution to if and only if . This isn’t hard to write in the end, but takes a little bit of thought to do it properly.

First of all, if and are relatively prime (that is, ) then we know from the last post that has an inverse modulo ; multiplying both sides by yields the solution .

OK, but what if ? In this case there might not even be any solutions. For example, has no solutions: any even number will be equivalent to or modulo , so there is no value of such that double it will be equivalent to . On the other hand, is OK: this will be true for any odd value of , that is, . In fact, it is easy to see that any common divisor of and must also divide in order to have any solutions. In case the GCD of and does divide , we can simply divide through by the GCD (including dividing the modulus !) and then solve the resulting equivalence.

-- solveMod a b m solves ax = b (mod m), returning a pair (y,k) (with
-- 0 <= y < k) such that x is a solution iff x = y (mod k).
solveMod :: Integer -> Integer -> Integer -> Maybe (Integer, Integer)
solveMod a b m
  | g == 1         = Just ((b * inverse m a) `mod` m, m)
  | b `mod` g == 0 = solveMod (a `div` g) (b `div` g) (m `div` g)
  | otherwise      = Nothing
  where
    g = gcd a m

Solving systems of congruences with CRT

In its most basic form, the Chinese remainder theorem (CRT) says that if we have a system of two modular equations

then as long as and are relatively prime, there is a unique solution for modulo the product ; that is, the system of two equations is equivalent to a single equation of the form

We first compute the Bézout coefficients and such that using egcd, and then compute the solution as . Indeed,

and hence ; similarly .

However, this is not quite general enough: we want to still be able to say something useful even if . I won’t go through the whole proof, but it turns out that there is a solution if and only if , and we can just divide everything through by , as we did for solving linear congruences. Here’s the code:

-- gcrt2 (a,n) (b,m) solves the pair of modular equations
--
--   x = a (mod n)
--   x = b (mod m)
--
-- It returns a pair (c, k) such that all solutions for x satisfy x =
-- c (mod k), that is, solutions are of the form x = kt + c for
-- integer t.
gcrt2 :: (Integer, Integer) -> (Integer, Integer) -> Maybe (Integer, Integer)
gcrt2 (a,n) (b,m)
  | a `mod` g == b `mod` g = Just (((a*v*m + b*u*n) `div` g) `mod` k, k)
  | otherwise              = Nothing
  where
    (g,u,v) = egcd n m
    k = (m*n) `div` g

From here we can bootstrap ourselves into solving systems of more than two equations, by iteratively combining two equations into one.

-- gcrt solves a system of modular equations.  Each equation x = a
-- (mod n) is given as a pair (a,n).  Returns a pair (z, k) such that
-- solutions for x satisfy x = z (mod k), that is, solutions are of
-- the form x = kt + z for integer t.
gcrt :: [(Integer, Integer)] -> Maybe (Integer, Integer)
gcrt []         = Nothing
gcrt [e]        = Just e
gcrt (e1:e2:es) = gcrt2 e1 e2 >>= \e -> gcrt (e:es)

Practice problems

And here are a bunch of problems for you to practice!

• Chinese Remainder
• Chinese Remainder Theorem (non-relatively prime moduli)
• Disgruntled Judge
• Astro
• Heliocentric
• Remainder Reminder
• DVD Screensaver