package for fast polynomial evaluation

[gmoroz]

The attached package provides conversion of univariate and multivariate polynomials into object that are optimized for fast evaluation on python object or low-levels c++ classes (see examples at the end).

It could enhanced the fast_callable function for several types, and also enhances in general the evaluation of polynomials on polynomials.

To test it, you can install it as a standard sage package with:

sage -i fast_polynomial-0.9.2.spkg​

Main features:

• handles univariate and multivariate polynomials
• specialized for several low-level types (mpfi, mpz, mpq, boost::interval)
• different evaluation layouts (horner, estrin, expanded, ...)
• easily extensible:
  • add new types (see fast_polynomial/interfaces/README)
  • add new layouts (see docstring of fast_polynomial.method)
• handles generic python/sage objects
• can be multi-threaded

Main limitations:

• only handles polynomial (no evaluation of trigonometric functions,...)
• polynomial needs to be converted to a fast callable object before evaluation (there is room for speed up on conversion time)

Examples and benchmarks:

from fast_polynomial import *
R.<x> = ZZ[x]
p = R.random_element(500,-100,100)

# evaluation of polynomials
q = python_polynomial(p, mode='horner')
r = python_polynomial(p, mode='estrin')
%timeit p(x+1) #5 loops, best of 3: 40.3 ms per loop
%timeit q(x+1) #5 loops, best of 3: 40.3 ms per loop
%timeit r(x+1) #125 loops, best of 3: 2.26 ms per loop
%timeit python_polynomial(p)(x+1) #125 loops, best of 3: 3.2 ms per loop

# evaluation of long integers
q = mpz_polynomial(p, num_threads=1)
r = mpz_polynomial(p, num_threads=2)
%timeit p(100) #625 loops, best of 3: 50.4 µs per loop
%timeit q(100) #625 loops, best of 3: 48.1 µs per loop
%timeit r(100) #625 loops, best of 3: 34.9 µs per loop


# evaluation of mpfi interval with precision 1000
q = mpfi_polynomial(p, 1000)
e = RealIntervalField(1000)(2^500, 2^500+1)
cmp(p(e),q(e)) #0
%timeit p(e)   #125 loops, best of 3: 2.71 ms per loop
%timeit q(e)   #625 loops, best of 3: 513 µs per loop
%timeit mpfi_polynomial(p)(e) #125 loops, best of 3: 1.15 ms per loop

# evaluation of boost interval (précision 53)
q = boost_polynomial(p, mode='horner')
r = boost_polynomial(p, mode='balanced', num_threads=2)
f = fast_callable(p, domain=float)
e = RIF(0.01)
%timeit p(e)    #125 loops, best of 3: 2.14 ms per loop
%timeit f(0.01) #625 loops, best of 3: 9.54 µs per loop
%timeit q(e)    #625 loops, best of 3: 13.4 µs per loop
%timeit r(e)    #625 loops, best of 3: 11.7 µs per loop
# Note that boost_polynomial evaluation offers more guarantees than raw float evaluation

# multivariate polynomials
R20 = PolynomialRing(QQ, 20,'x')
p = R20.random_element(5,100)
q = mpq_polynomial(p)
%timeit p((2/3,)*20) #125 loops, best of 3: 2.06 ms per loop
%timeit q((2/3,)*20) #625 loops, best of 3: 178 µs per loop
%timeit mpq_polynomial(p) #125 loops, best of 3: 1.91 ms per loop

[gmoroz]

fast_polynomial package compatible with sage >= 4.8

[gmoroz]

A minimal spkg (without boost dependency) to make the installation easier.

[gmoroz]

bug fix and add changelog.txt file

[vdelecroix]

Hello,

This would be a very nice addition in Sage!

First of all, I was not able to install your package on 6.2.rc0 with sage -i. On which version did you test it? It might come from the fact that the package structure has changed: did you have a look at the developer guide in Sage 6.2.rc0 and in particular the packaging section (it has been modified for the version 6.2.rc0 in the ticket #16048)?

The code you wrote looks like sage code but you wrote an external package. Was it on purpose that you did not write it directly inside Sage sources? It makes perfect sense to have an external package. But in that case, it should be relatively independent from Sage (I do not know if it is feasible, please tell me). There still might be some compilation options that depend on Sage (especially in the interfaces part).

Vincent

[gmoroz]

Hello,

Thanks for your interest. Indeed I only tested it for sage 5.9. I will look into the new package structure and update it.

About the package, I did write it as an external package such that it can easily be used directly within python only. All the code related to sage should in theory be in the interface directory only. The idea is that in order to use fast_polynomial with mpmath for example, it is only required to add a corresponding interface file in the interfaces directory (telling how to convert polynomials from mpmath to fast_polynomial) and tell in the setup.py file which interface to use. I must emphasize that this is in theory only, since I only wrote interface files for sage.

The other reason I wrote it as an external package is that some part also depends on the boost::interval library, that was not in sage at the time I started the project.

Guillaume

[vdelecroix]

Salut Guillaume,

It seems that your package is less independent of Sage than what you said: you import some components of the Sage library in fast_polynomial/generic/evaluation/graph.pyx and fast_polynomial/generic/polynomial.pyx.

To my mind, it would be better (for your work and for Sage) to distribute your library independently of Sage. It can be on your webpage, github or whatever. That would be a Python library with its own testing module. Once the library is ready and run within pure Python, it will be trivial to build a Sage spkg. I think that your library might interest some other projects such as ​Anaconda, ​GMPY and ​the Ipython notebook.

All best, Vincent

[jdemeyer]

Note: I personally don't care at all about this, but you should make a new-style package, see ​http://www.sagemath.org/doc/developer/#packaging-third-party-code