I am trying to implement a factorization algorithm for polynomials with integer coefficients, based on Kronecker's principle.
# Kronecker Factorization Algorithm
def kronecker_splitting(f):
"""
INPUT:
a polynomial f in ZZ[x]
OUTPUT:
a non trivial factor of f, or f if it is irreducible
"""
ZZX.<x> = ZZ[]
if f.parent() != ZZ[x]: raise ValueError("f must be in ZZ[x]")
if f.degree() == 1: return f
a = []
a.append(ZZ.random_element())
# Check if a0 is a root of f
if f.subs(a[0])==0: return x-a[0]
div = f.subs(a[0]).divisors() # compute positive divisors
M = [(d,) for d in div] # Turn m into a collection of tuples
for e in range(1,ceil(f.degree()/2)+1):
# Choose a not previously chosen integer
while(True):
aux = ZZ.random_element()
if aux not in a:
a.append(aux)
break
# Check that a[e] is not a root
if f.subs(a[e])==0: return x-a[e]
# compute divisors of a[e]
pos_div = f.subs(a[e]).divisors()
neg_div = [-d for d in pos_div]
div = pos_div + neg_div
Mp = M = [c_prev + (d,) for c_prev in M for d in div]
while len(Mp)>0:
c = Mp.pop()
# solve linear system
A = matrix.vandermonde(a).transpose()
b = A.solve_right(vector(c))
vector2poly = lambda h: sum([ci*x^i for i,ci in enumerate(h)])
g = vector2poly(b)
if b[e] != 0 and f % g == 0:
print("found non trivial factor {}".format(g))
print(g.parent())
return g.numerator()
return f # f is irreducible
def kronecker_factorization(f):
"""
INPUT:
f :- monic element of ZZ[x]
OUTPUT:
[f_1,...,f_s] :- monic irreducible factorization of f
"""
ZZX = ZZ[x]
factors = [f]; irreducible_factors = []
while len(factors) > 0:
g = factors.pop()
print("splitting {}".format(g))
h = kronecker_splitting(g)
# if g was already irreducible
if h == g: irreducible_factors.append(g)
else: factors.extend([h, ZZX(g/h)])
return irreducible_factors
ZZX.<x> = ZZ[x]
f = (x+1)*(x+2)*(x+3)
factors = kronecker_factorization(f)
print("{} factorizes as {}".format(f, factors))
assert(f == reduce(lambda a,b:a*b, factors))
for g in factors:
assert(g.parent() == ZZX)
assert(g.is_irreducible())
When I run the snippet above, I often get this output (try running it a few times to reproduce it, bcs of the randomness):
splitting x^3 + 6*x^2 + 11*x + 6
splitting x^2 + 3*x + 2
found non trivial factor 2*x + 2
Univariate Polynomial Ring in x over Rational Field
Error in lines 46-46
Traceback (most recent call last):
File "/cocalc/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 1188, in execute
flags=compile_flags) in namespace, locals
File "", line 1, in <module>
File "", line 15, in kronecker_factorization
File "sage/structure/parent.pyx", line 921, in sage.structure.parent.Parent.__call__ (build/cythonized/sage/structure/parent.c:9679)
return mor._call_(x)
File "sage/categories/map.pyx", line 790, in sage.categories.map.Map._call_ (build/cythonized/sage/categories/map.c:7325)
cpdef Element _call_(self, x):
File "/ext/sage/sage-8.4_1804/local/lib/python2.7/site-packages/sage/rings/fraction_field.py", line 1155, in _call_
raise TypeError("fraction must have unit denominator")
TypeError: fraction must have unit denominator
This confuses me, since if I run the following snippet:
QQX.<x> = QQ[x]
g = 3*x + 3
g.parent()
g.numerator()
I get the output: I expected, no errors in sight:
Univariate Polynomial Ring in x over Rational Field
3*x + 3
What is going on?
