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