mirror of
https://gitlab.com/scemama/qmcchem.git
synced 2024-09-06 10:04:45 +02:00
353 lines
11 KiB
FortranFixed
353 lines
11 KiB
FortranFixed
|
BEGIN_PROVIDER [ logical, compute_all_aos ]
|
||
|
implicit none
|
||
|
BEGIN_DOC
|
||
|
! If true, compute all AOs
|
||
|
END_DOC
|
||
|
|
||
|
compute_all_aos = .not. do_nucl_fitcusp
|
||
|
END_PROVIDER
|
||
|
|
||
|
|
||
|
BEGIN_PROVIDER [ integer, ao_elec ]
|
||
|
implicit none
|
||
|
BEGIN_DOC
|
||
|
! Current electron in AO evaluation
|
||
|
END_DOC
|
||
|
ao_elec = 1
|
||
|
END_PROVIDER
|
||
|
|
||
|
|
||
|
BEGIN_PROVIDER [ real, ao_radius, (ao_num_8) ]
|
||
|
&BEGIN_PROVIDER [ real, ao_radius_sorted, (ao_num_8) ]
|
||
|
&BEGIN_PROVIDER [ real, nucl_radius, (nucl_num_8) ]
|
||
|
&BEGIN_PROVIDER [ integer, ao_radius_order, (ao_num_8) ]
|
||
|
implicit none
|
||
|
BEGIN_DOC
|
||
|
! Radius of the AOs
|
||
|
!
|
||
|
! Radius of the nuclei
|
||
|
!
|
||
|
! Indices of the AOs per nucleus ordered by ao_radius in descending order
|
||
|
END_DOC
|
||
|
integer :: i,k, inucl
|
||
|
real :: to_sort(ao_num)
|
||
|
!DIR$ VECTOR ALIGNED
|
||
|
nucl_radius = 0.
|
||
|
!DIR$ VECTOR ALIGNED
|
||
|
do i=1,ao_num
|
||
|
ao_radius_order(i) = i
|
||
|
enddo
|
||
|
do i=1,ao_num
|
||
|
ao_radius(i) = 0.
|
||
|
do k=1,ao_prim_num(i)
|
||
|
ao_radius(i) = max(20./ao_expo_transp(k,i), ao_radius(i))
|
||
|
enddo
|
||
|
to_sort(i) = 1./ao_radius(i)
|
||
|
inucl = ao_nucl(i)
|
||
|
nucl_radius(inucl) = max(nucl_radius(inucl),ao_radius(i))
|
||
|
enddo
|
||
|
!DIR$ VECTOR ALIGNED
|
||
|
nucl_radius = sqrt(nucl_radius)
|
||
|
integer :: istart, iend
|
||
|
do inucl=1,nucl_num
|
||
|
istart = ao_nucl_idx(1,inucl)
|
||
|
iend = ao_nucl_idx(2,inucl)
|
||
|
call isort(to_sort(istart),ao_radius_order(istart),(iend-istart+1))
|
||
|
enddo
|
||
|
ao_radius_sorted = ao_radius
|
||
|
call set_order(ao_radius_sorted,ao_radius_order,ao_num)
|
||
|
|
||
|
END_PROVIDER
|
||
|
|
||
|
BEGIN_PROVIDER [ real, ao_oneD_block, (ao_num_8) ]
|
||
|
&BEGIN_PROVIDER [ integer, ao_oneD_prim_non_zero_idx, ((-simd_sp+1):ao_num_8) ]
|
||
|
&BEGIN_PROVIDER [ real, ao_oneD_lapl_block, (ao_num_8) ]
|
||
|
&BEGIN_PROVIDER [ real, ao_oneD_grad_block_z, (ao_num_8) ]
|
||
|
&BEGIN_PROVIDER [ real, ao_oneD_grad_block_y, (ao_num_8) ]
|
||
|
&BEGIN_PROVIDER [ real, ao_oneD_grad_block_x, (ao_num_8) ]
|
||
|
&BEGIN_PROVIDER [ integer, ao_block_num ]
|
||
|
&BEGIN_PROVIDER [ integer, ao_block_num_8 ]
|
||
|
implicit none
|
||
|
BEGIN_DOC
|
||
|
! Gradients and Laplacians of the primitive AOs
|
||
|
END_DOC
|
||
|
|
||
|
PROVIDE nucl_elec_dist
|
||
|
PROVIDE ao_power
|
||
|
PROVIDE ao_expo_unique
|
||
|
PROVIDE ao_expo_transp
|
||
|
PROVIDE ao_coef_transp
|
||
|
PROVIDE compute_all_aos
|
||
|
PROVIDE ao_prim_num
|
||
|
PROVIDE ao_expo
|
||
|
PROVIDE nucl_fitcusp_radius
|
||
|
|
||
|
integer, save :: ifirst = 0
|
||
|
if (ifirst == 0) then
|
||
|
ifirst = 1
|
||
|
!DIR$ VECTOR ALIGNED
|
||
|
ao_oneD_prim_non_zero_idx = 0.
|
||
|
ao_oneD_block = 0.
|
||
|
ao_oneD_lapl_block = 0.
|
||
|
ao_oneD_grad_block_x = 0.
|
||
|
ao_oneD_grad_block_y = 0.
|
||
|
ao_oneD_grad_block_z = 0.
|
||
|
endif
|
||
|
|
||
|
real :: ao_oneD_prim_block
|
||
|
real :: buffer2
|
||
|
|
||
|
real :: b1, b2, b3, b4
|
||
|
integer :: i, j, k, idx, l, m
|
||
|
real :: rtemp, r, r2
|
||
|
|
||
|
integer :: inucl, istart, iend, idx0, ii
|
||
|
real :: buffer(-simd_sp+1:max(ao_nucl_idx(2,nucl_num), ao_expo_unique_nucl_idx(2,nucl_num)) )
|
||
|
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: buffer
|
||
|
|
||
|
j = ao_elec
|
||
|
|
||
|
real :: d
|
||
|
real, parameter :: rdm = 20.
|
||
|
idx = 1
|
||
|
!DIR$ VECTOR ALIGNED
|
||
|
buffer = 0.
|
||
|
do inucl=1,nucl_num
|
||
|
r = nucl_elec_dist(inucl,j)
|
||
|
if (r > nucl_radius(inucl)) then
|
||
|
cycle
|
||
|
else
|
||
|
r2 = r*r
|
||
|
istart = ao_expo_unique_nucl_idx(1,inucl)
|
||
|
if (ao_expo_unique(istart)*r2 > rdm) then
|
||
|
cycle
|
||
|
endif
|
||
|
real :: dcomp
|
||
|
dcomp = rdm/r2
|
||
|
b1 = -nucl_elec_dist_vec(1,inucl,j)-nucl_elec_dist_vec(1,inucl,j)
|
||
|
b2 = -nucl_elec_dist_vec(2,inucl,j)-nucl_elec_dist_vec(2,inucl,j)
|
||
|
b3 = -nucl_elec_dist_vec(3,inucl,j)-nucl_elec_dist_vec(3,inucl,j)
|
||
|
b4 = 4.*r2
|
||
|
iend = ao_expo_unique_nucl_idx(2,inucl)
|
||
|
|
||
|
integer :: imax, istep
|
||
|
imax=istart+1
|
||
|
!DIR$ LOOP COUNT(20)
|
||
|
do while (imax<=iend .and. ao_expo_unique(imax) <= dcomp)
|
||
|
imax += 1
|
||
|
enddo
|
||
|
imax = imax-1
|
||
|
|
||
|
buffer(istart) = exp(-r2*ao_expo_unique(istart))
|
||
|
if (imax - istart == 0) then
|
||
|
continue
|
||
|
else if (imax - istart == 1) then
|
||
|
buffer(istart+1) = exp(-r2*ao_expo_unique(istart+1))
|
||
|
else if (imax - istart == 2) then
|
||
|
buffer(istart+1) = exp(-r2*ao_expo_unique(istart+1))
|
||
|
buffer(istart+2) = exp(-r2*ao_expo_unique(istart+2))
|
||
|
else
|
||
|
!DIR$ VECTOR UNALIGNED
|
||
|
!DIR$ LOOP COUNT MIN(3)
|
||
|
do i=istart+1,imax
|
||
|
buffer(i) = exp(-r2*ao_expo_unique(i))
|
||
|
enddo
|
||
|
endif
|
||
|
istart = ao_nucl_idx(1,inucl)
|
||
|
iend = ao_nucl_idx(2,inucl)
|
||
|
idx0 = idx
|
||
|
|
||
|
if ( r > nucl_fitcusp_radius(inucl) ) then
|
||
|
|
||
|
!DIR$ VECTOR UNALIGNED
|
||
|
do ii=istart,iend
|
||
|
if ( r2 > ao_radius_sorted(ii)) then
|
||
|
exit
|
||
|
endif
|
||
|
ao_oneD_prim_non_zero_idx(idx) = ao_radius_order(ii)
|
||
|
idx += 1
|
||
|
enddo
|
||
|
|
||
|
else if (.not.compute_all_aos) then
|
||
|
|
||
|
!DIR$ VECTOR UNALIGNED
|
||
|
do ii=istart,iend
|
||
|
if ( r2 > ao_radius_sorted(ii)) then
|
||
|
exit
|
||
|
endif
|
||
|
i = ao_radius_order(ii)
|
||
|
if ( .not.ao_power_is_zero(i) ) then
|
||
|
ao_oneD_prim_non_zero_idx(idx) = i
|
||
|
idx += 1
|
||
|
endif
|
||
|
enddo
|
||
|
|
||
|
else
|
||
|
|
||
|
!DIR$ VECTOR UNALIGNED
|
||
|
do ii=istart,iend
|
||
|
if ( r2 > ao_radius_sorted(ii)) then
|
||
|
exit
|
||
|
endif
|
||
|
ao_oneD_prim_non_zero_idx(idx) = ao_radius_order(ii)
|
||
|
idx += 1
|
||
|
enddo
|
||
|
|
||
|
endif
|
||
|
|
||
|
|
||
|
iend = idx-1
|
||
|
|
||
|
real :: f, g, h
|
||
|
integer :: lmax
|
||
|
i = ao_oneD_prim_non_zero_idx(idx0)
|
||
|
if (i==0) then
|
||
|
cycle
|
||
|
endif
|
||
|
lmax = ao_prim_num(i)
|
||
|
g = buffer(ao_expo_unique_idx_1(i))
|
||
|
f = ao_coef(i,1)*g*ao_expo(i,1)
|
||
|
ao_oneD_block(idx0) = ao_coef(i,1)*g
|
||
|
ao_oneD_lapl_block(idx0) = f*(b4 * ao_expo(i,1)-6.)
|
||
|
ao_oneD_grad_block_x(idx0) = b1*f
|
||
|
ao_oneD_grad_block_y(idx0) = b2*f
|
||
|
ao_oneD_grad_block_z(idx0) = b3*f
|
||
|
|
||
|
select case (iend-idx0)
|
||
|
case (0)
|
||
|
continue
|
||
|
case (1)
|
||
|
i = ao_oneD_prim_non_zero_idx(idx0+1)
|
||
|
g = buffer(ao_expo_unique_idx_1(i))
|
||
|
f = ao_coef(i,1)*g*ao_expo(i,1)
|
||
|
ao_oneD_block(idx0+1) = ao_coef(i,1)*g
|
||
|
ao_oneD_lapl_block(idx0+1) = f*(b4 * ao_expo(i,1)-6.)
|
||
|
ao_oneD_grad_block_x(idx0+1) = b1*f
|
||
|
ao_oneD_grad_block_y(idx0+1) = b2*f
|
||
|
ao_oneD_grad_block_z(idx0+1) = b3*f
|
||
|
lmax = max(ao_prim_num(i),lmax)
|
||
|
case (2)
|
||
|
i = ao_oneD_prim_non_zero_idx(idx0+1)
|
||
|
g = buffer(ao_expo_unique_idx_1(i))
|
||
|
f = ao_coef(i,1)*g*ao_expo(i,1)
|
||
|
ao_oneD_block(idx0+1) = ao_coef(i,1)*g
|
||
|
ao_oneD_lapl_block(idx0+1) = f*(b4 * ao_expo(i,1)-6.)
|
||
|
ao_oneD_grad_block_x(idx0+1) = b1*f
|
||
|
ao_oneD_grad_block_y(idx0+1) = b2*f
|
||
|
ao_oneD_grad_block_z(idx0+1) = b3*f
|
||
|
lmax = max(ao_prim_num(i),lmax)
|
||
|
|
||
|
i = ao_oneD_prim_non_zero_idx(idx0+2)
|
||
|
g = buffer(ao_expo_unique_idx_1(i))
|
||
|
f = ao_coef(i,1)*g*ao_expo(i,1)
|
||
|
ao_oneD_block(idx0+2) = ao_coef(i,1)*g
|
||
|
ao_oneD_lapl_block(idx0+2) = f*(b4 * ao_expo(i,1)-6.)
|
||
|
ao_oneD_grad_block_x(idx0+2) = b1*f
|
||
|
ao_oneD_grad_block_y(idx0+2) = b2*f
|
||
|
ao_oneD_grad_block_z(idx0+2) = b3*f
|
||
|
lmax = max(ao_prim_num(i),lmax)
|
||
|
case (3)
|
||
|
i = ao_oneD_prim_non_zero_idx(idx0+1)
|
||
|
g = buffer(ao_expo_unique_idx_1(i))
|
||
|
f = ao_coef(i,1)*g*ao_expo(i,1)
|
||
|
ao_oneD_block(idx0+1) = ao_coef(i,1)*g
|
||
|
ao_oneD_lapl_block(idx0+1) = f*(b4 * ao_expo(i,1)-6.)
|
||
|
ao_oneD_grad_block_x(idx0+1) = b1*f
|
||
|
ao_oneD_grad_block_y(idx0+1) = b2*f
|
||
|
ao_oneD_grad_block_z(idx0+1) = b3*f
|
||
|
lmax = max(ao_prim_num(i),lmax)
|
||
|
|
||
|
i = ao_oneD_prim_non_zero_idx(idx0+2)
|
||
|
g = buffer(ao_expo_unique_idx_1(i))
|
||
|
f = ao_coef(i,1)*g*ao_expo(i,1)
|
||
|
ao_oneD_block(idx0+2) = ao_coef(i,1)*g
|
||
|
ao_oneD_lapl_block(idx0+2) = f*(b4 * ao_expo(i,1)-6.)
|
||
|
ao_oneD_grad_block_x(idx0+2) = b1*f
|
||
|
ao_oneD_grad_block_y(idx0+2) = b2*f
|
||
|
ao_oneD_grad_block_z(idx0+2) = b3*f
|
||
|
lmax = max(ao_prim_num(i),lmax)
|
||
|
|
||
|
i = ao_oneD_prim_non_zero_idx(idx0+3)
|
||
|
g = buffer(ao_expo_unique_idx_1(i))
|
||
|
f = ao_coef(i,1)*g*ao_expo(i,1)
|
||
|
ao_oneD_block(idx0+3) = ao_coef(i,1)*g
|
||
|
ao_oneD_lapl_block(idx0+3) = f*(b4 * ao_expo(i,1)-6.)
|
||
|
ao_oneD_grad_block_x(idx0+3) = b1*f
|
||
|
ao_oneD_grad_block_y(idx0+3) = b2*f
|
||
|
ao_oneD_grad_block_z(idx0+3) = b3*f
|
||
|
lmax = max(ao_prim_num(i),lmax)
|
||
|
|
||
|
case default
|
||
|
!DIR$ VECTOR UNALIGNED
|
||
|
!DIR$ LOOP COUNT MIN(4)
|
||
|
do idx=idx0+1,iend
|
||
|
i = ao_oneD_prim_non_zero_idx(idx)
|
||
|
g = buffer(ao_expo_unique_idx_1(i))
|
||
|
f = ao_coef(i,1)*g*ao_expo(i,1)
|
||
|
ao_oneD_block(idx) = ao_coef(i,1)*g
|
||
|
ao_oneD_lapl_block(idx) = f*(b4 * ao_expo(i,1)-6.)
|
||
|
ao_oneD_grad_block_x(idx) = b1*f
|
||
|
ao_oneD_grad_block_y(idx) = b2*f
|
||
|
ao_oneD_grad_block_z(idx) = b3*f
|
||
|
lmax = max(ao_prim_num(i),lmax)
|
||
|
enddo
|
||
|
|
||
|
end select
|
||
|
|
||
|
if (lmax == 1) then
|
||
|
idx = iend+1
|
||
|
cycle
|
||
|
endif
|
||
|
|
||
|
do idx=idx0,iend
|
||
|
i = ao_oneD_prim_non_zero_idx(idx)
|
||
|
if (ao_prim_num(i) == 1) then
|
||
|
cycle
|
||
|
endif
|
||
|
l = ao_prim_num(i)
|
||
|
|
||
|
ao_oneD_block(idx) = 0.
|
||
|
ao_oneD_lapl_block(idx) = 0.
|
||
|
f = 0.
|
||
|
g = 0.
|
||
|
h = 0.
|
||
|
!DIR$ VECTOR ALIGNED
|
||
|
!DIR$ LOOP COUNT (8)
|
||
|
do k=1,l
|
||
|
ao_oneD_prim_block = buffer(ao_expo_unique_idx(k,i))*ao_coef_transp(k,i)
|
||
|
g = g + ao_oneD_prim_block
|
||
|
buffer2 = ao_expo_transp(k,i)*ao_oneD_prim_block
|
||
|
f = f+buffer2
|
||
|
h = h + buffer2 * (b4*ao_expo_transp(k,i)-6.)
|
||
|
enddo
|
||
|
ao_oneD_block(idx) = g
|
||
|
ao_oneD_lapl_block(idx) = h
|
||
|
ao_oneD_grad_block_x(idx) = b1*f
|
||
|
ao_oneD_grad_block_y(idx) = b2*f
|
||
|
ao_oneD_grad_block_z(idx) = b3*f
|
||
|
enddo
|
||
|
idx = iend+1
|
||
|
endif
|
||
|
enddo
|
||
|
|
||
|
ao_oneD_prim_non_zero_idx(0) = idx-1
|
||
|
if (ao_oneD_prim_non_zero_idx(0) > ao_block_num) then
|
||
|
ao_block_num = ao_oneD_prim_non_zero_idx(0)
|
||
|
integer :: mod_align
|
||
|
ao_block_num_8 = mod_align(ao_block_num)
|
||
|
endif
|
||
|
|
||
|
if (ao_oneD_prim_non_zero_idx(0) == 0) then
|
||
|
ao_oneD_prim_non_zero_idx(0) = 1
|
||
|
ao_oneD_prim_non_zero_idx(1) = 1
|
||
|
ao_oneD_block(idx) = 0.
|
||
|
ao_oneD_lapl_block(idx) = 0.
|
||
|
ao_oneD_grad_block_x(idx) = 0.
|
||
|
ao_oneD_grad_block_y(idx) = 0.
|
||
|
ao_oneD_grad_block_z(idx) = 0.
|
||
|
endif
|
||
|
|
||
|
END_PROVIDER
|
||
|
|
||
|
|