DISCLAIMER - LICENCING AND SUPPORT POLICY This documentation is dated as of May 25, 1991. ** PLEASE DO NOT SEND MONEY FOR THIS PACKAGE ** This is NOT shareware! ANUM is good old fashioned Public Domain Software. Source code, totally written with Turbo C ++ rel. 1.0 and Turbo C 2.0, is not included with this archive for the reasons developped in the next section. Turbo C 2.0 version is again delivered back with release 4 of ANUM. I had to remove TC 2.0 from my previous hard disk to gain space. But as I have also changed of disk since them (from 40 Mbyte to 334 Mbyte) I am able to support both versions. I have not forgotten the (not so far ago) time when I was a student not earning a single coin, and for whom the buying of a compiler meant months of restrictions. The TC 2.0 version is planned to be supported until the end of the year. While the provided functions have been passed to a minimal test, I will have no liability or responsibilty to the user or any other person or entity with respect to any liability, loss or damage caused, or alleged to be caused directly or indirectly by this library, including, but not limited to, any interruption of service, loss of business, anticipatory profits or consequential damages resulting from the use of this program. Furthermore, although this program has been placed into the public domain, I retain all copyrights to this library around the world. Pursuant to this, this library MAY NOT BE SOLD BY ANY PARTY unless specifically authorized by the author, Jean- Marie de Montarby, in writing prior to the first copy being sold. As well, this program MAY NOT BE INCLUDED IN OR USED TO PRODUCE ANY OTHER PACKAGE FOR SALE even if this program is indicated as being 'in the public domain'. All of the above applies to both the original as well as derived or modified copies of the original. Any modified copies of this library MUST NOT have this documentation -including its copyright notice- violated, changed or altered. Occasionally, if you want to use ANUM, or a port of it to some other computers, in a commercial product, take contact with me first. We will perhaps find a deal. Thus this library, as delivered here, may only be used for personal and/or research and/or educational purposes. Comments, bug reports and improvement suggestions are welcome and should be adressed only to me. I can be reached at the following adresses: Electronic mail: JMM…FRESTP11.BITNET (until... ?) Normal mail : M. Jean-Marie de Montarby 40 avenue des Tilleuls 75016 PARIS FRANCE Do not try to phone me at home. I am not registered in the phone listings and due to the possible jet lags I am not much happy when waked up in the middle of the night. However during working hours in France (UTC+1, UTC+2), I may be reached at my employer's: DATAID SystŠmes Techniques et Industriels Voice : +33 (1) 40.99.19.19. Better use FAX, however : +33 (1) 40.99.06.85. If you wish to do so, you may send me in case of a presumed bug, a copy of your code. However I cannot guarantee neither a quick response nor the return of your media if a postal shipment. Anum is a yet (see next section) non-profitable piece of software and I prefer to work harder on a paid job rather than on Anum. As I work since Nov 5 1990, besides my current very partial time job as teacher in the E.S.T.P, I cannot foresee yet any availability delays for the next release of Anum. GENERAL OVERVIEW I have since July 1990 the diploma of "Ing‚nieur G‚omŠtre" from the "Ecole Sp‚ciale des Travaux Publics, du Bƒtiment et de l'Industrie" (E.S.T.P.). I have planned to release in three years a major piece of software for micro computers which will be a commercial product. Do'nt ask me neither what it is nor why this software will be, according to me, somewhat revolutionnary. It could be compared to the introduction of Turbo Pascal on the micro computers programming languages market, either by its planned efficiency and financial impact. This major piece of software requires several different fields of the so-called "art of programming". One of its aspect is that it will extensively use mathematics and scientific number crunching algorithms. I am currently developping the corresponding functions on a IBM PC compatible under MS-DOS which allows low-cost developpement at home. For the reasons explained at the end of the disclaimer section I cannot spend much time thoroughly testing those functions. This is why I decided to make public the first functions already done. More functions should come in the future as this mathematical subroutines package is not yet complete. ANUM is a library for use with the Turbo C++ 1.0 (and Turbo C 2.0) compilers from Borland International. It provides several routines for use by any programmer interested in Numerical Analysis. These routines are thoroughly explained further in one of the following sections of this documentation. A few of these routines are only of interest for Anum's own routines but have been left as public as their use is shared by several functions in several modules. However they will remain marked as reserved in this file. Usage of functions with similar names can be of some problem, especially if the "warn duplicate symbol" toggle is not set active at link edit time. To use the ANUM.LIB file, you should either. - specify its name in one of the line of your project file. - append it at the end of the command line if you use TCC. Do not forget then to include the MATHL.LIB file and the coprocessor chip or emulator library file depending on your system configuration. Please refer to Turbo C manuals for more information on that topic. ANUM code is compiled with word data alignement and adefault signed char type. Consequently you should ensure your compiler settings will correspond to this. Otherwise the result could be undefined. ANUM code is compiled without the VROOM overlay facilities (for the Turbo C++ 1.0 version of course). If you want a special copy of ANUM with them, send me a notice. The ANUM.LIB file was created using Borland International library utility TLIB. Its version number is 3.0 for TC++ version and 2.0 for TC 2.0 version. It is highly recommended ANSI prototyping be used in your programs making use of Anum. This is why, two header files are also provided to include in your own routines. The SYSIO.H file contains declarations for the error and warning messages interface, whilst the ANUM.H file contains those relevant to numerical codification of errors and warnings and of the other routines. A special effort has been made to keep the interface for parameters passings quite close to the one used by Borland International Turbo Pascal Numerical Toolbox to help those of you who would like to port their already existing applications from Pascal to C. In fact, most of the high level functions of Anum use the same algorithm than the Numerical Toolbox. This is due to the fact that they offer quite a wide coverage of usual numerical analysis within quite a few routines, and it is quite important when one starts writing its own numerical toolbox. CONVENTIONS The following conventions are applied: All of the functions present in the library, which are documented here, follow the C parameter calling conventions and have been provided for use with the LARGE memory model exclusively. The only near functions provided here are of the static class and then totally hidden to the application programmer. No function returns any pointers allocated. Some functions may allocate space for temporary items. But this is hidden to the application programmer and the space is freed at the end of the routine. It is hereby up to the programmer to allocate the necessary amount of space by using the farcalloc function before calling any function of the present library. Rows and columns are numbered from 1 unless explicitely stated. Polynomial indexation strictly follows the following rule: The degree is always increasing i.e. *poly points to poly[0], *(poly+1) points to poly[1] and so on. Matrices are indexed as follows: the first index is the row index, the second is the column one just as the classical way C handles arrays. The typedef'ed COMPLEX type is required to use some of the provided routines dealing with complex arithmetic. It is defined in the ANUM.H header file. All real values in the Anum library are of the type double. All Anum error codes take negative values. A zero value indicates there is no error. PROVIDED ROUTINES void backwardsubst() Reserved. void beep (void) This function outputs an audible tone to the speaker. Parameters: none. Return code: none Possible error codes: none void bisect(double lb, double ub, double tol, int maxiter, double *root, double *value, int *iter, int *errcode, double f()) This function uses a dichotomic algorithm to find one root of the equation f(x)=0 on a given interval. Parameters: lb : lower bound of the search interval ub : upper bound of the search interval root : value of the root, if found. value : = f(root) iter : number of iterations taken the algorithm actually went through. errcode : possible error code. f : function of the equation. Returns: nothing Possible error codes: EYSAMESIGN ETOLLE0 EMAXITERLT0 EOVERMAXITER void cadd(COMPLEX z1, COMPLEX z2, COMPLEX *z3) This function performs the addition of two complex numbers. Parameters: z1 : First operand. z2 : Second operand. z3 : = z1+z2 Returns: nothing Possible error codes: none double carg(COMPLEX z) This function returns the characteristic angle of a complex number. Parameters: z : complex number Returns : the value of the angle. Possible error codes: none Notes : carg(0) = PI/2 (convention) void cdiv(COMPLEX z1, COMPLEX z2, COMPLEX *z3, int *errcode) This functions performs complex division. Parameters: z1 : complex number to be divided. z2 : complex number to divide by. z3 : = z1/z3. errcode : possible error code. Returns: nothing Possible error codes: E0DIVIDE checkslope() Reserved function void cmult(COMPLEX z1, COMPLEX z2, COMPLEX *z3) This function performs complex multiplication. Parameters: z1 : First complex operand. z2 : Second complex operand. z3 : = z2 x z3 Returns: nothing Possible error codes: none void conjugue(COMPLEX z1, COMPLEX *z2) This function calculates a-ib if z=a+ib Parameters: z1 : a + ib z2 : a - ib Returns: nothing Possible error codes: none void cpoly(int degree, COMPLEX *poly, COMPLEX z, COMPLEX *y, *yderiv, COMPLEX *yderiv2, int *errcode) This function computes the values of the complex plynomial poly, and of its first and second derivatives at a given point. Parameters: degree : degree of the polynomial poly : pointer to the coefficients of the polynomial z : value of the point at which omputation must be performed y : = poly(z) yderiv : = poly'(z) yderiv2 : = poly"(z) errcode : possible error code. Returns: nothing Possible error codes: ENOWCORE void cpoly2_solve(COMPLEX *poly, COMPLEX *root1, COMPLEX *root2) This function solves a second degree complex polynomial equation. Parameters: poly : complex second degree polynome. root1 : first root. root2 : second root. Returns: nothing Possible error codes: none void cpoly_one_root(int degree, COMPLEX *poly, COMPLEX z0, double tol, int maxiter, COMPLEX *zroot, COMPLEX *yroot, int *iter, int *errcode) This function approximates a root of a complex polynomial. Parameters: degree : degree of the polynomial poly : pointer to the coefficients of the polynomial z0 : value of the root guess tol : tolerance used for convergence test. maxiter : maximum number of allowed iterations zroot : value of the root if found yroot : = poly(zroot) iter : actual number of iterations the algorith went through. errcode : possible error code. Returns: nothing Possible error codes: ENOWCORE E0DIVIDE EOVERMAXITER void cpoly_reduce(int degree, COMPLEX *poly, COMPLEX zroot, int *errcode) This function reduces the degree of a complex polynomial by factoring out one of its given root. Parameters: degree : degree of the polynomial poly : pointer to the coefficients of the polynomial zroot : value of the root errcode : possible error code. Returns: nothing Possible error codes: ENOWCORE void cpow(COMPLEX z1, int n, COMPLEX *z2) This function performs a complex power evelation. Parameters: z1 : complex to elevate n : value of the power z2 : z1**n Returns: nothing Possible error codes: none void csqrt(COMPLEX z1, COMPLEX *z2) This function calculates a complex square root. Parameters: z1 : complex operand z2 : = z1^(1/2) Returns: nothing Possible error codes: none void csub(COMPLEX z1, COMPLEX z2, COMPLEX *z3) This function performs complex substraction. Parameters: z1 : First complex operand. z2 : Second complex operand. z3 : = z1 - z3. Returns: nothing Possible error codes: none ctestroot Reserved function. void determinant(int dim, double *mat0,double *det, int *errcode) This function calculates the determinant of a square matrix. Parameters: dim : order of the matrix. mat0 : matrix the determinant of which will be calculated. det : value of the calculated determinant. errcode : error code. Returns: nothing Possible error codes: EDIMLE0 ENOWCORE void disa_err(void) This function disables the generation of error messages. Parameters: none Returns: nothing Possible error codes: none void disa_m_beep(void) This function disables the generation of an audible tone while printing error or warning messages. Parameters: none Returns: nothing Possible error codes: none void disa_wrn(void) This function disables the generation of warning messages. Parameters: none Returns: nothing Possible error codes: none void eigen_power(int dim, double *mat, double *vect0, int maxiter, double tol, double *eigenval, double *eigenvect, int *iter, int *errcode) This function computes the dominant eigenvalue of a square matrix. Parameters: dim : Order of the matrix mat : Matrix vect0 : Approximate value of the eigenvector. maxiter : Maximum number of iterations to perform. tol : tolerance used for convergence test. eigenval : Computed eigenvalue. eigenvect : Computed eigenvector. iter : Actual number of iterations the algorith went through. errcode : Error code. Returns: nothing Possible error codes: EDIMLE0 ETOLLE0 EMAXITERLT0 E0DIVIDE ENOWCORE Note : This algoritm will not converge if vect0 happens to be orthogonal to the actual result. void eigen_vect(int dim, double *mat, int nbeval, double *eval, double *evect, int *errcode) This function computes the eigenvectors of a square matrix, the corresponding eigen values of which are known. Parameters: dim : Order of the matrix mat : Matrix nbeval : Number of given eigen values. eval : Known eigen values evect : Computed eigenvectors. errcode : Error code. Returns: nothing Possible error codes: EDIMLE0 E0DIVIDE ENOWCORE Note : This algoritm will not converge if vect0 happens to be orthogonal to the actual result. void ena_err(void) This function enables the generation of error messages. Parameters: none Returns: nothing Possible error codes: none void ena_m_beep(void) This function enables the generation of an audible tone while printing error or warning messages. Parameters: none Returns: nothing Possible error codes: none void ena_wrn(void) This function enables the generation of warning messages. Parameters: none Returns: nothing Possible error codes: none void errmsg(char *orgid, char *date, char *time, int code, char *msg FILE *fp) This function outputs user defined error messages. Parameters: orgid : origin identifier (usually the caller __FILE__ pseudo variable) date : date (usually the caller __DATE__ pseudo variable) time : time (usually the caller __TIME__ pseudo variable) code : error code msg : explanatory message to be printed fp : file strea, descriptor on which output must be performed Returns: nothing Possible error codes: none void gauss_elimination(int dim, double *mat_A, double *vect_B, double *vect_X int *errcode) This routine solves a linear system (A*X=B) using gaussian elimination. Parameters dim : order of the matrix. mat_A: matrix A vect_B : second member vector B vect_X : solution vector X. errcode : error code. Possible error codes: EDIMLE0 ENOSOL EMATSING void gauss_seidel(int dim, double *mat_A, double *vect_B, double tol, int maxiter, double *vect_X, int *iter, int *errcode) This routine solves a linear system (A*X=B) using the Gauss Seidel iterative algorithm. Parameters dim : order of the matrix to invert. mat_A : matrix A vect_B : second member vector B tol : tolerance used for convergence test maxiter : maximum number of allowed iterations vect_X : solution vector X. iter : actual number of performed iterations errcode : error code. Possible error codes: EDIMLE0 ENOWCORE EOVERMAXITER WNODIAGDOM ETOLLE0 EDIMLE0 EMAXITERLE0 E0MATDIAG ENOSOL EMATSING void householder_givens(int dim, double *mat, int nbeval, double *eval, int *errcode) This routine computes a given number of eigen values of a symmetric matrix using the the Householder Givens algorithm. (It is an adaptation of the version written by J. Ortega in the 1960's for the IBM 7094). Parameters dim : order of the matrix. mat : matrix nbeval : number of eigen values to be computed. eval : computed eigen values errcode : error code. Possible error codes: EDIMLE0 ENOWCORE EMATNOTSYM Note : If nbeval is lower than or equal to 0 then no computation is performed. If it is greater than dim then all eigenvalues are computed. void inverse(int dim, double *mat0, double *invmat, int *errcode) This routine inverts a square matrix Parameters dim : order of the matrix to invert. mat0 : matrix to invert invmat : inverted matrix. errcode : error code. Possible error codes: EDIMLE0 EMATSING ENOWCORE isnear0() Reserved function. isinfinite() Reserved function. bool ismatsym(int dim, double *mat, int *errcode) This function determinates if a matrix is symmetric or not. Parameters : dim : Order of the matrix. mat : matrix to examine. errcode : error code Returns : TRUE or FALSE Possible error codes : EDIMLE0 void jacobi(int dim, double *mat, int maxiter, double tol, double *eval, double *evect, int *iter, int *errcode) This function computes the eigen vectors and the eigen values of a symmetric matrix using the Jacobi's rotations algorithm. Parameters: dim : Order of the matrix. mat : Matrix. maxiter : maximum number of allowed iterations. tol : tolerance used in convergence tests. eval : computed eigen values. evect : computed eigen vectors. errcode : error code. Returns : nothing Possible error codes: EDIMLE0 ETOLLE0 EMAXITERLT0 EOVERMAXITER EMATNOTSYM void laguerre (int degree, COMPLEX *poly, COMPLEX z0 double tol, int maxiter, int *nbroot, COMPLEX *roots, COMPLEX *yroots, int *errcode) This function solves a complex polynome. Parameters: degree : degree of the polynome. poly : polynome to solve. z0 : initial guess of an actual root. tol : tolerance to use. maxiter : maximum number of allowed iterations. nbroot : number of found roots. roots : array containing the found roots. yroots : array containing the value of the polynome at the found roots. errcode : error code. Returns : nothing Possible error codes: ENOWCORE EDEGLT2 EMAXITERLT0 EOVERMAXITER void le_verrier(int dim, double *mat, double *poly, int *errcode) This function compute the characteristic polynome of a square matrix. Parameters: dim : Order of the matrix. mat : matrix poly : characteristic polynomial. errcode : error code. Returns : nothing Possible error codes: ENOWCORE void lsq(int nbpoints, double *xdata, double *ydata, int *nbtermes, double *vectsol, double *yfit, double *residus, double *stddev, double *variance, int *errcode, lsqfit fit) This function computes least square fit of points. Parameters: nbpoints : number of data points. xdata : x of each point. ydata : y of each point. nbterms : number of terms in the least square fit. vectsol : computed LS coefficients. yfit : values of each point. stddev : standard deviation. variance : variance. errcode : error code. lsq : least square fit type Returns : nothing Possible error codes: ENOWCORE E0DIVIDE ELOGLE0 ENDPTSLT2 ENTRMLT1 ENDTPTLTLTRMS ENOSOL char *lsqname(lsqfit fit) This function returns the name associated to a least square method. Parameters: lsq : least square fit type Returns : the name of the least square type fit void lu_decompose(int dim, double *mat_A, double *decomp, double *matperm, int *errcode) This function performs LU decomposition in order to solve several linear systems having the same matrix (AxX=B, AxX'=B', AxX"=B", ...) with the lu_solve function. Parameters: dim : order of the matrix. mat_A : matrix A. decomp, matperm : matrices to be passed as arguments to lu_solve.. errcode : error code. Possible error codes: EDIMLE0 ENOWCORE EMATSING void lu_solve(int dim, double *decomp, double *vect_B, double *matperm, double *vect_X, int *errcode) This function solves a linear system (AxX=B) the matrix of which has been previously processed with lu_decompose. Parameters: dim : order of the matrix A. decomp : matrix resulting from the call to lu_decompose. vect_B : Second member vector B. matperm : matrix resulting from the call to lu_decompose. vect_X : solution vector X. errcode : error code. Possible error codes: EDIMLE0 ENOWCORE double modulus(COMPLEX z) This function calculates the modulus of a complex number. Parameters: z : complex operand. Returns: value of the modulus Possible error codes: none void muller(COMPLEX z0, double tol, int maxiter, COMPLEX *root, COMPLEX *yroot, int *iter, int *errcode,COMPLEX f(COMPLEX z)) This function solves a general complex equation, using the Mller's algorithm. Parameters: z0 : initial approximate value of the root. tol : tolerance used for convergence tests. maxiter : maximum number of allowed iterations. root : value of the complex root, if found. yroot : = f(root). errcode : error code. f : function to solve. Returns: nothing Possible error codes: EOVERMAXITER ETOLLE0 EMAXITERLT0 ENOINTERSECT ENOFITPARAB void mult_mat_vect(int dim, double *mat, double *vect, double *res, int *errcode) This function multiplies a square matrix by vector. Parameters: dim : Order of the matrix (= vector dimension) mat : Matrix. vect : Vector. res : Resulting vector. Returns: nothing Possible error codes: EDIMLE0 void newton_horner(int degree0, double *poly0, double x0, double tol, int maxiter, int *degree, int *nbroot,s double *poly, double *roots, double *image, double *value, double *deriv, int *iter, int *errcode) This function solves a real polynomial equation using the Newton-Horner algoritm. Parameters: degree0 : Degree of the polynomial equation. poly0 : Polynome defining the equation. x0 : Approximate value of the solution. tol : Tolerance used for convergence tests. maxiter : Maximum number of allowed iterations. degree : Degree of the resulting polynome after treatment. nbroots : Number of found roots poly : Resulting deflated polynome. roots : Array containing the values of the real parts of the roots found so far. image : Array containing the values of the imaginary parts of the roots found so far. value : Array containing the value of poly0 at the found roots. deriv : Array containing the value of the first derivate of poly0 at the found roots. iter : Array containing the number of iterations the algorithm actually went through for each found root. errcode : Error code. Returns: nothing Possible error codes: ETOLLE0 EMAXITERLT0 EDEGLE0 EOVERMAXITER ENULLSLOPE void newton_raphson(double x0, double tol, int maxiter, double *root, double *yroot, double *deriv, int *iter, int *errcode, double f(double t), double fprime(double t)) This routine solves a real equation f(x)=0 using the general ewton algorithm (more often called Newton-Raphson). Parameters: x0 : Approximate value of the solution. tol : Tolerance used for convergence tests. maxiter : Maximum number of allowed iterations. roots : Value of the root, if found. yroot : = f(root) deriv : = fprime(root). iter : Number of iterations the algorithm actually went through. errcode : Error code. f : function defining the equation. fprime : first derivate of f. Returns: nothing Possible error codes: ETOLLE0 EMAXITERLT0 EOVERMAXITER ENULLSLOPE void partial_pivot(int dim, double *mat_A, double *vect_B, double *vect_X, int *errcode) This function solves a linear system (AxX=B) using the partial pivoting algorithm. Parameters: dim : order of matrix A. vect_B : second member vector B. vect_X : solution vector X. errcode : Error code. Returns: nothing Possible error codes: EDIMLE0 EMATSING ENOWCORE void rassign(double lambda, COMPLEX *z) This function assigns a real value to a complex number. Parameters: lambda : value to assign. z : = labda + 0i. Returns: nothing Possible error codes: none void rdiv_vect(int dim, double *vect, double lambda, int *errcode) This function allows to divide all the coordinates composing a vector by a real constant Parameters: none. Returns: nothing Possible error codes: none void rmult(double lambda, COMPLEX *z) This functions multiplies a complex by a real value. Parameters: lambda : real multiplier. z : complex which will be multiplied by lambda. Returns: nothing Possible error codes: none void romberg(double lb, double ub, double tol, int maxiter, double *I, int *iter, int *errcode, double f(double)) This function computes the integration of a real function on an interval using the Romberg algorithm. Parameters: lb : lower bound of the interval. ub : upper bound of the interval. tol : tolerance used for convergence test. maxiter : maximum number of allowed iterations. I : value of the computed integral. iter : actual number of iterations the algorithm went through. errcode : Error code. f : function to integrate Returns: nothing Possible error codes: ENOWCORE EOVERMAXITER EMAXITERLT0 ETOLLE0 void row_mat_div(int dim, double *mat, double lambda, int nrow, int *errcode) This function divides a given row of a square matrix by a given constant. Parameters: dim : Order of the matrix. mat : Points to the matrix. lambda : The real constant. nrow : Row index (1 to dim). errcode : Error code. Returns: nothing Possible error codes: E0DIVIDE rowmultadd() Reserved function. void rpoly2_solve(double a, double b, double c, double *root1r, double *root1im, double *root2r, double *root2im) This function solves the real polynomial equation ax^2+bx+c=0. Parameters: a, b, c : coefficient of the polynome. root1r : Real part of the first root. root1im : Imaginary part of the first root. root2r : Real part of the second root. root2im : Imaginary part of the second root. Returns: nothing Possible error codes: none void rpoly_one_root(int degree double *poly, double x, double tol, int maxiter, double *root, double *value, double *deriv, int *iter, int *errcode) This function finds one real root of a given real polynome. Parameters: degree : Degree of the polynome. poly : Polnome, one root of which must be found. x : Initial approximation of the root. tol : Tolerance used for convergence tests. maxiter : Maximum number of iterations allowed. root : Found root, if any. value : = poly(root) deriv : value of the first derivate of poly for x = root. iter : Actual number of iteration the algorithm went through. errcode : Error Code. Returns: nothing Possible error codes: EDEGLE0 EOVERMAXITER ENULLSLOPE void rpoly_reduce(int *degree, double *poly, double root, int *errcode) This functions reduce the degree of a given polynome, a root of which is known. Parameters: degree : Degree of the polynome. poly : Polynome. root : Known root. errcode : Error code. Returns: nothing Possible error codes: EDEGLE0 ENOWCORE void rpoly_value_deriv(int degree, double *poly, double x, double *y, double *deriv) This routines computes the value of a polynome and the one of its first derivate at a given point. Parameters: degree : Degree of the polynome. poly : Polynome. x : Real value. y : = poly(x). deriv : = poly'(x). Returns: nothing Possible error codes: none rtestroot() Reserved function. void secant(double x01, double x02, double tol, int maxiter, double *root, double *y, int *iter, int *errcode, double f(double t)) This function solves a real equation f(x)=0 using the secant algorithm. Parameters: x01 , x02 : Approximate values of the root. tol : Tolerance used for convergence tests. maxiter : Maximum number of allowed iterations. root : Found root, if any. iter : Actual number of iterations the algorithm went through. y : = f(root). errcode: Error code f : function defining the equation. Returns: nothing Possible error codes: EYSAMESIGN EMAXITERLT0 EOVERMAXITER ETOLLE0 void simpson(double lb, double ub, int nbinter, double *I, int *errcode, double f(double)) This function computes the integration of a real function on an interval using the Simpson method. Parameters: lb : lower bound of the interval. ub : upper bound of the interval. nbinter : number of sub-intervals used for the computation. I : value of the computed integral. errcode : Error code. f : function to integrate Returns: nothing Possible error codes: ENUMINTLE0 void slargest_vect(int dim, double *vect, double *res) This function finds the largest element in a vector. Parameters: dim : dimension of the vector. vect : vector. res : result. Returns: nothing void steffensen(double x0, double tol, int maxiter, double *root, double *value, int *iter, int *errcode, double f(double t)) This function finds the root of f(x)=0 using the steffensen algorithm. Parameters: x0 : Approximate value of the root. tol : tolerance to be used as stp criteria. root : value of the found root. value : =f(root). iter : number of iterations the algorith actually went through. f : function to solve. Possible error codes: E0DIVIDE EOVERMAXITER void swap_rows(int dim, double *mat, int nrow1, int nrow2, int *errcode) This routines exchanges two rows of a square matrix. Parameters: none. dim : Order of the matrix. mat : Matrix to be treated. nrow1, nrow2 : indexes of the given rows (1 to n) errcode : Error code. Returns: nothing Possible error codes: ENOWCORE void trans_size_system(int dim1, double *m1, int dim2, double *m2) This function converts a matrix of a given order to another one. Parameters: none. dim1 : order of the initial matrix. m1: initial Matrix to be treated. dim2 : order of the target matrix m2 : target matrix Note : if dim2>dim1 the rest of the m2 will be completed by zeroes, else the dim1-dim2th columns and rows of the matrix m1 will be lost in m2. void werrmsg(char *orgid, char *date, char *time, int code, FILE *fp); This function ouputs an error or warning message associated to a predefined code: Parameters: orgid : Message indicating the posistion of the call to werrmsg.(usually __FILE__ or __LINE__ or a combination of the two). date : Date of the last release of the file (usually __DATE__). time : time of the last release of the file (usually __TIME__). code : Error code returned by the last Anum called function. fp : Stream, the output should be performed to which. Returns: nothing Possible error codes: none Note : if code is zero, no message is output. You can then insert a call to werrmsg after any function call to anum returning an error code at debugging stage. Note : The SYSMSG macro defined in "sysio.h" helps programmers to automate calls to standard ANUM error and warning messages. void wrnmsg(char *orgid, char *date, char *time, int code, char *msg FILE *fp) This function outputs user defined warning messages. Parameters: orgid : origin identifier (usually the caller __FILE__ pseudo variable) date : date (usually the caller __DATE__ pseudo variable) time : time (usually the caller __TIME__ pseudo variable) code : error code msg : explanatory message to be printed fp : file strea, descriptor on which output must be performed Returns: nothing Possible error codes: none ERROR AND WARNING CODES Here are a few of the possible returnes error codes. A more complete list is given in the "anum.h" file. I promise I will include it in the next release of ANUM documentation :-). It is no joke working until 2 a.m. on a non paid job ;-). ETOLLE0 Class : Error Cause : The called routine was given a nil or negative tolerance value. User's Response : Correct the 'tol' parameter to a strictly positive value. EOVERMAXITER Class : Error Cause : An iterative function could not find the solution within 'maxiter' iterations. User's Response : Set the 'maxiter' parameter to a greater value or check the algorithm will converge with the given approximate values of the solution or assume there is no solution. ENULLSLOPE Class : Error Cause : The user function attempted to use a Newton algorithm with a point, the slope of the function is zero at which. User's Response : Choose another parameter value for approximate value of the root. WCPLXROOTS Class : Warning Cause : A function solving real polynomes could not find a number of real roots equal to its degree. This polynome could possibly have one or more complex roots. User's Response : Use a routine giving al so complex roots or set the 'maxiter' parameter to a greater value. EYSAMESIGN Class : Error Cause : The bisection algorithm was given two bounds the values of which are of the same sign. User's Response : Correct one of the bounds. Check if there is at least one root between the bounds or assume there is no solution. EMAXITERLT0 Class : Error Cause : An iterative routine was given a negative value for the maximum number of iterations as parameter. User's Response : Input a positive value as 'maxiter' parameter. A nil value will only be suitable if you want to treat only cases for which a trivial and direct solution exists. ENOPXINTERSECT Class : Error Cause : The computed parabola does not intersect the X axis. EDIMLT0 Class : Error Cause : A function attempted to use a negative length vector or a negative order square matrix. User's Response : Change the dim parameter or check your polynome is not totally deflated. EMATSING Class : Error Cause : A routine attempted to invert a singular matrix. EDIMLE0 Class : Error Cause : A function attempted to use a negative or zero length vector or a negative or zero order square matrix. User's Response : Change the dim parameter or check your polynome is not totally deflated. E0DIVIDE Class : Error Cause : A function attempted to do a division by zero. EDEGLT2 Class : Error Cause : The selected algorithm requires a polynome the degree of which is at least two. User's Response : Check your calling function for errors or try to use a direct solving method. ENOWCORE Class : Error Cause : There is not enough working core for one of the routine to allocate its working arrays. User's Response : Try to free some main storage by eliminating Terminate-and-Stay-Resident programs. Use overlays if your compilers allows them. Check your code for unused stuff. EDEGLE0 Class : Error Cause : A function attempted to use a negative or zero degree polynomial. User's Response : Check passed arguments or try to directly solve a trivial case (if degree is nil or polynomial equals 0). ENUMINTLE0 Class : Error Cause : The given number of arguments passed to a routine is lower than or equal to 0. HISTORY LOG - CHANGES TO COME From first to second release Fixed a bug concerning a hardwire coded error in the Newton- Raphson algorithm. Added : steffensen function laguerre function ena_beep function ena_warn function ena_err function ena_m_beep disa_m_beep The five last functions allow or disable sytem generated warning or messages and the audible tone associated to it. Not to cause some "vapourware" announces and as programming generates unexpected delays, no date is given for the next release of Anum. From second to third release Fixed severe bug in swap_rows function giving some random results with inverse and determinant functions. Added : a whole bunch of new functions. From third to fourth release Fixed severe bug in swap_rows function giving some random results with inverse and determinant functions. Added : (lsq, eigen_vector, householder_givens, jacobi, le_verrier). Added : C/C++ support, ANUM release 4 routines should link without problem with C++ modules. This is achieved by suppressing automatically type checking at link edit time (See the #ifdef compilation directives in ANUM.H and STDIO.H). In release 4.1, fixed a bug affecting the inverse function and the jacobi method due to a wrong parameter passing in row_mat_div. Solved an ambiguity in the eigen_vect routine. fixed a (very minor) bug in testlib.bat. In the next release of this package: Changes of the interface to support Turbo Pascal are considered. This would primarily affect the parameter calling conventions. Microsoft compilers are too slow, require too much space on the secondary storage, are too user unfriendly and too expansive for me to forecast any support yet. More functions are to come out, of course. Especially in the differential equations solving and interpolation fields. A better documentation will also arise, especially from your comments. Thanks to : CFGROP at WEIZMANN.BITNET ROSSBACH at DMRHRZ11.BITNET C3DOUTRE at FRECCL11.BITNET m_rives at paris.paris.ingr.com USERGB3W at UMICHUM.BITNET for support and/or suggestions. .D‚but Table M. DISCLAIMER - LICENCING AND SUPPORT POLICY 1 GENERAL OVERVIEW 3 CONVENTIONS 5 PROVIDED ROUTINES 6 void backwardsubst() 6 void beep (void) 7 void bisect(double lb, double ub, double tol, int maxiter, double *root, double *value, int *iter, int *errcode, double f()) 8 void cadd(COMPLEX z1, COMPLEX z2, COMPLEX *z3) 9 double carg(COMPLEX z) 10 void cdiv(COMPLEX z1, COMPLEX z2, COMPLEX *z3, int *errcode) 11 checkslope() 12 void cmult(COMPLEX z1, COMPLEX z2, COMPLEX *z3) 13 void conjugue(COMPLEX z1, COMPLEX *z2) 14 void cpoly(int degree, COMPLEX *poly, COMPLEX z, COMPLEX *y, *yderiv, COMPLEX *yderiv2, int *errcode) 15 void cpoly2_solve(COMPLEX *poly, COMPLEX *root1, COMPLEX *root2) 16 void cpoly_one_root(int degree, COMPLEX *poly, COMPLEX z0, double tol, int maxiter, COMPLEX *zroot, COMPLEX *yroot, int *iter, int *errcode) 17 void cpoly_reduce(int degree, COMPLEX *poly, COMPLEX zroot, int *errcode) 18 void cpow(COMPLEX z1, int n, COMPLEX *z2) 19 void csqrt(COMPLEX z1, COMPLEX *z2) 20 void csub(COMPLEX z1, COMPLEX z2, COMPLEX *z3) 21 ctestroot 22 void determinant(int dim, double *mat0,double *det, int *errcode) 23 void disa_err(void) 24 void disa_m_beep(void) 25 void disa_wrn(void) 26 void eigen_power(int dim, double *mat, double *vect0, int maxiter, double tol, double *eigenval, double *eigenvect, int *iter, int *errcode) 27 void eigen_vect(int dim, double *mat, int nbeval, double *eval, double *evect, int *errcode) 28 void ena_err(void) 29 void ena_m_beep(void) 30 void ena_wrn(void) 31 void errmsg(char *orgid, char *date, char *time, int code, char *msg FILE *fp) 32 void gauss_elimination(int dim, double *mat_A, double *vect_B, double *vect_X int *errcode) 33 void gauss_seidel(int dim, double *mat_A, double *vect_B, double tol, int maxiter, double *vect_X, int *iter, int *errcode) 34 void householder_givens(int dim, double *mat, int nbeval, double *eval, int *errcode) 35 void inverse(int dim, double *mat0, double *invmat, int *errcode) 36 isnear0() 37 isinfinite() 38 bool ismatsym(int dim, double *mat, int *errcode) 39 void jacobi(int dim, double *mat, int maxiter, double tol, double *eval, double *evect, int *iter, int *errcode) 40 void laguerre (int degree, COMPLEX *poly, COMPLEX z0 double tol, int maxiter, int *nbroot, COMPLEX *roots, COMPLEX *yroots, int *errcode) 41 void le_verrier(int dim, double *mat, double *poly, int *errcode) 42 void lsq(int nbpoints, double *xdata, double *ydata, int *nbtermes, double *vectsol, double *yfit, double *residus, double *stddev, double *variance, int *errcode, lsqfit fit) 43 char *lsqname(lsqfit fit) 44 void lu_decompose(int dim, double *mat_A, double *decomp, double *matperm, int *errcode) 45 void lu_solve(int dim, double *decomp, double *vect_B, double *matperm, double *vect_X, int *errcode) 46 double modulus(COMPLEX z) 47 void muller(COMPLEX z0, double tol, int maxiter, COMPLEX *root, COMPLEX *yroot, int *iter, int *errcode,COMPLEX f(COMPLEX z)) 48 void mult_mat_vect(int dim, double *mat, double *vect, double *res, int *errcode) 49 void newton_horner(int degree0, double *poly0, double x0, double tol, int maxiter, int *degree, int *nbroot,s double *poly, double *roots, double *image, double *value, double *deriv, int *iter, int *errcode) 50 void newton_raphson(double x0, double tol, int maxiter, double *root, double *yroot, double *deriv, int *iter, int *errcode, double f(double t), double fprime(double t)) 52 void partial_pivot(int dim, double *mat_A, double *vect_B, double *vect_X, int *errcode) 53 void rassign(double lambda, COMPLEX *z) 54 void rdiv_vect(int dim, double *vect, double lambda, int *errcode) 55 void rmult(double lambda, COMPLEX *z) 56 void romberg(double lb, double ub, double tol, int maxiter, double *I, int *iter, int *errcode, double f(double)) 57 void row_mat_div(int dim, double *mat, double lambda, int nrow, int *errcode) 58 rowmultadd() 59 void rpoly2_solve(double a, double b, double c, double *root1r, double *root1im, double *root2r, double *root2im) 60 void rpoly_one_root(int degree double *poly, double x, double tol, int maxiter, double *root, double *value, double *deriv, int *iter, int *errcode) 61 void rpoly_reduce(int *degree, double *poly, double root, int *errcode) 62 void rpoly_value_deriv(int degree, double *poly, double x, double *y, double *deriv) 63 rtestroot() 64 void secant(double x01, double x02, double tol, int maxiter, double *root, double *y, int *iter, int *errcode, double f(double t)) 65 void simpson(double lb, double ub, int nbinter, double *I, int *errcode, double f(double)) 66 void slargest_vect(int dim, double *vect, double *res) 67 void steffensen(double x0, double tol, int maxiter, double *root, double *value, int *iter, int *errcode, double f(double t)) 68 void swap_rows(int dim, double *mat, int nrow1, int nrow2, int *errcode) 69 void trans_size_system(int dim1, double *m1, int dim2, double *m2) 70 void werrmsg(char *orgid, char *date, char *time, int code, FILE *fp); 71 void wrnmsg(char *orgid, char *date, char *time, int code, char *msg FILE *fp) 72 ERROR AND WARNING CODES 73 ETOLLE0 73 EOVERMAXITER 73 ENULLSLOPE 73 WCPLXROOTS 73 EYSAMESIGN 74 EMAXITERLT0 74 ENOPXINTERSECT 74 EDIMLT0 74 EMATSING 74 EDIMLE0 74 E0DIVIDE 75 EDEGLT2 75 ENOWCORE 75 EDEGLE0 75 ENUMINTLE0 75 HISTORY LOG - CHANGES TO COME 76 .Fin Table M.