poly
poly copied to clipboard
Fast polynomial arithmetic in Haskell (dense and sparse, univariate and multivariate, usual and Laurent)
poly

Haskell library for univariate and multivariate polynomials, backed by Vectors.
> -- Univariate polynomials
> (X + 1) + (X - 1) :: VPoly Integer
2 * X
> (X + 1) * (X - 1) :: UPoly Int
1 * X^2 + (-1)
> -- Multivariate polynomials
> (X + Y) * (X - Y) :: VMultiPoly 2 Integer
1 * X^2 + (-1) * Y^2
> (X + Y + Z) ^ 2 :: UMultiPoly 3 Int
1 * X^2 + 2 * X * Y + 2 * X * Z + 1 * Y^2 + 2 * Y * Z + 1 * Z^2
> -- Laurent polynomials
> (X^-2 + 1) * (X - X^-1) :: VLaurent Integer
1 * X + (-1) * X^-3
> (X^-1 + Y) * (X + Y^-1) :: UMultiLaurent 2 Int
1 * X * Y + 2 + 1 * X^-1 * Y^-1
Vectors
Poly v a is polymorphic over a container v, implementing the Vector interface, and coefficients of type a. Usually v is either a boxed vector from Data.Vector or an unboxed vector from Data.Vector.Unboxed. Use unboxed vectors whenever possible, e. g., when the coefficients are Ints or Doubles.
There are handy type synonyms:
type VPoly a = Poly Data.Vector.Vector a
type UPoly a = Poly Data.Vector.Unboxed.Vector a
Construction
The simplest way to construct a polynomial is using the pattern X:
> X^2 - 3 * X + 2 :: UPoly Int
1 * X^2 + (-3) * X + 2
(Unfortunately, types are often ambiguous and must be given explicitly.)
While being convenient to read and write in REPL, X is relatively slow. The fastest approach is to use toPoly, providing it with a vector of coefficients (constant term first):
> toPoly (Data.Vector.Unboxed.fromList [2, -3, 1 :: Int])
1 * X^2 + (-3) * X + 2
Alternatively one can enable {-# LANGUAGE OverloadedLists #-} and simply write
> [2, -3, 1] :: UPoly Int
1 * X^2 + (-3) * X + 2
There is a shortcut to construct a monomial:
> monomial 2 3.5 :: UPoly Double
3.5 * X^2 + 0.0 * X + 0.0
Operations
Most operations are provided by means of instances, like Eq and Num. For example,
> (X^2 + 1) * (X^2 - 1) :: UPoly Int
1 * X^4 + 0 * X^3 + 0 * X^2 + 0 * X + (-1)
One can also find it convenient to scale by a monomial (cf. monomial above):
> scale 2 3.5 (X^2 + 1) :: UPoly Double
3.5 * X^4 + 0.0 * X^3 + 3.5 * X^2 + 0.0 * X + 0.0
While Poly cannot be made an instance of Integral (because there is no meaningful toInteger),
it is an instance of GcdDomain and Euclidean from the semirings package. These type classes
cover the main functionality of Integral, providing division with remainder and gcd / lcm:
> Data.Euclidean.gcd (X^2 + 7 * X + 6) (X^2 - 5 * X - 6) :: UPoly Int
1 * X + 1
> Data.Euclidean.quotRem (X^3 + 2) (X^2 - 1 :: UPoly Double)
(1.0 * X + 0.0,1.0 * X + 2.0)
Miscellaneous utilities include eval for evaluation at a given point,
and deriv / integral for taking the derivative and an indefinite integral, respectively:
> eval (X^2 + 1 :: UPoly Int) 3
10
> deriv (X^3 + 3 * X) :: UPoly Double
3.0 * X^2 + 0.0 * X + 3.0
> integral (3 * X^2 + 3) :: UPoly Double
1.0 * X^3 + 0.0 * X^2 + 3.0 * X + 0.0
Deconstruction
Use unPoly to deconstruct a polynomial to a vector of coefficients (constant term first):
> unPoly (X^2 - 3 * X + 2 :: UPoly Int)
[2,-3,1]
Further, leading is a shortcut to obtain the leading term of a non-zero polynomial,
expressed as a power and a coefficient:
> leading (X^2 - 3 * X + 2 :: UPoly Double)
Just (2,1.0)
Flavours
-
Data.Polyprovides dense univariate polynomials with aNum-based interface. This is a default choice for most users. -
Data.Poly.Semiringprovides dense univariate polynomials with aSemiring-based interface. -
Data.Poly.Laurentprovides dense univariate Laurent polynomials with aSemiring-based interface. -
Data.Poly.Sparseprovides sparse univariate polynomials with aNum-based interface. Besides that, you may find it easier to use in the REPL because of a more readableShowinstance, skipping zero coefficients. -
Data.Poly.Sparse.Semiringprovides sparse univariate polynomials with aSemiring-based interface. -
Data.Poly.Sparse.Laurentprovides sparse univariate Laurent polynomials with aSemiring-based interface. -
Data.Poly.Multiprovides sparse multivariate polynomials with aNum-based interface. -
Data.Poly.Multi.Semiringprovides sparse multivariate polynomials with aSemiring-based interface. -
Data.Poly.Multi.Laurentprovides sparse multivariate Laurent polynomials with aSemiring-based interface.
All flavours are available backed by boxed or unboxed vectors.
Performance
As a rough guide, poly is at least 20x-40x faster than the polynomial library.
Multiplication is implemented via the Karatsuba algorithm.
Here are a couple of benchmarks for UPoly Int:
| Benchmark | polynomial, μs | poly, μs | speedup |
|---|---|---|---|
| addition, 100 coeffs. | 45 | 2 | 22x |
| addition, 1000 coeffs. | 441 | 17 | 25x |
| addition, 10000 coeffs. | 6545 | 167 | 39x |
| multiplication, 100 coeffs. | 1733 | 33 | 52x |
| multiplication, 1000 coeffs. | 442000 | 1456 | 303x |