Title: | Multivariate Polynomials with Rational Coefficients |
---|---|
Description: | Symbolic calculation and evaluation of multivariate polynomials with rational coefficients. This package is strongly inspired by the 'spray' package. It provides a function to compute Gröbner bases (reference <doi:10.1007/978-3-319-16721-3>). It also includes some features for symmetric polynomials, such as the Hall inner product. The header file of the C++ code can be used by other packages. It provides the templated class 'Qspray' that can be used to represent and to deal with multivariate polynomials with another type of coefficients. |
Authors: | Stéphane Laurent [aut, cre], Robin Hankin [ctb, cph] (author of the 'spray' package, which strongly inspired this package) |
Maintainer: | Stéphane Laurent <[email protected]> |
License: | GPL-3 |
Version: | 3.1.0 |
Built: | 2024-10-27 04:04:47 UTC |
Source: | https://github.com/stla/qspray |
Coerces a qspray
polynomial into a function.
## S3 method for class 'qspray' as.function(x, N = FALSE, ...)
## S3 method for class 'qspray' as.function(x, N = FALSE, ...)
x |
object of class |
N |
Boolean, whether the function must numerically approximate the result |
... |
ignored |
A function having the same variables as the polynomial. If
N=FALSE
, it returns a string. If N=TRUE
, it returns a number
if the result does not contain any variable, otherwise it returns a
R expression.
library(qspray) P <- (qlone(1) + "1/2"*qlone(2))^2 + 5 f <- as.function(P) g <- as.function(P, N = TRUE) f(2, "3/7") g(2, "3/7") f("x", "y") g("x", "y") # the evaluation is performed by (R)yacas and complex numbers are # allowed; the imaginary unit is denoted by `I` f("2 + 2*I", "Sqrt(2)") g("2 + 2*I", "Sqrt(2)")
library(qspray) P <- (qlone(1) + "1/2"*qlone(2))^2 + 5 f <- as.function(P) g <- as.function(P, N = TRUE) f(2, "3/7") g(2, "3/7") f("x", "y") g("x", "y") # the evaluation is performed by (R)yacas and complex numbers are # allowed; the imaginary unit is denoted by `I` f("2 + 2*I", "Sqrt(2)") g("2 + 2*I", "Sqrt(2)")
Coercion to a 'qspray' object
## S4 method for signature 'character' as.qspray(x) ## S4 method for signature 'qspray' as.qspray(x) ## S4 method for signature 'numeric' as.qspray(x) ## S4 method for signature 'bigz' as.qspray(x) ## S4 method for signature 'bigq' as.qspray(x)
## S4 method for signature 'character' as.qspray(x) ## S4 method for signature 'qspray' as.qspray(x) ## S4 method for signature 'numeric' as.qspray(x) ## S4 method for signature 'bigz' as.qspray(x) ## S4 method for signature 'bigq' as.qspray(x)
x |
a |
A qspray
object.
as.qspray(2) as.qspray("1/3")
as.qspray(2) as.qspray("1/3")
Replaces the variables of a qspray
polynomial with
some qspray
polynomials. E.g. you have a polynomial
and you want the polynomial
. This is an alias of
composeQspray
.
## S4 method for signature 'qspray,list' changeVariables(x, listOfQsprays)
## S4 method for signature 'qspray,list' changeVariables(x, listOfQsprays)
x |
a |
listOfQsprays |
a list containing at least |
The qspray
polynomial obtained by replacing the variables of
the polynomial given in the x
argument with the polynomials given
in the listOfQsprays
argument.
library(qspray) f <- function(x, y) x*y/2 + 4*y x <- qlone(1) y <- qlone(2) P <- f(x, y) X <- x^2 Y <- x + y + 1 changeVariables(P, list(X, Y)) == f(X, Y) # should be TRUE
library(qspray) f <- function(x, y) x*y/2 + 4*y x <- qlone(1) y <- qlone(2) P <- f(x, y) X <- x^2 Y <- x + y + 1 changeVariables(P, list(X, Y)) == f(X, Y) # should be TRUE
Characteristic polynomial of a matrix.
characteristicPolynomial(A)
characteristicPolynomial(A)
A |
a square matrix with numeric, character, or |
A univariate qspray
polynomial.
set.seed(666) A <- matrix(rpois(9L, 10), nrow = 3, ncol = 3) ( P <- characteristicPolynomial(A) ) # check the roots are the eigen values: f <- as.function(P, N = TRUE) sapply(eigen(A)$values, f) # approx c(0, 0, 0)
set.seed(666) A <- matrix(rpois(9L, 10), nrow = 3, ncol = 3) ( P <- characteristicPolynomial(A) ) # check the roots are the eigen values: f <- as.function(P, N = TRUE) sapply(eigen(A)$values, f) # approx c(0, 0, 0)
Checks whether the polynomials represented by two qspray
objects are collinear, that is, whether they are equal up to a scalar
factor.
collinearQsprays(qspray1, qspray2)
collinearQsprays(qspray1, qspray2)
qspray1 , qspray2
|
two |
A Boolean value.
library(qspray) qspray1 <- qsprayMaker(string = "1/2 x^(1, 1) + 4 x^(0, 2) + 5") qspray2 <- "4/7" * qspray1 collinearQsprays(qspray1, qspray2)
library(qspray) qspray1 <- qsprayMaker(string = "1/2 x^(1, 1) + 4 x^(0, 2) + 5") qspray2 <- "4/7" * qspray1 collinearQsprays(qspray1, qspray2)
Prints a symmetric qspray polynomial as a linear combination of the monomial symmetric polynomials.
## S4 method for signature 'qspray,logical' compactSymmetricQspray(qspray, check) ## S4 method for signature 'qspray,ANY' compactSymmetricQspray(qspray)
## S4 method for signature 'qspray,logical' compactSymmetricQspray(qspray, check) ## S4 method for signature 'qspray,ANY' compactSymmetricQspray(qspray)
qspray |
a |
check |
Boolean, whether to check the symmetry (default |
A character string.
library(qspray) ( qspray <- PSFpoly(4, c(3, 1)) - ESFpoly(4, c(2, 2)) + 4L ) compactSymmetricQspray(qspray, check = TRUE)
library(qspray) ( qspray <- PSFpoly(4, c(3, 1)) - ESFpoly(4, c(2, 2)) + 4L ) compactSymmetricQspray(qspray, check = TRUE)
Substitutes the variables of a qspray
polynomial with
some qspray
polynomials. E.g. you have a polynomial
and you want the polynomial
(see example).
composeQspray(qspray, listOfQsprays)
composeQspray(qspray, listOfQsprays)
qspray |
a |
listOfQsprays |
a list containing at least |
The qspray
polynomial obtained by composing the polynomial
given in the qspray
argument with the polynomials given in the
listOfQsprays
argument.
library(qspray) x <- qlone(1) y <- qlone(2) P <- x*y/2 + 4*y X <- x^2 Y <- x + y + 1 composeQspray(P, list(X, Y)) # this is P(x^2, x+y+1)
library(qspray) x <- qlone(1) y <- qlone(2) P <- x*y/2 + 4*y X <- x^2 Y <- x + y + 1 composeQspray(P, list(X, Y)) # this is P(x^2, x+y+1)
Returns a complete homogeneous symmetric function as a
qspray
polynomial.
CSHFpoly(m, lambda)
CSHFpoly(m, lambda)
m |
integer, the number of variables |
lambda |
an integer partition, given as a vector of decreasing positive integers |
A qspray
object.
library(qspray) CSHFpoly(3, c(3, 1))
library(qspray) CSHFpoly(3, c(3, 1))
Partial derivative of a qspray
polynomial.
derivQspray(qspray, i, derivative = 1)
derivQspray(qspray, i, derivative = 1)
qspray |
object of class |
i |
integer, the dimension to differentiate with respect to, e.g.
|
derivative |
integer, how many times to differentiate |
A qspray
object.
library(qspray) x <- qlone(1) y <- qlone(2) qspray <- 2*x + 3*x*y derivQspray(qspray, 2) # derivative w.r.t. y
library(qspray) x <- qlone(1) y <- qlone(2) qspray <- 2*x + 3*x*y derivQspray(qspray, 2) # derivative w.r.t. y
Partial differentiation of a qspray
polynomial.
dQspray(qspray, orders)
dQspray(qspray, orders)
qspray |
object of class |
orders |
integer vector, the orders of the differentiation; e.g.
|
A qspray
object.
library(qspray) x <- qlone(1) y <- qlone(2) qspray <- x + 2*y + 3*x*y dQspray(qspray, c(1, 1)) derivQspray(derivQspray(qspray, 1), 2)
library(qspray) x <- qlone(1) y <- qlone(2) qspray <- x + 2*y + 3*x*y dQspray(qspray, c(1, 1)) derivQspray(derivQspray(qspray, 1), 2)
Returns an elementary symmetric function as a polynomial.
ESFpoly(m, lambda)
ESFpoly(m, lambda)
m |
integer, the number of variables |
lambda |
an integer partition, given as a vector of decreasing positive integers |
A qspray
object.
library(qspray) ESFpoly(3, c(3, 1))
library(qspray) ESFpoly(3, c(3, 1))
Evaluation of the multivariate polynomial represented by a
qspray
object.
evalQspray(qspray, values_re, values_im = NULL)
evalQspray(qspray, values_re, values_im = NULL)
qspray |
a |
values_re |
vector of the real parts of the values; each element of
|
values_im |
vector of the imaginary parts of the values; each element of
|
A bigq
number if values_im=NULL
, a pair of bigq
numbers otherwise: the real part and the imaginary part of the result.
x <- qlone(1); y <- qlone(2) P <- 2*x + "1/2"*y evalQspray(P, c("2", "5/2", "99999")) # "99999" will be ignored
x <- qlone(1); y <- qlone(2) P <- 2*x + "1/2"*y evalQspray(P, c("2", "5/2", "99999")) # "99999" will be ignored
Get the coefficient of the term with the given monomial.
## S4 method for signature 'qspray,numeric' getCoefficient(qspray, exponents)
## S4 method for signature 'qspray,numeric' getCoefficient(qspray, exponents)
qspray |
a |
exponents |
a vector of exponents, thereby defining a monomial; trailing zeros are ignored |
The coefficient as a bigq
number.
library(qspray) x <- qlone(1) y <- qlone(2) p <- 4*x^2 + 3*y - 5 getCoefficient(p, 2) # coefficient of x^2 getCoefficient(p, c(2, 0)) # same as getCoefficient(p, 2) getCoefficient(p, c(0, 1)) # coefficient of y because y=x^0*y^1 getCoefficient(p, 0) # the constant term getCoefficient(p, integer(0)) # the constant term getCoefficient(p, 3) # there's no x^3
library(qspray) x <- qlone(1) y <- qlone(2) p <- 4*x^2 + 3*y - 5 getCoefficient(p, 2) # coefficient of x^2 getCoefficient(p, c(2, 0)) # same as getCoefficient(p, 2) getCoefficient(p, c(0, 1)) # coefficient of y because y=x^0*y^1 getCoefficient(p, 0) # the constant term getCoefficient(p, integer(0)) # the constant term getCoefficient(p, 3) # there's no x^3
Get the constant term of a qspray
polynomial.
## S4 method for signature 'qspray' getConstantTerm(qspray)
## S4 method for signature 'qspray' getConstantTerm(qspray)
qspray |
a |
A bigq
number.
Returns a Gröbner basis following Buchberger's algorithm using the lexicographical order.
groebner(G, minimal = TRUE, reduced = TRUE)
groebner(G, minimal = TRUE, reduced = TRUE)
G |
a list of |
minimal |
Boolean, whether to return a minimal basis |
reduced |
Boolean, whether to return the reduced basis |
A Gröbner basis of the ideal generated by G
, given as a list
of qspray
polynomials.
Cox, Little & O'Shea. Ideals, Varieties, and Algorithms. An Introduction to Computational Algebraic Geometry and Commutative Algebra. Fourth edition, Springer 2015.
library(qspray) f <- qsprayMaker(string = "x^(3) - 2 x^(1,1)") g <- qsprayMaker(string = "x^(2,1) - 2 x^(0,2) + x^(1)") groebner(list(f, g), FALSE, FALSE) # other example x <- qlone(1); y <- qlone(2); z <- qlone(3) f1 <- x^2 + y + z^2 - 1 f2 <- x^2 + y + z - 1 f3 <- x + y^2 + z - 1 groebner(list(f1, f2, f3))
library(qspray) f <- qsprayMaker(string = "x^(3) - 2 x^(1,1)") g <- qsprayMaker(string = "x^(2,1) - 2 x^(0,2) + x^(1)") groebner(list(f, g), FALSE, FALSE) # other example x <- qlone(1); y <- qlone(2); z <- qlone(3) f1 <- x^2 + y + z^2 - 1 f2 <- x^2 + y + z - 1 f3 <- x + y^2 + z - 1 groebner(list(f1, f2, f3))
Hall inner product of two symmetric polynomials. It has a
parameter alpha
and the standard Hall inner product is the case
when alpha=1
. It is possible to get the Hall inner product with
a symbolic alpha
parameter.
HallInnerProduct(qspray1, qspray2, alpha = 1)
HallInnerProduct(qspray1, qspray2, alpha = 1)
qspray1 , qspray2
|
two symmetric |
alpha |
parameter equal to |
A bigq
number if alpha
is not NULL
, otherwise
a univariate qspray
polynomial.
Implicitization of a system of parametric equations (see example).
implicitization(nvariables, parameters, equations, relations)
implicitization(nvariables, parameters, equations, relations)
nvariables |
number of variables |
parameters |
character vector of the names of the parameters, or
|
equations |
named list of |
relations |
list of |
A list of qspray
polynomials.
library(qspray) # ellipse example #### # variables cost <- qlone(1) sint <- qlone(2) nvariables <- 2 # parameters a <- qlone(3) b <- qlone(4) parameters <- c("a", "b") # equations <- list( "x" = a * cost, "y" = b * sint ) relations <- list( cost^2 + sint^2 - 1 ) # implicitization(nvariables, parameters, equations, relations)
library(qspray) # ellipse example #### # variables cost <- qlone(1) sint <- qlone(2) nvariables <- 2 # parameters a <- qlone(3) b <- qlone(4) parameters <- c("a", "b") # equations <- list( "x" = a * cost, "y" = b * sint ) relations <- list( cost^2 + sint^2 - 1 ) # implicitization(nvariables, parameters, equations, relations)
Returns the exact value of the integral of a multivariate polynomial with rational coefficients over a simplex whose vertices have rational coordinates.
integratePolynomialOnSimplex(P, S)
integratePolynomialOnSimplex(P, S)
P |
a |
S |
the simplex, a |
A bigq
number, the exact value of the integral.
library(qspray) x <- qlone(1); y <- qlone(2) P <- x/2 + x*y S <- rbind(c("0", "0"), c("1", "0"), c("1", "1")) # a triangle integratePolynomialOnSimplex(P, S)
library(qspray) x <- qlone(1); y <- qlone(2) P <- x/2 + x*y S <- rbind(c("0", "0"), c("1", "0"), c("1", "1")) # a triangle integratePolynomialOnSimplex(P, S)
Variables involved in a qspray
object.
## S4 method for signature 'qspray' involvedVariables(x)
## S4 method for signature 'qspray' involvedVariables(x)
x |
a |
A vector of integers. Each integer represents the index of a
variable involved in x
.
x <- qlone(1); z <- qlone(3) involvedVariables(x^2 + x*z + 1) # should be c(1L, 3L)
x <- qlone(1); z <- qlone(3) involvedVariables(x^2 + x*z + 1) # should be c(1L, 3L)
Checks whether a qspray
object defines a constant
polynomial.
## S4 method for signature 'qspray' isConstant(x)
## S4 method for signature 'qspray' isConstant(x)
x |
a |
A Boolean value.
Checks whether the polynomial defined by a qspray
object is homogeneous, and also returns the degree if this is true.
isHomogeneousQspray(qspray)
isHomogeneousQspray(qspray)
qspray |
a |
A Boolean value indicating whether the polynomial defined by
qspray
is homogeneous. Moreover, if it is homogeneous, the degree
is given in the attribute "degree"
of the output.
lambda <- c(3, 2, 1) p <- PSFpoly(4, lambda) ( homogeneous <- isHomogeneousQspray(p) ) # should be TRUE attr(homogeneous, "degree") == sum(lambda) # should be TRUE
lambda <- c(3, 2, 1) p <- PSFpoly(4, lambda) ( homogeneous <- isHomogeneousQspray(p) ) # should be TRUE attr(homogeneous, "degree") == sum(lambda) # should be TRUE
Checks whether a qspray
polynomial can be written as
a polynomial of some given qspray
polynomials. If TRUE
,
this polynomial is returned.
isPolynomialOf(qspray, qsprays)
isPolynomialOf(qspray, qsprays)
qspray |
a |
qsprays |
a list of |
A Boolean value indicating whether the polynomial defined by
qspray
can be written as a polynomial of the polynomials defined
by the qspray
objects given in the qsprays
list. If this is
TRUE
, this polynomial is returned as an attribute named
"polynomial"
.
library(qspray) x <- qlone(1); y <- qlone(2); z <- qlone(3) q1 <- x + y q2 <- x*z^2 + 4 qspray <- q1^2*q2 + 2*q1 + 3 ( check <- isPolynomialOf(qspray, list(q1, q2)) ) POLYNOMIAL <- attr(check, "polynomial") changeVariables(POLYNOMIAL, list(q1, q2)) == qspray # should be TRUE
library(qspray) x <- qlone(1); y <- qlone(2); z <- qlone(3) q1 <- x + y q2 <- x*z^2 + 4 qspray <- q1^2*q2 + 2*q1 + 3 ( check <- isPolynomialOf(qspray, list(q1, q2)) ) POLYNOMIAL <- attr(check, "polynomial") changeVariables(POLYNOMIAL, list(q1, q2)) == qspray # should be TRUE
Checks whether a qspray
object defines the unit
polynomial.
## S4 method for signature 'qspray' isQone(qspray)
## S4 method for signature 'qspray' isQone(qspray)
qspray |
a |
A Boolean value.
Checks whether a qspray
object defines the zero
polynomial.
## S4 method for signature 'qspray' isQzero(qspray)
## S4 method for signature 'qspray' isQzero(qspray)
qspray |
a |
A Boolean value.
Check whether a qspray
polynomial is symmetric.
isSymmetricQspray(qspray)
isSymmetricQspray(qspray)
qspray |
a |
A Boolean value indicating whether the polynomial defined by
qspray
is symmetric.
MSPcombination
, compactSymmetricQspray
e1 <- ESFpoly(3, 1) e2 <- ESFpoly(3, 2) e3 <- ESFpoly(3, 3) q <- e1 + 2*e2 + 3*e3 + 4*e1*e3 isSymmetricQspray(q)
e1 <- ESFpoly(3, 1) e2 <- ESFpoly(3, 2) e3 <- ESFpoly(3, 3) q <- e1 + 2*e2 + 3*e3 + 4*e1*e3 isSymmetricQspray(q)
Checks whether a qspray
object defines a
univariate polynomial.
## S4 method for signature 'qspray' isUnivariate(x)
## S4 method for signature 'qspray' isUnivariate(x)
x |
a |
A Boolean value.
It is considered that a constant qspray
is univariate, and
that the qspray
object y^2+1
where y=qlone(2)
is
not univariate, although only one variable is present (see the note in
the documentation of numberOfVariables
).
Returns the leading coefficient of a qspray
polynomial.
leadingCoefficient(qspray)
leadingCoefficient(qspray)
qspray |
a |
The coefficient of the leading term of qspray
,
a bigq
rational number.
Returns the leading term of a qspray
polynomial.
leadingTerm(qspray)
leadingTerm(qspray)
qspray |
a |
A list providing the exponents of the leading term in the field
powers
, an integer vector, and the coefficient of the leading term
in the field coeff
, a bigq
rational number.
Returns a monomial symmetric function as a polynomial.
MSFpoly(m, lambda)
MSFpoly(m, lambda)
m |
integer, the number of variables |
lambda |
an integer partition, given as a vector of decreasing positive integers |
A qspray
object.
library(qspray) MSFpoly(3, c(3, 1))
library(qspray) MSFpoly(3, c(3, 1))
Expression of a symmetric polynomial as a linear combination of the monomial symmetric polynomials.
MSPcombination(qspray, check = TRUE)
MSPcombination(qspray, check = TRUE)
qspray |
a |
check |
Boolean, whether to check the symmetry |
A list defining the combination. Each element of this list is a
list with two elements: coeff
, a bigq
number, and
lambda
, an integer partition; then this list corresponds to the
term coeff * MSFpoly(n, lambda)
, where n
is the number of
variables in the symmetric polynomial.
qspray <- PSFpoly(4, c(3, 1)) + ESFpoly(4, c(2, 2)) + 4L MSPcombination(qspray)
qspray <- PSFpoly(4, c(3, 1)) + ESFpoly(4, c(2, 2)) + 4L MSPcombination(qspray)
Number of terms of the polynomial defined by a
qspray
object.
## S4 method for signature 'qspray' numberOfTerms(qspray)
## S4 method for signature 'qspray' numberOfTerms(qspray)
qspray |
a |
An integer.
Number of variables involved in a qspray
object (see
the note for the precise meaning).
## S4 method for signature 'qspray' numberOfVariables(x)
## S4 method for signature 'qspray' numberOfVariables(x)
x |
a |
An integer.
The number of variables in the qspray
object y^2+1
where
y=qlone(2)
is 2
, not 1
, although only one variable is
present. Strictly speaking, the function returns the maximal integer
d
such that the variable qlone(d)
occurs in the polynomial.
Reorders the terms of a qspray
object according to the
lexicographic order of the powers. This function is rather used
internally only but it is exported for internal usage in other packages.
orderedQspray(qspray)
orderedQspray(qspray)
qspray |
a |
A qspray
object. It defines the same polynomial as the
input qspray
object but it is ordered.
qspray <- rQspray() qspray == orderedQspray(qspray) # should be TRUE
qspray <- rQspray() qspray == orderedQspray(qspray) # should be TRUE
Permute the variables of a qspray
polynomial.
## S4 method for signature 'qspray,numeric' permuteVariables(x, permutation)
## S4 method for signature 'qspray,numeric' permuteVariables(x, permutation)
x |
a |
permutation |
a permutation |
A qspray
object.
library(qspray) f <- function(x, y, z) { x^2 + 5*y + z - 1 } x <- qlone(1) y <- qlone(2) z <- qlone(3) P <- f(x, y, z) permutation <- c(3, 1, 2) Q <- permuteVariables(P, permutation) Q == f(z, x, y) # should be TRUE
library(qspray) f <- function(x, y, z) { x^2 + 5*y + z - 1 } x <- qlone(1) y <- qlone(2) z <- qlone(3) P <- f(x, y, z) permutation <- c(3, 1, 2) Q <- permuteVariables(P, permutation) Q == f(z, x, y) # should be TRUE
Pretty form of a qspray
polynomial.
prettyQspray(qspray, vars = NULL)
prettyQspray(qspray, vars = NULL)
qspray |
a |
vars |
variable names; |
A character string.
library(qspray) P <- (qlone(1) + "1/2"*qlone(2))^2 + 5 prettyP <- prettyQspray(P, vars = c("x", "y")) prettyP cat(Ryacas::yac_str(sprintf("PrettyForm(%s)", prettyP))) Ryacas::yac_str(sprintf("TeXForm(%s)", prettyP))
library(qspray) P <- (qlone(1) + "1/2"*qlone(2))^2 + 5 prettyP <- prettyQspray(P, vars = c("x", "y")) prettyP cat(Ryacas::yac_str(sprintf("PrettyForm(%s)", prettyP))) Ryacas::yac_str(sprintf("TeXForm(%s)", prettyP))
Returns a power sum function as a polynomial.
PSFpoly(m, lambda)
PSFpoly(m, lambda)
m |
integer, the number of variables |
lambda |
an integer partition, given as a vector of decreasing positive integers |
A qspray
object.
library(qspray) PSFpoly(3, c(3, 1))
library(qspray) PSFpoly(3, c(3, 1))
Expression of a symmetric qspray
polynomial as a
linear combination of some power sum polynomials.
PSPcombination(qspray)
PSPcombination(qspray)
qspray |
a symmetric |
A list of pairs. Each pair is made of a bigq
number, the
coefficient of the term of the linear combination, and an integer
partition, corresponding to a power sum polynomial.
# take a symmetric polynomial ( qspray <- ESFpoly(4, c(2, 1)) + ESFpoly(4, c(2, 2)) ) # compute the power sum combination ( pspCombo <- PSPcombination(qspray) ) # then the polynomial can be reconstructed as follows: Reduce(`+`, lapply(pspCombo, function(term) { term[["coeff"]] * PSFpoly(4, term[["lambda"]]) }))
# take a symmetric polynomial ( qspray <- ESFpoly(4, c(2, 1)) + ESFpoly(4, c(2, 2)) ) # compute the power sum combination ( pspCombo <- PSPcombination(qspray) ) # then the polynomial can be reconstructed as follows: Reduce(`+`, lapply(pspCombo, function(term) { term[["coeff"]] * PSFpoly(4, term[["lambda"]]) }))
Expression of a symmetric qspray
polynomial as a
polynomial in the power sum polynomials.
PSPexpression(qspray)
PSPexpression(qspray)
qspray |
a symmetric |
A qspray
polynomial, say , such that
equals the input symmetric polynomial,
where
is the i-th power sum polynomial (
PSFpoly(n, i)
).
# take a symmetric polynomial ( qspray <- ESFpoly(4, c(2, 1)) + ESFpoly(4, c(2, 2)) ) # compute the power sum expression ( pspExpr <- PSPexpression(qspray) ) # take the involved power sum polynomials psPolys <- lapply(1:numberOfVariables(pspExpr), function(i) PSFpoly(4, i)) # then this should be TRUE: qspray == changeVariables(pspExpr, psPolys)
# take a symmetric polynomial ( qspray <- ESFpoly(4, c(2, 1)) + ESFpoly(4, c(2, 2)) ) # compute the power sum expression ( pspExpr <- PSPexpression(qspray) ) # take the involved power sum polynomials psPolys <- lapply(1:numberOfVariables(pspExpr), function(i) PSFpoly(4, i)) # then this should be TRUE: qspray == changeVariables(pspExpr, psPolys)
Division of a qspray polynomial by a list of qspray polynomials. See the reference for the definition.
qdivision(qspray, divisors)
qdivision(qspray, divisors)
qspray |
the dividend, a |
divisors |
the divisors, a list of |
The remainder of the division, a qspray
object.
Michael Weiss, 2010. Computing Gröbner Bases in Python with Buchberger’s Algorithm.
# a univariate example library(qspray) x <- qlone(1) f <- x^4 - 4*x^3 + 4*x^2 - x # 0 and 1 are trivial roots g <- x * (x - 1) qdivision(f, list(g)) # should be zero
# a univariate example library(qspray) x <- qlone(1) f <- x^4 - 4*x^3 + 4*x^2 - x # 0 and 1 are trivial roots g <- x * (x - 1) qdivision(f, list(g)) # should be zero
Creates a polynomial variable. Using this function is the main
way to build qspray
objects.
qlone(n)
qlone(n)
n |
positive integer, the index of the variable |
A qspray
object.
x <- qlone(1) y <- qlone(2) (x + y) * (x - y)
x <- qlone(1) y <- qlone(2) (x + y) * (x - y)
Returns the qspray
polynomial identically equal to 1.
qone()
qone()
A qspray
object.
This function is for internal usage. It is exported because it is also used for internal usage in others packages.
qspray_from_list(qspray_as_list)
qspray_from_list(qspray_as_list)
qspray_as_list |
list returned by the Rcpp function
|
A qspray
object.
Unary operators for qspray objects.
## S4 method for signature 'qspray,missing' e1 + e2 ## S4 method for signature 'qspray,missing' e1 - e2
## S4 method for signature 'qspray,missing' e1 + e2 ## S4 method for signature 'qspray,missing' e1 - e2
e1 |
object of class |
e2 |
nothing |
A qspray
object.
Division of two polynomials
qsprayDivision(qsprayA, qsprayB)
qsprayDivision(qsprayA, qsprayB)
qsprayA |
a |
qsprayB |
a |
A list with two qspray
objects, the quotient and the
remainder.
library(qspray) x <- qlone(1) y <- qlone(2) z <- qlone(3) B <- x*y^2 + z*x^2 + 1 A <- B * (x^2*y^2*z^2 - 3) + x*y divis <- qsprayDivision(A, B) B * divis[["Q"]] + divis[["R"]] == A # should be TRUE
library(qspray) x <- qlone(1) y <- qlone(2) z <- qlone(3) B <- x*y^2 + z*x^2 + 1 A <- B * (x^2*y^2*z^2 - 3) + x*y divis <- qsprayDivision(A, B) B * divis[["Q"]] + divis[["R"]] == A # should be TRUE
Make a qspray
object from a list of exponents and a
vector of coefficients.
qsprayMaker(powers, coeffs, string = NULL)
qsprayMaker(powers, coeffs, string = NULL)
powers |
list of positive integer vectors |
coeffs |
a vector such that each element of |
string |
if not |
A qspray
object.
powers <- list(c(1, 1), c(0, 2)) coeffs <- c("1/2", "4") qsprayMaker(powers, coeffs) qsprayMaker(string = "1/2 x^(1, 1) + 4 x^(0, 2)")
powers <- list(c(1, 1), c(0, 2)) coeffs <- c("1/2", "4") qsprayMaker(powers, coeffs) qsprayMaker(string = "1/2 x^(1, 1) + 4 x^(0, 2)")
Returns the qspray
polynomial identically equal to 0.
qzero()
qzero()
A qspray
object.
Generates a random qspray
object.
rQspray()
rQspray()
A qspray
object with at most 4 terms and at most 3 variables.
Prints a monomial like "x^(1, 0, 2)"
. This way of
showing a monomial was used by default in previous versions of this
package.
showMonomialOld(x = "x")
showMonomialOld(x = "x")
x |
a string, usually a letter such as |
A function which takes as argument a sequence of exponents and which prints the corresponding monomial.
showMonomialX1X2X3
, showMonomialXYZ
,
showQspray
, showQsprayOption<-
.
showMonomialOld("X")(c(1, 0, 2)) showMonomialOld("X")(NULL)
showMonomialOld("X")(c(1, 0, 2)) showMonomialOld("X")(NULL)
Prints a monomial in the style of "x1.x3^2"
.
showMonomialX1X2X3(x = "x", collapse = ".")
showMonomialX1X2X3(x = "x", collapse = ".")
x |
a string, usually a letter such as |
collapse |
a string to denote the symbol representing the
multiplication, e.g. |
A function which takes as argument a sequence of exponents and which prints the corresponding monomial.
The function returned by this function can be used as the option
"showMonomial"
in the showQsprayOption<-
function.
But if you are happy with the default collapse
argument, then you
can equivalently set the "x"
option instead, thereby typing less
code. See the example.
showQsprayX1X2X3
,
showMonomialXYZ
, showQsprayOption<-
.
showMonomialX1X2X3("X")(c(1, 0, 2)) showMonomialX1X2X3("X", collapse = "*")(c(1, 0, 2)) showMonomialX1X2X3("X")(c(1, 0, 2)) == showMonomialXYZ(c("X1", "X2", "X3"))(c(1, 0, 2)) showMonomialX1X2X3()(NULL) # setting a show option: set.seed(3141) ( qspray <- rQspray() ) showQsprayOption(qspray, "showMonomial") <- showMonomialX1X2X3("X") qspray # this is equivalent to: showQsprayOption(qspray, "showQspray") <- showQsprayX1X2X3("X") # and also equivalent to: showQsprayOption(qspray, "x") <- "X"
showMonomialX1X2X3("X")(c(1, 0, 2)) showMonomialX1X2X3("X", collapse = "*")(c(1, 0, 2)) showMonomialX1X2X3("X")(c(1, 0, 2)) == showMonomialXYZ(c("X1", "X2", "X3"))(c(1, 0, 2)) showMonomialX1X2X3()(NULL) # setting a show option: set.seed(3141) ( qspray <- rQspray() ) showQsprayOption(qspray, "showMonomial") <- showMonomialX1X2X3("X") qspray # this is equivalent to: showQsprayOption(qspray, "showQspray") <- showQsprayX1X2X3("X") # and also equivalent to: showQsprayOption(qspray, "x") <- "X"
Prints a monomial like "x.z^2"
if possible (see details).
showMonomialXYZ(letters = c("x", "y", "z"), collapse = ".")
showMonomialXYZ(letters = c("x", "y", "z"), collapse = ".")
letters |
a vector of strings, usually some letters such as |
collapse |
a string to denote the symbol representing the
multiplication, e.g. |
If the function returned by this function is applied to a vector
of exponents whose length is higher than the length of the letters
vector, then showMonomialX1X2X3(x=letters[1])
is applied
(see the last example).
A function which takes as argument a sequence of exponents and which prints the corresponding monomial.
The function returned by this function can be used as the option
"showMonomial"
in the showQsprayOption<-
function.
showQsprayXYZ
,
showMonomialX1X2X3
, showQsprayOption<-
.
showMonomialXYZ()(c(1, 0, 2)) showMonomialXYZ(collapse = "*")(c(1, 0, 2)) showMonomialXYZ()(NULL) # what happens if there are more exponents than letters: showMonomialXYZ(c("a", "b"), collapse = "*")(c(1, 2, 3)) # same as: showMonomialX1X2X3("a", collapse = "*")(c(1, 2, 3)) # setting a show option: set.seed(3141) ( qspray <- rQspray() ) showQsprayOption(qspray, "showMonomial") <- showMonomialXYZ(c("A", "B", "C")) qspray # this is equivalent to: showQsprayOption(qspray, "showQspray") <- showQsprayXYZ(c("A", "B", "C"))
showMonomialXYZ()(c(1, 0, 2)) showMonomialXYZ(collapse = "*")(c(1, 0, 2)) showMonomialXYZ()(NULL) # what happens if there are more exponents than letters: showMonomialXYZ(c("a", "b"), collapse = "*")(c(1, 2, 3)) # same as: showMonomialX1X2X3("a", collapse = "*")(c(1, 2, 3)) # setting a show option: set.seed(3141) ( qspray <- rQspray() ) showQsprayOption(qspray, "showMonomial") <- showMonomialXYZ(c("A", "B", "C")) qspray # this is equivalent to: showQsprayOption(qspray, "showQspray") <- showQsprayXYZ(c("A", "B", "C"))
Prints a qspray
object given a function which prints
the monomials.
showQspray(showMonomial, compact = FALSE, multiplication = "*")
showQspray(showMonomial, compact = FALSE, multiplication = "*")
showMonomial |
a function which takes as argument a sequence of exponents and which returns a string representing the corresponding monomial |
compact |
a Boolean value; if |
multiplication |
used to separate the coefficient and the monomial within a term |
A function which prints a qspray
object.
The function returned by this function can be used as the option
"showQspray"
in the showQsprayOption<-
function.
But one generally prefers to use showQsprayX1X2X3
or
showQsprayXYZ
instead, which are both built with
showQspray
.
showQsprayX1X2X3
, showQsprayXYZ
,
showQsprayOption<-
.
set.seed(3141) ( qspray <- rQspray() ) f <- showQspray(showMonomialX1X2X3("X"), compact = TRUE) f(qspray) # this is equivalent to: f <- showQsprayX1X2X3("X", compact = TRUE) f(qspray) # if you want to adopt this way to show a qspray, use # the setter function \code{\link{showQsprayOption<-}}: showQsprayOption(qspray, "showQspray") <- showQsprayX1X2X3("X", compact = TRUE) qspray # then this show option will be preserved by some operations on the qspray: qspray^2
set.seed(3141) ( qspray <- rQspray() ) f <- showQspray(showMonomialX1X2X3("X"), compact = TRUE) f(qspray) # this is equivalent to: f <- showQsprayX1X2X3("X", compact = TRUE) f(qspray) # if you want to adopt this way to show a qspray, use # the setter function \code{\link{showQsprayOption<-}}: showQsprayOption(qspray, "showQspray") <- showQsprayX1X2X3("X", compact = TRUE) qspray # then this show option will be preserved by some operations on the qspray: qspray^2
Set a show option to a qspray
object
showQsprayOption(x, which) <- value
showQsprayOption(x, which) <- value
x |
a |
which |
which option to set; this can be |
value |
the value of the option |
This returns the updated qspray
.
The interest of setting some show options to a 'qspray' is that these options are preserved by some operations. See the examples and the README.
set.seed(3141) ( qspray <- rQspray() ) showQsprayOption(qspray, "x") <- "a" qspray # this is identical to: showQsprayOption(qspray, "showMonomial") <- showMonomialX1X2X3("a") # and also identical to: showQsprayOption(qspray, "showQspray") <- showQsprayX1X2X3("a") # old show method: showQsprayOption(qspray, "showMonomial") <- showMonomialOld() qspray # show options are preserved by some operations: qspray^2 3*qspray derivQspray(qspray, 1) swapVariables(qspray, 1, 2) substituteQspray(qspray, c(NA, NA, "3/2")) # for the binary arithmetic operations, the show options of the first # operand are transferred to the result when possible: ( qspray2 <- rQspray() ) qspray + qspray2
set.seed(3141) ( qspray <- rQspray() ) showQsprayOption(qspray, "x") <- "a" qspray # this is identical to: showQsprayOption(qspray, "showMonomial") <- showMonomialX1X2X3("a") # and also identical to: showQsprayOption(qspray, "showQspray") <- showQsprayX1X2X3("a") # old show method: showQsprayOption(qspray, "showMonomial") <- showMonomialOld() qspray # show options are preserved by some operations: qspray^2 3*qspray derivQspray(qspray, 1) swapVariables(qspray, 1, 2) substituteQspray(qspray, c(NA, NA, "3/2")) # for the binary arithmetic operations, the show options of the first # operand are transferred to the result when possible: ( qspray2 <- rQspray() ) qspray + qspray2
Prints a qspray
object given a string for the variable.
showQsprayX1X2X3(x = "x", collapse = ".", ...)
showQsprayX1X2X3(x = "x", collapse = ".", ...)
x , collapse
|
|
... |
arguments passed to |
A function which prints a qspray
object.
The way qspray
objects are displayed can be controlled with the
help of the function showQsprayOption<-
, and
showQsprayX1X2X3()
is a possible option to pass in
showQsprayOption<-
.
showMonomialX1X2X3
, showQspray
,
showQsprayOption<-
.
set.seed(3141) ( qspray <- rQspray() ) showQsprayX1X2X3("X")(qspray) # setting a show option: showQsprayOption(qspray, "showQspray") <- showQsprayX1X2X3("A") qspray # this is equivalent to: showQsprayOption(qspray, "showMonomial") <- showMonomialX1X2X3("A") # and also equivalent to: showQsprayOption(qspray, "x") <- "A"
set.seed(3141) ( qspray <- rQspray() ) showQsprayX1X2X3("X")(qspray) # setting a show option: showQsprayOption(qspray, "showQspray") <- showQsprayX1X2X3("A") qspray # this is equivalent to: showQsprayOption(qspray, "showMonomial") <- showMonomialX1X2X3("A") # and also equivalent to: showQsprayOption(qspray, "x") <- "A"
Prints a polynomial by printing monomials like "x^2.yz"
.
showQsprayXYZ(letters = c("x", "y", "z"), collapse = ".", ...)
showQsprayXYZ(letters = c("x", "y", "z"), collapse = ".", ...)
letters , collapse
|
see |
... |
arguments passed to |
A function which prints a qspray
object. It is constructed
with showQspray
and showMonomialXYZ
.
The function returned by this function can be used as the option
"showQspray"
in the showQsprayOption<-
function.
showMonomialXYZ
, showQspray
,
showQsprayOption<-
.
set.seed(3141) ( qspray <- rQspray() ) showQsprayXYZ(c("X", "Y", "Z"))(qspray) showQsprayXYZ(c("X", "Y", "Z"))(qlone(1) + qlone(2) + qlone(3) + qlone(4)) # setting a show option: showQsprayOption(qspray, "showQspray") <- showQsprayXYZ(c("A", "B", "C")) qspray # this is equivalent to: showQsprayOption(qspray, "showMonomial") <- showMonomialXYZ(c("A", "B", "C"))
set.seed(3141) ( qspray <- rQspray() ) showQsprayXYZ(c("X", "Y", "Z"))(qspray) showQsprayXYZ(c("X", "Y", "Z"))(qlone(1) + qlone(2) + qlone(3) + qlone(4)) # setting a show option: showQsprayOption(qspray, "showQspray") <- showQsprayXYZ(c("A", "B", "C")) qspray # this is equivalent to: showQsprayOption(qspray, "showMonomial") <- showMonomialXYZ(c("A", "B", "C"))
Substitute some variables in a qspray
polynomial.
substituteQspray(qspray, values)
substituteQspray(qspray, values)
qspray |
a |
values |
the values to be substituted; this must be a vector whose
length equals the number of variables of |
A qspray
object.
library(qspray) x <- qlone(1) y <- qlone(2) z <- qlone(3) p <- x^2 + y^2 + x*y*z - 1 substituteQspray(p, c("2", NA, "3/2"))
library(qspray) x <- qlone(1) y <- qlone(2) z <- qlone(3) p <- x^2 + y^2 + x*y*z - 1 substituteQspray(p, c("2", NA, "3/2"))
Swap two variables in a qspray
polynomial.
## S4 method for signature 'qspray,numeric,numeric' swapVariables(x, i, j)
## S4 method for signature 'qspray,numeric,numeric' swapVariables(x, i, j)
x |
a |
i , j
|
indices of the variables to be swapped |
A qspray
object.
library(qspray) f <- function(x, y, z) { x^2 + 5*y + z - 1 } x <- qlone(1) y <- qlone(2) z <- qlone(3) P <- f(x, y, z) Q <- swapVariables(P, 2, 3) Q == f(x, z, y) # should be TRUE
library(qspray) f <- function(x, y, z) { x^2 + 5*y + z - 1 } x <- qlone(1) y <- qlone(2) z <- qlone(3) P <- f(x, y, z) Q <- swapVariables(P, 2, 3) Q == f(x, z, y) # should be TRUE