triangular_to_full Subroutine

public subroutine triangular_to_full(a, n, uplo)

Uses

  • proc~~triangular_to_full~~UsesGraph proc~triangular_to_full triangular_to_full module~messages messages proc~triangular_to_full->module~messages module~precision precision proc~triangular_to_full->module~precision module~messages->module~precision comm_IOFILE comm_IOFILE module~messages->comm_IOFILE comm_PAR comm_PAR module~messages->comm_PAR module~io_constants io_constants module~messages->module~io_constants iso_fortran_env iso_fortran_env module~precision->iso_fortran_env

@brief Fill the upper/lower triangle of the symmetric matrix @param[inout] a square NxN matrix in triangular form @param[in] n matrix dimension @param[in] uplo U if A is upper triangular, L if lower triangular

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(inout) :: a(n,*)
integer, intent(in) :: n
character(len=1), intent(in) :: uplo

Calls

proc~~triangular_to_full~~CallsGraph proc~triangular_to_full triangular_to_full interface~show_message show_message proc~triangular_to_full->interface~show_message

Called by

proc~~triangular_to_full~~CalledByGraph proc~triangular_to_full triangular_to_full proc~electric_moments electric_moments proc~electric_moments->proc~triangular_to_full proc~form_rohf_fock form_rohf_fock proc~form_rohf_fock->proc~triangular_to_full proc~oqp_tdhf_z_vector oqp_tdhf_z_vector proc~oqp_tdhf_z_vector->proc~triangular_to_full proc~tddft_fxc tddft_fxc proc~oqp_tdhf_z_vector->proc~tddft_fxc proc~tddft_gxc tddft_gxc proc~oqp_tdhf_z_vector->proc~tddft_gxc proc~tddft_fxc->proc~triangular_to_full proc~tddft_gxc->proc~triangular_to_full proc~utddft_fxc utddft_fxc proc~utddft_fxc->proc~triangular_to_full proc~scf_driver scf_driver proc~scf_driver->proc~form_rohf_fock proc~tdhf_energy tdhf_energy proc~tdhf_energy->proc~tddft_fxc proc~tdhf_z_vector_c tdhf_z_vector_C proc~tdhf_z_vector_c->proc~oqp_tdhf_z_vector proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver proc~tdhf_energy_c tdhf_energy_C proc~tdhf_energy_c->proc~tdhf_energy

Source Code

  subroutine triangular_to_full(a,n,uplo)
    use precision, only: dp
    use messages, only: show_message, with_abort

    implicit none

    real(kind=dp), intent(inout) :: a(n,*)
    integer, intent(in) :: n
    character(len=1), intent(in) :: uplo

    integer :: i

    if (uplo=='u' .or. uplo=='U') then
      do i = 1, n
        a(i+1:n,i) = a(i,i+1:n)
      end do
    else if (uplo=='l' .or. uplo=='L') then
      do i = 1, n
        a(i,i+1:n) = a(i+1:n,i)
      end do
    else
      call show_message('Invalid parameter UPLO='//uplo// &
              ' in `triangular_to_full`. Use either `L` or `U`.', with_abort)
    end if
  end subroutine triangular_to_full