[  Download   ]   [ Table of Contents ]   [ Key features  ]  [ Publications  ]     [  Contact   ]    [  About  ]

NAClab -- Numerical Algebraic Computing Toolbox for Matlab
Copyright © 2012 by Zhonggang Zeng.  All rights reserved.  Updated on August 18, 2018

Unique capabilities:  Solving hypersensitive problems from empirical data
                                   Intuitive WYSIWYG polynomial computation interface in Matlab
                                   Solving linear operator equation  L(z) = b directly from linear transformation L even if it is rank-deficient
                          New!  Rank-r Newton's iteration solving nonlinear systems  f(x) = 0 with non-isolated solutions and singular Jacobians
                                   solving nonlinear least squares problem with Gauss-Newtion iterations.

Featured modules:     Polynomial systems, Jordan Canonical Forms, multiplicities of nonlinear systems, 
                                                     polynomial GCD / factorizations, ....
Research toolbox:       Generic linear equation solver without matrices, Gauss-Newton iteration, rank-r Newton's iteration ....
Software building blocks  [Link]

NAClab a joint research project in development by Tien-Yien Li and Zhonggang Zeng  It is an expansion from its predecessor Apalab (and its Maple counterpart is  ApaTools).  .  Disclaimer:   NAClab is an on-going research project,  Bugs and shortcomings likely exist.  The authors are not responsible for any damage  or loss from using this software. Users should verify their results through whatever means.

Project leaders:   Tien-Yien Li  and  Zhonggang Zeng  (Contact via email:  zzeng at  neiu dot edu )
Colaborators/contributors:  Liping Chen, Tianran Chen, Wenrui Hao, Tsung-Lin Lee, Wenyuan Wu, Andrew Sommese

!!!New!!!Solving linear equations  L(z) = b  
                - directly from linear transformation L without representation matrix
                - including homogeneous equation  L(z) = 0  
                - even if  L is rank-deficient  

NAClab can directly solve a general linear equation   L(z) = b  for its general numerical solution  where L is a linear transformation between (products of) vector spaces. As a special case, the module LinearSolve solves L(z) = 0 for the numerical kernel.    [ A short paper on  LinearSolve ]

A general numerical solution within an error tolerance is in the form of 
                                   z = z0 + c1*z1+...+ck*zk  
where z0 is the minimum norm solution, c1,..,ck are constants and {z1,...,zk} forms a basis for the numerical kernel of L within the error tolerance

For a demo, NAClab outputs polynomials f and g  such that    p f + q g  = u   where  p, q  are given polynomials of degrees m and n  respectivly and  u = GCD(p,q) of degree k. In other words, NAClab solves the linear equation
                             L(f,g)  = u
where   L is the linear transformation defined by  L(f,g)  = p f + q g  with fixed parameters p, q.  Users do not need to provide representation matrices for the linear transformations unless they choose to do so.  NAClab generates representation matrices internally. 

The complete process for this example is as follows

Define an annonymous function for the linear transformation 

>> map2gcd = @(f,g,p,q) pplus(ptimes(p,f),ptimes(q,g));

Enter polynomials f, g, u as character strings, provide the domain and parameter for the linear transformation:

>> f = '5.99999 - 9*x + 3*x^3 - 4*x^4 + 6*x^5 - 2*x^7';
>> g = '-2 + 5*x - 3*x^2 - x^3 + x^4 + 2*x^5 - 3.00002*x^6 + x^8';
>> u = '2-2.99999*x+1.00001*x^3';
>> domain = {'1+x+x^2+x^3+x^4+x^5','1+x+x^2+x^3+x^4'};
>> parameter = {f,g};

Just like that, we are ready to call LinearSolve with error tolerance 1e-4

>> [Z,K,lcond,res] = LinearSolve({map2gcd,domain,parameter}, u, 1e-4);

The minimum norm solution:

>> Z

Z =

    '0.315213638005548 + 0.0325895733429298*x + 0.021703…'     '-0.054358027177057 + 0.0434095337991119*x + 0.10852…'

The numerical kernel within the error tolerance 1e-4:

>> K{:}

ans =

    '0.249999389461494 - 0.250001189336557*x - 0.2500021…'     '0.74999762066815 - 0.500002203566776*x^4'

!!!New!!!Rank-r Newton's iteration (preprint with examples, syntax)
                - solve nonlinear system  f(x) = 0 even if the solution is non-isolated
                - directly from smooth mapping f between products of vector spaces  
                - even if  the Jacobian is rank-deficient  

 Solving polynomial systems
             by the homotopy method based on HOM4PS

Polynomials systems can be solved  for numerical solutions in a straightforward and intuitive setting.  For example: To solve the following polynomial system with 20 zeros

     -x^5 + y^5 - 3*y - 1 = 0
     5*y^4-3                    = 0
      -20*x+y-z                = 0

simply call "psolve" with straitforward input [ more details ]

>> P = {'-x^5+y^5-3*y-1','5*y^4-3','-20*x+y-z'};  % define the polynomial system
>> [Solutions, variables] = psolve(P)             % call the polynomial system solver

   Solutions =

     Column 1

     0.778497746685646 + 0.893453081179308i
    -0.000000000000000 + 0.880111736793394i
   -15.569954933712914 -16.988949886792764i

     ... ...

     Column 20

     0.778497746685646 - 0.893453081179308i
     0.000000000000000 - 0.880111736793393i
   -15.569954933712925 +16.988949886792767i


   variables =

    'x'    'y'    'z'

 Direct polynomial manipulation:  

Polynomials can now be entered and shown as character strings, such as  
   >> f = '3*x^2 - (2-5i)*x^3*y^4 - 1e-3*z^5-6.8'
   >> g = '-2*y^3 - 5*x^2*z + 8.2'

Using polynomials strings as input, users can now perform common polynomial operations  such as addition, multiplication, power, evaluation, differentiation, factorization, extracting coefficients, finding greatest common divison (GCD), by calling NAClab functions, such as
   >> p = PolynomialPlus('2*x^5-3*y','4+x*y')     % add any number of polynomials
   >> q = PolynomialTimes(f,g,h)  % multiply any number of polynomials
   >> v = PolynomialEvaluate(f,{'x','z'},[3,4]) % evaluate f(x,y,z) for x=3, z=4
and a lot more.  For example, to compute a greatest common divisor:

   >> f = '10 - 5*x^2*y + 6*x*y^2 - 3*x^3*y^3';
   >> g = '30 + 10*y + 18*x*y^2 + 6*x*y^3'
   >> u = PolynomialGCD(f,g)
   u =
   33.5410196624968 + 20.1246117974981*x*y^2

 KeyKey features:

The functions for
are major development in both algorithm design and software implementations. For instance, the univariate factorization is based on an award winning paper on accurate computation of multiple roots

Table of contents


List of functions:

     Functions are in three categories:  Queck-access functions,  Programming tools for polynomial compuations,  and  Matrix compuation tools

Quick-access polynomial functions:

Functions in this category accept input polynomials as character strings, such as
   >> f = '3*x^2 - (2-5i)*x^3*y^4 - 1e-3*z^5-6.8'
   >> g = '-2*y^3 - 5*x^2*z + 8.2'
A list of quick access functions are as follows.  Use, say >> help GetVariables, to access documentations.  Many functions have shortened aliases. For example, pvar is an alias for GetVariables.  

Programming tools for polynomial computations:

Functions in this category are designed for experts for developing algorithms and implementations for numerical polynomial algebra. These functions requires polynomials represented as either a coefficient matrix or a coefficient vector.  

A univariate polynomial is represented by its coefficient vector.  For example,  f(x) = x3 + 2x2 + 3x + 4  is represented as
               >> f  =  [1, 2, 3, 4]   % or  its transpose

An  m-multivariate polynomial  f  of  n  terms is represented by a coefficient matrix  F  of size  (m+1) by  n  according to the following rules:
     1.  Every column of  F  represents one term (monomial)
     2.  F(i,j)    for i<m+1   is the exponent of the i-th variable on j-th term
     3.  F(m+1,j)  is the coefficient of the j-th term
For example, let    p(x, y, z)  =  8.5  +  (3-2i)x3 y  +  5 x2 z5   -  2 y3   + 6 y3 z2  .   Users can transform its string representation to coefficient matrix representation by calling the quick access function PolynString2CoefMat :  
   >> f = '8.5 + (3-2i)*x^3*y + 5*x^2*z^5 - 2*y^3 + 6*y^3*z^2';
   >> F = PolynString2CoefMat(f,{'x','y','z'})
   F =

        0     3.0000                     0          0     2.0000         
        0     1.0000                3.0000     3.0000          0         
        0          0                     0     2.0000     5.0000         
   8.5000     3.0000 - 2.0000i     -2.0000     6.0000     5.0000     

The following are a list of expert functions as software building blocks

Matrix computation functions

Relevant publications

               
               NAClab:  A Matlab toolbox for numerical algebraic computation,  Z. Zeng and T.Y, Li,  ACM Comm. in Computer Algebra, (47) 2013, pp. 170-173
               Intuitive interface for solving linear and nonlinear system of equations,  Z. Zeng,  Proc. of  ICMS 2018, 2018
               Computing multiple roots of inexact polynomials,  Z. Zeng,   Mathematics of Computation, 74(2005),  pp 869 - 903
                               Algorithm 835: MultRoot -- A Matlab package for computing polynomial roots and multiplicities, Z. Zeng,  ACM Transaction on Mathematical Software, 30(2004), pp 218-315
                    Multiple zeros of nonlinear systems  ,   B.H. Dayton, T.-Y. Li and Z. Zeng,  Mathematics of Computation, (80)2011, pp. 2143-2168
                    The numerical factorization of polynomials   W. Wu and Z. Zeng, Foundations of Computational Mathmatics, DOI 10.1007/s10208-015-9289-1          
                     Sensitivity and computation of a defective eigenvalue    Z. Zeng,  preprint,   [Resources / Matlab codes]  
                     Regularization and matrix computation in numerical polynomial algebra   Z. Zeng,
                Mixed volume  computation. A revisit   T.-L. Lee  and T..-Y. Li
                     HOM4PS-2.0:  A software package for solving polynomial systems by the polyhedral homotopy continuation method   T.-L. Lee,  T..-Y. Li and C.-H. Tsai  
                    Algorithm 931: An algorithm and software for computing multiplicity structure at zeros of nonlinear systems,  W. Hao, A.J. Sommese and Z. Zeng,  ACM TOMS
                    ApaTools:  A software toolbox for approximate polynomial algebra,  Zhonggang Zeng,  in Software for Algebraic Geometry,  IMA Volume 148, M.S. Stillman et al. eds.,  Springer,  pp149-167,  2008

Download and installation:  
    1. Download the zip file   (If your Matlab is older than 2013 edition, try this file )
    2. Unzip and extract the .p and .m  files in a folder, such as c:/NAClab
    3. Open Matlab and set path to the folder, e.g.  
              >> path(path,'c:/NAClab')

NAClab is then ready to go.


Modules with brief explanations and examples

LinearSolve

   Solve for the general numerical solution of
      (i)  matrix-vector equation   A*z = b  for z = z0+Kernel(A)
      (ii) linear operator equation    L(z) = b  for z = z0+Kernel(L)
   within an error tolerance tol, where z0 is the minimum norm solution
   and the kernel is the numerical kernle within the error tolerance tol
 
     Syntax for solving A*z = b (within error tolerance tol) 
     >> [z, K, cnd, res] = LinearSolve(A, b)
     >> [z, K, cnd, res] = LinearSolve(A, b, tol) 
 
     syntax for solving L(z) = b (within error tolerance tol)
     >> [z, K, cnd, res] = LinearSolve({L, Domain, para}, b)
     >> [z, K, cnd, res] = LinearSolve({L, Domain, para}, b, tol)
     >> [z, K, cnd, res] = LinearSolve({L, Domain, para, r}, b)
     >> [z, K, cnd, res] = LinearSolve({L, Domain, para, r}, b, tol)
 
     Input (case i: Matrix-vector equation A*z = b):
             A, b  -- coefficient matrix and right-side vector in A*z=b
              tol  -- error tolerance
 
     Input (case ii: Linear operator equation L(z) = b): 
              L  -- Matlab function handle for calculating the linear
                     transformation with syntax >> L(X1,...,Xm,Q1,...,Qk)
                     where X1,...,Xm and Q1,...,Qk are variables and
                     parameters of L respectively. Each one of X1,...,Xm is
                        (i) a matrix (including vector), or
                       (ii) a polynomial
                            to be transformed by L.
                     The remaining parameters Q1,...,Qk must be
                         provided as cell entries for para.
                     If the codomain is a product space, the output of the
                         function L must be a cell array
         Domain  -- domain of the linear transformation L represented by
                       (i) an mxn matrix of 0' and 1's representing
                           0 entries and variable entries respectively, or
                      (ii) a polynomial with variable coefficients being
                             nonzero, or
                     (iii) a cell array of matrices in (i) and/or (ii)
                     assuming the domain is a product of matrix spaces and
                     polynomial spaces
           para  -- (cell) parameters needed for running L
                     **** L must run with >> L(Domain{:},para{:}) ***
              r  -- (integer, optional) if provided, the rank-r TSVD
                     solution and kernel is computed (In other words,
                     substitute L and b with rank-r projections)
              b  -- (matrix, polynomial or cell) The right-hand side of
                         the equation L(z) = b, must be in the codomain
                         of L
             tol -- (numeric) error tolerance (< 1),
                       or r = tol > 1 indicating the rank-r projection of
                           the linear transformation is used
 
    Output:   z -- the minimum norm solution so that ||L(z)-b|| is minimized
              K -- an orthonormal basis for the kernel {y | L(y) = 0}
            cnd -- condition number of the representation matrix.
            res -- (vector) the residuals of ||L(Z)-b||,||L(K{1}||,...||L(K{n}||