Dense univariate polynomials over \RR, implemented using MPFR

class sage.rings.polynomial.polynomial_real_mpfr_dense.PolynomialRealDense
__call__()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-2, 0, 1])
sage: f(10)
98.0000000000000
sage: f(CC.0)
-3.00000000000000
sage: f(2.0000000000000000000000000000000000000000000)
2.00000000000000
sage: f(RealField(10)(2))
2.0
sage: f(pi)
pi^2 - 2.00000000000000


sage: f = PolynomialRealDense(RR['x'], range(5))
sage: f(1)
10.0000000000000
sage: f(-1)
2.00000000000000
sage: f(0)
0.000000000000000
sage: f = PolynomialRealDense(RR['x'])
sage: f(12)
0.000000000000000
__getitem__()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], range(5)); f
4.00000000000000*x^4 + 3.00000000000000*x^3 + 2.00000000000000*x^2 + 1.00000000000000*x
sage: f[0]
0.000000000000000
sage: f[3]
3.00000000000000
sage: f[5]
0.000000000000000

Test slices:

sage: R.<x> = RealField(10)[]
sage: f = (x+1)^5; f
1.0*x^5 + 5.0*x^4 + 10.*x^3 + 10.*x^2 + 5.0*x + 1.0
sage: f[:3]
10.*x^2 + 5.0*x + 1.0
sage: f[3:]
1.0*x^5 + 5.0*x^4 + 10.*x^3
sage: f[1:4]
10.*x^3 + 10.*x^2 + 5.0*x
__init__()
x.__init__(...) initializes x; see x.__class__.__doc__ for signature
__neg__()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-2,0,1])
sage: -f
-1.00000000000000*x^2 + 2.00000000000000
static __new__()
T.__new__(S, ...) -> a new object with type S, a subtype of T
__reduce__()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-2, 0, 1])
sage: loads(dumps(f)) == f
True
_add_()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-2,0,1]); f
1.00000000000000*x^2 - 2.00000000000000
sage: g = PolynomialRealDense(RR['x'], range(5)); g
4.00000000000000*x^4 + 3.00000000000000*x^3 + 2.00000000000000*x^2 + 1.00000000000000*x
sage: f+g
4.00000000000000*x^4 + 3.00000000000000*x^3 + 3.00000000000000*x^2 + 1.00000000000000*x - 2.00000000000000
sage: g + f == f + g
True
sage: f + (-f)
0
_derivative()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [pi, 0, 2, 1]);
sage: f.derivative()
3.00000000000000*x^2 + 4.00000000000000*x
_lmul_()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-5,0,0,1]); f
1.00000000000000*x^3 - 5.00000000000000
sage: 4.0 * f
4.00000000000000*x^3 - 20.0000000000000
sage: f * -0.2
-0.200000000000000*x^3 + 1.00000000000000
_mul_()

Here we use the naive O(n^2) algorithm, as asymptotically faster algorithms such as Karatsuba can have very inaccurate results due to intermediate rounding errors.

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [1e20, 1])
sage: g = PolynomialRealDense(RR['x'], [1e30, 1])
sage: f*g
1.00000000000000*x^2 + 1.00000000010000e30*x + 1.00000000000000e50
sage: f._mul_karatsuba(g)
1.00000000000000*x^2 + 1.00000000000000e50
sage: f = PolynomialRealDense(RR['x'], range(5))
sage: g = PolynomialRealDense(RR['x'], range(3))
sage: f*g
8.00000000000000*x^6 + 10.0000000000000*x^5 + 7.00000000000000*x^4 + 4.00000000000000*x^3 + 1.00000000000000*x^2
_rmul_()
_sub_()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-3,0,1]); f
1.00000000000000*x^2 - 3.00000000000000
sage: g = PolynomialRealDense(RR['x'], range(4)); g
3.00000000000000*x^3 + 2.00000000000000*x^2 + 1.00000000000000*x
sage: f-g
-3.00000000000000*x^3 - 1.00000000000000*x^2 - 1.00000000000000*x - 3.00000000000000
sage: (f-g) == -(g-f)
True
change_ring()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-2, 0, 1.5])
sage: f.change_ring(QQ)
3/2*x^2 - 2
sage: f.change_ring(RealField(10))
1.5*x^2 - 2.0
sage: f.change_ring(RealField(100))
1.5000000000000000000000000000*x^2 - 2.0000000000000000000000000000
degree()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [1, 2, 3]); f
3.00000000000000*x^2 + 2.00000000000000*x + 1.00000000000000
sage: f.degree()
2
gcd()

Returns the gcd of self and other as a monic polynomial. Due to the inherit instability of division in this inexact ring, the results may not be entirely stable.

EXAMPLES:

sage: R.<x> = RR[]
sage: (x^3).gcd(x^5+1)
1.00000000000000
sage: (x^3).gcd(x^5+x^2)
1.00000000000000*x^2
sage: f = (x+3)^2 * (x-1)
sage: g = (x+3)^5
sage: f.gcd(g)
1.00000000000000*x^2 + 6.00000000000000*x + 9.00000000000000

Unless the division is exact (i.e. no rounding occurs) the returned gcd is almost certain to be 1.

sage: f = (x+RR.pi())^2 * (x-1)
sage: g = (x+RR.pi())^5
sage: f.gcd(g)
1.00000000000000
integral()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [3, pi, 1])
sage: f.integral()
0.333333333333333*x^3 + 1.57079632679490*x^2 + 3.00000000000000*x
list()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [1, 0, -2]); f
-2.00000000000000*x^2 + 1.00000000000000
sage: f.list()
[1.00000000000000, 0.000000000000000, -2.00000000000000]
quo_rem()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-2, 0, 1])
sage: g = PolynomialRealDense(RR['x'], [5, 1])
sage: q, r = f.quo_rem(g)
sage: q
1.00000000000000*x - 5.00000000000000
sage: r
23.0000000000000
sage: q*g + r == f
True
sage: fg = f*g
sage: fg.quo_rem(f)
(1.00000000000000*x + 5.00000000000000, 0)
sage: fg.quo_rem(g)
(1.00000000000000*x^2 - 2.00000000000000, 0)

sage: f = PolynomialRealDense(RR['x'], range(5))
sage: g = PolynomialRealDense(RR['x'], [pi,3000,4])
sage: q, r = f.quo_rem(g)
sage: g*q + r == f
True
reverse()

Returns x^d f(1/x) where d is the degree of f.

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [-3, pi, 0, 1])
sage: f.reverse()
-3.00000000000000*x^3 + 3.14159265358979*x^2 + 1.00000000000000
shift()

Returns this polynomial multiplied by the power x^n. If n is negative, terms below x^n will be discarded. Does not change this polynomial.

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RR['x'], [1, 2, 3]); f
3.00000000000000*x^2 + 2.00000000000000*x + 1.00000000000000
sage: f.shift(10)
3.00000000000000*x^12 + 2.00000000000000*x^11 + 1.00000000000000*x^10
sage: f.shift(-1)
3.00000000000000*x + 2.00000000000000
sage: f.shift(-10)
0

TESTS:

sage: f = RR['x'](0)
sage: f.shift(3).is_zero()
True
sage: f.shift(-3).is_zero()
True
truncate()

Returns the polynomial of degree < n which is equivalent to self modulo x^n.

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RealField(10)['x'], [1, 2, 4, 8])
sage: f.truncate(3)
4.0*x^2 + 2.0*x + 1.0
sage: f.truncate(100)
8.0*x^3 + 4.0*x^2 + 2.0*x + 1.0
sage: f.truncate(1)
1.0
sage: f.truncate(0)
0
truncate_abs()

Truncate all high order coefficients below bound.

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import PolynomialRealDense
sage: f = PolynomialRealDense(RealField(10)['x'], [10^-k for k in range(10)])
sage: f
1.0e-9*x^9 + 1.0e-8*x^8 + 1.0e-7*x^7 + 1.0e-6*x^6 + 0.000010*x^5 + 0.00010*x^4 + 0.0010*x^3 + 0.010*x^2 + 0.10*x + 1.0
sage: f.truncate_abs(0.5e-6)
1.0e-6*x^6 + 0.000010*x^5 + 0.00010*x^4 + 0.0010*x^3 + 0.010*x^2 + 0.10*x + 1.0
sage: f.truncate_abs(10.0)
0
sage: f.truncate_abs(1e-100) == f
True
sage.rings.polynomial.polynomial_real_mpfr_dense.make_PolynomialRealDense()

EXAMPLES:

sage: from sage.rings.polynomial.polynomial_real_mpfr_dense import make_PolynomialRealDense
sage: make_PolynomialRealDense(RR['x'], [1,2,3])
3.00000000000000*x^2 + 2.00000000000000*x + 1.00000000000000

Previous topic

Dense univariate polynomials over \ZZ/n\ZZ, implemented using NTL.

Next topic

Polynomial Interfaces to Singular

This Page