schmd Subroutine

public subroutine schmd(v, m, n, ldv, x)

Uses

  • proc~~schmd~~UsesGraph proc~schmd schmd iso_fortran_env iso_fortran_env proc~schmd->iso_fortran_env module~messages messages proc~schmd->module~messages 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 module~precision precision module~messages->module~precision module~precision->iso_fortran_env

Arguments

Type IntentOptional Attributes Name
real(kind=real64), intent(inout) :: v(ldv,n)
integer, intent(in) :: m
integer, intent(in) :: n
integer, intent(in) :: ldv
real(kind=real64), intent(inout) :: x(n)

Calls

proc~~schmd~~CallsGraph proc~schmd schmd interface~show_message show_message proc~schmd->interface~show_message proc~oqp_dgeqrf_i64 oqp_dgeqrf_i64 proc~schmd->proc~oqp_dgeqrf_i64 proc~oqp_dorgqr_i64 oqp_dorgqr_i64 proc~schmd->proc~oqp_dorgqr_i64 proc~oqp_dgeqrf_i64->interface~show_message dgeqrf dgeqrf proc~oqp_dgeqrf_i64->dgeqrf proc~oqp_dorgqr_i64->interface~show_message dorgqr dorgqr proc~oqp_dorgqr_i64->dorgqr

Called by

proc~~schmd~~CalledByGraph proc~schmd schmd proc~corresponding_orbital_projection corresponding_orbital_projection proc~corresponding_orbital_projection->proc~schmd proc~huckel_guess huckel_guess proc~huckel_guess->proc~corresponding_orbital_projection proc~guess_huckel guess_huckel proc~guess_huckel->proc~huckel_guess proc~guess_huckel_c guess_huckel_C proc~guess_huckel_c->proc~guess_huckel

Source Code

  subroutine schmd(v, m, n, ldv, x)
    use, intrinsic :: iso_fortran_env, only: real64
    use messages, only: show_message, WITH_ABORT
    implicit none

    integer, intent(IN) :: ldv, m, n
    real(real64), intent(INOUT) :: v(ldv, n), x(n)

    real(real64), allocatable :: work(:)
    integer :: lwork
    integer :: info
    real(real64) :: wrksize(1)

    if (M > N) then
      call show_message("SCHMD: M > N", WITH_ABORT)
    end if
    if (N > LDV) then
      call show_message("SCHMD: N > LDV", WITH_ABORT)
    end if
!   Householder QR-based version using LAPACK:
    call dgeqrf(n, m, v, ldv, x, wrksize, -1, info)
    lwork = max(int(wrksize(1)), n)
    allocate (work(lwork))
    call dgeqrf(n, m, v, ldv, x, work, lwork, info)
    call dorgqr(n, n, m, v, ldv, x, work, lwork, info)
  end subroutine schmd