"Fossies" - the Fresh Open Source Software Archive

Member "NZMATH-1.2.0/nzmath/factor/find.py" (19 Nov 2012, 3740 Bytes) of package /linux/misc/old/NZMATH-1.2.0.tar.gz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) Python source code syntax highlighting (style: standard) with prefixed line numbers. Alternatively you can here view or download the uninterpreted source code file. For more information about "find.py" see the Fossies "Dox" file reference documentation.

    1 """
    2 All methods defined here return one of a factor of given integer.
    3 When 1 is returned, the method has failed to factor,
    4 but 1 is a factor anyway.
    5 
    6 'verbose' boolean flag can be specified for verbose reports.
    7 """
    8 
    9 import logging
   10 import nzmath.arith1 as arith1
   11 import nzmath.bigrandom as bigrandom
   12 import nzmath.gcd as gcd
   13 import nzmath.prime as prime
   14 
   15 _log = logging.getLogger('nzmath.factor.find')
   16 
   17 # Pollard's rho method
   18 def rhomethod(n, **options):
   19     """
   20     Find a non-trivial factor of n using Pollard's rho algorithm.
   21     The implementation refers the explanation in C.Pomerance's book.
   22     """
   23     # verbosity
   24     verbose = options.get('verbose', False)
   25     if not verbose:
   26         _silence()
   27 
   28     if n <= 3:
   29         return 1
   30 
   31     g = n
   32     while g == n:
   33         # x^2 + a is iterated. Starting value x = u.
   34         a = bigrandom.randrange(1, n-2)
   35         u = v = bigrandom.randrange(0, n-1)
   36         _log.info("%d %d" % (a, u))
   37         g = gcd.gcd((v**2 + v + a) % n - u, n)
   38         while g == 1:
   39             u = (u**2 + a) % n
   40             v = ((pow(v, 2, n) + a)**2 + a) % n
   41             g = gcd.gcd(v - u, n)
   42     if not verbose:
   43         _verbose()
   44     return g
   45 
   46 # p-1 method
   47 def pmom(n, **options):
   48     """
   49     This function tries to find a non-trivial factor of n using
   50     Algorithm 8.8.2 (p-1 first stage) of Cohen's book.
   51     In case of N = pow(2,i), this program will not terminate.
   52     """
   53     # verbosity
   54     verbose = options.get('verbose', False)
   55     if not verbose:
   56         _silence()
   57 
   58     # initialize
   59     x = y = 2
   60     primes = []
   61     if 'B' in options:
   62         B = options['B']
   63     else:
   64         B = 10000
   65 
   66     for q in prime.generator():
   67         primes.append(q)
   68         if q > B:
   69             if gcd.gcd(x-1, n) == 1:
   70                 if not verbose:
   71                     _verbose()
   72                 return 1
   73             x = y
   74             break
   75         q1 = q
   76         l = B // q
   77         while q1 <= l:
   78             q1 *= q
   79         x = pow(x, q1, n)
   80         if len(primes) >= 20:
   81             if gcd.gcd(x-1, n) == 1:
   82                 primes, y = [], x
   83             else:
   84                 x = y
   85                 break
   86 
   87     for q in primes:
   88         q1 = q
   89         while q1 <= B:
   90             x = pow(x, q, n)
   91             g = gcd.gcd(x-1, n)
   92             if g != 1:
   93                 if not verbose:
   94                     _verbose()
   95                 if g == n:
   96                     return 1
   97                 return g
   98             q1 *= q
   99 
  100 def trialDivision(n, **options):
  101     """
  102     Return a factor of given integer by trial division.
  103 
  104     options can be either:
  105     1) 'start' and 'stop' as range parameters.
  106     2) 'iterator' as an iterator of primes.
  107     If both options are not given, prime factor is searched from 2
  108     to the square root of the given integer.
  109     """
  110     # verbosity
  111     verbose = options.get('verbose', False)
  112     if not verbose:
  113         _silence()
  114 
  115     if 'start' in options and 'stop' in options:
  116         if 'step' in options:
  117             trials = range(options['start'], options['stop'], options['step'])
  118         else:
  119             trials = range(options['start'], options['stop'])
  120     elif 'iterator' in options:
  121         trials = options['iterator']
  122     elif n < 1000000:
  123         trials = prime.generator_eratosthenes(arith1.floorsqrt(n))
  124     else:
  125         trials = prime.generator()
  126 
  127     limit = arith1.floorsqrt(n)
  128     for p in trials:
  129         if limit < p:
  130             break
  131         if 0 == n % p:
  132             if not verbose:
  133                 _verbose()
  134             return p
  135     if not verbose:
  136         _verbose()
  137     return 1
  138 
  139 def _silence():
  140     """
  141     Stop verbose outputs.
  142     """
  143     _log.setLevel(logging.NOTSET)
  144 
  145 def _verbose():
  146     """
  147     Stop silencing.
  148     """
  149     _log.setLevel(logging.DEBUG)