28#include "../../base/base_uses.f90"
33 LOGICAL,
PARAMETER :: debug_this_module = .false.
36 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'al_system_dynamics'
53 SUBROUTINE al_particles(al, force_env, molecule_kind_set, molecule_set, &
54 particle_set, local_molecules, local_particles, group, vel)
63 REAL(kind=
dp),
INTENT(INOUT),
OPTIONAL :: vel(:, :)
65 CHARACTER(len=*),
PARAMETER :: routinen =
'al_particles'
68 LOGICAL :: my_shell_adiabatic
71 CALL timeset(routinen, handle)
72 my_shell_adiabatic = .false.
73 map_info => al%map_info
75 IF (debug_this_module) &
76 CALL dump_vel(molecule_kind_set, molecule_set, local_molecules, particle_set, vel,
"INIT")
78 IF (al%tau_nh <= 0.0_dp)
THEN
79 CALL al_ou_step(0.5_dp, al, force_env, map_info, molecule_kind_set, molecule_set, &
80 particle_set, local_molecules, local_particles, vel)
81 IF (debug_this_module) &
82 CALL dump_vel(molecule_kind_set, molecule_set, local_molecules, particle_set, vel,
"post OU")
85 CALL al_ou_step(0.25_dp, al, force_env, map_info, molecule_kind_set, molecule_set, &
86 particle_set, local_molecules, local_particles, vel)
87 IF (debug_this_module) &
88 CALL dump_vel(molecule_kind_set, molecule_set, local_molecules, particle_set, vel,
"post 1st OU")
92 local_molecules, molecule_set, group, vel=vel)
94 CALL al_nh_quarter_step(al, map_info, set_half_step_vel_factors=.true.)
98 local_molecules, my_shell_adiabatic, vel=vel)
101 local_molecules, molecule_set, group, vel=vel)
102 IF (debug_this_module) &
103 CALL dump_vel(molecule_kind_set, molecule_set, local_molecules, particle_set, vel,
"post rescale_vel")
106 CALL al_nh_quarter_step(al, map_info, set_half_step_vel_factors=.false.)
109 CALL al_ou_step(0.25_dp, al, force_env, map_info, molecule_kind_set, molecule_set, &
110 particle_set, local_molecules, local_particles, vel)
111 IF (debug_this_module) &
112 CALL dump_vel(molecule_kind_set, molecule_set, local_molecules, particle_set, vel,
"post 2nd OU")
117 local_molecules, molecule_set, group, vel=vel)
119 CALL timestop(handle)
131 SUBROUTINE dump_vel(molecule_kind_set, molecule_set, local_molecules, particle_set, vel, label)
136 REAL(
dp),
OPTIONAL :: vel(:, :)
137 CHARACTER(len=*) :: label
139 INTEGER :: first_atom, ikind, imol, imol_local, &
140 ipart, last_atom, nmol_local
143 DO ikind = 1,
SIZE(molecule_kind_set)
144 nmol_local = local_molecules%n_el(ikind)
145 DO imol_local = 1, nmol_local
146 imol = local_molecules%list(ikind)%array(imol_local)
147 molecule => molecule_set(imol)
148 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
149 DO ipart = first_atom, last_atom
150 IF (
PRESENT(vel))
THEN
151 WRITE (unit=*, fmt=
'("VEL ",A20," IPART ",I6," V ",3F20.10)') trim(label), ipart, vel(:, ipart)
153 WRITE (unit=*, fmt=
'("PARTICLE_SET%VEL ",A20," IPART ",I6," V ",3F20.10)') trim(label), &
154 ipart, particle_set(ipart)%v(:)
159 END SUBROUTINE dump_vel
174 SUBROUTINE al_ou_step(step, al, force_env, map_info, molecule_kind_set, molecule_set, &
175 particle_set, local_molecules, local_particles, vel)
176 REAL(
dp),
INTENT(in) :: step
184 REAL(kind=
dp),
INTENT(INOUT),
OPTIONAL :: vel(:, :)
186 INTEGER :: first_atom, i, ii, ikind, imap, imol, imol_local, ipart, iparticle_kind, &
187 iparticle_local, jj, last_atom, nmol_local, nparticle, nparticle_kind, nparticle_local
188 LOGICAL :: check, present_vel
189 REAL(kind=
dp) :: mass
190 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: w(:, :)
194 present_vel =
PRESENT(vel)
197 DO i = 1, al%loc_num_al
198 imap = map_info%map_index(i)
200 IF (al%tau_langevin > 0.0_dp)
THEN
201 map_info%v_scale(imap) = exp(-step*al%dt/al%tau_langevin)
202 map_info%s_kin(imap) = sqrt((al%nvt(i)%nkt/al%nvt(i)%degrees_of_freedom)*(1.0_dp - map_info%v_scale(imap)**2))
204 map_info%v_scale(imap) = 1.0_dp
205 map_info%s_kin(imap) = 0.0_dp
210 nparticle =
SIZE(particle_set)
211 nparticle_kind =
SIZE(local_particles%n_el)
212 ALLOCATE (w(3, nparticle))
214 check = (nparticle_kind <=
SIZE(local_particles%n_el) .AND. nparticle_kind <=
SIZE(local_particles%list))
216 check =
ASSOCIATED(local_particles%local_particle_set)
218 DO iparticle_kind = 1, nparticle_kind
219 nparticle_local = local_particles%n_el(iparticle_kind)
220 check = (nparticle_local <=
SIZE(local_particles%list(iparticle_kind)%array))
222 DO iparticle_local = 1, nparticle_local
223 ipart = local_particles%list(iparticle_kind)%array(iparticle_local)
224 w(1, ipart) = local_particles%local_particle_set(iparticle_kind)%rng(iparticle_local)%stream%next(variance=1.0_dp)
225 w(2, ipart) = local_particles%local_particle_set(iparticle_kind)%rng(iparticle_local)%stream%next(variance=1.0_dp)
226 w(3, ipart) = local_particles%local_particle_set(iparticle_kind)%rng(iparticle_local)%stream%next(variance=1.0_dp)
233 DO ikind = 1,
SIZE(molecule_kind_set)
234 nmol_local = local_molecules%n_el(ikind)
235 DO imol_local = 1, nmol_local
236 imol = local_molecules%list(ikind)%array(imol_local)
237 molecule => molecule_set(imol)
238 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
239 DO ipart = first_atom, last_atom
241 atomic_kind => particle_set(ipart)%atomic_kind
243 IF (present_vel)
THEN
245 vel(jj, ipart) = vel(jj, ipart)*map_info%p_scale(jj, ii)%point + &
246 map_info%p_kin(jj, ii)%point/sqrt(mass)*w(jj, ipart)
250 particle_set(ipart)%v(jj) = particle_set(ipart)%v(jj)*map_info%p_scale(jj, ii)%point + &
251 map_info%p_kin(jj, ii)%point/sqrt(mass)*w(jj, ipart)
260 END SUBROUTINE al_ou_step
269 SUBROUTINE al_nh_quarter_step(al, map_info, set_half_step_vel_factors)
272 LOGICAL,
INTENT(in) :: set_half_step_vel_factors
275 REAL(kind=
dp) :: decay, delta_k
279 DO i = 1, al%loc_num_al
280 IF (al%nvt(i)%mass > 0.0_dp)
THEN
281 imap = map_info%map_index(i)
282 delta_k = 0.5_dp*(map_info%s_kin(imap) - al%nvt(i)%nkt)
283 al%nvt(i)%chi = al%nvt(i)%chi + 0.5_dp*al%dt*delta_k/al%nvt(i)%mass
284 IF (set_half_step_vel_factors)
THEN
285 decay = exp(-0.5_dp*al%dt*al%nvt(i)%chi)
286 map_info%v_scale(imap) = decay
289 al%nvt(i)%chi = 0.0_dp
290 IF (set_half_step_vel_factors)
THEN
291 map_info%v_scale(imap) = 1.0_dp
296 END SUBROUTINE al_nh_quarter_step
subroutine, public al_particles(al, force_env, molecule_kind_set, molecule_set, particle_set, local_molecules, local_particles, group, vel)
...
Type for the canonical sampling through velocity rescaling.
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
subroutine, public fix_atom_control(force_env, w)
allows for fix atom constraints
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
Lumps all possible extended system variables into one type for easy access and passing.
Interface for the force calculations.
Defines the basic variable types.
integer, parameter, public dp
Interface to the message passing library MPI.
Define the molecule kind structure types and the corresponding functionality.
Define the data structure for the molecule information.
subroutine, public get_molecule(molecule, molecule_kind, lmi, lci, lg3x3, lg4x6, lcolv, first_atom, last_atom, first_shell, last_shell)
Get components from a molecule data set.
Define the data structure for the particle information.
Utilities for thermostats.
subroutine, public vel_rescale_particles(map_info, molecule_kind_set, molecule_set, particle_set, local_molecules, shell_adiabatic, shell_particle_set, core_particle_set, vel, shell_vel, core_vel)
...
subroutine, public ke_region_particles(map_info, particle_set, molecule_kind_set, local_molecules, molecule_set, group, vel)
...
Provides all information about an atomic kind.
structure to store local (to a processor) ordered lists of integers.
wrapper to abstract the force evaluation of the various methods