linalg 1.6.1
A linear algebra library that provides a user-friendly interface to several BLAS and LAPACK routines.
Loading...
Searching...
No Matches
linalg_core::cholesky_factor Interface Reference

Computes the Cholesky factorization of a symmetric, positive definite matrix. More...

Detailed Description

Computes the Cholesky factorization of a symmetric, positive definite matrix.

Syntax
subroutine cholesky_factor(real(real64) a(:,:), optional logical upper, optional class(errors) err)
subroutine cholesky_factor(complex(real64) a(:,:), optional logical upper, optional class(errors) err)
Parameters
[in,out]aOn input, the N-by-N matrix to factor. On output, the factored matrix is returned in either the upper or lower triangular portion of the matrix, dependent upon the value of upper.
[in]upperAn optional input that, if specified, provides control over whether the factorization is computed as \( A = U^T U \) (set to true), or as \( A = L L^T \) (set to false). The default value is true such that \( A = U^T U \).
[in,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.
  • LA_MATRIX_FORMAT_ERROR: Occurs if a is not positive definite.
Notes
This routine utilizes the LAPACK routine DPOTRF (ZPOTRF in the complex case).
Usage
The following example illustrates the solution of a positive-definite system of equations via Cholesky factorization.
program example
use iso_fortran_env, only : real64, int32
implicit none
! Variables
real(real64) :: a(3, 3), b(3), bu(3)
integer(int32) :: i
! Build the 3-by-3 positive-definite matrix A.
! | 4 12 -16 |
! A = | 12 37 -43 |
! |-16 -43 98 |
a = reshape([4.0d0, 12.0d0, -16.0d0, 12.0d0, 37.0d0, -43.0d0, -16.0d0, &
-43.0d0, 98.0d0], [3, 3])
! Build the 3-element array B
! | 5 |
! b = | 1 |
! | 3 |
b = [5.0d0, 1.0d0, 3.0d0]
! Make a copy of B for later use - not necessary, but just for example to
! illustrate the long or manual method of solving a Cholesky factored system
bu = b
! Compute the Cholesky factorization of A considering only the upper
! triangular portion of A (the default configuration).
! Compute the solution
call solve_cholesky(.true., a, b)
! Display the results
print '(A)', "Cholesky Solution: X = "
print '(F8.4)', (b(i), i = 1, size(b))
! The solution could also be computed manually noting the Cholesky
! factorization causes A = U**T * U. Then U**T * U * X = B.
! Step 1 would then be to solve the problem U**T * Y = B, for Y.
call solve_triangular_system(.true., .true., .true., a, bu)
! Now, solve the problem U * X = Y, for X
call solve_triangular_system(.true., .false., .true., a, bu)
! Display the results
print '(A)', "Cholesky Solution (Manual Approach): X = "
print '(F8.4)', (bu(i), i = 1, size(bu))
end program
Computes the Cholesky factorization of a symmetric, positive definite matrix.
Solves a system of Cholesky factored equations.
Solves a triangular system of equations.
Provides a set of common linear algebra routines.
Definition: linalg_core.f90:15
The above program produces the following output.
Cholesky Solution: X =
239.5833
-65.6667
10.3333
Cholesky Solution (Manual Approach): X =
239.5833
-65.6667
10.3333

Definition at line 1393 of file linalg_core.f90.


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