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

Solves a system of M QR-factored equations of N unknowns. More...

Private Member Functions

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

Detailed Description

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

Usage
The following example illustrates the solution of a system of equations using QR factorization.
program example
use iso_fortran_env, only : real64, int32
use linalg_solve, only : solve_qr
! Local Variables
real(real64) :: a(3,3), tau(3), b(3)
integer(int32) :: i, pvt(3)
! Build the 3-by-3 matrix A.
! | 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])
! Build the right-hand-side vector B.
! | -1 |
! b = | -2 |
! | -3 |
b = [-1.0d0, -2.0d0, -3.0d0]
! The solution is:
! | 1/3 |
! x = | -2/3 |
! | 0 |
! Compute the QR factorization, using pivoting
pvt = 0 ! Zero every entry in order not to lock any column in place
call qr_factor(a, tau, pvt)
! Compute the solution. The results overwrite b.
call solve_qr(a, tau, pvt, b)
! Display the results.
print '(A)', "QR Solution: X = "
print '(F8.4)', (b(i), i = 1, size(b))
! Notice, QR factorization without pivoting could be accomplished in the
! same manner. The only difference is to omit the PVT array (column pivot
! tracking array).
end program
The above program produces the following output.
QR Solution: X =
0.3333
-0.6667
0.0000
See Also

Definition at line 222 of file linalg_solve.f90.

Member Function/Subroutine Documentation

◆ solve_qr_no_pivot_mtx()

subroutine linalg_solve::solve_qr::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::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::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::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.


The documentation for this interface was generated from the following file: