Working on Slater dressing

This commit is contained in:
Anthony Scemama 2017-06-26 18:11:49 +02:00
parent ddcc1a5e51
commit 5f357ee68b
8 changed files with 101 additions and 83 deletions

View File

@ -14,7 +14,13 @@ interface: ezfio, provider
type: double precision
doc: Orthogonal AO basis
size: (ao_basis.ao_num,ao_basis.ao_num)
interface: ezfio, provider
interface: ezfio
[ao_orthoSlaOverlap]
type: double precision
doc: Orthogonal AO basis
size: (ao_basis.ao_num,nuclei.nucl_num)
interface: ezfio

View File

@ -14,7 +14,6 @@ BEGIN_PROVIDER [ double precision, cusp_A, (nucl_num, nucl_num) ]
! Projector
do mu=1,mo_tot_num
cusp_A(A,B) += AO_orthoSlaOverlap_matrix(mu,B) * ao_ortho_value_at_nucl(mu,A)
! cusp_A(A,B) += MOSlaOverlap_matrix(mu,B) * mo_value_at_nucl(mu,A)
enddo
enddo
enddo
@ -60,6 +59,5 @@ BEGIN_PROVIDER [ double precision, cusp_C, (nucl_num, mo_tot_num) ]
stop 'dgetrs failed'
endif
END_PROVIDER

View File

@ -86,8 +86,10 @@ subroutine run
enddo
mo_coef(1:ao_num,1:mo_tot_num) = cusp_corrected_mos(1:ao_num,1:mo_tot_num)
SOFT_TOUCH mo_coef slater_coef
! call ezfio_set_Hartree_Fock_SlaterDressed_slater_coef_ezfio(slater_coef)
! call save_mos
call ezfio_set_Hartree_Fock_SlaterDressed_slater_coef_ezfio(slater_coef)
call ezfio_set_Hartree_Fock_SlaterDressed_projector(ao_ortho_canonical_coef(1:ao_num,1:ao_num))
call ezfio_set_Hartree_Fock_SlaterDressed_ao_orthoSlaOverlap(AO_orthoSlaOverlap_matrix)
call save_mos
print *, 'ci'
print *, mo_coef(1:ao_num,1)
print *, 'cAi'

View File

@ -91,74 +91,82 @@ BEGIN_PROVIDER [ double precision, cusp_corrected_mos, (ao_num_align,mo_tot_num)
cusp_corrected_mos(1:ao_num,1:mo_tot_num) = mo_coef(1:ao_num,1:mo_tot_num)
slater_coef(1:nucl_num,1:mo_tot_num) = cusp_C(1:nucl_num,1:mo_tot_num)
return
else
do idx_dressing=1,mo_tot_num
if (idx_dressing>1) then
TOUCH idx_dressing
endif
do j=1,mo_tot_num
do i=1,mo_tot_num
F(i,j) = Fock_matrix_mo(i,j)
enddo
enddo
do j=1,mo_tot_num
do i=1,ao_num
M(i,j) = mo_coef(i,j)
enddo
enddo
integer :: it
do it=1,128
! print *, 'X', ao_ortho_canonical_coef(1:ao_num,1:ao_num)
! print *, 'C', mo_coef(1:ao_num,1:mo_tot_num)
! print *, 'Cp', mo_coef_in_ao_ortho_basis(1:ao_num,1:mo_tot_num)
! print *, 'cAi', cusp_C(1:nucl_num,1:mo_tot_num)
! print *, 'FmuA', AO_orthoSlaH_matrix(1:ao_num,1:nucl_num)
! print *, 'Fock:', Fock_matrix_ao(1:ao_num,1:ao_num)
! print *, 'Diag Dressing:', ao_ortho_mono_elec_integral_dressing(1:ao_num,1:ao_num)
! print *, 'Dressing:', ao_mono_elec_integral_dressing(1:ao_num,1:ao_num)
! print *, 'Dressed Fock:', Fock_matrix_ao(1:ao_num,1:ao_num) + ao_mono_elec_integral_dressing(1:ao_num,1:ao_num)
! print *, 'AO_orthoSlaOverlap_matrix', AO_orthoSlaOverlap_matrix(1:ao_num,1:nucl_num)
! print *, 'AO_orthoSlaH_matrix', AO_orthoSlaH_matrix(1:ao_num,1:nucl_num)
! print *, 'ao_ortho_mono_elec_integral', ao_ortho_mono_elec_integral(1:ao_num,1:ao_num)
! print *, 'Fock MO:', Fock_matrix_mo(1:mo_tot_num,1:mo_tot_num)
do j=1,mo_tot_num
do i=1,mo_tot_num
Fock_matrix_mo(i,j) += mo_mono_elec_integral_dressing(i,j)
enddo
enddo
do i=1,mo_tot_num
Fock_matrix_diag_mo(i) = Fock_matrix_mo(i,i)
enddo
! print *, 'Dressed Fock MO:', Fock_matrix_mo(1:mo_tot_num,1:mo_tot_num)
double precision :: conv
conv = 0.d0
do j=1,mo_tot_num
do i=1,mo_tot_num
if (i==j) cycle
conv = max(conv,Fock_matrix_mo(i,j))
enddo
enddo
TOUCH Fock_matrix_mo Fock_matrix_diag_mo
mo_coef(1:ao_num,1:mo_tot_num) = eigenvectors_fock_matrix_mo(1:ao_num,1:mo_tot_num)
TOUCH mo_coef
!print *, 'C', mo_coef(1:ao_num,1:mo_tot_num)
!print *, '-----'
print *, idx_dressing, it, real(mo_coef(1,idx_dressing)), real(conv)
if (conv < 1.d-5) exit
!stop
enddo
cusp_corrected_mos(1:ao_num,idx_dressing) = mo_coef(1:ao_num,idx_dressing)
slater_coef(1:nucl_num,idx_dressing) = cusp_C(1:nucl_num,idx_dressing)
enddo
idx_dressing = 1
mo_coef(1:ao_num,1:mo_tot_num) = M(1:ao_num,1:mo_tot_num)
soft_TOUCH mo_coef idx_dressing slater_coef
endif
do idx_dressing=1,mo_tot_num
if (idx_dressing>1) then
TOUCH idx_dressing
endif
do j=1,mo_tot_num
do i=1,mo_tot_num
F(i,j) = Fock_matrix_mo(i,j)
enddo
enddo
do j=1,mo_tot_num
do i=1,ao_num
M(i,j) = mo_coef(i,j)
enddo
enddo
integer :: it
do it=1,128
! print *, 'C', mo_coef(1:ao_num,1:mo_tot_num)
! print *, 'Cp', mo_coef_in_ao_ortho_basis(1:ao_num,1:mo_tot_num)
! print *, 'cAi', cusp_C(1:nucl_num,1:mo_tot_num)
! print *, 'FmuA', AO_orthoSlaH_matrix(1:ao_num,1:nucl_num)
! print *, 'Fock:', Fock_matrix_ao(1:ao_num,1:ao_num)
! print *, 'Diag Dressing:', ao_ortho_mono_elec_integral_dressing(1:ao_num,1:ao_num)
! print *, 'Dressing:', ao_mono_elec_integral_dressing(1:ao_num,1:ao_num)
! print *, 'Dressed Fock:', Fock_matrix_ao(1:ao_num,1:ao_num) + ao_mono_elec_integral_dressing(1:ao_num,1:ao_num)
! print *, 'AO_orthoSlaOverlap_matrix', AO_orthoSlaOverlap_matrix(1:ao_num,1:nucl_num)
! print *, 'AO_orthoSlaH_matrix', AO_orthoSlaH_matrix(1:ao_num,1:nucl_num)
! print *, 'ao_ortho_mono_elec_integral', ao_ortho_mono_elec_integral(1:ao_num,1:ao_num)
do j=1,mo_tot_num
do i=1,mo_tot_num
Fock_matrix_mo(i,j) += mo_mono_elec_integral_dressing(i,j)
enddo
enddo
do i=1,mo_tot_num
Fock_matrix_diag_mo(i) = Fock_matrix_mo(i,i)
enddo
double precision :: conv
conv = 0.d0
do j=1,mo_tot_num
do i=1,mo_tot_num
if (i==j) cycle
conv = max(conv,Fock_matrix_mo(i,j))
enddo
enddo
TOUCH Fock_matrix_mo Fock_matrix_diag_mo
mo_coef(1:ao_num,1:mo_tot_num) = eigenvectors_fock_matrix_mo(1:ao_num,1:mo_tot_num)
TOUCH mo_coef
! print *, 'C', mo_coef(1:ao_num,1:mo_tot_num)
! print *, '-----'
print *, idx_dressing, it, real(mo_coef(1,idx_dressing)), real(conv)
if (conv < 1.d-5) exit
enddo
cusp_corrected_mos(1:ao_num,idx_dressing) = mo_coef(1:ao_num,idx_dressing)
slater_coef(1:nucl_num,idx_dressing) = cusp_C(1:nucl_num,idx_dressing)
enddo
idx_dressing = 1
mo_coef(1:ao_num,1:mo_tot_num) = M(1:ao_num,1:mo_tot_num)
soft_TOUCH mo_coef idx_dressing slater_coef
END_PROVIDER

View File

@ -344,7 +344,7 @@ subroutine GauSlaNuclear(expGau,cGau,aGau,expSla,cSla,ZNuc,cNuc,result)
ss = k*ss
! Print result
write(*,*) ss
! write(*,*) ss
result = 0.d0
end
@ -520,6 +520,7 @@ BEGIN_PROVIDER [ double precision, AO_orthoSla$X_matrix, (ao_num, nucl_num) ]
END_PROVIDER
SUBST [ X ]
Overlap ;;
@ -616,9 +617,9 @@ BEGIN_PROVIDER [ double precision, AO_orthoSlaH_matrix, (ao_num, nucl_num) ]
BEGIN_DOC
! <AO ortho | Slater>
END_DOC
call dgemm('T','N',ao_num,nucl_num,ao_num,1.d0, &
ao_ortho_canonical_coef, size(ao_ortho_canonical_coef,1), &
GauSlaH_matrix, size(GauSlaH_matrix,1), &
0.d0, AO_orthoSlaH_matrix, size(AO_orthoSlaH_matrix,1))
call dgemm('T','N',ao_num,nucl_num,ao_num,1.d0, &
ao_ortho_canonical_coef, size(ao_ortho_canonical_coef,1), &
GauSlaH_matrix, size(GauSlaH_matrix,1), &
0.d0, AO_orthoSlaH_matrix, size(AO_orthoSlaH_matrix,1))
END_PROVIDER

View File

@ -8,7 +8,10 @@ BEGIN_PROVIDER [ double precision, slater_expo, (nucl_num) ]
if (exists) then
slater_expo(1:nucl_num) = slater_expo_ezfio(1:nucl_num)
else
slater_expo(1:nucl_num) = nucl_charge(1:nucl_num)
integer :: i
do i=1,nucl_num
slater_expo(i) = nucl_charge(i)
enddo
call ezfio_set_Hartree_Fock_SlaterDressed_slater_expo_ezfio(slater_expo)
endif
END_PROVIDER

View File

@ -106,9 +106,9 @@ END_PROVIDER
ao_ortho_canonical_coef(i,i) = 1.d0
enddo
!call ortho_lowdin(ao_overlap,size(ao_overlap,1),ao_num,ao_ortho_canonical_coef,size(ao_ortho_canonical_coef,1),ao_num)
!ao_ortho_canonical_num=ao_num
!return
call ortho_lowdin(ao_overlap,size(ao_overlap,1),ao_num,ao_ortho_canonical_coef,size(ao_ortho_canonical_coef,1),ao_num)
ao_ortho_canonical_num=ao_num
return
if (ao_cartesian) then

View File

@ -75,7 +75,7 @@ BEGIN_PROVIDER [ double precision, mo_coef_in_ao_ortho_basis, (ao_num_align, mo_
!
! C^(-1).C_mo
END_DOC
call dgemm('T','N',ao_num,mo_tot_num,ao_num,1.d0, &
call dgemm('N','N',ao_num,mo_tot_num,ao_num,1.d0, &
ao_ortho_canonical_coef_inv, size(ao_ortho_canonical_coef_inv,1),&
mo_coef, size(mo_coef,1), 0.d0, &
mo_coef_in_ao_ortho_basis, size(mo_coef_in_ao_ortho_basis,1))
@ -290,13 +290,13 @@ subroutine ao_ortho_cano_to_ao(A_ao,LDA_ao,A,LDA)
allocate ( T(ao_num_align,ao_num) )
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: T
call dgemm('N','N', ao_num, ao_num, ao_num, &
call dgemm('T','N', ao_num, ao_num, ao_num, &
1.d0, &
ao_ortho_canonical_coef_inv, size(ao_ortho_canonical_coef_inv,1), &
A_ao,LDA_ao, &
0.d0, T, ao_num_align)
call dgemm('N','T', ao_num, ao_num, ao_num, 1.d0, &
call dgemm('N','N', ao_num, ao_num, ao_num, 1.d0, &
T, size(T,1), &
ao_ortho_canonical_coef_inv,size(ao_ortho_canonical_coef_inv,1),&
0.d0, A, LDA)