The package defines the following public procedures (several exist as specialised procedures, see below):
Constructingmatricesandvectors::math::linearalgebra::mkVectorndimvalue
Create a vector with ndim elements, each with the value value.
integer ndim
Dimension of the vector (number of components)
double value
Uniform value to be used (default: 0.0)
::math::linearalgebra::mkUnitVectorndimndir
Create a unit vector in ndim-dimensional space, along the ndir-th direction.
integer ndim
Dimension of the vector (number of components)
integer ndir
Direction (0, ..., ndim-1)
::math::linearalgebra::mkMatrixnrowsncolsvalue
Create a matrix with nrows rows and ncols columns. All elements have the value value.
integer nrows
Number of rows
integer ncols
Number of columns
double value
Uniform value to be used (default: 0.0)
::math::linearalgebra::getrowmatrixrow ?imin? ?imax?
Returns a single row of a matrix as a list
list matrix
Matrix in question
integer row
Index of the row to return
integer imin
Minimum index of the column (default: 0)
integer imax
Maximum index of the column (default: ncols-1)
::math::linearalgebra::setrowmatrixrownewvalues ?imin? ?imax?
Set a single row of a matrix to new values (this list must have the same number of elements as the
number of columns in the matrix)
list matrixname of the matrix in question
integer row
Index of the row to update
list newvalues
List of new values for the row
integer imin
Minimum index of the column (default: 0)
integer imax
Maximum index of the column (default: ncols-1)
::math::linearalgebra::getcolmatrixcol ?imin? ?imax?
Returns a single column of a matrix as a list
list matrix
Matrix in question
integer col
Index of the column to return
integer imin
Minimum index of the row (default: 0)
integer imax
Maximum index of the row (default: nrows-1)
::math::linearalgebra::setcolmatrixcolnewvalues ?imin? ?imax?
Set a single column of a matrix to new values (this list must have the same number of elements as
the number of rows in the matrix)
list matrixname of the matrix in question
integer col
Index of the column to update
list newvalues
List of new values for the column
integer imin
Minimum index of the row (default: 0)
integer imax
Maximum index of the row (default: nrows-1)
::math::linearalgebra::getelemmatrixrowcol
Returns a single element of a matrix/vector
list matrix
Matrix or vector in question
integer row
Row of the element
integer col
Column of the element (not present for vectors)
::math::linearalgebra::setelemmatrixrow ?col? newvalue
Set a single element of a matrix (or vector) to a new value
list matrixname of the matrix in question
integer row
Row of the element
integer col
Column of the element (not present for vectors)
::math::linearalgebra::swaprowsmatrixirow1irow2 ?imin? ?imax?
Swap two rows in a matrix completely or only a selected part
list matrixname of the matrix in question
integer irow1
Index of first row
integer irow2
Index of second row
integer imin
Minimum column index (default: 0)
integer imin
Maximum column index (default: ncols-1)
::math::linearalgebra::swapcolsmatrixicol1icol2 ?imin? ?imax?
Swap two columns in a matrix completely or only a selected part
list matrixname of the matrix in question
integer irow1
Index of first column
integer irow2
Index of second column
integer imin
Minimum row index (default: 0)
integer imin
Maximum row index (default: nrows-1)
Queryingmatricesandvectors::math::linearalgebra::showobj ?format? ?rowsep? ?colsep?
Return a string representing the vector or matrix, for easy printing. (There is currently no way
to print fixed sets of columns)
list obj
Matrix or vector in question
string format
Format for printing the numbers (default: %6.4f)
string rowsep
String to use for separating rows (default: newline)
string colsep
String to use for separating columns (default: space)
::math::linearalgebra::dimobj
Returns the number of dimensions for the object (either 0 for a scalar, 1 for a vector and 2 for a
matrix)
any obj
Scalar, vector, or matrix
::math::linearalgebra::shapeobj
Returns the number of elements in each dimension for the object (either an empty list for a
scalar, a single number for a vector and a list of the number of rows and columns for a matrix)
any obj
Scalar, vector, or matrix
::math::linearalgebra::conformingtypeobj1obj2
Checks if two objects (vector or matrix) have conforming shapes, that is if they can be applied in
an operation like addition or matrix multiplication.
string type
Type of check:
• "shape" - the two objects have the same shape (for all element-wise operations)
• "rows" - the two objects have the same number of rows (for use as A and b in a
system of linear equations Ax=b
• "matmul" - the first object has the same number of columns as the number of rows of
the second object. Useful for matrix-matrix or matrix-vector multiplication.
list obj1
First vector or matrix (left operand)
list obj2
Second vector or matrix (right operand)
::math::linearalgebra::symmetricmatrix ?eps?
Checks if the given (square) matrix is symmetric. The argument eps is the tolerance.
list matrix
Matrix to be inspected
float eps
Tolerance for determining approximate equality (defaults to 1.0e-8)
Basicoperations::math::linearalgebra::normvectortype
Returns the norm of the given vector. The type argument can be: 1, 2, inf or max, respectively the
sum of absolute values, the ordinary Euclidean norm or the max norm.
list vector
Vector, list of coefficients
string type
Type of norm (default: 2, the Euclidean norm)
::math::linearalgebra::norm_onevector
Returns the L1 norm of the given vector, the sum of absolute values
list vector
Vector, list of coefficients
::math::linearalgebra::norm_twovector
Returns the L2 norm of the given vector, the ordinary Euclidean norm
list vector
Vector, list of coefficients
::math::linearalgebra::norm_maxvector ?index?
Returns the Linf norm of the given vector, the maximum absolute coefficient
list vector
Vector, list of coefficients
integer index
(optional) if non zero, returns a list made of the maximum value and the index where that
maximum was found. if zero, returns the maximum value.
::math::linearalgebra::normMatrixmatrixtype
Returns the norm of the given matrix. The type argument can be: 1, 2, inf or max, respectively the
sum of absolute values, the ordinary Euclidean norm or the max norm.
list matrix
Matrix, list of row vectors
string type
Type of norm (default: 2, the Euclidean norm)
::math::linearalgebra::dotproductvect1vect2
Determine the inproduct or dot product of two vectors. These must have the same shape (number of
dimensions)
list vect1
First vector, list of coefficients
list vect2
Second vector, list of coefficients
::math::linearalgebra::unitLengthVectorvector
Return a vector in the same direction with length 1.
list vector
Vector to be normalized
::math::linearalgebra::normalizeStatmv
Normalize the matrix or vector in a statistical sense: the mean of the elements of the columns of
the result is zero and the standard deviation is 1.
list mv
Vector or matrix to be normalized in the above sense
::math::linearalgebra::axpyscalemv1mv2
Return a vector or matrix that results from a "daxpy" operation, that is: compute a*x+y (a a
scalar and x and y both vectors or matrices of the same shape) and return the result.
Specialised variants are: axpy_vect and axpy_mat (slightly faster, but no check on the arguments)
double scale
The scale factor for the first vector/matrix (a)
list mv1
First vector or matrix (x)
list mv2
Second vector or matrix (y)
::math::linearalgebra::addmv1mv2
Return a vector or matrix that is the sum of the two arguments (x+y)
Specialised variants are: add_vect and add_mat (slightly faster, but no check on the arguments)
list mv1
First vector or matrix (x)
list mv2
Second vector or matrix (y)
::math::linearalgebra::submv1mv2
Return a vector or matrix that is the difference of the two arguments (x-y)
Specialised variants are: sub_vect and sub_mat (slightly faster, but no check on the arguments)
list mv1
First vector or matrix (x)
list mv2
Second vector or matrix (y)
::math::linearalgebra::scalescalemv
Scale a vector or matrix and return the result, that is: compute a*x.
Specialised variants are: scale_vect and scale_mat (slightly faster, but no check on the
arguments)
double scale
The scale factor for the vector/matrix (a)
list mv
Vector or matrix (x)
::math::linearalgebra::rotatecsvect1vect2
Apply a planar rotation to two vectors and return the result as a list of two vectors: c*x-s*y and
s*x+c*y. In algorithms you can often easily determine the cosine and sine of the angle, so it is
more efficient to pass that information directly.
double c
The cosine of the angle
double s
The sine of the angle
list vect1
First vector (x)
list vect2
Seocnd vector (x)
::math::linearalgebra::transposematrix
Transpose a matrix
list matrix
Matrix to be transposed
::math::linearalgebra::matmulmv1mv2
Multiply a vector/matrix with another vector/matrix. The result is a matrix, if both x and y are
matrices or both are vectors, in which case the "outer product" is computed. If one is a vector
and the other is a matrix, then the result is a vector.
list mv1
First vector/matrix (x)
list mv2
Second vector/matrix (y)
::math::linearalgebra::anglevect1vect2
Compute the angle between two vectors (in radians)
list vect1
First vector
list vect2
Second vector
::math::linearalgebra::crossproductvect1vect2
Compute the cross product of two (three-dimensional) vectors
list vect1
First vector
list vect2
Second vector
::math::linearalgebra::matmulmv1mv2
Multiply a vector/matrix with another vector/matrix. The result is a matrix, if both x and y are
matrices or both are vectors, in which case the "outer product" is computed. If one is a vector
and the other is a matrix, then the result is a vector.
list mv1
First vector/matrix (x)
list mv2
Second vector/matrix (y)
Commonmatricesandtestmatrices::math::linearalgebra::mkIdentitysize
Create an identity matrix of dimension size.
integer size
Dimension of the matrix
::math::linearalgebra::mkDiagonaldiag
Create a diagonal matrix whose diagonal elements are the elements of the vector diag.
list diag
Vector whose elements are used for the diagonal
::math::linearalgebra::mkRandomsize
Create a square matrix whose elements are uniformly distributed random numbers between 0 and 1 of
dimension size.
integer size
Dimension of the matrix
::math::linearalgebra::mkTriangularsize ?uplo? ?value?
Create a triangular matrix with non-zero elements in the upper or lower part, depending on
argument uplo.
integer size
Dimension of the matrix
string uplo
Fill the upper (U) or lower part (L)
double value
Value to fill the matrix with
::math::linearalgebra::mkHilbertsize
Create a Hilbert matrix of dimension size. Hilbert matrices are very ill-conditioned with respect
to eigenvalue/eigenvector problems. Therefore they are good candidates for testing the accuracy of
algorithms and implementations.
integer size
Dimension of the matrix
::math::linearalgebra::mkDingdongsize
Create a "dingdong" matrix of dimension size. Dingdong matrices are imprecisely represented, but
have the property of being very stable in such algorithms as Gauss elimination.
integer size
Dimension of the matrix
::math::linearalgebra::mkOnessize
Create a square matrix of dimension size whose entries are all 1.
integer size
Dimension of the matrix
::math::linearalgebra::mkMolersize
Create a Moler matrix of size size. (Moler matrices have a very simple Choleski decomposition. It
has one small eigenvalue and it can easily upset elimination methods for systems of linear
equations.)
integer size
Dimension of the matrix
::math::linearalgebra::mkFranksize
Create a Frank matrix of size size. (Frank matrices are fairly well-behaved matrices)
integer size
Dimension of the matrix
::math::linearalgebra::mkBordersize
Create a bordered matrix of size size. (Bordered matrices have a very low rank and can upset
certain specialised algorithms.)
integer size
Dimension of the matrix
::math::linearalgebra::mkWilkinsonW+size
Create a Wilkinson W+ of size size. This kind of matrix has pairs of eigenvalues that are very
close together. Usually the order (size) is odd.
integer size
Dimension of the matrix
::math::linearalgebra::mkWilkinsonW-size
Create a Wilkinson W- of size size. This kind of matrix has pairs of eigenvalues with opposite
signs, when the order (size) is odd.
integer size
Dimension of the matrix
Commonalgorithms::math::linearalgebra::solveGaussmatrixbvect
Solve a system of linear equations (Ax=b) using Gauss elimination. Returns the solution (x) as a
vector or matrix of the same shape as bvect.
list matrix
Square matrix (matrix A)
list bvect
Vector or matrix whose columns are the individual b-vectors
::math::linearalgebra::solvePGaussmatrixbvect
Solve a system of linear equations (Ax=b) using Gauss elimination with partial pivoting. Returns
the solution (x) as a vector or matrix of the same shape as bvect.
list matrix
Square matrix (matrix A)
list bvect
Vector or matrix whose columns are the individual b-vectors
::math::linearalgebra::solveTriangularmatrixbvect ?uplo?
Solve a system of linear equations (Ax=b) by backward substitution. The matrix is supposed to be
upper-triangular.
list matrix
Lower or upper-triangular matrix (matrix A)
list bvect
Vector or matrix whose columns are the individual b-vectors
string uplo
Indicates whether the matrix is lower-triangular (L) or upper-triangular (U). Defaults to
"U".
::math::linearalgebra::solveGaussBandmatrixbvect
Solve a system of linear equations (Ax=b) using Gauss elimination, where the matrix is stored as a
band matrix (cf.STORAGE). Returns the solution (x) as a vector or matrix of the same shape as
bvect.
list matrix
Square matrix (matrix A; in band form)
list bvect
Vector or matrix whose columns are the individual b-vectors
::math::linearalgebra::solveTriangularBandmatrixbvect
Solve a system of linear equations (Ax=b) by backward substitution. The matrix is supposed to be
upper-triangular and stored in band form.
list matrix
Upper-triangular matrix (matrix A)
list bvect
Vector or matrix whose columns are the individual b-vectors
::math::linearalgebra::determineSVDAeps
Determines the Singular Value Decomposition of a matrix: A = U S Vtrans. Returns a list with the
matrix U, the vector of singular values S and the matrix V.
list A Matrix to be decomposed
float eps
Tolerance (defaults to 2.3e-16)
::math::linearalgebra::eigenvectorsSVDAeps
Determines the eigenvectors and eigenvalues of a real symmetric matrix, using SVD. Returns a list
with the matrix of normalized eigenvectors and their eigenvalues.
list A Matrix whose eigenvalues must be determined
float eps
Tolerance (defaults to 2.3e-16)
::math::linearalgebra::leastSquaresSVDAyqmineps
Determines the solution to a least-sqaures problem Ax ~ y via singular value decomposition. The
result is the vector x.
Note that if you add a column of 1s to the matrix, then this column will represent a constant like
in: y = a*x1 + b*x2 + c. To force the intercept to be zero, simply leave it out.
list A Matrix of independent variables
list y List of observed values
float qmin
Minimum singular value to be considered (defaults to 0.0)
float eps
Tolerance (defaults to 2.3e-16)
::math::linearalgebra::choleskimatrix
Determine the Choleski decomposition of a symmetric positive semidefinite matrix (this condition
is not checked!). The result is the lower-triangular matrix L such that L Lt = matrix.
list matrix
Matrix to be decomposed
::math::linearalgebra::orthonormalizeColumnsmatrix
Use the modified Gram-Schmidt method to orthogonalize and normalize the columns of the given
matrix and return the result.
list matrix
Matrix whose columns must be orthonormalized
::math::linearalgebra::orthonormalizeRowsmatrix
Use the modified Gram-Schmidt method to orthogonalize and normalize the rows of the given matrix
and return the result.
list matrix
Matrix whose rows must be orthonormalized
::math::linearalgebra::dgermatrixalphaxy ?scope?
Perform the rank 1 operation A + alpha*x*y' inline (that is: the matrix A is adjusted). For
convenience the new matrix is also returned as the result.
list matrix
Matrix whose rows must be adjusted
double alpha
Scale factor
list x A column vector
list y A column vector
list scope
If not provided, the operation is performed on all rows/columns of A if provided, it is
expected to be the list {imin imax jmin jmax} where:
• imin Minimum row index
• imax Maximum row index
• jmin Minimum column index
• jmax Maximum column index
::math::linearalgebra::dgetrfmatrix
Computes an LU factorization of a general matrix, using partial, pivoting with row interchanges.
Returns the permutation vector.
The factorization has the form
P * A = L * U
where P is a permutation matrix, L is lower triangular with unit diagonal elements, and U is upper
triangular. Returns the permutation vector, as a list of length n-1. The last entry of the
permutation is not stored, since it is implicitely known, with value n (the last row is not
swapped with any other row). At index #i of the permutation is stored the index of the row #j
which is swapped with row #i at step #i. That means that each index of the permutation gives the
permutation at each step, not the cumulated permutation matrix, which is the product of
permutations.
list matrix
On entry, the matrix to be factored. On exit, the factors L and U from the factorization
P*A = L*U; the unit diagonal elements of L are not stored.
::math::linearalgebra::detmatrix
Returns the determinant of the given matrix, based on PA=LU decomposition, i.e. Gauss partial
pivotal.
list matrix
Square matrix (matrix A)
list ipiv
The pivots (optionnal). If the pivots are not provided, a PA=LU decomposition is
performed. If the pivots are provided, we assume that it contains the pivots and that the
matrix A contains the L and U factors, as provided by dgterf. b-vectors
::math::linearalgebra::largesteigenmatrixtolerancemaxiter
Returns a list made of the largest eigenvalue (in magnitude) and associated eigenvector. Uses
iterative Power Method as provided as algorithm #7.3.3 of Golub & Van Loan. This algorithm is
used here for a dense matrix (but is usually used for sparse matrices).
list matrix
Square matrix (matrix A)
double tolerance
The relative tolerance of the eigenvalue (default:1.e-8).
integer maxiter
The maximum number of iterations (default:10).
CompabilitywiththeLApackage Two procedures are provided for compatibility with Hume's LA package:
::math::linearalgebra::to_LAmv
Transforms a vector or matrix into the format used by the original LA package.
list mv
Matrix or vector
::math::linearalgebra::from_LAmv
Transforms a vector or matrix from the format used by the original LA package into the format used
by the present implementation.
list mv
Matrix or vector as used by the LA package