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