mirror of
https://github.com/QuantumPackage/qp2.git
synced 2025-04-25 17:54:44 +02:00
Accelerated PT2-cholesky with dgemm
This commit is contained in:
parent
eb6e1d4339
commit
70ad9f31b3
2
external/ezfio
vendored
2
external/ezfio
vendored
@ -1 +1 @@
|
||||
Subproject commit d02132ea79217c16fd24242e8f8b8a6c3ff68091
|
||||
Subproject commit dba01c4fe0ff7b84c5ecfb1c7c77ec68781311b3
|
||||
@ -1478,17 +1478,21 @@ subroutine get_d0(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)
|
||||
|
||||
integer, parameter :: bant=1
|
||||
|
||||
double precision, allocatable :: hij_cache1(:), hij_cache2(:)
|
||||
allocate (hij_cache1(mo_num),hij_cache2(mo_num))
|
||||
! double precision, allocatable :: hij_cache1(:), hij_cache2(:)
|
||||
double precision, allocatable :: hij_cache1(:,:), hij_cache2(:,:)
|
||||
! allocate (hij_cache1(mo_num),hij_cache2(mo_num))
|
||||
|
||||
PROVIDE mo_integrals_threshold
|
||||
|
||||
if(sp == 3) then ! AB
|
||||
|
||||
allocate(hij_cache1(mo_num,mo_num))
|
||||
|
||||
h1 = p(1,1)
|
||||
h2 = p(1,2)
|
||||
call get_mo_two_e_integrals_ij(h2,h1,mo_num,hij_cache1,mo_integrals_map)
|
||||
do p1=1, mo_num
|
||||
if(bannedOrb(p1, 1)) cycle
|
||||
call get_mo_two_e_integrals(p1,h2,h1,mo_num,hij_cache1,mo_integrals_map)
|
||||
do p2=1, mo_num
|
||||
if(bannedOrb(p2,2)) cycle
|
||||
if(banned(p1, p2, bant)) cycle ! rentable?
|
||||
@ -1497,7 +1501,7 @@ subroutine get_d0(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)
|
||||
call i_h_j(gen, det, N_int, hij)
|
||||
else
|
||||
phase = get_phase_bi(phasemask, 1, 2, h1, p1, h2, p2, N_int)
|
||||
hij = hij_cache1(p2) * phase
|
||||
hij = hij_cache1(p2,p1) * phase
|
||||
end if
|
||||
if (dabs(hij) < mo_integrals_threshold) cycle
|
||||
!DIR$ LOOP COUNT AVG(4)
|
||||
@ -1507,13 +1511,18 @@ subroutine get_d0(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)
|
||||
end do
|
||||
end do
|
||||
|
||||
deallocate(hij_cache1)
|
||||
|
||||
else ! AA BB
|
||||
|
||||
allocate(hij_cache1(mo_num,mo_num),hij_cache2(mo_num,mo_num))
|
||||
|
||||
p1 = p(1,sp)
|
||||
p2 = p(2,sp)
|
||||
call get_mo_two_e_integrals_ij(p2,p1,mo_num,hij_cache1,mo_integrals_map)
|
||||
call get_mo_two_e_integrals_ij(p1,p2,mo_num,hij_cache2,mo_integrals_map)
|
||||
do puti=1, mo_num
|
||||
if (bannedOrb(puti, sp)) cycle
|
||||
call get_mo_two_e_integrals(puti,p2,p1,mo_num,hij_cache1,mo_integrals_map)
|
||||
call get_mo_two_e_integrals(puti,p1,p2,mo_num,hij_cache2,mo_integrals_map)
|
||||
do putj=puti+1, mo_num
|
||||
if(bannedOrb(putj, sp)) cycle
|
||||
if(banned(puti, putj, bant)) cycle ! rentable?
|
||||
@ -1522,7 +1531,7 @@ subroutine get_d0(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)
|
||||
call i_h_j(gen, det, N_int, hij)
|
||||
if (dabs(hij) < mo_integrals_threshold) cycle
|
||||
else
|
||||
hij = hij_cache1(putj) - hij_cache2(putj)
|
||||
hij = hij_cache1(putj,puti) - hij_cache2(putj,puti)
|
||||
if (dabs(hij) < mo_integrals_threshold) cycle
|
||||
hij = hij * get_phase_bi(phasemask, sp, sp, puti, p1 , putj, p2, N_int)
|
||||
end if
|
||||
@ -1532,9 +1541,12 @@ subroutine get_d0(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)
|
||||
enddo
|
||||
end do
|
||||
end do
|
||||
end if
|
||||
|
||||
deallocate(hij_cache1,hij_cache2)
|
||||
|
||||
end if
|
||||
|
||||
! deallocate(hij_cache1,hij_cache2)
|
||||
end
|
||||
|
||||
|
||||
|
||||
@ -262,8 +262,14 @@ subroutine get_mo_two_e_integrals(j,k,l,sze,out_val,map)
|
||||
out_val, 1)
|
||||
else
|
||||
integer :: isplit
|
||||
out_val = 0.d0
|
||||
do isplit=1,4
|
||||
call sgemv('T', cholesky_mo_num_split(isplit+1) - cholesky_mo_num_split(isplit), &
|
||||
mo_integrals_cache_min-1, 1., &
|
||||
cholesky_mo_transp_sp(cholesky_mo_num_split(isplit),1,k), cholesky_mo_num, &
|
||||
cholesky_mo_transp_sp(cholesky_mo_num_split(isplit),j,l), 1, 0., &
|
||||
out_val_sp, 1)
|
||||
out_val(1:mo_integrals_cache_min-1) += out_val_sp(1:mo_integrals_cache_min-1)
|
||||
out_val(1:mo_integrals_cache_min-1) = out_val_sp(1:mo_integrals_cache_min-1)
|
||||
do isplit=2,4
|
||||
call sgemv('T', cholesky_mo_num_split(isplit+1) - cholesky_mo_num_split(isplit), &
|
||||
mo_integrals_cache_min-1, 1., &
|
||||
cholesky_mo_transp_sp(cholesky_mo_num_split(isplit),1,k), cholesky_mo_num, &
|
||||
@ -442,16 +448,36 @@ subroutine get_mo_two_e_integrals_ij(k,l,sze,out_array,map)
|
||||
double precision, intent(out) :: out_array(sze,sze)
|
||||
type(map_type), intent(inout) :: map
|
||||
integer :: j
|
||||
real(integral_kind), allocatable :: tmp_val(:)
|
||||
|
||||
if ( (mo_integrals_cache_min>1).or.(mo_integrals_cache_max<mo_num) ) then
|
||||
|
||||
if (do_mo_cholesky) then
|
||||
|
||||
if (mo_cholesky_double) then
|
||||
call dgemm('T', 'N', mo_num, mo_num, cholesky_mo_num, 1.d0, &
|
||||
cholesky_mo_transp(1,1,k), cholesky_mo_num, &
|
||||
cholesky_mo_transp(1,1,l), cholesky_mo_num, 0.d0, &
|
||||
out_array, sze)
|
||||
else
|
||||
integer :: isplit
|
||||
double precision, allocatable :: out_array_sp(:,:)
|
||||
allocate(out_array_sp(sze,sze))
|
||||
call sgemm('T', 'N', mo_num, mo_num, &
|
||||
cholesky_mo_num_split(isplit+1) - cholesky_mo_num_split(isplit), 1.d0, &
|
||||
cholesky_mo_transp_sp(cholesky_mo_num_split(isplit),1,k), cholesky_mo_num, &
|
||||
cholesky_mo_transp_sp(cholesky_mo_num_split(isplit),1,l), cholesky_mo_num, 0.d0, &
|
||||
out_array_sp, sze)
|
||||
out_array(1:sze,1:sze) = out_array_sp(1:sze,1:sze)
|
||||
do isplit=2,4
|
||||
call sgemm('T', 'N', mo_num, mo_num, &
|
||||
cholesky_mo_num_split(isplit+1) - cholesky_mo_num_split(isplit), 1.d0, &
|
||||
cholesky_mo_transp_sp(cholesky_mo_num_split(isplit),1,k), cholesky_mo_num, &
|
||||
cholesky_mo_transp_sp(cholesky_mo_num_split(isplit),1,l), cholesky_mo_num, 0.d0, &
|
||||
out_array_sp, sze)
|
||||
out_array(1:sze,1:sze) = out_array(1:sze,1:sze) + out_array_sp(1:sze,1:sze)
|
||||
enddo
|
||||
deallocate(out_array_sp)
|
||||
endif
|
||||
|
||||
else
|
||||
|
||||
|
||||
Loading…
x
Reference in New Issue
Block a user