Chinese Remainder Theorem | Diophantine Equation | Modular Division

blockchain/chinese_remainder_theorem.py

@@ -0,0 +1,91 @@
+# Chinese Remainder Theorem:
+# GCD ( Greatest Common Divisor ) or HCF ( Highest Common Factor )
+
+# If GCD(a,b) = 1, then for any remainder ra modulo a and any remainder rb modulo b there exists integer n,
+# such that n = ra (mod a) and n = ra(mod b).  If n1 and n2 are two such integers, then n1=n2(mod ab)
+
+# Algorithm :
+
+# 1. Use extended euclid algorithm to find x,y such that a*x + b*y = 1
+# 2. Take n = ra*by + rb*ax
+
+
+# Extended Euclid
+def extended_euclid(a, b):
+    """
+    >>> extended_euclid(10, 6)
+    (-1, 2)
+
+    >>> extended_euclid(7, 5)
+    (-2, 3)
+
+    """
+    if b == 0:
+        return (1, 0)
+    (x, y) = extended_euclid(b, a % b)
+    k = a // b
+    return (y, x - k * y)
+
+
+# Uses ExtendedEuclid to find inverses
+def chinese_remainder_theorem(n1, r1, n2, r2):
+    """
+    >>> chinese_remainder_theorem(5,1,7,3)
+    31
+
+    Explanation : 31 is the smallest number such that
+                (i)  When we divide it by 5, we get remainder 1
+                (ii) When we divide it by 7, we get remainder 3
+
+    >>> chinese_remainder_theorem(6,1,4,3)
+    14
+
+    """
+    (x, y) = extended_euclid(n1, n2)
+    m = n1 * n2
+    n = r2 * x * n1 + r1 * y * n2
+    return ((n % m + m) % m)
+
+
+# ----------SAME SOLUTION USING InvertModulo instead ExtendedEuclid----------------
+
+# This function find the inverses of a i.e., a^(-1)
+def invert_modulo(a, n):
+    """
+    >>> invert_modulo(2, 5)
+    3
+
+    >>> invert_modulo(8,7)
+    1
+
+    """
+    (b, x) = extended_euclid(a, n)
+    if b < 0:
+        b = (b % n + n) % n
+    return b
+
+
+# Same a above using InvertingModulo
+def chinese_remainder_theorem2(n1, r1, n2, r2):
+    """
+    >>> chinese_remainder_theorem2(5,1,7,3)
+    31
+
+    >>> chinese_remainder_theorem2(6,1,4,3)
+    14
+
+    """
+    x, y = invert_modulo(n1, n2), invert_modulo(n2, n1)
+    m = n1 * n2
+    n = r2 * x * n1 + r1 * y * n2
+    return (n % m + m) % m
+
+
+# import testmod for testing our function
+from doctest import testmod
+
+if __name__ == '__main__':
+    testmod(name='chinese_remainder_theorem', verbose=True)
+    testmod(name='chinese_remainder_theorem2', verbose=True)
+    testmod(name='invert_modulo', verbose=True)
+    testmod(name='extended_euclid', verbose=True)

blockchain/diophantine_equation.py

@@ -0,0 +1,124 @@
+# Diophantine Equation : Given integers a,b,c ( at least one of a and b != 0), the diophantine equation
+# a*x + b*y = c has a solution (where x and y are integers) iff gcd(a,b) divides c.
+
+# GCD ( Greatest Common Divisor ) or HCF ( Highest Common Factor )
+
+
+def diophantine(a, b, c):
+    """
+    >>> diophantine(10,6,14)
+    (-7.0, 14.0)
+
+    >>> diophantine(391,299,-69)
+    (9.0, -12.0)
+
+    But above equation has one more solution i.e., x = -4, y = 5.
+    That's why we need diophantine all solution function.
+
+    """
+
+    assert c % greatest_common_divisor(a, b) == 0  # greatest_common_divisor(a,b) function implemented below
+    (d, x, y) = extended_gcd(a, b)  # extended_gcd(a,b) function implemented below
+    r = c / d
+    return (r * x, r * y)
+
+
+# Lemma : if n|ab and gcd(a,n) = 1, then n|b.
+
+# Finding All solutions of Diophantine Equations:
+
+# Theorem : Let gcd(a,b) = d, a = d*p, b = d*q. If (x0,y0) is a solution of Diophantine Equation a*x + b*y = c.
+# a*x0 + b*y0 = c, then all the solutions have the form a(x0 + t*q) + b(y0 - t*p) = c, where t is an arbitrary integer.
+
+# n is the number of solution you want, n = 2 by default
+
+def diophantine_all_soln(a, b, c, n=2):
+    """
+    >>> diophantine_all_soln(10, 6, 14)
+    -7.0 14.0
+    -4.0 9.0
+
+    >>> diophantine_all_soln(10, 6, 14, 4)
+    -7.0 14.0
+    -4.0 9.0
+    -1.0 4.0
+    2.0 -1.0
+
+    >>> diophantine_all_soln(391, 299, -69, n = 4)
+    9.0 -12.0
+    22.0 -29.0
+    35.0 -46.0
+    48.0 -63.0
+
+    """
+    (x0, y0) = diophantine(a, b, c)  # Initial value
+    d = greatest_common_divisor(a, b)
+    p = a // d
+    q = b // d
+
+    for i in range(n):
+        x = x0 + i * q
+        y = y0 - i * p
+        print(x, y)
+
+
+# Euclid's Lemma :  d divides a and b, if and only if d divides a-b and b
+
+# Euclid's Algorithm
+
+def greatest_common_divisor(a, b):
+    """
+    >>> greatest_common_divisor(7,5)
+    1
+
+    Note : In number theory, two integers a and b are said to be relatively prime, mutually prime, or co-prime
+           if the only positive integer (factor) that divides both of them is 1  i.e., gcd(a,b) = 1.
+
+    >>> greatest_common_divisor(121, 11)
+    11
+
+    """
+    if a < b:
+        a, b = b, a
+
+    while a % b != 0:
+        a, b = b, a % b
+
+    return b
+
+
+# Extended Euclid's Algorithm : If d divides a and b and d = a*x + b*y for integers x and y, then d = gcd(a,b)
+
+
+def extended_gcd(a, b):
+    """
+    >>> extended_gcd(10, 6)
+    (2, -1, 2)
+
+    >>> extended_gcd(7, 5)
+    (1, -2, 3)
+
+    """
+    assert a >= 0 and b >= 0
+
+    if b == 0:
+        d, x, y = a, 1, 0
+    else:
+        (d, p, q) = extended_gcd(b, a % b)
+        x = q
+        y = p - q * (a // b)
+
+    assert a % d == 0 and b % d == 0
+    assert d == a * x + b * y
+
+    return (d, x, y)
+
+
+# import testmod for testing our function
+from doctest import testmod
+
+if __name__ == '__main__':
+    testmod(name='diophantine', verbose=True)
+    testmod(name='diophantine_all_soln', verbose=True)
+    testmod(name='extended_gcd', verbose=True)
+    testmod(name='greatest_common_divisor', verbose=True)

blockchain/modular_division.py

@@ -0,0 +1,149 @@
+# Modular Division :
+# An efficient algorithm for dividing b by a modulo n.
+
+# GCD ( Greatest Common Divisor ) or HCF ( Highest Common Factor )
+
+# Given three integers a, b, and n, such that gcd(a,n)=1 and n>1, the algorithm should return an integer x such that
+#        0≤x≤n−1, and  b/a=x(modn) (that is, b=ax(modn)).
+
+# Theorem:
+# a has a multiplicative inverse modulo n iff gcd(a,n) = 1
+
+
+# This find x = b*a^(-1) mod n
+# Uses ExtendedEuclid to find the inverse of a
+
+
+def modular_division(a, b, n):
+    """
+    >>> modular_division(4,8,5)
+    2
+
+    >>> modular_division(3,8,5)
+    1
+
+    >>> modular_division(4, 11, 5)
+    4
+
+    """
+    assert n > 1 and a > 0 and greatest_common_divisor(a, n) == 1
+    (d, t, s) = extended_gcd(n, a)  # Implemented below
+    x = (b * s) % n
+    return x
+
+
+# This function find the inverses of a i.e., a^(-1)
+def invert_modulo(a, n):
+    """
+    >>> invert_modulo(2, 5)
+    3
+
+    >>> invert_modulo(8,7)
+    1
+
+    """
+    (b, x) = extended_euclid(a, n)  # Implemented below
+    if b < 0:
+        b = (b % n + n) % n
+    return b
+
+
+# ------------------ Finding Modular division using invert_modulo -------------------
+
+# This function used the above inversion of a to find x = (b*a^(-1))mod n
+def modular_division2(a, b, n):
+    """
+    >>> modular_division2(4,8,5)
+    2
+
+    >>> modular_division2(3,8,5)
+    1
+
+    >>> modular_division2(4, 11, 5)
+    4
+
+    """
+    s = invert_modulo(a, n)
+    x = (b * s) % n
+    return x
+
+
+# Extended Euclid's Algorithm : If d divides a and b and d = a*x + b*y for integers x and y, then d = gcd(a,b)
+
+def extended_gcd(a, b):
+    """
+    >>> extended_gcd(10, 6)
+    (2, -1, 2)
+
+    >>> extended_gcd(7, 5)
+    (1, -2, 3)
+
+   ** extended_gcd function is used when d = gcd(a,b) is required in output
+
+    """
+    assert a >= 0 and b >= 0
+
+    if b == 0:
+        d, x, y = a, 1, 0
+    else:
+        (d, p, q) = extended_gcd(b, a % b)
+        x = q
+        y = p - q * (a // b)
+
+    assert a % d == 0 and b % d == 0
+    assert d == a * x + b * y
+
+    return (d, x, y)
+
+
+# Extended Euclid
+def extended_euclid(a, b):
+    """
+    >>> extended_euclid(10, 6)
+    (-1, 2)
+
+    >>> extended_euclid(7, 5)
+    (-2, 3)
+
+    """
+    if b == 0:
+        return (1, 0)
+    (x, y) = extended_euclid(b, a % b)
+    k = a // b
+    return (y, x - k * y)
+
+
+# Euclid's Lemma :  d divides a and b, if and only if d divides a-b and b
+# Euclid's Algorithm
+
+def greatest_common_divisor(a, b):
+    """
+    >>> greatest_common_divisor(7,5)
+    1
+
+    Note : In number theory, two integers a and b are said to be relatively prime, mutually prime, or co-prime
+           if the only positive integer (factor) that divides both of them is 1  i.e., gcd(a,b) = 1.
+
+    >>> greatest_common_divisor(121, 11)
+    11
+
+    """
+    if a < b:
+        a, b = b, a
+
+    while a % b != 0:
+        a, b = b, a % b
+
+    return b
+
+
+# Import testmod for testing our function
+from doctest import testmod
+
+if __name__ == '__main__':
+    testmod(name='modular_division', verbose=True)
+    testmod(name='modular_division2', verbose=True)
+    testmod(name='invert_modulo', verbose=True)
+    testmod(name='extended_gcd', verbose=True)
+    testmod(name='extended_euclid', verbose=True)
+    testmod(name='greatest_common_divisor', verbose=True)