linalg  1.4.3
A linear algebra library that provides a user-friendly interface to several BLAS and LAPACK routines.
linalg_solve Module Reference

linalg_solve More...

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...
 

Detailed Description

linalg_solve

Purpose
Provides a set of routines for solving systems of linear equations.

Function/Subroutine Documentation

◆ mtx_inverse()

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.

Parameters
[in,out]aOn input, the N-by-N matrix to invert. On output, the inverted matrix.
[out]iworkAn optional N-element integer workspace array.
[out]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if a is not square. Will also occur if incorrectly sized workspace arrays are provided.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_SINGULAR_MATRIX_ERROR: Occurs if the input matrix is singular.
Usage
The following example illustrates the inversion of a 3-by-3 matrix.
program example
use iso_fortran_env, only : real64, int32
use linalg_solve, only : mtx_inverse
implicit none
! Variables
real(real64) :: a(3,3), ai(3,3), c(3,3)
integer(int32) :: i
! Construct the 3-by-3 matrix A to invert
! | 1 2 3 |
! A = | 4 5 6 |
! | 7 8 0 |
a = reshape([1.0d0, 4.0d0, 7.0d0, 2.0d0, 5.0d0, 8.0d0, 3.0d0, 6.0d0, &
0.0d0], [3, 3])
! Compute the inverse of A. Notice, the original matrix is overwritten
! with it's inverse.
ai = a
call mtx_inverse(ai)
! Show that A * inv(A) = I
c = matmul(a, ai)
! Display the inverse
print '(A)', "Inverse:"
do i = 1, size(ai, 1)
print *, ai(i,:)
end do
! Display the result of A * inv(A)
print '(A)', "A * A**-1:"
do i = 1, size(c, 1)
print *, c(i,:)
end do
end program
The above program produces the following output.
Inverse:
-1.7777777777777777 0.88888888888888884 -0.11111111111111110
1.5555555555555556 -0.77777777777777779 0.22222222222222221
-0.11111111111111119 0.22222222222222227 -0.11111111111111112
A * A**-1:
0.99999999999999989 5.5511151231257827E-017 -4.1633363423443370E-017
5.5511151231257827E-017 1.0000000000000000 -8.3266726846886741E-017
1.7763568394002505E-015 -8.8817841970012523E-016 1.0000000000000000
Notes
This routine utilizes the LAPACK routines DGETRF to perform an LU factorization of the matrix, and DGETRI to invert the LU factored matrix.
See Also

Definition at line 1598 of file linalg_solve.f90.

◆ mtx_pinverse()

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.

Parameters
[in,out]aOn input, the M-by-N matrix to invert. The matrix is overwritten on output.
[out]ainvThe N-by-M matrix where the pseudo-inverse of a will be written.
[in]tolAn 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]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_CONVERGENCE_ERROR: Occurs as a warning if the QR iteration process could not converge to a zero value.
Usage
The following example illustrates how to compute the Moore-Penrose pseudo-inverse of a matrix.
program example
use iso_fortran_env, only : int32, real64
use linalg_solve, only : mtx_pinverse
implicit none
! Variables
real(real64) :: a(3,2), ai(2,3), ao(3,2), c(2,2)
integer(int32) :: i
! Create the 3-by-2 matrix A
! | 1 0 |
! A = | 0 1 |
! | 0 1 |
a = reshape([1.0d0, 0.0d0, 0.0d0, 0.0d0, 1.0d0, 1.0d0], [3, 2])
ao = a ! Just making a copy for later as mtx_pinverse will destroy the
! contents of the original matrix
! The Moore-Penrose pseudo-inverse of this matrix is:
! | 1 0 0 |
! A**-1 = | |
! | 0 1/2 1/2 |
call mtx_pinverse(a, ai)
! Notice, A**-1 * A is an identity matrix.
c = matmul(ai, ao)
! Display the inverse
print '(A)', "Inverse:"
do i = 1, size(ai, 1)
print *, ai(i,:)
end do
! Display the result of inv(A) * A
print '(A)', "A**-1 * A:"
do i = 1, size(c, 1)
print *, c(i,:)
end do
end program
The above program produces the following output.
Inverse:
1.0000000000000000 0.0000000000000000 0.0000000000000000
0.0000000000000000 0.49999999999999978 0.49999999999999989
A**-1 * A:
1.0000000000000000 0.0000000000000000
0.0000000000000000 0.99999999999999967
See Also

Definition at line 1789 of file linalg_solve.f90.

◆ solve_cholesky_mtx()

subroutine linalg_solve::solve_cholesky_mtx ( logical, intent(in)  upper,
real(real64), dimension(:,:), intent(in)  a,
real(real64), dimension(:,:), intent(inout)  b,
class(errors), intent(inout), optional, target  err 
)
private

Solves a system of Cholesky factored equations.

Parameters
[in]upperSet 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]aThe N-by-N Cholesky factored matrix.
[in,out]bOn input, the N-by-NRHS right-hand-side matrix B. On output, the solution matrix X.
[out]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input array sizes are incorrect.
Notes
This routine utilizes the LAPACK routine DPOTRS.

Definition at line 1400 of file linalg_solve.f90.

◆ solve_cholesky_vec()

subroutine linalg_solve::solve_cholesky_vec ( logical, intent(in)  upper,
real(real64), dimension(:,:), intent(in)  a,
real(real64), dimension(:), intent(inout)  b,
class(errors), intent(inout), optional, target  err 
)
private

Solves a system of Cholesky factored equations.

Parameters
[in]upperSet 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]aThe N-by-N Cholesky factored matrix.
[in,out]bOn input, the N-element right-hand-side vector B. On output, the solution vector X.
[out]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input array sizes are incorrect.
Notes
This routine utilizes the LAPACK routine DPOTRS.

Definition at line 1466 of file linalg_solve.f90.

◆ solve_least_squares_mtx()

subroutine linalg_solve::solve_least_squares_mtx ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:,:), intent(inout)  b,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
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.

Parameters
[in,out]aOn 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]bIf 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]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_INVALID_OPERATION_ERROR: Occurs if a is not of full rank.
Notes
This routine utilizes the LAPACK routine DGELS.

Definition at line 1966 of file linalg_solve.f90.

◆ solve_least_squares_mtx_pvt()

subroutine linalg_solve::solve_least_squares_mtx_pvt ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:,:), intent(inout)  b,
integer(int32), dimension(:), intent(inout), optional, target  ipvt,
integer(int32), intent(out), optional  arnk,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a complete orthogonal factorization of matrix A.

Parameters
[in,out]aOn input, the M-by-N matrix A. On output, the matrix is overwritten by the details of its complete orthogonal factorization.
[in,out]bIf 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]ipvtAn 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]arnkAn optional output, that if provided, will return the rank of a.
[out]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
Notes
This routine utilizes the LAPACK routine DGELSY.

Definition at line 2184 of file linalg_solve.f90.

◆ solve_least_squares_mtx_svd()

subroutine linalg_solve::solve_least_squares_mtx_svd ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:,:), intent(inout)  b,
real(real64), dimension(:), intent(out), optional, target  s,
integer(int32), intent(out), optional  arnk,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a singular value decomposition of matrix A.

Parameters
[in,out]aOn input, the M-by-N matrix A. On output, the matrix is overwritten by the details of its complete orthogonal factorization.
[in,out]bIf 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]arnkAn optional output, that if provided, will return the rank of a.
[out]sAn 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]arnkAn optional output, that if provided, will return the rank of a.
[out]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_CONVERGENCE_ERROR: Occurs as a warning if the QR iteration process could not converge to a zero value.
Notes
This routine utilizes the LAPACK routine DGELSS.

Definition at line 2472 of file linalg_solve.f90.

◆ solve_least_squares_vec()

subroutine linalg_solve::solve_least_squares_vec ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:), intent(inout)  b,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
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.

Parameters
[in,out]aOn 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]bIf 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]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_INVALID_OPERATION_ERROR: Occurs if a is not of full rank.
Notes
This routine utilizes the LAPACK routine DGELS.

Definition at line 2072 of file linalg_solve.f90.

◆ solve_least_squares_vec_pvt()

subroutine linalg_solve::solve_least_squares_vec_pvt ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:), intent(inout)  b,
integer(int32), dimension(:), intent(inout), optional, target  ipvt,
integer(int32), intent(out), optional  arnk,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a complete orthogonal factorization of matrix A.

Parameters
[in,out]aOn input, the M-by-N matrix A. On output, the matrix is overwritten by the details of its complete orthogonal factorization.
[in,out]bIf 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]ipvtAn 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]arnkAn optional output, that if provided, will return the rank of a.
[out]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
Notes
This routine utilizes the LAPACK routine DGELSY.

Definition at line 2328 of file linalg_solve.f90.

◆ solve_least_squares_vec_svd()

subroutine linalg_solve::solve_least_squares_vec_svd ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:), intent(inout)  b,
real(real64), dimension(:), intent(out), optional, target  s,
integer(int32), intent(out), optional  arnk,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Solves the overdetermined or underdetermined system (A*X = B) of M equations of N unknowns using a singular value decomposition of matrix A.

Parameters
[in,out]aOn input, the M-by-N matrix A. On output, the matrix is overwritten by the details of its complete orthogonal factorization.
[in,out]bIf 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]sAn 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]arnkAn optional output, that if provided, will return the rank of a.
[out]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
  • LA_CONVERGENCE_ERROR: Occurs as a warning if the QR iteration process could not converge to a zero value.
Notes
This routine utilizes the LAPACK routine DGELSS.

Definition at line 2619 of file linalg_solve.f90.

◆ solve_lu_mtx()

subroutine linalg_solve::solve_lu_mtx ( real(real64), dimension(:,:), intent(in)  a,
integer(int32), dimension(:), intent(in)  ipvt,
real(real64), dimension(:,:), intent(inout)  b,
class(errors), intent(inout), optional, target  err 
)
private

Solves a system of LU-factored equations.

Parameters
[in]aThe N-by-N LU factored matrix as output by lu_factor.
[in]ipvtThe N-element pivot array as output by lu_factor.
[in,out]bOn input, the N-by-NRHS right-hand-side matrix. On output, the N-by-NRHS solution matrix.
[out]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input array sizes are incorrect.
Notes
The routine is based upon the LAPACK routine DGETRS.

Definition at line 678 of file linalg_solve.f90.

◆ solve_lu_vec()

subroutine linalg_solve::solve_lu_vec ( real(real64), dimension(:,:), intent(in)  a,
integer(int32), dimension(:), intent(in)  ipvt,
real(real64), dimension(:), intent(inout)  b,
class(errors), intent(inout), optional, target  err 
)
private

Solves a system of LU-factored equations.

Parameters
[in]aThe N-by-N LU factored matrix as output by lu_factor.
[in]ipvtThe N-element pivot array as output by lu_factor.
[in,out]bOn input, the N-element right-hand-side array. On output, the N-element solution array.
[out]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input array sizes are incorrect.
Notes
The routine is based upon the LAPACK routine DGETRS.

Definition at line 739 of file linalg_solve.f90.

◆ solve_qr_no_pivot_mtx()

subroutine linalg_solve::solve_qr_no_pivot_mtx ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:), intent(in)  tau,
real(real64), dimension(:,:), intent(inout)  b,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Solves a system of M QR-factored equations of N unknowns where M >= N.

Parameters
[in]aOn 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]tauA MIN(M, N)-element array containing the scalar factors of the elementary reflectors as returned by qr_factor.
[in]bOn input, the M-by-NRHS right-hand-side matrix. On output, the first N columns are overwritten by the solution matrix X.
[out]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
Notes
This routine is based upon a subset of the LAPACK routine DGELS.

Definition at line 814 of file linalg_solve.f90.

◆ solve_qr_no_pivot_vec()

subroutine linalg_solve::solve_qr_no_pivot_vec ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:), intent(in)  tau,
real(real64), dimension(:), intent(inout)  b,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Solves a system of M QR-factored equations of N unknowns where M >= N.

Parameters
[in]aOn 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]tauA MIN(M, N)-element array containing the scalar factors of the elementary reflectors as returned by qr_factor.
[in]bOn input, the M-element right-hand-side vector. On output, the first N elements are overwritten by the solution vector X.
[out]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
Notes
This routine is based upon a subset of the LAPACK routine DGELS.

Definition at line 929 of file linalg_solve.f90.

◆ solve_qr_pivot_mtx()

subroutine linalg_solve::solve_qr_pivot_mtx ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:), intent(in)  tau,
integer(int32), dimension(:), intent(in)  jpvt,
real(real64), dimension(:,:), intent(inout)  b,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Solves a system of M QR-factored equations of N unknowns where the QR factorization made use of column pivoting.

Parameters
[in]aOn input, the M-by-N QR factored matrix as returned by qr_factor. On output, the contents of this matrix are altered.
[in]tauA MIN(M, N)-element array containing the scalar factors of the elementary reflectors as returned by qr_factor.
[in]jpvtAn N-element array, as output by qr_factor, used to track the column pivots.
[in]bOn 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]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
Notes
This routine is based upon a subset of the LAPACK routine DGELSY.

Definition at line 1042 of file linalg_solve.f90.

◆ solve_qr_pivot_vec()

subroutine linalg_solve::solve_qr_pivot_vec ( real(real64), dimension(:,:), intent(inout)  a,
real(real64), dimension(:), intent(in)  tau,
integer(int32), dimension(:), intent(in)  jpvt,
real(real64), dimension(:), intent(inout)  b,
real(real64), dimension(:), intent(out), optional, target  work,
integer(int32), intent(out), optional  olwork,
class(errors), intent(inout), optional, target  err 
)
private

Solves a system of M QR-factored equations of N unknowns where the QR factorization made use of column pivoting.

Parameters
[in]aOn input, the M-by-N QR factored matrix as returned by qr_factor. On output, the contents of this matrix are altered.
[in]tauA MIN(M, N)-element array containing the scalar factors of the elementary reflectors as returned by qr_factor.
[in]jpvtAn N-element array, as output by qr_factor, used to track the column pivots.
[in]bOn 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]workAn 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]olworkAn 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]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if any of the input arrays are not sized appropriately.
  • LA_OUT_OF_MEMORY_ERROR: Occurs if local memory must be allocated, and there is insufficient memory available.
Notes
This routine is based upon a subset of the LAPACK routine DGELSY.

Definition at line 1229 of file linalg_solve.f90.

◆ solve_tri_mtx()

subroutine linalg_solve::solve_tri_mtx ( logical, intent(in)  lside,
logical, intent(in)  upper,
logical, intent(in)  trans,
logical, intent(in)  nounit,
real(real64), intent(in)  alpha,
real(real64), dimension(:,:), intent(in)  a,
real(real64), dimension(:,:), intent(inout)  b,
class(errors), intent(inout), optional, target  err 
)
private

Solves one of the matrix equations: op(A) * X = alpha * B, or X * op(A) = alpha * B, where A is a triangular matrix.

Parameters
[in]lsideSet to true to solve op(A) * X = alpha * B; else, set to false to solve X * op(A) = alpha * B.
[in]upperSet to true if A is an upper triangular matrix; else, set to false if A is a lower triangular matrix.
[in]transSet to true if op(A) = A**T; else, set to false if op(A) = A.
[in]nounitSet 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]alphaThe scalar multiplier to B.
[in]aIf 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]bOn input, the M-by-N right-hand-side. On output, the M-by-N solution.
[out]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if a is not square, or if the sizes of a and b are not compatible.
Notes
This routine is based upon the BLAS routine DTRSM.

Definition at line 496 of file linalg_solve.f90.

◆ solve_tri_vec()

subroutine linalg_solve::solve_tri_vec ( logical, intent(in)  upper,
logical, intent(in)  trans,
logical, intent(in)  nounit,
real(real64), dimension(:,:), intent(in)  a,
real(real64), dimension(:), intent(inout)  x,
class(errors), intent(inout), optional, target  err 
)
private

Solves the system of equations: op(A) * X = B, where A is a triangular matrix.

Parameters
[in]upperSet to true if A is an upper triangular matrix; else, set to false if A is a lower triangular matrix.
[in]transSet to true if op(A) = A**T; else, set to false if op(A) = A.
[in]nounitSet 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]aThe N-by-N triangular matrix.
[in,out]xOn input, the N-element right-hand-side array. On output, the N-element solution array.
[out]errAn 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.
  • LA_ARRAY_SIZE_ERROR: Occurs if a is not square, or if the sizes of a and b are not compatible.
Usage
To solve a triangular system of N equations of N unknowns A*X = B, where A is an N-by-N upper triangular matrix, and B and X are N-element arrays, the following code will suffice.
! Solve the system: A*X = B, where A is an upper triangular N-by-N
! matrix, and B and X are N-elements in size.
! Variables
integer(int32) :: info
real(real64), dimension(n, n) :: a
real(real64), dimension(n) :: b
! Initialize A and B...
! Solve A*X = B for X - Note: X overwrites B.
call solve_triangular_system(.true., .false., a, b)
Notes
This routine is based upon the BLAS routine DTRSV.

Definition at line 602 of file linalg_solve.f90.