NZMATH  1.2.0
About: NZMATH is a Python based number theory oriented calculation system.
  Fossies Dox: NZMATH-1.2.0.tar.gz  ("inofficial" and yet experimental doxygen-generated source code documentation)  

mpqs.py
Go to the documentation of this file.
1 import math
2 import time
3 import logging
4 import nzmath.arith1 as arith1
5 import nzmath.gcd as gcd
6 import nzmath.prime as prime
7 import nzmath.factor.find as find
8 
9 # TODO: add psyco to config entry
10 
15 
16 _log = logging.getLogger('nzmath.factor.mpqs')
17 
18 class QS(object):
19  def __init__(self, n, sieverange, factorbase):
20  self.number = n
21  self.sqrt_n = int(math.sqrt(n))
22  for i in [2,3,5,7,11,17,19]:
23  if n % i == 0:
24  raise ValueError("This number is divided by %d" % i)
25 
26  self.digit = arith1.log(self.number, 10) + 1
27  self.Srange = sieverange
28  self.FBN = factorbase
29 
30  self.move_range = range(self.sqrt_n-self.Srange, self.sqrt_n+self.Srange+1)
31 
32  i = 0
33  k = 0
34  factor_base = [-1]
35  FB_log = [0]
36  while True:
37  ii = primes_table[i]
38  if arith1.legendre(self.number, ii) == 1:
39  factor_base.append(ii)
40  FB_log.append(primes_log_table[i])
41  k += 1
42  i += 1
43  if k == self.FBN:
44  break
45  else:
46  i += 1
47 
48  self.FB = factor_base
49  self.FB_log = FB_log
50  self.maxFB = factor_base[-1]
51  N_sqrt_list = []
52  for i in self.FB:
53  if i != 2 and i != -1:
54  e = int(math.log(2*self.Srange, i))
55  N_sqrt_modp = sqroot_power(self.number, i, e)
56  N_sqrt_list.append(N_sqrt_modp)
57  self.solution = N_sqrt_list #This is square roots of N on Z/pZ, p in factor base.
58 
59  poly_table = []
60  log_poly = []
61  minus_val = []
62  for j in self.move_range:
63  jj = (j**2)-self.number
64  if jj < 0:
65  jj = -jj
66  minus_val.append(j-self.sqrt_n+self.Srange)
67  elif jj == 0:
68  jj = 1
69  lj = int((math.log(jj)*30)*0.97) # 0.97 is an erroe
70  poly_table.append(jj)
71  log_poly.append(lj)
72  self.poly_table = poly_table # This is Q(x) value , x in [-M+sqrt_n,M+sqrt_n].
73  self.log_poly = log_poly # This is log(Q(x)) value.
74  self.minus_check = minus_val # This is "x" that Q(x) is minus value.
75 
76  def run_sieve(self):
77  T = time.time()
78  M = self.Srange
79  start_location = []
80  logp = [0]*(2*M+1)
81  j = 2
82  for i in self.solution:
83  k = 0
84  start_p = []
85  while k < len(i):
86  q = int((self.sqrt_n)/(self.FB[j]**(k+1)))
87  s_1 = q*(self.FB[j]**(k+1))+i[k][0]
88  s_2 = q*(self.FB[j]**(k+1))+i[k][1]
89  while True:
90  if s_1 < self.sqrt_n-M:
91  s_1 = s_1+(self.FB[j]**(k+1))
92  break
93  else:
94  s_1 = s_1-(self.FB[j]**(k+1))
95  while True:
96  if s_2 < self.sqrt_n-M:
97  s_2 = s_2+(self.FB[j]**(k+1))
98  break
99  else:
100  s_2 = s_2-(self.FB[j]**(k+1))
101  start_p.append([s_1-self.sqrt_n+M,s_2-self.sqrt_n+M])
102 
103  k += 1
104  start_location.append(start_p)
105  j += 1
106  self.start_location = start_location
107 
108  if self.poly_table[0] & 1 == 0:
109  i = 0
110  while i <= 2*M:
111  j = 1
112  while True:
113  if self.poly_table[i] % (2**(j+1)) == 0:
114  j+=1
115  else:
116  break
117  logp[i] += self.FB_log[1]*j
118  i += 2
119  else:
120  i = 1
121  while i <= 2*M:
122  j = 1
123  while True:
124  if self.poly_table[i] % (2**(j+1)) == 0:
125  j += 1
126  else:
127  break
128  logp[i] += self.FB_log[1]*j
129  i += 2
130  L = 2
131  for j in self.start_location:
132  k = 0
133  while k < len(j):
134  s_1 = j[k][0]
135  s_2 = j[k][1]
136  h_1 = 0
137  h_2 = 0
138  while s_1+h_1 <= 2*M:
139  logp[s_1+h_1] += self.FB_log[L]
140  h_1 += self.FB[L]**(k+1)
141  while s_2+h_2 <= 2*M:
142  logp[s_2+h_2] += self.FB_log[L]
143  h_2 += self.FB[L]**(k+1)
144  k += 1
145  L += 1
146 
147  self.logp = logp
148  smooth = []
149  for t in range(2*M+1):
150  if logp[t] >= self.log_poly[t]:
151  poly_val = self.poly_table[t]
152  index_vector = []
153  for p in self.FB:
154  if p == -1:
155  if t in self.minus_check:
156  index_vector.append(1)
157  else:
158  index_vector.append(0)
159  else:
160  r = 0
161  while poly_val % (p**(r+1)) == 0:
162  r += 1
163  v = r & 1
164  index_vector.append(v)
165  smooth.append([index_vector, (poly_val, t+self.sqrt_n-M)])
166  _log.info(" Sieve time is %f sec" % (time.time()-T))
167  _log.info(" Found smooth numbers are %d / %d" % (len(smooth), len(self.FB)))
168  self.smooth = smooth
169  return smooth
170 
171 
172 class MPQS(object):
173  def __init__(self, n, sieverange=0, factorbase=0, multiplier=0):
174  self.number = n
175  _log.info("The number is %d MPQS starting" % n)
176 
177  if prime.primeq(self.number):
178  raise ValueError("This number is Prime Number")
179  for i in [2,3,5,7,11,13]:
180  if n % i == 0:
181  raise ValueError("This number is divided by %d" % i)
182 
183  self.sievingtime = 0
185  self.d_list = []
186  self.a_list = []
187  self.b_list = []
188 
189  #Decide prameters for each digits
190  self.digit = arith1.log(self.number, 10) + 1
191 
192  if sieverange != 0:
193  self.Srange = sieverange
194  if factorbase != 0:
195  self.FBN = factorbase
196  elif self.digit < 9:
197  self.FBN = parameters_for_mpqs[0][1]
198  else:
199  self.FBN = parameters_for_mpqs[self.digit-9][1]
200  elif factorbase != 0:
201  self.FBN = factorbase
202  if self.digit < 9:
203  self.Srange = parameters_for_mpqs[0][0]
204  else:
205  self.Srange = parameters_for_mpqs[self.digit-9][0]
206  elif self.digit < 9:
207  self.Srange = parameters_for_mpqs[0][0]
208  self.FBN = parameters_for_mpqs[0][1]
209  elif self.digit > 53:
210  self.Srange = parameters_for_mpqs[44][0]
211  self.FBN = parameters_for_mpqs[44][1]
212  else:
213  self.Srange = parameters_for_mpqs[self.digit-9][0]
214  self.FBN = parameters_for_mpqs[self.digit-9][1]
215 
216  self.move_range = range(-self.Srange, self.Srange+1)
217 
218  # Decide k such that k*n = 1 (mod4) and k*n has many factor base
219  if multiplier == 0:
220  self.sqrt_state = []
221  for i in [3,5,7,11,13]:
222  s = arith1.legendre(self.number, i)
223  self.sqrt_state.append(s)
224 
225  if self.number % 8 == 1 and self.sqrt_state == [1,1,1,1,1]:
226  k=1
227  else:
228  index8 = (self.number & 7) >> 1
229  j = 0
230  while self.sqrt_state != prime_8[index8][j][1]:
231  j += 1
232  k = prime_8[index8][j][0]
233  else:
234  if n & 3 == 1:
235  k = 1
236  else:
237  if multiplier == 1:
238  raise ValueError("This number is 3 mod 4 ")
239  else:
240  k = multiplier
241  self.number = k*self.number
242  self.multiplier = k
243 
244  _log.info("%d - digits Number" % self.digit)
245  _log.info("Multiplier is %d" % self.multiplier)
246 
247  # Table of (log p) , p in FB
248  i = 0
249  k = 0
250  factor_base = [-1]
251  FB_log = [0]
252  while k < self.FBN:
253  ii = primes_table[i]
254  if arith1.legendre(self.number,ii) == 1:
255  factor_base.append(ii)
256  FB_log.append(primes_log_table[i])
257  k += 1
258  i += 1
259 
260  self.FB = factor_base
261  self.FB_log = FB_log
262  self.maxFB = factor_base[-1]
263 
264  # Solve x^2 = n (mod p^e)
265  N_sqrt_list = []
266  for i in self.FB:
267  if i != 2 and i != -1:
268  e = int(math.log(2*self.Srange, i))
269  N_sqrt_modp = sqroot_power(self.number, i, e)
270  N_sqrt_list.append(N_sqrt_modp)
271  self.Nsqrt = N_sqrt_list
272 
273  def make_poly(self):
274  """
275  Make coefficients of f(x)= ax^2+b*x+c
276  """
277  T = time.time()
278  if self.d_list == []:
279  d = int(math.sqrt((math.sqrt(self.number)/(math.sqrt(2)*self.Srange))))
280  if d& 1 == 0:
281  if (d+1)& 3 == 1: #case d=0 mod4
282  d += 3
283  else:
284  d += 1 #case d=2 mod4
285  elif d& 3 == 1: #case d=1 mod4
286  d += 2
287  #case d=3 mod4
288  else:
289  d = self.d_list[-1]
290 
291  while d in self.d_list or not prime.primeq(d) or \
292  arith1.legendre(self.number,d) != 1 or d in self.FB:
293  d += 4
294  a = d**2
295  h_0 = pow(self.number, (d-3) >> 2, d)
296  h_1 = (h_0*self.number) % d
297  h_2 = ((arith1.inverse(2,d)*h_0*(self.number - h_1**2))//d) % d
298  b = (h_1 + h_2*d) % a
299  if b& 1 == 0:
300  b = b - a
301 
302  self.d_list.append(d)
303  self.a_list.append(a)
304  self.b_list.append(b)
305 
306  # Get solution of F(x) = 0 (mod p^i)
307  solution = []
308  i = 0
309  for s in self.Nsqrt:
310  k = 0
311  p_solution = []
312  ppow = 1
313  while k < len(s):
314  ppow *= self.FB[i+2]
315  a_inverse = arith1.inverse(2*self.a_list[-1], ppow)
316  x_1 = ((-b + s[k][0])*a_inverse) % ppow
317  x_2 = ((-b + s[k][1])*a_inverse) % ppow
318  p_solution.append([x_1, x_2])
319  k += 1
320  i += 1
321  solution.append(p_solution)
322  self.solution = solution
323  self.coefficienttime += time.time() - T
324 
325  def run_sieve(self):
326  self.make_poly()
327  T = time.time()
328  M = self.Srange
329  a = self.a_list[-1] #
330  b = self.b_list[-1] # These are coefficients of F(x)
331  c = (b**2-self.number)//(4*a) #
332  d = self.d_list[-1] #
333 
334  self.poly_table = [] # This is F(x) value , x in [-M,M].
335  self.log_poly = [] # This is log(F(x)) value.
336  self.minus_check = [] # This is "x" that F(x) is minus value.
337  for j in self.move_range:
338  jj = (a * j + b) * j + c
339  if jj < 0:
340  jj = -jj
341  self.minus_check.append(j + M)
342  elif jj == 0:
343  jj = 1
344  lj = int((math.log(jj)*30)*0.95) # 0.95 is an erroe
345  self.poly_table.append(jj)
346  self.log_poly.append(lj)
347 
348  y = arith1.inverse(2*d, self.number)
349  start_location = []
350  logp = [0]*(2*M+1)
351  j = 2
352  for i in self.solution:
353  start_p = []
354  ppow = 1
355  for k in range(len(i)):
356  ppow *= self.FB[j]
357  q = -M // ppow
358  s_1 = (q + 1) * ppow + i[k][0]
359  s_2 = (q + 1) * ppow + i[k][1]
360  while s_1 + M >= ppow:
361  s_1 = s_1 - ppow
362  while s_2 + M >= ppow:
363  s_2 = s_2 - ppow
364  start_p.append([s_1+M, s_2+M])
365  start_location.append(start_p)
366  j += 1
367  self.start_location = start_location
368 
369  i = self.poly_table[0] & 1
370  while i <= 2*M:
371  j = 1
372  while self.poly_table[i] % (2**(j+1)) == 0:
373  j += 1
374  logp[i] += self.FB_log[1]*j
375  i += 2
376  L = 2
377  for plocation in self.start_location:
378  for k in range(len(plocation)):
379  s_1 = plocation[k][0]
380  s_2 = plocation[k][1]
381  ppow = self.FB[L]**(k+1)
382  while s_1 <= 2*M:
383  logp[s_1] += self.FB_log[L]
384  s_1 += ppow
385  while s_2 <= 2*M:
386  logp[s_2] += self.FB_log[L]
387  s_2 += ppow
388  L += 1
389 
390  self.logp = logp
391  smooth = []
392  for t in range(2*M+1):
393  if logp[t] >= self.log_poly[t]:
394  poly_val = self.poly_table[t]
395  index_vector = []
396  H = (y*(2*a*(t-self.Srange)+b))%self.number
397  for p in self.FB:
398  if p == -1:
399  if t in self.minus_check:
400  index_vector.append(1)
401  else:
402  index_vector.append(0)
403  else:
404  r = 0
405  while poly_val % (p**(r+1)) == 0:
406  r += 1
407  v = r & 1
408  index_vector.append(v)
409  smooth.append([index_vector, (poly_val, H)])
410  self.sievingtime += time.time() - T
411  return smooth
412 
413  def get_vector(self):
414  T = time.time()
415  P = len(self.FB)
416  if P < 100:
417  V = -5
418  else:
419  V = 0
420  smooth = []
421  i = 0
422  while P*1 > V:
423  n = self.run_sieve()
424  V += len(n)
425  smooth += n
426  i += 1
427  if i % 50 == 0:
428  if i == 50:
429  _log.info("*/ Per 50 times changing poly report /* ")
430  _log.info("Time of deciding coefficient = %f sec" % self.coefficienttime)
431  _log.info("Sieving Time = %f sec" % self.sievingtime)
432  _log.info("Find smooth numbers are %d / %d" % (V, P))
433  if P < 100:
434  V += 5
435  self.smooth = smooth
436  _log.info("*/ Total %d times changing poly report /*" % i)
437  _log.info("Time of deciding coefficient = %f sec" % self.coefficienttime)
438  _log.info("Sieving Time = %f sec" % self.sievingtime)
439  _log.info("Found smooth numbers are %d / %d" % (V, P))
440  _log.info("Total time of getting enough smooth numbers = %f sec" % (time.time()-T))
441  return smooth
442 
443 
444 class Elimination(object):
445  def __init__(self, smooth):
446  self.vector = []
447  self.history = []
448  i = 0
449  for vec in smooth:
450  self.vector.append(vec[0])
451  self.history.append({i:1})
452  i += 1
453  self.FB_number = len(self.vector[0])
454  self.row_size = len(self.vector)
455  self.historytime = 0
456 
457  def vector_add(self, i, j):
458  V_i = self.vector[i]
459  V_j = self.vector[j]
460  k = 0
461  while k < len(V_i):
462  if V_i[k] == 1:
463  if V_j[k] == 1:
464  V_j[k] = 0
465  else:
466  V_j[k] = 1
467  k += 1
468 
469  def transpose(self):
470  Transe_vector = []
471  i = 0
472  while i < self.FB_number:
473  j = 0
474  vector = []
475  while j < self.row_size:
476  vector.append(self.vector[j][i])
477  j += 1
478  Transe_vector.append(vector)
479  i += 1
480  self.Transe = Transe_vector
481 
482  def history_add(self, i, j):
483  T = time.time()
484  H_i = self.history[i].keys()
485  H_j = self.history[j].keys()
486  for k in H_i:
487  if k in H_j:
488  del self.history[j][k]
489  else:
490  self.history[j][k] = 1
491  self.historytime += time.time() - T
492 
493  def gaussian(self):
494  T = time.time()
495  pivot = []
496  FBnum = self.FB_number
497  Smooth = len(self.vector)
498  for j in range(self.FB_number):
499  for k in range(Smooth):
500  if k in pivot or not self.vector[k][j]:
501  continue
502  pivot.append(k)
503  V_k = self.vector[k]
504  for h in range(Smooth):
505  if h in pivot or not self.vector[h][j]:
506  continue
507  self.history_add(k, h)
508  V_h = self.vector[h]
509  for q in range(j, FBnum):
510  if V_k[q]:
511  V_h[q] = not V_h[q]
512  break
513  self.pivot = pivot
514 
515  zero_vector = []
516  for check in range(Smooth):
517  if check not in pivot:
518  g = 0
519  while g < FBnum:
520  if self.vector[check][g] == 1:
521  break
522  g += 1
523  if g == FBnum:
524  zero_vector.append(check)
525  _log.info("Time of Gaussian Elimination = %f sec" % (time.time() - T))
526  return zero_vector
527 
528 def qs(n, s, f):
529  """
530  This is main function of QS
531  Arguments are (composite_number, sieve_range, factorbase_size)
532  You must input these 3 arguments.
533  """
534  a = time.time()
535  Q = QS(n, s, f)
536  _log.info("Sieve range is [ %d , %d ] , Factorbase size = %d , Max Factorbase %d" % (Q.move_range[0], Q.move_range[-1], len(Q.FB), Q.maxFB))
537  Q.run_sieve()
538  V = Elimination(Q.smooth)
539  A = V.gaussian()
540  _log.info("Found %d linearly dependent relations" % len(A))
541  answerX_Y = []
542  N_factors = []
543  for i in A:
544  B = V.history[i].keys()
545  X = 1
546  Y = 1
547  for j in B:
548  X *= Q.smooth[j][1][0]
549  Y *= Q.smooth[j][1][1]
550  Y = Y % Q.number
551  X = sqrt_modn(X, Q.number)
552  answerX_Y.append(X-Y)
553  for k in answerX_Y:
554  if k != 0:
555  factor = gcd.gcd(k, Q.number)
556  if factor not in N_factors and factor != 1 and \
557  factor != Q.number and prime.primeq(factor) == 1:
558  N_factors.append(factor)
559  N_factors.sort()
560  _log.info("Total time = %f sec" % (time.time()-a))
561  _log.info(str(N_factors))
562 
563 def mpqs(n, s=0, f=0, m=0):
564  """
565  This is main function of MPQS.
566  Arguments are (composite_number, sieve_range, factorbase_size, multiplier)
567  You must input composite_number at least.
568  """
569  T = time.time()
570  M = MPQS(n, s, f, m)
571  _log.info("Sieve range is [ %d , %d ] , Factorbase size = %d , Max Factorbase %d" % (M.move_range[0], M.move_range[-1], len(M.FB), M.maxFB))
572  M.get_vector()
573  N = M.number // M.multiplier
574  V = Elimination(M.smooth)
575  A = V.gaussian()
576  _log.info("Found %d linerly dependent relations" % len(A))
577  answerX_Y = []
578  N_prime_factors = []
579  N_factors = []
580  output = []
581  for i in A:
582  B = V.history[i].keys()
583  X = 1
584  Y = 1
585  for j in B:
586  X *= M.smooth[j][1][0]
587  Y *= M.smooth[j][1][1]
588  Y = Y % M.number
589  X = sqrt_modn(X, M.number)
590  if X != Y:
591  answerX_Y.append(X-Y)
592  NN = 1
593  for k in answerX_Y:
594  factor = gcd.gcd(k, N)
595  if factor not in N_factors and factor != 1 and factor != N \
596  and factor not in N_prime_factors:
597  if prime.primeq(factor):
598  NN = NN*factor
599  N_prime_factors.append(factor)
600  else:
601  N_factors.append(factor)
602 
603  _log.info("Total time = %f sec" % (time.time() - T))
604 
605  if NN == N:
606  _log.debug("Factored completely!")
607  N_prime_factors.sort()
608  for p in N_prime_factors:
609  N = N // p
610  i = arith1.vp(N, p, 1)[0]
611  output.append((p, i))
612  return output
613  elif NN != 1:
614  f = N // NN
615  if prime.primeq(f):
616  N_prime_factors.append(f)
617  _log.debug("Factored completely !")
618  N_prime_factors.sort()
619  for p in N_prime_factors:
620  N = N // p
621  i = arith1.vp(N, p, 1)[0]
622  output.append((p, i))
623  return output
624 
625  for F in N_factors:
626  for FF in N_factors:
627  if F != FF:
628  Q = gcd.gcd(F, FF)
629  if prime.primeq(Q) and Q not in N_prime_factors:
630  N_prime_factors.append(Q)
631  NN = NN*Q
632 
633  N_prime_factors.sort()
634  for P in N_prime_factors:
635  i, N = arith1.vp(N, P)
636  output.append((P, i))
637 
638  if N == 1:
639  _log.debug("Factored completely!! ")
640  return output
641 
642  for F in N_factors:
643  g = gcd.gcd(N, F)
644  if prime.primeq(g):
645  N_prime_factors.append(g)
646  N = N // g
647  i = arith1.vp(N, g, 1)[0]
648  output.append((g, i))
649  if N == 1:
650  _log.debug("Factored completely !! ")
651  return output
652  elif prime.primeq(N):
653  output.append((N, 1))
654  _log.debug("Factored completely!!! ")
655  return output
656  else:
657  N_factors.sort()
658  _log.error("Sorry, not factored completely")
659  return output, N_factors
660 
661 
662 
663 
668 
670  if n < 2:
671  raise ValueError("Input value must be bigger than 1")
672  else:
673  primes = [2]
674  sieves = [1] * (((n+1) >> 1) - 1)
675  k = 3
676  i = 0
677  sieves_len = len(sieves)
678  while k <= math.sqrt(n):
679  if sieves[i] == 1:
680  j = 1
681  while i+j*k <= sieves_len-1:
682  sieves[i+j*k] = 0
683  j = j+1
684  k, i = k+2, i+1
685  else:
686  k, i = k+2, i+1
687  y = 3
688  for x in sieves:
689  if x == 1:
690  primes.append(y)
691  y = y+2
692  return primes
693 
694 def prime_mod8(n):
695  """
696  This is table for choosing multiplier which makes N to have
697  factorbase(2,3,5,7,11,13)
698  """
699  primes = eratosthenes(n)
700  PrimeList = {1:[], 3:[], 5:[], 7:[]}
701  LegendreList = {1:[], 3:[], 5:[], 7:[]}
702  sp = [2, 3, 5, 7, 11, 13]
703  for p in primes:
704  if p not in sp:
705  leg = [arith1.legendre(p, q) for q in sp[1:]]
706  if leg not in PrimeList[p & 7]:
707  LegendreList[p & 7].append(leg)
708  PrimeList[p & 7].append([p, leg])
709  return [PrimeList[1], PrimeList[3], PrimeList[5], PrimeList[7]]
710 
712  primes = eratosthenes(n)
713  primes_log = []
714  for i in primes:
715  l = int(math.log(i)*30) #30 is scale
716  primes_log.append(l)
717  return primes_log
718 
719 def sqrt_modn(n, modulo):
720  import nzmath.factor.methods as methods
721  factorOfN = methods.trialDivision(n)
722  prod = 1
723  for p, e in factorOfN:
724  prod = (prod * pow(p, e >> 1, modulo)) % modulo
725  return prod
726 
727 def sqroot_power(a, p, n):
728  """
729  return squareroot of a mod p^k for k = 2,3,...,n
730  """
731  x = arith1.modsqrt(a, p)
732  answer = [[x, p - x]]
733  ppower = p
734  inverse = arith1.inverse(x << 1, p)
735  for i in range(n - 1):
736  x += (a - x ** 2) // ppower * inverse % p * ppower
737  ppower *= p
738  answer.append([x, ppower - x])
739  return answer
740 
741 
744 primes_table = eratosthenes(10**5)
745 primes_log_table = eratosthenes_log(10**5)
746 prime_8 = prime_mod8(8090)
747 parameters_for_mpqs = [[100,20], #9 -digits
748  [100,21], #10
749  [100,22], #11
750  [100,24], #12
751  [100,26], #13
752  [100,29], #14
753  [100,32], #15
754  [200,35], #16
755  [300,40], #17
756  [300,60], #18
757  [300,80], #19
758  [300,100], #20
759  [300,100], #21
760  [300,120], #22
761  [300,140], #23
762  [600,160], #24
763  [900,180], #25
764  [1200,200], #26
765  [1000,220], #27
766  [2000,240], #28
767  [2000,260], #29
768  [2000,325], #30
769  [2000,355], #31
770  [2000,375], #32
771  [3000,400], #33
772  [2000,425], #34
773  [2000,550], #35
774  [3000,650], #36
775  [5000,750], #37
776  [4000,850], #38
777  [4000,950], #39
778  [5000,1000], #40
779  [14000,1150], #41
780  [15000,1300], #42
781  [15000,1600], #43
782  [15000,1900], #44
783  [15000,2200], #45
784  [20000,2500], #46
785  [25000,2500], #47
786  [27500,2700], #48
787  [30000,2800], #49
788  [35000,2900], #50
789  [40000,3000], #51
790  [50000,3200], #52
791  [50000,3500]] #53
792 
793 
794 
797 
798 def mpqsfind(n, s=0, f=0, m=0, verbose=False):
799  """
800  This is main function of MPQS.
801  Arguments are (composite_number, sieve_range, factorbase_size, multiplier)
802  You must input composite_number at least.
803  """
804  # verbosity
805  if verbose:
806  _log.setLevel(logging.DEBUG)
807  _log.debug("verbose")
808  else:
809  _log.setLevel(logging.NOTSET)
810 
811  starttime = time.time()
812  M = MPQS(n, s, f, m)
813  _log.info("Sieve range is [%d, %d]" % (M.move_range[0], M.move_range[-1]))
814  _log.info("Factorbase size = %d, Max Factorbase %d" % (len(M.FB), M.maxFB))
815  M.get_vector()
816  N = M.number // M.multiplier
817  V = Elimination(M.smooth)
818  A = V.gaussian()
819  _log.info("Found %d linearly dependent relations" % len(A))
820  differences = []
821  for i in A:
822  B = V.history[i].keys()
823  X = 1
824  Y = 1
825  for j in B:
826  X *= M.smooth[j][1][0]
827  Y *= M.smooth[j][1][1]
828  Y = Y % M.number
829  X = arith1.floorsqrt(X) % M.number
830  if X != Y:
831  differences.append(X-Y)
832 
833  for diff in differences:
834  divisor = gcd.gcd(diff, N)
835  if 1 < divisor < N:
836  _log.info("Total time = %f sec" % (time.time() - starttime))
837  return divisor
nzmath.factor.mpqs.MPQS.Nsqrt
Nsqrt
Definition: mpqs.py:271
nzmath.bigrange.range
def range(start, stop=None, step=None)
Definition: bigrange.py:19
nzmath.factor.find
Definition: find.py:1
nzmath.factor.mpqs.QS.run_sieve
def run_sieve(self)
Definition: mpqs.py:76
nzmath.factor.mpqs.prime_mod8
def prime_mod8(n)
Definition: mpqs.py:694
nzmath.factor.mpqs.MPQS.logp
logp
Definition: mpqs.py:390
nzmath.factor.mpqs.MPQS.digit
digit
Definition: mpqs.py:190
nzmath.factor.mpqs.MPQS.sqrt_state
sqrt_state
Definition: mpqs.py:220
nzmath.factor.mpqs.qs
def qs(n, s, f)
Definition: mpqs.py:528
nzmath.factor.mpqs.MPQS.solution
solution
Definition: mpqs.py:322
nzmath.factor.mpqs.MPQS.run_sieve
def run_sieve(self)
Definition: mpqs.py:325
nzmath.factor.mpqs.QS.log_poly
log_poly
Definition: mpqs.py:73
nzmath.factor.mpqs.Elimination.vector
vector
Definition: mpqs.py:446
nzmath.factor.mpqs.Elimination.__init__
def __init__(self, smooth)
Definition: mpqs.py:445
nzmath.gcd
Definition: gcd.py:1
nzmath.factor.mpqs.MPQS.move_range
move_range
Definition: mpqs.py:216
nzmath.factor.mpqs.QS.Srange
Srange
Definition: mpqs.py:27
nzmath.factor.mpqs.QS.move_range
move_range
Definition: mpqs.py:30
nzmath.factor.mpqs.mpqs
def mpqs(n, s=0, f=0, m=0)
Definition: mpqs.py:563
nzmath.factor.mpqs.Elimination.history
history
Definition: mpqs.py:447
nzmath.factor.mpqs.MPQS.b_list
b_list
Definition: mpqs.py:187
nzmath.factor.mpqs.QS.logp
logp
Definition: mpqs.py:147
nzmath.factor.mpqs.QS.maxFB
maxFB
Definition: mpqs.py:50
nzmath.factor.mpqs.MPQS.minus_check
minus_check
Definition: mpqs.py:336
nzmath.factor.mpqs.Elimination.historytime
historytime
Definition: mpqs.py:455
nzmath.factor.mpqs.Elimination
Definition: mpqs.py:444
nzmath.factor.mpqs.MPQS.maxFB
maxFB
Definition: mpqs.py:262
nzmath.factor.mpqs.QS.FB
FB
Definition: mpqs.py:48
nzmath.factor.mpqs.MPQS.FBN
FBN
Definition: mpqs.py:195
nzmath.factor.mpqs.QS.digit
digit
Definition: mpqs.py:26
nzmath.factor.mpqs.MPQS.__init__
def __init__(self, n, sieverange=0, factorbase=0, multiplier=0)
Definition: mpqs.py:173
nzmath.factor.mpqs.QS.minus_check
minus_check
Definition: mpqs.py:74
nzmath.factor.mpqs.MPQS.get_vector
def get_vector(self)
Definition: mpqs.py:413
nzmath.factor.mpqs.Elimination.history_add
def history_add(self, i, j)
Definition: mpqs.py:482
nzmath.factor.mpqs.MPQS.poly_table
poly_table
Definition: mpqs.py:334
nzmath.factor.mpqs.QS.start_location
start_location
Definition: mpqs.py:106
nzmath.factor.mpqs.QS.poly_table
poly_table
Definition: mpqs.py:72
nzmath.factor.mpqs.mpqsfind
def mpqsfind(n, s=0, f=0, m=0, verbose=False)
only find a factor
Definition: mpqs.py:798
nzmath.factor.mpqs.MPQS.log_poly
log_poly
Definition: mpqs.py:335
nzmath.factor.mpqs.Elimination.gaussian
def gaussian(self)
Definition: mpqs.py:493
nzmath.factor.mpqs.Elimination.FB_number
FB_number
Definition: mpqs.py:453
nzmath.factor.mpqs.Elimination.pivot
pivot
Definition: mpqs.py:513
nzmath.factor.mpqs.MPQS.start_location
start_location
Definition: mpqs.py:367
nzmath.factor.mpqs.Elimination.row_size
row_size
Definition: mpqs.py:454
nzmath.arith1
Definition: arith1.py:1
nzmath.factor.mpqs.MPQS.number
number
Definition: mpqs.py:174
nzmath.factor.mpqs.MPQS.FB_log
FB_log
Definition: mpqs.py:261
nzmath.factor.mpqs.QS.sqrt_n
sqrt_n
Definition: mpqs.py:21
nzmath.factor.mpqs.MPQS.smooth
smooth
Definition: mpqs.py:435
nzmath.factor.mpqs.sqrt_modn
def sqrt_modn(n, modulo)
Definition: mpqs.py:719
nzmath.factor.mpqs.QS.number
number
Definition: mpqs.py:20
nzmath.prime
Definition: prime.py:1
nzmath.factor.mpqs.MPQS.sievingtime
sievingtime
Definition: mpqs.py:183
nzmath.factor.mpqs.eratosthenes
def eratosthenes(n)
Definition: mpqs.py:669
nzmath.factor.mpqs.QS.solution
solution
Definition: mpqs.py:57
nzmath.factor.mpqs.MPQS.multiplier
multiplier
Definition: mpqs.py:242
nzmath.factor.mpqs.MPQS.a_list
a_list
Definition: mpqs.py:186
nzmath.factor.mpqs.MPQS.make_poly
def make_poly(self)
Definition: mpqs.py:273
nzmath.factor.mpqs.MPQS.coefficienttime
coefficienttime
Definition: mpqs.py:184
nzmath.factor.mpqs.eratosthenes_log
def eratosthenes_log(n)
Definition: mpqs.py:711
nzmath.factor.mpqs.sqroot_power
def sqroot_power(a, p, n)
Definition: mpqs.py:727
nzmath.factor.methods
Definition: methods.py:1
nzmath.factor.mpqs.QS.smooth
smooth
Definition: mpqs.py:168
nzmath.factor.mpqs.MPQS.d_list
d_list
Definition: mpqs.py:185
nzmath.factor.mpqs.Elimination.transpose
def transpose(self)
Definition: mpqs.py:469
nzmath.factor.mpqs.Elimination.Transe
Transe
Definition: mpqs.py:480
nzmath.factor.mpqs.MPQS.Srange
Srange
Definition: mpqs.py:193
nzmath.factor.mpqs.MPQS.FB
FB
Definition: mpqs.py:260
nzmath.factor.mpqs.QS
Definition: mpqs.py:18
nzmath.factor.mpqs.QS.__init__
def __init__(self, n, sieverange, factorbase)
Definition: mpqs.py:19
nzmath.factor.mpqs.MPQS
Definition: mpqs.py:172
nzmath.factor.mpqs.QS.FB_log
FB_log
Definition: mpqs.py:49
nzmath.factor.mpqs.Elimination.vector_add
def vector_add(self, i, j)
Definition: mpqs.py:457
nzmath.factor.mpqs.QS.FBN
FBN
Definition: mpqs.py:28