From 3a4b41257f63398c26527b3e0fbb5dff008fe031 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Wed, 2 Dec 2015 00:32:04 +0100 Subject: [PATCH] Transposed arrays for pseudos --- src/Integrals_Monoelec/pot_ao_ints.irp.f | 56 +++++++------ .../pot_ao_pseudo_ints.irp.f | 80 ++++++++++++------- 2 files changed, 76 insertions(+), 60 deletions(-) diff --git a/src/Integrals_Monoelec/pot_ao_ints.irp.f b/src/Integrals_Monoelec/pot_ao_ints.irp.f index 09b67425..eadc0b72 100644 --- a/src/Integrals_Monoelec/pot_ao_ints.irp.f +++ b/src/Integrals_Monoelec/pot_ao_ints.irp.f @@ -26,7 +26,7 @@ n_pt_in = n_pt_max_integrals - !$OMP DO SCHEDULE (guided) + !$OMP DO SCHEDULE (dynamic) do j = 1, ao_num num_A = ao_nucl(j) @@ -81,23 +81,17 @@ integer :: power_A(3),power_B(3) integer :: i,j,k,l,n_pt_in,m double precision ::overlap_x,overlap_y,overlap_z,overlap,dx,NAI_pol_mult - ! Important for OpenMP ao_nucl_elec_integral_per_atom = 0.d0 - - do k = 1, nucl_num - C_center(1) = nucl_coord(k,1) - C_center(2) = nucl_coord(k,2) - C_center(3) = nucl_coord(k,3) !$OMP PARALLEL & !$OMP DEFAULT (NONE) & - !$OMP PRIVATE (i,j,l,m,alpha,beta,A_center,B_center,power_A,power_B, & - !$OMP num_A,num_B,c,n_pt_in) & - !$OMP SHARED (k,ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp, & - !$OMP n_pt_max_integrals,ao_nucl_elec_integral_per_atom,nucl_num,C_center) + !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,power_A,power_B, & + !$OMP num_A,num_B,c,n_pt_in,C_center) & + !$OMP SHARED (ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp, & + !$OMP n_pt_max_integrals,ao_nucl_elec_integral_per_atom,nucl_num) n_pt_in = n_pt_max_integrals - !$OMP DO SCHEDULE (guided) + !$OMP DO SCHEDULE (dynamic) double precision :: c do j = 1, ao_num @@ -108,29 +102,33 @@ A_center(1) = nucl_coord(num_A,1) A_center(2) = nucl_coord(num_A,2) A_center(3) = nucl_coord(num_A,3) - do i = 1, ao_num - power_B(1)= ao_power(i,1) - power_B(2)= ao_power(i,2) - power_B(3)= ao_power(i,3) - num_B = ao_nucl(i) - B_center(1) = nucl_coord(num_B,1) - B_center(2) = nucl_coord(num_B,2) - B_center(3) = nucl_coord(num_B,3) - c = 0.d0 - do l=1,ao_prim_num(j) - alpha = ao_expo_ordered_transp(l,j) - do m=1,ao_prim_num(i) - beta = ao_expo_ordered_transp(m,i) - c = c + NAI_pol_mult(A_center,B_center,power_A,power_B,alpha,beta,C_center,n_pt_in) & - * ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i) + do k = 1, nucl_num + C_center(1) = nucl_coord(k,1) + C_center(2) = nucl_coord(k,2) + C_center(3) = nucl_coord(k,3) + do i = 1, ao_num + power_B(1)= ao_power(i,1) + power_B(2)= ao_power(i,2) + power_B(3)= ao_power(i,3) + num_B = ao_nucl(i) + B_center(1) = nucl_coord(num_B,1) + B_center(2) = nucl_coord(num_B,2) + B_center(3) = nucl_coord(num_B,3) + c = 0.d0 + do l=1,ao_prim_num(j) + alpha = ao_expo_ordered_transp(l,j) + do m=1,ao_prim_num(i) + beta = ao_expo_ordered_transp(m,i) + c = c + NAI_pol_mult(A_center,B_center,power_A,power_B,alpha,beta,C_center,n_pt_in) & + * ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i) + enddo enddo + ao_nucl_elec_integral_per_atom(i,j,k) = -c enddo - ao_nucl_elec_integral_per_atom(i,j,k) = -c enddo enddo !$OMP END DO !$OMP END PARALLEL - enddo END_PROVIDER diff --git a/src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f b/src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f index e18bc006..2856807b 100644 --- a/src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f +++ b/src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f @@ -28,12 +28,6 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu ao_pseudo_integral_local = 0.d0 - !! Dump array - integer, allocatable :: n_k_dump(:) - double precision, allocatable :: v_k_dump(:), dz_k_dump(:) - - allocate(n_k_dump(1:pseudo_klocmax), v_k_dump(1:pseudo_klocmax), dz_k_dump(1:pseudo_klocmax)) - print*, 'Providing the nuclear electron pseudo integrals (local)' call wall_time(wall_1) @@ -44,11 +38,10 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu !$OMP DEFAULT (NONE) & !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B,& !$OMP num_A,num_B,Z,c,n_pt_in, & - !$OMP v_k_dump,n_k_dump, dz_k_dump, & !$OMP wall_0,wall_2,thread_num) & !$OMP SHARED (ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp,& !$OMP ao_pseudo_integral_local,nucl_num,nucl_charge, & - !$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_v_k,pseudo_n_k, pseudo_dz_k,& + !$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_v_k_transp,pseudo_n_k_transp, pseudo_dz_k_transp,& !$OMP wall_1) !$ thread_num = omp_get_thread_num() @@ -84,11 +77,10 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu C_center(1:3) = nucl_coord(k,1:3) - v_k_dump = pseudo_v_k(k,1:pseudo_klocmax) - n_k_dump = pseudo_n_k(k,1:pseudo_klocmax) - dz_k_dump = pseudo_dz_k(k,1:pseudo_klocmax) - - c = c + Vloc(pseudo_klocmax, v_k_dump,n_k_dump, dz_k_dump,& + c = c + Vloc(pseudo_klocmax, & + pseudo_v_k_transp (1,k), & + pseudo_n_k_transp (1,k), & + pseudo_dz_k_transp(1,k), & A_center,power_A,alpha,B_center,power_B,beta,C_center) enddo @@ -112,8 +104,6 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu !$OMP END PARALLEL - deallocate(n_k_dump,v_k_dump, dz_k_dump) - END_PROVIDER @@ -135,12 +125,6 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu ao_pseudo_integral_non_local = 0.d0 - !! Dump array - integer, allocatable :: n_kl_dump(:,:) - double precision, allocatable :: v_kl_dump(:,:), dz_kl_dump(:,:) - - allocate(n_kl_dump(pseudo_kmax,0:pseudo_lmax), v_kl_dump(pseudo_kmax,0:pseudo_lmax), dz_kl_dump(pseudo_kmax,0:pseudo_lmax)) - print*, 'Providing the nuclear electron pseudo integrals (non-local)' call wall_time(wall_1) @@ -150,14 +134,14 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu !$OMP DEFAULT (NONE) & !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B,& !$OMP num_A,num_B,Z,c,n_pt_in, & - !$OMP n_kl_dump, v_kl_dump, dz_kl_dump, & !$OMP wall_0,wall_2,thread_num) & !$OMP SHARED (ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp,& !$OMP ao_pseudo_integral_non_local,nucl_num,nucl_charge,& - !$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_n_kl, pseudo_v_kl, pseudo_dz_kl,& + !$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_n_kl_transp, pseudo_v_kl_transp, pseudo_dz_kl_transp,& !$OMP wall_1) !$ thread_num = omp_get_thread_num() + !$OMP DO SCHEDULE (guided) do j = 1, ao_num @@ -191,12 +175,11 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu C_center(1:3) = nucl_coord(k,1:3) - n_kl_dump = pseudo_n_kl(k,1:pseudo_kmax,0:pseudo_lmax) - v_kl_dump = pseudo_v_kl(k,1:pseudo_kmax,0:pseudo_lmax) - dz_kl_dump = pseudo_dz_kl(k,1:pseudo_kmax,0:pseudo_lmax) - - c = c + Vpseudo(pseudo_lmax,pseudo_kmax,v_kl_dump,n_kl_dump,dz_kl_dump,A_center,power_A,alpha,B_center,power_B,beta,C_center) - + c = c + Vpseudo(pseudo_lmax,pseudo_kmax, & + pseudo_v_kl_transp(1,0,k), & + pseudo_n_kl_transp(1,0,k), & + pseudo_dz_kl_transp(1,0,k), & + A_center,power_A,alpha,B_center,power_B,beta,C_center) enddo ao_pseudo_integral_non_local(i,j) = ao_pseudo_integral_non_local(i,j) +& ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i)*c @@ -215,13 +198,48 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu enddo !$OMP END DO + !$OMP END PARALLEL - deallocate(n_kl_dump,v_kl_dump, dz_kl_dump) - END_PROVIDER + BEGIN_PROVIDER [ double precision, pseudo_v_k_transp, (pseudo_klocmax,nucl_num) ] +&BEGIN_PROVIDER [ integer , pseudo_n_k_transp, (pseudo_klocmax,nucl_num) ] +&BEGIN_PROVIDER [ double precision, pseudo_dz_k_transp, (pseudo_klocmax,nucl_num)] + implicit none + BEGIN_DOC + ! Transposed arrays for pseudopotentials + END_DOC + integer :: i,j + do j=1,nucl_num + do i=1,pseudo_klocmax + pseudo_v_k_transp (i,j) = pseudo_v_k (j,i) + pseudo_n_k_transp (i,j) = pseudo_n_k (j,i) + pseudo_dz_k_transp(i,j) = pseudo_dz_k(j,i) + enddo + enddo +END_PROVIDER + + BEGIN_PROVIDER [ double precision, pseudo_v_kl_transp, (pseudo_kmax,0:pseudo_lmax,nucl_num) ] +&BEGIN_PROVIDER [ integer , pseudo_n_kl_transp, (pseudo_kmax,0:pseudo_lmax,nucl_num) ] +&BEGIN_PROVIDER [ double precision, pseudo_dz_kl_transp, (pseudo_kmax,0:pseudo_lmax,nucl_num)] + implicit none + BEGIN_DOC + ! Transposed arrays for pseudopotentials + END_DOC + + integer :: i,j,l + do j=1,nucl_num + do l=0,pseudo_lmax + do i=1,pseudo_kmax + pseudo_v_kl_transp (i,l,j) = pseudo_v_kl (j,i,l) + pseudo_n_kl_transp (i,l,j) = pseudo_n_kl (j,i,l) + pseudo_dz_kl_transp(i,l,j) = pseudo_dz_kl(j,i,l) + enddo + enddo + enddo +END_PROVIDER