LCPQ/qmcchem
10
1
Fork 1
mirror of https://gitlab.com/scemama/qmcchem.git synced 2024-12-21 11:53:30 +01:00
Code Releases Activity

DIR -> OMP

Browse Source
This commit is contained in:
Anthony Scemama 2021-05-31 14:01:33 +02:00
parent 216ddc4351
commit 248545e558
3 changed files with 257 additions and 257 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
configure.sh
View File

@ -39,11 +39,11 @@ export PATH="\${QMCCHEM_PATH}/bin:\${PATH}"
export LD_LIBRARY_PATH="\${QMCCHEM_PATH}/lib:\${LD_LIBRARY_PATH}" export LD_LIBRARY_PATH="\${QMCCHEM_PATH}/lib:\${LD_LIBRARY_PATH}"
export LIBRARY_PATH="\${QMCCHEM_PATH}/lib:\${LIBRARY_PATH}" export LIBRARY_PATH="\${QMCCHEM_PATH}/lib:\${LIBRARY_PATH}"
export QMCCHEM_MPIRUN="mpirun" export QMCCHEM_MPIRUN="mpirun"
export QMCCHEM_MPIRUN_FLAGS="--bind-to-core" export QMCCHEM_MPIRUN_FLAGS=""
export C_INCLUDE_PATH="\${QMCCHEM_PATH}/include:\${C_INCLUDE_PATH}" export C_INCLUDE_PATH="\${QMCCHEM_PATH}/include:\${C_INCLUDE_PATH}"
#export QMCCHEM_NIC=ib0 #export QMCCHEM_NIC=ib0
source \${QMCCHEM_PATH}/irpf90/bin/irpman source \${QMCCHEM_PATH}/irpf90/bin/irpman
#source \${QMCCHEM_PATH}/EZFIO/Bash/ezfio.sh #source \${QMCCHEM_PATH}/EZFIO/Bash/ezfio.sh
eval $(opam env) eval \$(opam env)
EOF EOF

340
src/det.irp.f
View File

File diff suppressed because it is too large Load Diff

170
src/mo.irp.f
View File

@ -5,12 +5,12 @@
! Number of Molecular orbitals ! Number of Molecular orbitals
END_DOC END_DOC
integer, external :: mod_align integer, external :: mod_align
mo_num = maxval(present_mos) mo_num = maxval(present_mos)
call iinfo(irp_here,'mo_num',mo_num) call iinfo(irp_here,'mo_num',mo_num)
mo_num_8 = mod_align(mo_num) mo_num_8 = mod_align(mo_num)
END_PROVIDER END_PROVIDER
@ -22,7 +22,7 @@ BEGIN_PROVIDER [ real, mo_coef_input, (ao_num_8,mo_tot_num) ]
integer :: i, j integer :: i, j
real,allocatable :: buffer(:,:) real,allocatable :: buffer(:,:)
allocate (buffer(ao_num,mo_tot_num)) allocate (buffer(ao_num,mo_tot_num))
buffer = 0. buffer = 0.
call get_mo_basis_mo_coef(buffer) call get_mo_basis_mo_coef(buffer)
do i=1,mo_tot_num do i=1,mo_tot_num
@ -35,7 +35,7 @@ BEGIN_PROVIDER [ real, mo_coef_input, (ao_num_8,mo_tot_num) ]
call set_order(mo_coef_input(1,i),ao_nucl_sort_idx,ao_num) call set_order(mo_coef_input(1,i),ao_nucl_sort_idx,ao_num)
enddo enddo
deallocate(buffer) deallocate(buffer)
END_PROVIDER END_PROVIDER
@ -56,7 +56,7 @@ BEGIN_PROVIDER [ real, mo_coef, (ao_num_8,mo_num_8) ]
! Molecular orbital coefficients ! Molecular orbital coefficients
END_DOC END_DOC
integer :: i, j integer :: i, j
do j=1,mo_num do j=1,mo_num
do i=1,ao_num_8 do i=1,ao_num_8
mo_coef(i,j) = mo_coef_input(i,j) mo_coef(i,j) = mo_coef_input(i,j)
@ -68,10 +68,10 @@ BEGIN_PROVIDER [ real, mo_coef, (ao_num_8,mo_num_8) ]
mo_coef(i,j) = 0. mo_coef(i,j) = 0.
enddo enddo
enddo enddo
! Input MOs are not needed any more ! Input MOs are not needed any more
FREE mo_coef_input FREE mo_coef_input
real :: f real :: f
f = 1./mo_scale f = 1./mo_scale
do j=1,mo_num do j=1,mo_num
@ -80,7 +80,7 @@ BEGIN_PROVIDER [ real, mo_coef, (ao_num_8,mo_num_8) ]
mo_coef(i,j) *= f mo_coef(i,j) *= f
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
@ -90,7 +90,7 @@ BEGIN_PROVIDER [ real, mo_coef_transp, (mo_num_8,ao_num_8) ]
! Transpose of the Molecular orbital coefficients ! Transpose of the Molecular orbital coefficients
END_DOC END_DOC
call transpose(mo_coef,ao_num_8,mo_coef_transp,mo_num_8,ao_num_8,mo_num_8) call transpose(mo_coef,ao_num_8,mo_coef_transp,mo_num_8,ao_num_8,mo_num_8)
END_PROVIDER END_PROVIDER
@ -102,7 +102,7 @@ END_PROVIDER
! orbital coefficients ! orbital coefficients
END_DOC END_DOC
integer :: i, j integer :: i, j
integer :: idx integer :: idx
mo_coef_transp_sparsity = 0. mo_coef_transp_sparsity = 0.
do j=1,ao_num do j=1,ao_num
@ -117,7 +117,7 @@ END_PROVIDER
mo_coef_transp_sparsity += float(idx) mo_coef_transp_sparsity += float(idx)
enddo enddo
mo_coef_transp_sparsity *= 1./(mo_num*ao_num) mo_coef_transp_sparsity *= 1./(mo_num*ao_num)
END_PROVIDER END_PROVIDER
@ -133,7 +133,7 @@ BEGIN_PROVIDER [ real, mo_coef_transp_present, (num_present_mos_8,ao_num_8) ]
mo_coef_transp_present(j,i) = mo_coef_transp(present_mos(j),i) mo_coef_transp_present(j,i) = mo_coef_transp(present_mos(j),i)
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
@ -143,24 +143,24 @@ END_PROVIDER
&BEGIN_PROVIDER [ real, mo_grad_transp_z, (mo_num_8,elec_num) ] &BEGIN_PROVIDER [ real, mo_grad_transp_z, (mo_num_8,elec_num) ]
&BEGIN_PROVIDER [ real, mo_lapl_transp, (mo_num_8,elec_num) ] &BEGIN_PROVIDER [ real, mo_lapl_transp, (mo_num_8,elec_num) ]
implicit none implicit none
BEGIN_DOC BEGIN_DOC
! Values, gradients, laplacians of the molecular orbitals ! Values, gradients, laplacians of the molecular orbitals
! !
! Arrays are padded for efficiency ! Arrays are padded for efficiency
END_DOC END_DOC
integer :: i, j, k, l, m integer :: i, j, k, l, m
PROVIDE primitives_reduced PROVIDE primitives_reduced
if (do_nucl_fitcusp) then if (do_nucl_fitcusp) then
PROVIDE nucl_fitcusp_param PROVIDE nucl_fitcusp_param
PROVIDE nucl_elec_dist_vec PROVIDE nucl_elec_dist_vec
PROVIDE nucl_elec_dist_inv PROVIDE nucl_elec_dist_inv
endif endif
do i=1,elec_num do i=1,elec_num
if (i>1) then if (i>1) then
ao_elec = i ao_elec = i
@ -180,7 +180,7 @@ END_PROVIDER
mo_grad_transp_z(1,i), & mo_grad_transp_z(1,i), &
mo_lapl_transp(1,i), & mo_lapl_transp(1,i), &
ao_num) ao_num)
else else
call sparse_full_mv(mo_coef_transp_present,num_present_mos_8, & call sparse_full_mv(mo_coef_transp_present,num_present_mos_8, &
ao_value_block(1),ao_num_8, & ao_value_block(1),ao_num_8, &
@ -195,7 +195,7 @@ END_PROVIDER
mo_grad_transp_z(1,i), & mo_grad_transp_z(1,i), &
mo_lapl_transp(1,i), & mo_lapl_transp(1,i), &
ao_num) ao_num)
do j=num_present_mos,1,-1 do j=num_present_mos,1,-1
mo_value_transp (present_mos(j),i) = mo_value_transp (j,i) mo_value_transp (present_mos(j),i) = mo_value_transp (j,i)
mo_grad_transp_x(present_mos(j),i) = mo_grad_transp_x(j,i) mo_grad_transp_x(present_mos(j),i) = mo_grad_transp_x(j,i)
@ -207,10 +207,10 @@ END_PROVIDER
endif endif
enddo enddo
endif endif
if (do_nucl_fitcusp) then if (do_nucl_fitcusp) then
real :: r, r2, r_inv, d, expzr, Z, Z2, a, b, c, phi, rx, ry, rz real :: r, r2, r_inv, d, expzr, Z, Z2, a, b, c, phi, rx, ry, rz
do k=1,nucl_num do k=1,nucl_num
r = nucl_elec_dist(k,i) r = nucl_elec_dist(k,i)
if (r > nucl_fitcusp_radius(k)) then if (r > nucl_fitcusp_radius(k)) then
@ -235,14 +235,14 @@ END_PROVIDER
enddo enddo
exit exit
enddo ! k enddo ! k
endif endif
enddo ! i enddo ! i
ao_elec = 1 ao_elec = 1
SOFT_TOUCH ao_elec SOFT_TOUCH ao_elec
if (do_prepare) then if (do_prepare) then
real :: lambda, t real :: lambda, t
! Scale off-diagonal elements ! Scale off-diagonal elements
@ -268,7 +268,7 @@ END_PROVIDER
enddo enddo
enddo enddo
endif endif
do i=1,mo_num do i=1,mo_num
do j=1,elec_num do j=1,elec_num
mo_value_transp(i,j) *= mo_cusp_rescale(i) mo_value_transp(i,j) *= mo_cusp_rescale(i)
@ -285,10 +285,10 @@ BEGIN_PROVIDER [ real, mo_value, (elec_num_8,mo_num) ]
BEGIN_DOC BEGIN_DOC
! Values of the molecular orbitals ! Values of the molecular orbitals
END_DOC END_DOC
integer :: i,j integer :: i,j
integer, save :: ifirst = 0 integer, save :: ifirst = 0
if (ifirst == 0) then if (ifirst == 0) then
ifirst = 1 ifirst = 1
PROVIDE primitives_reduced PROVIDE primitives_reduced
@ -296,7 +296,7 @@ BEGIN_PROVIDER [ real, mo_value, (elec_num_8,mo_num) ]
mo_value = 0. mo_value = 0.
endif endif
call transpose(mo_value_transp(1,1),mo_num_8,mo_value,elec_num_8,mo_num,elec_num) call transpose(mo_value_transp(1,1),mo_num_8,mo_value,elec_num_8,mo_num,elec_num)
END_PROVIDER END_PROVIDER
@ -304,14 +304,14 @@ END_PROVIDER
&BEGIN_PROVIDER [ double precision, mo_grad_y, (elec_num_8,mo_num) ] &BEGIN_PROVIDER [ double precision, mo_grad_y, (elec_num_8,mo_num) ]
&BEGIN_PROVIDER [ double precision, mo_grad_z, (elec_num_8,mo_num) ] &BEGIN_PROVIDER [ double precision, mo_grad_z, (elec_num_8,mo_num) ]
implicit none implicit none
BEGIN_DOC BEGIN_DOC
! Gradients of the molecular orbitals ! Gradients of the molecular orbitals
END_DOC END_DOC
integer :: i,j integer :: i,j
integer, save :: ifirst = 0 integer, save :: ifirst = 0
if (ifirst == 0) then if (ifirst == 0) then
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
mo_grad_x = 0.d0 mo_grad_x = 0.d0
@ -327,7 +327,7 @@ END_PROVIDER
call transpose_to_dp(mo_grad_transp_z(1,1),mo_num_8,mo_grad_z(1,1),elec_num_8,mo_num,elec_num) call transpose_to_dp(mo_grad_transp_z(1,1),mo_num_8,mo_grad_z(1,1),elec_num_8,mo_num,elec_num)
call transpose_to_dp(mo_grad_transp_x(1,1),mo_num_8,mo_grad_x(1,1),elec_num_8,mo_num,elec_num) call transpose_to_dp(mo_grad_transp_x(1,1),mo_num_8,mo_grad_x(1,1),elec_num_8,mo_num,elec_num)
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ double precision, mo_lapl, (elec_num_8,mo_num) ] BEGIN_PROVIDER [ double precision, mo_lapl, (elec_num_8,mo_num) ]
@ -335,7 +335,7 @@ BEGIN_PROVIDER [ double precision, mo_lapl, (elec_num_8,mo_num) ]
BEGIN_DOC BEGIN_DOC
! Laplacians of the molecular orbitals ! Laplacians of the molecular orbitals
END_DOC END_DOC
integer :: i,j integer :: i,j
integer, save :: ifirst = 0 integer, save :: ifirst = 0
if (ifirst == 0) then if (ifirst == 0) then
@ -345,7 +345,7 @@ BEGIN_PROVIDER [ double precision, mo_lapl, (elec_num_8,mo_num) ]
mo_lapl = 0.d0 mo_lapl = 0.d0
endif endif
call transpose_to_dp(mo_lapl_transp(1,1),mo_num_8,mo_lapl,elec_num_8,mo_num,elec_num) call transpose_to_dp(mo_lapl_transp(1,1),mo_num_8,mo_lapl,elec_num_8,mo_num,elec_num)
END_PROVIDER END_PROVIDER
@ -358,21 +358,21 @@ END_PROVIDER
integer :: i,j integer :: i,j
do j=1,mo_num do j=1,mo_num
do i=1,elec_alpha_num do i=1,elec_alpha_num
mo_grad_lapl_alpha(1,i,j) = mo_grad_transp_x(j,i) mo_grad_lapl_alpha(1,i,j) = mo_grad_transp_x(j,i)
mo_grad_lapl_alpha(2,i,j) = mo_grad_transp_y(j,i) mo_grad_lapl_alpha(2,i,j) = mo_grad_transp_y(j,i)
mo_grad_lapl_alpha(3,i,j) = mo_grad_transp_z(j,i) mo_grad_lapl_alpha(3,i,j) = mo_grad_transp_z(j,i)
mo_grad_lapl_alpha(4,i,j) = mo_lapl_transp (j,i) mo_grad_lapl_alpha(4,i,j) = mo_lapl_transp (j,i)
enddo enddo
enddo enddo
do j=1,mo_num do j=1,mo_num
do i=elec_alpha_num+1,elec_num do i=elec_alpha_num+1,elec_num
mo_grad_lapl_beta(1,i,j) = mo_grad_transp_x(j,i) mo_grad_lapl_beta(1,i,j) = mo_grad_transp_x(j,i)
mo_grad_lapl_beta(2,i,j) = mo_grad_transp_y(j,i) mo_grad_lapl_beta(2,i,j) = mo_grad_transp_y(j,i)
mo_grad_lapl_beta(3,i,j) = mo_grad_transp_z(j,i) mo_grad_lapl_beta(3,i,j) = mo_grad_transp_z(j,i)
mo_grad_lapl_beta(4,i,j) = mo_lapl_transp (j,i) mo_grad_lapl_beta(4,i,j) = mo_lapl_transp (j,i)
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ double precision, mo_grad_lapl_transp, (4,mo_num,elec_num) ] BEGIN_PROVIDER [ double precision, mo_grad_lapl_transp, (4,mo_num,elec_num) ]
@ -389,7 +389,7 @@ BEGIN_PROVIDER [ double precision, mo_grad_lapl_transp, (4,mo_num,elec_num) ]
mo_grad_lapl_transp(4,j,i) = mo_lapl_transp (j,i) mo_grad_lapl_transp(4,j,i) = mo_lapl_transp (j,i)
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
@ -404,18 +404,18 @@ END_PROVIDER
BEGIN_PROVIDER [ integer, mo_tot_num ] BEGIN_PROVIDER [ integer, mo_tot_num ]
BEGIN_DOC BEGIN_DOC
! Total number of MOs in the EZFIO file ! Total number of MOs in the EZFIO file
END_DOC END_DOC
mo_tot_num = -1 mo_tot_num = -1
call get_mo_basis_mo_num(mo_tot_num) call get_mo_basis_mo_num(mo_tot_num)
if (mo_tot_num <= 0) then if (mo_tot_num <= 0) then
call abrt(irp_here,'Total number of MOs can''t be <0') call abrt(irp_here,'Total number of MOs can''t be <0')
endif endif
call iinfo(irp_here,'mo_tot_num',mo_tot_num) call iinfo(irp_here,'mo_tot_num',mo_tot_num)
END_PROVIDER END_PROVIDER
@ -432,16 +432,16 @@ BEGIN_PROVIDER [ double precision , mo_value_at_nucl, (mo_num_8,nucl_num) ]
END_DOC END_DOC
integer :: i, j, k, l integer :: i, j, k, l
real :: ao_value_at_nucl_no_S(ao_num) real :: ao_value_at_nucl_no_S(ao_num)
PROVIDE mo_fitcusp_normalization_before PROVIDE mo_fitcusp_normalization_before
do k=1,nucl_num do k=1,nucl_num
point(1) = nucl_coord(k,1) point(1) = nucl_coord(k,1)
point(2) = nucl_coord(k,2) point(2) = nucl_coord(k,2)
point(3) = nucl_coord(k,3) point(3) = nucl_coord(k,3)
TOUCH point TOUCH point
PROVIDE ao_value_p PROVIDE ao_value_p
do i=1,ao_num do i=1,ao_num
if (ao_nucl(i) /= k) then if (ao_nucl(i) /= k) then
ao_value_at_nucl_no_S(i) = ao_value_p(i) ao_value_at_nucl_no_S(i) = ao_value_p(i)
@ -449,7 +449,7 @@ BEGIN_PROVIDER [ double precision , mo_value_at_nucl, (mo_num_8,nucl_num) ]
ao_value_at_nucl_no_S(i) = 0. ao_value_at_nucl_no_S(i) = 0.
endif endif
enddo enddo
integer :: jj integer :: jj
do jj=1,num_present_mos do jj=1,num_present_mos
j = present_mos(jj) j = present_mos(jj)
@ -459,11 +459,11 @@ BEGIN_PROVIDER [ double precision , mo_value_at_nucl, (mo_num_8,nucl_num) ]
mo_value_at_nucl(j,k) = mo_value_at_nucl(j,k) + mo_coef(i,j)*ao_value_at_nucl_no_S(i) mo_value_at_nucl(j,k) = mo_value_at_nucl(j,k) + mo_coef(i,j)*ao_value_at_nucl_no_S(i)
enddo enddo
enddo enddo
enddo enddo
FREE ao_value_p ao_grad_p ao_lapl_p ao_axis_grad_p ao_oned_grad_p ao_oned_prim_grad_p ao_oned_lapl_p ao_axis_lapl_p ao_oned_prim_lapl_p ao_oned_p ao_oned_prim_p ao_axis_p ao_axis_power_p FREE ao_value_p ao_grad_p ao_lapl_p ao_axis_grad_p ao_oned_grad_p ao_oned_prim_grad_p ao_oned_lapl_p ao_axis_lapl_p ao_oned_prim_lapl_p ao_oned_p ao_oned_prim_p ao_axis_p ao_axis_power_p
SOFT_TOUCH point SOFT_TOUCH point
END_PROVIDER END_PROVIDER
@ -474,15 +474,15 @@ END_PROVIDER
BEGIN_DOC BEGIN_DOC
! Values of the atomic orbitals without S components on atoms ! Values of the atomic orbitals without S components on atoms
END_DOC END_DOC
integer :: i, j, k integer :: i, j, k
do k=1,nucl_num do k=1,nucl_num
point(1) = nucl_coord(k,1) point(1) = nucl_coord(k,1)
point(2) = nucl_coord(k,2) point(2) = nucl_coord(k,2)
point(3) = nucl_coord(k,3)+ nucl_fitcusp_radius(k) point(3) = nucl_coord(k,3)+ nucl_fitcusp_radius(k)
TOUCH point TOUCH point
do j=1,ao_num do j=1,ao_num
ao_value_at_fitcusp_radius(j,k) = ao_value_p(j) ao_value_at_fitcusp_radius(j,k) = ao_value_p(j)
ao_grad_at_fitcusp_radius(j,k) = ao_grad_p(j,3) ao_grad_at_fitcusp_radius(j,k) = ao_grad_p(j,3)
@ -496,7 +496,7 @@ END_PROVIDER
enddo enddo
FREE ao_value_p ao_grad_p ao_lapl_p ao_axis_grad_p ao_oned_grad_p ao_oned_prim_grad_p ao_oned_lapl_p ao_axis_lapl_p ao_oned_prim_lapl_p ao_oned_p ao_oned_prim_p ao_axis_p ao_axis_power_p FREE ao_value_p ao_grad_p ao_lapl_p ao_axis_grad_p ao_oned_grad_p ao_oned_prim_grad_p ao_oned_lapl_p ao_axis_lapl_p ao_oned_prim_lapl_p ao_oned_p ao_oned_prim_p ao_axis_p ao_axis_power_p
SOFT_TOUCH point SOFT_TOUCH point
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ double precision, mo_fitcusp_normalization_before, (mo_tot_num) ] BEGIN_PROVIDER [ double precision, mo_fitcusp_normalization_before, (mo_tot_num) ]
@ -560,7 +560,7 @@ BEGIN_PROVIDER [ double precision, mo_fitcusp_normalization_after, (mo_tot_num)
t = 0.d0 t = 0.d0
do j=1,ao_num do j=1,ao_num
if ( (ao_nucl(j) /= k).or.(ao_power(j,4) > 0) ) then if ( (ao_nucl(j) /= k).or.(ao_power(j,4) > 0) ) then
t = t + mo_coef(j,i) * ao_value_p(j) t = t + mo_coef(j,i) * ao_value_p(j)
endif endif
enddo enddo
t = t + nucl_fitcusp_param(1,i,k) + & t = t + nucl_fitcusp_param(1,i,k) + &
@ -601,7 +601,7 @@ END_PROVIDER
! Values of the molecular orbitals without S components on atoms ! Values of the molecular orbitals without S components on atoms
END_DOC END_DOC
integer :: i, j, k, l integer :: i, j, k, l
do k=1,nucl_num do k=1,nucl_num
do j=1,mo_num do j=1,mo_num
mo_value_at_fitcusp_radius(j,k) = 0.d0 mo_value_at_fitcusp_radius(j,k) = 0.d0
@ -615,7 +615,7 @@ END_PROVIDER
enddo enddo
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
@ -626,7 +626,7 @@ BEGIN_PROVIDER [ real, nucl_fitcusp_param, (4,mo_num,nucl_num) ]
END_DOC END_DOC
integer :: i,k, niter integer :: i,k, niter
character*(80) :: message character*(80) :: message
nucl_fitcusp_param = 0.d0 nucl_fitcusp_param = 0.d0
do k=1,nucl_num do k=1,nucl_num
double precision :: r, Z double precision :: r, Z
@ -642,14 +642,14 @@ BEGIN_PROVIDER [ real, nucl_fitcusp_param, (4,mo_num,nucl_num) ]
grad_phi = mo_grad_at_fitcusp_radius(i,k) grad_phi = mo_grad_at_fitcusp_radius(i,k)
phi = mo_value_at_fitcusp_radius(i,k) phi = mo_value_at_fitcusp_radius(i,k)
eta = mo_value_at_nucl(i,k) eta = mo_value_at_nucl(i,k)
nucl_fitcusp_param(1,i,k) = -(R*(2.d0*eta*Z-6.d0*grad_phi)+lap_phi*R*R+6.d0*phi)/(2.d0*R*Z-6.d0) nucl_fitcusp_param(1,i,k) = -(R*(2.d0*eta*Z-6.d0*grad_phi)+lap_phi*R*R+6.d0*phi)/(2.d0*R*Z-6.d0)
nucl_fitcusp_param(2,i,k) = (lap_phi*R*R*Z-6.d0*grad_phi*R*Z+6.d0*phi*Z+6.d0*eta*Z)/(2.d0*R*Z-6.d0) nucl_fitcusp_param(2,i,k) = (lap_phi*R*R*Z-6.d0*grad_phi*R*Z+6.d0*phi*Z+6.d0*eta*Z)/(2.d0*R*Z-6.d0)
nucl_fitcusp_param(3,i,k) = -(R*(-5.d0*grad_phi*Z-1.5d0*lap_phi)+lap_phi*R*R*Z+3.d0*phi*Z+& nucl_fitcusp_param(3,i,k) = -(R*(-5.d0*grad_phi*Z-1.5d0*lap_phi)+lap_phi*R*R*Z+3.d0*phi*Z+&
3.d0*eta*Z+6.d0*grad_phi)/(R*R*Z-3.d0*R) 3.d0*eta*Z+6.d0*grad_phi)/(R*R*Z-3.d0*R)
nucl_fitcusp_param(4,i,k) = (R*(-2.d0*grad_phi*Z-lap_phi)+0.5d0*lap_phi*R*R*Z+phi*Z+& nucl_fitcusp_param(4,i,k) = (R*(-2.d0*grad_phi*Z-lap_phi)+0.5d0*lap_phi*R*R*Z+phi*Z+&
eta*Z+3.d0*grad_phi)/(R*R*R*Z-3.d0*R*R) eta*Z+3.d0*grad_phi)/(R*R*R*Z-3.d0*R*R)
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
@ -689,7 +689,7 @@ subroutine sparse_full_mv(A,LDA, &
!DIR$ ASSUME_ALIGNED C3 : $IRP_ALIGN !DIR$ ASSUME_ALIGNED C3 : $IRP_ALIGN
!DIR$ ASSUME_ALIGNED C4 : $IRP_ALIGN !DIR$ ASSUME_ALIGNED C4 : $IRP_ALIGN
!DIR$ ASSUME_ALIGNED C5 : $IRP_ALIGN !DIR$ ASSUME_ALIGNED C5 : $IRP_ALIGN
integer :: kao, kmax, kmax2, kmax3 integer :: kao, kmax, kmax2, kmax3
integer :: i,j,k integer :: i,j,k
integer :: k_vec(8) integer :: k_vec(8)
@ -698,9 +698,9 @@ subroutine sparse_full_mv(A,LDA, &
real :: d21, d22, d23, d24, d25 real :: d21, d22, d23, d24, d25
real :: d31, d32, d33, d34, d35 real :: d31, d32, d33, d34, d35
real :: d41, d42, d43, d44, d45 real :: d41, d42, d43, d44, d45
! LDC and LDA have to be factors of simd_sp ! LDC and LDA have to be factors of simd_sp
! IRP_IF NO_PREFETCH ! IRP_IF NO_PREFETCH
! IRP_ELSE ! IRP_ELSE
! call MM_PREFETCH (A(1,indices(1)),3) ! call MM_PREFETCH (A(1,indices(1)),3)
@ -709,7 +709,7 @@ subroutine sparse_full_mv(A,LDA, &
! call MM_PREFETCH (A(1,indices(4)),3) ! call MM_PREFETCH (A(1,indices(4)),3)
! IRP_ENDIF ! IRP_ENDIF
!DIR$ SIMD !OMP$ SIMD
do j=1,LDC do j=1,LDC
C1(j) = 0. C1(j) = 0.
C2(j) = 0. C2(j) = 0.
@ -727,35 +727,35 @@ subroutine sparse_full_mv(A,LDA, &
k_vec(2) = indices(kao+1) k_vec(2) = indices(kao+1)
k_vec(3) = indices(kao+2) k_vec(3) = indices(kao+2)
k_vec(4) = indices(kao+3) k_vec(4) = indices(kao+3)
d11 = B1(kao ) d11 = B1(kao )
d21 = B1(kao+1) d21 = B1(kao+1)
d31 = B1(kao+2) d31 = B1(kao+2)
d41 = B1(kao+3) d41 = B1(kao+3)
d12 = B2(kao ) d12 = B2(kao )
d22 = B2(kao+1) d22 = B2(kao+1)
d32 = B2(kao+2) d32 = B2(kao+2)
d42 = B2(kao+3) d42 = B2(kao+3)
d13 = B3(kao ) d13 = B3(kao )
d23 = B3(kao+1) d23 = B3(kao+1)
d33 = B3(kao+2) d33 = B3(kao+2)
d43 = B3(kao+3) d43 = B3(kao+3)
d14 = B4(kao ) d14 = B4(kao )
d24 = B4(kao+1) d24 = B4(kao+1)
d34 = B4(kao+2) d34 = B4(kao+2)
d44 = B4(kao+3) d44 = B4(kao+3)
d15 = B5(kao ) d15 = B5(kao )
d25 = B5(kao+1) d25 = B5(kao+1)
d35 = B5(kao+2) d35 = B5(kao+2)
d45 = B5(kao+3) d45 = B5(kao+3)
do k=0,LDA-1,$IRP_ALIGN/4 do k=0,LDA-1,$IRP_ALIGN/4
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
!DIR$ SIMD FIRSTPRIVATE(d11,d21,d31,d41) !OMP$ SIMD FIRSTPRIVATE(d11,d21,d31,d41)
do j=1,$IRP_ALIGN/4 do j=1,$IRP_ALIGN/4
! IRP_IF NO_PREFETCH ! IRP_IF NO_PREFETCH
! IRP_ELSE ! IRP_ELSE
@ -767,18 +767,18 @@ subroutine sparse_full_mv(A,LDA, &
C1(j+k) = C1(j+k) + A(j+k,k_vec(1))*d11 + A(j+k,k_vec(2))*d21& C1(j+k) = C1(j+k) + A(j+k,k_vec(1))*d11 + A(j+k,k_vec(2))*d21&
+ A(j+k,k_vec(3))*d31 + A(j+k,k_vec(4))*d41 + A(j+k,k_vec(3))*d31 + A(j+k,k_vec(4))*d41
enddo enddo
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
!DIR$ SIMD FIRSTPRIVATE(d12,d22,d32,d42,d13,d23,d33,d43) !OMP$ SIMD FIRSTPRIVATE(d12,d22,d32,d42,d13,d23,d33,d43)
do j=1,$IRP_ALIGN/4 do j=1,$IRP_ALIGN/4
C2(j+k) = C2(j+k) + A(j+k,k_vec(1))*d12 + A(j+k,k_vec(2))*d22& C2(j+k) = C2(j+k) + A(j+k,k_vec(1))*d12 + A(j+k,k_vec(2))*d22&
+ A(j+k,k_vec(3))*d32 + A(j+k,k_vec(4))*d42 + A(j+k,k_vec(3))*d32 + A(j+k,k_vec(4))*d42
C3(j+k) = C3(j+k) + A(j+k,k_vec(1))*d13 + A(j+k,k_vec(2))*d23& C3(j+k) = C3(j+k) + A(j+k,k_vec(1))*d13 + A(j+k,k_vec(2))*d23&
+ A(j+k,k_vec(3))*d33 + A(j+k,k_vec(4))*d43 + A(j+k,k_vec(3))*d33 + A(j+k,k_vec(4))*d43
enddo enddo
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
!DIR$ SIMD FIRSTPRIVATE(d14,d24,d34,d44,d15,d25,d35,d45) !OMP$ SIMD FIRSTPRIVATE(d14,d24,d34,d44,d15,d25,d35,d45)
do j=1,$IRP_ALIGN/4 do j=1,$IRP_ALIGN/4
C4(j+k) = C4(j+k) + A(j+k,k_vec(1))*d14 + A(j+k,k_vec(2))*d24& C4(j+k) = C4(j+k) + A(j+k,k_vec(1))*d14 + A(j+k,k_vec(2))*d24&
+ A(j+k,k_vec(3))*d34 + A(j+k,k_vec(4))*d44 + A(j+k,k_vec(3))*d34 + A(j+k,k_vec(4))*d44
@ -786,9 +786,9 @@ subroutine sparse_full_mv(A,LDA, &
+ A(j+k,k_vec(3))*d35 + A(j+k,k_vec(4))*d45 + A(j+k,k_vec(3))*d35 + A(j+k,k_vec(4))*d45
enddo enddo
enddo enddo
enddo enddo
do kao = kmax2+1, kmax3 do kao = kmax2+1, kmax3
k_vec(1) = indices(kao) k_vec(1) = indices(kao)
d11 = B1(kao) d11 = B1(kao)
@ -799,7 +799,7 @@ subroutine sparse_full_mv(A,LDA, &
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
do k=0,LDA-1,$IRP_ALIGN/4 do k=0,LDA-1,$IRP_ALIGN/4
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
!DIR$ SIMD FIRSTPRIVATE(d11,d12,d13,d14,d15) !OMP$ SIMD FIRSTPRIVATE(d11,d12,d13,d14,d15)
do j=1,$IRP_ALIGN/4 do j=1,$IRP_ALIGN/4
C1(j+k) = C1(j+k) + A(j+k,k_vec(1))*d11 C1(j+k) = C1(j+k) + A(j+k,k_vec(1))*d11
C2(j+k) = C2(j+k) + A(j+k,k_vec(1))*d12 C2(j+k) = C2(j+k) + A(j+k,k_vec(1))*d12
@ -809,7 +809,7 @@ subroutine sparse_full_mv(A,LDA, &
enddo enddo
enddo enddo
enddo enddo
end end