Sharedwww / wiki / ant07(2f)projects(2f)shumow_raw.htmlOpen in CoCalc
Author: William A. Stein
ant07/projects/shumow raw
 [top] [TitleIndex] [WordIndex]

ant07/projects/shumow raw

Shumow -- GNFS -- raw notebook code

To use this in Sage, click edit in this page and edit in a sage worksheet and copy the resulting plain text into the sage worksheet.

GNFS Term Project system:sage

id=0| 
# 
# This is a Toy implementation of the General Number Field Sieve For Factoring
# 
# Author:
# 
# Math 581F Fall 2007
# Professor William Stein
# Final Project 
#
# "You haven't lived till you've sieved."
#

#  The following are the sources I used to learn about the General Number Field Sieve

#  Matthew E Briggs Master's Thesis from Virginia Polytech (This is probably the most clear and detailed)
#  The Lenstra et al book "The development of the number field sieve"
#  Pomerance's "Tale of Two Sieves"
#  Henri Cohn's "A Course in Computational Algebraic Number Theory"
#  The Math 581F Textbook / Lecture Notes (William Stein)
#  
#  Personal conversations with Peter Montgomery

id=61| 
#
# The General Number Field Sieve works like this
# 
# 1) Decide to use GNFS to factor N, determine algorithm parameters based on N
#    such as, the size of the rational factor base, size of the algebraic factor base
#    number of quadratic characters, degree d, and polynomial f to use
#
# 2) calculate the rational prime factor base (primes up to bound)
# 3) calculate the algebraic factor base (first degree prime ideals up to bound)
# 4) calculate the quadratic characters (first few 1st degree prime ideals beyond the algebraic factor base.)
# 5) use sieving to determine integer pairs (a,b) that have values
#    have a+b*m (m a root of f mod N) and N(a+b*theta) (theta complex root of f)
#    smooth with respect to both the algebraic and rational prime factor bases.
#    while sieving, also determine the degree mod 2 over each prime (rational and algebraic)
#    and save this information.
#    Also, record the power (mod 2) of -1 determining the legendre symbol (a+b*s / q) of the
#    first degree prime ideal (s,q) in the quadratic character base.
# 6) Once there are enough such primes more than the sum of the lengths of the bases
#    use linear algebra to determine the kernel of the transformation.
# 7) Use the result of step 6 to calculate a product of (a+b*theta)'s =  Beta^2 in the number field
# 8) Use the CRT to determine Beta from Beta^2 in the number field
# 9) replace theta in Beta^2 with m to obtain v^2 = phi(m)
# 10) calculate u^2 = product of (a+b*m)'s mod N, calculate u from (a+b*m)
# 11) u^2 - v^2 = 0 mod N, so use difference of squares trick to try to factor N.
# 12) repeat until N is factored.
#

id=67| 
#
# Major Weaknesses of this project:
# The major weaknesses of this project are that I relied a lot on sage
# to do math.  Although it served its purpose in teaching me about the GNFS.
#
# WARNING!!!
# On My Laptop This Notebook can not handle in the mid to high 8 digits.
#
# Other than that:
# 1) To calculate square roots, I should use the CRT, but instead
#    I let sage functions do the work for me.
#    I did this because, the CRT is a technical detail at the end of the
#    algorithm, and also, I understand that math already.
#    I wanted to focus on learning the stuff I didn't understand.
# 2) I rely on sage for all linear algebra, this wouldn't work at limits where
#    10s or 100s of thousands of elements may be in the factor base.
# 3) This program uses is_prime_power() to try and determine if an N is a prime power
#    (because the GNFS does not work well on prime powers.)
#    This function is actually the limit on sizes of N that can be factored.
# 4) If a number is a prime power, I use the ECM to factor it, in actuality
#    an algorithm would be selected heuristically.
#

id=62| 
#
# Global Flags: These can be set to change the behavior of the sieving algorithm
# Used for sort of anecdotal observation and experimentation of the algorithm
#

# if this value is set to false, the algorithm will not use quadratic characters
use_quadratic_characters = True;

id=1| 
#
# This cell implements simple heuristics for the
# parameters used in the GNFS based on given N
#

rational_factor_base_bound_ratio = 5;    # this factor seemed to work well
algebraic_factor_base_bound_ratio = 10;  # so did this
quadratic_character_base_ratio = 0.10;   # as a percentage of the algebraic factor base

degree = 3;     # hard coded in this toy implementation

a_bound_ratio = 0.004;

def determine_GNFS_bounds(N):
    B1 = max(10, ceil(rational_factor_base_bound_ratio*log(N)));
    B2 = max(20, ceil(algebraic_factor_base_bound_ratio*log(N)));

    if (use_quadratic_characters):
        Quad_character_count = max(10, ceil((1.0+quadratic_character_base_ratio)*B2))
    else:
        Quad_character_count = 0

    B3 = B2 + Quad_character_count;

    a_bound = min(max(50, ceil(a_bound_ratio * N)), (2^16 - 1) );

    print 'Parameter Selection for GNFS:'

    print 'upper bound on rational prime base (all primes under):', B1
    print 'upper bound on algebraic prime base (all first degree primes under):', B2
    print 'upper bound on quadratic character prime base (all subsequent primes under):', B3

    print 'search bounds for unit term: [', -a_bound, ',', a_bound, ']'

    return (B1, B2, B3, degree, a_bound);

id=65| 
#
# Note about the representation of prime base elements
# rational primes are represented as (m mod p, p) p prime in QQ
# algebraic primes are represented as (r, p) p prime in QQ, r root of f mod p.
# storing in this form helps sieving and allows common code to be used
# for both
#

id=2| 
# This cell defines the function that
# generates the rational factor base
# (rational primes for sieving)

def generate_rational_factor_base(m, Bound):

    rational_factor_base = [];

    for p in xrange(1,Bound):
        if(is_prime(p)):
            rational_factor_base.append((m % p, p));

    return rational_factor_base;

id=3| 
# This cell generate the algebraic factor base
# (first degree prime ideals of order Z[theta] in number field)

# helper functions:

# uses sage to determine roots of f mod p
def get_roots_mod_p(f, p):
    return f.change_ring(GF(p)).roots(None, False);

# generates the algebraic factor base
def generate_algebraic_factor_base(f, Bound):

    algebraic_factor_base = [];

    for p in xrange(Bound):
        if(is_prime(p)):
            roots_fmodp = get_roots_mod_p(f, p);

            for j in xrange(len(roots_fmodp)):
                algebraic_factor_base.append((roots_fmodp[j].lift(),p))
                
    return algebraic_factor_base;
# generates quadratic characters
def generate_quadratic_character_base(f, Lower_Bound, Upper_Bound):

    quadratic_character_base = [];

    for p in xrange(Lower_Bound, Upper_Bound):
        if(is_prime(p)):
            roots_fmodp = get_roots_mod_p(f, p);

            for j in xrange(len(roots_fmodp)):
                quadratic_character_base.append((roots_fmodp[j].lift(),p))
                
    return quadratic_character_base;

id=4| 
#
# this cell defines helper functions for polynomials
# 
#

#
# determines m-arry expansion of N
# and the coefficients of powers of m
def get_polynomial_coefficients(N, m):

    coefficient_list = [];

    temp = N;

    while(0 != temp):
        r = temp % m;
        temp = (temp - r) / m;
        coefficient_list.append(r);    

    return coefficient_list;
    

#
# function that tries to generate a
# monic, irreducible polynomial f with root m mod N
# with degree d.
# 
def generate_polynomial(polynomial_ring, N, d):
    # determine and lower and upper bounds for m

    q = 0;

    while (1 != q):

        u = floor(pow(N, 1.0/d));
        l = ceil(pow(N, (2.0*d - 3.0)/(2.0*d*(d-1))));
    
        m = floor(round(uniform(l, u)));

        (q,r) = (N).quo_rem(m^d);
    

    f = polynomial_ring(get_polynomial_coefficients(N, m));
   
    return (f, m)

#
# checks if a polynomial is monic and irreducible
#
def is_good_polynomial(f):
    return (f.is_monic() and f.is_irreducible())

#
# retries a bunch of times to generate a suitable random polynomial
#
def get_polynomial(polynomial_ring, N, d, retry_count):

    found_polynomial = False;

    for j in xrange(retry_count):
        (f, m) = generate_polynomial(polynomial_ring, N, d);
        if is_good_polynomial(f):
            found_polynomial = True;
            break;

    if (found_polynomial):
        return (f, m)
    else:
        return (0,0)

id=18| 
#
#    This cell contains the sieving code
#

quo_rem = Integer.quo_rem

#
#    Common code that will take a factor base and sieve it
#
def sieve_interval(interval_to_sieve, a_lowerbound, b, factor_base):
    
    factor_base_count = len(factor_base);
    
    exponent_lists = [ [0 for j in xrange(factor_base_count+1)] for k in xrange(len(interval_to_sieve))]

    #
    #    determine the sign exponents of the elements in each interval
    #    

    for j in xrange(len(interval_to_sieve)):
        if (interval_to_sieve[j] < 0):
            exponent_lists[j][factor_base_count] = 1;
            interval_to_sieve[j] = abs(interval_to_sieve[j]);


    #
    #    do the factor base sieving
    #

    for j in xrange(factor_base_count):
        p = factor_base[j][1];
        bt = -b*factor_base[j][0] % p;
        
        #
        #    Much of the power of the sieving algorithms
        #    come from the regular access pattern and
        #    advancing through the array p elements at a time
        #

        # in this case we have our factor bases stored as a pair (t,p)
        # such that the value in the sieve in that position
        # (for the rational base this is a+b*m and for the algebraic base it is N(a+b*theta)
        # is only divisible by p if a = -b*t + k*p, as such we
        # use the lower bound of the a range to determine the first
        # position in the array to start sieving at advance p positions from there

        for k in xrange((bt  - a_lowerbound) % p, len(interval_to_sieve), p):

            #
            # This loop goes through and calculates the largest integer e
            # such that p^e divides the array entry
            #
            
            e = 0;
            r = 0;

            while (0 == r):
                if (interval_to_sieve[k] != 0):
                    (q,r) = quo_rem(interval_to_sieve[k] , p);
                    if(0 == r):
                        r = 0;
                        e = e + 1;
                        interval_to_sieve[k] = q;
                    else:
                        r = 1;
                else:
                    r = 1;

            exponent_lists[k][j] = e;

    return (interval_to_sieve, exponent_lists)

#
# determines the quadratic character exponent vector
#
def determine_quad_char_expons(a,b, quad_char_base):

    qc_expons = [0 for j in xrange(len(quad_char_base))];

    for j in xrange(len(quad_char_base)):

        q = quad_char_base[j][1];
        s = quad_char_base[j][0];        

        qc_expons[j] = Integer((1 - legendre_symbol(a+s*b, q))/2);

    return qc_expons;

#
#   function: search_for_smooth_pairs
#   input: 6 parameters:
#   lower bound of a's
#   b
#   list of a+bm, sieved with respect to rational factor base
#   The exponent vector with respect to the rational factor base
#   list of N(a+b*theta) sieved with respect to the algebraic factor base
#   The exponent vector with respect to the algebraic factor base

#   output: a list of triplets of:
#   (a,b) pairs of integers that are smooth with respect to
#         both the algebraic and rational factor bases.
#   the exponent list with respect to the rational factor base
#   the exponent list with respect to the algebraic factor base
#

def search_for_smooth_pairs(a_lowerbound, b, rat_sieve, rat_exp_list, alg_sieve, alg_exp_list, quad_char_base):

    smooth_list = [];

    for j in xrange(len(rat_sieve)):
        if ( (1 == rat_sieve[j]) and (1 == alg_sieve[j]) ):
            
            qc_expons = determine_quad_char_expons(a_lowerbound+j, b, quad_char_base);

            smooth_list.append(((a_lowerbound+j,b),rat_exp_list[j],alg_exp_list[j], qc_expons));



    return smooth_list;


#
# Do the actual sieving
# returns a list of (a,b) pairs along with the rational, algebraic, and quadratic exponents
#

def perform_sieving(f, m, a_bound, rat_factor_base, alg_factor_base, quad_char_base):

    # determine the dimension of the space we are working over.    
    l1 = len(rat_factor_base);
    l2 = len(alg_factor_base);
    l3 = len(quad_char_base);

    l = l1 + l2 + l3 + 2; # add two for rational and algebraic sign exponents;

    d = f.degree()

    smooth_element_list = [];

    # loop until we have more smooth elements
    # than we have dimension

    b = 1;

    while(len(smooth_element_list) <= l):

        rat_interval = [];
        alg_interval = [];
    
        for j in xrange(2*a_bound + 1):

            a = -a_bound + j;

            rat_interval.append(a + b*m); # sieve the value a + b*m
            alg_interval.append(Integer(((-b)^d)*f(-a/b))); # sieve the value N(a+b*theta) = (-b)^d * f(-a/b)
        
        (rat_interval, rat_exp_list) = sieve_interval(rat_interval, -a_bound, b, rat_factor_base);
        
        (alg_interval, alg_exp_list) = sieve_interval(alg_interval, -a_bound, b, alg_factor_base);

        smooth_pairs_found = search_for_smooth_pairs(-a_bound, b, rat_interval, rat_exp_list, alg_interval, alg_exp_list, quad_char_base);
 
        smooth_element_list.extend(smooth_pairs_found);

        print 'found', len(smooth_element_list), 'simultaneously smooth pairs.'

        b = b+1;

    return smooth_element_list;

id=28| 

id=38| 
#
#    This cell contains a function that computes the product of the integer pairs
#    as an integer polynomial, and reduces it module the ideal generated by f
#

#
# computes product of linear polynomials mod Ideal generated by f
# (equivalent to doing math in number field
#
def compute_numberfield_product_of_pairs(polynomial_ring, f, integer_pairs, vector):

    I = polynomial_ring.ideal(f);

    product_polynomial = polynomial_ring([1]);

    for j in xrange(len(vector)):

        if (1 == vector[j]):

            linear_poly = polynomial_ring([integer_pairs[j][0], integer_pairs[j][1]]);
            product_polynomial = product_polynomial * linear_poly;
            product_polynomial.mod(I);

    return product_polynomial;

#
# computes the product of (a + b*m)'s mod N
# of the (a + b*m)'s corresponding to 1s in vector
#
def compute_integer_product_of_pairs(polynomial_ring, f, m, N, integer_pairs, vector):

    prod  = 1;

    for j in xrange(len(vector)):
        if (1 == vector[j]):
            prod  = prod*(integer_pairs[j][0] + m*integer_pairs[j][1]) % N;

    return prod;
    
#
# computes the difference of squares u^2 - v^2
# and returns (u+v) and (u-v)
#
def compute_difference_of_squares(polynomial_ring, f, m, N, integer_pairs, vector):

    found_squares = False;

    u_plus_v = 0; 
    u_minus_v = 0;

    ints_mod_N = ZZ.quo(ZZ.ideal(N));

    vsquared = ints_mod_N(compute_integer_product_of_pairs(polynomial_ring, f, m, N, integer_pairs, vector));

    if ( is_square(vsquared)):

        beta_squared = compute_numberfield_product_of_pairs(polynomial_ring, f, integer_pairs, vector);
        beta = find_square(f, beta_squared, polynomial_ring);
    
        if (False != beta):
        
            u = (beta(m)) % N;
            v = vsquared.sqrt().lift();
            u_plus_v = (u + v) % N; 
            u_minus_v = (u - v) % N;

            found_squares = True;

        else:
            print 'failed to find a square root in number field.'
    else:
        print 'integer was not square'

    return (found_squares, u_plus_v, u_minus_v);

id=7| 
#
#    Cheater Helper Functions
#
#    These are functions that I was too lame to implement myself
#    And Hence relied on sage's built in functionality.
#

ecm_instance = ECM();

def factor_prime_power(N):
    # The GNFS doesn't work one prime powers
    # so we use ecm to factor,
    # in reality, 
    # we would apply a heuristic on the best way to determine the prime and exponent    
    
    ecm_factors = ecm_instance.factor(N);

    prime_power_factorization = [(ecm_factors[0], len(ecm_factors))]

    return Factorization(prime_power_factorization);

#
#    This simple function does all the linear algebra (thank you sage!)
#
def find_kernel_vectors(M):
    return M.kernel().basis_matrix().rows();

#
# In reality this would be done with the CRT.
#
def find_square(f, beta_squared, polynomial_ring):
    # computes the square root of beta_squared
    # (a polynomial over ZZ of degree less than d, that represents the element in the order ZZ[theta])
    # and returns the polynomial representing the square root.
    
    K.<a> = NumberField(f);
    beta2 = beta_squared(a);
    if(not beta2.is_square()):
        print 'The supposedly found square is not square!'
        return False;
    return polynomial_ring(beta2.sqrt().list());

id=9| 
#
# This is the factoring algorithm
# It factors N, and returns a sage factorization object of N
#

def GNFS_Factoring_Algorithm(N):

    # determine if this is a good number to use the GNFS on
    # if it is not (i.e. it is a prime power)
    # call the specialized prime power factorization function

    if(is_prime_power(N)):

        print N, 'is a prime power, reducing to using ECM to factor'

        return factor_prime_power(N);

    # determine the algorithm bounds
    (B1, B2, B3, d, a_bound) = determine_GNFS_bounds(N);

    R.<x> = ZZ[];

    (f, m) = get_polynomial(R, N, d, 20);

    if((0 == f) and (0 == m)):
    
        # if the polynomial selection code did not find a degree d
        # polynomial in 20 tries, just use the built-in factoring
        # algorithm, because N is too small for the GNFS to work
        # at that point

        return factor(N);

    # get the rational factor base
    rational_factor_base = generate_rational_factor_base(m, B1);

    # get the algebraic factor base
    algebraic_factor_base = generate_algebraic_factor_base(f, B2);

    # get the quadratic character base
    quadratic_character_base = generate_quadratic_character_base(f, B2, B3);

    #
    # get a list of smooth elements
    #

    smooth_element_list = perform_sieving(f, m, a_bound, rational_factor_base, algebraic_factor_base, quadratic_character_base);

    #
    #    parse out the list of smooth elements into a list of integer pairs,
    #

    integer_pairs = [smooth_element_list[k][0] for k in xrange(len(smooth_element_list))];
    exponent_vectors = [];

    #
    # use the smooth element list vectors to make a matrix
    #

    for k in xrange(len(smooth_element_list)):
        exponent_vectors.extend(smooth_element_list[k][1]);
        exponent_vectors.extend(smooth_element_list[k][2]);
        exponent_vectors.extend(smooth_element_list[k][3]);

    num_expons = len(smooth_element_list[0][1]) + len(smooth_element_list[0][2]) + len(smooth_element_list[0][3]);

    M = matrix(GF(2), len(smooth_element_list), exponent_vectors )
    
    solutions = find_kernel_vectors(M);

    nontrivial_first_factor = False;
    nontrivial_second_factor = False;

    k = 0

    #
    #    Loop over the vectors that sum to zero to try and find one
    #    That leads to a non trivial factorization
    #
    while ((k < len(solutions)) and (not nontrivial_first_factor) and (not nontrivial_second_factor)):
        
        (found_squares, u_plus_v, u_minus_v) = compute_difference_of_squares(R, f, m, N, integer_pairs, solutions[k]);

        if(found_squares):
                
            first_factor_of_N = gcd(u_plus_v, N);

            second_factor_of_N = gcd(u_minus_v, Integer(N/first_factor_of_N));

            if((1 != first_factor_of_N) and (N != first_factor_of_N)):
                nontrivial_first_factor = True;
            if((1 != second_factor_of_N) and (N != second_factor_of_N)):
                nontrivial_second_factor = True;
            else:
                print 'failed to find a nontrivial factor.'
    
        else:
            print 'failed to find squares.'

        k = k + 1;


    factors_of_N = [];

    N_remaining = N;

    if(nontrivial_first_factor):
        factors_of_N.append(first_factor_of_N);
        N_remaining = Integer(N_remaining/first_factor_of_N)

    if(nontrivial_second_factor):
        remaining_second_factor = gcd(second_factor_of_N, N_remaining);
        factors_of_N.append(remaining_second_factor);
        N_remaining = Integer(N_remaining/second_factor_of_N);

    if(1 != N_remaining):
        factors_of_N.append(N_remaining);

    factorization_of_N = factor(1);

    #
    #    Recursively factor non trivial factors
    #    Or N, if we completely failed to find a factor
    #

    for j in xrange(len(factors_of_N)):
        factorization_of_N = factorization_of_N * GNFS_Factoring_Algorithm(factors_of_N[j]);

    return factorization_of_N;

id=68| 
#
# Some test cases follow:
#

id=29| 
time GNFS_Factoring_Algorithm(45113);
///
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       54
upper bound on algebraic prime base (all first degree primes under):                                       108
upper bound on quadratic character prime base (all subsequent primes under):                                       227
search bounds for unit term: [ -181 , 181 ]
found 25 simultaneously smooth pairs.
found 33 simultaneously smooth pairs.
found 60 simultaneously smooth pairs.
found 65 simultaneously smooth pairs.
229 is a prime power, reducing to using ECM to factor
197 is a prime power, reducing to using ECM to factor
197 * 229
CPU time: 0.60 s,  Wall time: 1.24 s

id=69| 
197 * 229
///
45113

id=50| 
GNFS_Factoring_Algorithm(99999);
///
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       58
upper bound on algebraic prime base (all first degree primes under):                                       116
upper bound on quadratic character prime base (all subsequent primes under):                                       244
search bounds for unit term: [ -400 , 400 ]
found 34 simultaneously smooth pairs.
found 45 simultaneously smooth pairs.
found 82 simultaneously smooth pairs.
failed to find a nontrivial factor.
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       30
upper bound on algebraic prime base (all first degree primes under):                                       60
upper bound on quadratic character prime base (all subsequent primes under):                                       126
search bounds for unit term: [ -50 , 50 ]
found 18 simultaneously smooth pairs.
found 27 simultaneously smooth pairs.
found 45 simultaneously smooth pairs.
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       25
upper bound on algebraic prime base (all first degree primes under):                                       49
upper bound on quadratic character prime base (all subsequent primes under):                                       103
search bounds for unit term: [ -50 , 50 ]
found 17 simultaneously smooth pairs.
found 24 simultaneously smooth pairs.
found 44 simultaneously smooth pairs.
failed to find a nontrivial factor.
failed to find a nontrivial factor.
failed to find a nontrivial factor.
3 is a prime power, reducing to using ECM to factor
41 is a prime power, reducing to using ECM to factor
3 is a prime power, reducing to using ECM to factor
271 is a prime power, reducing to using ECM to factor
3^2 * 41 * 271

id=40| 
GNFS_Factoring_Algorithm(784352)
///
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       68
upper bound on algebraic prime base (all first degree primes under):                                       136
upper bound on quadratic character prime base (all subsequent primes under):                                       286
search bounds for unit term: [ -3138 , 3138 ]
found 29 simultaneously smooth pairs.
found 68 simultaneously smooth pairs.
found 83 simultaneously smooth pairs.
found 123 simultaneously smooth pairs.
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       44
upper bound on algebraic prime base (all first degree primes under):                                       88
upper bound on quadratic character prime base (all subsequent primes under):                                       185
search bounds for unit term: [ -50 , 50 ]
found 20 simultaneously smooth pairs.
found 47 simultaneously smooth pairs.
found 58 simultaneously smooth pairs.
16 is a prime power, reducing to using ECM to factor
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       30
upper bound on algebraic prime base (all first degree primes under):                                       60
upper bound on quadratic character prime base (all subsequent primes under):                                       126
search bounds for unit term: [ -50 , 50 ]
found 19 simultaneously smooth pairs.
found 45 simultaneously smooth pairs.
failed to find a nontrivial factor.
2 is a prime power, reducing to using ECM to factor
193 is a prime power, reducing to using ECM to factor
127 is a prime power, reducing to using ECM to factor
2^5 * 127 * 193

id=57| 

id=58| 
#
#    Factoring test code
#    This function will randomly generate an integer and factor it with the GNFS
#
#    It fails sometimes, due to stuff like low memory conditions, and bugs
#

def test_GNFS():

    N = Integer(ceil(uniform(1, 1000000)));

    print 'Factoring N =', N

    GNFS_Factorization = GNFS_Factoring_Algorithm(N);    
    
    if (GNFS_Factorization.value() != N):
        print 'Test Case N =', N, ' FAILED!!!'
        print 'says factorization is:', GNFS_Factorization
        return False

    print GNFS_Factorization

    print 'Test Case PASSED!'

    return True

id=56| 
test_GNFS()
///
Factoring N = 606063
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       67
upper bound on algebraic prime base (all first degree primes under):                                       134
upper bound on quadratic character prime base (all subsequent primes under):                                       282
search bounds for unit term: [ -2425 , 2425 ]
found 31 simultaneously smooth pairs.
found 43 simultaneously smooth pairs.
found 85 simultaneously smooth pairs.
3 is a prime power, reducing to using ECM to factor
202021 is a prime power, reducing to using ECM to factor
3 * 202021
Test Case PASSED!
True

id=60| 
#
#    At first, when reading up on the General Number Field Sieve algorithm
#    I had a hard time understanding the Quadratic character portion of
#    The algorithm, as such, I set this flag to test just what happens if
#    it is not used
#    observe that when the flag is set to false, the algorithm takes longer
#    and has many more errors about failing to find a square in the number field
#
use_quadratic_characters = False

GNFS_Factoring_Algorithm(784352)
///
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       68
upper bound on algebraic prime base (all first degree primes under):                                       136
upper bound on quadratic character prime base (all subsequent primes under):                                       136
search bounds for unit term: [ -3138 , 3138 ]
found 48 simultaneously smooth pairs.
found 101 simultaneously smooth pairs.
The supposedly found square is not square!
failed to find a square root in number field.
failed to find squares.
failed to find a nontrivial factor.
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       42
upper bound on algebraic prime base (all first degree primes under):                                       84
upper bound on quadratic character prime base (all subsequent primes under):                                       84
search bounds for unit term: [ -50 , 50 ]
found 40 simultaneously smooth pairs.
The supposedly found square is not square!
failed to find a square root in number field.
failed to find squares.
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       35
upper bound on algebraic prime base (all first degree primes under):                                       70
upper bound on quadratic character prime base (all subsequent primes under):                                       70
search bounds for unit term: [ -50 , 50 ]
found 24 simultaneously smooth pairs.
found 53 simultaneously smooth pairs.
Parameter Selection for GNFS:
upper bound on rational prime base (all primes under):                                       32
upper bound on algebraic prime base (all first degree primes under):                                       63
upper bound on quadratic character prime base (all subsequent primes under):                                       63
search bounds for unit term: [ -50 , 50 ]
found 23 simultaneously smooth pairs.
found 50 simultaneously smooth pairs.
The supposedly found square is not square!
failed to find a square root in number field.
failed to find squares.
4 is a prime power, reducing to using ECM to factor
127 is a prime power, reducing to using ECM to factor
2 is a prime power, reducing to using ECM to factor
4 is a prime power, reducing to using ECM to factor
193 is a prime power, reducing to using ECM to factor
2^5 * 127 * 193

id=63| 
use_quadratic_characters = True

id=64| 


2013-05-11 18:33