QR_decomposition Subroutine

public subroutine QR_decomposition(A, method, Q, R)

QR decomposition of a matrix A using various methods This subroutine performs QR decomposition of a given matrix A using the specified method (Householder, Givens, Classical Gram-Schmidt, or Modified Gram-Schmidt). The output matrices Q is an orthogonal matrix and R is an upper triangular matrix.

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(in), dimension(:, :) :: A
character(len=*), intent(in), optional :: method
real(kind=dp), intent(out), dimension(size(A, 1), size(A, 2)) :: Q
real(kind=dp), intent(out), dimension(size(A, 1), size(A, 2)) :: R

Calls

proc~~qr_decomposition~~CallsGraph proc~qr_decomposition QR_decomposition proc~qr_givens_decomposition QR_Givens_decomposition proc~qr_decomposition->proc~qr_givens_decomposition proc~qr_gram_schmidt_classical_decomposition QR_Gram_Schmidt_Classical_decomposition proc~qr_decomposition->proc~qr_gram_schmidt_classical_decomposition proc~qr_gram_schmidt_modified_decomposition QR_Gram_Schmidt_Modified_decomposition proc~qr_decomposition->proc~qr_gram_schmidt_modified_decomposition proc~qr_householder_decomposition QR_Householder_decomposition proc~qr_decomposition->proc~qr_householder_decomposition proc~identity_n Identity_n proc~qr_givens_decomposition->proc~identity_n proc~rotation_matrix rotation_matrix proc~qr_givens_decomposition->proc~rotation_matrix proc~qr_householder_decomposition->proc~identity_n proc~rotation_matrix->proc~identity_n

Called by

proc~~qr_decomposition~~CalledByGraph proc~qr_decomposition QR_decomposition proc~eigen Eigen proc~eigen->proc~qr_decomposition proc~is_spd is_SPD proc~is_spd->proc~eigen

Source Code

    subroutine QR_decomposition(A, method, Q, R)

        real(dp), dimension(:, :), intent(in) :: A
        character(LEN=*), optional, intent(in) :: method
        real(dp), dimension(size(A, 1), size(A, 2)), intent(out) :: Q, R

        if (method == "QR_Householder") then
            call QR_Householder_decomposition(A, Q, R)
        else if (method == "QR_Givens") then
            call QR_Givens_decomposition(A, Q, R)
        else if (method == "QR_Gram_Schmidt_Classical") then
            call QR_Gram_Schmidt_Classical_decomposition(A, Q, R)
        else if (method == "QR_Gram_Schmidt_Modified") then
            call QR_Gram_Schmidt_Modified_decomposition(A, Q, R)
        end if

    end subroutine QR_decomposition