library(jack)
library(microbenchmark)
Schur polynomials have applications in combinatorics and zonal polynomials have applications in multivariate statistics. They are particular cases of Jack polynomials. This package allows to evaluate these polynomials. It can also compute their symbolic form.
The functions JackPol
, ZonalPol
,
ZonalQPol
and SchurPol
respectively return the
Jack polynomial, the zonal polynomial, the quaternionic zonal
polynomial, and the Schur polynomial.
Each of these polynomials corresponds is given by a positive integer,
the number of variables, and an integer partition, the
lambda
argument; the Jack polynomial has one more
parameter, the alpha
argument, a positive number.
To get an exact symbolic polynomial with JackPol
, you
have to supply a bigq
rational number for the parameter
alpha
:
<- JackPol(2, lambda = c(3, 1), alpha = gmp::as.bigq("2/5"))
jpol
jpol## 98/25*x^(3, 1) + 98/25*x^(1, 3) + 28/5*x^(2, 2)
This is a qspray
object, from the qspray
package. Here is how you can evaluate this polynomial:
::evalQspray(jpol, c("2", "3/2"))
qspray## Big Rational ('bigq') :
## [1] 1239/10
By default, ZonalPol
, ZonalQPol
and
SchurPol
return exact symbolic polynomials.
<- ZonalPol(2, lambda = c(3, 1))
zpol
zpol## 24/7*x^(3, 1) + 24/7*x^(1, 3) + 16/7*x^(2, 2)
It is also possible to convert a qspray
polynomial to a
function whose evaluation is performed by the Ryacas
package:
<- as.function(zpol) zyacas
You can provide the values of the variables of this function as numbers or character strings:
zyacas(2, "3/2")
## [1] "594/7"
You can even pass a variable name to this function:
zyacas("x", "x")
## [1] "(64*x^4)/7"
If you want to substitute a variable with a complex number, use a
character string which represents this number, with I
denoting the imaginary unit:
zyacas("2 + 2*I", "2/3")
## [1] "Complex((-2176)/63,2944/63)"
As of version 2.0.0, the Jack polynomials can be calculated with Julia. The speed is amazing:
<- Jack_julia()
julia ## Starting Julia ...
<- c(1/2, 2/3, 1, 2/3, -1, -2, 1)
x <- c(5, 3, 2, 2, 1)
lambda <- 3
alpha print(
microbenchmark(
R = Jack(x, lambda, alpha),
Julia = julia$Jack(x, lambda, alpha),
times = 6L,
unit = "seconds"
),signif = 6L
)## Unit: seconds
## expr min lq mean median uq max neval cld
## R 7.73008000 7.8626500 7.903070 7.9140600 7.9534400 8.04412 6 b
## Julia 0.00338613 0.0459389 0.216757 0.0697467 0.0948473 1.01688 6 a
Jack_julia()
returns a list of functions.
ZonalPol
, ZonalQPol
and SchurPol
always return an exact expression of the polynomial, i.e. with rational
coefficients (integers for SchurPol
). If you want an exact
expression with JackPol
, you have to give a rational number
for the argument alpha
, as a character string:
<- julia$JackPol(m = 2, lambda = c(3, 1), alpha = "2/5")
JP
JP## 98/25*x^(1, 3) + 28/5*x^(2, 2) + 98/25*x^(3, 1)
Again, Julia is faster:
<- 5
n <- c(4, 3, 3)
lambda <- "2/3"
alpha <- gmp::as.bigq(alpha)
alphaq microbenchmark(
R = JackPol(n, lambda, alphaq),
Julia = julia$JackPol(n, lambda, alpha),
times = 6L
)## Unit: milliseconds
## expr min lq mean median uq max neval cld
## R 1170.8920 1188.7068 1225.4087 1211.0391 1236.0073 1334.7680 6 b
## Julia 464.3045 491.7447 549.3006 496.9847 546.5387 799.2461 6 a
::stopJulia() JuliaConnectoR