TABLE OF CONTENTS


ABINIT/m_hmc [ Modules ]

[ Top ] [ Modules ]

NAME

  m_hmc

FUNCTION

  Auxiliary hmc functions

COPYRIGHT

  Copyright (C) 2018-2022 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

SOURCE

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

ABINIT/m_hmc/compute_kinetic_energy [ Functions ]

[ Top ] [ Functions ]

NAME

  comute_kintic_energy

FUNCTION

  Computes kintetic energy

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

64 subroutine compute_kinetic_energy(ab_mover,vel,ekin)
65     
66 !Arguments ------------------------------------
67  type(abimover),intent(in)   :: ab_mover              
68  real(dp),      intent(in)   :: vel(3,ab_mover%natom) ! velocities
69  real(dp),      intent(out)  :: ekin                  ! output kinetic energy     
70 
71 !Local variables-------------------------------
72  integer :: ii,jj                                
73 !character(len=500) :: msg                   
74  
75 ! *************************************************************************
76 
77  ekin=0.0
78  do ii=1,ab_mover%natom
79    do jj=1,3
80      ekin=ekin+half*ab_mover%amass(ii)*vel(jj,ii)**2
81    end do
82  end do
83 
84 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

SOURCE

109 subroutine generate_random_velocities(ab_mover,kbtemp,seed,vel,ekin)
110     
111 !Arguments ------------------------------------
112  type(abimover),intent(in)   :: ab_mover              
113  integer,       intent(inout):: seed
114  real(dp),      intent(in)   :: kbtemp
115  real(dp),      intent(inout):: vel(3,ab_mover%natom) ! velocities
116  real(dp),      intent(out)  :: ekin                  ! output kinetic energy     
117 !Local variables-------------------------------
118  integer :: ii,jj
119  real(dp):: mtot,mvtot(3),mv2tot,factor                                
120 !character(len=500) :: msg                   
121  
122 ! *************************************************************************
123 
124 
125  mtot=sum(ab_mover%amass(:))         ! total mass to eventually get rid of total center of mass (CoM) momentum
126  !generate velocities from normal distribution with zero mean and correct standard deviation
127  do ii=1,ab_mover%natom
128    do jj=1,3
129      vel(jj,ii)=sqrt(kbtemp/ab_mover%amass(ii))*cos(two_pi*uniformrandom(seed))
130      vel(jj,ii)=vel(jj,ii)*sqrt(-2.0*log(uniformrandom(seed)))
131    end do
132  end do
133  !since number of atoms is most probably not big enough to obtain overall zero CoM momentum, shift the velocities
134  !and then renormalize
135  mvtot=zero ! total momentum
136  do ii=1,ab_mover%natom
137    do jj=1,3
138      mvtot(jj)=mvtot(jj)+ab_mover%amass(ii)*vel(jj,ii)
139    end do
140  end do
141  do ii=1,ab_mover%natom
142    do jj=1,3
143      vel(jj,ii)=vel(jj,ii)-(mvtot(jj)/mtot)
144    end do
145  end do
146  !now the total cell momentum is zero
147  mv2tot=0.0
148  do ii=1,ab_mover%natom
149    do jj=1,3
150      mv2tot=mv2tot+ab_mover%amass(ii)*vel(jj,ii)**2
151    end do
152  end do
153  factor = mv2tot/(dble(3*ab_mover%natom))
154  factor = sqrt(kbtemp/factor)
155  vel(:,:)=vel(:,:)*factor
156 
157  call compute_kinetic_energy(ab_mover,vel,ekin)
158 
159 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

SOURCE

182 subroutine metropolis_check(seed,de,kbtemp,iacc)
183     
184 !Arguments ------------------------------------
185  integer,       intent(inout):: seed
186  real(dp),      intent(in)   :: de         
187  real(dp),      intent(in)   :: kbtemp         
188  integer,       intent(inout):: iacc    
189 
190 !Local variables-------------------------------
191  real(dp)   :: rnd
192 !character(len=500) :: msg                   
193  
194 ! *************************************************************************
195 
196  iacc=0
197  rnd=uniformrandom(seed)
198  if(de<0)then
199    iacc=1
200  else
201    if(exp(-de/kbtemp)>rnd)then
202       iacc=1
203    end if
204  end if
205 
206 end subroutine metropolis_check