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|