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

Solves a triangular system of equations. More...

Private Member Functions

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

Detailed Description

Solves a triangular system of equations.

Usage
The following example illustrates the solution of two triangular systems to solve a system of LU factored equations.
program example
use iso_fortran_env, only : real64, int32
use linalg_factor, only : lu_factor, form_lu
implicit none
! Variables
real(real64) :: a(3,3), b(3), u(3,3), p(3,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 LU factorization
call lu_factor(a, pvt)
! Extract the L and U matrices. A is overwritten with L.
call form_lu(a, pvt, u, p)
! Solve the lower triangular system L * Y = P * B for Y, but first compute
! P * B, and store the results in B
b = matmul(p, b)
! Now, compute the solution to the lower triangular system. Store the
! result in B. Remember, L is unit diagonal (ones on its diagonal)
call solve_triangular_system(.false., .false., .false., a, b)
! Solve the upper triangular system U * X = Y for X.
call solve_triangular_system(.true., .false., .true., u, b)
! Display the results.
print '(A)', "LU Solution: X = "
print '(F8.4)', (b(i), i = 1, size(b))
end program
The above program produces the following output.
LU Solution: X =
0.3333
-0.6667
0.0000

Definition at line 92 of file linalg_solve.f90.

Member Function/Subroutine Documentation

◆ solve_tri_mtx()

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


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