10
0
mirror of https://github.com/LCPQ/quantum_package synced 2024-12-25 05:43:47 +01:00

Fixed ao_overlap ezfio problem

This commit is contained in:
Anthony Scemama 2018-10-17 16:16:05 +02:00
parent 6566529ec4
commit 75790e6972
2 changed files with 9 additions and 18 deletions

View File

@ -18,7 +18,7 @@ END_PROVIDER
pt2_n_tasks_max = 1 + min((e*(e-1))/2, int(dsqrt(dble(N_det_generators)))/10) pt2_n_tasks_max = 1 + min((e*(e-1))/2, int(dsqrt(dble(N_det_generators)))/10)
pt2_F(:) = 1 pt2_F(:) = 1
do i=1,min(10000,N_det_generators) do i=1,min(10000,N_det_generators)
pt2_F(i) = 1 + dble(pt2_n_tasks_max)*maxval(dsqrt(dabs(psi_coef_sorted_gen(i,1:N_states)))) pt2_F(i) = 1 + int(dble(pt2_n_tasks_max)*maxval(dsqrt(dabs(psi_coef_sorted_gen(i,1:N_states)))))
enddo enddo
if(N_det_generators < 128) then if(N_det_generators < 128) then

View File

@ -1,4 +1,4 @@
BEGIN_PROVIDER [ double precision, ao_overlap,(ao_num,ao_num) ] BEGIN_PROVIDER [ double precision, ao_overlap_matrix,(ao_num,ao_num) ]
&BEGIN_PROVIDER [ double precision, ao_overlap_x,(ao_num,ao_num) ] &BEGIN_PROVIDER [ double precision, ao_overlap_x,(ao_num,ao_num) ]
&BEGIN_PROVIDER [ double precision, ao_overlap_y,(ao_num,ao_num) ] &BEGIN_PROVIDER [ double precision, ao_overlap_y,(ao_num,ao_num) ]
&BEGIN_PROVIDER [ double precision, ao_overlap_z,(ao_num,ao_num) ] &BEGIN_PROVIDER [ double precision, ao_overlap_z,(ao_num,ao_num) ]
@ -14,10 +14,6 @@
double precision :: alpha, beta, c double precision :: alpha, beta, c
double precision :: A_center(3), B_center(3) double precision :: A_center(3), B_center(3)
integer :: power_A(3), power_B(3) integer :: power_A(3), power_B(3)
if (read_ao_one_integrals) then
call ezfio_get_ao_basis_integral_overlap(ao_overlap(1:ao_num, 1:ao_num))
print *, 'AO overlap integrals read from disk'
else
dim1=100 dim1=100
!$OMP PARALLEL DO SCHEDULE(GUIDED) & !$OMP PARALLEL DO SCHEDULE(GUIDED) &
!$OMP DEFAULT(NONE) & !$OMP DEFAULT(NONE) &
@ -25,7 +21,7 @@
!$OMP overlap_x,overlap_y, overlap_z, overlap, & !$OMP overlap_x,overlap_y, overlap_z, overlap, &
!$OMP alpha, beta,i,j,c) & !$OMP alpha, beta,i,j,c) &
!$OMP SHARED(nucl_coord,ao_power,ao_prim_num, & !$OMP SHARED(nucl_coord,ao_power,ao_prim_num, &
!$OMP ao_overlap_x,ao_overlap_y,ao_overlap_z,ao_overlap,ao_num,ao_coef_normalized_ordered_transp,ao_nucl, & !$OMP ao_overlap_x,ao_overlap_y,ao_overlap_z,ao_overlap_matrix,ao_num,ao_coef_normalized_ordered_transp,ao_nucl, &
!$OMP ao_expo_ordered_transp,dim1) !$OMP ao_expo_ordered_transp,dim1)
do j=1,ao_num do j=1,ao_num
A_center(1) = nucl_coord( ao_nucl(j), 1 ) A_center(1) = nucl_coord( ao_nucl(j), 1 )
@ -35,7 +31,7 @@
power_A(2) = ao_power( j, 2 ) power_A(2) = ao_power( j, 2 )
power_A(3) = ao_power( j, 3 ) power_A(3) = ao_power( j, 3 )
do i= 1,ao_num do i= 1,ao_num
ao_overlap(i,j)= 0.d0 ao_overlap_matrix(i,j)= 0.d0
ao_overlap_x(i,j)= 0.d0 ao_overlap_x(i,j)= 0.d0
ao_overlap_y(i,j)= 0.d0 ao_overlap_y(i,j)= 0.d0
ao_overlap_z(i,j)= 0.d0 ao_overlap_z(i,j)= 0.d0
@ -51,7 +47,7 @@
beta = ao_expo_ordered_transp(l,i) beta = ao_expo_ordered_transp(l,i)
call overlap_gaussian_xyz(A_center,B_center,alpha,beta,power_A,power_B,overlap_x,overlap_y,overlap_z,overlap,dim1) call overlap_gaussian_xyz(A_center,B_center,alpha,beta,power_A,power_B,overlap_x,overlap_y,overlap_z,overlap,dim1)
c = ao_coef_normalized_ordered_transp(n,j) * ao_coef_normalized_ordered_transp(l,i) c = ao_coef_normalized_ordered_transp(n,j) * ao_coef_normalized_ordered_transp(l,i)
ao_overlap(i,j) += c * overlap ao_overlap_matrix(i,j) += c * overlap
ao_overlap_x(i,j) += c * overlap_x ao_overlap_x(i,j) += c * overlap_x
ao_overlap_y(i,j) += c * overlap_y ao_overlap_y(i,j) += c * overlap_y
ao_overlap_z(i,j) += c * overlap_z ao_overlap_z(i,j) += c * overlap_z
@ -60,11 +56,6 @@
enddo enddo
enddo enddo
!$OMP END PARALLEL DO !$OMP END PARALLEL DO
endif
if (write_ao_one_integrals) then
call ezfio_set_ao_basis_integral_overlap(ao_overlap(1:ao_num, 1:ao_num))
print *, 'AO overlap integrals written to disk'
endif
END_PROVIDER END_PROVIDER
@ -128,7 +119,7 @@ BEGIN_PROVIDER [ double precision, S_inv,(ao_num,ao_num) ]
BEGIN_DOC BEGIN_DOC
! S^-1 ! S^-1
END_DOC END_DOC
call get_pseudo_inverse(ao_overlap,size(ao_overlap,1),ao_num,ao_num,S_inv,size(S_inv,1)) call get_pseudo_inverse(ao_overlap_matrix,size(ao_overlap_matrix,1),ao_num,ao_num,S_inv,size(S_inv,1))
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ double precision, S_half_inv, (AO_num,AO_num) ] BEGIN_PROVIDER [ double precision, S_half_inv, (AO_num,AO_num) ]
@ -145,7 +136,7 @@ BEGIN_PROVIDER [ double precision, S_half_inv, (AO_num,AO_num) ]
integer :: info, i, j, k integer :: info, i, j, k
double precision, parameter :: threshold_overlap_AO_eigenvalues = 1.d-6 double precision, parameter :: threshold_overlap_AO_eigenvalues = 1.d-6
LDA = size(AO_overlap,1) LDA = size(AO_overlap_matrix,1)
LDC = size(S_half_inv,1) LDC = size(S_half_inv,1)
allocate( & allocate( &
@ -154,7 +145,7 @@ BEGIN_PROVIDER [ double precision, S_half_inv, (AO_num,AO_num) ]
D(AO_num)) D(AO_num))
call svd( & call svd( &
AO_overlap,LDA, & AO_overlap_matrix,LDA, &
U,LDC, & U,LDC, &
D, & D, &
Vt,LDA, & Vt,LDA, &
@ -203,7 +194,7 @@ BEGIN_PROVIDER [ double precision, S_half, (ao_num,ao_num) ]
allocate(U(ao_num,ao_num),Vt(ao_num,ao_num),D(ao_num)) allocate(U(ao_num,ao_num),Vt(ao_num,ao_num),D(ao_num))
call svd(ao_overlap,size(ao_overlap,1),U,size(U,1),D,Vt,size(Vt,1),ao_num,ao_num) call svd(ao_overlap_matrix,size(ao_overlap_matrix,1),U,size(U,1),D,Vt,size(Vt,1),ao_num,ao_num)
do i=1,ao_num do i=1,ao_num
D(i) = dsqrt(D(i)) D(i) = dsqrt(D(i))