Line data Source code
1 : !!****m* ABINIT/m_lattice_langevin_mover
2 : !! NAME
3 : !! m_lattice_langevin_mover
4 : !!
5 : !! FUNCTION
6 : !! This module contains the langevin (NVT) lattice mover.
7 : !! It is a translation from the ASE (GPL licenced) Langevin mover python code to fortran.
8 : !! The original code can be found at
9 : !! https://gitlab.com/ase/ase/blob/master/ase/md/langevin.py
10 : !! The method is described in
11 : !! E. V.-Eijnden, and G. Ciccotti, Chem. Phys. Lett. 429, 310 (2006)
12 : !! https://doi.org/10.1016/j.cplett.2006.07.086
13 : !!
14 : !!
15 : !! Datatypes:
16 : !!
17 : !! * lattice_langevin_mover_t: defines the lattice movers
18 : !!
19 : !! Subroutines:
20 : !! TODO: add this when F2003 doc style is determined.
21 : !!
22 : !!
23 : !! COPYRIGHT
24 : !! Copyright (C) 2001-2026 ABINIT group (hexu)
25 : !! This file is distributed under the terms of the
26 : !! GNU General Public License, see ~abinit/COPYING
27 : !! or http://www.gnu.org/copyleft/gpl.txt .
28 : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
29 : !!
30 : !! SOURCE
31 :
32 :
33 :
34 : #if defined HAVE_CONFIG_H
35 : #include "config.h"
36 : #endif
37 :
38 : #include "abi_common.h"
39 :
40 : module m_lattice_langevin_mover
41 : use defs_basis
42 : use m_abicore
43 : use m_errors
44 :
45 : use m_multibinit_dataset, only: multibinit_dtset_type
46 : use m_abstract_potential, only: abstract_potential_t
47 : use m_abstract_mover, only: abstract_mover_t
48 : use m_lattice_mover, only: lattice_mover_t
49 : use m_multibinit_cell, only: mbcell_t, mbsupercell_t
50 : use m_random_xoroshiro128plus, only: rng_t
51 : use m_hashtable_strval, only: hash_table_t
52 : !!***
53 :
54 : implicit none
55 :
56 : private
57 :
58 : !-------------------------------------------------------------------!
59 : ! Lattice_langevin_mover_t
60 : !
61 : !-------------------------------------------------------------------!
62 : type, public, extends(lattice_mover_t) :: lattice_langevin_mover_t
63 : ! c1 to c5 are a constants (temperature and mass dependent)
64 : ! c3 c4 c5 has dimension of natom.
65 : real(dp) :: c1, c2
66 : real(dp), allocatable :: c3(:), c4(:), c5(:)
67 : real(dp) :: fr =1e-4 ! friction, usually 1e-2~1e-4
68 : ! xi and eta: random numbers of dimension (3, natom)
69 : real(dp), allocatable :: xi(:,:), eta(:,:)
70 : contains
71 : procedure :: initialize
72 : procedure :: finalize
73 : procedure :: run_one_step
74 : procedure :: update_vars
75 : end type lattice_langevin_mover_t
76 :
77 : contains
78 :
79 : !-------------------------------------------------------------------!
80 : ! Initialize:
81 : ! read parameters
82 : ! allocate memory and call update_vars()
83 : !-------------------------------------------------------------------!
84 1 : subroutine initialize(self,params, supercell, rng)
85 : class(lattice_langevin_mover_t), intent(inout) :: self
86 : type(multibinit_dtset_type), target, intent(in):: params
87 : type(mbsupercell_t), target, intent(in) :: supercell
88 : type(rng_t), target, intent(in) :: rng
89 1 : call self%lattice_mover_t%initialize(params, supercell, rng)
90 :
91 : ! TODO: add friction
92 1 : self%fr = params%latt_friction
93 :
94 3 : ABI_MALLOC(self%c3, (self%natom))
95 3 : ABI_MALLOC(self%c4, (self%natom))
96 3 : ABI_MALLOC(self%c5, (self%natom))
97 3 : ABI_MALLOC(self%xi, (3,self%natom))
98 3 : ABI_MALLOC(self%eta, (3,self%natom))
99 :
100 1 : call self%update_vars()
101 1 : end subroutine initialize
102 :
103 :
104 : !-------------------------------------------------------------------!
105 : ! Finalize
106 : !-------------------------------------------------------------------!
107 1 : subroutine finalize(self)
108 : class(lattice_langevin_mover_t), intent(inout) :: self
109 1 : ABI_FREE(self%c3)
110 1 : ABI_FREE(self%c4)
111 1 : ABI_FREE(self%c5)
112 1 : ABI_FREE(self%xi)
113 1 : ABI_FREE(self%eta)
114 1 : call self%lattice_mover_t%finalize()
115 1 : end subroutine finalize
116 :
117 : !-------------------------------------------------------------------!
118 : ! update_vars:
119 : ! calculate c1 to c5 from dt, masses and temperature
120 : ! It is called by the initialization function.
121 : !-------------------------------------------------------------------!
122 1 : subroutine update_vars(self)
123 : class(lattice_langevin_mover_t), intent(inout) :: self
124 2 : real(dp) :: dt, T, fr, sigma(self%natom)
125 :
126 1 : dt=self%dt
127 1 : T= self%temperature
128 1 : fr=self%fr
129 321 : sigma(:) = sqrt(2.0*T*fr/self%masses(:))
130 1 : self%c1 = dt / 2.0 - dt * dt * fr / 8.0
131 1 : self%c2 = dt * fr / 2.0 - dt * dt * fr * fr / 8.0
132 322 : self%c3 = sqrt(dt) * sigma / 2.0 - dt**1.5 * fr * sigma / 8.0
133 322 : self%c5 = dt**1.5 * sigma / (2.0 * sqrt(3.0))
134 322 : self%c4 = fr / 2. * self%c5
135 1 : end subroutine update_vars
136 :
137 : !-------------------------------------------------------------------!
138 : !> @brief: Run_one_step using a Langevin heat bath.
139 : !> effpot: the potential (which do the calculation of E and dE/dvar)
140 : !> param[in]: effpot
141 : !> param[in]: (optional) displacement
142 : !> param[in]: (optional) strain
143 : !> param[in]: (optional) spin
144 : !> param[in]: (optional) lwf
145 : ! NOTE: No need to pass the variable already saved in the mover.
146 : ! e.g. For spin mover, do NOT pass the spin to it.
147 : ! The other variables are only required if there is coupling with
148 : ! the mover variable.
149 : !-------------------------------------------------------------------!
150 10 : subroutine run_one_step(self, effpot,displacement, strain, spin, lwf , energy_table)
151 : class(lattice_langevin_mover_t), intent(inout) :: self
152 : class(abstract_potential_t), intent(inout) :: effpot
153 : real(dp), optional, intent(inout) :: displacement(:,:), strain(:,:), spin(:,:), lwf(:)
154 : type(hash_table_t), optional, intent(inout) :: energy_table
155 : integer :: i
156 : character(len=40) :: key
157 :
158 : ! do not use displacement and strain because they are stored in the mover.
159 10 : ABI_UNUSED_A(displacement)
160 10 : ABI_UNUSED_A(strain)
161 :
162 10 : self%energy = 0.0
163 12810 : self%forces(:,:) =0.0
164 : call effpot%calculate( displacement=self%displacement, strain=self%strain, &
165 : & spin=spin, lwf=lwf, force=self%forces, stress=self%stress, &
166 30 : & energy=self%energy, energy_table=energy_table)
167 10 : call self%rng%rand_normal_array(self%xi, 3*self%natom)
168 10 : call self%rng%rand_normal_array(self%eta, 3*self%natom)
169 :
170 :
171 :
172 : ! First half of velocity update
173 3210 : do i =1, self%natom
174 : self%current_vcart(:,i) = self%current_vcart(:,i) + &
175 : & self%c1 * self%forces(:,i) / self%masses(i) - &
176 : & self%c2 * self%current_vcart(:,i) + &
177 12800 : & self%c3(i) * self%xi(:, i) - self%c4(i) * self%eta(:,i)
178 :
179 : self%displacement(:, i) = self%displacement(:, i) &
180 12810 : & + self%dt * self%current_vcart(:, i) + self%c5( i) *self%eta(:,i)
181 : end do
182 :
183 : ! second half, update the velocity but not the displacement.
184 10 : self%energy=0.0
185 12810 : self%forces=0.0
186 : call effpot%calculate( displacement=self%displacement, strain=self%strain, &
187 : & spin=spin, lwf=lwf, force=self%forces, stress=self%stress, &
188 30 : & energy=self%energy, energy_table=energy_table)
189 3210 : do i =1, self%natom
190 : self%current_vcart(:,i) = self%current_vcart(:,i) + &
191 : &self%c1 * self%forces(:,i) / self%masses(i) - &
192 : & self%c2 * self%current_vcart(:,i) + &
193 12810 : & self%c3(i) * self%xi(:, i) - self%c4(i) * self%eta(:,i)
194 : end do
195 :
196 10 : call self%get_T_and_Ek()
197 10 : if (present(energy_table)) then
198 10 : key = 'Lattice kinetic energy'
199 10 : call energy_table%put(key, self%Ek)
200 : end if
201 :
202 10 : end subroutine run_one_step
203 :
204 4 : end module m_lattice_langevin_mover
205 :
|
