mirror of
https://gitlab.com/scemama/qmcchem.git
synced 2024-07-27 13:17:25 +02:00
230 lines
4.9 KiB
Org Mode
230 lines
4.9 KiB
Org Mode
|
#+TITLE: Asynchronous DMC
|
||
|
#+AUTHOR: Anthony Scemama
|
||
|
#+EMAIL: scemama@irsamc.ups-tlse.fr
|
||
|
|
||
|
|
||
|
#+PROPERTY: header-args :tangle no :noweb yes
|
||
|
|
||
|
|
||
|
* Main program
|
||
|
|
||
|
** Declarations
|
||
|
|
||
|
#+NAME: declarations
|
||
|
#+begin_src f90
|
||
|
include '../types.F'
|
||
|
|
||
|
integer :: iter
|
||
|
integer :: k_full ! Index of walkers in elec_coord_full
|
||
|
integer :: k_new ! Index of walkers in elec_coord_new
|
||
|
integer :: iw ! Number of copies in branching
|
||
|
integer :: l
|
||
|
|
||
|
double precision :: w, E_out, w_sum
|
||
|
|
||
|
double precision, allocatable :: elec_coord_new(:,:,:)
|
||
|
|
||
|
double precision, external :: qmc_ranf
|
||
|
|
||
|
allocate(elec_coord_new(elec_num+1,3,walk_num))
|
||
|
#+end_src
|
||
|
|
||
|
** Main flow
|
||
|
|
||
|
- Fetch ~walk_num~ electron coordinates in ~elec_coord_full~
|
||
|
- For each set of coordinates,
|
||
|
- Make a PDMC trajectory, and output the weight ~w~
|
||
|
- Perform branching depending on the value of the weight
|
||
|
- Store the new sets of coordinates in ~elec_coord_new~
|
||
|
- When ~elec_coord_new~ is full, send it to the server
|
||
|
|
||
|
#+begin_src f90 :tangle "admc.irp.f"
|
||
|
program admc
|
||
|
call run
|
||
|
call ezfio_finish
|
||
|
end program admc
|
||
|
|
||
|
subroutine run
|
||
|
implicit none
|
||
|
|
||
|
<<declarations>>
|
||
|
|
||
|
! Initialization
|
||
|
if (vmc_algo /= t_Brownian) then
|
||
|
call abrt(irp_here,'DMC should run with Brownian algorithm')
|
||
|
endif
|
||
|
|
||
|
do iter=1,2
|
||
|
call read_coords()
|
||
|
|
||
|
k_new = 1
|
||
|
do while (k_new < walk_num)
|
||
|
call pdmc_trajectory(k_full, w, E_out, w_sum)
|
||
|
write(*,*) 'E', E_out, w_sum
|
||
|
k_full = k_full+1
|
||
|
if (k_full > walk_num) k_full = 1
|
||
|
|
||
|
<<branching>>
|
||
|
|
||
|
end do
|
||
|
|
||
|
elec_coord_full(1:elec_num+1,1:3,1:walk_num) = &
|
||
|
elec_coord_new(1:elec_num+1,1:3,1:walk_num)
|
||
|
|
||
|
call write_coords()
|
||
|
end do
|
||
|
|
||
|
end subroutine run
|
||
|
|
||
|
<<read_coords>>
|
||
|
<<write_coords>>
|
||
|
<<pdmc_trajectory>>
|
||
|
#+end_src
|
||
|
|
||
|
** Branching
|
||
|
|
||
|
#+NAME: branching
|
||
|
#+begin_src f90
|
||
|
! Find number of copies
|
||
|
iw = int(w)
|
||
|
w = w - int(w)
|
||
|
if (qmc_ranf() < w) then
|
||
|
iw = iw+1
|
||
|
end if
|
||
|
|
||
|
! Duplicate walker
|
||
|
do l=1,iw
|
||
|
elec_coord_new(1:elec_num+1,1:3,k_new) = &
|
||
|
elec_coord_full(1:elec_num+1,1:3,k_full)
|
||
|
k_new = k_new+1
|
||
|
if (k_new > walk_num) exit
|
||
|
end do
|
||
|
|
||
|
#+end_src
|
||
|
|
||
|
* Read/write coordinates
|
||
|
|
||
|
** Read
|
||
|
|
||
|
Fetch a new set of coordinates for ~walk_num~ walkers from the pool of coordinates.
|
||
|
|
||
|
#+NAME: read_coords
|
||
|
#+begin_src f90
|
||
|
subroutine read_coords()
|
||
|
implicit none
|
||
|
|
||
|
integer :: i, k
|
||
|
|
||
|
do k=1,walk_num
|
||
|
do i=1,elec_num
|
||
|
read(*,*) elec_coord_full(i,1:3,k)
|
||
|
end do
|
||
|
end do
|
||
|
|
||
|
SOFT_TOUCH elec_coord_full
|
||
|
|
||
|
end subroutine read_coords
|
||
|
#+end_src
|
||
|
|
||
|
** Write
|
||
|
|
||
|
Send the current set of coordinates for ~walk_num~ walkers to the pool of coordinates.
|
||
|
|
||
|
#+NAME: write_coords
|
||
|
#+begin_src f90
|
||
|
subroutine write_coords()
|
||
|
implicit none
|
||
|
|
||
|
integer :: i, k
|
||
|
|
||
|
do k=1,walk_num
|
||
|
do i=1,elec_num
|
||
|
write(*,*) 'C', elec_coord_full(i,1:3,k)
|
||
|
end do
|
||
|
end do
|
||
|
|
||
|
end subroutine write_coords
|
||
|
#+end_src
|
||
|
|
||
|
* PDMC trajectory
|
||
|
|
||
|
Computes a PDMC trajectory until the weight ~w~ is $1/2 < w < 3/2$.
|
||
|
The energy of the trajectory is computed as
|
||
|
\[
|
||
|
E = \frac{\sum_i w_i E(R_i)}{\sum_i w_i}
|
||
|
\]
|
||
|
|
||
|
The function returns:
|
||
|
- ~pdmc_weight~: the last of all $w_i$
|
||
|
- ~E_out~: The average energy $E$ of the trajectory
|
||
|
- ~w_sum~: The sum of the weights
|
||
|
|
||
|
#+NAME: declarations_pdmc
|
||
|
#+begin_src f90
|
||
|
integer :: i,j,l
|
||
|
double precision :: delta
|
||
|
|
||
|
! If true, continue to make more steps
|
||
|
logical :: loop
|
||
|
|
||
|
! Brownian step variables
|
||
|
double precision :: p,q
|
||
|
real :: delta_x
|
||
|
logical :: accepted
|
||
|
|
||
|
! Local energies from the past
|
||
|
double precision :: E_loc_save(4)
|
||
|
|
||
|
#+end_src
|
||
|
|
||
|
#+NAME: pdmc_trajectory
|
||
|
#+begin_src f90
|
||
|
subroutine pdmc_trajectory(k_full, pdmc_weight, E_out, w_sum)
|
||
|
implicit none
|
||
|
|
||
|
integer, intent(in) :: k_full
|
||
|
double precision, intent(out) :: pdmc_weight, E_out, w_sum
|
||
|
|
||
|
<<declarations_pdmc>>
|
||
|
|
||
|
elec_coord(1:elec_num+1,1:3) = elec_coord_full(1:elec_num+1,1:3,k_full)
|
||
|
TOUCH elec_coord
|
||
|
|
||
|
E_out = 0.d0
|
||
|
w_sum = 0.d0
|
||
|
|
||
|
E_loc_save(1:4) = E_loc
|
||
|
|
||
|
pdmc_weight = 1.d0
|
||
|
loop = .True.
|
||
|
|
||
|
do while (loop)
|
||
|
|
||
|
call brownian_step(p,q,accepted,delta_x)
|
||
|
|
||
|
delta = (9.d0*E_loc+19.d0*E_loc_save(1)-5.d0*E_loc_save(2)+E_loc_save(3))/24.d0
|
||
|
delta = (delta - E_ref)*p
|
||
|
|
||
|
if (delta >= 0.d0) then
|
||
|
pdmc_weight = dexp(-dtime_step*delta)
|
||
|
else
|
||
|
pdmc_weight = 2.d0-dexp(dtime_step*delta)
|
||
|
endif
|
||
|
elec_coord(elec_num+1,1) += p*time_step
|
||
|
elec_coord(elec_num+1,2) = E_loc
|
||
|
elec_coord(elec_num+1,3) = pdmc_weight
|
||
|
|
||
|
w_sum = w_sum + pdmc_weight
|
||
|
E_out = E_out + pdmc_weight * E_loc
|
||
|
|
||
|
loop = pdmc_weight > 0.5d0 .or. pdmc_weight < 1.5d0
|
||
|
|
||
|
end do
|
||
|
|
||
|
elec_coord_full(1:elec_num+1,1:3,k_full) = elec_coord(1:elec_num+1,1:3)
|
||
|
SOFT_TOUCH elec_coord_full
|
||
|
|
||
|
end subroutine pdmc_trajectory
|
||
|
#+end_src
|
||
|
|