Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Define the molecule kind structure types and the corresponding
10 : !> functionality
11 : !> \par History
12 : !> Teodoro Laino [tlaino] 12.2008 - Preparing for VIRTUAL SITE constraints
13 : !> (patch by Marcel Baer)
14 : !> \author Matthias Krack (22.08.2003)
15 : ! **************************************************************************************************
16 : MODULE molecule_kind_types
17 : USE atomic_kind_types, ONLY: atomic_kind_type,&
18 : get_atomic_kind
19 : USE cell_types, ONLY: periodicity_string,&
20 : use_perd_x,&
21 : use_perd_xy,&
22 : use_perd_xyz,&
23 : use_perd_xz,&
24 : use_perd_y,&
25 : use_perd_yz,&
26 : use_perd_z
27 : USE colvar_types, ONLY: &
28 : Wc_colvar_id, acid_hyd_dist_colvar_id, acid_hyd_shell_colvar_id, angle_colvar_id, &
29 : colvar_counters, combine_colvar_id, coord_colvar_id, dfunct_colvar_id, dist_colvar_id, &
30 : distance_from_path_colvar_id, gyration_colvar_id, hbp_colvar_id, hydronium_dist_colvar_id, &
31 : hydronium_shell_colvar_id, mindist_colvar_id, no_colvar_id, plane_distance_colvar_id, &
32 : plane_plane_angle_colvar_id, population_colvar_id, qparm_colvar_id, &
33 : reaction_path_colvar_id, ring_puckering_colvar_id, rmsd_colvar_id, rotation_colvar_id, &
34 : torsion_colvar_id, u_colvar_id, voronoiipz_colvar_id, xyz_diag_colvar_id, &
35 : xyz_outerdiag_colvar_id
36 : USE cp_log_handling, ONLY: cp_get_default_logger,&
37 : cp_logger_type
38 : USE cp_output_handling, ONLY: cp_print_key_finished_output,&
39 : cp_print_key_unit_nr
40 : USE cp_units, ONLY: cp_unit_from_cp2k
41 : USE force_field_kind_types, ONLY: &
42 : bend_kind_type, bond_kind_type, do_ff_undef, impr_kind_dealloc_ref, impr_kind_type, &
43 : opbend_kind_type, torsion_kind_dealloc_ref, torsion_kind_type, ub_kind_dealloc_ref, &
44 : ub_kind_type
45 : USE input_section_types, ONLY: section_vals_type
46 : USE kinds, ONLY: default_string_length,&
47 : dp
48 : USE shell_potential_types, ONLY: shell_kind_type
49 : #include "../base/base_uses.f90"
50 :
51 : IMPLICIT NONE
52 :
53 : PRIVATE
54 :
55 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'molecule_kind_types'
56 :
57 : ! Define the derived structure types
58 :
59 : TYPE atom_type
60 : TYPE(atomic_kind_type), POINTER :: atomic_kind => NULL()
61 : INTEGER :: id_name = 0
62 : END TYPE atom_type
63 :
64 : TYPE shell_type
65 : INTEGER :: a = 0
66 : CHARACTER(LEN=default_string_length) :: name = ""
67 : TYPE(shell_kind_type), POINTER :: shell_kind => NULL()
68 : END TYPE shell_type
69 :
70 : TYPE bond_type
71 : INTEGER :: a = 0, b = 0
72 : INTEGER :: id_type = do_ff_undef, itype = 0
73 : TYPE(bond_kind_type), POINTER :: bond_kind => NULL()
74 : END TYPE bond_type
75 :
76 : TYPE bend_type
77 : INTEGER :: a = 0, b = 0, c = 0
78 : INTEGER :: id_type = do_ff_undef, itype = 0
79 : TYPE(bend_kind_type), POINTER :: bend_kind => NULL()
80 : END TYPE bend_type
81 :
82 : TYPE ub_type
83 : INTEGER :: a = 0, b = 0, c = 0
84 : INTEGER :: id_type = do_ff_undef, itype = 0
85 : TYPE(ub_kind_type), POINTER :: ub_kind => NULL()
86 : END TYPE ub_type
87 :
88 : TYPE torsion_type
89 : INTEGER :: a = 0, b = 0, c = 0, d = 0
90 : INTEGER :: id_type = do_ff_undef, itype = 0
91 : TYPE(torsion_kind_type), POINTER :: torsion_kind => NULL()
92 : END TYPE torsion_type
93 :
94 : TYPE impr_type
95 : INTEGER :: a = 0, b = 0, c = 0, d = 0
96 : INTEGER :: id_type = do_ff_undef, itype = 0
97 : TYPE(impr_kind_type), POINTER :: impr_kind => NULL()
98 : END TYPE impr_type
99 :
100 : TYPE opbend_type
101 : INTEGER :: a = 0, b = 0, c = 0, d = 0
102 : INTEGER :: id_type = do_ff_undef, itype = 0
103 : TYPE(opbend_kind_type), POINTER :: opbend_kind => NULL()
104 : END TYPE opbend_type
105 :
106 : TYPE restraint_type
107 : LOGICAL :: active = .FALSE.
108 : REAL(KIND=dp) :: k0 = 0.0_dp
109 : END TYPE restraint_type
110 :
111 : ! Constraint types
112 : TYPE colvar_constraint_type
113 : INTEGER :: type_id = no_colvar_id
114 : INTEGER :: inp_seq_num = 0
115 : LOGICAL :: use_points = .FALSE.
116 : REAL(KIND=dp) :: expected_value = 0.0_dp
117 : REAL(KIND=dp) :: expected_value_growth_speed = 0.0_dp
118 : REAL(KIND=dp) :: target_velocity = 0.0_dp
119 : INTEGER, POINTER, DIMENSION(:) :: i_atoms => NULL()
120 : TYPE(restraint_type) :: restraint = restraint_type()
121 : END TYPE colvar_constraint_type
122 :
123 : TYPE g3x3_constraint_type
124 : INTEGER :: a = 0, b = 0, c = 0
125 : REAL(KIND=dp) :: dab = 0.0_dp, dac = 0.0_dp, dbc = 0.0_dp
126 : TYPE(restraint_type) :: restraint = restraint_type()
127 : END TYPE g3x3_constraint_type
128 :
129 : TYPE g4x6_constraint_type
130 : INTEGER :: a = 0, b = 0, c = 0, d = 0
131 : REAL(KIND=dp) :: dab = 0.0_dp, dac = 0.0_dp, dbc = 0.0_dp, &
132 : dad = 0.0_dp, dbd = 0.0_dp, dcd = 0.0_dp
133 : TYPE(restraint_type) :: restraint = restraint_type()
134 : END TYPE g4x6_constraint_type
135 :
136 : TYPE vsite_constraint_type
137 : INTEGER :: a = 0, b = 0, c = 0, d = 0
138 : REAL(KIND=dp) :: wbc = 0.0_dp, wdc = 0.0_dp
139 : TYPE(restraint_type) :: restraint = restraint_type()
140 : END TYPE vsite_constraint_type
141 :
142 : TYPE fixd_constraint_type
143 : TYPE(restraint_type) :: restraint = restraint_type()
144 : INTEGER :: fixd = 0, itype = 0
145 : REAL(KIND=dp), DIMENSION(3) :: coord = 0.0_dp
146 : END TYPE fixd_constraint_type
147 :
148 : TYPE local_fixd_constraint_type
149 : INTEGER :: ifixd_index = 0, ikind = 0
150 : END TYPE local_fixd_constraint_type
151 :
152 : ! Molecule kind type
153 : TYPE molecule_kind_type
154 : TYPE(atom_type), DIMENSION(:), POINTER :: atom_list => NULL()
155 : TYPE(bond_kind_type), DIMENSION(:), POINTER :: bond_kind_set => NULL()
156 : TYPE(bond_type), DIMENSION(:), POINTER :: bond_list => NULL()
157 : TYPE(bend_kind_type), DIMENSION(:), POINTER :: bend_kind_set => NULL()
158 : TYPE(bend_type), DIMENSION(:), POINTER :: bend_list => NULL()
159 : TYPE(ub_kind_type), DIMENSION(:), POINTER :: ub_kind_set => NULL()
160 : TYPE(ub_type), DIMENSION(:), POINTER :: ub_list => NULL()
161 : TYPE(torsion_kind_type), DIMENSION(:), POINTER :: torsion_kind_set => NULL()
162 : TYPE(torsion_type), DIMENSION(:), POINTER :: torsion_list => NULL()
163 : TYPE(impr_kind_type), DIMENSION(:), POINTER :: impr_kind_set => NULL()
164 : TYPE(impr_type), DIMENSION(:), POINTER :: impr_list => NULL()
165 : TYPE(opbend_kind_type), DIMENSION(:), POINTER :: opbend_kind_set => NULL()
166 : TYPE(opbend_type), DIMENSION(:), POINTER :: opbend_list => NULL()
167 : TYPE(colvar_constraint_type), DIMENSION(:), &
168 : POINTER :: colv_list => NULL()
169 : TYPE(g3x3_constraint_type), DIMENSION(:), POINTER :: g3x3_list => NULL()
170 : TYPE(g4x6_constraint_type), DIMENSION(:), POINTER :: g4x6_list => NULL()
171 : TYPE(vsite_constraint_type), DIMENSION(:), POINTER :: vsite_list => NULL()
172 : TYPE(fixd_constraint_type), DIMENSION(:), POINTER :: fixd_list => NULL()
173 : TYPE(shell_type), DIMENSION(:), POINTER :: shell_list => NULL()
174 : CHARACTER(LEN=default_string_length) :: name = ""
175 : REAL(KIND=dp) :: charge = 0.0_dp, &
176 : mass = 0.0_dp
177 : INTEGER :: kind_number = 0, &
178 : natom = 0, &
179 : nbond = 0, &
180 : nbend = 0, &
181 : nimpr = 0, &
182 : nopbend = 0, &
183 : ntorsion = 0, &
184 : nub = 0, &
185 : ng3x3 = 0, &
186 : ng3x3_restraint = 0, &
187 : ng4x6 = 0, &
188 : ng4x6_restraint = 0, &
189 : nvsite = 0, &
190 : nvsite_restraint = 0, &
191 : nfixd = 0, &
192 : nfixd_restraint = 0, &
193 : nmolecule = 0, &
194 : nshell = 0
195 : TYPE(colvar_counters) :: ncolv = colvar_counters()
196 : INTEGER :: nsgf = 0, &
197 : nelectron = 0, &
198 : nelectron_alpha = 0, &
199 : nelectron_beta = 0
200 : INTEGER, DIMENSION(:), POINTER :: molecule_list => NULL()
201 : LOGICAL :: molname_generated = .FALSE.
202 : END TYPE molecule_kind_type
203 :
204 : ! Public subroutines
205 : PUBLIC :: allocate_molecule_kind_set, &
206 : deallocate_molecule_kind_set, &
207 : get_molecule_kind, &
208 : get_molecule_kind_set, &
209 : set_molecule_kind, &
210 : write_molecule_kind_set, &
211 : setup_colvar_counters, &
212 : write_colvar_constraint, &
213 : write_fixd_constraint, &
214 : write_g3x3_constraint, &
215 : write_g4x6_constraint, &
216 : write_vsite_constraint
217 :
218 : ! Public data types
219 : PUBLIC :: atom_type, &
220 : bend_type, &
221 : bond_type, &
222 : ub_type, &
223 : torsion_type, &
224 : impr_type, &
225 : opbend_type, &
226 : colvar_constraint_type, &
227 : g3x3_constraint_type, &
228 : g4x6_constraint_type, &
229 : vsite_constraint_type, &
230 : fixd_constraint_type, &
231 : local_fixd_constraint_type, &
232 : molecule_kind_type, &
233 : shell_type
234 :
235 : CONTAINS
236 :
237 : ! **************************************************************************************************
238 : !> \brief ...
239 : !> \param colv_list ...
240 : !> \param ncolv ...
241 : ! **************************************************************************************************
242 162568 : SUBROUTINE setup_colvar_counters(colv_list, ncolv)
243 : TYPE(colvar_constraint_type), DIMENSION(:), &
244 : POINTER :: colv_list
245 : TYPE(colvar_counters), INTENT(OUT) :: ncolv
246 :
247 : INTEGER :: k
248 :
249 162568 : IF (ASSOCIATED(colv_list)) THEN
250 1120 : DO k = 1, SIZE(colv_list)
251 458 : IF (colv_list(k)%restraint%active) ncolv%nrestraint = ncolv%nrestraint + 1
252 662 : SELECT CASE (colv_list(k)%type_id)
253 : CASE (angle_colvar_id)
254 50 : ncolv%nangle = ncolv%nangle + 1
255 : CASE (coord_colvar_id)
256 2 : ncolv%ncoord = ncolv%ncoord + 1
257 : CASE (population_colvar_id)
258 0 : ncolv%npopulation = ncolv%npopulation + 1
259 : CASE (gyration_colvar_id)
260 0 : ncolv%ngyration = ncolv%ngyration + 1
261 : CASE (rotation_colvar_id)
262 0 : ncolv%nrot = ncolv%nrot + 1
263 : CASE (dist_colvar_id)
264 342 : ncolv%ndist = ncolv%ndist + 1
265 : CASE (dfunct_colvar_id)
266 4 : ncolv%ndfunct = ncolv%ndfunct + 1
267 : CASE (plane_distance_colvar_id)
268 0 : ncolv%nplane_dist = ncolv%nplane_dist + 1
269 : CASE (plane_plane_angle_colvar_id)
270 4 : ncolv%nplane_angle = ncolv%nplane_angle + 1
271 : CASE (torsion_colvar_id)
272 38 : ncolv%ntorsion = ncolv%ntorsion + 1
273 : CASE (qparm_colvar_id)
274 0 : ncolv%nqparm = ncolv%nqparm + 1
275 : CASE (xyz_diag_colvar_id)
276 6 : ncolv%nxyz_diag = ncolv%nxyz_diag + 1
277 : CASE (xyz_outerdiag_colvar_id)
278 6 : ncolv%nxyz_outerdiag = ncolv%nxyz_outerdiag + 1
279 : CASE (hydronium_shell_colvar_id)
280 0 : ncolv%nhydronium_shell = ncolv%nhydronium_shell + 1
281 : CASE (hydronium_dist_colvar_id)
282 0 : ncolv%nhydronium_dist = ncolv%nhydronium_dist + 1
283 : CASE (acid_hyd_dist_colvar_id)
284 0 : ncolv%nacid_hyd_dist = ncolv%nacid_hyd_dist + 1
285 : CASE (acid_hyd_shell_colvar_id)
286 0 : ncolv%nacid_hyd_shell = ncolv%nacid_hyd_shell + 1
287 : CASE (reaction_path_colvar_id)
288 2 : ncolv%nreactionpath = ncolv%nreactionpath + 1
289 : CASE (combine_colvar_id)
290 2 : ncolv%ncombinecvs = ncolv%ncombinecvs + 1
291 : CASE (voronoiipz_colvar_id)
292 2 : ncolv%nvoronoiipz = ncolv%nvoronoiipz + 1
293 : CASE DEFAULT
294 458 : CPABORT("Unknown colvar type")
295 : END SELECT
296 : END DO
297 : END IF
298 : ncolv%ntot = ncolv%ndist + &
299 : ncolv%nangle + &
300 : ncolv%ntorsion + &
301 : ncolv%ncoord + &
302 : ncolv%nplane_dist + &
303 : ncolv%nplane_angle + &
304 : ncolv%ndfunct + &
305 : ncolv%nrot + &
306 : ncolv%nqparm + &
307 : ncolv%nxyz_diag + &
308 : ncolv%nxyz_outerdiag + &
309 : ncolv%nhydronium_shell + &
310 : ncolv%nhydronium_dist + &
311 : ncolv%nacid_hyd_dist + &
312 : ncolv%nacid_hyd_shell + &
313 : ncolv%nreactionpath + &
314 : ncolv%ncombinecvs + &
315 : ncolv%npopulation + &
316 : ncolv%ngyration + &
317 162568 : ncolv%nvoronoiipz
318 :
319 162568 : END SUBROUTINE setup_colvar_counters
320 :
321 : ! **************************************************************************************************
322 : !> \brief Allocate and initialize a molecule kind set.
323 : !> \param molecule_kind_set ...
324 : !> \param nmolecule_kind ...
325 : !> \date 22.08.2003
326 : !> \author Matthias Krack
327 : !> \version 1.0
328 : ! **************************************************************************************************
329 12139 : SUBROUTINE allocate_molecule_kind_set(molecule_kind_set, nmolecule_kind)
330 : TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
331 : INTEGER, INTENT(IN) :: nmolecule_kind
332 :
333 : INTEGER :: imolecule_kind
334 :
335 12139 : IF (ASSOCIATED(molecule_kind_set)) THEN
336 0 : CALL deallocate_molecule_kind_set(molecule_kind_set)
337 : END IF
338 :
339 186736 : ALLOCATE (molecule_kind_set(nmolecule_kind))
340 :
341 162458 : DO imolecule_kind = 1, nmolecule_kind
342 150319 : molecule_kind_set(imolecule_kind)%kind_number = imolecule_kind
343 : CALL setup_colvar_counters(molecule_kind_set(imolecule_kind)%colv_list, &
344 162458 : molecule_kind_set(imolecule_kind)%ncolv)
345 : END DO
346 :
347 12139 : END SUBROUTINE allocate_molecule_kind_set
348 :
349 : ! **************************************************************************************************
350 : !> \brief Deallocate a molecule kind set.
351 : !> \param molecule_kind_set ...
352 : !> \date 22.08.2003
353 : !> \author Matthias Krack
354 : !> \version 1.0
355 : ! **************************************************************************************************
356 12138 : SUBROUTINE deallocate_molecule_kind_set(molecule_kind_set)
357 :
358 : TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
359 :
360 : INTEGER :: i, imolecule_kind, j, nmolecule_kind
361 :
362 12138 : IF (ASSOCIATED(molecule_kind_set)) THEN
363 :
364 12138 : nmolecule_kind = SIZE(molecule_kind_set)
365 :
366 162455 : DO imolecule_kind = 1, nmolecule_kind
367 :
368 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%atom_list)) THEN
369 150317 : DEALLOCATE (molecule_kind_set(imolecule_kind)%atom_list)
370 : END IF
371 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bend_kind_set)) THEN
372 123045 : DO i = 1, SIZE(molecule_kind_set(imolecule_kind)%bend_kind_set)
373 123045 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bend_kind_set(i)%legendre%coeffs)) THEN
374 2079 : DEALLOCATE (molecule_kind_set(imolecule_kind)%bend_kind_set(i)%legendre%coeffs)
375 : END IF
376 : END DO
377 29127 : DEALLOCATE (molecule_kind_set(imolecule_kind)%bend_kind_set)
378 : END IF
379 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bend_list)) THEN
380 150317 : DEALLOCATE (molecule_kind_set(imolecule_kind)%bend_list)
381 : END IF
382 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%ub_list)) THEN
383 150317 : DEALLOCATE (molecule_kind_set(imolecule_kind)%ub_list)
384 : END IF
385 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%ub_kind_set)) THEN
386 29113 : CALL ub_kind_dealloc_ref(molecule_kind_set(imolecule_kind)%ub_kind_set)
387 : END IF
388 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%impr_list)) THEN
389 150317 : DEALLOCATE (molecule_kind_set(imolecule_kind)%impr_list)
390 : END IF
391 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%impr_kind_set)) THEN
392 4888 : DO i = 1, SIZE(molecule_kind_set(imolecule_kind)%impr_kind_set)
393 4888 : CALL impr_kind_dealloc_ref() !This Subroutine doesn't deallocate anything, maybe needs to be implemented
394 : END DO
395 1678 : DEALLOCATE (molecule_kind_set(imolecule_kind)%impr_kind_set)
396 : END IF
397 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%opbend_list)) THEN
398 150317 : DEALLOCATE (molecule_kind_set(imolecule_kind)%opbend_list)
399 : END IF
400 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%opbend_kind_set)) THEN
401 1678 : DEALLOCATE (molecule_kind_set(imolecule_kind)%opbend_kind_set)
402 : END IF
403 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bond_kind_set)) THEN
404 29457 : DEALLOCATE (molecule_kind_set(imolecule_kind)%bond_kind_set)
405 : END IF
406 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bond_list)) THEN
407 150317 : DEALLOCATE (molecule_kind_set(imolecule_kind)%bond_list)
408 : END IF
409 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%colv_list)) THEN
410 990 : DO j = 1, SIZE(molecule_kind_set(imolecule_kind)%colv_list)
411 990 : DEALLOCATE (molecule_kind_set(imolecule_kind)%colv_list(j)%i_atoms)
412 : END DO
413 608 : DEALLOCATE (molecule_kind_set(imolecule_kind)%colv_list)
414 : END IF
415 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%g3x3_list)) THEN
416 268 : DEALLOCATE (molecule_kind_set(imolecule_kind)%g3x3_list)
417 : END IF
418 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%g4x6_list)) THEN
419 20 : DEALLOCATE (molecule_kind_set(imolecule_kind)%g4x6_list)
420 : END IF
421 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%vsite_list)) THEN
422 10 : DEALLOCATE (molecule_kind_set(imolecule_kind)%vsite_list)
423 : END IF
424 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%fixd_list)) THEN
425 4908 : DEALLOCATE (molecule_kind_set(imolecule_kind)%fixd_list)
426 : END IF
427 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%torsion_kind_set)) THEN
428 84275 : DO i = 1, SIZE(molecule_kind_set(imolecule_kind)%torsion_kind_set)
429 84275 : CALL torsion_kind_dealloc_ref(molecule_kind_set(imolecule_kind)%torsion_kind_set(i))
430 : END DO
431 5544 : DEALLOCATE (molecule_kind_set(imolecule_kind)%torsion_kind_set)
432 : END IF
433 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%shell_list)) THEN
434 10766 : DEALLOCATE (molecule_kind_set(imolecule_kind)%shell_list)
435 : END IF
436 150317 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%torsion_list)) THEN
437 150317 : DEALLOCATE (molecule_kind_set(imolecule_kind)%torsion_list)
438 : END IF
439 162455 : IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%molecule_list)) THEN
440 150317 : DEALLOCATE (molecule_kind_set(imolecule_kind)%molecule_list)
441 : END IF
442 : END DO
443 :
444 12138 : DEALLOCATE (molecule_kind_set)
445 : END IF
446 12138 : NULLIFY (molecule_kind_set)
447 :
448 12138 : END SUBROUTINE deallocate_molecule_kind_set
449 :
450 : ! **************************************************************************************************
451 : !> \brief Get informations about a molecule kind.
452 : !> \param molecule_kind ...
453 : !> \param atom_list ...
454 : !> \param bond_list ...
455 : !> \param bend_list ...
456 : !> \param ub_list ...
457 : !> \param impr_list ...
458 : !> \param opbend_list ...
459 : !> \param colv_list ...
460 : !> \param fixd_list ...
461 : !> \param g3x3_list ...
462 : !> \param g4x6_list ...
463 : !> \param vsite_list ...
464 : !> \param torsion_list ...
465 : !> \param shell_list ...
466 : !> \param name ...
467 : !> \param mass ...
468 : !> \param charge ...
469 : !> \param kind_number ...
470 : !> \param natom ...
471 : !> \param nbend ...
472 : !> \param nbond ...
473 : !> \param nub ...
474 : !> \param nimpr ...
475 : !> \param nopbend ...
476 : !> \param nconstraint ...
477 : !> \param nconstraint_fixd ...
478 : !> \param nfixd ...
479 : !> \param ncolv ...
480 : !> \param ng3x3 ...
481 : !> \param ng4x6 ...
482 : !> \param nvsite ...
483 : !> \param nfixd_restraint ...
484 : !> \param ng3x3_restraint ...
485 : !> \param ng4x6_restraint ...
486 : !> \param nvsite_restraint ...
487 : !> \param nrestraints ...
488 : !> \param nmolecule ...
489 : !> \param nsgf ...
490 : !> \param nshell ...
491 : !> \param ntorsion ...
492 : !> \param molecule_list ...
493 : !> \param nelectron ...
494 : !> \param nelectron_alpha ...
495 : !> \param nelectron_beta ...
496 : !> \param bond_kind_set ...
497 : !> \param bend_kind_set ...
498 : !> \param ub_kind_set ...
499 : !> \param impr_kind_set ...
500 : !> \param opbend_kind_set ...
501 : !> \param torsion_kind_set ...
502 : !> \param molname_generated ...
503 : !> \date 27.08.2003
504 : !> \author Matthias Krack
505 : !> \version 1.0
506 : ! **************************************************************************************************
507 15825386 : SUBROUTINE get_molecule_kind(molecule_kind, atom_list, bond_list, bend_list, &
508 : ub_list, impr_list, opbend_list, colv_list, fixd_list, &
509 : g3x3_list, g4x6_list, vsite_list, torsion_list, shell_list, &
510 : name, mass, charge, kind_number, natom, nbend, nbond, nub, &
511 : nimpr, nopbend, nconstraint, nconstraint_fixd, nfixd, ncolv, ng3x3, ng4x6, &
512 : nvsite, nfixd_restraint, ng3x3_restraint, ng4x6_restraint, &
513 : nvsite_restraint, nrestraints, nmolecule, nsgf, nshell, ntorsion, &
514 : molecule_list, nelectron, nelectron_alpha, nelectron_beta, &
515 : bond_kind_set, bend_kind_set, &
516 : ub_kind_set, impr_kind_set, opbend_kind_set, torsion_kind_set, &
517 : molname_generated)
518 :
519 : TYPE(molecule_kind_type), INTENT(IN) :: molecule_kind
520 : TYPE(atom_type), DIMENSION(:), OPTIONAL, POINTER :: atom_list
521 : TYPE(bond_type), DIMENSION(:), OPTIONAL, POINTER :: bond_list
522 : TYPE(bend_type), DIMENSION(:), OPTIONAL, POINTER :: bend_list
523 : TYPE(ub_type), DIMENSION(:), OPTIONAL, POINTER :: ub_list
524 : TYPE(impr_type), DIMENSION(:), OPTIONAL, POINTER :: impr_list
525 : TYPE(opbend_type), DIMENSION(:), OPTIONAL, POINTER :: opbend_list
526 : TYPE(colvar_constraint_type), DIMENSION(:), &
527 : OPTIONAL, POINTER :: colv_list
528 : TYPE(fixd_constraint_type), DIMENSION(:), &
529 : OPTIONAL, POINTER :: fixd_list
530 : TYPE(g3x3_constraint_type), DIMENSION(:), &
531 : OPTIONAL, POINTER :: g3x3_list
532 : TYPE(g4x6_constraint_type), DIMENSION(:), &
533 : OPTIONAL, POINTER :: g4x6_list
534 : TYPE(vsite_constraint_type), DIMENSION(:), &
535 : OPTIONAL, POINTER :: vsite_list
536 : TYPE(torsion_type), DIMENSION(:), OPTIONAL, &
537 : POINTER :: torsion_list
538 : TYPE(shell_type), DIMENSION(:), OPTIONAL, POINTER :: shell_list
539 : CHARACTER(LEN=default_string_length), &
540 : INTENT(OUT), OPTIONAL :: name
541 : REAL(KIND=dp), OPTIONAL :: mass, charge
542 : INTEGER, INTENT(OUT), OPTIONAL :: kind_number, natom, nbend, nbond, nub, &
543 : nimpr, nopbend, nconstraint, &
544 : nconstraint_fixd, nfixd
545 : TYPE(colvar_counters), INTENT(out), OPTIONAL :: ncolv
546 : INTEGER, INTENT(OUT), OPTIONAL :: ng3x3, ng4x6, nvsite, nfixd_restraint, ng3x3_restraint, &
547 : ng4x6_restraint, nvsite_restraint, nrestraints, nmolecule, nsgf, nshell, ntorsion
548 : INTEGER, DIMENSION(:), OPTIONAL, POINTER :: molecule_list
549 : INTEGER, INTENT(OUT), OPTIONAL :: nelectron, nelectron_alpha, &
550 : nelectron_beta
551 : TYPE(bond_kind_type), DIMENSION(:), OPTIONAL, &
552 : POINTER :: bond_kind_set
553 : TYPE(bend_kind_type), DIMENSION(:), OPTIONAL, &
554 : POINTER :: bend_kind_set
555 : TYPE(ub_kind_type), DIMENSION(:), OPTIONAL, &
556 : POINTER :: ub_kind_set
557 : TYPE(impr_kind_type), DIMENSION(:), OPTIONAL, &
558 : POINTER :: impr_kind_set
559 : TYPE(opbend_kind_type), DIMENSION(:), OPTIONAL, &
560 : POINTER :: opbend_kind_set
561 : TYPE(torsion_kind_type), DIMENSION(:), OPTIONAL, &
562 : POINTER :: torsion_kind_set
563 : LOGICAL, INTENT(OUT), OPTIONAL :: molname_generated
564 :
565 : INTEGER :: i
566 :
567 15825386 : IF (PRESENT(atom_list)) atom_list => molecule_kind%atom_list
568 15825386 : IF (PRESENT(bend_list)) bend_list => molecule_kind%bend_list
569 15825386 : IF (PRESENT(bond_list)) bond_list => molecule_kind%bond_list
570 15825386 : IF (PRESENT(impr_list)) impr_list => molecule_kind%impr_list
571 15825386 : IF (PRESENT(opbend_list)) opbend_list => molecule_kind%opbend_list
572 15825386 : IF (PRESENT(ub_list)) ub_list => molecule_kind%ub_list
573 15825386 : IF (PRESENT(bond_kind_set)) bond_kind_set => molecule_kind%bond_kind_set
574 15825386 : IF (PRESENT(bend_kind_set)) bend_kind_set => molecule_kind%bend_kind_set
575 15825386 : IF (PRESENT(ub_kind_set)) ub_kind_set => molecule_kind%ub_kind_set
576 15825386 : IF (PRESENT(impr_kind_set)) impr_kind_set => molecule_kind%impr_kind_set
577 15825386 : IF (PRESENT(opbend_kind_set)) opbend_kind_set => molecule_kind%opbend_kind_set
578 15825386 : IF (PRESENT(torsion_kind_set)) torsion_kind_set => molecule_kind%torsion_kind_set
579 15825386 : IF (PRESENT(colv_list)) colv_list => molecule_kind%colv_list
580 15825386 : IF (PRESENT(g3x3_list)) g3x3_list => molecule_kind%g3x3_list
581 15825386 : IF (PRESENT(g4x6_list)) g4x6_list => molecule_kind%g4x6_list
582 15825386 : IF (PRESENT(vsite_list)) vsite_list => molecule_kind%vsite_list
583 15825386 : IF (PRESENT(fixd_list)) fixd_list => molecule_kind%fixd_list
584 15825386 : IF (PRESENT(torsion_list)) torsion_list => molecule_kind%torsion_list
585 15825386 : IF (PRESENT(shell_list)) shell_list => molecule_kind%shell_list
586 15825386 : IF (PRESENT(name)) name = molecule_kind%name
587 15825386 : IF (PRESENT(molname_generated)) molname_generated = molecule_kind%molname_generated
588 15825386 : IF (PRESENT(mass)) mass = molecule_kind%mass
589 15825386 : IF (PRESENT(charge)) charge = molecule_kind%charge
590 15825386 : IF (PRESENT(kind_number)) kind_number = molecule_kind%kind_number
591 15825386 : IF (PRESENT(natom)) natom = molecule_kind%natom
592 15825386 : IF (PRESENT(nbend)) nbend = molecule_kind%nbend
593 15825386 : IF (PRESENT(nbond)) nbond = molecule_kind%nbond
594 15825386 : IF (PRESENT(nub)) nub = molecule_kind%nub
595 15825386 : IF (PRESENT(nimpr)) nimpr = molecule_kind%nimpr
596 15825386 : IF (PRESENT(nopbend)) nopbend = molecule_kind%nopbend
597 15825386 : IF (PRESENT(nconstraint)) nconstraint = (molecule_kind%ncolv%ntot - molecule_kind%ncolv%nrestraint) + &
598 : 3*(molecule_kind%ng3x3 - molecule_kind%ng3x3_restraint) + &
599 : 6*(molecule_kind%ng4x6 - molecule_kind%ng4x6_restraint) + &
600 3145531 : 3*(molecule_kind%nvsite - molecule_kind%nvsite_restraint)
601 15825386 : IF (PRESENT(ncolv)) ncolv = molecule_kind%ncolv
602 15825386 : IF (PRESENT(ng3x3)) ng3x3 = molecule_kind%ng3x3
603 15825386 : IF (PRESENT(ng4x6)) ng4x6 = molecule_kind%ng4x6
604 15825386 : IF (PRESENT(nvsite)) nvsite = molecule_kind%nvsite
605 : ! Number of atoms that have one or more components fixed
606 15825386 : IF (PRESENT(nfixd)) nfixd = molecule_kind%nfixd
607 : ! Number of degrees of freedom fixed
608 15825386 : IF (PRESENT(nconstraint_fixd)) THEN
609 295070 : nconstraint_fixd = 0
610 295070 : IF (molecule_kind%nfixd /= 0) THEN
611 171892 : DO i = 1, SIZE(molecule_kind%fixd_list)
612 170007 : IF (molecule_kind%fixd_list(i)%restraint%active) CYCLE
613 1885 : SELECT CASE (molecule_kind%fixd_list(i)%itype)
614 : CASE (use_perd_x, use_perd_y, use_perd_z)
615 62976 : nconstraint_fixd = nconstraint_fixd + 1
616 : CASE (use_perd_xy, use_perd_xz, use_perd_yz)
617 20992 : nconstraint_fixd = nconstraint_fixd + 2
618 : CASE (use_perd_xyz)
619 169579 : nconstraint_fixd = nconstraint_fixd + 3
620 : END SELECT
621 : END DO
622 : END IF
623 : END IF
624 15825386 : IF (PRESENT(ng3x3_restraint)) ng3x3_restraint = molecule_kind%ng3x3_restraint
625 15825386 : IF (PRESENT(ng4x6_restraint)) ng4x6_restraint = molecule_kind%ng4x6_restraint
626 15825386 : IF (PRESENT(nvsite_restraint)) nvsite_restraint = molecule_kind%nvsite_restraint
627 15825386 : IF (PRESENT(nfixd_restraint)) nfixd_restraint = molecule_kind%nfixd_restraint
628 15825386 : IF (PRESENT(nrestraints)) nrestraints = molecule_kind%ncolv%nrestraint + &
629 : molecule_kind%ng3x3_restraint + &
630 : molecule_kind%ng4x6_restraint + &
631 280070 : molecule_kind%nvsite_restraint
632 15825386 : IF (PRESENT(nmolecule)) nmolecule = molecule_kind%nmolecule
633 15825386 : IF (PRESENT(nshell)) nshell = molecule_kind%nshell
634 15825386 : IF (PRESENT(ntorsion)) ntorsion = molecule_kind%ntorsion
635 15825386 : IF (PRESENT(nsgf)) nsgf = molecule_kind%nsgf
636 15825386 : IF (PRESENT(nelectron)) nelectron = molecule_kind%nelectron
637 15825386 : IF (PRESENT(nelectron_alpha)) nelectron_alpha = molecule_kind%nelectron_beta
638 15825386 : IF (PRESENT(nelectron_beta)) nelectron_beta = molecule_kind%nelectron_alpha
639 15825386 : IF (PRESENT(molecule_list)) molecule_list => molecule_kind%molecule_list
640 :
641 15825386 : END SUBROUTINE get_molecule_kind
642 :
643 : ! **************************************************************************************************
644 : !> \brief Get informations about a molecule kind set.
645 : !> \param molecule_kind_set ...
646 : !> \param maxatom ...
647 : !> \param natom ...
648 : !> \param nbond ...
649 : !> \param nbend ...
650 : !> \param nub ...
651 : !> \param ntorsion ...
652 : !> \param nimpr ...
653 : !> \param nopbend ...
654 : !> \param nconstraint ...
655 : !> \param nconstraint_fixd ...
656 : !> \param nmolecule ...
657 : !> \param nrestraints ...
658 : !> \date 27.08.2003
659 : !> \author Matthias Krack
660 : !> \version 1.0
661 : ! **************************************************************************************************
662 51730 : SUBROUTINE get_molecule_kind_set(molecule_kind_set, maxatom, natom, &
663 : nbond, nbend, nub, ntorsion, nimpr, nopbend, &
664 : nconstraint, nconstraint_fixd, nmolecule, &
665 : nrestraints)
666 :
667 : TYPE(molecule_kind_type), DIMENSION(:), INTENT(IN) :: molecule_kind_set
668 : INTEGER, INTENT(OUT), OPTIONAL :: maxatom, natom, nbond, nbend, nub, &
669 : ntorsion, nimpr, nopbend, nconstraint, &
670 : nconstraint_fixd, nmolecule, &
671 : nrestraints
672 :
673 : INTEGER :: ibend, ibond, iimpr, imolecule_kind, iopbend, itorsion, iub, na, nc, nc_fixd, &
674 : nfixd_restraint, nm, nmolecule_kind, nrestraints_tot
675 :
676 51730 : IF (PRESENT(maxatom)) maxatom = 0
677 51730 : IF (PRESENT(natom)) natom = 0
678 51730 : IF (PRESENT(nbond)) nbond = 0
679 51730 : IF (PRESENT(nbend)) nbend = 0
680 51730 : IF (PRESENT(nub)) nub = 0
681 51730 : IF (PRESENT(ntorsion)) ntorsion = 0
682 51730 : IF (PRESENT(nimpr)) nimpr = 0
683 51730 : IF (PRESENT(nopbend)) nopbend = 0
684 51730 : IF (PRESENT(nconstraint)) nconstraint = 0
685 51730 : IF (PRESENT(nconstraint_fixd)) nconstraint_fixd = 0
686 51730 : IF (PRESENT(nrestraints)) nrestraints = 0
687 51730 : IF (PRESENT(nmolecule)) nmolecule = 0
688 :
689 51730 : nmolecule_kind = SIZE(molecule_kind_set)
690 :
691 331800 : DO imolecule_kind = 1, nmolecule_kind
692 51730 : ASSOCIATE (molecule_kind => molecule_kind_set(imolecule_kind))
693 :
694 : CALL get_molecule_kind(molecule_kind=molecule_kind, &
695 : natom=na, &
696 : nbond=ibond, &
697 : nbend=ibend, &
698 : nub=iub, &
699 : ntorsion=itorsion, &
700 : nimpr=iimpr, &
701 : nopbend=iopbend, &
702 : nconstraint=nc, &
703 : nconstraint_fixd=nc_fixd, &
704 : nfixd_restraint=nfixd_restraint, &
705 : nrestraints=nrestraints_tot, &
706 280070 : nmolecule=nm)
707 280070 : IF (PRESENT(maxatom)) maxatom = MAX(maxatom, na)
708 280070 : IF (PRESENT(natom)) natom = natom + na*nm
709 280070 : IF (PRESENT(nbond)) nbond = nbond + ibond*nm
710 280070 : IF (PRESENT(nbend)) nbend = nbend + ibend*nm
711 280070 : IF (PRESENT(nub)) nub = nub + iub*nm
712 280070 : IF (PRESENT(ntorsion)) ntorsion = ntorsion + itorsion*nm
713 280070 : IF (PRESENT(nimpr)) nimpr = nimpr + iimpr*nm
714 280070 : IF (PRESENT(nopbend)) nopbend = nopbend + iopbend*nm
715 280070 : IF (PRESENT(nconstraint)) nconstraint = nconstraint + nc*nm + nc_fixd
716 280070 : IF (PRESENT(nconstraint_fixd)) nconstraint_fixd = nconstraint_fixd + nc_fixd
717 280070 : IF (PRESENT(nmolecule)) nmolecule = nmolecule + nm
718 560140 : IF (PRESENT(nrestraints)) nrestraints = nrestraints + nm*nrestraints_tot + nfixd_restraint
719 :
720 : END ASSOCIATE
721 : END DO
722 :
723 51730 : END SUBROUTINE get_molecule_kind_set
724 :
725 : ! **************************************************************************************************
726 : !> \brief Set the components of a molecule kind.
727 : !> \param molecule_kind ...
728 : !> \param name ...
729 : !> \param mass ...
730 : !> \param charge ...
731 : !> \param kind_number ...
732 : !> \param molecule_list ...
733 : !> \param atom_list ...
734 : !> \param nbond ...
735 : !> \param bond_list ...
736 : !> \param nbend ...
737 : !> \param bend_list ...
738 : !> \param nub ...
739 : !> \param ub_list ...
740 : !> \param nimpr ...
741 : !> \param impr_list ...
742 : !> \param nopbend ...
743 : !> \param opbend_list ...
744 : !> \param ntorsion ...
745 : !> \param torsion_list ...
746 : !> \param fixd_list ...
747 : !> \param ncolv ...
748 : !> \param colv_list ...
749 : !> \param ng3x3 ...
750 : !> \param g3x3_list ...
751 : !> \param ng4x6 ...
752 : !> \param nfixd ...
753 : !> \param g4x6_list ...
754 : !> \param nvsite ...
755 : !> \param vsite_list ...
756 : !> \param ng3x3_restraint ...
757 : !> \param ng4x6_restraint ...
758 : !> \param nfixd_restraint ...
759 : !> \param nshell ...
760 : !> \param shell_list ...
761 : !> \param nvsite_restraint ...
762 : !> \param bond_kind_set ...
763 : !> \param bend_kind_set ...
764 : !> \param ub_kind_set ...
765 : !> \param torsion_kind_set ...
766 : !> \param impr_kind_set ...
767 : !> \param opbend_kind_set ...
768 : !> \param nelectron ...
769 : !> \param nsgf ...
770 : !> \param molname_generated ...
771 : !> \date 27.08.2003
772 : !> \author Matthias Krack
773 : !> \version 1.0
774 : ! **************************************************************************************************
775 2198931 : SUBROUTINE set_molecule_kind(molecule_kind, name, mass, charge, kind_number, &
776 : molecule_list, atom_list, nbond, bond_list, &
777 : nbend, bend_list, nub, ub_list, nimpr, impr_list, &
778 : nopbend, opbend_list, ntorsion, &
779 : torsion_list, fixd_list, ncolv, colv_list, ng3x3, &
780 : g3x3_list, ng4x6, nfixd, g4x6_list, nvsite, &
781 : vsite_list, ng3x3_restraint, ng4x6_restraint, &
782 : nfixd_restraint, nshell, shell_list, &
783 : nvsite_restraint, bond_kind_set, bend_kind_set, &
784 : ub_kind_set, torsion_kind_set, impr_kind_set, &
785 : opbend_kind_set, nelectron, nsgf, &
786 : molname_generated)
787 :
788 : TYPE(molecule_kind_type), INTENT(INOUT) :: molecule_kind
789 : CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: name
790 : REAL(KIND=dp), OPTIONAL :: mass, charge
791 : INTEGER, INTENT(IN), OPTIONAL :: kind_number
792 : INTEGER, DIMENSION(:), OPTIONAL, POINTER :: molecule_list
793 : TYPE(atom_type), DIMENSION(:), OPTIONAL, POINTER :: atom_list
794 : INTEGER, INTENT(IN), OPTIONAL :: nbond
795 : TYPE(bond_type), DIMENSION(:), OPTIONAL, POINTER :: bond_list
796 : INTEGER, INTENT(IN), OPTIONAL :: nbend
797 : TYPE(bend_type), DIMENSION(:), OPTIONAL, POINTER :: bend_list
798 : INTEGER, INTENT(IN), OPTIONAL :: nub
799 : TYPE(ub_type), DIMENSION(:), OPTIONAL, POINTER :: ub_list
800 : INTEGER, INTENT(IN), OPTIONAL :: nimpr
801 : TYPE(impr_type), DIMENSION(:), OPTIONAL, POINTER :: impr_list
802 : INTEGER, INTENT(IN), OPTIONAL :: nopbend
803 : TYPE(opbend_type), DIMENSION(:), OPTIONAL, POINTER :: opbend_list
804 : INTEGER, INTENT(IN), OPTIONAL :: ntorsion
805 : TYPE(torsion_type), DIMENSION(:), OPTIONAL, &
806 : POINTER :: torsion_list
807 : TYPE(fixd_constraint_type), DIMENSION(:), &
808 : OPTIONAL, POINTER :: fixd_list
809 : TYPE(colvar_counters), INTENT(IN), OPTIONAL :: ncolv
810 : TYPE(colvar_constraint_type), DIMENSION(:), &
811 : OPTIONAL, POINTER :: colv_list
812 : INTEGER, INTENT(IN), OPTIONAL :: ng3x3
813 : TYPE(g3x3_constraint_type), DIMENSION(:), &
814 : OPTIONAL, POINTER :: g3x3_list
815 : INTEGER, INTENT(IN), OPTIONAL :: ng4x6, nfixd
816 : TYPE(g4x6_constraint_type), DIMENSION(:), &
817 : OPTIONAL, POINTER :: g4x6_list
818 : INTEGER, INTENT(IN), OPTIONAL :: nvsite
819 : TYPE(vsite_constraint_type), DIMENSION(:), &
820 : OPTIONAL, POINTER :: vsite_list
821 : INTEGER, INTENT(IN), OPTIONAL :: ng3x3_restraint, ng4x6_restraint, &
822 : nfixd_restraint, nshell
823 : TYPE(shell_type), DIMENSION(:), OPTIONAL, POINTER :: shell_list
824 : INTEGER, INTENT(IN), OPTIONAL :: nvsite_restraint
825 : TYPE(bond_kind_type), DIMENSION(:), OPTIONAL, &
826 : POINTER :: bond_kind_set
827 : TYPE(bend_kind_type), DIMENSION(:), OPTIONAL, &
828 : POINTER :: bend_kind_set
829 : TYPE(ub_kind_type), DIMENSION(:), OPTIONAL, &
830 : POINTER :: ub_kind_set
831 : TYPE(torsion_kind_type), DIMENSION(:), OPTIONAL, &
832 : POINTER :: torsion_kind_set
833 : TYPE(impr_kind_type), DIMENSION(:), OPTIONAL, &
834 : POINTER :: impr_kind_set
835 : TYPE(opbend_kind_type), DIMENSION(:), OPTIONAL, &
836 : POINTER :: opbend_kind_set
837 : INTEGER, INTENT(IN), OPTIONAL :: nelectron, nsgf
838 : LOGICAL, INTENT(IN), OPTIONAL :: molname_generated
839 :
840 : INTEGER :: n
841 :
842 2198931 : IF (PRESENT(atom_list)) THEN
843 300638 : n = SIZE(atom_list)
844 300638 : molecule_kind%natom = n
845 300638 : molecule_kind%atom_list => atom_list
846 : END IF
847 2198931 : IF (PRESENT(molname_generated)) molecule_kind%molname_generated = molname_generated
848 2198931 : IF (PRESENT(name)) molecule_kind%name = name
849 2198931 : IF (PRESENT(mass)) molecule_kind%mass = mass
850 2198931 : IF (PRESENT(charge)) molecule_kind%charge = charge
851 2198931 : IF (PRESENT(kind_number)) molecule_kind%kind_number = kind_number
852 2198931 : IF (PRESENT(nbond)) molecule_kind%nbond = nbond
853 2198931 : IF (PRESENT(bond_list)) molecule_kind%bond_list => bond_list
854 2198931 : IF (PRESENT(nbend)) molecule_kind%nbend = nbend
855 2198931 : IF (PRESENT(nelectron)) molecule_kind%nelectron = nelectron
856 2198931 : IF (PRESENT(nsgf)) molecule_kind%nsgf = nsgf
857 2198931 : IF (PRESENT(bend_list)) molecule_kind%bend_list => bend_list
858 2198931 : IF (PRESENT(nub)) molecule_kind%nub = nub
859 2198931 : IF (PRESENT(ub_list)) molecule_kind%ub_list => ub_list
860 2198931 : IF (PRESENT(ntorsion)) molecule_kind%ntorsion = ntorsion
861 2198931 : IF (PRESENT(torsion_list)) molecule_kind%torsion_list => torsion_list
862 2198931 : IF (PRESENT(nimpr)) molecule_kind%nimpr = nimpr
863 2198931 : IF (PRESENT(impr_list)) molecule_kind%impr_list => impr_list
864 2198931 : IF (PRESENT(nopbend)) molecule_kind%nopbend = nopbend
865 2198931 : IF (PRESENT(opbend_list)) molecule_kind%opbend_list => opbend_list
866 2198931 : IF (PRESENT(ncolv)) molecule_kind%ncolv = ncolv
867 2198931 : IF (PRESENT(colv_list)) molecule_kind%colv_list => colv_list
868 2198931 : IF (PRESENT(ng3x3)) molecule_kind%ng3x3 = ng3x3
869 2198931 : IF (PRESENT(g3x3_list)) molecule_kind%g3x3_list => g3x3_list
870 2198931 : IF (PRESENT(ng4x6)) molecule_kind%ng4x6 = ng4x6
871 2198931 : IF (PRESENT(nvsite)) molecule_kind%nvsite = nvsite
872 2198931 : IF (PRESENT(nfixd)) molecule_kind%nfixd = nfixd
873 2198931 : IF (PRESENT(nfixd_restraint)) molecule_kind%nfixd_restraint = nfixd_restraint
874 2198931 : IF (PRESENT(ng3x3_restraint)) molecule_kind%ng3x3_restraint = ng3x3_restraint
875 2198931 : IF (PRESENT(ng4x6_restraint)) molecule_kind%ng4x6_restraint = ng4x6_restraint
876 2198931 : IF (PRESENT(nvsite_restraint)) molecule_kind%nvsite_restraint = nvsite_restraint
877 2198931 : IF (PRESENT(g4x6_list)) molecule_kind%g4x6_list => g4x6_list
878 2198931 : IF (PRESENT(vsite_list)) molecule_kind%vsite_list => vsite_list
879 2198931 : IF (PRESENT(fixd_list)) molecule_kind%fixd_list => fixd_list
880 2198931 : IF (PRESENT(bond_kind_set)) molecule_kind%bond_kind_set => bond_kind_set
881 2198931 : IF (PRESENT(bend_kind_set)) molecule_kind%bend_kind_set => bend_kind_set
882 2198931 : IF (PRESENT(ub_kind_set)) molecule_kind%ub_kind_set => ub_kind_set
883 2198931 : IF (PRESENT(torsion_kind_set)) molecule_kind%torsion_kind_set => torsion_kind_set
884 2198931 : IF (PRESENT(impr_kind_set)) molecule_kind%impr_kind_set => impr_kind_set
885 2198931 : IF (PRESENT(opbend_kind_set)) molecule_kind%opbend_kind_set => opbend_kind_set
886 2198931 : IF (PRESENT(nshell)) molecule_kind%nshell = nshell
887 2198931 : IF (PRESENT(shell_list)) molecule_kind%shell_list => shell_list
888 2198931 : IF (PRESENT(molecule_list)) THEN
889 150319 : n = SIZE(molecule_list)
890 150319 : molecule_kind%nmolecule = n
891 150319 : molecule_kind%molecule_list => molecule_list
892 : END IF
893 2198931 : END SUBROUTINE set_molecule_kind
894 :
895 : ! **************************************************************************************************
896 : !> \brief Write a molecule kind data set to the output unit.
897 : !> \param molecule_kind ...
898 : !> \param output_unit ...
899 : !> \date 24.09.2003
900 : !> \author Matthias Krack
901 : !> \version 1.0
902 : ! **************************************************************************************************
903 2253 : SUBROUTINE write_molecule_kind(molecule_kind, output_unit)
904 : TYPE(molecule_kind_type), INTENT(IN) :: molecule_kind
905 : INTEGER, INTENT(in) :: output_unit
906 :
907 : CHARACTER(LEN=default_string_length) :: name
908 : INTEGER :: iatom, imolecule, natom, nmolecule
909 : TYPE(atomic_kind_type), POINTER :: atomic_kind
910 :
911 2253 : IF (output_unit > 0) THEN
912 2253 : natom = SIZE(molecule_kind%atom_list)
913 2253 : nmolecule = SIZE(molecule_kind%molecule_list)
914 :
915 2253 : IF (natom == 1) THEN
916 211 : atomic_kind => molecule_kind%atom_list(1)%atomic_kind
917 211 : CALL get_atomic_kind(atomic_kind=atomic_kind, name=name)
918 : WRITE (UNIT=output_unit, FMT="(/,T2,I5,A,T36,A,A,T64,A)") &
919 211 : molecule_kind%kind_number, &
920 211 : ". Molecule kind: "//TRIM(molecule_kind%name), &
921 422 : "Atomic kind name: ", TRIM(name)
922 : WRITE (UNIT=output_unit, FMT="(T9,A,L1,T55,A,T75,I6)") &
923 211 : "Automatic name: ", molecule_kind%molname_generated, &
924 422 : "Number of molecules:", nmolecule
925 : ELSE
926 : WRITE (UNIT=output_unit, FMT="(/,T2,I5,A,T50,A,T75,I6,/,T22,A)") &
927 2042 : molecule_kind%kind_number, &
928 2042 : ". Molecule kind: "//TRIM(molecule_kind%name), &
929 2042 : "Number of atoms: ", natom, &
930 4084 : "Atom Atomic kind name"
931 17246 : DO iatom = 1, natom
932 15204 : atomic_kind => molecule_kind%atom_list(iatom)%atomic_kind
933 15204 : CALL get_atomic_kind(atomic_kind=atomic_kind, name=name)
934 : WRITE (UNIT=output_unit, FMT="(T20,I6,(7X,A18))") &
935 17246 : iatom, TRIM(name)
936 : END DO
937 : WRITE (UNIT=output_unit, FMT="(/,T9,A,L1)") &
938 2042 : "The name was automatically generated: ", &
939 4084 : molecule_kind%molname_generated
940 : WRITE (UNIT=output_unit, FMT="(T9,A,I6,/,T9,A,(T30,5I10))") &
941 2042 : "Number of molecules: ", nmolecule, "Molecule list:", &
942 33860 : (molecule_kind%molecule_list(imolecule), imolecule=1, nmolecule)
943 2042 : IF (molecule_kind%nbond > 0) THEN
944 : WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
945 1767 : "Number of bonds: ", molecule_kind%nbond
946 : END IF
947 2042 : IF (molecule_kind%nbend > 0) THEN
948 : WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
949 1621 : "Number of bends: ", molecule_kind%nbend
950 : END IF
951 2042 : IF (molecule_kind%nub > 0) THEN
952 : WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
953 266 : "Number of Urey-Bradley:", molecule_kind%nub
954 : END IF
955 2042 : IF (molecule_kind%ntorsion > 0) THEN
956 : WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
957 1125 : "Number of torsions: ", molecule_kind%ntorsion
958 : END IF
959 2042 : IF (molecule_kind%nimpr > 0) THEN
960 : WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
961 179 : "Number of improper: ", molecule_kind%nimpr
962 : END IF
963 2042 : IF (molecule_kind%nopbend > 0) THEN
964 : WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
965 4 : "Number of out opbends: ", molecule_kind%nopbend
966 : END IF
967 : END IF
968 : END IF
969 2253 : END SUBROUTINE write_molecule_kind
970 :
971 : ! **************************************************************************************************
972 : !> \brief Write a moleculeatomic kind set data set to the output unit.
973 : !> \param molecule_kind_set ...
974 : !> \param subsys_section ...
975 : !> \date 24.09.2003
976 : !> \author Matthias Krack
977 : !> \version 1.0
978 : ! **************************************************************************************************
979 12104 : SUBROUTINE write_molecule_kind_set(molecule_kind_set, subsys_section)
980 : TYPE(molecule_kind_type), DIMENSION(:), INTENT(IN) :: molecule_kind_set
981 : TYPE(section_vals_type), INTENT(IN) :: subsys_section
982 :
983 : CHARACTER(len=*), PARAMETER :: routineN = 'write_molecule_kind_set'
984 :
985 : INTEGER :: handle, imolecule_kind, natom, nbend, &
986 : nbond, nimpr, nmolecule, &
987 : nmolecule_kind, nopbend, ntors, &
988 : ntotal, nub, output_unit
989 : LOGICAL :: all_single_atoms
990 : TYPE(cp_logger_type), POINTER :: logger
991 :
992 12104 : CALL timeset(routineN, handle)
993 :
994 12104 : NULLIFY (logger)
995 12104 : logger => cp_get_default_logger()
996 : output_unit = cp_print_key_unit_nr(logger, subsys_section, &
997 12104 : "PRINT%MOLECULES", extension=".Log")
998 12104 : IF (output_unit > 0) THEN
999 2856 : WRITE (UNIT=output_unit, FMT="(/,/,T2,A)") "MOLECULE KIND INFORMATION"
1000 :
1001 2856 : nmolecule_kind = SIZE(molecule_kind_set)
1002 :
1003 2856 : all_single_atoms = .TRUE.
1004 33119 : DO imolecule_kind = 1, nmolecule_kind
1005 30263 : natom = SIZE(molecule_kind_set(imolecule_kind)%atom_list)
1006 30263 : nmolecule = SIZE(molecule_kind_set(imolecule_kind)%molecule_list)
1007 33119 : IF (natom*nmolecule > 1) all_single_atoms = .FALSE.
1008 : END DO
1009 :
1010 2856 : IF (all_single_atoms) THEN
1011 : WRITE (UNIT=output_unit, FMT="(/,/,T2,A)") &
1012 2158 : "All atoms are their own molecule, skipping detailed information"
1013 : ELSE
1014 2951 : DO imolecule_kind = 1, nmolecule_kind
1015 2951 : CALL write_molecule_kind(molecule_kind_set(imolecule_kind), output_unit)
1016 : END DO
1017 : END IF
1018 :
1019 : CALL get_molecule_kind_set(molecule_kind_set=molecule_kind_set, &
1020 : nbond=nbond, &
1021 : nbend=nbend, &
1022 : nub=nub, &
1023 : ntorsion=ntors, &
1024 : nimpr=nimpr, &
1025 2856 : nopbend=nopbend)
1026 2856 : ntotal = nbond + nbend + nub + ntors + nimpr + nopbend
1027 2856 : IF (ntotal > 0) THEN
1028 : WRITE (UNIT=output_unit, FMT="(/,/,T2,A,T45,A30,I6)") &
1029 596 : "MOLECULE KIND SET INFORMATION", &
1030 1192 : "Total Number of bonds: ", nbond
1031 : WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
1032 596 : "Total Number of bends: ", nbend
1033 : WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
1034 596 : "Total Number of Urey-Bradley:", nub
1035 : WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
1036 596 : "Total Number of torsions: ", ntors
1037 : WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
1038 596 : "Total Number of improper: ", nimpr
1039 : WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
1040 596 : "Total Number of opbends: ", nopbend
1041 : END IF
1042 : END IF
1043 : CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
1044 12104 : "PRINT%MOLECULES")
1045 :
1046 12104 : CALL timestop(handle)
1047 :
1048 12104 : END SUBROUTINE write_molecule_kind_set
1049 :
1050 : ! **************************************************************************************************
1051 : !> \brief Write collective variable constraint information to output unit
1052 : !> \param colvar_constraint Data set of the collective variable constraint
1053 : !> \param icolv Collective variable number (index)
1054 : !> \param iw Logical unit number of the output unit
1055 : !> \author Matthias Krack (25.11.2025)
1056 : ! **************************************************************************************************
1057 30 : SUBROUTINE write_colvar_constraint(colvar_constraint, icolv, iw)
1058 :
1059 : TYPE(colvar_constraint_type), INTENT(IN), POINTER :: colvar_constraint
1060 : INTEGER, INTENT(IN) :: icolv, iw
1061 :
1062 : CHARACTER(LEN=30) :: type_string
1063 :
1064 30 : IF (iw > 0) THEN
1065 30 : CPASSERT(ASSOCIATED(colvar_constraint))
1066 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1067 30 : "COLVAR| Number", icolv
1068 30 : SELECT CASE (colvar_constraint%type_id)
1069 : CASE (no_colvar_id)
1070 0 : type_string = "Undefined"
1071 : CASE (dist_colvar_id)
1072 23 : type_string = "Distance"
1073 : CASE (coord_colvar_id)
1074 1 : type_string = "Coordination number"
1075 : CASE (torsion_colvar_id)
1076 0 : type_string = "Torsion"
1077 : CASE (angle_colvar_id)
1078 2 : type_string = "Angle"
1079 : CASE (plane_distance_colvar_id)
1080 0 : type_string = "Plane distance"
1081 : CASE (rotation_colvar_id)
1082 0 : type_string = "Rotation"
1083 : CASE (dfunct_colvar_id)
1084 2 : type_string = "Distance function"
1085 : CASE (qparm_colvar_id)
1086 0 : type_string = "Q parameter"
1087 : CASE (hydronium_shell_colvar_id)
1088 0 : type_string = "Hydronium shell"
1089 : CASE (reaction_path_colvar_id)
1090 0 : type_string = "Reaction path"
1091 : CASE (combine_colvar_id)
1092 0 : type_string = "Combine"
1093 : CASE (population_colvar_id)
1094 0 : type_string = "Population"
1095 : CASE (plane_plane_angle_colvar_id)
1096 2 : type_string = "Angle plane-plane"
1097 : CASE (gyration_colvar_id)
1098 0 : type_string = "Gyration radius"
1099 : CASE (rmsd_colvar_id)
1100 0 : type_string = "RMSD"
1101 : CASE (distance_from_path_colvar_id)
1102 0 : type_string = "Distance from path"
1103 : CASE (xyz_diag_colvar_id)
1104 0 : type_string = "XYZ diag"
1105 : CASE (xyz_outerdiag_colvar_id)
1106 0 : type_string = "XYZ outerdiag"
1107 : CASE (u_colvar_id)
1108 0 : type_string = "U"
1109 : CASE (Wc_colvar_id)
1110 0 : type_string = "WC"
1111 : CASE (HBP_colvar_id)
1112 0 : type_string = "HBP"
1113 : CASE (ring_puckering_colvar_id)
1114 0 : type_string = "Ring puckering"
1115 : CASE (mindist_colvar_id)
1116 0 : type_string = "Distance point-plane"
1117 : CASE (acid_hyd_dist_colvar_id)
1118 0 : type_string = "Acid hydronium distance"
1119 : CASE (acid_hyd_shell_colvar_id)
1120 0 : type_string = "Acid hydronium shell"
1121 : CASE (hydronium_dist_colvar_id)
1122 0 : type_string = "Hydronium distance"
1123 : CASE (voronoiipz_colvar_id)
1124 0 : type_string = "Voronoi ion position"
1125 : CASE DEFAULT
1126 30 : CPABORT("Invalid collective variable ID specified. Check the code!")
1127 : END SELECT
1128 30 : IF (colvar_constraint%restraint%active) THEN
1129 : WRITE (UNIT=iw, FMT="(T2,A,T51,A30)") &
1130 7 : "COLVAR| Restraint type", ADJUSTR(TRIM(type_string))
1131 : WRITE (UNIT=iw, FMT="(T2,A,T66,ES15.6)") &
1132 7 : "COLVAR| Restraint constant k [a.u.]", colvar_constraint%restraint%k0
1133 : ELSE
1134 : WRITE (UNIT=iw, FMT="(T2,A,T51,A30)") &
1135 23 : "COLVAR| Constraint type", ADJUSTR(TRIM(type_string))
1136 : END IF
1137 : WRITE (UNIT=iw, FMT="(T2,A,T66,ES15.6)") &
1138 30 : "COLVAR| Target value", colvar_constraint%expected_value, &
1139 60 : "COLVAR| Target value growth speed", colvar_constraint%expected_value_growth_speed
1140 30 : IF (colvar_constraint%use_points) THEN
1141 4 : WRITE (UNIT=iw, FMT="(T2,A,T78,A3)") "COLVAR| Use points", "Yes"
1142 : ELSE
1143 26 : WRITE (UNIT=iw, FMT="(T2,A,T79,A2)") "COLVAR| Use points", "No"
1144 : END IF
1145 : END IF
1146 :
1147 30 : END SUBROUTINE write_colvar_constraint
1148 :
1149 : ! **************************************************************************************************
1150 : !> \brief Write fix atom constraint information to output unit
1151 : !> \param fixd_constraint Data set of the fix atom constraint
1152 : !> \param ifixd Fix atom constraint/restraint number (index)
1153 : !> \param iw Logical unit number of the output unit
1154 : !> \author Matthias Krack (26.11.2025)
1155 : ! **************************************************************************************************
1156 2 : SUBROUTINE write_fixd_constraint(fixd_constraint, ifixd, iw)
1157 :
1158 : TYPE(fixd_constraint_type), INTENT(IN), POINTER :: fixd_constraint
1159 : INTEGER, INTENT(IN) :: ifixd, iw
1160 :
1161 2 : IF (iw > 0) THEN
1162 2 : CPASSERT(ASSOCIATED(fixd_constraint))
1163 2 : IF (fixd_constraint%restraint%active) THEN
1164 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1165 2 : "FIX_ATOM| Number (restraint)", ifixd
1166 : WRITE (UNIT=iw, FMT="(T2,A,T66,ES15.6)") &
1167 2 : "FIX_ATOM| Restraint constant k [a.u.]", fixd_constraint%restraint%k0
1168 : ELSE
1169 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1170 0 : "FIX_ATOM| Number (constraint)", ifixd
1171 : END IF
1172 : WRITE (UNIT=iw, FMT="(T2,A,T71,I10)") &
1173 2 : "FIX_ATOM| Atom index", fixd_constraint%fixd
1174 : WRITE (UNIT=iw, FMT="(T2,A,T78,A3)") &
1175 2 : "FIX_ATOM| Fixed Cartesian components", periodicity_string(fixd_constraint%itype)
1176 2 : IF (INDEX(periodicity_string(fixd_constraint%itype), "X") > 0) THEN
1177 : WRITE (UNIT=iw, FMT="(T2,A,T66,F15.8)") &
1178 2 : "FIX_ATOM| X coordinate [Angstrom]", cp_unit_from_cp2k(fixd_constraint%coord(1), "Angstrom")
1179 : END IF
1180 2 : IF (INDEX(periodicity_string(fixd_constraint%itype), "Y") > 0) THEN
1181 : WRITE (UNIT=iw, FMT="(T2,A,T66,F15.8)") &
1182 2 : "FIX_ATOM| Y coordinate [Angstrom]", cp_unit_from_cp2k(fixd_constraint%coord(2), "Angstrom")
1183 : END IF
1184 2 : IF (INDEX(periodicity_string(fixd_constraint%itype), "Z") > 0) THEN
1185 : WRITE (UNIT=iw, FMT="(T2,A,T66,F15.8)") &
1186 2 : "FIX_ATOM| Z coordinate [Angstrom]", cp_unit_from_cp2k(fixd_constraint%coord(3), "Angstrom")
1187 : END IF
1188 : END IF
1189 :
1190 2 : END SUBROUTINE write_fixd_constraint
1191 :
1192 : ! **************************************************************************************************
1193 : !> \brief Write G3x3 constraint information to output unit
1194 : !> \param g3x3_constraint Data set of the g3x3 constraint
1195 : !> \param ig3x3 G3x3 constraint/restraint number (index)
1196 : !> \param iw Logical unit number of the output unit
1197 : !> \author Matthias Krack (26.11.2025)
1198 : ! **************************************************************************************************
1199 2 : SUBROUTINE write_g3x3_constraint(g3x3_constraint, ig3x3, iw)
1200 :
1201 : TYPE(g3x3_constraint_type), INTENT(IN), POINTER :: g3x3_constraint
1202 : INTEGER, INTENT(IN) :: ig3x3, iw
1203 :
1204 2 : IF (iw > 0) THEN
1205 2 : CPASSERT(ASSOCIATED(g3x3_constraint))
1206 2 : IF (g3x3_constraint%restraint%active) THEN
1207 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1208 1 : "G3X3| Number (restraint)", ig3x3
1209 : WRITE (UNIT=iw, FMT="(T2,A,T66,ES15.6)") &
1210 1 : "G3X3| Restraint constant k [a.u.]", g3x3_constraint%restraint%k0
1211 : ELSE
1212 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1213 1 : "G3X3| Number (constraint)", ig3x3
1214 : END IF
1215 : WRITE (UNIT=iw, FMT="(T2,A,T71,I10)") &
1216 2 : "G3X3| Atom index a", g3x3_constraint%a, &
1217 2 : "G3X3| Atom index b", g3x3_constraint%b, &
1218 4 : "G3X3| Atom index c", g3x3_constraint%c
1219 : WRITE (UNIT=iw, FMT="(T2,A,T66,F15.8)") &
1220 2 : "G3X3| Distance (a,b) [Angstrom]", cp_unit_from_cp2k(g3x3_constraint%dab, "Angstrom"), &
1221 2 : "G3X3| Distance (a,c) [Angstrom]", cp_unit_from_cp2k(g3x3_constraint%dac, "Angstrom"), &
1222 4 : "G3X3| Distance (b,c) [Angstrom]", cp_unit_from_cp2k(g3x3_constraint%dbc, "Angstrom")
1223 : END IF
1224 :
1225 2 : END SUBROUTINE write_g3x3_constraint
1226 :
1227 : ! **************************************************************************************************
1228 : !> \brief Write G4x6 constraint information to output unit
1229 : !> \param g4x6_constraint Data set of the g4x6 constraint
1230 : !> \param ig4x6 G4x6 constraint/restraint number (index)
1231 : !> \param iw Logical unit number of the output unit
1232 : !> \author Matthias Krack (26.11.2025)
1233 : ! **************************************************************************************************
1234 2 : SUBROUTINE write_g4x6_constraint(g4x6_constraint, ig4x6, iw)
1235 :
1236 : TYPE(g4x6_constraint_type), INTENT(IN), POINTER :: g4x6_constraint
1237 : INTEGER, INTENT(IN) :: ig4x6, iw
1238 :
1239 2 : IF (iw > 0) THEN
1240 2 : CPASSERT(ASSOCIATED(g4x6_constraint))
1241 2 : IF (g4x6_constraint%restraint%active) THEN
1242 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1243 1 : "G4X6| Number (restraint)", ig4x6
1244 : WRITE (UNIT=iw, FMT="(T2,A,T66,ES15.6)") &
1245 1 : "G4X6| Restraint constant k [a.u.]", g4x6_constraint%restraint%k0
1246 : ELSE
1247 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1248 1 : "G4X6| Number (constraint)", ig4x6
1249 : END IF
1250 : WRITE (UNIT=iw, FMT="(T2,A,T71,I10)") &
1251 2 : "G4X6| Atom index a", g4x6_constraint%a, &
1252 2 : "G4X6| Atom index b", g4x6_constraint%b, &
1253 2 : "G4X6| Atom index c", g4x6_constraint%c, &
1254 4 : "G4X6| Atom index d", g4x6_constraint%d
1255 : WRITE (UNIT=iw, FMT="(T2,A,T66,F15.8)") &
1256 2 : "G4X6| Distance (a,b) [Angstrom]", cp_unit_from_cp2k(g4x6_constraint%dab, "Angstrom"), &
1257 2 : "G4X6| Distance (a,c) [Angstrom]", cp_unit_from_cp2k(g4x6_constraint%dac, "Angstrom"), &
1258 2 : "G4X6| Distance (a,d) [Angstrom]", cp_unit_from_cp2k(g4x6_constraint%dad, "Angstrom"), &
1259 2 : "G4X6| Distance (b,c) [Angstrom]", cp_unit_from_cp2k(g4x6_constraint%dbc, "Angstrom"), &
1260 2 : "G4X6| Distance (b,d) [Angstrom]", cp_unit_from_cp2k(g4x6_constraint%dbd, "Angstrom"), &
1261 4 : "G4X6| Distance (c,d) [Angstrom]", cp_unit_from_cp2k(g4x6_constraint%dcd, "Angstrom")
1262 : END IF
1263 :
1264 2 : END SUBROUTINE write_g4x6_constraint
1265 :
1266 : ! **************************************************************************************************
1267 : !> \brief Write virtual site constraint information to output unit
1268 : !> \param vsite_constraint Data set of the vsite constraint
1269 : !> \param ivsite Virtual site constraint/restraint number (index)
1270 : !> \param iw Logical unit number of the output unit
1271 : !> \author Matthias Krack (01.12.2025)
1272 : ! **************************************************************************************************
1273 0 : SUBROUTINE write_vsite_constraint(vsite_constraint, ivsite, iw)
1274 :
1275 : TYPE(vsite_constraint_type), INTENT(IN), POINTER :: vsite_constraint
1276 : INTEGER, INTENT(IN) :: ivsite, iw
1277 :
1278 0 : IF (iw > 0) THEN
1279 0 : CPASSERT(ASSOCIATED(vsite_constraint))
1280 0 : IF (vsite_constraint%restraint%active) THEN
1281 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1282 0 : "VSITE| Number (restraint)", ivsite
1283 : WRITE (UNIT=iw, FMT="(T2,A,T66,ES15.6)") &
1284 0 : "VSITE| Restraint constant k [a.u.]", vsite_constraint%restraint%k0
1285 : ELSE
1286 : WRITE (UNIT=iw, FMT="(/,T2,A,T71,I10)") &
1287 0 : "VSITE| Number (constraint)", ivsite
1288 : END IF
1289 : WRITE (UNIT=iw, FMT="(T2,A,T71,I10)") &
1290 0 : "VSITE| Atom index of virtual site", vsite_constraint%a, &
1291 0 : "VSITE| Atom index b", vsite_constraint%b, &
1292 0 : "VSITE| Atom index c", vsite_constraint%c, &
1293 0 : "VSITE| Atom index d", vsite_constraint%d
1294 : WRITE (UNIT=iw, FMT="(T2,A,T66,F15.8)") &
1295 0 : "VSITE| Distance (b,c) [Angstrom]", cp_unit_from_cp2k(vsite_constraint%wbc, "Angstrom"), &
1296 0 : "VSITE| Distance (d,c) [Angstrom]", cp_unit_from_cp2k(vsite_constraint%wdc, "Angstrom")
1297 : END IF
1298 :
1299 0 : END SUBROUTINE write_vsite_constraint
1300 :
1301 0 : END MODULE molecule_kind_types
|
