LCPQ/qp2
9
1
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-12-22 03:23:29 +01:00
Code Releases Activity

non_h_ints_mu: combined

Browse Source
This commit is contained in:
Abdallah Ammar 2023-03-04 02:19:15 +01:00
parent 10b461f5a2
commit 5b58b062d9
5 changed files with 579 additions and 241 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

75
src/non_h_ints_mu/grad_squared.irp.f
View File

@ -299,7 +299,6 @@ BEGIN_PROVIDER [ double precision, u12sq_j1bsq, (ao_num, ao_num, n_points_final_
END_PROVIDER END_PROVIDER
! ---
! --- ! ---
BEGIN_PROVIDER [ double precision, u12_grad1_u12_j1b_grad1_j1b, (ao_num, ao_num, n_points_final_grid) ] BEGIN_PROVIDER [ double precision, u12_grad1_u12_j1b_grad1_j1b, (ao_num, ao_num, n_points_final_grid) ]
@ -364,12 +363,28 @@ BEGIN_PROVIDER [double precision, tc_grad_square_ao, (ao_num, ao_num, ao_num, ao
integer :: ipoint, i, j, k, l integer :: ipoint, i, j, k, l
double precision :: weight1, ao_ik_r, ao_i_r double precision :: weight1, ao_ik_r, ao_i_r
double precision :: time0, time1 double precision :: time0, time1
double precision, allocatable :: ac_mat(:,:,:,:), b_mat(:,:,:), tmp(:,:,:) double precision, allocatable :: b_mat(:,:,:), tmp(:,:,:)
print*, ' providing tc_grad_square_ao ...' print*, ' providing tc_grad_square_ao ...'
call wall_time(time0) call wall_time(time0)
allocate(ac_mat(ao_num,ao_num,ao_num,ao_num), b_mat(n_points_final_grid,ao_num,ao_num), tmp(ao_num,ao_num,n_points_final_grid)) if(read_tc_integ) then
open(unit=11, form="unformatted", file='tc_grad_square_ao', action="read")
do i = 1, ao_num
do j = 1, ao_num
do k = 1, ao_num
do l = 1, ao_num
read(11) tc_grad_square_ao(l,k,j,i)
enddo
enddo
enddo
enddo
close(11)
else
allocate(b_mat(n_points_final_grid,ao_num,ao_num), tmp(ao_num,ao_num,n_points_final_grid))
b_mat = 0.d0 b_mat = 0.d0
!$OMP PARALLEL & !$OMP PARALLEL &
@ -403,31 +418,45 @@ BEGIN_PROVIDER [double precision, tc_grad_square_ao, (ao_num, ao_num, ao_num, ao
!$OMP END DO !$OMP END DO
!$OMP END PARALLEL !$OMP END PARALLEL
tc_grad_square_ao = 0.d0
ac_mat = 0.d0
call dgemm( "N", "N", ao_num*ao_num, ao_num*ao_num, n_points_final_grid, 1.d0 & call dgemm( "N", "N", ao_num*ao_num, ao_num*ao_num, n_points_final_grid, 1.d0 &
, tmp(1,1,1), ao_num*ao_num, b_mat(1,1,1), n_points_final_grid & , tmp(1,1,1), ao_num*ao_num, b_mat(1,1,1), n_points_final_grid &
, 1.d0, ac_mat, ao_num*ao_num) , 1.d0, tc_grad_square_ao, ao_num*ao_num)
deallocate(tmp, b_mat) deallocate(tmp, b_mat)
!$OMP PARALLEL & call sum_A_At(tc_grad_square_ao(1,1,1,1), ao_num*ao_num)
!$OMP DEFAULT (NONE) &
!$OMP PRIVATE (i, j, k, l) &
!$OMP SHARED (ac_mat, tc_grad_square_ao, ao_num)
!$OMP DO SCHEDULE (static)
do j = 1, ao_num
do l = 1, ao_num
do i = 1, ao_num
do k = 1, ao_num
tc_grad_square_ao(k,i,l,j) = ac_mat(k,i,l,j) + ac_mat(l,j,k,i)
enddo
enddo
enddo
enddo
!$OMP END DO
!$OMP END PARALLEL
deallocate(ac_mat) !!$OMP PARALLEL &
!!$OMP DEFAULT (NONE) &
!!$OMP PRIVATE (i, j, k, l) &
!!$OMP SHARED (ac_mat, tc_grad_square_ao, ao_num)
!!$OMP DO SCHEDULE (static)
! do j = 1, ao_num
! do l = 1, ao_num
! do i = 1, ao_num
! do k = 1, ao_num
! tc_grad_square_ao(k,i,l,j) = ac_mat(k,i,l,j) + ac_mat(l,j,k,i)
! enddo
! enddo
! enddo
! enddo
!!$OMP END DO
!!$OMP END PARALLEL
endif
if(write_tc_integ) then
open(unit=11, form="unformatted", file='tc_grad_square_ao', action="write")
do i = 1, ao_num
do j = 1, ao_num
do k = 1, ao_num
do l = 1, ao_num
write(11) tc_grad_square_ao(l,k,j,i)
enddo
enddo
enddo
enddo
close(11)
endif
call wall_time(time1) call wall_time(time1)
print*, ' Wall time for tc_grad_square_ao = ', time1 - time0 print*, ' Wall time for tc_grad_square_ao = ', time1 - time0

174
src/non_h_ints_mu/grad_squared_manu.irp.f
View File

@ -11,11 +11,177 @@ BEGIN_PROVIDER [double precision, tc_grad_square_ao_test, (ao_num, ao_num, ao_nu
integer :: ipoint, i, j, k, l integer :: ipoint, i, j, k, l
double precision :: weight1, ao_ik_r, ao_i_r,contrib,contrib2 double precision :: weight1, ao_ik_r, ao_i_r,contrib,contrib2
double precision :: time0, time1 double precision :: time0, time1
double precision, allocatable :: ac_mat(:,:,:,:), b_mat(:,:,:), tmp(:,:,:) double precision, allocatable :: b_mat(:,:,:), tmp(:,:,:)
print*, ' providing tc_grad_square_ao_test ...' print*, ' providing tc_grad_square_ao_test ...'
call wall_time(time0) call wall_time(time0)
if(read_tc_integ) then
open(unit=11, form="unformatted", file='tc_grad_square_ao_test', action="read")
do i = 1, ao_num
do j = 1, ao_num
do k = 1, ao_num
do l = 1, ao_num
read(11) tc_grad_square_ao_test(l,k,j,i)
enddo
enddo
enddo
enddo
close(11)
else
provide u12sq_j1bsq_test u12_grad1_u12_j1b_grad1_j1b_test grad12_j12_test
allocate(b_mat(n_points_final_grid,ao_num,ao_num), tmp(ao_num,ao_num,n_points_final_grid))
b_mat = 0.d0
!$OMP PARALLEL &
!$OMP DEFAULT (NONE) &
!$OMP PRIVATE (i, k, ipoint) &
!$OMP SHARED (aos_in_r_array_transp, b_mat, ao_num, n_points_final_grid, final_weight_at_r_vector)
!$OMP DO SCHEDULE (static)
do i = 1, ao_num
do k = 1, ao_num
do ipoint = 1, n_points_final_grid
b_mat(ipoint,k,i) = final_weight_at_r_vector(ipoint) * aos_in_r_array_transp(ipoint,i) * aos_in_r_array_transp(ipoint,k)
enddo
enddo
enddo
!$OMP END DO
!$OMP END PARALLEL
tmp = 0.d0
!$OMP PARALLEL &
!$OMP DEFAULT (NONE) &
!$OMP PRIVATE (j, l, ipoint) &
!$OMP SHARED (tmp, ao_num, n_points_final_grid, u12sq_j1bsq_test, u12_grad1_u12_j1b_grad1_j1b_test, grad12_j12_test)
!$OMP DO SCHEDULE (static)
do ipoint = 1, n_points_final_grid
do j = 1, ao_num
do l = 1, ao_num
tmp(l,j,ipoint) = u12sq_j1bsq_test(l,j,ipoint) + u12_grad1_u12_j1b_grad1_j1b_test(l,j,ipoint) + 0.5d0 * grad12_j12_test(l,j,ipoint)
enddo
enddo
enddo
!$OMP END DO
!$OMP END PARALLEL
tc_grad_square_ao_test = 0.d0
call dgemm( "N", "N", ao_num*ao_num, ao_num*ao_num, n_points_final_grid, 1.d0 &
, tmp(1,1,1), ao_num*ao_num, b_mat(1,1,1), n_points_final_grid &
, 1.d0, tc_grad_square_ao_test, ao_num*ao_num)
deallocate(tmp, b_mat)
call sum_A_At(tc_grad_square_ao_test(1,1,1,1), ao_num*ao_num)
!do i = 1, ao_num
! do j = 1, ao_num
! do k = i, ao_num
! do l = max(j,k), ao_num
! tc_grad_square_ao_test(i,j,k,l) = 0.5d0 * (tc_grad_square_ao_test(i,j,k,l) + tc_grad_square_ao_test(k,l,i,j))
! tc_grad_square_ao_test(k,l,i,j) = tc_grad_square_ao_test(i,j,k,l)
! end do
! !if (j.eq.k) then
! ! do l = j+1, ao_num
! ! tc_grad_square_ao_test(i,j,k,l) = 0.5d0 * (tc_grad_square_ao_test(i,j,k,l) + tc_grad_square_ao_test(k,l,i,j))
! ! tc_grad_square_ao_test(k,l,i,j) = tc_grad_square_ao_test(i,j,k,l)
! ! end do
! !else
! ! do l = j, ao_num
! ! tc_grad_square_ao_test(i,j,k,l) = 0.5d0 * (tc_grad_square_ao_test(i,j,k,l) + tc_grad_square_ao_test(k,l,i,j))
! ! tc_grad_square_ao_test(k,l,i,j) = tc_grad_square_ao_test(i,j,k,l)
! ! enddo
! !endif
! enddo
! enddo
!enddo
!tc_grad_square_ao_test = 2.d0 * tc_grad_square_ao_test
! !$OMP PARALLEL &
! !$OMP DEFAULT (NONE) &
! !$OMP PRIVATE (i, j, k, l) &
! !$OMP SHARED (tc_grad_square_ao_test, ao_num)
! !$OMP DO SCHEDULE (static)
! integer :: ii
! ii = 0
! do j = 1, ao_num
! do l = 1, ao_num
! do i = 1, ao_num
! do k = 1, ao_num
! if((i.lt.j) .and. (k.lt.l)) cycle
! ii = ii + 1
! tc_grad_square_ao_test(k,i,l,j) = tc_grad_square_ao_test(k,i,l,j) + tc_grad_square_ao_test(l,j,k,i)
! enddo
! enddo
! enddo
! enddo
! print *, ' ii =', ii
! !$OMP END DO
! !$OMP END PARALLEL
! !$OMP PARALLEL &
! !$OMP DEFAULT (NONE) &
! !$OMP PRIVATE (i, j, k, l) &
! !$OMP SHARED (tc_grad_square_ao_test, ao_num)
! !$OMP DO SCHEDULE (static)
! do j = 1, ao_num
! do l = 1, ao_num
! do i = 1, j-1
! do k = 1, l-1
! ii = ii + 1
! tc_grad_square_ao_test(k,i,l,j) = tc_grad_square_ao_test(l,j,k,i)
! enddo
! enddo
! enddo
! enddo
! print *, ' ii =', ii
! print *, ao_num * ao_num * ao_num * ao_num
! !$OMP END DO
! !$OMP END PARALLEL
endif
if(write_tc_integ) then
open(unit=11, form="unformatted", file='tc_grad_square_ao_test', action="write")
do i = 1, ao_num
do j = 1, ao_num
do k = 1, ao_num
do l = 1, ao_num
write(11) tc_grad_square_ao_test(l,k,j,i)
enddo
enddo
enddo
enddo
close(11)
endif
call wall_time(time1)
print*, ' Wall time for tc_grad_square_ao_test = ', time1 - time0
END_PROVIDER
! ---
BEGIN_PROVIDER [double precision, tc_grad_square_ao_test_ref, (ao_num, ao_num, ao_num, ao_num)]
BEGIN_DOC
!
! tc_grad_square_ao_test_ref(k,i,l,j) = -1/2 <kl | |\grad_1 u(r1,r2)|^2 + |\grad_1 u(r1,r2)|^2 | ij>
!
END_DOC
implicit none
integer :: ipoint, i, j, k, l
double precision :: weight1, ao_ik_r, ao_i_r,contrib,contrib2
double precision :: time0, time1
double precision, allocatable :: ac_mat(:,:,:,:), b_mat(:,:,:), tmp(:,:,:)
print*, ' providing tc_grad_square_ao_test_ref ...'
call wall_time(time0)
provide u12sq_j1bsq_test u12_grad1_u12_j1b_grad1_j1b_test grad12_j12_test provide u12sq_j1bsq_test u12_grad1_u12_j1b_grad1_j1b_test grad12_j12_test
allocate(ac_mat(ao_num,ao_num,ao_num,ao_num), b_mat(n_points_final_grid,ao_num,ao_num), tmp(ao_num,ao_num,n_points_final_grid)) allocate(ac_mat(ao_num,ao_num,ao_num,ao_num), b_mat(n_points_final_grid,ao_num,ao_num), tmp(ao_num,ao_num,n_points_final_grid))
@ -61,13 +227,13 @@ BEGIN_PROVIDER [double precision, tc_grad_square_ao_test, (ao_num, ao_num, ao_nu
!$OMP PARALLEL & !$OMP PARALLEL &
!$OMP DEFAULT (NONE) & !$OMP DEFAULT (NONE) &
!$OMP PRIVATE (i, j, k, l) & !$OMP PRIVATE (i, j, k, l) &
!$OMP SHARED (ac_mat, tc_grad_square_ao_test, ao_num) !$OMP SHARED (ac_mat, tc_grad_square_ao_test_ref, ao_num)
!$OMP DO SCHEDULE (static) !$OMP DO SCHEDULE (static)
do j = 1, ao_num do j = 1, ao_num
do l = 1, ao_num do l = 1, ao_num
do i = 1, ao_num do i = 1, ao_num
do k = 1, ao_num do k = 1, ao_num
tc_grad_square_ao_test(k,i,l,j) = ac_mat(k,i,l,j) + ac_mat(l,j,k,i) tc_grad_square_ao_test_ref(k,i,l,j) = ac_mat(k,i,l,j) + ac_mat(l,j,k,i)
enddo enddo
enddo enddo
enddo enddo
@ -78,7 +244,7 @@ BEGIN_PROVIDER [double precision, tc_grad_square_ao_test, (ao_num, ao_num, ao_nu
deallocate(ac_mat) deallocate(ac_mat)
call wall_time(time1) call wall_time(time1)
print*, ' Wall time for tc_grad_square_ao_test = ', time1 - time0 print*, ' Wall time for tc_grad_square_ao_test_ref = ', time1 - time0
END_PROVIDER END_PROVIDER

116
src/non_h_ints_mu/new_grad_tc.irp.f
View File

@ -25,7 +25,7 @@ BEGIN_PROVIDER [ double precision, int2_grad1_u12_ao, (ao_num, ao_num, n_points_
END_DOC END_DOC
implicit none implicit none
integer :: ipoint, i, j integer :: ipoint, i, j, m
double precision :: time0, time1 double precision :: time0, time1
double precision :: x, y, z, tmp_x, tmp_y, tmp_z, tmp0, tmp1, tmp2 double precision :: x, y, z, tmp_x, tmp_y, tmp_z, tmp0, tmp1, tmp2
@ -34,38 +34,46 @@ BEGIN_PROVIDER [ double precision, int2_grad1_u12_ao, (ao_num, ao_num, n_points_
PROVIDE j1b_type PROVIDE j1b_type
if(j1b_type .eq. 3) then if(read_tc_integ) then
open(unit=11, form="unformatted", file='int2_grad1_u12_ao', action="read")
do m = 1, 3
do ipoint = 1, n_points_final_grid
do j = 1, ao_num
do i = 1, ao_num
read(11) int2_grad1_u12_ao(i,j,ipoint,m)
enddo
enddo
enddo
enddo
close(11)
else
if(j1b_type .eq. 3) then
do ipoint = 1, n_points_final_grid do ipoint = 1, n_points_final_grid
x = final_grid_points(1,ipoint) x = final_grid_points(1,ipoint)
y = final_grid_points(2,ipoint) y = final_grid_points(2,ipoint)
z = final_grid_points(3,ipoint) z = final_grid_points(3,ipoint)
tmp0 = 0.5d0 * v_1b(ipoint) tmp0 = 0.5d0 * v_1b(ipoint)
tmp_x = v_1b_grad(1,ipoint) tmp_x = v_1b_grad(1,ipoint)
tmp_y = v_1b_grad(2,ipoint) tmp_y = v_1b_grad(2,ipoint)
tmp_z = v_1b_grad(3,ipoint) tmp_z = v_1b_grad(3,ipoint)
do j = 1, ao_num do j = 1, ao_num
do i = 1, ao_num do i = 1, ao_num
tmp1 = tmp0 * v_ij_erf_rk_cst_mu_j1b(i,j,ipoint) tmp1 = tmp0 * v_ij_erf_rk_cst_mu_j1b(i,j,ipoint)
tmp2 = v_ij_u_cst_mu_j1b(i,j,ipoint) tmp2 = v_ij_u_cst_mu_j1b(i,j,ipoint)
int2_grad1_u12_ao(i,j,ipoint,1) = tmp1 * x - tmp0 * x_v_ij_erf_rk_cst_mu_j1b(i,j,ipoint,1) - tmp2 * tmp_x int2_grad1_u12_ao(i,j,ipoint,1) = tmp1 * x - tmp0 * x_v_ij_erf_rk_cst_mu_j1b(i,j,ipoint,1) - tmp2 * tmp_x
int2_grad1_u12_ao(i,j,ipoint,2) = tmp1 * y - tmp0 * x_v_ij_erf_rk_cst_mu_j1b(i,j,ipoint,2) - tmp2 * tmp_y int2_grad1_u12_ao(i,j,ipoint,2) = tmp1 * y - tmp0 * x_v_ij_erf_rk_cst_mu_j1b(i,j,ipoint,2) - tmp2 * tmp_y
int2_grad1_u12_ao(i,j,ipoint,3) = tmp1 * z - tmp0 * x_v_ij_erf_rk_cst_mu_j1b(i,j,ipoint,3) - tmp2 * tmp_z int2_grad1_u12_ao(i,j,ipoint,3) = tmp1 * z - tmp0 * x_v_ij_erf_rk_cst_mu_j1b(i,j,ipoint,3) - tmp2 * tmp_z
enddo enddo
enddo enddo
enddo enddo
else else
do ipoint = 1, n_points_final_grid do ipoint = 1, n_points_final_grid
x = final_grid_points(1,ipoint) x = final_grid_points(1,ipoint)
y = final_grid_points(2,ipoint) y = final_grid_points(2,ipoint)
z = final_grid_points(3,ipoint) z = final_grid_points(3,ipoint)
do j = 1, ao_num do j = 1, ao_num
do i = 1, ao_num do i = 1, ao_num
tmp1 = v_ij_erf_rk_cst_mu(i,j,ipoint) tmp1 = v_ij_erf_rk_cst_mu(i,j,ipoint)
@ -76,11 +84,25 @@ BEGIN_PROVIDER [ double precision, int2_grad1_u12_ao, (ao_num, ao_num, n_points_
enddo enddo
enddo enddo
enddo enddo
int2_grad1_u12_ao *= 0.5d0 int2_grad1_u12_ao *= 0.5d0
endif
endif endif
if(write_tc_integ) then
open(unit=11, form="unformatted", file='int2_grad1_u12_ao', action="write")
do m = 1, 3
do ipoint = 1, n_points_final_grid
do j = 1, ao_num
do i = 1, ao_num
write(11) int2_grad1_u12_ao(i,j,ipoint,m)
enddo
enddo
enddo
enddo
close(11)
endif
call wall_time(time1) call wall_time(time1)
print*, ' Wall time for int2_grad1_u12_ao = ', time1 - time0 print*, ' Wall time for int2_grad1_u12_ao = ', time1 - time0
@ -290,12 +312,28 @@ BEGIN_PROVIDER [double precision, tc_grad_and_lapl_ao, (ao_num, ao_num, ao_num,
integer :: ipoint, i, j, k, l, m integer :: ipoint, i, j, k, l, m
double precision :: weight1, ao_k_r, ao_i_r double precision :: weight1, ao_k_r, ao_i_r
double precision :: time0, time1 double precision :: time0, time1
double precision, allocatable :: ac_mat(:,:,:,:), b_mat(:,:,:,:) double precision, allocatable :: b_mat(:,:,:,:)
print*, ' providing tc_grad_and_lapl_ao ...' print*, ' providing tc_grad_and_lapl_ao ...'
call wall_time(time0) call wall_time(time0)
allocate(b_mat(n_points_final_grid,ao_num,ao_num,3), ac_mat(ao_num,ao_num,ao_num,ao_num)) if(read_tc_integ) then
open(unit=11, form="unformatted", file='tc_grad_and_lapl_ao', action="read")
do i = 1, ao_num
do j = 1, ao_num
do k = 1, ao_num
do l = 1, ao_num
read(11) tc_grad_and_lapl_ao(l,k,j,i)
enddo
enddo
enddo
enddo
close(11)
else
allocate(b_mat(n_points_final_grid,ao_num,ao_num,3))
b_mat = 0.d0 b_mat = 0.d0
!$OMP PARALLEL & !$OMP PARALLEL &
@ -321,34 +359,48 @@ BEGIN_PROVIDER [double precision, tc_grad_and_lapl_ao, (ao_num, ao_num, ao_num,
!$OMP END DO !$OMP END DO
!$OMP END PARALLEL !$OMP END PARALLEL
ac_mat = 0.d0 tc_grad_and_lapl_ao = 0.d0
do m = 1, 3 do m = 1, 3
call dgemm( "N", "N", ao_num*ao_num, ao_num*ao_num, n_points_final_grid, 1.d0 & call dgemm( "N", "N", ao_num*ao_num, ao_num*ao_num, n_points_final_grid, 1.d0 &
, int2_grad1_u12_ao(1,1,1,m), ao_num*ao_num, b_mat(1,1,1,m), n_points_final_grid & , int2_grad1_u12_ao(1,1,1,m), ao_num*ao_num, b_mat(1,1,1,m), n_points_final_grid &
, 1.d0, ac_mat, ao_num*ao_num) , 1.d0, tc_grad_and_lapl_ao, ao_num*ao_num)
enddo enddo
deallocate(b_mat) deallocate(b_mat)
!$OMP PARALLEL & call sum_A_At(tc_grad_and_lapl_ao(1,1,1,1), ao_num*ao_num)
!$OMP DEFAULT (NONE) & ! !$OMP PARALLEL &
!$OMP PRIVATE (i, j, k, l) & ! !$OMP DEFAULT (NONE) &
!$OMP SHARED (ac_mat, tc_grad_and_lapl_ao, ao_num) ! !$OMP PRIVATE (i, j, k, l) &
!$OMP DO SCHEDULE (static) ! !$OMP SHARED (ac_mat, tc_grad_and_lapl_ao, ao_num)
do j = 1, ao_num ! !$OMP DO SCHEDULE (static)
do l = 1, ao_num ! do j = 1, ao_num
do i = 1, ao_num ! do l = 1, ao_num
do k = 1, ao_num ! do i = 1, ao_num
tc_grad_and_lapl_ao(k,i,l,j) = ac_mat(k,i,l,j) + ac_mat(l,j,k,i) ! do k = 1, ao_num
!tc_grad_and_lapl_ao(k,i,l,j) = ac_mat(k,i,l,j) ! tc_grad_and_lapl_ao(k,i,l,j) = ac_mat(k,i,l,j) + ac_mat(l,j,k,i)
enddo ! enddo
enddo ! enddo
enddo ! enddo
enddo ! enddo
!$OMP END DO ! !$OMP END DO
!$OMP END PARALLEL ! !$OMP END PARALLEL
deallocate(ac_mat) endif
if(write_tc_integ) then
open(unit=11, form="unformatted", file='tc_grad_and_lapl_ao', action="write")
do i = 1, ao_num
do j = 1, ao_num
do k = 1, ao_num
do l = 1, ao_num
write(11) tc_grad_and_lapl_ao(l,k,j,i)
enddo
enddo
enddo
enddo
close(11)
endif
call wall_time(time1) call wall_time(time1)
print*, ' Wall time for tc_grad_and_lapl_ao = ', time1 - time0 print*, ' Wall time for tc_grad_and_lapl_ao = ', time1 - time0

75
src/non_h_ints_mu/new_grad_tc_manu.irp.f
View File

@ -24,7 +24,7 @@ BEGIN_PROVIDER [ double precision, int2_grad1_u12_ao_test, (ao_num, ao_num, n_po
END_DOC END_DOC
implicit none implicit none
integer :: ipoint, i, j integer :: ipoint, i, j, m
double precision :: time0, time1 double precision :: time0, time1
double precision :: x, y, z, tmp_x, tmp_y, tmp_z, tmp0, tmp1, tmp2 double precision :: x, y, z, tmp_x, tmp_y, tmp_z, tmp0, tmp1, tmp2
@ -33,53 +33,74 @@ BEGIN_PROVIDER [ double precision, int2_grad1_u12_ao_test, (ao_num, ao_num, n_po
PROVIDE j1b_type PROVIDE j1b_type
if(j1b_type .eq. 3) then if(read_tc_integ) then
open(unit=11, form="unformatted", file='int2_grad1_u12_ao_test', action="read")
do m = 1, 3
do ipoint = 1, n_points_final_grid
do j = 1, ao_num
do i = 1, ao_num
read(11) int2_grad1_u12_ao_test(i,j,ipoint,m)
enddo
enddo
enddo
enddo
close(11)
else
if(j1b_type .eq. 3) then
do ipoint = 1, n_points_final_grid do ipoint = 1, n_points_final_grid
x = final_grid_points(1,ipoint) x = final_grid_points(1,ipoint)
y = final_grid_points(2,ipoint) y = final_grid_points(2,ipoint)
z = final_grid_points(3,ipoint) z = final_grid_points(3,ipoint)
tmp0 = 0.5d0 * v_1b(ipoint) tmp0 = 0.5d0 * v_1b(ipoint)
tmp_x = v_1b_grad(1,ipoint) tmp_x = v_1b_grad(1,ipoint)
tmp_y = v_1b_grad(2,ipoint) tmp_y = v_1b_grad(2,ipoint)
tmp_z = v_1b_grad(3,ipoint) tmp_z = v_1b_grad(3,ipoint)
do j = 1, ao_num do j = 1, ao_num
do i = 1, ao_num do i = 1, ao_num
tmp1 = tmp0 * v_ij_erf_rk_cst_mu_j1b_test(i,j,ipoint) tmp1 = tmp0 * v_ij_erf_rk_cst_mu_j1b_test(i,j,ipoint)
tmp2 = v_ij_u_cst_mu_j1b_test(i,j,ipoint) tmp2 = v_ij_u_cst_mu_j1b_test(i,j,ipoint)
int2_grad1_u12_ao_test(i,j,ipoint,1) = tmp1 * x - tmp0 * x_v_ij_erf_rk_cst_mu_j1b_test(i,j,ipoint,1) - tmp2 * tmp_x int2_grad1_u12_ao_test(i,j,ipoint,1) = tmp1 * x - tmp0 * x_v_ij_erf_rk_cst_mu_j1b_test(i,j,ipoint,1) - tmp2 * tmp_x
int2_grad1_u12_ao_test(i,j,ipoint,2) = tmp1 * y - tmp0 * x_v_ij_erf_rk_cst_mu_j1b_test(i,j,ipoint,2) - tmp2 * tmp_y int2_grad1_u12_ao_test(i,j,ipoint,2) = tmp1 * y - tmp0 * x_v_ij_erf_rk_cst_mu_j1b_test(i,j,ipoint,2) - tmp2 * tmp_y
int2_grad1_u12_ao_test(i,j,ipoint,3) = tmp1 * z - tmp0 * x_v_ij_erf_rk_cst_mu_j1b_test(i,j,ipoint,3) - tmp2 * tmp_z int2_grad1_u12_ao_test(i,j,ipoint,3) = tmp1 * z - tmp0 * x_v_ij_erf_rk_cst_mu_j1b_test(i,j,ipoint,3) - tmp2 * tmp_z
enddo enddo
enddo enddo
enddo enddo
else else
do ipoint = 1, n_points_final_grid do ipoint = 1, n_points_final_grid
x = final_grid_points(1,ipoint) x = final_grid_points(1,ipoint)
y = final_grid_points(2,ipoint) y = final_grid_points(2,ipoint)
z = final_grid_points(3,ipoint) z = final_grid_points(3,ipoint)
do j = 1, ao_num do j = 1, ao_num
do i = 1, ao_num do i = 1, ao_num
tmp1 = v_ij_erf_rk_cst_mu(i,j,ipoint) tmp1 = v_ij_erf_rk_cst_mu(i,j,ipoint)
int2_grad1_u12_ao_test(i,j,ipoint,1) = tmp1 * x - x_v_ij_erf_rk_cst_mu_tmp(i,j,ipoint,1) int2_grad1_u12_ao_test(i,j,ipoint,1) = tmp1 * x - x_v_ij_erf_rk_cst_mu_tmp(i,j,ipoint,1)
int2_grad1_u12_ao_test(i,j,ipoint,2) = tmp1 * y - x_v_ij_erf_rk_cst_mu_tmp(i,j,ipoint,2) int2_grad1_u12_ao_test(i,j,ipoint,2) = tmp1 * y - x_v_ij_erf_rk_cst_mu_tmp(i,j,ipoint,2)
int2_grad1_u12_ao_test(i,j,ipoint,3) = tmp1 * z - x_v_ij_erf_rk_cst_mu_tmp(i,j,ipoint,3) int2_grad1_u12_ao_test(i,j,ipoint,3) = tmp1 * z - x_v_ij_erf_rk_cst_mu_tmp(i,j,ipoint,3)
enddo enddo
enddo enddo
enddo enddo
int2_grad1_u12_ao_test *= 0.5d0 int2_grad1_u12_ao_test *= 0.5d0
endif
endif endif
if(write_tc_integ) then
open(unit=11, form="unformatted", file='int2_grad1_u12_ao_test', action="write")
do m = 1, 3
do ipoint = 1, n_points_final_grid
do j = 1, ao_num
do i = 1, ao_num
write(11) int2_grad1_u12_ao_test(i,j,ipoint,m)
enddo
enddo
enddo
enddo
close(11)
endif
call wall_time(time1) call wall_time(time1)
print*, ' Wall time for int2_grad1_u12_ao_test = ', time1 - time0 print*, ' Wall time for int2_grad1_u12_ao_test = ', time1 - time0
@ -109,6 +130,22 @@ BEGIN_PROVIDER [double precision, tc_grad_and_lapl_ao_test, (ao_num, ao_num, ao_
print*, ' providing tc_grad_and_lapl_ao_test ...' print*, ' providing tc_grad_and_lapl_ao_test ...'
call wall_time(time0) call wall_time(time0)
if(read_tc_integ) then
open(unit=11, form="unformatted", file='tc_grad_and_lapl_ao_test', action="read")
do i = 1, ao_num
do j = 1, ao_num
do k = 1, ao_num
do l = 1, ao_num
read(11) tc_grad_and_lapl_ao_test(l,k,j,i)
enddo
enddo
enddo
enddo
close(11)
else
provide int2_grad1_u12_ao_test provide int2_grad1_u12_ao_test
allocate(b_mat(n_points_final_grid,ao_num,ao_num,3), ac_mat(ao_num,ao_num,ao_num,ao_num)) allocate(b_mat(n_points_final_grid,ao_num,ao_num,3), ac_mat(ao_num,ao_num,ao_num,ao_num))
@ -165,6 +202,22 @@ BEGIN_PROVIDER [double precision, tc_grad_and_lapl_ao_test, (ao_num, ao_num, ao_
deallocate(ac_mat) deallocate(ac_mat)
endif
if(write_tc_integ) then
open(unit=11, form="unformatted", file='tc_grad_and_lapl_ao_test', action="write")
do i = 1, ao_num
do j = 1, ao_num
do k = 1, ao_num
do l = 1, ao_num
write(11) tc_grad_and_lapl_ao_test(l,k,j,i)
enddo
enddo
enddo
enddo
close(11)
endif
call wall_time(time1) call wall_time(time1)
print*, ' Wall time for tc_grad_and_lapl_ao_test = ', time1 - time0 print*, ' Wall time for tc_grad_and_lapl_ao_test = ', time1 - time0

38
src/non_h_ints_mu/total_tc_int.irp.f
View File

@ -1,6 +1,44 @@
! --- ! ---
BEGIN_PROVIDER [double precision, ao_vartc_int_chemist, (ao_num, ao_num, ao_num, ao_num)]
implicit none
integer :: i, j, k, l
double precision :: wall1, wall0
print *, ' providing ao_vartc_int_chemist ...'
call wall_time(wall0)
if(test_cycle_tc) then
do j = 1, ao_num
do l = 1, ao_num
do i = 1, ao_num
do k = 1, ao_num
ao_vartc_int_chemist(k,i,l,j) = tc_grad_square_ao_test(k,i,l,j) + ao_two_e_coul(k,i,l,j)
enddo
enddo
enddo
enddo
else
do j = 1, ao_num
do l = 1, ao_num
do i = 1, ao_num
do k = 1, ao_num
ao_vartc_int_chemist(k,i,l,j) = tc_grad_square_ao(k,i,l,j) + ao_two_e_coul(k,i,l,j)
enddo
enddo
enddo
enddo
endif
call wall_time(wall1)
print *, ' wall time for ao_vartc_int_chemist ', wall1 - wall0
END_PROVIDER
! ---
BEGIN_PROVIDER [double precision, ao_tc_int_chemist, (ao_num, ao_num, ao_num, ao_num)] BEGIN_PROVIDER [double precision, ao_tc_int_chemist, (ao_num, ao_num, ao_num, ao_num)]
implicit none implicit none