Increased integrals cache to 128 MOs

This commit is contained in:
Anthony Scemama 2017-08-20 12:27:08 +02:00
parent d1c9f179bc
commit f2cb73f006
6 changed files with 31 additions and 137 deletions

View File

@ -543,7 +543,7 @@ def create_ocaml_input(dict_ezfio_cfg, module_lower):
template += ["open Qptypes;;",
"open Qputils;;",
"open Core.Std;;",
"open Core;;",
"",
"module {0} : sig".format(module_lower.capitalize())]

View File

@ -352,7 +352,7 @@ class EZFIO_ocaml(object):
l_template = ['open Qputils;;',
'open Qptypes;;',
'open Core.Std;;',
'open Core;;',
'']
for m in self.l_module_lower:

View File

@ -9,14 +9,11 @@ BEGIN_PROVIDER [ integer, ao_num_align ]
ao_num_align = align_double(ao_num)
END_PROVIDER
BEGIN_PROVIDER [ integer, ao_prim_num_max ]
&BEGIN_PROVIDER [ integer, ao_prim_num_max_align ]
BEGIN_PROVIDER [ integer, ao_prim_num_max ]
implicit none
ao_prim_num_max = 0
PROVIDE ezfio_filename
call ezfio_get_ao_basis_ao_prim_num_max(ao_prim_num_max)
integer :: align_double
ao_prim_num_max_align = align_double(ao_prim_num_max)
END_PROVIDER
BEGIN_PROVIDER [ double precision, ao_coef_normalized, (ao_num_align,ao_prim_num_max) ]

View File

@ -97,7 +97,7 @@ type: double precision
size: (determinants.n_det)
[expected_s2]
interface: ezfio,provider
interface: ezfio
doc: Expected value of S^2
type: double precision

View File

@ -133,115 +133,6 @@ BEGIN_PROVIDER [ integer(bit_kind), psi_det, (N_int,2,psi_det_size) ]
END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), psi_occ_pattern, (N_int,2,psi_det_size) ]
&BEGIN_PROVIDER [ integer, N_occ_pattern ]
implicit none
BEGIN_DOC
! array of the occ_pattern present in the wf
! psi_occ_pattern(:,1,j) = jth occ_pattern of the wave function : represent all the single occupation
! psi_occ_pattern(:,2,j) = jth occ_pattern of the wave function : represent all the double occupation
END_DOC
integer :: i,j,k
! create
do i = 1, N_det
do k = 1, N_int
psi_occ_pattern(k,1,i) = ieor(psi_det(k,1,i),psi_det(k,2,i))
psi_occ_pattern(k,2,i) = iand(psi_det(k,1,i),psi_det(k,2,i))
enddo
enddo
! Sort
integer, allocatable :: iorder(:)
integer*8, allocatable :: bit_tmp(:)
integer*8, external :: occ_pattern_search_key
integer(bit_kind), allocatable :: tmp_array(:,:,:)
logical,allocatable :: duplicate(:)
allocate ( iorder(N_det), duplicate(N_det), bit_tmp(N_det), tmp_array(N_int,2,psi_det_size) )
do i=1,N_det
iorder(i) = i
!$DIR FORCEINLINE
bit_tmp(i) = occ_pattern_search_key(psi_occ_pattern(1,1,i),N_int)
enddo
call i8sort(bit_tmp,iorder,N_det)
!DIR$ IVDEP
do i=1,N_det
do k=1,N_int
tmp_array(k,1,i) = psi_occ_pattern(k,1,iorder(i))
tmp_array(k,2,i) = psi_occ_pattern(k,2,iorder(i))
enddo
duplicate(i) = .False.
enddo
i=1
integer (bit_kind) :: occ_pattern_tmp
do i=1,N_det
duplicate(i) = .False.
enddo
do i=1,N_det-1
if (duplicate(i)) then
cycle
endif
j = i+1
do while (bit_tmp(j)==bit_tmp(i))
if (duplicate(j)) then
j+=1
cycle
endif
duplicate(j) = .True.
do k=1,N_int
if ( (tmp_array(k,1,i) /= tmp_array(k,1,j)) &
.or. (tmp_array(k,2,i) /= tmp_array(k,2,j)) ) then
duplicate(j) = .False.
exit
endif
enddo
j+=1
if (j>N_det) then
exit
endif
enddo
enddo
N_occ_pattern=0
do i=1,N_det
if (duplicate(i)) then
cycle
endif
N_occ_pattern += 1
do k=1,N_int
psi_occ_pattern(k,1,N_occ_pattern) = tmp_array(k,1,i)
psi_occ_pattern(k,2,N_occ_pattern) = tmp_array(k,2,i)
enddo
enddo
deallocate(iorder,duplicate,bit_tmp,tmp_array)
! !TODO DEBUG
! integer :: s
! do i=1,N_occ_pattern
! do j=i+1,N_occ_pattern
! s = 0
! do k=1,N_int
! if((psi_occ_pattern(k,1,j) /= psi_occ_pattern(k,1,i)).or. &
! (psi_occ_pattern(k,2,j) /= psi_occ_pattern(k,2,i))) then
! s=1
! exit
! endif
! enddo
! if ( s == 0 ) then
! print *, 'Error : occ ', j, 'already in wf'
! call debug_det(psi_occ_pattern(1,1,j),N_int)
! stop
! endif
! enddo
! enddo
! !TODO DEBUG
END_PROVIDER
BEGIN_PROVIDER [ double precision, psi_coef, (psi_det_size,N_states) ]
implicit none

View File

@ -327,41 +327,45 @@ subroutine insert_into_mo_integrals_map(n_integrals, &
call map_update(mo_integrals_map, buffer_i, buffer_values, n_integrals, thr)
end
BEGIN_PROVIDER [ integer, mo_integrals_cache_min ]
&BEGIN_PROVIDER [ integer, mo_integrals_cache_max ]
BEGIN_PROVIDER [ integer*4, mo_integrals_cache_min ]
&BEGIN_PROVIDER [ integer*4, mo_integrals_cache_max ]
&BEGIN_PROVIDER [ integer*8, mo_integrals_cache_min_8 ]
&BEGIN_PROVIDER [ integer*8, mo_integrals_cache_max_8 ]
implicit none
BEGIN_DOC
! Min and max values of the MOs for which the integrals are in the cache
END_DOC
mo_integrals_cache_min = max(1,elec_alpha_num - 31)
mo_integrals_cache_max = min(mo_tot_num,mo_integrals_cache_min+63)
mo_integrals_cache_min_8 = max(1_8,elec_alpha_num - 63_8)
mo_integrals_cache_max_8 = min(int(mo_tot_num,8),mo_integrals_cache_min+127_8)
mo_integrals_cache_min = mo_integrals_cache_min_8
mo_integrals_cache_max = mo_integrals_cache_max_8
END_PROVIDER
BEGIN_PROVIDER [ double precision, mo_integrals_cache, (0:64*64*64*64) ]
BEGIN_PROVIDER [ double precision, mo_integrals_cache, (0_8:128_8*128_8*128_8*128_8) ]
implicit none
BEGIN_DOC
! Cache of MO integrals for fast access
END_DOC
PROVIDE mo_bielec_integrals_in_map
integer :: i,j,k,l
integer :: ii
integer*8 :: i,j,k,l
integer*8 :: ii
integer(key_kind) :: idx
real(integral_kind) :: integral
FREE ao_integrals_cache
!$OMP PARALLEL DO PRIVATE (i,j,k,l,idx,ii,integral)
do l=mo_integrals_cache_min,mo_integrals_cache_max
do k=mo_integrals_cache_min,mo_integrals_cache_max
do j=mo_integrals_cache_min,mo_integrals_cache_max
do i=mo_integrals_cache_min,mo_integrals_cache_max
do l=mo_integrals_cache_min_8,mo_integrals_cache_max_8
do k=mo_integrals_cache_min_8,mo_integrals_cache_max_8
do j=mo_integrals_cache_min_8,mo_integrals_cache_max_8
do i=mo_integrals_cache_min_8,mo_integrals_cache_max_8
!DIR$ FORCEINLINE
call bielec_integrals_index(i,j,k,l,idx)
!DIR$ FORCEINLINE
call map_get(mo_integrals_map,idx,integral)
ii = l-mo_integrals_cache_min
ii = ior( ishft(ii,6), k-mo_integrals_cache_min)
ii = ior( ishft(ii,6), j-mo_integrals_cache_min)
ii = ior( ishft(ii,6), i-mo_integrals_cache_min)
ii = l-mo_integrals_cache_min_8
ii = ior( ishft(ii,7), k-mo_integrals_cache_min_8)
ii = ior( ishft(ii,7), j-mo_integrals_cache_min_8)
ii = ior( ishft(ii,7), i-mo_integrals_cache_min_8)
mo_integrals_cache(ii) = integral
enddo
enddo
@ -381,6 +385,7 @@ double precision function get_mo_bielec_integral(i,j,k,l,map)
integer, intent(in) :: i,j,k,l
integer(key_kind) :: idx
integer :: ii
integer*8 :: ii_8
type(map_type), intent(inout) :: map
real(integral_kind) :: tmp
PROVIDE mo_bielec_integrals_in_map mo_integrals_cache
@ -388,18 +393,19 @@ double precision function get_mo_bielec_integral(i,j,k,l,map)
ii = ior(ii, k-mo_integrals_cache_min)
ii = ior(ii, j-mo_integrals_cache_min)
ii = ior(ii, i-mo_integrals_cache_min)
if (iand(ii, -64) /= 0) then
if (iand(ii, -128) /= 0) then
! if (.True.) then
!DIR$ FORCEINLINE
call bielec_integrals_index(i,j,k,l,idx)
!DIR$ FORCEINLINE
call map_get(map,idx,tmp)
get_mo_bielec_integral = dble(tmp)
else
ii = l-mo_integrals_cache_min
ii = ior( ishft(ii,6), k-mo_integrals_cache_min)
ii = ior( ishft(ii,6), j-mo_integrals_cache_min)
ii = ior( ishft(ii,6), i-mo_integrals_cache_min)
get_mo_bielec_integral = mo_integrals_cache(ii)
ii_8 = l-mo_integrals_cache_min_8
ii_8 = ior( ishft(ii_8,7), k-mo_integrals_cache_min_8)
ii_8 = ior( ishft(ii_8,7), j-mo_integrals_cache_min_8)
ii_8 = ior( ishft(ii_8,7), i-mo_integrals_cache_min_8)
get_mo_bielec_integral = mo_integrals_cache(ii_8)
endif
end