Incomplete_Cholesky_decomposition Subroutine

public subroutine Incomplete_Cholesky_decomposition(A, L, level)

Incomplete Cholesky decomposition of a matrix A This subroutine performs incomplete Cholesky decomposition of a given matrix A, where L is a lower triangular matrix and U is an upper triangular matrix.

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(in), dimension(:, :) :: A
real(kind=dp), intent(out), dimension(size(A, 1), size(A, 1)) :: L
integer, intent(in), optional :: level

Called by

proc~~incomplete_cholesky_decomposition~~CalledByGraph proc~incomplete_cholesky_decomposition Incomplete_Cholesky_decomposition proc~calculate_icf_preconditioner Calculate_ICF_preconditioner proc~calculate_icf_preconditioner->proc~incomplete_cholesky_decomposition

Source Code

    subroutine Incomplete_Cholesky_decomposition(A, L, level)

        real(dp), dimension(:, :), intent(in) :: A
        integer, optional, intent(in) :: level
        real(dp), dimension(size(A, 1), size(A, 1)), intent(out) :: L
        logical, dimension(size(A, 1), size(A, 1)) :: S
        integer :: N, i, j
        integer, dimension(size(A, 1), size(A, 1)) :: fill_level

        N = size(A, 1)

        L = 0.d0

        if (present(level)) then
            call compute_fill_pattern_IC(A, fill_level, level, N)
            S = (fill_level <= level)
        else
            S = A /= 0
        end if

        do i = 1, N
            do j = 1, i - 1
                if (S(i, j)) L(i, j) = (A(i, j) - &
                                        dot_product(L(i, 1:j - 1), L(j, 1:j - 1))) / &
                                       L(j, j)
            end do

            if (S(i, i)) L(i, i) = sqrt(A(i, i) - &
                                        dot_product(L(i, 1:i - 1), L(i, 1:i - 1)))
        end do

    end subroutine Incomplete_Cholesky_decomposition