mirror of
https://github.com/QuantumPackage/qp2.git
synced 2025-02-17 14:14:03 +01:00
135 lines
3.7 KiB
Fortran
135 lines
3.7 KiB
Fortran
BEGIN_PROVIDER [ logical, do_mo_cholesky ]
|
|
implicit none
|
|
BEGIN_DOC
|
|
! If True, use Cholesky vectors for MO integrals
|
|
END_DOC
|
|
do_mo_cholesky = do_ao_cholesky
|
|
! do_mo_cholesky = .False.
|
|
END_PROVIDER
|
|
|
|
BEGIN_PROVIDER [ integer, cholesky_mo_num ]
|
|
implicit none
|
|
BEGIN_DOC
|
|
! Number of Cholesky vectors in MO basis
|
|
END_DOC
|
|
integer, external :: getUnitAndOpen
|
|
integer :: iunit
|
|
if (read_mo_cholesky) then
|
|
iunit = getUnitAndOpen(trim(ezfio_work_dir)//'cholesky_mo_transp', 'R')
|
|
read(iunit) cholesky_mo_num
|
|
close(iunit)
|
|
else
|
|
cholesky_mo_num = cholesky_ao_num
|
|
endif
|
|
END_PROVIDER
|
|
|
|
BEGIN_PROVIDER [ double precision, cholesky_mo, (mo_num, mo_num, cholesky_mo_num) ]
|
|
implicit none
|
|
BEGIN_DOC
|
|
! Cholesky vectors in MO basis
|
|
END_DOC
|
|
|
|
integer :: k, i, j
|
|
|
|
call set_multiple_levels_omp(.False.)
|
|
!$OMP PARALLEL DO PRIVATE(k)
|
|
do k=1,cholesky_mo_num
|
|
do j=1,mo_num
|
|
do i=1,mo_num
|
|
cholesky_mo(i,j,k) = cholesky_mo_transp(k,i,j)
|
|
enddo
|
|
enddo
|
|
enddo
|
|
!$OMP END PARALLEL DO
|
|
|
|
END_PROVIDER
|
|
|
|
BEGIN_PROVIDER [ double precision, cholesky_mo_transp, (cholesky_mo_num, mo_num, mo_num) ]
|
|
implicit none
|
|
BEGIN_DOC
|
|
! Cholesky vectors in MO basis. Warning: it is transposed wrt cholesky_ao:
|
|
!
|
|
! - cholesky_ao is (ao_num^2 x cholesky_ao_num)
|
|
!
|
|
! - cholesky_mo_transp is (cholesky_mo_num x mo_num^2)
|
|
END_DOC
|
|
|
|
double precision, allocatable :: X(:,:,:)
|
|
double precision :: wall0, wall1
|
|
integer, external :: getUnitAndOpen
|
|
integer :: iunit, ierr, rank
|
|
|
|
if (read_mo_cholesky) then
|
|
print *, 'Reading Cholesky MO vectors from disk...'
|
|
iunit = getUnitAndOpen(trim(ezfio_work_dir)//'cholesky_mo_transp', 'R')
|
|
read(iunit) rank
|
|
if (cholesky_mo_num /= rank) then
|
|
stop 'inconsistent rank'
|
|
endif
|
|
read(iunit) cholesky_mo_transp
|
|
close(iunit)
|
|
else
|
|
print *, ''
|
|
print *, 'AO->MO Transformation of Cholesky vectors'
|
|
print *, '-----------------------------------------'
|
|
print *, ''
|
|
|
|
call wall_time(wall0)
|
|
|
|
allocate(X(mo_num,cholesky_mo_num,ao_num), stat=ierr)
|
|
if (ierr /= 0) then
|
|
print *, irp_here, ': Allocation failed'
|
|
endif
|
|
call dgemm('T','N', ao_num*cholesky_mo_num, mo_num, ao_num, 1.d0, &
|
|
cholesky_ao, ao_num, mo_coef, ao_num, 0.d0, X, ao_num*cholesky_mo_num)
|
|
call dgemm('T','N', cholesky_mo_num*mo_num, mo_num, ao_num, 1.d0, &
|
|
X, ao_num, mo_coef, ao_num, 0.d0, cholesky_mo_transp, cholesky_mo_num*mo_num)
|
|
deallocate(X)
|
|
call wall_time(wall1)
|
|
print*,'Time to provide MO cholesky vectors = ',(wall1-wall0)/60.d0, ' min'
|
|
|
|
|
|
if (write_mo_cholesky) then
|
|
print *, 'Writing Cholesky MO vectors to disk...'
|
|
iunit = getUnitAndOpen(trim(ezfio_work_dir)//'cholesky_mo_transp', 'W')
|
|
write(iunit) cholesky_mo_num
|
|
write(iunit) cholesky_mo_transp
|
|
close(iunit)
|
|
call ezfio_set_mo_two_e_ints_io_mo_cholesky('Read')
|
|
endif
|
|
endif
|
|
|
|
END_PROVIDER
|
|
|
|
|
|
BEGIN_PROVIDER [ double precision, cholesky_semi_mo_transp_simple, (cholesky_mo_num, ao_num, mo_num) ]
|
|
implicit none
|
|
BEGIN_DOC
|
|
! Cholesky vectors in MO basis
|
|
END_DOC
|
|
|
|
double precision, allocatable :: X(:,:,:)
|
|
double precision :: wall0, wall1
|
|
integer :: ierr
|
|
print *, 'Semi AO->MO Transformation of Cholesky vectors'
|
|
call wall_time(wall0)
|
|
|
|
allocate(X(mo_num,cholesky_mo_num,ao_num), stat=ierr)
|
|
if (ierr /= 0) then
|
|
print *, irp_here, ': Allocation failed'
|
|
endif
|
|
integer :: i_chol, i_mo, j_mo, i_ao
|
|
cholesky_semi_mo_transp_simple = 0.d0
|
|
do i_mo = 1, mo_num
|
|
do i_ao = 1, ao_num
|
|
do j_mo = 1, mo_num
|
|
do i_chol = 1, cholesky_mo_num
|
|
cholesky_semi_mo_transp_simple(i_chol, i_ao,i_mo) += cholesky_mo_transp(i_chol,j_mo,i_mo) * mo_coef_transp(j_mo,i_ao)
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
|
|
END_PROVIDER
|
|
|