Polynomial Numerical Linear Algebra package for Matlab©/Octave©
Multi-functional package for solving all kinds of problems with multivariate polynomials in double precision.
Author: Kim Batselier
- Basics
Only 1 monomial ordering is supported: the degree negative lex ordering. This graded ordering is defined for two n-variate monomials x^a and x^b as
a > b
if
|a| = \sum_i^n a_i > |b| = \sum_i^n b_i or
|a|=|b| and the leftmost nonzero entry of a-b is negative.
In this definition x^a stands for
x^a = x_1^a_1 * x_2^a_2 * ... * x_n^a_n,
where a is an n-tuple of nonnegative integers.
A system of s multivariate polynomials is represented by an s-by-2 cell. The first column elements are vectors containing the coefficients. The second column elements are matrices containing the corresponding exponents. For example, the system
f_1 = 5.3 x_1^2 + 9 x_2 x_3 -1,
f_2 = 2 x_1^3 + .5 x_2^2 - 7.89 x_3 - 94,
f_3 = x_1 - 2.13,
is represented by
polysys{1,1}=[5.3,9,-1]
polysys{1,2}=[2 0 0;0 1 1;0 0 0]
polysys{2,1}=[2,.5,-7.89,-94]
polysys{2,2}=[3 0 0;0 2 0;0 0 1;0 0 0]
polysys{3,1}=[1,-2.13]
polysys{3,2}=[1 0 0;0 0 0]
- Available functions
This is an incomplete list of all available functions. Most functions have more input and output arguments than specified here. Type 'help X' in MATLAB/OCTAVE to get a full description of the function X.
- index = feti(exponent)
Converts an exponent to its corresponding linear index with respect to the degree negative lex ordering. If the exponent argument is a matrix of exponents, then feti() is applied to each row of the matrix.
- exponent = fite(index,n)
Converts a linear index with respect to the degree negative lex ordering to the corresponding n-variate exponent.
- d0 = getD0(polysys)
Returns the maximal total degree of the given polynomial system.
- dorig = getDorig(polysys)
dorig(i) (i=1:s) contains the maximal total degree of the multivariate polynomial corresponding with polysys[i,1],polysys[i,2].
- [p,q] = getMDim(polysys,d)
Returns the number of columns p and the number of rows q of the Macaulay matrix of degree d.
- [M nzmax] = getM(polysys,d,sparseM)
Returns the Macaulay matrix of degree d for the given polynomial system polysys. sparseM is an optional boolean variable, set to 1 to use a column compressed data structure to store the matrix, default sparseM is 0.
- Mex = getMex(polysys,d,dmin,sparseM)
Returns the extra rows that need to be added to M(dmin) in order to construct M(d). Suitable for recursive updating of M.
- [q r] = polyDiv(psys,polysys)
Polynomial division of a multivariate polynomial psys by a set of multivariate divisors polysys, computed using oblique projections. Returns the coefficient vectors of the quotient q and remainder r.
- [a b R] = candecomp(polysys,d)
Computes the canonical decomposition of the vector space of n-variate polynomials up to degree d for a given polynomial system polysys. Vectors a and b contain the indices of the linearly independent and dependent monomials (degree negative lex ordered). Each row of R is a coefficient vector of a reduced polynomial.
- [a b R] = rcandecomp(polysys,d)
Computes the reduced canonical decomposition of the vector space of n-variate polynomials up to degree d for a given polynomial system polysys. Vectors a and b contain the indices of the reduced linearly independent and dependent monomials (degree negative lex ordered).
- [p d] = elim(polysys,var)
Returns the coefficient vector of a polynomial p that lies in the ideal spanned by the polynomials in polysys and from which all variables var are eliminated. Var is a a vector of indices, e.g. suppose n=5 and var=[1,2,5] then p will be a polynomial in the variables x3 and x4.
- puni = punivar(polysys,var)
Returns the univariate (in the variable var) polynomial puni that lies in the polynomial ideal generated by polysys.
- [dp dm] = aln(polysys,dmax)
Analyzes the left null space of the Macaulay matrix of polysys up to a degree dmax. The variables dp and dm contains the degrees of the binomial coefficients that generate the Hilbert Function. Each entry k of dp contributes a term nchoosek(d-k+n,n) to the Hilbert Function and likewise -nchoosek(d-k+n,n) for an entry of dm.
- [lcm h] = getLCM(fsys,gsys,tol)
Returns the approximate least common multiple (lcm) of two given polynomials fsys and gsys. Also the factor h such that lcm = g*h is calculated. The optional variable tol is to decide when a principal angle is numerically zero.
- [gcd e] = getGCD(psys,qsys,tol)
Returns the approximate greatest common divisor (gcd) of two given polynomials psys and qsys. The optional variable tol is to decide when a principal angle is numerically zero.
- [Nnew tol] = updateN(N,M,sparse)
Updates the numerical kernel N of M(d) using the recursive orthogonalization theorem. M contains the rows that need to be added to M(d) to form M(d+1).
- [Unew Nnew tol] = updateOrth(U,N,M,sparse)
Updates the numerical basis for the row space U and kernel N of M(d) using the recursive orthogonalization theorem. M contains the rows that need to be added to M(d) to form M(d+1).
- [root d c ns check cr digits] = qdsparf(polysys)
Quick & Dirty Sparse Affine Root-Finding algorithm. Uses SuiteSparseQR. Computes all affine roots of a given polynomial system polysys.
- [root d c ns check cr digits] = sparf(polysys)
Sparse affine root-finding algorithm. Try the quick & dirty algorithm first, it will be faster.