Simple multivariate polynomials in Haskell. This package deals with multivariate polynomials over a commutative ring, fractions of multivariate polynomials over a commutative field, and multivariate polynomials with symbolic parameters in their coefficients.
The main type provided by this package is Spray a
.
An object of type Spray a
represents a multivariate polynomial whose
coefficients are represented by the objects of type a
. For example:
import Math.Algebra.Hspray
x = lone 1 :: Spray Double
y = lone 2 :: Spray Double
z = lone 3 :: Spray Double
poly = (2 *^ (x^**^3 ^*^ y ^*^ z) ^+^ x^**^2) ^*^ (4 *^ (x ^*^ y ^*^ z))
putStrLn $ prettyNumSpray poly
-- 8.0*x^4.y^2.z^2 + 4.0*x^3.y.z
This is the easiest way to construct a spray: first introduce the polynomial
variables with the lone
function, and then combine them with arithmetic
operations.
There are numerous functions to print a spray. If you don't like the letters
x
, y
, z
in the output of prettyNumSpray
, you can use prettyNumSprayXYZ
to change them to whatever you want:
putStrLn $ prettyNumSprayXYZ ["A","B","C"] poly
-- 8.0*A^4.B^2.C^2 + 4.0*A^3.B.C
Note that this function does not throw an error if you don't provide enough
letters; in such a situation, it takes the first given letter and it appends
it with the digit i
to denote the i
-th variable:
putStrLn $ prettyNumSprayXYZ ["A","B"] poly
-- 8.0*A1^4.A2^2.A3^2 + 4.0*A1^3.A2.A3
This is the same output as the one of prettyNumSprayX1X2X3 "A" poly
.
More generally, one can use the type Spray a
as long as the type a
has
the instances Eq
and Algebra.Ring
(defined in the numeric-prelude
library). For example a = Rational
:
import Math.Algebra.Hspray
import Data.Ratio
x = lone 1 :: QSpray -- QSpray = Spray Rational
y = lone 2 :: QSpray
z = lone 3 :: QSpray
poly = ((2%3) *^ (x^**^3 ^*^ y ^*^ z) ^-^ x^**^2) ^*^ ((7%4) *^ (x ^*^ y ^*^ z))
putStrLn $ prettyQSpray poly
-- (7/6)*x^4.y^2.z^2 - (7/4)*x^3.y.z
Or a = Spray Double
:
import Math.Algebra.Hspray
alpha = lone 1 :: Spray Double
x = lone 1 :: Spray (Spray Double)
y = lone 2 :: Spray (Spray Double)
poly = ((alpha *^ x) ^+^ (alpha *^ y))^**^2
showSprayXYZ' (prettyNumSprayXYZ ["alpha"]) ["x","y"] poly
-- (alpha^2)*x^2 + (2.0*alpha^2)*x.y + (alpha^2)*y^2
We will come back to these sprays of type Spray (Spray a)
. They can
be used to represent parametric polynomials.
import Math.Algebra.Hspray
x = lone 1 :: Spray Double
y = lone 2 :: Spray Double
z = lone 3 :: Spray Double
spray = 2 *^ (x ^*^ y ^*^ z)
-- evaluate spray at x=2, y=1, z=2
evalSpray spray [2, 1, 2]
-- 8.0
import Math.Algebra.Hspray
import Data.Ratio
x1 = lone 1 :: Spray Rational
x2 = lone 2 :: Spray Rational
x3 = lone 3 :: Spray Rational
spray = x1^**^2 ^+^ x2 ^+^ x3 ^-^ unitSpray
putStrLn $ prettyQSprayX1X2X3 "x" spray
-- x1^2 + x2 + x3 - 1
--
-- substitute x1 -> 2 and x3 -> 3, and don't substitute x2
spray' = substituteSpray [Just 2, Nothing, Just 3] spray
putStrLn $ prettyQSprayX1X2X3 "x" spray'
-- x2 + 6
import Math.Algebra.Hspray
x = lone 1 :: Spray Double
y = lone 2 :: Spray Double
z = lone 3 :: Spray Double
spray = 2 *^ (x ^*^ y ^*^ z) ^+^ (3 *^ x^**^2)
putStrLn $ prettyNumSpray spray
-- 3.0*x^2 + 2.0*x.y.z
--
-- derivative with respect to x
putStrLn $ prettyNumSpray $ derivative 1 spray
-- 6.0*x + 2.0*y.z"
As of version 2.0.0, it is possible to compute a Gröbner basis of the ideal generated by a list of spray polynomials.
One of the numerous applications of Gröbner bases is the implicitization of a system of parametric equations. That means that they allow to get an implicit equation equivalent to a given set of parametric equations, for example an implicit equation of a parametrically defined surface. Let us give an example of implicitization for an Enneper surface. Here are the parametric equations of this surface:
import Math.Algebra.Hspray
import Data.Ratio ( (%) )
u = qlone 1
v = qlone 2
x = 3*^u ^+^ 3*^(u ^*^ v^**^2) ^-^ u^**^3
y = 3*^v ^+^ 3*^(u^**^2 ^*^ v) ^-^ v^**^3
z = 3*^u^**^2 ^-^ 3*^v^**^2
The first step of the implicitization is to compute a Gröbner basis of the following list of polynomials:
generators = [x ^-^ qlone 3, y ^-^ qlone 4, z ^-^ qlone 5]
Once the Gröbner basis is obtained, one has to identify which of its elements
do not involve the first two variables (u
and v
):
gbasis = groebnerBasis generators True
isFreeOfUV :: QSpray -> Bool
isFreeOfUV p = not (involvesVariable p 1) && not (involvesVariable p 2)
freeOfUV = filter isFreeOfUV gbasis
Then the implicit equations (there can be several ones) are obtained by removing the first two variables from these elements:
results = map (dropVariables 2) freeOfUV
Here we find only one implicit equation (i.e. length results
is 1
).
We only partially display it because it is long:
implicitEquation = results !! 0
putStrLn $ prettyQSpray implicitEquation
-- x^6 - 3*x^4.y^2 + (5/9)*x^4.z^3 + 6*x^4.z^2 - 3*x^4.z + 3*x^2.y^4 + ...
Let us check it is correct. We take a point on the Enneper surface, for
example the point corresponding to
xyz = map (evaluateAt [1%4, 2%3]) [x, y, z]
If the implicitization is correct, then the polynomial implicitEquation
should be zero at any point on the Enneper surface. This is true for our point:
evaluateAt xyz implicitEquation == 0
-- True
In the unit tests,
you will find a more complex example of implicitization: the implicitization
of the ellipse parametrically defined by
To construct a spray using the ordinary symbols +
, -
, *
and ^
,
one can hide these operators from Prelude and import them from the
numeric-prelude library; constructing a spray in this context is easier:
import Prelude hiding ((+), (-), (*), (^), (*>), (<*))
import qualified Prelude as P
import Algebra.Additive
import Algebra.Module
import Algebra.Ring
import Math.Algebra.Hspray
import Data.Ratio
x = lone 1 :: QSpray
y = lone 2 :: QSpray
z = lone 3 :: QSpray
spray = ((2%3) *^ (x^**^3 ^*^ y ^*^ z) ^-^ x^**^2) ^*^ ((7%4) *^ (x ^*^ y ^*^ z))
spray' = ((2%3) *^ (x^3 * y * z) - x^2) * ((7%4) *^ (x * y * z))
spray == spray'
-- True
Note that *>
could be used instead of *^
but running lambda *> spray
possibly throws an "ambiguous type" error regarding the type of lambda
.
Maybe better (I didn't try yet), follow the "Usage" section on the Hackage page of numeric-prelude.
Assume you have the polynomial a * (x² + y²) + 2b/3 * z
,
where a
and b
are symbolic rational numbers. You can represent this
polynomial by a Spray (Spray Rational)
spray as follows:
import Prelude hiding ((*), (+), (-), (^))
import qualified Prelude as P
import Algebra.Additive
import Algebra.Ring
import Math.Algebra.Hspray
x = lone 1 :: Spray (Spray Rational)
y = lone 2 :: Spray (Spray Rational)
z = lone 3 :: Spray (Spray Rational)
a = lone 1 :: Spray Rational
b = lone 2 :: Spray Rational
spray = a *^ (x^2 + y^2) + ((2 *^ b) /^ 3) *^ z
putStrLn $
showSprayXYZ' (prettyQSprayXYZ ["a","b"]) ["X","Y","Z"] spray
-- (a)*X^2 + (a)*Y^2 + ((2/3)*b)*Z
You can extract the powers and the coefficients as follows:
l = toList spray
map fst l
-- [[0,0,1],[2],[0,2]]
map toList $ map snd l
-- [[([0,1],2 % 3)],[([1],1 % 1)],[([1],1 % 1)]]
These Spray (Spray a)
sprays can be very useful. They represent polynomials
whose coefficients polynomially depend on some parameters.
Actually there is a type alias of Spray (Spray a)
in hspray, namely
SimpleParametricSpray a
, and there are some convenient functions to deal
with sprays of this type. There is also a type alias of
SimpleParametricSpray Rational
, namely SimpleParametricQSpray
.
For example we can print our SimpleParametricQSpray
spray spray
as follows:
putStrLn $
prettySimpleParametricQSprayABCXYZ ["a","b"] ["X","Y","Z"] spray
-- { a }*X^2 + { a }*Y^2 + { (2/3)*b }*Z
The
Gegenbauer polynomials
are a real-life example of polynomials that can be represented by
SimpleParametricQSpray
sprays. They are univariate polynomials whose
coefficients polynomially depend on a parameter
gegenbauerPolynomial :: Int -> SimpleParametricQSpray
gegenbauerPolynomial n
| n == 0 = unitSpray
| n == 1 = (2.^a) *^ x
| otherwise =
(2.^(n'' ^+^ a) /^ n') *^ (x ^*^ gegenbauerPolynomial (n - 1))
^-^ ((n'' ^+^ 2.^a ^-^ unitSpray) /^ n') *^ gegenbauerPolynomial (n - 2)
where
x = lone 1 :: SimpleParametricQSpray
a = lone 1 :: QSpray
n' = toRational n
n'' = constantSpray (n' - 1)
Let's try it:
n = 3
g = gegenbauerPolynomial n
putStrLn $
prettySimpleParametricQSprayABCXYZ ["alpha"] ["X"] g
-- { (4/3)*alpha^3 + 4*alpha^2 + (8/3)*alpha }*X^3 + { -2*alpha^2 - 2*alpha }*X
Let's check the differential equation given in the Wikipedia article:
g' = derivative 1 g
g'' = derivative 1 g'
alpha = lone 1 :: QSpray
x = lone 1 :: SimpleParametricQSpray
nAsSpray = constantSpray (toRational n)
shouldBeZero =
(unitSpray ^-^ x^**^2) ^*^ g''
^-^ (2.^alpha ^+^ unitSpray) *^ (x ^*^ g')
^+^ n.^(nAsSpray ^+^ 2.^alpha) *^ g
putStrLn $ prettySpray shouldBeZero
-- 0
Now, how to substitute a value to the parameter substituteParameters
to
perform this task:
import Data.Ratio (%)
putStrLn $
prettyQSpray'' $ substituteParameters g [1%2]
-- (5/2)*X^3 - (3/2)*X
This is a Spray Rational
spray.
The Wikipedia article also provides the value at evalParametricSpray
:
putStrLn $
prettyQSprayXYZ ["alpha"] $ evalParametricSpray g [1]
-- (4/3)*alpha^3 + 2*alpha^2 + (2/3)*alpha
This is also a Spray Rational
spray.
Since you have just seen that the type Spray (Spray a)
is named
SimpleParametricSpray a
, you probably guessed there is also a more general
type named ParametricSpray a
. Yes, and this is an alias of
Spray (RatioOfSprays a)
, where the type RatioOfSprays a
has not been
discussed yet. The objects of this type represent fractions of multivariate
polynomials and so this type is a considerable enlargment of the Spray a
type. Thus the Spray (RatioOfSprays a)
sprays can represent multivariate
polynomials whose coefficients depend on some parameters, with a dependence
described by a fraction of polynomials in these parameters. Let's start with
a short presentation of the ratios of sprays.
The type RatioOfSprays a
, whose objects represent ratios of sprays, has
been introduced in version 0.2.7.0.
To construct a ratio of sprays, apply %//%
between its numerator and
its denominator:
import Math.Algebra.Hspray
x = qlone 1 -- shortcut for lone 1 :: Spray Rational
y = qlone 2
rOS = (x ^-^ y) %//% (x^**^2 ^-^ y^**^2)
putStrLn $ prettyRatioOfQSprays rOS
-- [ 1 ] %//% [ x + y ]
The %//%
operator always returns an irreducible fraction. If you are
sure that your numerator and your denominator are coprime, you can use
the %:%
instead, to gain some efficiency. But if they are not coprime, this
can have unfortunate consequences.
The RatioOfSprays a
type makes sense when a
has a field instance, and then
it has a field instance too. To use the field operations, import the necessary
modules from numeric-prelude, and hide these operations from the Prelude
module; then you can also use the numeric-prelude operations for sprays,
instead of using ^+^
, ^-^
, ^*^
, ^**^
:
import Prelude hiding ((+), (-), (*), (/), (^), (*>), (<*))
import qualified Prelude as P
import Algebra.Additive
import Algebra.Module
import Algebra.RightModule
import Algebra.Ring
import Algebra.Field
import Math.Algebra.Hspray
x = qlone 1
y = qlone 2
p = x^2 - 3*^(x * y) + y^3
q = x - y
rOS1 = p^2 %//% q
rOS2 = rOS1 + unitRatioOfSprays
rOS = rOS1^2 + rOS1*rOS2 - rOS1/rOS2 + rOS2 -- slow!
(rOS1 + rOS2) * (rOS1 - rOS2) == rOS1^2 - rOS2^2
-- True
rOS / rOS == unitRatioOfSprays
-- True
Actually, as of version 0.5.0.0, it is possible to apply the operations
^+^
, ^-^
, ^*^
, ^**^
and *^
to the ratios of sprays.
The RatioOfSprays a
type also has left and right module instances over a
and over Spray a
as well. That means you can multiply a ratio of sprays by
a scalar and by a spray, by using, depending on the side, either *>
or <*
:
import Data.Ratio ( (%) )
rOS' = (3%4::Rational) *> rOS^2 + p *> rOS
You can also divide a ratio of sprays by a spray with %/%
:
p *> (rOS' %/% p) == rOS'
-- True
rOS1 %/% p == p %//% q
-- True
When a
has a field instance, both a Spray a
spray and a RatioOfSprays a
ratio of sprays can be divided by a scalar with the />
operator:
k = 3 :: Rational
(p /> k) *> rOS == p *> (rOS /> k)
-- True
Use evalRatioOfSprays
to evaluate a ratio of sprays:
import Data.Ratio ( (%) )
f :: Algebra.Field.C a => a -> a -> a
f u v = u^2 + u*v - u/v + v
rOS == f rOS1 rOS2
-- True
values = [2%3, 7%4]
r1 = evalRatioOfSprays rOS1 values
r2 = evalRatioOfSprays rOS2 values
evalRatioOfSprays rOS values == f r1 r2
-- True
Recall that SimpleParametricSpray a = Spray (Spray a)
and
ParametricSpray a = Spray (RatioOfSprays a)
, and we have the aliases
SimpleParametricQSpray = SimpleParametricSpray Rational
and
ParametricQSpray = ParametricSpray Rational
.
The functions substituteParameters
and evalParametricSpray
, that we
previously applied to a SimpleParametricSpray a
spray, are also applicable
to a ParametricSpray a
spray. We didn't mention the function
changeParameters
yet, which is also applicable to these two types of
parametric sprays. This function performs some polynomial transformations of
the parameters of a parametric spray. For example, consider the
Jacobi polynomials.
They are univariate polynomials with two parameters ParametricQSpray
sprays. In fact
it seems that the coefficients of the Jacobi polynomials polynomially
depend on SimpleParametricQSpray
sprays. I will come back to this point later. The
recurrence relation defining the Jacobi polynomials involves a division which
makes the type ParametricQSpray
necessary anyway.
The changeParameters
function is useful to derive the Gegenbauer polynomials
from the Jacobi polynomials. Indeed, as asserted in the Wikipedia article,
the Gegenbauer polynomials coincide, up to a factor, with the Jacobi
polynomials with parameters changeParameters
function to get this special case of Jacobi
polynomials:
import Data.Ratio ( (%) )
j = jacobiPolynomial 3
alpha = qlone 1
alpha' = alpha ^-^ constantSpray (1%2)
j' = changeParameters j [alpha', alpha']
Now let's come back to the conjecture claiming that the coefficients of the
Jacobi polynomials polynomially depend on SimpleParametricQSpray
sprays.
Maybe this can be deduced from a formula given in the Wikipedia article, I
didn't spend some time on this problem. I made this conjecture because I
observed this fact for some small values of canCoerceToSimpleParametricSpray
for other values, which always returned
True
. One can apply the function asSimpleParametricSpray
to perform the
coercion.
There is a third type of parametric sprays in the package, namely the
OneParameterSpray
sprays. The objects of this type represent
multivariate polynomials whose coefficients are fractions
of polynomials in only one variable (the parameter). So they are less
general than the ParametricSpray
sprays.
These sprays are no longer very useful. They have been introduced in version
0.2.5.0 and this is the first type of parametric sprays that has been provided
by the package. When the more general ParametricSpray
sprays have been
introduced, I continued to develop the OneParameterSpray
sprays because they
were more efficient than the univariate ParametricSpray
sprays. But as of
version 0.4.0.0, this is no longer the case. This is what I concluded from
some benchmarks on the Jack polynomials, implemented in the
jackpolynomials package.
These are sprays of type Spray (RatioOfPolynomials a)
, where the type
RatioOfPolynomials a
deals with objects that represent fractions of
univariate polynomials. The type of these univariate polynomials is
Polynomial a
.
The functions we have seen for the simple parametric sprays and the parametric
sprays, substituteParameters
, evalParametricSpray
, and changeParameters
,
are also applicable to the one-parameter sprays.
The OneParameterSpray
sprays were used in the
jackpolynomials package to
represent the Jack polynomials with a symbolic Jack parameter but they have
been replaced with the ParametricSpray
sprays.
Other features offered by the package include: resultant and subresultants of two polynomials, Sturm-Habicht sequence of a polynomial, number of real roots of a univariate polynomial, and greatest common divisor of two polynomials with coefficients in a field.