LCPQ/qmcchem
10
1
Fork 1
mirror of https://gitlab.com/scemama/qmcchem.git synced 2024-12-22 04:13:31 +01:00
Code Releases Activity

Avoid pseudo explosion

Browse Source
This commit is contained in:
Anthony Scemama 2020-05-26 01:46:14 +02:00
parent 9b191eaa7b
commit b837d4f346
1 changed files with 30 additions and 30 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

60
src/PROPERTIES/properties_energy.irp.f
View File

@ -13,7 +13,7 @@ BEGIN_PROVIDER [ double precision, single_det_E_kin ]
do i=1,elec_num do i=1,elec_num
single_det_E_kin -= 0.5d0*single_det_lapl(i)/single_det_value single_det_E_kin -= 0.5d0*single_det_lapl(i)/single_det_value
enddo enddo
END_PROVIDER END_PROVIDER
@ -22,7 +22,7 @@ BEGIN_PROVIDER [ double precision, single_det_E_loc ]
BEGIN_DOC BEGIN_DOC
! Local energy : single_det_E_kin + E_pot + E_nucl ! Local energy : single_det_E_kin + E_pot + E_nucl
END_DOC END_DOC
single_det_E_loc = single_det_E_kin + E_pot + E_nucl single_det_E_loc = single_det_E_kin + E_pot + E_nucl
END_PROVIDER END_PROVIDER
@ -32,7 +32,7 @@ BEGIN_PROVIDER [ double precision, E_pot_grad, (elec_num,3) ]
BEGIN_DOC BEGIN_DOC
! Gradient of the Electronic Potential energy ! Gradient of the Electronic Potential energy
END_DOC END_DOC
integer :: i,j integer :: i,j
double precision :: dinv double precision :: dinv
do i=1,elec_num do i=1,elec_num
@ -64,7 +64,7 @@ BEGIN_PROVIDER [ double precision, E_pot_grad, (elec_num,3) ]
E_pot_grad(i,3) += nucl_elec_dist_vec(3,j,i)*dinv E_pot_grad(i,3) += nucl_elec_dist_vec(3,j,i)*dinv
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
@ -73,7 +73,7 @@ BEGIN_PROVIDER [ double precision, E_pot_elec, (elec_num) ]
BEGIN_DOC BEGIN_DOC
! Electronic Potential energy ! Electronic Potential energy
END_DOC END_DOC
integer :: i, j integer :: i, j
if (do_pseudo) then if (do_pseudo) then
do i=1,elec_num do i=1,elec_num
@ -89,7 +89,7 @@ BEGIN_PROVIDER [ double precision, E_pot_elec, (elec_num) ]
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
!DIR$ LOOP COUNT(50) !DIR$ LOOP COUNT(50)
do j=1,elec_num do j=1,elec_num
E_pot_elec(i) = E_pot_elec(i) + 0.5d0*elec_dist_inv(j,i) E_pot_elec(i) = E_pot_elec(i) + 0.5d0*elec_dist_inv(j,i)
enddo enddo
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
!DIR$ LOOP COUNT(50) !DIR$ LOOP COUNT(50)
@ -97,7 +97,7 @@ BEGIN_PROVIDER [ double precision, E_pot_elec, (elec_num) ]
E_pot_elec(i) = E_pot_elec(i) - nucl_charge(j)*nucl_elec_dist_inv(j,i) E_pot_elec(i) = E_pot_elec(i) - nucl_charge(j)*nucl_elec_dist_inv(j,i)
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
@ -106,7 +106,7 @@ BEGIN_PROVIDER [ double precision, E_pot_elec_one, (elec_num) ]
BEGIN_DOC BEGIN_DOC
! Electronic Potential energy ! Electronic Potential energy
END_DOC END_DOC
integer :: i, j integer :: i, j
do i=1,elec_num do i=1,elec_num
E_pot_elec_one(i) = 0.d0 E_pot_elec_one(i) = 0.d0
@ -116,7 +116,7 @@ BEGIN_PROVIDER [ double precision, E_pot_elec_one, (elec_num) ]
E_pot_elec_one(i) -= nucl_charge(j)*nucl_elec_dist_inv(j,i) E_pot_elec_one(i) -= nucl_charge(j)*nucl_elec_dist_inv(j,i)
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
@ -125,7 +125,7 @@ BEGIN_PROVIDER [ double precision, E_pot_elec_two, (elec_num) ]
BEGIN_DOC BEGIN_DOC
! Electronic Potential energy ! Electronic Potential energy
END_DOC END_DOC
integer :: i, j integer :: i, j
do i=1,elec_num do i=1,elec_num
E_pot_elec_two(i) = 0.d0 E_pot_elec_two(i) = 0.d0
@ -138,22 +138,22 @@ BEGIN_PROVIDER [ double precision, E_pot_elec_two, (elec_num) ]
E_pot_elec_two(i) += 0.5d0*elec_dist_inv(j,i) E_pot_elec_two(i) += 0.5d0*elec_dist_inv(j,i)
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ double precision, E_kin_elec, (elec_num) ] BEGIN_PROVIDER [ double precision, E_kin_elec, (elec_num) ]
implicit none implicit none
BEGIN_DOC BEGIN_DOC
! Electronic Kinetic energy : -1/2 (Lapl.Psi)/Psi ! Electronic Kinetic energy : -1/2 (Lapl.Psi)/Psi
END_DOC END_DOC
integer :: i integer :: i
do i=1,elec_num do i=1,elec_num
E_kin_elec(i) = -0.5d0*psi_lapl_psi_inv(i) E_kin_elec(i) = -0.5d0*psi_lapl_psi_inv(i)
enddo enddo
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ double precision, dmc_zv_weight ] BEGIN_PROVIDER [ double precision, dmc_zv_weight ]
@ -182,7 +182,7 @@ BEGIN_PROVIDER [ double precision, E_nucl ]
BEGIN_DOC BEGIN_DOC
! Nuclear potential energy ! Nuclear potential energy
END_DOC END_DOC
E_nucl = 0.d0 E_nucl = 0.d0
integer :: i, j integer :: i, j
do i=1,nucl_num do i=1,nucl_num
@ -190,7 +190,7 @@ BEGIN_PROVIDER [ double precision, E_nucl ]
E_nucl += nucl_charge(i)*nucl_charge(j)/nucl_dist(j,i) E_nucl += nucl_charge(i)*nucl_charge(j)/nucl_dist(j,i)
enddo enddo
enddo enddo
E_nucl_min = min(E_nucl,E_nucl_min) E_nucl_min = min(E_nucl,E_nucl_min)
E_nucl_max = max(E_nucl,E_nucl_max) E_nucl_max = max(E_nucl,E_nucl_max)
SOFT_TOUCH E_nucl_min E_nucl_max SOFT_TOUCH E_nucl_min E_nucl_max
@ -202,13 +202,13 @@ BEGIN_PROVIDER [ double precision, E_pot ]
BEGIN_DOC BEGIN_DOC
! Electronic Potential energy ! Electronic Potential energy
END_DOC END_DOC
E_pot = 0.d0 E_pot = 0.d0
integer :: i, j integer :: i, j
do i=1,elec_num do i=1,elec_num
E_pot += E_pot_elec(i) E_pot += E_pot_elec(i)
enddo enddo
E_pot_min = min(E_pot,E_pot_min) E_pot_min = min(E_pot,E_pot_min)
E_pot_max = max(E_pot,E_pot_max) E_pot_max = max(E_pot,E_pot_max)
SOFT_TOUCH E_pot_min E_pot_max SOFT_TOUCH E_pot_min E_pot_max
@ -220,16 +220,16 @@ BEGIN_PROVIDER [ double precision, E_kin ]
BEGIN_DOC BEGIN_DOC
! Electronic Kinetic energy : -1/2 (Lapl.Psi)/Psi ! Electronic Kinetic energy : -1/2 (Lapl.Psi)/Psi
END_DOC END_DOC
E_kin = 0.d0 E_kin = 0.d0
integer :: i integer :: i
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
!DIR$ LOOP COUNT(200) !DIR$ LOOP COUNT(200)
do i=1,elec_num do i=1,elec_num
E_kin -= 0.5d0*psi_lapl_psi_inv(i) E_kin -= 0.5d0*psi_lapl_psi_inv(i)
enddo enddo
E_kin_min = min(E_kin,E_kin_min) E_kin_min = min(E_kin,E_kin_min)
E_kin_max = max(E_kin,E_kin_max) E_kin_max = max(E_kin,E_kin_max)
SOFT_TOUCH E_kin_min E_kin_max SOFT_TOUCH E_kin_min E_kin_max
@ -242,7 +242,7 @@ BEGIN_PROVIDER [ double precision, E_loc ]
BEGIN_DOC BEGIN_DOC
! Local energy : E_kin + E_pot + E_nucl ! Local energy : E_kin + E_pot + E_nucl
END_DOC END_DOC
integer :: i integer :: i
E_loc = E_nucl E_loc = E_nucl
!DIR$ VECTOR ALIGNED !DIR$ VECTOR ALIGNED
@ -250,14 +250,14 @@ BEGIN_PROVIDER [ double precision, E_loc ]
do i=1,elec_num do i=1,elec_num
E_loc += E_kin_elec(i) + E_pot_elec(i) E_loc += E_kin_elec(i) + E_pot_elec(i)
enddo enddo
! ! Avoid divergence of E_loc ! Avoid divergence of E_loc and population explosion
! if (qmc_method == t_DMC) then if (do_pseudo) then
! double precision :: delta_e double precision :: delta_e
! delta_e = E_loc-E_ref ! delta_e = E_loc-E_ref
! E_loc = E_ref + erf(1.d0/(time_step*delta_e*time_step*delta_e)) * delta_e ! E_loc = E_ref + erf(delta_e*time_step_sq)/time_step_sq
! endif E_loc = max(2.d0*E_ref, E_loc)
endif
E_loc_min = min(E_loc,E_loc_min) E_loc_min = min(E_loc,E_loc_min)
E_loc_max = max(E_loc,E_loc_max) E_loc_max = max(E_loc,E_loc_max)
SOFT_TOUCH E_loc_min E_loc_max SOFT_TOUCH E_loc_min E_loc_max
@ -270,7 +270,7 @@ BEGIN_PROVIDER [ double precision, E_loc_zv ]
BEGIN_DOC BEGIN_DOC
! Zero-variance parameter on E_loc ! Zero-variance parameter on E_loc
END_DOC END_DOC
E_loc_zv = E_loc E_loc_zv = E_loc
E_loc_zv += (E_trial-E_loc) * dmc_zv_weight E_loc_zv += (E_trial-E_loc) * dmc_zv_weight
! E_loc_zv += - time_step*(E_trial**2 + 1.44341217940434 - E_loc**2)*dmc_zv_weight ! E_loc_zv += - time_step*(E_trial**2 + 1.44341217940434 - E_loc**2)*dmc_zv_weight
! E_loc_zv(3) = dmc_zv_weight_half ! E_loc_zv(3) = dmc_zv_weight_half