sfdmat Subroutine

public subroutine sfdmat(bvec, abxc, mo_a, ta, tb, noca, nocb)

Uses

  • proc~~sfdmat~~UsesGraph proc~sfdmat sfdmat module~mathlib mathlib proc~sfdmat->module~mathlib module~precision precision proc~sfdmat->module~precision module~tdhf_lib tdhf_lib proc~sfdmat->module~tdhf_lib module~mathlib->module~precision module~oqp_linalg oqp_linalg module~mathlib->module~oqp_linalg iso_fortran_env iso_fortran_env module~precision->iso_fortran_env module~tdhf_lib->module~precision module~basis_tools basis_tools module~tdhf_lib->module~basis_tools module~int2_compute int2_compute module~tdhf_lib->module~int2_compute module~tdhf_lib->module~oqp_linalg module~basis_tools->module~precision module~basis_tools->iso_fortran_env module~atomic_structure_m atomic_structure_m module~basis_tools->module~atomic_structure_m module~constants constants module~basis_tools->module~constants module~io_constants io_constants module~basis_tools->module~io_constants module~parallel parallel module~basis_tools->module~parallel module~int2_compute->module~precision module~int2_compute->module~basis_tools module~int2_compute->module~atomic_structure_m module~int2_pairs int2_pairs module~int2_compute->module~int2_pairs module~int2e_libint int2e_libint module~int2_compute->module~int2e_libint module~int2e_rys int2e_rys module~int2_compute->module~int2e_rys module~messages messages module~int2_compute->module~messages module~int2_compute->module~parallel module~blas_wrap blas_wrap module~oqp_linalg->module~blas_wrap module~lapack_wrap lapack_wrap module~oqp_linalg->module~lapack_wrap iso_c_binding iso_c_binding module~atomic_structure_m->iso_c_binding module~blas_wrap->module~precision module~blas_wrap->module~messages module~mathlib_types mathlib_types module~blas_wrap->module~mathlib_types module~constants->module~precision module~int2_pairs->module~precision module~int2e_libint->module~precision module~int2e_libint->module~constants module~int2e_libint->module~int2_pairs module~int2e_libint->iso_c_binding module~libint_f libint_f module~int2e_libint->module~libint_f module~int2e_rys->module~precision module~int2e_rys->module~basis_tools module~int2e_rys->module~constants module~lapack_wrap->module~precision module~lapack_wrap->module~messages module~lapack_wrap->module~mathlib_types module~messages->module~precision module~messages->module~io_constants comm_IOFILE comm_IOFILE module~messages->comm_IOFILE comm_PAR comm_PAR module~messages->comm_PAR module~parallel->module~precision module~parallel->iso_fortran_env module~parallel->iso_c_binding mpi mpi module~parallel->mpi module~libint_f->iso_c_binding

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(in), dimension(:) :: bvec
real(kind=dp), intent(inout), dimension(:,:) :: abxc
real(kind=dp), intent(in), dimension(:,:) :: mo_a
real(kind=dp), intent(out), dimension(:) :: ta
real(kind=dp), intent(out), dimension(:) :: tb
integer, intent(in) :: noca
integer, intent(in) :: nocb

Calls

proc~~sfdmat~~CallsGraph proc~sfdmat sfdmat interface~pack_matrix pack_matrix proc~sfdmat->interface~pack_matrix proc~iatogen iatogen proc~sfdmat->proc~iatogen proc~oqp_dgemm_i64 oqp_dgemm_i64 proc~sfdmat->proc~oqp_dgemm_i64 proc~orthogonal_transform orthogonal_transform proc~sfdmat->proc~orthogonal_transform proc~pack_f90 PACK_F90 interface~pack_matrix->proc~pack_f90 dgemm dgemm proc~oqp_dgemm_i64->dgemm interface~show_message show_message proc~oqp_dgemm_i64->interface~show_message proc~orthogonal_transform->proc~oqp_dgemm_i64 proc~orthogonal_transform->interface~show_message proc~pack_f90->interface~show_message proc~oqp_dtrttp_i64 oqp_dtrttp_i64 proc~pack_f90->proc~oqp_dtrttp_i64 proc~oqp_dtrttp_i64->interface~show_message dtrttp dtrttp proc~oqp_dtrttp_i64->dtrttp

Called by

proc~~sfdmat~~CalledByGraph proc~sfdmat sfdmat proc~tdhf_sf_energy tdhf_sf_energy proc~tdhf_sf_energy->proc~sfdmat proc~tdhf_sf_energy_c tdhf_sf_energy_C proc~tdhf_sf_energy_c->proc~tdhf_sf_energy

Source Code

  subroutine sfdmat(bvec,abxc,mo_a,ta,tb, &
                    noca,nocb)
    use precision, only: dp
    use tdhf_lib, only: iatogen
    use mathlib, only: pack_matrix
    use mathlib, only: orthogonal_transform

    implicit none

    real(kind=dp), intent(in), dimension(:) :: bvec
    real(kind=dp), intent(in), dimension(:,:) :: mo_a
    real(kind=dp), intent(inout), dimension(:,:) :: abxc
    real(kind=dp), intent(out), dimension(:) :: ta, tb
    integer, intent(in) :: noca, nocb

    integer :: nvirb, nbf, nbf_tri, xvec_dim
    real(kind=dp), allocatable, dimension(:,:) :: scr1, scr2

    nbf = ubound(mo_a,1)
    nbf_tri = ubound(ta,1)
    xvec_dim = ubound(bvec,1)
    allocate(scr1(nbf,nbf), &
             scr2(nbf,nbf), &
             source=0.0_dp)

  ! MO(I+,A-) -> AO(M,N)
    nvirb = nbf-nocb

    call iatogen(bvec,scr1,noca,nocb)
    call orthogonal_transform('t', nbf, mo_a, scr1, abxc, scr2)

  ! Unrelaxed difference density matrix -----

  ! OCC(Alpha)-OCC(Alpha)
    call dgemm('n','t',noca,noca,nvirb, &
              -1.0_dp,bvec,noca, &
                      bvec,noca, &
               0.0_dp,scr1,noca)

  ! MO(I+,J+) -> AO(M,N)
    call dgemm('n','n',nbf,noca,noca, &
               1.0_dp,mo_a,nbf, &
                      scr1,noca, &
               0.0_dp,scr2,nbf)
    call dgemm('n','t',nbf,nbf,noca, &
               1.0_dp,scr2,nbf, &
                      mo_a,nbf, &
               0.0_dp,scr1,nbf)
    call pack_matrix(scr1,ta)

    call dgemm('t','n',nvirb,nvirb,noca, &
               1.0_dp,bvec,noca, &
                      bvec,noca, &
               0.0_dp,scr1,nvirb)

  ! MO(A-,B-) -> AO(M,N)
    call dgemm('n','n',nbf,nvirb,nvirb, &
               1.0_dp,mo_a(:,nocb+1:),nbf, &
                      scr1,nvirb, &
               0.0_dp,scr2,nbf)
    call dgemm('n','t',nbf,nbf,nvirb, &
               1.0_dp,scr2,nbf, &
                      mo_a(:,nocb+1:),nbf, &
               0.0_dp,scr1,nbf)
    call pack_matrix(scr1,tb)

    deallocate(scr1,scr2)
  end subroutine sfdmat