Next Previous Contents

8. gslmatrix: A Collection of Matrix-Oriented GSL functions

The S-lang interpreter has wide-spread native support for manipulating arrays and matrices. The gslmatrix supplements this by adding some linear algebra routines such LU decomposition as well as routines for dealing with eigenvalues and eigenvectors.

GSL has separate functions for complex numbers. Rather than creating separater wrappers for each of these functions, the complex-valued routines have been incorporated into single wrappers that support for both real and complex numbers. In this way the interface is polymorphic.

8.1 Linear Algebra and Matrix-Oriented Routines

linalg_LU_decomp

Synopsis

Factorize a square matrix into its LU decomposition

Usage

(LU,p) = linalg_LU_decomp (A [,&signum])

Description

This routines returns the LU decomposition of the square matrix A such that P#A == LU. See the corresponding GSL documentation for how L and U are stored in LU, and how the permutation matrix P is defined. For many applications, it is unnecessary to unpack the matrix LU into its separate components.

If the optional argument &signum is given, upon return signum will be set to the sign of the permutation that relates P to the identity matrix.

See Also

linalg_LU_det, linalg_LU_invert, linalg_LU_solve

linalg_LU_det

Synopsis

Compute the determinant of a matrix from its LU decomposition

Usage

det = linalg_LU_det (LU, signum)

Description

This function computes the determinant of a matrix from its LU decomposition. In the LU form, determinant is given by the product of the diagonal elements with the sign of the permutation.

     require ("gslmatrix");
     define determinant (A)
     {
        variable LU, sig;
        (LU,) = linalg_LU_decomp (A, &sig);
        return linalg_LU_det (LU,sig);
     }

See Also

linalg_LU_lndet, linalg_LU_decomp, linalg_LU_invert, linalg_LU_solve

linalg_LU_lndet

Synopsis

Compute the log of a determinant using LU decomposition

Usage

det = linalg_LU_lndet (LU)

Description

This function computes the natural logarithm of the determinant of a matrix from its LU decomposition. In the LU form, determinant is given by the product of the diagonal elements with the sign of the permutation. This function is useful for cases where the product of the diagonal elements would overflow.

See Also

linalg_LU_det, linalg_LU_decomp, linalg_LU_solve, linalg_LU_invert

linalg_LU_invert

Synopsis

Compute the inverse of a matrix via its LU decomposition

Usage

inv = linalg_LU_invert (LU, p)

Description

This function may be used to compute the inverse of a matrix from its LU decomposition. For the purposes of inverting a set of linear equations, it is preferable to use the linalg_LU_solve function rather than inverting the equations via the inverse.

    define matrix_inverse (A)
    {
       return linalg_LU_invert (linalg_LU_decomp (A));
    }

See Also

linalg_LU_decomp, linalg_LU_solve, linalg_LU_det

linalg_LU_solve

Synopsis

Solve a set of linear equations using LU decomposition

Usage

x = linalg_LU_solve (LU, p, b)

Description

This function solves the square linear system of equations A#x=b for the vector x via the LU decomposition of A.

   define solve_equations (A, b)
   {
      return linalg_LU_solve (linalg_LU_decomp (A), b);
   }

See Also

linalg_LU_decomp, linalg_LU_det, linalg_LU_invert

linalg_QR_decomp

Synopsis

Factor a matrix into its QR form

Usage

(QR, tau) = linalg_QR_decomp(A)

Description

This function may be used to decompose a rectangular matrix into its so-called QR such that A=Q#R where Q is a square orthogonal matrix and R is a rectangular right-triangular matrix.

The factor R encoded in the diagonal and upper-triangular elements of the first return value QR. The matrix Q is encoded in the lower triangular part of QR and the vector tau via Householder vectors and coefficients. See the corresponding GSL documentation for the details of the encoding. For most uses encoding details are not required.

See Also

linalg_QR_solve

linalg_QR_solve

Synopsis

Solve a system of linear equations using QR decomposition

Usage

x = linalg_QR_solve(QR, tau, b [,&residual])

Description

This function may be used to solve the linear system A#x=b using the QR decomposition of A.

If the optional fourth argument is present (&residual), or if QR is not a square matrix, then the linear system will be solved in the least-squares sense by minimizing the (Euclidean) norm of A#x-b. Upon return, the value of the variable residual is set to the the norm of A#x-b.

Notes

GSL has a separate function called gsl_linalg_QR_lssolve for computing this least-squares solution. The linalg_QR_solve combines both gsl_linalg_QR_lssolve and gsl_linalg_QR_solve into a single routine.

See Also

linalg_QR_decomp

linalg_SV_decomp

Synopsis

Perform a singular-value decomposition on a matrix

Usage

(U,S,V) = linalg_SV_decomp(A)

Description

This function factors a MxN (M>=N) rectangular matrix A into three factors such that A = U#S#transpose(V), where S is diagonal matrix containing the singular values of A and V is a square orthogonal matrix. Since S is diagonal, it is returned as a 1-d array.

See Also

linalg_SV_solve

linalg_SV_solve

Synopsis

Solve a linear system using Singular-Value Decomposition

Usage

x = linalg_SV_solve (U,V,S,b)

Description

This function ``solves'' the linear system A#x=b using the SVD form of A.

Example

   define svd_solve (A, b)
   {
      variable U, V, S;
      (U,V,S) = linalg_SV_decomp (A);
      return linalg_SV_solve (U,V,S,b);
   }

See Also

linalg_SV_decomp, linalg_QR_solve, linalg_LU_solve

eigen_symmv

Synopsis

Compute the eigenvalues and eigenvectors of a Hermitian matrix

Usage

(eigvecs, eigvals)=eigen_symmv(A)

Description

This function computes the eigenvalues and eigenvectors of a Hermitian (or real-symmetric) square matrix A. The eigenvalues are returned sorted on their absolute value (or norm) in descending order.

See Also

eigen_nonsymmv

eigen_nonsymmv

Synopsis

Compute the eigenvalues and eigenvectors of a matrix

Usage

(eigvecs, eigvals)=eigen_nonsymmv(A)

Description

This function returns the eigenvalues and eigenvectors of a real non-symmetric matrix A. As such quantities are in general complex, complex-valued arrays will be returned. The eigenvalues are returned in descending order sorted upon norm.

See Also

eigen_symmv


Next Previous Contents