(git:58e3e09)
extended_system_mapping.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 !> \par History
10 !> CJM, 20-Feb-01
11 !> JGH (10-Mar-2001)
12 !> CJM (10-Apr-2001)
13 !> \author CJM
14 ! **************************************************************************************************
16 
17  USE distribution_1d_types, ONLY: distribution_1d_type
19  lnhc_parameters_type,&
20  map_info_type
21  USE input_constants, ONLY: &
26  USE kinds, ONLY: dp
27  USE message_passing, ONLY: mp_para_env_type
28  USE molecule_kind_types, ONLY: molecule_kind_type
29  USE molecule_types, ONLY: global_constraint_type,&
30  molecule_type
31  USE simpar_types, ONLY: simpar_type
35  USE thermostat_types, ONLY: thermostat_info_type
36 #include "../../base/base_uses.f90"
37 
38  IMPLICIT NONE
39 
40  PRIVATE
41 
42  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'extended_system_mapping'
43 
47 
48 CONTAINS
49 
50 ! **************************************************************************************************
51 !> \brief Creates the thermostatting for the barostat
52 !> \param simpar ...
53 !> \param nhc ...
54 !> \par History
55 !> CJM, 20-Feb-01 : nhc structure allocated to zero when not in use
56 !> JGH (10-Mar-2001) : set nhc variables to zero when not in use
57 !> \author CJM
58 ! **************************************************************************************************
59  SUBROUTINE nhc_to_barostat_mapping(simpar, nhc)
60 
61  TYPE(simpar_type), POINTER :: simpar
62  TYPE(lnhc_parameters_type), POINTER :: nhc
63 
64  CHARACTER(LEN=*), PARAMETER :: routinen = 'nhc_to_barostat_mapping'
65 
66  INTEGER :: handle, i, number
67  TYPE(map_info_type), POINTER :: map_info
68 
69  CALL timeset(routinen, handle)
70 
71  SELECT CASE (simpar%ensemble)
72  CASE DEFAULT
73  cpabort('Never reach this point!')
75  map_info => nhc%map_info
76  map_info%dis_type = do_thermo_only_master
77 
78  ! Counting the total number of thermostats ( 1 for NPT_I, NPT_IA, and NPT_F )
79  nhc%loc_num_nhc = 1
80  nhc%glob_num_nhc = 1
81  IF (simpar%ensemble == npt_f_ensemble) THEN
82  number = 9
83  ELSE
84  number = 1
85  END IF
86 
87  CALL init_baro_map_info(map_info, number, nhc%loc_num_nhc)
88 
89  ALLOCATE (nhc%nvt(nhc%nhc_len, nhc%loc_num_nhc))
90  ! Now that we know how many there are stick this into nhc % nkt
91  ! (number of degrees of freedom times k_B T )
92  DO i = 1, nhc%loc_num_nhc
93  nhc%nvt(1, i)%nkt = simpar%temp_ext*number
94  nhc%nvt(1, i)%degrees_of_freedom = number
95  IF (debug_isotropic_limit) THEN
96  nhc%nvt(1, i)%nkt = simpar%temp_ext
97  END IF
98  END DO
99 
100  ! getting the number of degrees of freedom times k_B T for the rest of the chain
101  DO i = 2, nhc%nhc_len
102  nhc%nvt(i, :)%nkt = simpar%temp_ext
103  END DO
104 
105  ! Let's clean the arrays
106  map_info%s_kin = 0.0_dp
107  map_info%v_scale = 0.0_dp
108  END SELECT
109 
110  CALL timestop(handle)
111 
112  END SUBROUTINE nhc_to_barostat_mapping
113 
114 ! **************************************************************************************************
115 !> \brief Creates the thermostatting maps
116 !> \param thermostat_info ...
117 !> \param simpar ...
118 !> \param local_molecules ...
119 !> \param molecule_set ...
120 !> \param molecule_kind_set ...
121 !> \param nhc ...
122 !> \param para_env ...
123 !> \param gci ...
124 !> \par History
125 !> 29-Nov-00 (JGH) correct counting of DOF if constraints are off
126 !> CJM, 20-Feb-01 : nhc structure allocated to zero when not in use
127 !> JGH (10-Mar-2001) : set nhc variables to zero when not in use
128 !> CJM(10-NOV-2001) : New parallelization with new molecule structures
129 !> Teodoro Laino 09.2007 [tlaino] - University of Zurich - cleaning and updating
130 !> \author CJM
131 ! **************************************************************************************************
132  SUBROUTINE nhc_to_particle_mapping(thermostat_info, simpar, local_molecules, &
133  molecule_set, molecule_kind_set, nhc, para_env, gci)
134 
135  TYPE(thermostat_info_type), POINTER :: thermostat_info
136  TYPE(simpar_type), POINTER :: simpar
137  TYPE(distribution_1d_type), POINTER :: local_molecules
138  TYPE(molecule_type), POINTER :: molecule_set(:)
139  TYPE(molecule_kind_type), POINTER :: molecule_kind_set(:)
140  TYPE(lnhc_parameters_type), POINTER :: nhc
141  TYPE(mp_para_env_type), POINTER :: para_env
142  TYPE(global_constraint_type), POINTER :: gci
143 
144  CHARACTER(LEN=*), PARAMETER :: routinen = 'nhc_to_particle_mapping'
145 
146  INTEGER :: handle, i, imap, j, natoms_local, &
147  sum_of_thermostats
148  INTEGER, DIMENSION(:), POINTER :: deg_of_freedom, massive_atom_list
149  REAL(kind=dp) :: fac
150  TYPE(map_info_type), POINTER :: map_info
151 
152  CALL timeset(routinen, handle)
153 
154  NULLIFY (massive_atom_list, deg_of_freedom)
155 
156  SELECT CASE (simpar%ensemble)
157  CASE DEFAULT
158  cpabort('Unknown ensemble!')
161  cpabort('Never reach this point!')
163 
164  CALL setup_nhc_thermostat(nhc, thermostat_info, deg_of_freedom, massive_atom_list, &
165  molecule_kind_set, local_molecules, molecule_set, para_env, natoms_local, &
166  simpar, sum_of_thermostats, gci)
167 
168  ! Sum up the number of degrees of freedom on each thermostat.
169  ! first: initialize the target
170  map_info => nhc%map_info
171  map_info%s_kin = 0.0_dp
172  DO i = 1, 3
173  DO j = 1, natoms_local
174  map_info%p_kin(i, j)%point = map_info%p_kin(i, j)%point + 1
175  END DO
176  END DO
177 
178  ! if thermostats are replicated but molecules distributed, we have to
179  ! sum s_kin over all processors
180  IF (map_info%dis_type == do_thermo_communication) CALL para_env%sum(map_info%s_kin)
181 
182  ! We know the total number of system thermostats.
183  IF ((sum_of_thermostats == 1) .AND. (map_info%dis_type /= do_thermo_no_communication)) THEN
184  fac = map_info%s_kin(1) - deg_of_freedom(1) - simpar%nfree_rot_transl
185  IF (fac == 0.0_dp) THEN
186  cpabort('Zero degrees of freedom. Nothing to thermalize!')
187  END IF
188  nhc%nvt(1, 1)%nkt = simpar%temp_ext*fac
189  nhc%nvt(1, 1)%degrees_of_freedom = floor(fac)
190  ELSE
191  DO i = 1, nhc%loc_num_nhc
192  imap = map_info%map_index(i)
193  fac = (map_info%s_kin(imap) - deg_of_freedom(i))
194  nhc%nvt(1, i)%nkt = simpar%temp_ext*fac
195  nhc%nvt(1, i)%degrees_of_freedom = floor(fac)
196  END DO
197  END IF
198 
199  ! Getting the number of degrees of freedom times k_B T for the rest
200  ! of the chain
201  DO i = 2, nhc%nhc_len
202  nhc%nvt(i, :)%nkt = simpar%temp_ext
203  nhc%nvt(i, :)%degrees_of_freedom = 1
204  END DO
205  DEALLOCATE (deg_of_freedom)
206  DEALLOCATE (massive_atom_list)
207 
208  ! Let's clean the arrays
209  map_info%s_kin = 0.0_dp
210  map_info%v_scale = 0.0_dp
211  END SELECT
212 
213  CALL timestop(handle)
214 
215  END SUBROUTINE nhc_to_particle_mapping
216 
217 ! **************************************************************************************************
218 !> \brief Main general setup for Adiabatic Nose-Hoover thermostats
219 !> \param nhc ...
220 !> \param thermostat_info ...
221 !> \param deg_of_freedom ...
222 !> \param massive_atom_list ...
223 !> \param molecule_kind_set ...
224 !> \param local_molecules ...
225 !> \param molecule_set ...
226 !> \param para_env ...
227 !> \param natoms_local ...
228 !> \param simpar ...
229 !> \param sum_of_thermostats ...
230 !> \param gci ...
231 !> \param shell ...
232 !> \author CJM -PNNL -2011
233 ! **************************************************************************************************
234  SUBROUTINE setup_adiabatic_thermostat(nhc, thermostat_info, deg_of_freedom, &
235  massive_atom_list, molecule_kind_set, local_molecules, molecule_set, &
236  para_env, natoms_local, simpar, sum_of_thermostats, gci, shell)
237 
238  TYPE(lnhc_parameters_type), POINTER :: nhc
239  TYPE(thermostat_info_type), POINTER :: thermostat_info
240  INTEGER, DIMENSION(:), POINTER :: deg_of_freedom, massive_atom_list
241  TYPE(molecule_kind_type), POINTER :: molecule_kind_set(:)
242  TYPE(distribution_1d_type), POINTER :: local_molecules
243  TYPE(molecule_type), POINTER :: molecule_set(:)
244  TYPE(mp_para_env_type), POINTER :: para_env
245  INTEGER, INTENT(OUT) :: natoms_local
246  TYPE(simpar_type), POINTER :: simpar
247  INTEGER, INTENT(OUT) :: sum_of_thermostats
248  TYPE(global_constraint_type), POINTER :: gci
249  LOGICAL, INTENT(IN), OPTIONAL :: shell
250 
251  CHARACTER(LEN=*), PARAMETER :: routinen = 'setup_adiabatic_thermostat'
252 
253  INTEGER :: handle, nkind, number, region
254  LOGICAL :: do_shell
255  TYPE(map_info_type), POINTER :: map_info
256 
257  CALL timeset(routinen, handle)
258 
259  do_shell = .false.
260  IF (PRESENT(shell)) do_shell = shell
261  map_info => nhc%map_info
262 
263  nkind = SIZE(molecule_kind_set)
264  sum_of_thermostats = thermostat_info%sum_of_thermostats
265  map_info%dis_type = thermostat_info%dis_type
266  number = thermostat_info%number_of_thermostats
267  region = nhc%region
268 
269  CALL adiabatic_mapping_region(map_info, deg_of_freedom, massive_atom_list, &
270  molecule_kind_set, local_molecules, molecule_set, para_env, natoms_local, &
271  simpar, number, region, gci, do_shell, thermostat_info%map_loc_thermo_gen, &
272  sum_of_thermostats)
273  ALLOCATE (nhc%nvt(nhc%nhc_len, number))
274 
275  ! Now that we know how many there are stick this into nhc%nkt
276  ! (number of degrees of freedom times k_B T for the first thermostat
277  ! on the chain)
278  nhc%loc_num_nhc = number
279  nhc%glob_num_nhc = sum_of_thermostats
280 
281  CALL timestop(handle)
282 
283  END SUBROUTINE setup_adiabatic_thermostat
284 
285 ! **************************************************************************************************
286 !> \brief Creates the thermostatting maps
287 !> \param thermostat_info ...
288 !> \param simpar ...
289 !> \param local_molecules ...
290 !> \param molecule_set ...
291 !> \param molecule_kind_set ...
292 !> \param nhc ...
293 !> \param para_env ...
294 !> \param gci ...
295 !> \par History
296 !> \author CJM
297 ! **************************************************************************************************
298  SUBROUTINE nhc_to_particle_mapping_slow(thermostat_info, simpar, local_molecules, &
299  molecule_set, molecule_kind_set, nhc, para_env, gci)
300 
301  TYPE(thermostat_info_type), POINTER :: thermostat_info
302  TYPE(simpar_type), POINTER :: simpar
303  TYPE(distribution_1d_type), POINTER :: local_molecules
304  TYPE(molecule_type), POINTER :: molecule_set(:)
305  TYPE(molecule_kind_type), POINTER :: molecule_kind_set(:)
306  TYPE(lnhc_parameters_type), POINTER :: nhc
307  TYPE(mp_para_env_type), POINTER :: para_env
308  TYPE(global_constraint_type), POINTER :: gci
309 
310  CHARACTER(LEN=*), PARAMETER :: routinen = 'nhc_to_particle_mapping_slow'
311 
312  INTEGER :: handle, i, imap, j, natoms_local, &
313  sum_of_thermostats
314  INTEGER, DIMENSION(:), POINTER :: deg_of_freedom, massive_atom_list
315  REAL(kind=dp) :: fac
316  TYPE(map_info_type), POINTER :: map_info
317 
318  CALL timeset(routinen, handle)
319 
320  NULLIFY (massive_atom_list, deg_of_freedom)
321 
322  SELECT CASE (simpar%ensemble)
323  CASE DEFAULT
324  cpabort('Unknown ensemble!')
326  CALL setup_adiabatic_thermostat(nhc, thermostat_info, deg_of_freedom, massive_atom_list, &
327  molecule_kind_set, local_molecules, molecule_set, para_env, natoms_local, &
328  simpar, sum_of_thermostats, gci)
329 
330  ! Sum up the number of degrees of freedom on each thermostat.
331  ! first: initialize the target
332  map_info => nhc%map_info
333  map_info%s_kin = 0.0_dp
334  DO i = 1, 3
335  DO j = 1, natoms_local
336  IF (ASSOCIATED(map_info%p_kin(i, j)%point)) &
337  map_info%p_kin(i, j)%point = map_info%p_kin(i, j)%point + 1
338  END DO
339  END DO
340 
341  ! if thermostats are replicated but molecules distributed, we have to
342  ! sum s_kin over all processors
343  IF (map_info%dis_type == do_thermo_communication) CALL para_env%sum(map_info%s_kin)
344 
345  ! We know the total number of system thermostats.
346  IF ((sum_of_thermostats == 1) .AND. (map_info%dis_type /= do_thermo_no_communication)) THEN
347  fac = map_info%s_kin(1) - deg_of_freedom(1) - simpar%nfree_rot_transl
348  IF (fac == 0.0_dp) THEN
349  cpabort('Zero degrees of freedom. Nothing to thermalize!')
350  END IF
351  nhc%nvt(1, 1)%nkt = simpar%temp_slow*fac
352  nhc%nvt(1, 1)%degrees_of_freedom = floor(fac)
353  ELSE
354  DO i = 1, nhc%loc_num_nhc
355  imap = map_info%map_index(i)
356  fac = (map_info%s_kin(imap) - deg_of_freedom(i))
357  nhc%nvt(1, i)%nkt = simpar%temp_slow*fac
358  nhc%nvt(1, i)%degrees_of_freedom = floor(fac)
359  END DO
360  END IF
361 
362  ! Getting the number of degrees of freedom times k_B T for the rest
363  ! of the chain
364  DO i = 2, nhc%nhc_len
365  nhc%nvt(i, :)%nkt = simpar%temp_slow
366  nhc%nvt(i, :)%degrees_of_freedom = 1
367  END DO
368  DEALLOCATE (deg_of_freedom)
369  DEALLOCATE (massive_atom_list)
370 
371  ! Let's clean the arrays
372  map_info%s_kin = 0.0_dp
373  map_info%v_scale = 0.0_dp
374  END SELECT
375 
376  CALL timestop(handle)
377 
378  END SUBROUTINE nhc_to_particle_mapping_slow
379 
380 ! **************************************************************************************************
381 !> \brief Creates the thermostatting maps
382 !> \param thermostat_info ...
383 !> \param simpar ...
384 !> \param local_molecules ...
385 !> \param molecule_set ...
386 !> \param molecule_kind_set ...
387 !> \param nhc ...
388 !> \param para_env ...
389 !> \param gci ...
390 !> \par History
391 !> \author CJM
392 ! **************************************************************************************************
393  SUBROUTINE nhc_to_particle_mapping_fast(thermostat_info, simpar, local_molecules, &
394  molecule_set, molecule_kind_set, nhc, para_env, gci)
395 
396  TYPE(thermostat_info_type), POINTER :: thermostat_info
397  TYPE(simpar_type), POINTER :: simpar
398  TYPE(distribution_1d_type), POINTER :: local_molecules
399  TYPE(molecule_type), POINTER :: molecule_set(:)
400  TYPE(molecule_kind_type), POINTER :: molecule_kind_set(:)
401  TYPE(lnhc_parameters_type), POINTER :: nhc
402  TYPE(mp_para_env_type), POINTER :: para_env
403  TYPE(global_constraint_type), POINTER :: gci
404 
405  CHARACTER(LEN=*), PARAMETER :: routinen = 'nhc_to_particle_mapping_fast'
406 
407  INTEGER :: handle, i, imap, j, natoms_local, &
408  sum_of_thermostats
409  INTEGER, DIMENSION(:), POINTER :: deg_of_freedom, massive_atom_list
410  REAL(kind=dp) :: fac
411  TYPE(map_info_type), POINTER :: map_info
412 
413  CALL timeset(routinen, handle)
414 
415  NULLIFY (massive_atom_list, deg_of_freedom)
416 
417  SELECT CASE (simpar%ensemble)
418  CASE DEFAULT
419  cpabort('Unknown ensemble!')
421  CALL setup_adiabatic_thermostat(nhc, thermostat_info, deg_of_freedom, massive_atom_list, &
422  molecule_kind_set, local_molecules, molecule_set, para_env, natoms_local, &
423  simpar, sum_of_thermostats, gci)
424 
425  ! Sum up the number of degrees of freedom on each thermostat.
426  ! first: initialize the target
427  map_info => nhc%map_info
428  map_info%s_kin = 0.0_dp
429  DO i = 1, 3
430  DO j = 1, natoms_local
431  IF (ASSOCIATED(map_info%p_kin(i, j)%point)) &
432  map_info%p_kin(i, j)%point = map_info%p_kin(i, j)%point + 1
433  END DO
434  END DO
435 
436  ! if thermostats are replicated but molecules distributed, we have to
437  ! sum s_kin over all processors
438  IF (map_info%dis_type == do_thermo_communication) CALL para_env%sum(map_info%s_kin)
439 
440  ! We know the total number of system thermostats.
441  IF ((sum_of_thermostats == 1) .AND. (map_info%dis_type /= do_thermo_no_communication)) THEN
442  fac = map_info%s_kin(1) - deg_of_freedom(1) - simpar%nfree_rot_transl
443  IF (fac == 0.0_dp) THEN
444  cpabort('Zero degrees of freedom. Nothing to thermalize!')
445  END IF
446  nhc%nvt(1, 1)%nkt = simpar%temp_fast*fac
447  nhc%nvt(1, 1)%degrees_of_freedom = floor(fac)
448  ELSE
449  DO i = 1, nhc%loc_num_nhc
450  imap = map_info%map_index(i)
451  fac = (map_info%s_kin(imap) - deg_of_freedom(i))
452  nhc%nvt(1, i)%nkt = simpar%temp_fast*fac
453  nhc%nvt(1, i)%degrees_of_freedom = floor(fac)
454  END DO
455  END IF
456 
457  ! Getting the number of degrees of freedom times k_B T for the rest
458  ! of the chain
459  DO i = 2, nhc%nhc_len
460  nhc%nvt(i, :)%nkt = simpar%temp_fast
461  nhc%nvt(i, :)%degrees_of_freedom = 1
462  END DO
463  DEALLOCATE (deg_of_freedom)
464  DEALLOCATE (massive_atom_list)
465 
466  ! Let's clean the arrays
467  map_info%s_kin = 0.0_dp
468  map_info%v_scale = 0.0_dp
469  END SELECT
470 
471  CALL timestop(handle)
472 
473  END SUBROUTINE nhc_to_particle_mapping_fast
474 
475 ! **************************************************************************************************
476 !> \brief Main general setup for Nose-Hoover thermostats
477 !> \param nhc ...
478 !> \param thermostat_info ...
479 !> \param deg_of_freedom ...
480 !> \param massive_atom_list ...
481 !> \param molecule_kind_set ...
482 !> \param local_molecules ...
483 !> \param molecule_set ...
484 !> \param para_env ...
485 !> \param natoms_local ...
486 !> \param simpar ...
487 !> \param sum_of_thermostats ...
488 !> \param gci ...
489 !> \param shell ...
490 !> \author Teodoro Laino [tlaino] - University of Zurich - 10.2007
491 ! **************************************************************************************************
492  SUBROUTINE setup_nhc_thermostat(nhc, thermostat_info, deg_of_freedom, &
493  massive_atom_list, molecule_kind_set, local_molecules, molecule_set, &
494  para_env, natoms_local, simpar, sum_of_thermostats, gci, shell)
495 
496  TYPE(lnhc_parameters_type), POINTER :: nhc
497  TYPE(thermostat_info_type), POINTER :: thermostat_info
498  INTEGER, DIMENSION(:), POINTER :: deg_of_freedom, massive_atom_list
499  TYPE(molecule_kind_type), POINTER :: molecule_kind_set(:)
500  TYPE(distribution_1d_type), POINTER :: local_molecules
501  TYPE(molecule_type), POINTER :: molecule_set(:)
502  TYPE(mp_para_env_type), POINTER :: para_env
503  INTEGER, INTENT(OUT) :: natoms_local
504  TYPE(simpar_type), POINTER :: simpar
505  INTEGER, INTENT(OUT) :: sum_of_thermostats
506  TYPE(global_constraint_type), POINTER :: gci
507  LOGICAL, INTENT(IN), OPTIONAL :: shell
508 
509  CHARACTER(LEN=*), PARAMETER :: routinen = 'setup_nhc_thermostat'
510 
511  INTEGER :: handle, nkind, number, region
512  LOGICAL :: do_shell
513  TYPE(map_info_type), POINTER :: map_info
514 
515  CALL timeset(routinen, handle)
516 
517  do_shell = .false.
518  IF (PRESENT(shell)) do_shell = shell
519  map_info => nhc%map_info
520 
521  nkind = SIZE(molecule_kind_set)
522  sum_of_thermostats = thermostat_info%sum_of_thermostats
523  map_info%dis_type = thermostat_info%dis_type
524  number = thermostat_info%number_of_thermostats
525  region = nhc%region
526 
527  CALL thermostat_mapping_region(map_info, deg_of_freedom, massive_atom_list, &
528  molecule_kind_set, local_molecules, molecule_set, para_env, natoms_local, &
529  simpar, number, region, gci, do_shell, thermostat_info%map_loc_thermo_gen, &
530  sum_of_thermostats)
531 
532  ALLOCATE (nhc%nvt(nhc%nhc_len, number))
533 
534  ! Now that we know how many there are stick this into nhc%nkt
535  ! (number of degrees of freedom times k_B T for the first thermostat
536  ! on the chain)
537  nhc%loc_num_nhc = number
538  nhc%glob_num_nhc = sum_of_thermostats
539 
540  CALL timestop(handle)
541 
542  END SUBROUTINE setup_nhc_thermostat
543 
544 ! **************************************************************************************************
545 !> \brief ...
546 !> \param thermostat_info ...
547 !> \param simpar ...
548 !> \param local_molecules ...
549 !> \param molecule_set ...
550 !> \param molecule_kind_set ...
551 !> \param nhc ...
552 !> \param para_env ...
553 !> \param gci ...
554 ! **************************************************************************************************
555  SUBROUTINE nhc_to_shell_mapping(thermostat_info, simpar, local_molecules, &
556  molecule_set, molecule_kind_set, nhc, para_env, gci)
557 
558  TYPE(thermostat_info_type), POINTER :: thermostat_info
559  TYPE(simpar_type), POINTER :: simpar
560  TYPE(distribution_1d_type), POINTER :: local_molecules
561  TYPE(molecule_type), POINTER :: molecule_set(:)
562  TYPE(molecule_kind_type), POINTER :: molecule_kind_set(:)
563  TYPE(lnhc_parameters_type), POINTER :: nhc
564  TYPE(mp_para_env_type), POINTER :: para_env
565  TYPE(global_constraint_type), POINTER :: gci
566 
567  CHARACTER(LEN=*), PARAMETER :: routinen = 'nhc_to_shell_mapping'
568 
569  INTEGER :: handle, i, imap, j, nshell_local, &
570  sum_of_thermostats
571  INTEGER, DIMENSION(:), POINTER :: deg_of_freedom, massive_shell_list
572  TYPE(map_info_type), POINTER :: map_info
573 
574  CALL timeset(routinen, handle)
575 
576  NULLIFY (massive_shell_list, deg_of_freedom)
577 
578  SELECT CASE (simpar%ensemble)
579  CASE DEFAULT
580  cpabort('Unknown ensemble!')
583  cpabort('Never reach this point!')
586 
587  CALL setup_nhc_thermostat(nhc, thermostat_info, deg_of_freedom, massive_shell_list, &
588  molecule_kind_set, local_molecules, molecule_set, para_env, nshell_local, &
589  simpar, sum_of_thermostats, gci, shell=.true.)
590 
591  map_info => nhc%map_info
592  ! Sum up the number of degrees of freedom on each thermostat.
593  ! first: initialize the target, via p_kin init s_kin
594  map_info%s_kin = 0.0_dp
595  DO j = 1, nshell_local
596  DO i = 1, 3
597  map_info%p_kin(i, j)%point = map_info%p_kin(i, j)%point + 1
598  END DO
599  END DO
600 
601  ! If thermostats are replicated but molecules distributed, we have to
602  ! sum s_kin over all processors
603  IF (map_info%dis_type == do_thermo_communication) CALL para_env%sum(map_info%s_kin)
604 
605  ! Now that we know how many there are stick this into nhc%nkt
606  ! (number of degrees of freedom times k_B T )
607  DO i = 1, nhc%loc_num_nhc
608  imap = map_info%map_index(i)
609  nhc%nvt(1, i)%nkt = simpar%temp_sh_ext*map_info%s_kin(imap)
610  nhc%nvt(1, i)%degrees_of_freedom = int(map_info%s_kin(imap))
611  END DO
612 
613  ! Getting the number of degrees of freedom times k_B T for the rest of the chain
614  DO i = 2, nhc%nhc_len
615  nhc%nvt(i, :)%nkt = simpar%temp_sh_ext
616  nhc%nvt(i, :)%degrees_of_freedom = 1
617  END DO
618  DEALLOCATE (deg_of_freedom)
619  DEALLOCATE (massive_shell_list)
620 
621  ! Let's clean the arrays
622  map_info%s_kin = 0.0_dp
623  map_info%v_scale = 0.0_dp
624  END SELECT
625 
626  CALL timestop(handle)
627 
628  END SUBROUTINE nhc_to_shell_mapping
629 
630 END MODULE extended_system_mapping
static GRID_HOST_DEVICE double fac(const int i)
Factorial function, e.g. fac(5) = 5! = 120.
Definition: grid_common.h:48
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
subroutine, public nhc_to_particle_mapping_fast(thermostat_info, simpar, local_molecules, molecule_set, molecule_kind_set, nhc, para_env, gci)
Creates the thermostatting maps.
subroutine, public nhc_to_particle_mapping_slow(thermostat_info, simpar, local_molecules, molecule_set, molecule_kind_set, nhc, para_env, gci)
Creates the thermostatting maps.
subroutine, public nhc_to_barostat_mapping(simpar, nhc)
Creates the thermostatting for the barostat.
subroutine, public nhc_to_shell_mapping(thermostat_info, simpar, local_molecules, molecule_set, molecule_kind_set, nhc, para_env, gci)
...
subroutine, public nhc_to_particle_mapping(thermostat_info, simpar, local_molecules, molecule_set, molecule_kind_set, nhc, para_env, gci)
Creates the thermostatting maps.
Lumps all possible extended system variables into one type for easy access and passing.
logical, parameter, public debug_isotropic_limit
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public do_thermo_no_communication
integer, parameter, public do_thermo_only_master
integer, parameter, public nvt_adiabatic_ensemble
integer, parameter, public nph_uniaxial_ensemble
integer, parameter, public npt_i_ensemble
integer, parameter, public isokin_ensemble
integer, parameter, public nph_uniaxial_damped_ensemble
integer, parameter, public npe_f_ensemble
integer, parameter, public langevin_ensemble
integer, parameter, public npe_i_ensemble
integer, parameter, public npt_ia_ensemble
integer, parameter, public nve_ensemble
integer, parameter, public npt_f_ensemble
integer, parameter, public reftraj_ensemble
integer, parameter, public nvt_ensemble
integer, parameter, public do_thermo_communication
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.
Type for storing MD parameters.
Definition: simpar_types.F:14
subroutine, public adiabatic_mapping_region(map_info, deg_of_freedom, massive_atom_list, molecule_kind_set, local_molecules, molecule_set, para_env, natoms_local, simpar, number, region, gci, shell, map_loc_thermo_gen, sum_of_thermostats)
Main general setup for adiabatic thermostat regions (Nose only)
subroutine, public thermostat_mapping_region(map_info, deg_of_freedom, massive_atom_list, molecule_kind_set, local_molecules, molecule_set, para_env, natoms_local, simpar, number, region, gci, shell, map_loc_thermo_gen, sum_of_thermostats)
Main general setup thermostat regions (thermostat independent)
subroutine, public init_baro_map_info(map_info, ndeg, num_thermo)
Initialize the map_info for barostat thermostat.
Thermostat structure: module containing thermostat available for MD.