# Solving systems of congruences with Gaussian integers

I am writing my own function to solve systems of congruences in arbitrary euclidean domains, but I can't quite make it work for gaussians.

# Extended Euclides Algorithm
def extended_euclides(a,b,quo=lambda a,b:a//b):
r0 = a; r1 = b
s0 = 1; s1 = 0
t0 = 0; t1 = 1

while r1 != 0:
q = quo(r0, r1)
r0, r1 = r1, r0 - q * r1
s0, s1 = s1, s0 - q * s1
t0, t1 = t1, t0 - q * t1

return r0, s0, t0

# Chinese Remainder Theorem
def chinese_remainder(remainders, modules, quo=lambda a,b:a//b):
"""
INPUT:
A list of remainders v_i and pairwise coprime modules m_i
representing a system of congruences {x == v_i mod m_i},
where v_i, m_i belong to an euclidean domain ED
OUTPUT:
A solution of the system of congruences in the domain ED
"""
m = reduce(lambda x,y:x*y, modules)
c = 0
for v_i, m_i in zip(remainders, modules):
m_div_m_i = quo(m, m_i)
_, s_i, _ = extended_euclides(m_div_m_i, m_i, quo)
# c_i = (v_i * s_i) % m_i
# we need to do it convolutedly bcs gaussian integer quotient is not implemented in SageMath
a = v_i * s_i; b = m_i
c_i = a - quo(a,b)*b

c += c_i * m_div_m_i
return c

# EXAMPLE ON Z
remainders = [ZZ(2),ZZ(7)]
modules = [ZZ(11), ZZ(13)]
print("System to be solved:")
for r,m in zip(remainders, modules):
print("x == {} (mod {})".format(r,m))
solution = chinese_remainder(remainders, modules)
print("A solution of the system is {}".format(solution))
#check that the solution is correct
for r, m in zip(remainders, modules):
assert(solution % m == r)

# EXAMPLE ON Z[i]
gaussian_quo = lambda a,b : ZI(real(a/b).round() + I*imag(a/b).round())
gaussian_rem=lambda a,b: a - gaussian_quo(a,b)*b
remainders = [ZI(2+3*I),ZI(7+5*I)]
modules = [ZI(11), ZI(13-I)]
print("System to be solved:")
for r,m in zip(remainders, modules):
print("x == {} (mod {})".format(r,m))
solution = chinese_remainder(remainders, modules, gaussian_quo)
print("A solution of the system is {}".format(solution))
#check that the solution is correct
for r, m in zip(remainders, modules):
print(gaussian_rem(solution, m) / r)
assert(gaussian_rem(solution, m) == r)


This produces the following ouptut:

System to be solved:
x == 2 (mod 11)
x == 7 (mod 13)
A solution of the system is 46
System to be solved:
x == 3*I + 2 (mod 11)
x == 5*I + 7 (mod -I + 13)
A solution of the system is -41*I - 75
1
36/37*I - 6/37
Error in lines 31-33
Traceback (most recent call last):
File "/cocalc/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 1188, in execute
flags=compile_flags) in namespace, locals
File "", line 3, in <module>
AssertionError


So in the integer example we find a correct solution, but in the Gaussian example we find a solution that satisfies the first incongruence relation, but not the second.

I have checked with wolfram alpha, and if I did not mess up it seems like the problem is in the Chinese Remainder Theorem algorithm rather than on the assertion tests.

What is going on?

edit retag close merge delete

1

Note that CRT is a theorem for finitely generated abelian groups. That the additive group Z[i] also has a (euclidean) ring structure is not required for CRT. All you need is a list of finite index subgroups that are "pairwise coprime" (that any two of then generate the entire additive group). In particular, it works with any ring of integer and a list of coprime ideals. It doesn't need the ideals to be principal. See idealchinesein pari/gp: https://pari.math.u-bordeaux.fr/docht...

( 2018-11-15 18:03:41 -0500 )edit

Sort by » oldest newest most voted

The code and output is correct. The remainder after dividing -41*I - 75 by -I + 13 is 6*I - 6, which is congruent to 5*I + 7 modulo -I + 13. Your final comparison

gaussian_rem(solution, m) == r


is not the right thing to test, because it implicitly assumes that the right-hand sides of your system of equations were already in normal form. Instead, you want to test

gaussian_rem(solution, m) == gaussian_rem(r, m)


or more efficiently

gaussian_rem(solution - r, m) == 0

more