linalg
1.4.3
A linear algebra library that provides a user-friendly interface to several BLAS and LAPACK routines.
|
Data Types | |
interface | solve_cholesky |
Solves a system of Cholesky factored equations. More... | |
interface | solve_least_squares |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns. More... | |
interface | solve_least_squares_full |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns, but uses a full orthogonal factorization of the system. More... | |
interface | solve_least_squares_svd |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a singular value decomposition of matrix A. More... | |
interface | solve_lu |
Solves a system of LU-factored equations. More... | |
interface | solve_qr |
Solves a system of M QR-factored equations of N unknowns. More... | |
interface | solve_triangular_system |
Solves a triangular system of equations. More... | |
Functions/Subroutines | |
subroutine | solve_tri_mtx (lside, upper, trans, nounit, alpha, a, b, err) |
Solves one of the matrix equations: op(A) * X = alpha * B, or X * op(A) = alpha * B, where A is a triangular matrix. More... | |
subroutine | solve_tri_vec (upper, trans, nounit, a, x, err) |
Solves the system of equations: op(A) * X = B, where A is a triangular matrix. More... | |
subroutine | solve_lu_mtx (a, ipvt, b, err) |
Solves a system of LU-factored equations. More... | |
subroutine | solve_lu_vec (a, ipvt, b, err) |
Solves a system of LU-factored equations. More... | |
subroutine | solve_qr_no_pivot_mtx (a, tau, b, work, olwork, err) |
Solves a system of M QR-factored equations of N unknowns where M >= N. More... | |
subroutine | solve_qr_no_pivot_vec (a, tau, b, work, olwork, err) |
Solves a system of M QR-factored equations of N unknowns where M >= N. More... | |
subroutine | solve_qr_pivot_mtx (a, tau, jpvt, b, work, olwork, err) |
Solves a system of M QR-factored equations of N unknowns where the QR factorization made use of column pivoting. More... | |
subroutine | solve_qr_pivot_vec (a, tau, jpvt, b, work, olwork, err) |
Solves a system of M QR-factored equations of N unknowns where the QR factorization made use of column pivoting. More... | |
subroutine | solve_cholesky_mtx (upper, a, b, err) |
Solves a system of Cholesky factored equations. More... | |
subroutine | solve_cholesky_vec (upper, a, b, err) |
Solves a system of Cholesky factored equations. More... | |
subroutine, public | mtx_inverse (a, iwork, work, olwork, err) |
Computes the inverse of a square matrix. More... | |
subroutine, public | mtx_pinverse (a, ainv, tol, work, olwork, err) |
Computes the Moore-Penrose pseudo-inverse of a M-by-N matrix using the singular value decomposition of the matrix. More... | |
subroutine | solve_least_squares_mtx (a, b, work, olwork, err) |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a QR or LQ factorization of the matrix A. Notice, it is assumed that matrix A has full rank. More... | |
subroutine | solve_least_squares_vec (a, b, work, olwork, err) |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a QR or LQ factorization of the matrix A. Notice, it is assumed that matrix A has full rank. More... | |
subroutine | solve_least_squares_mtx_pvt (a, b, ipvt, arnk, work, olwork, err) |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a complete orthogonal factorization of matrix A. More... | |
subroutine | solve_least_squares_vec_pvt (a, b, ipvt, arnk, work, olwork, err) |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a complete orthogonal factorization of matrix A. More... | |
subroutine | solve_least_squares_mtx_svd (a, b, s, arnk, work, olwork, err) |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a singular value decomposition of matrix A. More... | |
subroutine | solve_least_squares_vec_svd (a, b, s, arnk, work, olwork, err) |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a singular value decomposition of matrix A. More... | |
subroutine, public linalg_solve::mtx_inverse | ( | real(real64), dimension(:,:), intent(inout) | a, |
integer(int32), dimension(:), intent(out), optional, target | iwork, | ||
real(real64), dimension(:), intent(out), optional, target | work, | ||
integer(int32), intent(out), optional | olwork, | ||
class(errors), intent(inout), optional, target | err | ||
) |
Computes the inverse of a square matrix.
[in,out] | a | On input, the N-by-N matrix to invert. On output, the inverted matrix. |
[out] | iwork | An optional N-element integer workspace array. |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 1598 of file linalg_solve.f90.
subroutine, public linalg_solve::mtx_pinverse | ( | real(real64), dimension(:,:), intent(inout) | a, |
real(real64), dimension(:,:), intent(out) | ainv, | ||
real(real64), intent(in), optional | tol, | ||
real(real64), dimension(:), intent(out), optional, target | work, | ||
integer(int32), intent(out), optional | olwork, | ||
class(errors), intent(inout), optional, target | err | ||
) |
Computes the Moore-Penrose pseudo-inverse of a M-by-N matrix using the singular value decomposition of the matrix.
[in,out] | a | On input, the M-by-N matrix to invert. The matrix is overwritten on output. |
[out] | ainv | The N-by-M matrix where the pseudo-inverse of a will be written. |
[in] | tol | An optional input, that if supplied, overrides the default tolerance on singular values such that singular values less than this tolerance are forced to have a reciprocal of zero, as opposed to 1/S(I). The default tolerance is: MAX(M, N) * EPS * MAX(S). If the supplied value is less than a value that causes an overflow, the tolerance reverts back to its default value, and the operation continues; however, a warning message is issued. |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 1789 of file linalg_solve.f90.
|
private |
Solves a system of Cholesky factored equations.
[in] | upper | Set to true if the original matrix A was factored such that A = U**T * U; else, set to false if the factorization of A was A = L**T * L. |
[in] | a | The N-by-N Cholesky factored matrix. |
[in,out] | b | On input, the N-by-NRHS right-hand-side matrix B. On output, the solution matrix X. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 1400 of file linalg_solve.f90.
|
private |
Solves a system of Cholesky factored equations.
[in] | upper | Set to true if the original matrix A was factored such that A = U**T * U; else, set to false if the factorization of A was A = L**T * L. |
[in] | a | The N-by-N Cholesky factored matrix. |
[in,out] | b | On input, the N-element right-hand-side vector B. On output, the solution vector X. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 1466 of file linalg_solve.f90.
|
private |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a QR or LQ factorization of the matrix A. Notice, it is assumed that matrix A has full rank.
[in,out] | a | On input, the M-by-N matrix A. On output, if M >= N, the QR factorization of A in the form as output by qr_factor; else, if M < N, the LQ factorization of A. |
[in,out] | b | If M >= N, the M-by-NRHS matrix B. On output, the first N rows contain the N-by-NRHS solution matrix X. If M < N, an N-by-NRHS matrix with the first M rows containing the matrix B. On output, the N-by-NRHS solution matrix X. |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 1966 of file linalg_solve.f90.
|
private |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a complete orthogonal factorization of matrix A.
[in,out] | a | On input, the M-by-N matrix A. On output, the matrix is overwritten by the details of its complete orthogonal factorization. |
[in,out] | b | If M >= N, the M-by-NRHS matrix B. On output, the first N rows contain the N-by-NRHS solution matrix X. If M < N, an N-by-NRHS matrix with the first M rows containing the matrix B. On output, the N-by-NRHS solution matrix X. |
[out] | ipvt | An optional input that on input, an N-element array that if IPVT(I) .ne. 0, the I-th column of A is permuted to the front of A * P; if IPVT(I) = 0, the I-th column of A is a free column. On output, if IPVT(I) = K, then the I-th column of A * P was the K-th column of A. If not supplied, memory is allocated internally, and IPVT is set to all zeros such that all columns are treated as free. |
[out] | arnk | An optional output, that if provided, will return the rank of a . |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 2184 of file linalg_solve.f90.
|
private |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a singular value decomposition of matrix A.
[in,out] | a | On input, the M-by-N matrix A. On output, the matrix is overwritten by the details of its complete orthogonal factorization. |
[in,out] | b | If M >= N, the M-by-NRHS matrix B. On output, the first N rows contain the N-by-NRHS solution matrix X. If M < N, an N-by-NRHS matrix with the first M rows containing the matrix B. On output, the N-by-NRHS solution matrix X. |
[out] | arnk | An optional output, that if provided, will return the rank of a . |
[out] | s | An optional MIN(M, N)-element array that on output contains the singular values of a in descending order. Notice, the condition number of a can be determined by S(1) / S(MIN(M, N)). |
[out] | arnk | An optional output, that if provided, will return the rank of a . |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 2472 of file linalg_solve.f90.
|
private |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a QR or LQ factorization of the matrix A. Notice, it is assumed that matrix A has full rank.
[in,out] | a | On input, the M-by-N matrix A. On output, if M >= N, the QR factorization of A in the form as output by qr_factor; else, if M < N, the LQ factorization of A. |
[in,out] | b | If M >= N, the M-element array B. On output, the first N elements contain the N-element solution array X. If M < N, an N-element array with the first M elements containing the array B. On output, the N-element solution array X. |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 2072 of file linalg_solve.f90.
|
private |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a complete orthogonal factorization of matrix A.
[in,out] | a | On input, the M-by-N matrix A. On output, the matrix is overwritten by the details of its complete orthogonal factorization. |
[in,out] | b | If M >= N, the M-element array B. On output, the first N elements contain the N-element solution array X. If M < N, an N-element array with the first M elements containing the array B. On output, the N-element solution array X. |
[out] | ipvt | An optional input that on input, an N-element array that if IPVT(I) .ne. 0, the I-th column of A is permuted to the front of A * P; if IPVT(I) = 0, the I-th column of A is a free column. On output, if IPVT(I) = K, then the I-th column of A * P was the K-th column of A. If not supplied, memory is allocated internally, and IPVT is set to all zeros such that all columns are treated as free. |
[out] | arnk | An optional output, that if provided, will return the rank of a . |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 2328 of file linalg_solve.f90.
|
private |
Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a singular value decomposition of matrix A.
[in,out] | a | On input, the M-by-N matrix A. On output, the matrix is overwritten by the details of its complete orthogonal factorization. |
[in,out] | b | If M >= N, the M-by-NRHS matrix B. On output, the first N rows contain the N-by-NRHS solution matrix X. If M < N, an N-by-NRHS matrix with the first M rows containing the matrix B. On output, the N-by-NRHS solution matrix X. |
[out] | s | An optional MIN(M, N)-element array that on output contains the singular values of a in descending order. Notice, the condition number of a can be determined by S(1) / S(MIN(M, N)). |
[out] | arnk | An optional output, that if provided, will return the rank of a . |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 2619 of file linalg_solve.f90.
|
private |
Solves a system of LU-factored equations.
[in] | a | The N-by-N LU factored matrix as output by lu_factor. |
[in] | ipvt | The N-element pivot array as output by lu_factor. |
[in,out] | b | On input, the N-by-NRHS right-hand-side matrix. On output, the N-by-NRHS solution matrix. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 678 of file linalg_solve.f90.
|
private |
Solves a system of LU-factored equations.
[in] | a | The N-by-N LU factored matrix as output by lu_factor. |
[in] | ipvt | The N-element pivot array as output by lu_factor. |
[in,out] | b | On input, the N-element right-hand-side array. On output, the N-element solution array. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 739 of file linalg_solve.f90.
|
private |
Solves a system of M QR-factored equations of N unknowns where M >= N.
[in] | a | On input, the M-by-N QR factored matrix as returned by qr_factor. On output, the contents of this matrix are restored. Notice, M must be greater than or equal to N. |
[in] | tau | A MIN(M, N)-element array containing the scalar factors of the elementary reflectors as returned by qr_factor. |
[in] | b | On input, the M-by-NRHS right-hand-side matrix. On output, the first N columns are overwritten by the solution matrix X. |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 814 of file linalg_solve.f90.
|
private |
Solves a system of M QR-factored equations of N unknowns where M >= N.
[in] | a | On input, the M-by-N QR factored matrix as returned by qr_factor. On output, the contents of this matrix are restored. Notice, M must be greater than or equal to N. |
[in] | tau | A MIN(M, N)-element array containing the scalar factors of the elementary reflectors as returned by qr_factor. |
[in] | b | On input, the M-element right-hand-side vector. On output, the first N elements are overwritten by the solution vector X. |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 929 of file linalg_solve.f90.
|
private |
Solves a system of M QR-factored equations of N unknowns where the QR factorization made use of column pivoting.
[in] | a | On input, the M-by-N QR factored matrix as returned by qr_factor. On output, the contents of this matrix are altered. |
[in] | tau | A MIN(M, N)-element array containing the scalar factors of the elementary reflectors as returned by qr_factor. |
[in] | jpvt | An N-element array, as output by qr_factor, used to track the column pivots. |
[in] | b | On input, the MAX(M, N)-by-NRHS matrix where the first M rows contain the right-hand-side matrix B. On output, the first N rows are overwritten by the solution matrix X. |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 1042 of file linalg_solve.f90.
|
private |
Solves a system of M QR-factored equations of N unknowns where the QR factorization made use of column pivoting.
[in] | a | On input, the M-by-N QR factored matrix as returned by qr_factor. On output, the contents of this matrix are altered. |
[in] | tau | A MIN(M, N)-element array containing the scalar factors of the elementary reflectors as returned by qr_factor. |
[in] | jpvt | An N-element array, as output by qr_factor, used to track the column pivots. |
[in] | b | On input, the MAX(M, N)-element array where the first M elements contain the right-hand-side vector B. On output, the first N elements are overwritten by the solution vector X. |
[out] | work | An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork . |
[out] | olwork | An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work , and returns without performing any actual calculations. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 1229 of file linalg_solve.f90.
|
private |
Solves one of the matrix equations: op(A) * X = alpha * B, or X * op(A) = alpha * B, where A is a triangular matrix.
[in] | lside | Set to true to solve op(A) * X = alpha * B; else, set to false to solve X * op(A) = alpha * B. |
[in] | upper | Set to true if A is an upper triangular matrix; else, set to false if A is a lower triangular matrix. |
[in] | trans | Set to true if op(A) = A**T; else, set to false if op(A) = A. |
[in] | nounit | Set to true if A is not a unit-diagonal matrix (ones on every diagonal element); else, set to false if A is a unit-diagonal matrix. |
[in] | alpha | The scalar multiplier to B. |
[in] | a | If lside is true, the M-by-M triangular matrix on which to operate; else, if lside is false, the N-by-N triangular matrix on which to operate. |
[in,out] | b | On input, the M-by-N right-hand-side. On output, the M-by-N solution. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 496 of file linalg_solve.f90.
|
private |
Solves the system of equations: op(A) * X = B, where A is a triangular matrix.
[in] | upper | Set to true if A is an upper triangular matrix; else, set to false if A is a lower triangular matrix. |
[in] | trans | Set to true if op(A) = A**T; else, set to false if op(A) = A. |
[in] | nounit | Set to true if A is not a unit-diagonal matrix (ones on every diagonal element); else, set to false if A is a unit-diagonal matrix. |
[in] | a | The N-by-N triangular matrix. |
[in,out] | x | On input, the N-element right-hand-side array. On output, the N-element solution array. |
[out] | err | An optional errors-based object that if provided can be used to retrieve information relating to any errors encountered during execution. If not provided, a default implementation of the errors class is used internally to provide error handling. Possible errors and warning messages that may be encountered are as follows.
|
Definition at line 602 of file linalg_solve.f90.