(git:0de0cc2)
csvr_system_dynamics.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \author 10.2007 [tlaino] - Teodoro Laino - University of Zurich
10 ! **************************************************************************************************
12 
13  USE atomic_kind_types, ONLY: atomic_kind_type
14  USE csvr_system_types, ONLY: csvr_system_type
16  USE distribution_1d_types, ONLY: distribution_1d_type
17  USE extended_system_types, ONLY: map_info_type,&
18  npt_info_type
19  USE kinds, ONLY: dp
20  USE message_passing, ONLY: mp_comm_type
21  USE molecule_kind_types, ONLY: molecule_kind_type
22  USE molecule_types, ONLY: molecule_type
23  USE parallel_rng_types, ONLY: rng_stream_type
24  USE particle_types, ONLY: particle_type
25  USE thermostat_utils, ONLY: ke_region_baro,&
31 #include "../../base/base_uses.f90"
32 
33  IMPLICIT NONE
34 
35  PRIVATE
36  LOGICAL, PARAMETER :: debug_this_module = .false.
37  PUBLIC :: csvr_particles, &
38  csvr_barostat, &
40 
41  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'csvr_system_dynamics'
42 
43 CONTAINS
44 
45 ! **************************************************************************************************
46 !> \brief ...
47 !> \param csvr ...
48 !> \param npt ...
49 !> \param group ...
50 !> \author 10.2007 [tlaino] - Teodoro Laino - University of Zurich
51 ! **************************************************************************************************
52  SUBROUTINE csvr_barostat(csvr, npt, group)
53 
54  TYPE(csvr_system_type), POINTER :: csvr
55  TYPE(npt_info_type), DIMENSION(:, :), &
56  INTENT(INOUT) :: npt
57  TYPE(mp_comm_type), INTENT(IN) :: group
58 
59  CHARACTER(len=*), PARAMETER :: routinen = 'csvr_barostat'
60 
61  INTEGER :: handle
62  TYPE(map_info_type), POINTER :: map_info
63 
64  CALL timeset(routinen, handle)
65  map_info => csvr%map_info
66 
67  ! Compute the kinetic energy of the barostat
68  CALL ke_region_baro(map_info, npt, group)
69 
70  ! Apply the Canonical Sampling through Velocity Rescaling
71  CALL do_csvr(csvr, map_info)
72 
73  ! Now scale the particle velocities
74  CALL vel_rescale_baro(map_info, npt)
75 
76  ! Re-Compute the kinetic energy of the barostat
77  CALL ke_region_baro(map_info, npt, group)
78 
79  ! Compute thermostat energy
80  CALL do_csvr_eval_energy(csvr, map_info)
81 
82  CALL timestop(handle)
83  END SUBROUTINE csvr_barostat
84 
85 ! **************************************************************************************************
86 !> \brief ...
87 !> \param csvr ...
88 !> \param molecule_kind_set ...
89 !> \param molecule_set ...
90 !> \param particle_set ...
91 !> \param local_molecules ...
92 !> \param group ...
93 !> \param shell_adiabatic ...
94 !> \param shell_particle_set ...
95 !> \param core_particle_set ...
96 !> \param vel ...
97 !> \param shell_vel ...
98 !> \param core_vel ...
99 !> \author 10.2007 [tlaino] - Teodoro Laino - University of Zurich
100 ! **************************************************************************************************
101  SUBROUTINE csvr_particles(csvr, molecule_kind_set, molecule_set, &
102  particle_set, local_molecules, group, shell_adiabatic, &
103  shell_particle_set, core_particle_set, vel, shell_vel, core_vel)
104 
105  TYPE(csvr_system_type), POINTER :: csvr
106  TYPE(molecule_kind_type), POINTER :: molecule_kind_set(:)
107  TYPE(molecule_type), POINTER :: molecule_set(:)
108  TYPE(particle_type), POINTER :: particle_set(:)
109  TYPE(distribution_1d_type), POINTER :: local_molecules
110  TYPE(mp_comm_type), INTENT(IN) :: group
111  LOGICAL, INTENT(IN), OPTIONAL :: shell_adiabatic
112  TYPE(particle_type), OPTIONAL, POINTER :: shell_particle_set(:), &
113  core_particle_set(:)
114  REAL(kind=dp), INTENT(INOUT), OPTIONAL :: vel(:, :), shell_vel(:, :), &
115  core_vel(:, :)
116 
117  CHARACTER(len=*), PARAMETER :: routinen = 'csvr_particles'
118 
119  INTEGER :: handle
120  LOGICAL :: my_shell_adiabatic
121  TYPE(map_info_type), POINTER :: map_info
122 
123  CALL timeset(routinen, handle)
124  my_shell_adiabatic = .false.
125  IF (PRESENT(shell_adiabatic)) my_shell_adiabatic = shell_adiabatic
126  map_info => csvr%map_info
127 
128  ! Compute the kinetic energy for the region to thermostat
129  CALL ke_region_particles(map_info, particle_set, molecule_kind_set, &
130  local_molecules, molecule_set, group, vel)
131 
132  ! Apply the Canonical Sampling through Velocity Rescaling
133  CALL do_csvr(csvr, map_info)
134 
135  ! Now scale the particle velocities
136  CALL vel_rescale_particles(map_info, molecule_kind_set, molecule_set, particle_set, &
137  local_molecules, my_shell_adiabatic, shell_particle_set, core_particle_set, &
138  vel, shell_vel, core_vel)
139 
140  ! Re-Compute the kinetic energy for the region to thermostat
141  CALL ke_region_particles(map_info, particle_set, molecule_kind_set, &
142  local_molecules, molecule_set, group, vel)
143 
144  ! Compute thermostat energy
145  CALL do_csvr_eval_energy(csvr, map_info)
146 
147  CALL timestop(handle)
148  END SUBROUTINE csvr_particles
149 
150 ! **************************************************************************************************
151 !> \brief ...
152 !> \param csvr ...
153 !> \param atomic_kind_set ...
154 !> \param particle_set ...
155 !> \param local_particles ...
156 !> \param group ...
157 !> \param shell_particle_set ...
158 !> \param core_particle_set ...
159 !> \param vel ...
160 !> \param shell_vel ...
161 !> \param core_vel ...
162 !> \author 10.2007 [tlaino] - Teodoro Laino - University of Zurich
163 ! **************************************************************************************************
164  SUBROUTINE csvr_shells(csvr, atomic_kind_set, particle_set, local_particles, &
165  group, shell_particle_set, core_particle_set, vel, shell_vel, core_vel)
166 
167  TYPE(csvr_system_type), POINTER :: csvr
168  TYPE(atomic_kind_type), POINTER :: atomic_kind_set(:)
169  TYPE(particle_type), POINTER :: particle_set(:)
170  TYPE(distribution_1d_type), POINTER :: local_particles
171  TYPE(mp_comm_type), INTENT(IN) :: group
172  TYPE(particle_type), OPTIONAL, POINTER :: shell_particle_set(:), &
173  core_particle_set(:)
174  REAL(kind=dp), INTENT(INOUT), OPTIONAL :: vel(:, :), shell_vel(:, :), &
175  core_vel(:, :)
176 
177  CHARACTER(len=*), PARAMETER :: routinen = 'csvr_shells'
178 
179  INTEGER :: handle
180  TYPE(map_info_type), POINTER :: map_info
181 
182  CALL timeset(routinen, handle)
183  map_info => csvr%map_info
184 
185  ! Compute the kinetic energy of the region to thermostat
186  CALL ke_region_shells(map_info, particle_set, atomic_kind_set, local_particles, &
187  group, core_particle_set, shell_particle_set, core_vel, shell_vel)
188 
189  ! Apply the Canonical Sampling through Velocity Rescaling
190  CALL do_csvr(csvr, map_info)
191 
192  ! Now scale the particle velocities
193  CALL vel_rescale_shells(map_info, atomic_kind_set, particle_set, local_particles, &
194  shell_particle_set, core_particle_set, shell_vel, core_vel, vel)
195 
196  ! Re-Compute the kinetic energy of the region to thermostat
197  CALL ke_region_shells(map_info, particle_set, atomic_kind_set, local_particles, &
198  group, core_particle_set, shell_particle_set, core_vel, shell_vel)
199 
200  ! Compute thermostat energy
201  CALL do_csvr_eval_energy(csvr, map_info)
202 
203  CALL timestop(handle)
204  END SUBROUTINE csvr_shells
205 
206 ! **************************************************************************************************
207 !> \brief ...
208 !> \param csvr ...
209 !> \param map_info ...
210 !> \author 10.2007 [tlaino] - Teodoro Laino - University of Zurich
211 ! **************************************************************************************************
212  SUBROUTINE do_csvr(csvr, map_info)
213  TYPE(csvr_system_type), POINTER :: csvr
214  TYPE(map_info_type), POINTER :: map_info
215 
216  INTEGER :: i, imap, ndeg
217  REAL(kind=dp) :: kin_energy, kin_target, taut
218  TYPE(rng_stream_type), POINTER :: rng_stream
219 
220  DO i = 1, csvr%loc_num_csvr
221  imap = map_info%map_index(i)
222  kin_energy = map_info%s_kin(imap)
223  csvr%nvt(i)%region_kin_energy = kin_energy
224  kin_energy = kin_energy*0.5_dp
225  kin_target = csvr%nvt(i)%nkt*0.5_dp
226  ndeg = csvr%nvt(i)%degrees_of_freedom
227  taut = csvr%tau_csvr/csvr%dt_fact
228  rng_stream => csvr%nvt(i)%gaussian_rng_stream
229  map_info%v_scale(imap) = rescaling_factor(kin_energy, kin_target, ndeg, taut, &
230  rng_stream)
231  END DO
232 
233  END SUBROUTINE do_csvr
234 
235 ! **************************************************************************************************
236 !> \brief ...
237 !> \param csvr ...
238 !> \param map_info ...
239 !> \author 10.2007 [tlaino] - Teodoro Laino - University of Zurich
240 ! **************************************************************************************************
241  SUBROUTINE do_csvr_eval_energy(csvr, map_info)
242  TYPE(csvr_system_type), POINTER :: csvr
243  TYPE(map_info_type), POINTER :: map_info
244 
245  INTEGER :: i, imap
246  REAL(kind=dp) :: kin_energy_ar, kin_energy_br
247 
248  DO i = 1, csvr%loc_num_csvr
249  imap = map_info%map_index(i)
250  kin_energy_br = csvr%nvt(i)%region_kin_energy
251  kin_energy_ar = map_info%s_kin(imap)
252  csvr%nvt(i)%thermostat_energy = csvr%nvt(i)%thermostat_energy + &
253  0.5_dp*(kin_energy_br - kin_energy_ar)
254  END DO
255 
256  END SUBROUTINE do_csvr_eval_energy
257 
258 END MODULE csvr_system_dynamics
Define the atomic kind types and their sub types.
subroutine, public csvr_shells(csvr, atomic_kind_set, particle_set, local_particles, group, shell_particle_set, core_particle_set, vel, shell_vel, core_vel)
...
subroutine, public csvr_particles(csvr, molecule_kind_set, molecule_set, particle_set, local_molecules, group, shell_adiabatic, shell_particle_set, core_particle_set, vel, shell_vel, core_vel)
...
subroutine, public csvr_barostat(csvr, npt, group)
...
Type for the canonical sampling through velocity rescaling.
real(kind=dp) function, public rescaling_factor(kk, sigma, ndeg, taut, rng_stream)
Stochastic velocity rescale, as described in Bussi, Donadio and Parrinello, J. Chem....
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.
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public dp
Definition: kinds.F:34
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.
Parallel (pseudo)random number generator (RNG) for multiple streams and substreams of random numbers.
Define the data structure for the particle information.
Utilities for thermostats.
subroutine, public vel_rescale_shells(map_info, atomic_kind_set, particle_set, local_particles, shell_particle_set, core_particle_set, shell_vel, core_vel, vel)
...
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_shells(map_info, particle_set, atomic_kind_set, local_particles, group, core_particle_set, shell_particle_set, core_vel, shell_vel)
...
subroutine, public ke_region_baro(map_info, npt, group)
...
subroutine, public vel_rescale_baro(map_info, npt)
...
subroutine, public ke_region_particles(map_info, particle_set, molecule_kind_set, local_molecules, molecule_set, group, vel)
...