iframe-proxy | Sunbelt Computer Software

History

152 lines (106 loc) · 3.16 KB

Raw

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

# 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)

>>> modular_division(3,8,5)

>>> modular_division(4, 11, 5)

"""

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)

>>> invert_modulo(8,7)

"""

(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)

>>> modular_division2(3,8,5)

>>> modular_division2(4, 11, 5)

"""

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)

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)

"""

if a < b:

a, b = b, a

while a % b != 0:

a, b = b, a % b

return b

if __name__ == "__main__":

from doctest import testmod

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)

Sunbelt Computer Software

PL/B Language Development and Support

Navigation Menu

Sunbelt Computer Software

PL/B Language Development and Support

FilesExpand file tree

modular_division.py

Latest commit

History

modular_division.py

File metadata and controls