14.2 Univariate Polynomials

Module: sage.rings.polynomial_element

Author Log:

Module-level Functions

is_Polynomial( f)

Class: Polynomial

class Polynomial
Polynomial base class.
Polynomial( self, parent, [is_gen=False], [construct=False])

The following examples illustrate creation of elements of polynomial rings, and some basic arithmetic.

First we make a polynomial over the integers and do some arithmetic:

sage: x = PolynomialRing(IntegerRing()).gen()
sage: f = x^5 + 2*x^2 + (-1); f
x^5 + 2*x^2 - 1
sage: f^2
x^10 + 4*x^7 - 2*x^5 + 4*x^4 - 4*x^2 + 1

Next we do arithmetic in a sparse polynomial ring over the integers:

sage: R = PolynomialRing(IntegerRing(), "x"); x = R.gen(); R
Univariate Polynomial Ring in x over Integer Ring
sage: S = PolynomialRing(R, "Z"); Z = S.gen(); S
Univariate Polynomial Ring in Z over Univariate Polynomial Ring in x over
Integer Ring
sage: f = Z^3 + (x^2-2*x+1)*Z - 3; f
Z^3 + (x^2 - 2*x + 1)*Z + -3
sage: f*f
Z^6 + (2*x^2 - 4*x + 2)*Z^4 + (-6)*Z^3 + (x^4 - 4*x^3 + 6*x^2 - 4*x +
1)*Z^2 + (-6*x^2 + 12*x - 6)*Z + 9
sage: f^3 == f*f*f
True

To have the element print as 'y', give 'y' as the second argument to the PolynomialRing constructor.

sage: y = PolynomialRing(IntegerRing(), 'y').gen()
sage: y^3 - 2*y
y^3 - 2*y

Functions: base_extend,$  $ base_ring,$  $ coeffs,$  $ constant_coefficient,$  $ copy,$  $ degree,$  $ denominator,$  $ derivative,$  $ dict,$  $ factor,$  $ integral,$  $ inverse_of_unit,$  $ is_constant,$  $ is_gen,$  $ is_irreducible,$  $ is_monic,$  $ is_unit,$  $ is_zero,$  $ leading_coefficient,$  $ list,$  $ monic,$  $ name,$  $ newton_raphson,$  $ newton_slopes,$  $ polynomial,$  $ resultant,$  $ reverse,$  $ roots,$  $ truncate,$  $ valuation

base_extend( self, R)

Return a copy of this polynomial but with coefficients in R.

base_ring( self)

Return the base ring of the parent of self.

sage: x = PolynomialRing(ZZ).gen()
sage: x.base_ring()
Integer Ring
sage: (2*x+3).base_ring()
Integer Ring

coeffs( self)

Returns self.list().

(It potentially slightly faster better to use self.list() directly.)

sage: x = Q['x'].0
sage: f = 10*x^3 + 5*x + 2/17
sage: f.coeffs()
[2/17, 5, 0, 10]

copy( self)

Return a copy of self.

We create the polynomial $ f=x+3$ , then set $ g=f$ , and change the coefficient of $ x$ in $ g$ , which also changes the coefficient of $ x$ in $ f$ . If we instead copy $ f$ , then changing the coefficient of $ x$ of $ g$ does not change $ f$ .

sage: x = PolynomialRing(IntegerRing()).gen()
sage: f = x+3
sage: g = f
sage: g[1]=3
sage: f
3*x + 3
sage: g = f.copy()
sage: g[1]=5
sage: f
3*x + 3
sage: g
5*x + 3

degree( self)

Return the degree of this polynomial. The zero polynomial has degree -1.

sage: x = ZZ['x'].0
sage: f = x^93 + 2*x + 1
sage: f.degree()
93
sage: x = PolynomialRing(QQ, sparse=True).gen()
sage: f = x^100000 
sage: f.degree()
100000

sage: x = QQ['x'].0
sage: f = 2006*x^2006 - x^2 + 3
sage: f.degree()
2006
sage: f = 0*x
sage: f.degree()
-1
sage: f = x + 33
sage: f.degree()
1

Author: Naqi Jaffery (2006-01-24): examples

denominator( self)

Return the least common multiple of the denominators of the entries of self, when this makes sense, i.e., when the coefficients have a denominator function.

WARNING: This is not the denominator of the rational function defined by self, which would always be 1 since self is a polynomial.

First we compute the denominator of a polynomial with integer coefficients, which is of course 1.

sage: x = PolynomialRing(IntegerRing()).gen()
sage: f = x^3 + 17*x + 1
sage: f.denominator()
1

Next we compute the denominator of a polynomial with rational coefficients.

sage: Q = RationalField()
sage: x = PolynomialRing(Q).gen()
sage: f = Q('1/17')*x^19 - Q('2/3')*x + Q('1/3'); f
1/17*x^19 - 2/3*x + 1/3
sage: f.denominator()
51

Finally, we try to compute the denominator of a polynomial with coefficients in the real numbers, which is a ring whose elements do not have a denominator method.

sage: R = RealField()
sage: x = PolynomialRing(R).gen()
sage: f = x + R('0.3'); f
1.0000000000000000*x + 0.29999999999999999
sage: f.denominator()
Traceback (most recent call last):
...
AttributeError: 'mpfr.RealNumber' object has no attribute 'denominator'

factor( self)

Return the factorization of self.

INPUT:
    a polynomial

OUTPUT:
    Factorization -- the factorization of self, which is
    a product of a unit with a product of powers of irreducible
    factors.

Over a field the irreducible factors are all monic.

We factor some polynomials over $ \mathbf{Q}$ .

sage: x = QQ['x'].0
sage: f = (x^3 - 1)^2
sage: f.factor()
(x - 1)^2 * (x^2 + x + 1)^2

Notice that over the field $ \mathbf{Q}$ the irreducible factors are monic.

sage: f = 10*x^5 - 1
sage: f.factor()
(10) * (x^5 - 1/10)
sage: f = 10*x^5 - 10
sage: f.factor()
(10) * (x - 1) * (x^4 + x^3 + x^2 + x + 1)

Over $ \mathbf{Z}$ the irreducible factors need not be monic:

sage: x = ZZ['x'].0
sage: f = 10*x^5 - 1
sage: f.factor() 
(10*x^5 - 1)

We factor a non-monic polynomial over the finite field $ F_{25}$ .

sage: k, a = GF(25,'a').objgen()
sage: R, x = PolynomialRing(k).objgen()
sage: f = 2*x^10 + 2*x + 2*a
sage: F = f.factor(); F
(2) * (x + a + 2) * (x^2 + (a + 1)*x + a + 2) * (x^2 + 3*x + 4*a + 4) *
(x^5 + (3*a + 4)*x^4 + (3*a + 3)*x^3 + 2*a*x^2 + (3*a + 1)*x + 3*a + 1)

Notice that the unit factor is included when we multiply $ F$ back out.

sage: F.mul()
2*x^10 + 2*x + 2*a

Factorization also works even if the variable of the finite field is nefariously labeled "x".

sage: R, x = PolynomialRing(GF(3^2, 'x')).objgen()
sage: f = x^10 +7*x -13
sage: f.factor()
(x + 2*x + 1) * (x + x) * (x^4 + 2*x*x^3 + (x + 1)*x + 2) * (x^4 + (x +
2)*x^3 + (2*x + 2)*x + 2)
sage: f.parent().base_ring().assign_names(['a'])
sage: f.factor()
(x + 2*a + 1) * (x + a) * (x^4 + 2*a*x^3 + (a + 1)*x + 2) * (x^4 + (a +
2)*x^3 + (2*a + 2)*x + 2)

sage: k, a = GF(9,'x').objgen()
sage: x = PolynomialRing(k,'x0').gen()
sage: f = x^3 + x + 1
sage: f.factor()
(x0 + 2*x + 1) * (x0 + x) * (x0 + 2)

sage: f = 0*x
sage: f.factor()
Traceback (most recent call last):
...
ValueError: factorization of 0 not defined

sage: f = x^0
sage: f.factor()
1

is_monic( self)

Returns True if this polynomial is monic. The zero polynomial is by definition not monic.

sage: x = QQ['x'].0
sage: f = x + 33
sage: f.is_monic()
True
sage: f = 0*x
sage: f.is_monic()
False
sage: f = 3*x^3 + x^4 + x^2
sage: f.is_monic()
True
sage: f = 2*x^2 + x^3 + 56*x^5
sage: f.is_monic()
False

Author: Naqi Jaffery (2006-01-24): examples

monic( self)

Return this polynomial divided by its leading coefficient. Does not change this polynomial.

sage: x = QQ['x'].0
sage: f = 2*x^2 + x^3 + 56*x^5
sage: f.monic()
x^5 + 1/56*x^3 + 1/28*x^2
sage: f = (1/4)*x^2 + 3*x + 1
sage: f.monic()
x^2 + 12*x + 4

The following happens because $ f = 0$ cannot be made into a monic polynomial

sage: f = 0*x
sage: f.monic()
Traceback (most recent call last):
...
ZeroDivisionError: rational division by zero

Notice that the monic version of a polynomial over the integers is defined over the rationals.

sage: x = ZZ['x'].0
sage: f = 3*x^19 + x^2 - 37
sage: g = f.monic(); g
x^19 + 1/3*x^2 - 37/3
sage: g.parent()
Univariate Polynomial Ring in x over Rational Field

Author: Naqi Jaffery (2006-01-24): examples

newton_raphson( self, n, x0)

Return a list of n iterative approximations to a root of this polynomial, computed using the Newton-Raphson method.

The Newton-Raphson method is an iterative root-finding algorithm. For f(x) a polynomial, as is the case here, this is essentially the same as Horner's method.

INPUT:
   n -- an integer (=the number of iterations),
   x0 -- an initial guess x0.
   
OUTPUT:
   A list of numbers hopefully approximating a root of f(x)=0.

** If one of the iterates is a critical point of f then a ZeroDivisionError exception is raised.

sage: x = PolynomialRing(RealField(), 'x').gen()
sage: f = x^2 - 2
sage: f.newton_raphson(4, 1)
[1.5000000000000000, 1.4166666666666667, 1.4142156862745099,
1.4142135623746899]

Author: David Joyner and William Stein (2005-11-28)

newton_slopes( self, p)

Return the $ p$ -adic slopes of the Newton polygon of self, when this makes sense.

OUTPUT: - list of rational numbers

sage: x = QQ['x'].0
sage: f = x^3 + 2
sage: f.newton_slopes(2)
[1/3, 1/3, 1/3]

ALGORITHM: Uses PARI.

roots( self)

Return all roots of this polynomial.

sage: x = PolynomialRing(RationalField()).gen()
sage: f = x^3 - 1
sage: f.roots()
[(1, 1)]
sage: f = (x^3 - 1)^2
sage: f.roots()
[(1, 2)]

sage: f = -19*x + 884736
sage: f.roots()
[(884736/19, 1)]
sage: (f^20).roots()
[(884736/19, 20)]

truncate( self, n)

Replace this polynomial by $ \sum a_m x^m$ where the sum is over $ m < n$ . The resulting polynomial is equivalent to self modulo $ x^n$ .

valuation( self)

If $ f = a_r x^r + a_{r+1}x^{r+1} + \cdots$ , with $ a_r$ nonzero, then the valuation of $ f$ is $ r$ . The valuation of the zero polynomial is $ \infty$ .

Special Functions: __call__,$  $ __div__,$  $ __float__,$  $ __floordiv__,$  $ __getitem__,$  $ __int__,$  $ __invert__,$  $ __iter__,$  $ __long__,$  $ __mod__,$  $ __setitem__,$  $ _add_,$  $ _cmp_,$  $ _factor_pari_helper,$  $ _im_gens_,$  $ _integer_,$  $ _is_atomic,$  $ _latex_,$  $ _lcm,$  $ _mul_,$  $ _mul_fateman,$  $ _mul_generic,$  $ _mul_karatsuba,$  $ _pari_,$  $ _pari_init_,$  $ _pow,$  $ _repr,$  $ _repr_,$  $ _xgcd

__call__( self)

Evaluate polynomial at x=a using Horner's rule

INPUT:
    a -- ring element a; need not be in the coefficient
         ring of the polynomial.
         
OUTPUT:
    the value of f at a.

sage: x = Q['x'].gen()
sage: f = x/2 - 5
sage: f(3)
-7/2
sage: x = Z['x'].gen()
sage: f = (x-1)^5
sage: f(2/3)
-1/243

Author Log:

__div__( self, right)

sage: x = QQ['x'].gen()
sage: f = (x^3 + 5)/3; f
1/3*x^3 + 5/3
sage: f.parent()
Univariate Polynomial Ring in x over Rational Field

If we do the same over $ \mathbf{Z}$ the result has to lie in the fraction field.

sage: x  = ZZ['x'].gen()
sage: f = (x^3 + 5)/3; f
(x^3 + 5)/3
sage: f.parent()
Fraction Field of Univariate Polynomial Ring in x over Integer Ring

Note that / is a constructor for elements of the fraction field in all cases as long as both arguments have the same parent. This agrees with the behavior for integers and rational numbers.

sage: R, x = QQ['x'].objgen()
sage: f = x^3 + 5
sage: g = R(3)
sage: h = f/g; h
1/3*x^3 + 5/3
sage: h.parent()
Fraction Field of Univariate Polynomial Ring in x over Rational Field

This is another example over a non-prime finite field (submited by a student of Jon Hanke). It illustrates cancellation between the numerator and denominator over a non-prime finite field.

sage: R.<x> = PolynomialRing(GF(5^2), 'x')
sage: f = x^3 + 4*x
sage: f / (x - 1)
x^2 + x

__floordiv__( self, right)

Quotient of division of self by other. This is denoted //.

__mod__( self, other)

Remainder of division of self by other.

sage: x = PolynomialRing(IntegerRing()).gen()
sage: x % (x+1)
-1
sage: (x^3 + x - 1) % (x^2 - 1)
2*x - 1

_cmp_( self, other)

sage: x = QQ['x'].0
sage: 3*x^3  + 5 > 10*x^2 + 19
True
sage: f = x^2 - 2*x + 1; g= x^2 - 1
sage: f < g
True
sage: f > g
False
sage: g < f
False
sage: g > f
True

_im_gens_( self, codomain, im_gens)

sage: R, x = PolynomialRing(ZZ).objgen()
sage: H = Hom(R, QQ); H
Set of Homomorphisms from Univariate Polynomial Ring in x over Integer Ring
to Rational Field
sage: f = H([5]); f
Ring morphism:
  From: Univariate Polynomial Ring in x over Integer Ring
  To:   Rational Field
  Defn: x |--> 5
sage: f(x)
5
sage: f(x^2 + 3)
28

_lcm( self, other)

Let f and g be two polynomials. Then this function returns the monic least common multiple of f and g.

_mul_( self, right)

sage: x = PolynomialRing(IntegerRing()).gen()
sage: (x - 4)*(x^2 - 8*x + 16)
x^3 - 12*x^2 + 48*x - 64

_mul_fateman( self, right)

Returns the product of two polynomials using Kronecker's trick to do the multiplication. This could be used used over a generic base ring.

NOTES:

INPUT:
   self -- Polynomial
   right -- Polynomial (over same base ring as self)

OUTPUT: Polynomial
   The product self*right.

ALGORITHM: Based on a paper by R. Fateman

http://www.cs.berkeley.edu/ fateman/papers/polysbyGMP.pdf

The idea is to encode dense univariate polynomials as big integers, instead of sequences of coefficients. The paper argues that because integer multiplication is so cheap, that encoding 2 polynomials to big numbers and then decoding the result might be faster than popular multiplication algorithms. This seems true when the degree is larger than 200.

sage: S.<y> = PolynomialRing(RR)
sage: f = y^10 - 1.393493*y + 0.3
sage: f._mul_karatsuba(f)
1.0000000000000000*y^20 - 2.7869860000000002*y^11 +
0.60000000000000031*y^10 + 0.00000000000000011102230246251565*y^8 -
0.00000000000000011102230246251565*y^6 -
0.00000000000000011102230246251565*y^3 + 1.9418227410490003*y^2 -
0.83609580000000017*y + 0.089999999999999997
sage: f._mul_fateman(f)
1.0000000000000000*y^20 - 2.7869860000000002*y^11 +
0.59999999999999998*y^10 + 1.9418227410490003*y^2 - 0.83609580000000006*y +
0.089999999999999997

Advantages:

Drawbacks:

Author: Didier Deshommes (2006-05-25)

_mul_karatsuba( self, right)

Returns the product of two polynomials using the Karatsuba divide and conquer multiplication algorithm. This is only used over a generic base ring. (Special libraries like NTL are used, e.g., for the integers and rationals, which are much faster.)

INPUT:
   self: Polynomial
   right: Polynomial (over same base ring as self)

OUTPUT: Polynomial
   The product self*right.

ALGORITHM: The basic idea is to use that

$\displaystyle (aX + b) (cX + d) = acX^2 + ((a+b)(c+d)-ac-bd)X + bd
$

where ac=a*c and bd=b*d, which requires three multiplications instead of the naive four. (In my examples, strangely just doing the above with four multiplications does tend to speed things up noticeably.) Given f and g of arbitrary degree bigger than one, let e be min(deg(f),deg(g))/2. Write

$\displaystyle f = a X^e + b$    and $\displaystyle g = c X^e + d
$

and use the identity

$\displaystyle (aX^e + b) (cX^e + d) = ac X^{2e} +((a+b)(c+d) - ac - bd)X^e + bd
$

to recursively compute $ fg$ .

TIMINGS: On a Pentium M 1.8Ghz laptop: f=R.random(1000,bound=100) g=R.random(1000,bound=100) time h=f._mul_karatsuba(g) Time: 0.42 seconds The naive multiplication algorithm takes 14.58 seconds. In contrast, MAGMA does this sort of product almost instantly, and can easily deal with degree 5000. Basically MAGMA is 100 times faster at polynomial multiplication.

Over Z using NTL, multiplying two polynomials constructed using R.random(10000,bound=100) takes 0.10 seconds. Using MAGMA V2.11-10 the same takes 0.14 seconds. So in this case NTL is somewhat faster than MAGMA.

Over Q using PARI, multiplying two polynomials constructed using R.random(10000,bound=100) takes 1.23 seconds. Not good! TODO: use NTL polynomials over Z with a denominator instead of PARI.

NOTES: * Karatsuba multiplication of polynomials is also implemented in PARI in src/basemath/polarit3.c * The MAGMA documentation appears to give no information about how polynomial multiplication is implemented.

_pari_( self, [variable=None])

Return polynomial as a PARI object. Note that the variable will be "x" unless you explicitly specify otherwise, no matter what the polynomial indeterminate is.

sage: f = PolynomialRing(QQ)([0,1,2/3,3])
sage: pari(f)
3*x^3 + 2/3*x^2 + x

_xgcd( self, other)

Extended gcd of self and polynomial other.

Returns g, u, and v such that g = u*self + v*other.

sage: P, x = PolynomialRing(QQ).objgen()
sage: F = (x^2 + 2)*x^3; G = (x^2+2)*(x-3)
sage: g, u, v = F.xgcd(G)
sage: g, u, v
(27*x^2 + 54, 1, -x^2 - 3*x - 9)
sage: u*F + v*G
27*x^2 + 54
sage: x.xgcd(P(0))
(1, 0, x)
sage: f = P(0)
sage: f.xgcd(x)
(x, 0, 1)

Class: Polynomial_dense_mod_n

class Polynomial_dense_mod_n
A dense polynomial over the integers modulo n, where n is composite.

sage: R, x = PolynomialRing(Integers(16)).objgen()
sage: f = x^3 - x + 17
sage: loads(f.dumps()) == f
True
Polynomial_dense_mod_n( self, parent, [x=False], [check=False], [is_gen=True], [construct=None])

Functions: copy,$  $ degree,$  $ int_list,$  $ is_irreducible,$  $ list,$  $ ntl_set_directly,$  $ ntl_ZZ_pX,$  $ quo_rem

degree( self)

Return the degree of this polynomial. The zero polynomial has degree -1.

list( self)

sage: x = PolynomialRing(Integers(100)).gen()
sage: f = x^3 + 3*x - 17
sage: f.list()
[83, 3, 0, 1]

ntl_set_directly( self, v)

Set the value of this polynomial directly from a vector or string.

Polynomials over the integers modulo n are stored internally using NTL's ZZ_pX class. Use this function to set the value of this polynomial using the NTL constructor, which is potentially very fast. The input v is either a vector of ints or a string of the form '[ n1 n2 n3 ... ]' where the ni are integers and there are no commas between them. The optimal input format is the string format, since that's what NTL uses by default.

sage: R = PolynomialRing(Integers(100))
sage: R([1,-2,3])
3*x^2 + 98*x + 1
sage: f = R(0)
sage: f.ntl_set_directly([1,-2,3])
sage: f
3*x^2 + 98*x + 1
sage: f.ntl_set_directly('[1 -2 3 4]')
sage: f
4*x^3 + 3*x^2 + 98*x + 1

ntl_ZZ_pX( self)

Return underlying NTL representation of this polynomial. Additional ``bonus'' functionality is available through this function.

quo_rem( self, right)

Returns a tuple (quotient, remainder) where self = quotient*other + remainder.

Special Functions: __getitem__,$  $ __getslice__,$  $ __reduce__,$  $ __setitem__,$  $ _add_,$  $ _mul_,$  $ _pari_,$  $ _pow,$  $ _sub_

_mul_( self, right)

sage: x = PolynomialRing(Integers(100)).gen()
sage: (x - 2)*(x^2 - 8*x + 16)
x^3 + 90*x^2 + 32*x + 68

Class: Polynomial_dense_mod_p

class Polynomial_dense_mod_p
A dense polynomial over the integers modulo p, where p is prime.

Functions: discriminant,$  $ resultant

discriminant( self)

sage: x = PolynomialRing(GF(19)).gen()        
sage: f = x^3 + 3*x - 17
sage: f.discriminant()
12

resultant( self, other)

Returns the resultant of self and other, which must lie in the same polynomial ring.

INPUT:
    other -- a polynomial
OUTPUT:
    an element of the base ring of the polynomial ring

sage: x = PolynomialRing(GF(19)).gen()
sage: f = x^3 + x + 1;  g = x^3 - x - 1
sage: f.resultant(g)
11

Special Functions: __reduce__,$  $ _gcd,$  $ _xgcd

_gcd( self, right)

Return the GCD of self and other, as a monic polynomial.

_xgcd( self, right)

Return $ g, u, v$ such that g = u*self + v*right.

Class: Polynomial_generic_dense

class Polynomial_generic_dense
A generic dense polynomial.

sage: R, x = PolynomialRing(PolynomialRing(Q)).objgen()
sage: f = x^3 - x + 17
sage: type(f)
<class 'sage.rings.polynomial_element.Polynomial_generic_dense'>
sage: loads(f.dumps()) == f
True
Polynomial_generic_dense( self, parent, [x=False], [check=False], [is_gen=True], [construct=None])

Functions: degree,$  $ list

Special Functions: __getitem__,$  $ __getslice__,$  $ __setitem__,$  $ _Polynomial_generic_dense__normalize

Class: Polynomial_generic_dense_field

class Polynomial_generic_dense_field
Polynomial_generic_dense_field( self, parent, [x=False], [check=False], [is_gen=True], [construct=None])

Class: Polynomial_generic_field

class Polynomial_generic_field
Polynomial_generic_field( self, parent, [is_gen=False], [construct=False])

Functions: quo_rem

quo_rem( self, other)

Returns a tuple (quotient, remainder) where self = quotient*other + remainder.

sage: K = NumberField(x**ZZ(2)-ZZ(2),'t')
sage: P, x = PolynomialRing(K).objgen()
sage: x.quo_rem(K(1))
(x, 0)
sage: x.xgcd(K(1)) 
(1, 0, 1)

Special Functions: _gcd

_gcd( self, other)

Return the GCD of self and other, as a monic polynomial.

Class: Polynomial_generic_sparse

class Polynomial_generic_sparse
A generic sparse polynomial.

sage: R, x = PolynomialRing(PolynomialRing(Q), sparse=True).objgen()
sage: f = x^3 - x + 17
sage: type(f)
<class 'sage.rings.polynomial_element.Polynomial_generic_sparse'>
sage: loads(f.dumps()) == f
True
Polynomial_generic_sparse( self, parent, [x=False], [check=False], [is_gen=True], [construct=None])

Functions: degree,$  $ list

Special Functions: __getitem__,$  $ __getslice__,$  $ __setitem__,$  $ _Polynomial_generic_sparse__normalize

Class: Polynomial_generic_sparse_field

class Polynomial_generic_sparse_field

sage: R, x = PolynomialRing(RealField(), sparse=True).objgen()
sage: f = x^3 - x + 17
sage: type(f)
<class 'sage.rings.polynomial_element.Polynomial_generic_sparse_field'>
sage: loads(f.dumps()) == f
True
Polynomial_generic_sparse_field( self, parent, [x=False], [check=False], [is_gen=True], [construct=None])

Class: Polynomial_integer_dense

class Polynomial_integer_dense
A dense polynomial over the integers.
Polynomial_integer_dense( self, parent, [x=False], [check=False], [is_gen=True], [construct=None])

Functions: complex_roots,$  $ content,$  $ copy,$  $ degree,$  $ discriminant,$  $ factor_mod,$  $ factor_padic,$  $ gcd,$  $ lcm,$  $ list,$  $ ntl_set_directly,$  $ ntl_ZZX,$  $ quo_rem,$  $ resultant,$  $ xgcd

complex_roots( self, [flag=0])

Returns the complex roots of this polynomial.

INPUT:
    flag -- optional, and can be
            0: (default), uses Schonhage's method modified by Gourdon,
            1: uses a modified Newton method.
OUTPUT:
    list of complex roots of this polynomial, counted with multiplicities.

NOTE: Calls the pari function polroots.

We compute the roots of the characteristic polynomial of some Salem numbers:

sage: R = PolynomialRing(ZZ); x = R.gen()
sage: f = 1 - x^2 - x^3 - x^4 + x^6
sage: f.complex_roots()[0]    # todo: known bug in PARI 2.2.10 !!
0.71363917353690087

content( self)

Return the greatest common divisor of the coefficients of this polynomial.

degree( self)

Return the degree of this polynomial. The zero polynomial has degree -1.

discriminant( self)

sage: x = PolynomialRing(ZZ).gen()        
sage: f = x^3 + 3*x - 17
sage: f.discriminant()
-7911

factor_mod( self, p)

Return the factorization of self modulo the prime p.

INPUT:
    p -- prime

OUTPUT:
    factorization of self reduced modulo p.

sage: x = Z['x'].0
sage: f = -3*x*(x-2)*(x-9) + x
sage: f.factor_mod(3)
x
sage: f = -3*x*(x-2)*(x-9)
sage: f.factor_mod(3)
Traceback (most recent call last):
...
ValueError: factorization of 0 not defined

sage: f = 2*x*(x-2)*(x-9)
sage: f.factor_mod(7)
(2) * x * (x + 5)^2

factor_padic( self, p, [prec=10])

Return p-adic factorization of self to given precision.

INPUT:
    p -- prime
    prec -- integer; the precision

OUTPUT:
    factorization of self reduced modulo p.

gcd( self, right)

Return the GCD of self and other. The leading coefficient need not be 1.

lcm( self, right)

Return the LCM of self and other, as a monic polynomial.

list( self)

sage: x = PolynomialRing(ZZ).gen()
sage: f = x^3 + 3*x - 17
sage: f.list()
[-17, 3, 0, 1]

ntl_set_directly( self, v)

Set the value of this polynomial directly from a vector or string.

Polynomials over the integers are stored internally using NTL's ZZX class. Use this function to set the value of this polynomial using the NTL constructor, which is potentially quicker. The input v is either a vector of ints or a string of the form '[ n1 n2 n3 ... ]' where the ni are integers and there are no commas between them. The optimal input format is the string format, since that's what NTL uses.

sage: R = PolynomialRing(ZZ)
sage: R([1,2,3])
3*x^2 + 2*x + 1
sage: f = R(0)
sage: f.ntl_set_directly([1,2,3])
sage: f
3*x^2 + 2*x + 1
sage: f.ntl_set_directly('[1 2 3 4]')
sage: f
4*x^3 + 3*x^2 + 2*x + 1

ntl_ZZX( self)

Return underlying NTL representation of this polynomial. Additional ``bonus'' functionality may be available through this function.

quo_rem( self, right)

Returns a tuple (quotient, remainder) where self = quotient*other + remainder.

resultant( self, other)

Returns the resultant of self and other, which must lie in the same polynomial ring.

INPUT:
    other -- a polynomial
OUTPUT:
    an element of the base ring of the polynomial ring

NOTES: Implemented using NTL's polresultant function.

sage: x = PolynomialRing(ZZ).gen()
sage: f = x^3 + x + 1;  g = x^3 - x - 1
sage: f.resultant(g)
-8

xgcd( self, right)

Return $ g, u, v$ such that g = u*self + v*right.

If self and right are coprime as polynomials over the rationals, then $ g$ is guaranteed to be the resultant of self and right, as a constant polynomial.

sage: P, x = PolynomialRing(ZZ).objgen()
sage: F = (x^2 + 2)*x^3; G = (x^2+2)*(x-3)
sage: g, u, v = F.xgcd(G)
sage: g, u, v
(27*x^2 + 54, 1, -x^2 - 3*x - 9)
sage: u*F + v*G
27*x^2 + 54
sage: x.xgcd(P(0))
(1, 0, x)
sage: f = P(0)
sage: f.xgcd(x)
(x, 0, 1)
sage: F = (x-3)^3; G = (x-15)^2
sage: g, u, v = F.xgcd(G)
sage: g, u, v
(2985984, -432*x + 8208, 432*x^2 + 864*x + 14256)
sage: u*F + v*G
2985984

Special Functions: __getitem__,$  $ __getslice__,$  $ __reduce__,$  $ __setitem__,$  $ _add_,$  $ _mul_,$  $ _pari_,$  $ _pow,$  $ _sub_

_mul_( self, right)

sage: x = PolynomialRing(ZZ).gen()
sage: (x - 2)*(x^2 - 8*x + 16)
x^3 - 10*x^2 + 32*x - 32

Class: Polynomial_rational_dense

class Polynomial_rational_dense
A dense polynomial over the rational numbers.
Polynomial_rational_dense( self, parent, [x=False], [check=False], [is_gen=True], [construct=None])

Functions: complex_roots,$  $ copy,$  $ degree,$  $ disc,$  $ discriminant,$  $ factor_mod,$  $ factor_padic,$  $ galois_group,$  $ hensel_lift,$  $ is_irreducible,$  $ list,$  $ quo_rem,$  $ rescale,$  $ resultant

complex_roots( self, [flag=0])

Returns the complex roots of this polynomial.

INPUT:
    flag -- optional, and can be
            0: (default), uses Schonhage's method modified by Gourdon,
            1: uses a modified Newton method.
OUTPUT:
    list of complex roots of this polynomial, counted with multiplicities.

NOTE: Calls the pari function polroots.

We compute the roots of the characteristic polynomial of some Salem numbers:

sage: R = PolynomialRing(QQ); x = R.gen()
sage: f = 1 - x^2 - x^3 - x^4 + x^6
sage: f.complex_roots()[0]
0.71363917353690087

degree( self)

Return the degree of this polynomial. The zero polynomial has degree -1.

disc( self)

Same as discriminant().

discriminant( self)

sage: x = PolynomialRing(QQ).gen()        
sage: f = x^3 + 3*x - 17
sage: f.discriminant()
-7911

factor_mod( self, p)

Return the factorization of self modulo the prime p.

INPUT:
    p -- prime

OUTPUT:
    factorization of self reduced modulo p.

factor_padic( self, p, [prec=10])

Return p-adic factorization of self to given precision.

INPUT:
    p -- prime
    prec -- integer; the precision

OUTPUT:
    factorization of self reduced modulo p.

galois_group( self, [pari_group=False], [use_kash=False])

Return the Galois group of f as a permutation group.

INPUT:
    self -- an irreducible polynomial

pari_group - bool (default: False); if True instead return the Galois group as a PARI group. This has a useful label in it, and may be slightly faster since it doesn't require looking up a group in Gap. To get a permutation group from a PARI group P, type PermutationGroup(P).

use_kash - bool (default: False); if True use KASH's Galois command instead of using the PARI C library. An attempt is always made to use KASH if the degree of the polynomial is >= 12.

ALGORITHM: The Galois group is computed using PARI in C library mode, or possibly kash if available.

Note: The PARI documentation contains the following warning: The method used is that of resolvent polynomials and is sensitive to the current precision. The precision is updated internally but, in very rare cases, a wrong result may be returned if the initial precision was not sufficient.

sage: f = x^4 - 17*x^3 - 2*x + 1
sage: G = f.galois_group(); G            # uses optional database_gap package
Transitive group number 5 of degree 4
sage: G.gens()                           # uses optional database_gap package
((1,2,3,4), (1,2))
sage: G.order()                          # uses optional database_gap package
24

It is potentially useful to instead obtain the corresponding PARI group, which is little more than a $ 4$ -tuple. See the PARI manual for the exact details. (Note that the third entry in the tuple is in the new standard ordering.)

sage: f = x^4 - 17*x^3 - 2*x + 1
sage: G = f.galois_group(pari_group=True); G
PARI group [24, -1, 5, "S4"] of degree 4
sage: PermutationGroup(G)                # uses optional database_gap package
Transitive group number 5 of degree 4

You can use KASH to compute Galois groups as well. The avantage is that KASH can compute Galois groups of fields up to degree 23, whereas PARI only goes to degree 11. (In my not-so-thorough experiments PARI is faster than KASH.)

sage: f = x^4 - 17*x^3 - 2*x + 1
sage: f.galois_group(use_kash=true)      # requires optional KASH
Transitive group number 5 of degree 4

hensel_lift( self, p, e)

Assuming that self factors modulo $ p$ into distinct factors, computes the Hensel lifts of these factors modulo $ p^e$ . We assume that $ p$ has integer coefficients.

list( self)

sage: x = PolynomialRing(QQ).gen()
sage: f = x^3 + 3*x - QQ('17/13')
sage: f.list()
[-17/13, 3, 0, 1]

quo_rem( self, right)

Returns a tuple (quotient, remainder) where self = quotient*other + remainder.

rescale( self, a)

Return f(a*X).

resultant( self, other)

Returns the resultant of self and other, which must lie in the same polynomial ring.

INPUT:
    other -- a polynomial
OUTPUT:
    an element of the base ring of the polynomial ring

NOTES: Implemented using pari's polresultant function.

sage: x = PolynomialRing(RationalField()).gen()
sage: f = x^3 + x + 1;  g = x^3 - x - 1
sage: f.resultant(g)
-8

Special Functions: __getitem__,$  $ __getslice__,$  $ __reduce__,$  $ __setitem__,$  $ _add_,$  $ _mul_,$  $ _pow,$  $ _repr,$  $ _repr_,$  $ _sub_

_mul_( self, right)

sage: x = PolynomialRing(QQ).gen()
sage: (x - QQ('2/3'))*(x^2 - 8*x + 16)
x^3 - 26/3*x^2 + 64/3*x - 32/3

See About this document... for information on suggesting changes.