(git:b279b6b)
input_restart_force_eval.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 
9 
10  USE atomic_kind_list_types, ONLY: atomic_kind_list_type
12  USE cell_types, ONLY: cell_type,&
14  USE cp_control_types, ONLY: dft_control_type
15  USE cp_files, ONLY: close_file,&
16  open_file
19  cp_sll_val_type
20  USE cp_subsys_types, ONLY: cp_subsys_get,&
21  cp_subsys_type
22  USE cp_units, ONLY: cp_unit_from_cp2k
23  USE distribution_1d_types, ONLY: distribution_1d_type
25  USE force_env_types, ONLY: force_env_get,&
26  force_env_type,&
30  use_qmmm,&
32  USE input_constants, ONLY: do_coord_cp2k,&
33  ehrenfest,&
34  mol_dyn_run,&
35  mon_car_run,&
36  pint_run,&
40  USE input_section_types, ONLY: &
44  USE input_val_types, ONLY: val_create,&
45  val_release,&
46  val_type
47  USE kinds, ONLY: default_string_length,&
48  dp
49  USE message_passing, ONLY: mp_para_env_type
51  USE molecule_list_types, ONLY: molecule_list_type
52  USE molecule_types, ONLY: get_molecule,&
53  molecule_type
54  USE multipole_types, ONLY: multipole_type
57  USE particle_list_types, ONLY: particle_list_type
61  USE virial_types, ONLY: virial_type
62 #include "./base/base_uses.f90"
63 
64  IMPLICIT NONE
65 
66  PRIVATE
67 
68  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'input_restart_force_eval'
69 
71 
72 CONTAINS
73 
74 ! **************************************************************************************************
75 !> \brief Updates the force_eval section of the input file
76 !> \param force_env ...
77 !> \param root_section ...
78 !> \param write_binary_restart_file ...
79 !> \param respa ...
80 !> \par History
81 !> 01.2006 created [teo]
82 !> \author Teodoro Laino
83 ! **************************************************************************************************
84  SUBROUTINE update_force_eval(force_env, root_section, &
85  write_binary_restart_file, respa)
86 
87  TYPE(force_env_type), POINTER :: force_env
88  TYPE(section_vals_type), POINTER :: root_section
89  LOGICAL, INTENT(IN) :: write_binary_restart_file
90  LOGICAL, OPTIONAL :: respa
91 
92  INTEGER :: i, i_subsys, iforce_eval, myid, &
93  nforce_eval
94  INTEGER, DIMENSION(:), POINTER :: i_force_eval
95  LOGICAL :: is_present, multiple_subsys, &
96  skip_vel_section
97  REAL(kind=dp), DIMENSION(:), POINTER :: work
98  TYPE(cell_type), POINTER :: cell
99  TYPE(cp_subsys_type), POINTER :: subsys
100  TYPE(dft_control_type), POINTER :: dft_control
101  TYPE(section_vals_type), POINTER :: cell_section, dft_section, &
102  force_env_sections, qmmm_section, &
103  rng_section, subsys_section, &
104  tmp_section
105  TYPE(virial_type), POINTER :: virial
106 
107  NULLIFY (rng_section, subsys_section, cell_section, virial, subsys, cell, dft_control, work, tmp_section)
108  ! If it's not a dynamical run don't update the velocity section
109  CALL section_vals_val_get(root_section, "GLOBAL%RUN_TYPE", i_val=myid)
110  skip_vel_section = ((myid /= mol_dyn_run) .AND. &
111  (myid /= mon_car_run) .AND. &
112  (myid /= pint_run) .AND. &
113  (myid /= ehrenfest))
114 
115  ! Go on updatig the force_env_section
116  force_env_sections => section_vals_get_subs_vals(root_section, "FORCE_EVAL")
117  CALL multiple_fe_list(force_env_sections, root_section, i_force_eval, nforce_eval)
118  ! The update of the input MUST be realized only on the main force_eval
119  ! All the others will be left not updated because there is no real need to update them...
120  iforce_eval = 1
121  subsys_section => section_vals_get_subs_vals3(force_env_sections, "SUBSYS", &
122  i_rep_section=i_force_eval(iforce_eval))
123  CALL update_subsys(subsys_section, force_env, skip_vel_section, &
124  write_binary_restart_file)
125 
126  ! For RESPA we need to update all subsystems
127  IF (PRESENT(respa)) THEN
128  IF (respa) THEN
129  DO i_subsys = 1, 2
130  iforce_eval = i_subsys
131  subsys_section => section_vals_get_subs_vals3(force_env_sections, "SUBSYS", &
132  i_rep_section=i_force_eval(iforce_eval))
133  CALL update_subsys(subsys_section, force_env, skip_vel_section, &
134  write_binary_restart_file)
135  END DO
136  END IF
137  END IF
138 
139  rng_section => section_vals_get_subs_vals(subsys_section, "RNG_INIT")
140  CALL update_rng_particle(rng_section, force_env)
141 
142  qmmm_section => section_vals_get_subs_vals3(force_env_sections, "QMMM", &
143  i_rep_section=i_force_eval(iforce_eval))
144  CALL update_qmmm(qmmm_section, force_env)
145 
146  ! Non-mixed CDFT calculation: update constraint Lagrangian and counter
147  IF (nforce_eval == 1) THEN
148  IF (ASSOCIATED(force_env%qs_env)) THEN
149  CALL get_qs_env(force_env%qs_env, dft_control=dft_control)
150  IF (dft_control%qs_control%cdft) THEN
151  dft_section => section_vals_get_subs_vals3(force_env_sections, "DFT", &
152  i_rep_section=i_force_eval(iforce_eval))
153  CALL section_vals_val_set(dft_section, "QS%CDFT%COUNTER", &
154  i_val=dft_control%qs_control%cdft_control%ienergy)
155  ALLOCATE (work(SIZE(dft_control%qs_control%cdft_control%strength)))
156  work = dft_control%qs_control%cdft_control%strength
157  CALL section_vals_val_set(dft_section, "QS%CDFT%STRENGTH", &
158  r_vals_ptr=work)
159  END IF
160  END IF
161  END IF
162  ! And now update only the cells of all other force_eval
163  ! This is to make consistent for cell variable runs..
164  IF (nforce_eval > 1) THEN
165  CALL force_env_get(force_env, subsys=subsys, cell=cell)
166  CALL cp_subsys_get(subsys, virial=virial)
167  CALL section_vals_val_get(root_section, "MULTIPLE_FORCE_EVALS%MULTIPLE_SUBSYS", l_val=multiple_subsys)
168  IF (virial%pv_availability .AND. multiple_subsys) THEN
169  DO iforce_eval = 2, nforce_eval
170  subsys_section => section_vals_get_subs_vals3(force_env_sections, "SUBSYS", &
171  i_rep_section=i_force_eval(iforce_eval))
172  cell_section => section_vals_get_subs_vals(subsys_section, "CELL")
173  CALL update_cell_section(cell, cell_section)
174  END DO
175  END IF
176  ! With mixed CDFT, update value of constraint Lagrangian
177  IF (ASSOCIATED(force_env%mixed_env)) THEN
178  IF (force_env%mixed_env%do_mixed_cdft) THEN
179  DO iforce_eval = 2, nforce_eval
180  dft_section => section_vals_get_subs_vals3(force_env_sections, "DFT", &
181  i_rep_section=i_force_eval(iforce_eval))
182  ALLOCATE (work(SIZE(force_env%mixed_env%strength, 2)))
183  work = force_env%mixed_env%strength(iforce_eval - 1, :)
184  CALL section_vals_val_set(dft_section, "QS%CDFT%STRENGTH", &
185  r_vals_ptr=work)
186  CALL section_vals_val_set(dft_section, "QS%CDFT%COUNTER", &
187  i_val=force_env%mixed_env%cdft_control%sim_step)
188  END DO
189  END IF
190  END IF
191  END IF
192 
193  IF (ASSOCIATED(force_env%qs_env)) THEN
194  CALL get_qs_env(force_env%qs_env, dft_control=dft_control)
195  dft_section => section_vals_get_subs_vals(force_env_sections, "DFT", &
196  i_rep_section=i_force_eval(iforce_eval))
197  tmp_section => section_vals_get_subs_vals(dft_section, "REAL_TIME_PROPAGATION")
198  CALL section_vals_get(tmp_section, explicit=is_present)
199  IF (is_present) THEN
200  CALL section_vals_val_set(root_section, "FORCE_EVAL%DFT%REAL_TIME_PROPAGATION%INITIAL_WFN", &
201  i_val=use_rt_restart)
202  IF (ASSOCIATED(dft_control%efield_fields)) THEN
203  tmp_section => section_vals_get_subs_vals(dft_section, "EFIELD")
204  ALLOCATE (work(3))
205  DO i = 1, SIZE(dft_control%efield_fields)
206  work = dft_control%efield_fields(1)%efield%vec_pot_initial
207  CALL section_vals_val_set(tmp_section, "VEC_POT_INITIAL", i_rep_section=i, &
208  r_vals_ptr=work)
209  END DO
210  END IF
211  END IF
212  END IF
213  DEALLOCATE (i_force_eval)
214 
215  END SUBROUTINE update_force_eval
216 
217 ! **************************************************************************************************
218 !> \brief Updates the qmmm section if needed
219 !> \param qmmm_section ...
220 !> \param force_env ...
221 !> \par History
222 !> 08.2007 created [teo]
223 !> \author Teodoro Laino
224 ! **************************************************************************************************
225  SUBROUTINE update_qmmm(qmmm_section, force_env)
226 
227  TYPE(section_vals_type), POINTER :: qmmm_section
228  TYPE(force_env_type), POINTER :: force_env
229 
230  LOGICAL :: explicit
231  REAL(kind=dp), DIMENSION(:), POINTER :: work
232 
233  SELECT CASE (force_env%in_use)
234  CASE (use_qmmm)
235  CALL section_vals_get(qmmm_section, explicit=explicit)
236  cpassert(explicit)
237 
238  ALLOCATE (work(3))
239  work = force_env%qmmm_env%qm%transl_v
240  CALL section_vals_val_set(qmmm_section, "INITIAL_TRANSLATION_VECTOR", r_vals_ptr=work)
241  END SELECT
242 
243  END SUBROUTINE update_qmmm
244 
245 ! **************************************************************************************************
246 !> \brief Updates the rng section of the input file
247 !> Write current status of the parallel random number generator (RNG)
248 !> \param rng_section ...
249 !> \param force_env ...
250 !> \par History
251 !> 01.2006 created [teo]
252 !> \author Teodoro Laino
253 ! **************************************************************************************************
254  SUBROUTINE update_rng_particle(rng_section, force_env)
255 
256  TYPE(section_vals_type), POINTER :: rng_section
257  TYPE(force_env_type), POINTER :: force_env
258 
259  CHARACTER(LEN=rng_record_length) :: rng_record
260  INTEGER :: iparticle, iparticle_kind, &
261  iparticle_local, nparticle, &
262  nparticle_kind, nparticle_local
263  INTEGER, ALLOCATABLE, DIMENSION(:, :) :: ascii
264  TYPE(atomic_kind_list_type), POINTER :: atomic_kinds
265  TYPE(cp_subsys_type), POINTER :: subsys
266  TYPE(distribution_1d_type), POINTER :: local_particles
267  TYPE(mp_para_env_type), POINTER :: para_env
268  TYPE(particle_list_type), POINTER :: particles
269 
270  CALL force_env_get(force_env, subsys=subsys, para_env=para_env)
271  CALL cp_subsys_get(subsys, atomic_kinds=atomic_kinds, local_particles=local_particles, &
272  particles=particles)
273 
274  IF (ASSOCIATED(local_particles%local_particle_set)) THEN
275  nparticle_kind = atomic_kinds%n_els
276  nparticle = particles%n_els
277 
278  ALLOCATE (ascii(rng_record_length, nparticle))
279  ascii = 0
280 
281  DO iparticle = 1, nparticle
282  DO iparticle_kind = 1, nparticle_kind
283  nparticle_local = local_particles%n_el(iparticle_kind)
284  DO iparticle_local = 1, nparticle_local
285  IF (iparticle == local_particles%list(iparticle_kind)%array(iparticle_local)) THEN
286  CALL local_particles%local_particle_set(iparticle_kind)% &
287  rng(iparticle_local)%stream%dump(rng_record=rng_record)
288  CALL string_to_ascii(rng_record, ascii(:, iparticle))
289  END IF
290  END DO
291  END DO
292  END DO
293 
294  CALL para_env%sum(ascii)
295 
296  CALL section_rng_val_set(rng_section=rng_section, nsize=nparticle, ascii=ascii)
297 
298  DEALLOCATE (ascii)
299  END IF
300 
301  END SUBROUTINE update_rng_particle
302 
303 ! **************************************************************************************************
304 !> \brief Updates the subsys section of the input file
305 !> \param subsys_section ...
306 !> \param force_env ...
307 !> \param skip_vel_section ...
308 !> \param write_binary_restart_file ...
309 !> \par History
310 !> 01.2006 created [teo]
311 !> \author Teodoro Laino
312 ! **************************************************************************************************
313  SUBROUTINE update_subsys(subsys_section, force_env, skip_vel_section, &
314  write_binary_restart_file)
315 
316  TYPE(section_vals_type), POINTER :: subsys_section
317  TYPE(force_env_type), POINTER :: force_env
318  LOGICAL, INTENT(IN) :: skip_vel_section, &
319  write_binary_restart_file
320 
321  CHARACTER(LEN=*), PARAMETER :: routinen = 'update_subsys'
322 
323  CHARACTER(LEN=default_string_length) :: coord_file_name, unit_str
324  INTEGER :: coord_file_format, handle, output_unit
325  INTEGER, DIMENSION(:), POINTER :: multiple_unit_cell
326  LOGICAL :: scale, use_ref_cell
327  REAL(kind=dp) :: conv_factor
328  TYPE(cell_type), POINTER :: cell
329  TYPE(cp_subsys_type), POINTER :: subsys
330  TYPE(molecule_list_type), POINTER :: molecules
331  TYPE(mp_para_env_type), POINTER :: para_env
332  TYPE(multipole_type), POINTER :: multipoles
333  TYPE(particle_list_type), POINTER :: core_particles, particles, &
334  shell_particles
335  TYPE(section_vals_type), POINTER :: work_section
336 
337  NULLIFY (work_section, core_particles, particles, shell_particles, &
338  subsys, cell, para_env, multipoles)
339  CALL timeset(routinen, handle)
340  CALL force_env_get(force_env, subsys=subsys, cell=cell, para_env=para_env)
341 
342  CALL cp_subsys_get(subsys, particles=particles, molecules=molecules, &
343  shell_particles=shell_particles, core_particles=core_particles, &
344  multipoles=multipoles)
345 
346  ! Remove the multiple_unit_cell from the input structure.. at this point we have
347  ! already all the information we need..
348  ALLOCATE (multiple_unit_cell(3))
349  multiple_unit_cell = 1
350  CALL section_vals_val_set(subsys_section, "TOPOLOGY%MULTIPLE_UNIT_CELL", &
351  i_vals_ptr=multiple_unit_cell)
352 
353  ! Coordinates and Velocities
354  work_section => section_vals_get_subs_vals(subsys_section, "COORD")
355  CALL section_vals_val_get(work_section, "UNIT", c_val=unit_str)
356  CALL section_vals_val_get(work_section, "SCALED", l_val=scale)
357  conv_factor = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
358  IF (.NOT. write_binary_restart_file) THEN
359  CALL section_vals_val_get(subsys_section, "TOPOLOGY%COORD_FILE_FORMAT", &
360  i_val=coord_file_format)
361  IF (coord_file_format == do_coord_cp2k) THEN
362  CALL section_vals_val_get(subsys_section, "TOPOLOGY%COORD_FILE_NAME", &
363  c_val=coord_file_name)
364  output_unit = 0
365  IF (para_env%is_source()) THEN
366  CALL open_file(file_name=trim(adjustl(coord_file_name)), &
367  file_action="READWRITE", &
368  file_form="FORMATTED", &
369  file_position="REWIND", &
370  file_status="UNKNOWN", &
371  unit_number=output_unit)
372  CALL dump_coordinates_cp2k(particles, molecules, cell, conv_factor, &
373  output_unit=output_unit, &
374  core_or_shell=.false., &
375  scaled_coordinates=scale)
376  CALL close_file(unit_number=output_unit)
377  END IF
378  ELSE
379  CALL section_coord_val_set(work_section, particles, molecules, conv_factor, scale, &
380  cell)
381  END IF
382  END IF
383  CALL section_vals_val_set(subsys_section, "TOPOLOGY%NUMBER_OF_ATOMS", i_val=particles%n_els)
384  work_section => section_vals_get_subs_vals(subsys_section, "VELOCITY")
385  IF (.NOT. skip_vel_section) THEN
386  IF (.NOT. write_binary_restart_file) THEN
387  CALL section_velocity_val_set(work_section, particles, conv_factor=1.0_dp)
388  END IF
389  ELSE
390  CALL section_vals_remove_values(work_section)
391  END IF
392 
393  ! Write restart input for core-shell model
394  IF (.NOT. write_binary_restart_file) THEN
395  IF (ASSOCIATED(shell_particles)) THEN
396  work_section => section_vals_get_subs_vals(subsys_section, "SHELL_COORD")
397  CALL section_vals_val_get(work_section, "UNIT", c_val=unit_str)
398  CALL section_vals_val_get(work_section, "SCALED", l_val=scale)
399  conv_factor = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
400  CALL section_coord_val_set(work_section, shell_particles, molecules, &
401  conv_factor, scale, cell, shell=.true.)
402  IF (.NOT. skip_vel_section) THEN
403  work_section => section_vals_get_subs_vals(subsys_section, "SHELL_VELOCITY")
404  CALL section_velocity_val_set(work_section, shell_particles, conv_factor=1.0_dp)
405  END IF
406  END IF
407  IF (ASSOCIATED(core_particles)) THEN
408  work_section => section_vals_get_subs_vals(subsys_section, "CORE_COORD")
409  CALL section_vals_val_get(work_section, "UNIT", c_val=unit_str)
410  CALL section_vals_val_get(work_section, "SCALED", l_val=scale)
411  conv_factor = cp_unit_from_cp2k(1.0_dp, trim(unit_str))
412  CALL section_coord_val_set(work_section, core_particles, molecules, &
413  conv_factor, scale, cell, shell=.true.)
414  IF (.NOT. skip_vel_section) THEN
415  work_section => section_vals_get_subs_vals(subsys_section, "CORE_VELOCITY")
416  CALL section_velocity_val_set(work_section, core_particles, conv_factor=1.0_dp)
417  END IF
418  END IF
419  END IF
420 
421  ! Updating cell info
422  CALL force_env_get(force_env, cell=cell)
423  work_section => section_vals_get_subs_vals(subsys_section, "CELL")
424  CALL update_cell_section(cell, cell_section=work_section)
425 
426  ! Updating cell_ref info
427  use_ref_cell = .false.
428  SELECT CASE (force_env%in_use)
429  CASE (use_qs_force)
430  CALL get_qs_env(force_env%qs_env, cell_ref=cell, use_ref_cell=use_ref_cell)
431  CASE (use_eip_force)
432  CALL eip_env_get(force_env%eip_env, cell_ref=cell, use_ref_cell=use_ref_cell)
433  CASE (use_nnp_force)
434  CALL nnp_env_get(force_env%nnp_env, cell_ref=cell, use_ref_cell=use_ref_cell)
435  END SELECT
436  IF (use_ref_cell) THEN
437  work_section => section_vals_get_subs_vals(subsys_section, "CELL%CELL_REF")
438  CALL update_cell_section(cell, cell_section=work_section)
439  END IF
440 
441  ! Updating multipoles
442  IF (ASSOCIATED(multipoles)) THEN
443  work_section => section_vals_get_subs_vals(subsys_section, "MULTIPOLES")
444  DO
445  IF (SIZE(work_section%values, 2) == 1) EXIT
446  CALL section_vals_add_values(work_section)
447  END DO
448  IF (ASSOCIATED(multipoles%dipoles)) THEN
449  work_section => section_vals_get_subs_vals(subsys_section, "MULTIPOLES%DIPOLES")
450  CALL update_dipoles_section(multipoles%dipoles, work_section)
451  END IF
452  IF (ASSOCIATED(multipoles%quadrupoles)) THEN
453  work_section => section_vals_get_subs_vals(subsys_section, "MULTIPOLES%QUADRUPOLES")
454  CALL update_quadrupoles_section(multipoles%quadrupoles, work_section)
455  END IF
456  END IF
457 
458  CALL timestop(handle)
459 
460  END SUBROUTINE update_subsys
461 
462 ! **************************************************************************************************
463 !> \brief Routine to update a cell section
464 !> \param cell ...
465 !> \param cell_section ...
466 !> \author Ole Schuett
467 ! **************************************************************************************************
468  SUBROUTINE update_cell_section(cell, cell_section)
469 
470  TYPE(cell_type), POINTER :: cell
471  TYPE(section_vals_type), POINTER :: cell_section
472 
473  CHARACTER(LEN=*), PARAMETER :: routinen = 'update_cell_section'
474 
475  INTEGER :: handle
476  INTEGER, DIMENSION(:), POINTER :: iwork
477  REAL(kind=dp), DIMENSION(:), POINTER :: work
478 
479  NULLIFY (work, iwork)
480  CALL timeset(routinen, handle)
481 
482  ! CELL VECTORS - A
483  ALLOCATE (work(3))
484  work(1:3) = cell%hmat(1:3, 1)
485  CALL section_vals_val_set(cell_section, "A", r_vals_ptr=work)
486 
487  ! CELL VECTORS - B
488  ALLOCATE (work(3))
489  work(1:3) = cell%hmat(1:3, 2)
490  CALL section_vals_val_set(cell_section, "B", r_vals_ptr=work)
491 
492  ! CELL VECTORS - C
493  ALLOCATE (work(3))
494  work(1:3) = cell%hmat(1:3, 3)
495  CALL section_vals_val_set(cell_section, "C", r_vals_ptr=work)
496 
497  ! MULTIPLE_UNIT_CELL
498  ALLOCATE (iwork(3))
499  iwork = 1
500  CALL section_vals_val_set(cell_section, "MULTIPLE_UNIT_CELL", i_vals_ptr=iwork)
501 
502  ! Unset unused or misleading information
503  CALL section_vals_val_unset(cell_section, "ABC")
504  CALL section_vals_val_unset(cell_section, "ALPHA_BETA_GAMMA")
505 
506  CALL timestop(handle)
507 
508  END SUBROUTINE update_cell_section
509 
510 ! **************************************************************************************************
511 !> \brief routine to dump coordinates.. fast implementation
512 !> \param coord_section ...
513 !> \param particles ...
514 !> \param molecules ...
515 !> \param conv_factor ...
516 !> \param scale ...
517 !> \param cell ...
518 !> \param shell ...
519 !> \par History
520 !> 02.2006 created [teo]
521 !> \author Teodoro Laino
522 ! **************************************************************************************************
523  SUBROUTINE section_coord_val_set(coord_section, particles, molecules, conv_factor, &
524  scale, cell, shell)
525 
526  TYPE(section_vals_type), POINTER :: coord_section
527  TYPE(particle_list_type), POINTER :: particles
528  TYPE(molecule_list_type), POINTER :: molecules
529  REAL(kind=dp) :: conv_factor
530  LOGICAL, INTENT(IN) :: scale
531  TYPE(cell_type), POINTER :: cell
532  LOGICAL, INTENT(IN), OPTIONAL :: shell
533 
534  CHARACTER(LEN=*), PARAMETER :: routinen = 'section_coord_val_set'
535 
536  CHARACTER(LEN=2*default_string_length) :: line
537  CHARACTER(LEN=default_string_length) :: my_tag, name
538  INTEGER :: handle, ik, imol, irk, last_atom, nlist
539  LOGICAL :: ldum, molname_generated, my_shell
540  REAL(kind=dp), DIMENSION(3) :: r, s
541  TYPE(cp_sll_val_type), POINTER :: new_pos, vals
542  TYPE(molecule_type), POINTER :: molecule_now
543  TYPE(section_type), POINTER :: section
544  TYPE(val_type), POINTER :: my_val, old_val
545 
546  CALL timeset(routinen, handle)
547 
548  NULLIFY (my_val, old_val, section, vals)
549  my_shell = .false.
550  IF (PRESENT(shell)) my_shell = shell
551  cpassert(ASSOCIATED(coord_section))
552  cpassert(coord_section%ref_count > 0)
553  section => coord_section%section
554  ik = section_get_keyword_index(section, "_DEFAULT_KEYWORD_")
555  IF (ik == -2) &
556  CALL cp_abort(__location__, &
557  "section "//trim(section%name)//" does not contain keyword "// &
558  "_DEFAULT_KEYWORD_")
559 
560  DO
561  IF (SIZE(coord_section%values, 2) == 1) EXIT
562  CALL section_vals_add_values(coord_section)
563  END DO
564  vals => coord_section%values(ik, 1)%list
565  nlist = 0
566  IF (ASSOCIATED(vals)) THEN
567  nlist = cp_sll_val_get_length(vals)
568  END IF
569  imol = 0
570  last_atom = 0
571  DO irk = 1, particles%n_els
572  CALL get_atomic_kind(particles%els(irk)%atomic_kind, name=name)
573  IF (my_shell) THEN
574  s = particles%els(irk)%r
575  IF (scale) THEN
576  r = s
577  CALL real_to_scaled(s, r, cell)
578  ELSE
579  s = s*conv_factor
580  END IF
581  WRITE (unit=line, fmt="(T7,A,3ES25.16E3,1X,I0)") &
582  trim(adjustl(name)), s(1:3), particles%els(irk)%atom_index
583  CALL val_create(my_val, lc_val=line)
584  IF (nlist /= 0) THEN
585  IF (irk == 1) THEN
586  new_pos => vals
587  ELSE
588  new_pos => new_pos%rest
589  END IF
590  old_val => new_pos%first_el
591  CALL val_release(old_val)
592  new_pos%first_el => my_val
593  ELSE
594  IF (irk == 1) THEN
595  NULLIFY (new_pos)
596  CALL cp_sll_val_create(new_pos, first_el=my_val)
597  vals => new_pos
598  ELSE
599  NULLIFY (new_pos%rest)
600  CALL cp_sll_val_create(new_pos%rest, first_el=my_val)
601  new_pos => new_pos%rest
602  END IF
603  END IF
604  NULLIFY (my_val)
605  ELSE
606  IF (last_atom < irk) THEN
607  imol = imol + 1
608  molecule_now => molecules%els(imol)
609  CALL get_molecule(molecule_now, last_atom=last_atom)
610  CALL get_molecule_kind(molecule_now%molecule_kind, molname_generated=molname_generated, &
611  name=my_tag)
612  IF (molname_generated) my_tag = ""
613  END IF
614  ldum = qmmm_ff_precond_only_qm(my_tag)
615  ldum = qmmm_ff_precond_only_qm(name)
616  s = particles%els(irk)%r
617  IF (scale) THEN
618  r = s
619  CALL real_to_scaled(s, r, cell)
620  ELSE
621  s = s*conv_factor
622  END IF
623  IF (len_trim(my_tag) > 0) THEN
624  WRITE (unit=line, fmt="(T7,A,3ES25.16E3,1X,A,1X,I0)") &
625  trim(adjustl(name)), s(1:3), trim(my_tag), imol
626  ELSE
627  WRITE (unit=line, fmt="(T7,A,3ES25.16E3)") &
628  trim(adjustl(name)), s(1:3)
629  END IF
630  CALL val_create(my_val, lc_val=line)
631 
632  IF (nlist /= 0) THEN
633  IF (irk == 1) THEN
634  new_pos => vals
635  ELSE
636  new_pos => new_pos%rest
637  END IF
638  old_val => new_pos%first_el
639  CALL val_release(old_val)
640  new_pos%first_el => my_val
641  ELSE
642  IF (irk == 1) THEN
643  NULLIFY (new_pos)
644  CALL cp_sll_val_create(new_pos, first_el=my_val)
645  vals => new_pos
646  ELSE
647  NULLIFY (new_pos%rest)
648  CALL cp_sll_val_create(new_pos%rest, first_el=my_val)
649  new_pos => new_pos%rest
650  END IF
651  END IF
652  NULLIFY (my_val)
653  END IF
654  END DO
655 
656  coord_section%values(ik, 1)%list => vals
657 
658  CALL timestop(handle)
659 
660  END SUBROUTINE section_coord_val_set
661 
662 ! **************************************************************************************************
663 !> \brief routine to dump dipoles.. fast implementation
664 !> \param dipoles ...
665 !> \param dipoles_section ...
666 !> \par History
667 !> 12.2007 created [teo]
668 !> \author Teodoro Laino
669 ! **************************************************************************************************
670  SUBROUTINE update_dipoles_section(dipoles, dipoles_section)
671 
672  REAL(kind=dp), DIMENSION(:, :), POINTER :: dipoles
673  TYPE(section_vals_type), POINTER :: dipoles_section
674 
675  CHARACTER(LEN=*), PARAMETER :: routinen = 'update_dipoles_section'
676 
677  INTEGER :: handle, ik, irk, nlist, nloop
678  REAL(kind=dp), DIMENSION(:), POINTER :: work
679  TYPE(cp_sll_val_type), POINTER :: new_pos, vals
680  TYPE(section_type), POINTER :: section
681  TYPE(val_type), POINTER :: my_val, old_val
682 
683  CALL timeset(routinen, handle)
684  NULLIFY (my_val, old_val, section, vals)
685  cpassert(ASSOCIATED(dipoles_section))
686  cpassert(dipoles_section%ref_count > 0)
687  section => dipoles_section%section
688  ik = section_get_keyword_index(section, "_DEFAULT_KEYWORD_")
689  IF (ik == -2) &
690  CALL cp_abort(__location__, &
691  "section "//trim(section%name)//" does not contain keyword "// &
692  "_DEFAULT_KEYWORD_")
693 
694  ! At least one of the two arguments must be present..
695  nloop = SIZE(dipoles, 2)
696  DO
697  IF (SIZE(dipoles_section%values, 2) == 1) EXIT
698  CALL section_vals_add_values(dipoles_section)
699  END DO
700  vals => dipoles_section%values(ik, 1)%list
701  nlist = 0
702  IF (ASSOCIATED(vals)) THEN
703  nlist = cp_sll_val_get_length(vals)
704  END IF
705  DO irk = 1, nloop
706  ALLOCATE (work(3))
707  ! Always stored in A.U.
708  work = dipoles(1:3, irk)
709  CALL val_create(my_val, r_vals_ptr=work)
710 
711  IF (nlist /= 0) THEN
712  IF (irk == 1) THEN
713  new_pos => vals
714  ELSE
715  new_pos => new_pos%rest
716  END IF
717  old_val => new_pos%first_el
718  CALL val_release(old_val)
719  new_pos%first_el => my_val
720  ELSE
721  IF (irk == 1) THEN
722  NULLIFY (new_pos)
723  CALL cp_sll_val_create(new_pos, first_el=my_val)
724  vals => new_pos
725  ELSE
726  NULLIFY (new_pos%rest)
727  CALL cp_sll_val_create(new_pos%rest, first_el=my_val)
728  new_pos => new_pos%rest
729  END IF
730  END IF
731  NULLIFY (my_val)
732  END DO
733  dipoles_section%values(ik, 1)%list => vals
734 
735  CALL timestop(handle)
736 
737  END SUBROUTINE update_dipoles_section
738 
739 ! **************************************************************************************************
740 !> \brief routine to dump quadrupoles.. fast implementation
741 !> \param quadrupoles ...
742 !> \param quadrupoles_section ...
743 !> \par History
744 !> 12.2007 created [teo]
745 !> \author Teodoro Laino
746 ! **************************************************************************************************
747  SUBROUTINE update_quadrupoles_section(quadrupoles, quadrupoles_section)
748 
749  REAL(kind=dp), DIMENSION(:, :, :), POINTER :: quadrupoles
750  TYPE(section_vals_type), POINTER :: quadrupoles_section
751 
752  CHARACTER(LEN=*), PARAMETER :: routinen = 'update_quadrupoles_section'
753 
754  INTEGER :: handle, i, ik, ind, irk, j, nlist, nloop
755  REAL(kind=dp), DIMENSION(:), POINTER :: work
756  TYPE(cp_sll_val_type), POINTER :: new_pos, vals
757  TYPE(section_type), POINTER :: section
758  TYPE(val_type), POINTER :: my_val, old_val
759 
760  CALL timeset(routinen, handle)
761  NULLIFY (my_val, old_val, section, vals)
762  cpassert(ASSOCIATED(quadrupoles_section))
763  cpassert(quadrupoles_section%ref_count > 0)
764  section => quadrupoles_section%section
765  ik = section_get_keyword_index(section, "_DEFAULT_KEYWORD_")
766  IF (ik == -2) &
767  CALL cp_abort(__location__, &
768  "section "//trim(section%name)//" does not contain keyword "// &
769  "_DEFAULT_KEYWORD_")
770 
771  ! At least one of the two arguments must be present..
772  nloop = SIZE(quadrupoles, 2)
773  DO
774  IF (SIZE(quadrupoles_section%values, 2) == 1) EXIT
775  CALL section_vals_add_values(quadrupoles_section)
776  END DO
777  vals => quadrupoles_section%values(ik, 1)%list
778  nlist = 0
779  IF (ASSOCIATED(vals)) THEN
780  nlist = cp_sll_val_get_length(vals)
781  END IF
782  DO irk = 1, nloop
783  ALLOCATE (work(6))
784  ! Always stored in A.U.
785  ind = 0
786  DO i = 1, 3
787  DO j = i, 3
788  ind = ind + 1
789  work(ind) = quadrupoles(j, i, irk)
790  END DO
791  END DO
792  CALL val_create(my_val, r_vals_ptr=work)
793 
794  IF (nlist /= 0) THEN
795  IF (irk == 1) THEN
796  new_pos => vals
797  ELSE
798  new_pos => new_pos%rest
799  END IF
800  old_val => new_pos%first_el
801  CALL val_release(old_val)
802  new_pos%first_el => my_val
803  ELSE
804  IF (irk == 1) THEN
805  NULLIFY (new_pos)
806  CALL cp_sll_val_create(new_pos, first_el=my_val)
807  vals => new_pos
808  ELSE
809  NULLIFY (new_pos%rest)
810  CALL cp_sll_val_create(new_pos%rest, first_el=my_val)
811  new_pos => new_pos%rest
812  END IF
813  END IF
814  NULLIFY (my_val)
815  END DO
816 
817  quadrupoles_section%values(ik, 1)%list => vals
818 
819  CALL timestop(handle)
820 
821  END SUBROUTINE update_quadrupoles_section
822 
823 ! **************************************************************************************************
824 !> \brief Dump atomic, core, or shell coordinates to a file in CP2K &COORD
825 !> section format
826 !> \param particles ...
827 !> \param molecules ...
828 !> \param cell ...
829 !> \param conv_factor ...
830 !> \param output_unit ...
831 !> \param core_or_shell ...
832 !> \param scaled_coordinates ...
833 !> \par History
834 !> 07.02.2011 (Creation, MK)
835 !> \author Matthias Krack (MK)
836 !> \version 1.0
837 ! **************************************************************************************************
838  SUBROUTINE dump_coordinates_cp2k(particles, molecules, cell, conv_factor, &
839  output_unit, core_or_shell, &
840  scaled_coordinates)
841 
842  TYPE(particle_list_type), POINTER :: particles
843  TYPE(molecule_list_type), POINTER :: molecules
844  TYPE(cell_type), POINTER :: cell
845  REAL(kind=dp), INTENT(IN) :: conv_factor
846  INTEGER, INTENT(IN) :: output_unit
847  LOGICAL, INTENT(IN) :: core_or_shell, scaled_coordinates
848 
849  CHARACTER(LEN=*), PARAMETER :: routinen = 'dump_coordinates_cp2k'
850 
851  CHARACTER(LEN=default_string_length) :: kind_name, tag
852  INTEGER :: handle, imolecule, iparticle, last_atom
853  LOGICAL :: ldum, molname_generated
854  REAL(kind=dp), DIMENSION(3) :: r, s
855  TYPE(molecule_type), POINTER :: molecule
856 
857  CALL timeset(routinen, handle)
858 
859  cpassert(ASSOCIATED(particles))
860  IF (.NOT. core_or_shell) THEN
861  cpassert(ASSOCIATED(molecules))
862  END IF
863  IF (scaled_coordinates) THEN
864  cpassert(ASSOCIATED(cell))
865  END IF
866 
867  kind_name = ""
868  tag = ""
869  imolecule = 0
870  last_atom = 0
871  DO iparticle = 1, particles%n_els
872  CALL get_atomic_kind(particles%els(iparticle)%atomic_kind, name=kind_name)
873  IF (.NOT. core_or_shell) THEN
874  IF (iparticle > last_atom) THEN
875  imolecule = imolecule + 1
876  molecule => molecules%els(imolecule)
877  CALL get_molecule(molecule, last_atom=last_atom)
878  CALL get_molecule_kind(molecule%molecule_kind, &
879  molname_generated=molname_generated, &
880  name=tag)
881  IF (molname_generated) tag = ""
882  END IF
883  ldum = qmmm_ff_precond_only_qm(tag)
884  ldum = qmmm_ff_precond_only_qm(kind_name)
885  END IF
886  IF (scaled_coordinates) THEN
887  CALL real_to_scaled(s, particles%els(iparticle)%r, cell)
888  r(1:3) = s(1:3)
889  ELSE
890  r(1:3) = particles%els(iparticle)%r(1:3)*conv_factor
891  END IF
892  IF (core_or_shell) THEN
893  WRITE (unit=output_unit, fmt="(A,3ES25.16E3,1X,I0)") &
894  trim(adjustl(kind_name)), r(1:3), particles%els(iparticle)%atom_index
895  ELSE
896  IF (len_trim(tag) > 0) THEN
897  WRITE (unit=output_unit, fmt="(A,3ES25.16E3,1X,A,1X,I0)") &
898  trim(adjustl(kind_name)), r(1:3), trim(tag), imolecule
899  ELSE
900  WRITE (unit=output_unit, fmt="(A,3ES25.16E3)") &
901  trim(adjustl(kind_name)), r(1:3)
902  END IF
903  END IF
904  END DO
905 
906  CALL timestop(handle)
907 
908  END SUBROUTINE dump_coordinates_cp2k
909 
910 END MODULE input_restart_force_eval
represent a simple array based list of the given type
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.
Handles all functions related to the CELL.
Definition: cell_types.F:15
subroutine, public real_to_scaled(s, r, cell)
Transform real to scaled cell coordinates. s=h_inv*r.
Definition: cell_types.F:486
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Utility routines to open and close files. Tracking of preconnections.
Definition: cp_files.F:16
subroutine, public open_file(file_name, file_status, file_form, file_action, file_position, file_pad, unit_number, debug, skip_get_unit_number, file_access)
Opens the requested file using a free unit number.
Definition: cp_files.F:308
subroutine, public close_file(unit_number, file_status, keep_preconnection)
Close an open file given by its logical unit number. Optionally, keep the file and unit preconnected.
Definition: cp_files.F:119
integer function, public cp_sll_val_get_length(sll)
returns the length of the list
subroutine, public cp_sll_val_create(sll, first_el, rest)
allocates and initializes a single linked list
types that represent a subsys, i.e. a part of the system
subroutine, public cp_subsys_get(subsys, ref_count, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell)
returns information about various attributes of the given subsys
unit conversion facility
Definition: cp_units.F:30
real(kind=dp) function, public cp_unit_from_cp2k(value, unit_str, defaults, power)
converts from the internal cp2k units to the given unit
Definition: cp_units.F:1179
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
The environment for the empirical interatomic potential methods.
subroutine, public eip_env_get(eip_env, eip_model, eip_energy, eip_energy_var, eip_forces, coord_avg, coord_var, count, subsys, atomic_kind_set, particle_set, local_particles, molecule_kind_set, molecule_set, local_molecules, eip_input, force_env_input, cell, cell_ref, use_ref_cell, eip_kinetic_energy, eip_potential_energy, virial)
Returns various attributes of the eip environment.
Interface for the force calculations.
integer, parameter, public use_qmmm
subroutine, public multiple_fe_list(force_env_sections, root_section, i_force_eval, nforce_eval)
returns the order of the multiple force_env
integer, parameter, public use_eip_force
recursive subroutine, public force_env_get(force_env, in_use, fist_env, qs_env, meta_env, fp_env, subsys, para_env, potential_energy, additional_potential, kinetic_energy, harmonic_shell, kinetic_shell, cell, sub_force_env, qmmm_env, qmmmx_env, eip_env, pwdft_env, globenv, input, force_env_section, method_name_id, root_section, mixed_env, nnp_env, embed_env)
returns various attributes about the force environment
integer, parameter, public use_qs_force
integer, parameter, public use_nnp_force
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public ehrenfest
integer, parameter, public mon_car_run
integer, parameter, public mol_dyn_run
integer, parameter, public use_rt_restart
integer, parameter, public pint_run
integer, parameter, public do_coord_cp2k
subroutine, public section_velocity_val_set(velocity_section, particles, velocity, conv_factor)
routine to dump velocities.. fast implementation
subroutine, public update_subsys(subsys_section, force_env, skip_vel_section, write_binary_restart_file)
Updates the subsys section of the input file.
subroutine, public update_force_eval(force_env, root_section, write_binary_restart_file, respa)
Updates the force_eval section of the input file.
subroutine, public section_rng_val_set(rng_section, nsize, ascii)
routine to dump rngs.. fast implementation
objects that represent the structure of input sections and the data contained in an input section
subroutine, public section_vals_val_unset(section_vals, keyword_name, i_rep_section, i_rep_val)
unsets (removes) the requested value (if it is a keyword repetitions removes the repetition,...
subroutine, public section_vals_val_set(section_vals, keyword_name, i_rep_section, i_rep_val, val, l_val, i_val, r_val, c_val, l_vals_ptr, i_vals_ptr, r_vals_ptr, c_vals_ptr)
sets the requested value
integer function, public section_get_keyword_index(section, keyword_name)
returns the index of the requested keyword (or -2 if not found)
subroutine, public section_vals_remove_values(section_vals)
removes the values of a repetition of the section
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
type(section_vals_type) function, pointer, public section_vals_get_subs_vals3(section_vals, subsection_name, i_rep_section)
returns the values of the n-th non default subsection (null if no such section exists (not so many no...
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
subroutine, public section_vals_add_values(section_vals)
adds the place to store the values of a repetition of the section
a wrapper for basic fortran types.
subroutine, public val_create(val, l_val, l_vals, l_vals_ptr, i_val, i_vals, i_vals_ptr, r_val, r_vals, r_vals_ptr, c_val, c_vals, c_vals_ptr, lc_val, lc_vals, lc_vals_ptr, enum)
creates a keyword value
subroutine, public val_release(val)
releases the given val
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public dp
Definition: kinds.F:34
integer, parameter, public default_string_length
Definition: kinds.F:57
Interface to the message passing library MPI.
Define the molecule kind structure types and the corresponding functionality.
subroutine, public get_molecule_kind(molecule_kind, atom_list, bond_list, bend_list, ub_list, impr_list, opbend_list, colv_list, fixd_list, g3x3_list, g4x6_list, vsite_list, torsion_list, shell_list, name, mass, charge, kind_number, natom, nbend, nbond, nub, nimpr, nopbend, nconstraint, nconstraint_fixd, nfixd, ncolv, ng3x3, ng4x6, nvsite, nfixd_restraint, ng3x3_restraint, ng4x6_restraint, nvsite_restraint, nrestraints, nmolecule, nsgf, nshell, ntorsion, molecule_list, nelectron, nelectron_alpha, nelectron_beta, bond_kind_set, bend_kind_set, ub_kind_set, impr_kind_set, opbend_kind_set, torsion_kind_set, molname_generated)
Get informations about a molecule kind.
represent a simple array based list of the given type
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.
Multipole structure: for multipole (fixed and induced) in FF based MD.
Data types for neural network potentials.
subroutine, public nnp_env_get(nnp_env, nnp_forces, subsys, atomic_kind_set, particle_set, local_particles, molecule_kind_set, molecule_set, local_molecules, nnp_input, force_env_input, cell, cell_ref, use_ref_cell, nnp_potential_energy, virial)
Returns various attributes of the nnp environment.
Parallel (pseudo)random number generator (RNG) for multiple streams and substreams of random numbers.
integer, parameter, public rng_record_length
represent a simple array based list of the given type
logical function, public qmmm_ff_precond_only_qm(id1, id2, id3, id4, is_link)
This function handles the atom names and modifies the "_QM_" prefix, in order to find the parameters ...
Definition: qmmm_ff_fist.F:39
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_RI_aux_kp, matrix_s, matrix_s_RI_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, ecoul_1c, rho0_s_rs, rho0_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, WannierCentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Get the QUICKSTEP environment.
Utilities for string manipulations.
subroutine, public string_to_ascii(string, nascii)
Convert a string to sequence of integer numbers.