Skip to content
This repository has been archived by the owner on Mar 28, 2019. It is now read-only.

Commit

Permalink
Upgrade: v1.1.2
Browse files Browse the repository at this point in the history
* Renew Bézout Equation algorithm
* Fix Eratosthenes Sieve when given large upper bound
* Update Congruence solution with Bézout in Chinese Remainder Theorem
  • Loading branch information
JarryShaw committed Jul 15, 2017
1 parent 9ce44e8 commit 423e4da
Show file tree
Hide file tree
Showing 5 changed files with 88 additions and 80 deletions.
61 changes: 32 additions & 29 deletions jsntlib/NTLArchive/NTLBezoutEquation.py
Original file line number Diff line number Diff line change
Expand Up @@ -14,50 +14,53 @@

'''Usage sample:
print('%d*-179 + %d*-367 = (-179,-367)' % bezout(-179,-367))
print('%d*-179 + %d*-367 = (-179,-367)' % bezout(-179, -367))
'''


def bezoutEquation(a, b):
int_check(a, b)

exflag = pn_a = pn_b = 0
if a < 0: a *= -1; pn_a = 1
if b < 0: b *= -1; pn_b = 1
if a < b: a, b = b, a; pn_a, pn_b = pn_b, pn_a; exflag = 1 # 交換a與b的次序,使得a≥b
exflag = pn_a = pn_b = False
if a < 0: a *= -1; pn_a = True
if b < 0: b *= -1; pn_b = True
if a < b: a, b = b, a; pn_a, pn_b = pn_b, pn_a; exflag = True # 交換a與b的次序,使得a≥b

q = [0] + euclideanAlgorithm(a, b) # 廣義歐幾里德除法(輾轉相除法),求不完全商數組q
s = coefficient_s(q, 1, 0, 0) # 求係數s
t = coefficient_t(q, 0, 1, 0) # 求係數t
s, t = coefficients(q) # 求係數s和t

if pn_a ^ pn_b:
if exflag:
tmp = s
s = -t if pn_b else t
t = -tmp if pn_a else tmp

else:
if pn_b: s = -s
if pn_a: t = -t

elif pn_a and pn_b:
if exflag: s, t = -t, -s
else: s, t = -s, -t

if exflag:
s *= -1 if pn_b else 1; t *= -1 if pn_a else 1
else:
s *= -1 if pn_a else 1; t *= -1 if pn_b else 1
if exflag: s, t = t, s

return s, t


def coefficient_s(q_j, s_j1, s_j2, ctr):
try:
s = -1 * q_j[ctr] * s_j1 + s_j2 # s_j = (-q_j) * s_j-1 + s_j-2
except IndexError:
return s_j1
def coefficients(q):
s_j1 = 1; s_j2 = 0
t_j1 = 0; t_j2 = 1

s_j2 = s_j1
s_j1 = s
ctr += 1
return coefficient_t(q_j, s_j1, s_j2, ctr)
for q_j in q:
s_j = -q_j * s_j1 + s_j2 # s_j = (-q_j) * s_j-1 + s_j-2
t_j = -q_j * t_j1 + t_j2 # t_j = (-q_j) * t_j-1 + t_j-2

s_j2 = s_j1; s_j1 = s_j
t_j2 = t_j1; t_j1 = t_j
else:
return s_j, t_j

def coefficient_t(q_j, t_j1, t_j2, ctr):
try:
t = -1 * q_j[ctr] * t_j1 + t_j2 # t_j = (-q_j) * t_j-1 + t_j-2
except IndexError:
return t_j1

t_j2 = t_j1
t_j1 = t
ctr += 1
return coefficient_t(q_j, t_j1, t_j2, ctr)
return (0, 1)
17 changes: 4 additions & 13 deletions jsntlib/NTLArchive/NTLChineseRemainderTheorem.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,7 @@
# 求基本同餘式組的通解


from .NTLBezoutEquation import bezoutEquation
from .NTLExceptions import DefinitionError
from .NTLUtilities import jsrange
from .NTLValidations import int_check, list_check, tuple_check
Expand Down Expand Up @@ -50,25 +51,15 @@ def CHNRemainderTheorem(*args):
modulo *= tmpMod1 # M(original modulo) = ∏m_i

bList = []
for tmpMod2 in mod:
M = modulo // tmpMod2 # M_i = M / m_i
t = solve(M, tmpMod2) # t_i * M_i ≡ 1 (mod m_i)
for tmpMod in mod:
M = modulo // tmpMod # M_i = M / m_i
t = bezoutEquation(M, tmpMod)[0] # t_i * M_i ≡ 1 (mod m_i)
bList.append(t * M) # b_i = t_i * M_i

remainder = iterCalc(rmd, bList, modulo) # x_j = Σ(b_i * r_i) (mod M)
return sorted(remainder)


# 求解M_i^-1 (mod m_i)
def solve(variable, modulo):
polyCgc = '%d*x - 1' % variable # 將係數與指數數組生成多項式

r = lambda x: eval(polyCgc) # 用於計算多項式的取值
for x in jsrange(modulo): # 逐一驗算,如模為0則加入結果數組
if r(x) % modulo == 0:
return x


# 對rmd多維數組(層,號)中的數進行全排列並計算結果
def iterCalc(ognList, coeList, modulo):
ptrList = [] # 寄存指向每一數組層的號
Expand Down
13 changes: 2 additions & 11 deletions jsntlib/NTLArchive/NTLCongruence.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@
# 由多項式衍生的同餘式


from .NTLBezoutEquation import bezoutEquation
from .NTLExceptions import DefinitionError, PCError, PolyError, SolutionError
from .NTLGreatestCommonDivisor import greatestCommonDivisor
from .NTLPolynomial import Polynomial
Expand Down Expand Up @@ -235,22 +236,12 @@ def _primePro(self, rem, exp):

# 中國剩餘定理
def _CTR(self, rem, mod):

# 求解M_i^-1 (mod m_i)
def solve(variable, modulo):
polyCgc = '%d*x - 1' % variable # 將係數與指數數組生成多項式

r = lambda x: eval(polyCgc) # 用於計算多項式的取值
for x in jsrange(modulo): # 逐一驗算,如模為0則加入結果數組
if r(x) % modulo == 0:
return x

modulo = self._modulo # M(original modulo) = ∏m_i

bList = []
for tmpMod in mod:
M = modulo // tmpMod # M_i = M / m_i
t = solve(M, tmpMod) # t_i * M_i ≡ 1 (mod m_i)
t = bezoutEquation(M, tmpMod)[0] # t_i * M_i ≡ 1 (mod m_i)
bList.append(t * M) # b_i = t_i * M_i

_rem = self._iterCalc(rem, bList) # x_j = Σ(b_i * r_i) (mod M)
Expand Down
73 changes: 48 additions & 25 deletions jsntlib/NTLArchive/NTLEratosthenesSieve.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,8 @@
from .NTLUtilities import jsrange
from .NTLValidations import int_check, pos_check

# from NTLTrivialDivision import trivialDivision


__all__ = ['eratosthenesSieve']
nickname = 'primelist'
Expand All @@ -24,7 +26,7 @@


# 前1280個素數表,簡化小素數的求取
PRIME_LIST = [
_PRIME_LIST = [
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97],
[101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199],
[211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293],
Expand Down Expand Up @@ -142,42 +144,63 @@ def eratosthenesSieve(upper, lower=None):
if upper <= 10100:
_llevel = lower // 100; _ulevel = upper // 100

for ptr_1 in jsrange(len(PRIME_LIST[_llevel])):
if PRIME_LIST[_llevel][ptr_1] >= lower:
for ptr_1 in jsrange(len(_PRIME_LIST[_llevel])):
if _PRIME_LIST[_llevel][ptr_1] >= lower:
llevel = ptr_1; break
else:
_llevel += 1; llevel = 0

for ptr_2 in jsrange(len(PRIME_LIST[_ulevel])):
if PRIME_LIST[_ulevel][ptr_2] >= upper:
for ptr_2 in jsrange(len(_PRIME_LIST[_ulevel])):
if _PRIME_LIST[_ulevel][ptr_2] >= upper:
ulevel = ptr_2; break
else:
ulevel = None

if _llevel != _ulevel:
rst = PRIME_LIST[_llevel][llevel:]
rst = _PRIME_LIST[_llevel][llevel:]
for _level in jsrange(_llevel+1, _ulevel):
rst += PRIME_LIST[_level]
rst += PRIME_LIST[_ulevel][:ulevel]
rst += _PRIME_LIST[_level]
rst += _PRIME_LIST[_ulevel][:ulevel]
else:
rst = PRIME_LIST[_llevel][llevel:ulevel]
rst = _PRIME_LIST[_llevel][llevel:ulevel]

else:
# 用於存儲upper個正整數的表格/狀態;其中,0表示篩去,1表示保留
table = [1]*(upper+1)

# 篩法(平凡除法)
for index in jsrange(2, int(math.sqrt(upper))+1):
tmp = index * 2
if table[index] == 1:
while tmp <= upper:
# 將index的倍數篩去
table[tmp] = 0; tmp += index

# 獲取結果
rst = []
for ptr in jsrange(lower, upper+1):
if table[ptr] == 1:
rst.append(ptr)
from NTLTrivialDivision import trivialDivision

def _eratosthenesSieve(upper, lower):
# 用於存儲upper個正整數的表格/狀態;其中,0表示篩去,1表示保留
table = [1]*(upper+1)

# 篩法(平凡除法)
for index in jsrange(2, int(math.sqrt(upper))+1):
tmp = index * 2
if table[index] == 1:
while tmp <= upper:
# 將index的倍數篩去
table[tmp] = 0; tmp += index

# 獲取結果
rst = []
for ptr in jsrange(lower, upper+1):
if table[ptr] == 1:
rst.append(ptr)

return rst

if math.log(upper, 10) > 4.0:
if math.log(lower, 10) > 4.0:
rst = []
start = lower if lower % 2 == 1 else lower + 1

else:
rst = _eratosthenesSieve(10000, lower)
start = 10001

for _int in jsrange(start, upper, 2):
if trivialDivision(_int):
rst.append(_int)

else:
rst = _eratosthenesSieve(upper, lower)

return rst
4 changes: 2 additions & 2 deletions setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@

name = 'jsntlib',

version = '1.1.1',
version = '1.1.2',

description =
'A number theory library adapted from mathematic fundamentals of information security homework codes.',
Expand All @@ -32,7 +32,7 @@

url = 'https://github.com/JarryShaw/jsntlib/tree/release',

download_url = 'https://github.com/JarryShaw/jsntlib/archive/1.1.1.tar.gz',
download_url = 'https://github.com/JarryShaw/jsntlib/archive/1.1.2.tar.gz',

packages = [
'jsntlib',
Expand Down

0 comments on commit 423e4da

Please sign in to comment.