TABLE OF CONTENTS


ABINIT/m_hmc [ Modules ]

[ Top ] [ Modules ]

NAME

  m_hmc

FUNCTION

  Auxiliary hmc functions

COPYRIGHT

  Copyright (C) 2018 ABINIT group (SPr)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

NOTES

PARENTS

CHILDREN

SOURCE

23 #if defined HAVE_CONFIG_H
24 #include "config.h"
25 #endif
26 
27 #include "abi_common.h"
28 
29 module m_hmc
30     
31  use defs_basis
32  use m_abicore
33  use m_errors
34  use m_abimover
35  use m_io_tools
36 
37 !use m_geometry,       only : xred2xcart
38  use m_numeric_tools,  only : uniformrandom
39 
40  implicit none
41 
42  private
43 
44 ! *************************************************************************
45  public :: compute_kinetic_energy 
46  public :: generate_random_velocities
47  public :: metropolis_check
48 
49 contains 

ABINIT/m_hmc/compute_kinetic_energy [ Functions ]

[ Top ] [ Functions ]

NAME

  comute_kintic_energy

FUNCTION

  Computes kintetic energy

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

CHILDREN

SOURCE

 73 subroutine compute_kinetic_energy(ab_mover,vel,ekin)
 74     
 75  use defs_basis
 76 
 77 !This section has been created automatically by the script Abilint (TD).
 78 !Do not modify the following lines by hand.
 79 #undef ABI_FUNC
 80 #define ABI_FUNC 'compute_kinetic_energy'
 81 !End of the abilint section
 82 
 83  implicit none
 84 
 85 !Arguments ------------------------------------
 86  type(abimover),intent(in)   :: ab_mover              
 87  real(dp),      intent(in)   :: vel(3,ab_mover%natom) ! velocities
 88  real(dp),      intent(out)  :: ekin                  ! output kinetic energy     
 89 
 90 !Local variables-------------------------------
 91  integer :: ii,jj                                
 92 !character(len=500) :: msg                   
 93  
 94 ! *************************************************************************
 95 
 96  ekin=0.0
 97  do ii=1,ab_mover%natom
 98    do jj=1,3
 99      ekin=ekin+half*ab_mover%amass(ii)*vel(jj,ii)**2
100    end do
101  end do
102 
103 end subroutine compute_kinetic_energy

ABINIT/m_hmc/generate_random_velocities [ Functions ]

[ Top ] [ Functions ]

NAME

  generate_random_velocities

FUNCTION

  Generate normally distributed random velocities

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

CHILDREN

SOURCE

132 subroutine generate_random_velocities(ab_mover,kbtemp,seed,vel,ekin)
133     
134  use defs_basis
135 
136 !This section has been created automatically by the script Abilint (TD).
137 !Do not modify the following lines by hand.
138 #undef ABI_FUNC
139 #define ABI_FUNC 'generate_random_velocities'
140 !End of the abilint section
141 
142  implicit none
143 
144 !Arguments ------------------------------------
145  type(abimover),intent(in)   :: ab_mover              
146  integer,       intent(inout):: seed
147  real(dp),      intent(in)   :: kbtemp
148  real(dp),      intent(inout):: vel(3,ab_mover%natom) ! velocities
149  real(dp),      intent(out)  :: ekin                  ! output kinetic energy     
150 !Local variables-------------------------------
151  integer :: ii,jj
152  real(dp):: mtot,mvtot(3),mv2tot,factor                                
153 !character(len=500) :: msg                   
154  
155 ! *************************************************************************
156 
157 
158  mtot=sum(ab_mover%amass(:))         ! total mass to eventually get rid of total center of mass (CoM) momentum
159  !generate velocities from normal distribution with zero mean and correct standard deviation
160  do ii=1,ab_mover%natom
161    do jj=1,3
162      vel(jj,ii)=sqrt(kbtemp/ab_mover%amass(ii))*cos(two_pi*uniformrandom(seed))
163      vel(jj,ii)=vel(jj,ii)*sqrt(-2.0*log(uniformrandom(seed)))
164    end do
165  end do
166  !since number of atoms is most probably not big enough to obtain overall zero CoM momentum, shift the velocities
167  !and then renormalize
168  mvtot=zero ! total momentum
169  do ii=1,ab_mover%natom
170    do jj=1,3
171      mvtot(jj)=mvtot(jj)+ab_mover%amass(ii)*vel(jj,ii)
172    end do
173  end do
174  do ii=1,ab_mover%natom
175    do jj=1,3
176      vel(jj,ii)=vel(jj,ii)-(mvtot(jj)/mtot)
177    end do
178  end do
179  !now the total cell momentum is zero
180  mv2tot=0.0
181  do ii=1,ab_mover%natom
182    do jj=1,3
183      mv2tot=mv2tot+ab_mover%amass(ii)*vel(jj,ii)**2
184    end do
185  end do
186  factor = mv2tot/(dble(3*ab_mover%natom))
187  factor = sqrt(kbtemp/factor)
188  vel(:,:)=vel(:,:)*factor
189 
190  call compute_kinetic_energy(ab_mover,vel,ekin)
191 
192 end subroutine generate_random_velocities

ABINIT/m_hmc/metropolis_check [ Functions ]

[ Top ] [ Functions ]

NAME

  metropolis_check      

FUNCTION

  Make an acceptance decision based on the energy differences

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

CHILDREN

SOURCE

219 subroutine metropolis_check(seed,de,kbtemp,iacc)
220     
221  use defs_basis
222 
223 !This section has been created automatically by the script Abilint (TD).
224 !Do not modify the following lines by hand.
225 #undef ABI_FUNC
226 #define ABI_FUNC 'metropolis_check'
227 !End of the abilint section
228 
229  implicit none
230 
231 !Arguments ------------------------------------
232  integer,       intent(inout):: seed
233  real(dp),      intent(in)   :: de         
234  real(dp),      intent(in)   :: kbtemp         
235  integer,       intent(inout):: iacc    
236 
237 !Local variables-------------------------------
238  real(dp)   :: rnd
239 !character(len=500) :: msg                   
240  
241 ! *************************************************************************
242 
243  iacc=0
244  rnd=uniformrandom(seed)
245  if(de<0)then
246    iacc=1
247  else
248    if(exp(-de/kbtemp)>rnd)then
249       iacc=1
250    end if
251  end if
252 
253 end subroutine metropolis_check