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

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

Private Member Functions

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

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

Usage
The following example illustrates the least squares solution of an overdetermined system of linear equations.
program example
use iso_fortran_env, only : real64, int32
implicit none
! Local Variables
real(real64) :: a(3,2), b(3)
integer(int32) :: i
! Build the 3-by-2 matrix A
! | 2 1 |
! A = |-3 1 |
! |-1 1 |
a = reshape([2.0d0, -3.0d0, -1.0d0, 1.0d0, 1.0d0, 1.0d0], [3, 2])
! Build the right-hand-side vector B.
! |-1 |
! b = |-2 |
! | 1 |
b = [-1.0d0, -2.0d0, 1.0d0]
! The solution is:
! x = [0.13158, -0.57895]**T
! Compute the solution via a least-squares approach. The results overwrite
! the first 2 elements in b.
! Display the results
print '(A)', "Least Squares Solution: X = "
print '(F9.5)', (b(i), i = 1, size(a, 2))
end program
The above program produces the following output.
Least Squares Solution: X =
0.13158
-0.57895

Definition at line 457 of file linalg_solve.f90.

Member Function/Subroutine Documentation

◆ solve_least_squares_mtx_svd()

subroutine linalg_solve::solve_least_squares_svd::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_svd()

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


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