Search code examples
pythonpython-3.xsympypolynomial-mathgmpy

Sympy polynomials with `mpfr` coefficients?


I want to use Sympy's polynomials, but I also want to use higher-precision coefficients.

Just Doing It seems to give me polynomials with sympy.core.numbers.float coefficients.

import sympy
from sympy import Poly
from sympy.abc import x
from gmpy2 import mpfr, get_context

get_context().precision = 150

#float64 can't tell this from 1.0
one_and_change = mpfr('1.0000000000000000000000000000000000001')
#mpfr('1.0000000000000000000000000000000000001000000005',150)

p = [one_and_change]
px = Poly(p, x)

print(px)
# Poly(1.0, x, domain='RR')
print(px.is_one)
# True
print(type(px.all_coeffs()[0]))
# <class 'sympy.core.numbers.Float'>

I've also tried sympy.mpmath.mpf, with the same results.

This also didn't work:[1]

domain = sympy.polys.domains.realfield.RealField(150)
px = Poly(p, x, domain=domain)
print(type(px.all_coeffs()[0]))
# <class 'sympy.core.numbers.Float'>

Solution

  • There are a few obstacles:

    • gmpy.mpfr has no ._sympy_ method, so it will be converted to float in an intermediate step.
    • sympy.Poly, by default, uses sympy.polys.domain.RR for floating point coefficients, and will use it to convert.
    • RR, when it loads, uses 53 as its precision, and ignores mpmath.mp.precision.
    • When converting, RR uses RR.precision and ignores the argument's precision.

    Solution:

    1. Coefficients must be a Sympy type (e.g. sympy.Float, which has extended precision).
    2. For the domain, one of the following:
      • Set sympy.polys.domain.RR._context.prec = my_precision.
      • Pass domain='EX', which doesn't do conversion.
      • Pass in a custom domain (such as sympy.polys.domains.realfield.RealField(150)).