mo_reorder Subroutine

public subroutine mo_reorder(infos, Va, Ea, Vb, Eb, Sq)

Uses

  • proc~~mo_reorder~~UsesGraph proc~mo_reorder mo_reorder module~io_constants io_constants proc~mo_reorder->module~io_constants module~precision precision proc~mo_reorder->module~precision module~types types proc~mo_reorder->module~types iso_fortran_env iso_fortran_env module~precision->iso_fortran_env module~types->module~precision iso_c_binding iso_c_binding module~types->iso_c_binding module~atomic_structure_m atomic_structure_m module~types->module~atomic_structure_m module~basis_tools basis_tools module~types->module~basis_tools module~functionals functionals module~types->module~functionals module~parallel parallel module~types->module~parallel tagarray tagarray module~types->tagarray module~atomic_structure_m->iso_c_binding module~basis_tools->module~io_constants module~basis_tools->module~precision module~basis_tools->iso_fortran_env module~basis_tools->module~atomic_structure_m module~basis_tools->module~parallel module~constants constants module~basis_tools->module~constants module~functionals->module~precision module~functionals->iso_c_binding xc_f03_lib_m xc_f03_lib_m module~functionals->xc_f03_lib_m module~parallel->module~precision module~parallel->iso_c_binding module~parallel->iso_fortran_env mpi mpi module~parallel->mpi module~constants->module~precision

Arguments

Type IntentOptional Attributes Name
type(information), intent(inout) :: infos
real(kind=dp), intent(inout), dimension(:,:) :: Va
real(kind=dp), intent(inout), dimension(:) :: Ea
real(kind=dp), intent(inout), dimension(:,:) :: Vb
real(kind=dp), intent(inout), dimension(:) :: Eb
real(kind=dp), intent(in), dimension(:,:) :: Sq

Calls

proc~~mo_reorder~~CallsGraph proc~mo_reorder mo_reorder dgemm dgemm proc~mo_reorder->dgemm proc~reordermos reordermos proc~mo_reorder->proc~reordermos dswap dswap proc~reordermos->dswap

Called by

proc~~mo_reorder~~CalledByGraph proc~mo_reorder mo_reorder proc~scf_driver scf_driver proc~scf_driver->proc~mo_reorder proc~hf_energy hf_energy proc~hf_energy->proc~scf_driver

Source Code

 subroutine mo_reorder(infos, Va, Ea, Vb, Eb, Sq)
!  Va, Ea: MO and MO energies of Previous Point
!  Vb, Eb: MO and MO energies of Current Point
!  The Vb and Eb are reordered based on the orders of Va and Ea.
!  Sq: Square Form of Overlap Matrix in AO
   use precision, only: dp
   use io_constants, only: iw
   use types, only: information

   implicit none

   type(information), intent(inout) :: infos
   real(kind=dp), intent(inout), dimension(:,:) :: Va, Vb ! (nbf, nbf)
   real(kind=dp), intent(inout), dimension(:) :: Ea, Eb    ! (nbf)
   real(kind=dp), intent(in), dimension(:,:) :: Sq    ! (nbf, nbf)

   integer :: i, loc
   real(kind=dp) :: tmp, tmp2, tmp_abs
   real(kind=dp), allocatable, dimension(:,:) :: scr, smo
   logical, allocatable, dimension(:) :: pr
   integer, allocatable, dimension(:) :: locs
   integer :: na, nbf

   write(iw, fmt='(/," MOM is on: Reordering the New MO and MO energies....")')

   nbf = infos%basis%nbf
   na = infos%mol_prop%nelec_a

   allocate(scr(nbf,nbf), &
            smo(nbf,nbf), &
            source=0.0_dp)
   allocate(pr(nbf), source=.false.)
   allocate(locs(nbf), source=0)

   call dgemm('t', 'n', nbf, nbf, nbf, &
               1.0_dp, Va, nbf, &
                       Sq, nbf, &
               0.0_dp, scr, nbf)
   call dgemm('n', 'n', nbf, nbf, nbf, &
               1.0_dp, scr, nbf, &
                       Vb, nbf, &
               0.0_dp, smo, nbf)
!  Now Smo is overlap matrix in MO, row and column corresponds to Old and New MOs.
   do i = 1, nbf
     smo(:,i) = smo(:,i) / norm2(smo(:,i))
   end do
!  Finding out the location of maximum value in column i
   do i = 1, nbf
     loc = maxloc(abs(smo(:nbf,i)), dim=1)
     locs(i) = loc
     pr(loc) = .true.
   enddo
!  Checking printouts
   if (.not.all(pr)) then
     write(iw, fmt='(/,"   Warning")')
     write(iw, fmt='(" Some orbitals are missing in reordering")')
     write(iw, advance='no', fmt='(" Their indices are:")')
     do i = 1, nbf
       if (.not.pr(i)) write(iw, advance='no', fmt='(I4,",")') i
     end do
     write(iw,*)
     write(iw, fmt='(/,"   Error. Stop")')
     write(iw, fmt='(" Some orbitals are missing in reordering")')
     write(iw,*)
   end if

   write(iw,fmt='(x,a)') 'Old MO Index ---> New MO Index'
   do i = 1, nbf
     tmp_abs = maxval(abs(smo(:nbf,i)), dim=1)
     loc = maxloc(abs(smo(:nbf,i)), dim=1)
     tmp = smo(loc,i)
     tmp2 = smo(i,i)

     if ((loc.ne.i).and.(i.le.na+1)) then
       write(iw, advance='no', &
        fmt='(5x,i4,14x,i4)') i, loc
       if (i == na-1) write(iw, advance='no',fmt='(2x,a)') ' HOMO'
       if (i == na) write(iw, advance='no',fmt='(2x,a)') ' LUMO'
       if (i /= loc .and. tmp_abs < 0.9d+00) then
         write(iw, fmt='(2x,a)') &
         ' rearranged, WARNING'
       elseif (i == loc .and. tmp_abs < 0.9d+00) then
         write(iw, fmt='(2x,a)')  &
         ' WARNING'
       elseif (i /= loc .and. tmp_abs > 0.9d+00) then
         write(iw, fmt='(2x,a)') &
         ' rearranged'
       else
         write(iw,*)
       end if
     endif
   enddo

!  Reordering MO and MO energies of Current point.
   call reorderMOs(Vb, Eb, Smo, nbf, nbf, 1, na+1)

 end subroutine mo_reorder