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.zThis 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.CNote 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.A3This 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.zOr 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^2We 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.0import 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 + 6import 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^**^2The 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 gbasisThen 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 -- TrueIn 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' -- TrueNote 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)*ZYou 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 }*ZThe 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 }*XLet'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 -- 0Now, 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)*XThis 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)*alphaThis 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 -- TrueActually, 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 *> rOSYou can also divide a ratio of sprays by a spray with %/%:
p *> (rOS' %/% p) == rOS' -- True rOS1 %/% p == p %//% q -- TrueWhen 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) -- TrueUse 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 -- TrueRecall 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.