(git:374b731)
Loading...
Searching...
No Matches
topology_constraint_util.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Collection of subroutine needed for topology related things
10!> \par History
11!> jgh (23-05-2004) Last atom of molecule information added
12! **************************************************************************************************
13MODULE topology_constraint_util
14 USE atomic_kind_types, ONLY: atomic_kind_type,&
15 get_atomic_kind,&
16 is_hydrogen
17 USE cell_types, ONLY: use_perd_x,&
18 use_perd_xy,&
19 use_perd_xyz,&
20 use_perd_xz,&
21 use_perd_y,&
22 use_perd_yz,&
23 use_perd_z
24 USE colvar_methods, ONLY: colvar_eval_mol_f
25 USE colvar_types, ONLY: &
26 colvar_clone, colvar_counters, colvar_create, colvar_p_reallocate, colvar_release, &
27 colvar_setup, colvar_type, dist_colvar_id, torsion_colvar_id, xyz_diag_colvar_id, &
28 xyz_outerdiag_colvar_id
29 USE colvar_utils, ONLY: post_process_colvar
30 USE cp_log_handling, ONLY: cp_get_default_logger,&
31 cp_logger_type,&
32 cp_to_string
33 USE cp_output_handling, ONLY: cp_print_key_finished_output,&
34 cp_print_key_unit_nr
35 USE input_constants, ONLY: do_constr_atomic,&
36 do_constr_molec
37 USE input_section_types, ONLY: section_vals_get,&
38 section_vals_get_subs_vals,&
39 section_vals_type,&
40 section_vals_val_get,&
41 section_vals_val_set
42 USE kinds, ONLY: default_string_length,&
43 dp
44 USE memory_utilities, ONLY: reallocate
45 USE molecule_kind_types, ONLY: &
46 atom_type, bond_type, colvar_constraint_type, fixd_constraint_type, g3x3_constraint_type, &
47 g4x6_constraint_type, get_molecule_kind, molecule_kind_type, set_molecule_kind, &
48 setup_colvar_counters, vsite_constraint_type
49 USE molecule_types, ONLY: get_molecule,&
50 global_constraint_type,&
51 local_colvar_constraint_type,&
52 local_constraint_type,&
53 local_g3x3_constraint_type,&
54 local_g4x6_constraint_type,&
55 molecule_type,&
56 set_molecule
57 USE particle_types, ONLY: particle_type
58 USE qmmm_ff_fist, ONLY: qmmm_ff_precond_only_qm
59 USE qmmm_types_low, ONLY: qmmm_env_mm_type
60 USE topology_types, ONLY: constr_list_type,&
61 constraint_info_type,&
62 topology_parameters_type
63#include "./base/base_uses.f90"
64
65 IMPLICIT NONE
66
67 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'topology_constraint_util'
68
69 PRIVATE
70 PUBLIC :: topology_constraint_pack
71
72CONTAINS
73
74! **************************************************************************************************
75!> \brief Pack in all the information needed for the constraints
76!> \param molecule_kind_set ...
77!> \param molecule_set ...
78!> \param topology ...
79!> \param qmmm_env ...
80!> \param particle_set ...
81!> \param input_file ...
82!> \param subsys_section ...
83!> \param gci ...
84! **************************************************************************************************
85 SUBROUTINE topology_constraint_pack(molecule_kind_set, molecule_set, &
86 topology, qmmm_env, particle_set, input_file, subsys_section, gci)
87 TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
88 TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
89 TYPE(topology_parameters_type), INTENT(INOUT) :: topology
90 TYPE(qmmm_env_mm_type), OPTIONAL, POINTER :: qmmm_env
91 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
92 TYPE(section_vals_type), POINTER :: input_file, subsys_section
93 TYPE(global_constraint_type), POINTER :: gci
94
95 CHARACTER(len=*), PARAMETER :: routinen = 'topology_constraint_pack'
96
97 CHARACTER(LEN=2) :: element_symbol
98 CHARACTER(LEN=default_string_length) :: molname, name
99 CHARACTER(LEN=default_string_length), &
100 DIMENSION(:), POINTER :: atom_typeh, cnds
101 INTEGER :: cind, first, first_atom, gind, handle, handle2, i, ii, itype, iw, j, k, k1loc, &
102 k2loc, kk, last, last_atom, m, n_start_colv, natom, nbond, ncolv_glob, ncolv_mol, &
103 nfixd_list_gci, nfixd_restart, nfixd_restraint, nfixed_atoms, ng3x3, ng3x3_restraint, &
104 ng4x6, ng4x6_restraint, nhdist, nmolecule, nrep, nvsite, nvsite_restraint, offset
105 INTEGER, DIMENSION(:), POINTER :: constr_x_glob, inds, molecule_list
106 LOGICAL :: exclude_mm, exclude_qm, fix_atom_mm, fix_atom_molname, fix_atom_qm, &
107 fix_atom_qmmm, fix_fixed_atom, found_molname, is_qm, ishbond, ldummy, &
108 restart_restraint_clv, restart_restraint_pos, use_clv_info
109 LOGICAL, ALLOCATABLE, DIMENSION(:) :: missed_molname
110 REAL(kind=dp) :: rmod, rvec(3)
111 REAL(kind=dp), DIMENSION(:), POINTER :: hdist, r
112 TYPE(atom_type), DIMENSION(:), POINTER :: atom_list
113 TYPE(atomic_kind_type), POINTER :: atomic_kind
114 TYPE(bond_type), DIMENSION(:), POINTER :: bond_list
115 TYPE(colvar_constraint_type), DIMENSION(:), &
116 POINTER :: colv_list
117 TYPE(colvar_counters) :: ncolv
118 TYPE(constr_list_type), DIMENSION(:), POINTER :: constr_x_mol
119 TYPE(constraint_info_type), POINTER :: cons_info
120 TYPE(cp_logger_type), POINTER :: logger
121 TYPE(fixd_constraint_type), DIMENSION(:), POINTER :: fixd_list, fixd_list_gci
122 TYPE(g3x3_constraint_type), DIMENSION(:), POINTER :: g3x3_list
123 TYPE(g4x6_constraint_type), DIMENSION(:), POINTER :: g4x6_list
124 TYPE(local_colvar_constraint_type), DIMENSION(:), &
125 POINTER :: lcolv
126 TYPE(local_constraint_type), POINTER :: lci
127 TYPE(local_g3x3_constraint_type), DIMENSION(:), &
128 POINTER :: lg3x3
129 TYPE(local_g4x6_constraint_type), DIMENSION(:), &
130 POINTER :: lg4x6
131 TYPE(molecule_kind_type), POINTER :: molecule_kind
132 TYPE(molecule_type), POINTER :: molecule
133 TYPE(section_vals_type), POINTER :: colvar_func_info, colvar_rest, &
134 fixd_restr_rest, hbonds_section
135 TYPE(vsite_constraint_type), DIMENSION(:), POINTER :: vsite_list
136
137 NULLIFY (logger, constr_x_mol, constr_x_glob)
138 logger => cp_get_default_logger()
139 iw = cp_print_key_unit_nr(logger, subsys_section, "PRINT%TOPOLOGY_INFO/UTIL_INFO", &
140 extension=".subsysLog")
141 CALL timeset(routinen, handle)
142 CALL timeset(routinen//"_1", handle2)
143
144 cons_info => topology%cons_info
145 hbonds_section => section_vals_get_subs_vals(input_file, &
146 "MOTION%CONSTRAINT%HBONDS")
147 fixd_restr_rest => section_vals_get_subs_vals(input_file, &
148 "MOTION%CONSTRAINT%FIX_ATOM_RESTART")
149 CALL section_vals_get(fixd_restr_rest, explicit=restart_restraint_pos)
150 colvar_rest => section_vals_get_subs_vals(input_file, &
151 "MOTION%CONSTRAINT%COLVAR_RESTART")
152 CALL section_vals_get(colvar_rest, explicit=restart_restraint_clv)
153 colvar_func_info => section_vals_get_subs_vals(subsys_section, &
154 "COLVAR%COLVAR_FUNC_INFO")
155 CALL section_vals_get(colvar_func_info, explicit=use_clv_info)
156 !-----------------------------------------------------------------------------
157 !-----------------------------------------------------------------------------
158 ! 1. NULLIFY the molecule_set(imol)%lci via set_molecule_set
159 !-----------------------------------------------------------------------------
160 DO i = 1, topology%nmol
161 molecule => molecule_set(i)
162 NULLIFY (lci)
163 ! only allocate the lci if constraints are active. Can this stuff be distributed ?
164 IF (topology%const_atom .OR. topology%const_hydr .OR. &
165 topology%const_33 .OR. topology%const_46 .OR. &
166 topology%const_colv .OR. topology%const_vsite) THEN
167 ALLOCATE (lci)
168 NULLIFY (lci%lcolv)
169 NULLIFY (lci%lg3x3)
170 NULLIFY (lci%lg4x6)
171 END IF
172 CALL set_molecule(molecule, lci=lci)
173 END DO
174 ALLOCATE (gci)
175 NULLIFY (gci%lcolv, &
176 gci%lg3x3, &
177 gci%lg4x6, &
178 gci%fixd_list, &
179 gci%colv_list, &
180 gci%g3x3_list, &
181 gci%g4x6_list, &
182 gci%vsite_list)
183 gci%ntot = 0
184 gci%ng3x3 = 0
185 gci%ng4x6 = 0
186 gci%nvsite = 0
187 gci%ng3x3_restraint = 0
188 gci%ng4x6_restraint = 0
189 gci%nvsite_restraint = 0
190 CALL setup_colvar_counters(gci%colv_list, gci%ncolv)
191 gci%nrestraint = gci%ng3x3_restraint + &
192 gci%ng4x6_restraint + &
193 gci%nvsite_restraint + &
194 gci%ncolv%nrestraint
195 CALL timestop(handle2)
196 CALL timeset(routinen//"_2", handle2)
197 !-----------------------------------------------------------------------------
198 !-----------------------------------------------------------------------------
199 ! 2. Add more stuff to COLVAR constraint if constraint hydrogen is on
200 !-----------------------------------------------------------------------------
201 IF (topology%const_hydr) THEN
202 topology%const_colv = .true.
203 NULLIFY (atom_typeh, hdist)
204 ALLOCATE (constr_x_mol(SIZE(molecule_kind_set)))
205 DO i = 1, SIZE(molecule_kind_set)
206 ALLOCATE (constr_x_mol(i)%constr(1))
207 constr_x_mol(i)%constr(1) = 1
208 END DO
209 CALL section_vals_val_get(hbonds_section, "MOLECULE", n_rep_val=nrep)
210 IF (nrep /= 0) THEN
211 NULLIFY (inds)
212 DO i = 1, SIZE(molecule_kind_set)
213 constr_x_mol(i)%constr(1) = 0
214 END DO
215 CALL section_vals_val_get(hbonds_section, "MOLECULE", i_vals=inds)
216 DO i = 1, SIZE(inds)
217 constr_x_mol(inds(i))%constr(1) = 1
218 END DO
219 ELSE
220 CALL section_vals_val_get(hbonds_section, "MOLNAME", n_rep_val=nrep)
221 IF (nrep /= 0) THEN
222 NULLIFY (cnds)
223 DO i = 1, SIZE(molecule_kind_set)
224 constr_x_mol(i)%constr(1) = 0
225 END DO
226 CALL section_vals_val_get(hbonds_section, "MOLNAME", c_vals=cnds)
227 DO i = 1, SIZE(cnds)
228 found_molname = .false.
229 DO k = 1, SIZE(molecule_kind_set)
230 molecule_kind => molecule_kind_set(k)
231 name = molecule_kind%name
232 ldummy = qmmm_ff_precond_only_qm(id1=name)
233 IF (cnds(i) == name) THEN
234 constr_x_mol(k)%constr(1) = 1
235 found_molname = .true.
236 END IF
237 END DO
238 CALL print_warning_molname(found_molname, cnds(i))
239 END DO
240 END IF
241 END IF
242 CALL section_vals_val_get(hbonds_section, "ATOM_TYPE", n_rep_val=nrep)
243 IF (nrep /= 0) &
244 CALL section_vals_val_get(hbonds_section, "ATOM_TYPE", c_vals=atom_typeh)
245 CALL section_vals_val_get(hbonds_section, "TARGETS", n_rep_val=nrep)
246 IF (nrep /= 0) &
247 CALL section_vals_val_get(hbonds_section, "TARGETS", r_vals=hdist)
248 IF (ASSOCIATED(hdist)) THEN
249 cpassert(SIZE(hdist) == SIZE(atom_typeh))
250 END IF
251 CALL section_vals_val_get(hbonds_section, "exclude_qm", l_val=exclude_qm)
252 CALL section_vals_val_get(hbonds_section, "exclude_mm", l_val=exclude_mm)
253 nhdist = 0
254 DO i = 1, SIZE(molecule_kind_set)
255 molecule_kind => molecule_kind_set(i)
256 IF (constr_x_mol(i)%constr(1) == 0) cycle
257 CALL get_molecule_kind(molecule_kind=molecule_kind, &
258 bond_list=bond_list, nbond=nbond, atom_list=atom_list, &
259 molecule_list=molecule_list)
260 ! Let's tag all requested atoms involving Hydrogen
261 ! on the first molecule of this kind
262 molecule => molecule_set(molecule_list(1))
263 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
264 natom = last_atom - first_atom + 1
265 DO k = 1, nbond
266 ishbond = .false.
267 j = bond_list(k)%a
268 IF (j < 1 .OR. j > natom) cycle
269 atomic_kind => atom_list(j)%atomic_kind
270 CALL get_atomic_kind(atomic_kind=atomic_kind, name=name)
271 is_qm = qmmm_ff_precond_only_qm(id1=name)
272 IF ((name(1:1) == "H") .OR. is_hydrogen(atomic_kind)) ishbond = .true.
273 IF (is_qm .AND. exclude_qm) ishbond = .false.
274 IF (.NOT. (is_qm) .AND. exclude_mm) ishbond = .false.
275 IF (.NOT. ishbond) THEN
276 j = bond_list(k)%b
277 IF (j < 1 .OR. j > natom) cycle
278 atomic_kind => atom_list(j)%atomic_kind
279 CALL get_atomic_kind(atomic_kind=atomic_kind, name=name)
280 is_qm = qmmm_ff_precond_only_qm(id1=name)
281 IF ((name(1:1) == "H") .OR. is_hydrogen(atomic_kind)) ishbond = .true.
282 IF (is_qm .AND. exclude_qm) ishbond = .false.
283 IF (.NOT. (is_qm) .AND. exclude_mm) ishbond = .false.
284 END IF
285 IF (ishbond) THEN
286 nhdist = nhdist + 1
287 END IF
288 END DO
289 END DO
290 n_start_colv = cons_info%nconst_colv
291 cons_info%nconst_colv = nhdist + n_start_colv
292 CALL reallocate(cons_info%const_colv_mol, 1, cons_info%nconst_colv)
293 CALL reallocate(cons_info%const_colv_molname, 1, cons_info%nconst_colv)
294 CALL reallocate(cons_info%const_colv_target, 1, cons_info%nconst_colv)
295 CALL reallocate(cons_info%const_colv_target_growth, 1, cons_info%nconst_colv)
296 CALL colvar_p_reallocate(cons_info%colvar_set, 1, cons_info%nconst_colv)
297 ! Fill in Restraints info
298 CALL reallocate(cons_info%colv_intermolecular, 1, cons_info%nconst_colv)
299 CALL reallocate(cons_info%colv_restraint, 1, cons_info%nconst_colv)
300 CALL reallocate(cons_info%colv_k0, 1, cons_info%nconst_colv)
301 CALL reallocate(cons_info%colv_exclude_qm, 1, cons_info%nconst_colv)
302 CALL reallocate(cons_info%colv_exclude_mm, 1, cons_info%nconst_colv)
303 ! Bonds involving hydrogens are by their nature only intramolecular
304 cons_info%colv_intermolecular(n_start_colv + 1:cons_info%nconst_colv) = .false.
305 cons_info%colv_exclude_qm(n_start_colv + 1:cons_info%nconst_colv) = .false.
306 cons_info%colv_exclude_mm(n_start_colv + 1:cons_info%nconst_colv) = .false.
307 cons_info%colv_restraint(n_start_colv + 1:cons_info%nconst_colv) = cons_info%hbonds_restraint
308 cons_info%colv_k0(n_start_colv + 1:cons_info%nconst_colv) = cons_info%hbonds_k0
309 !
310 nhdist = 0
311 DO i = 1, SIZE(molecule_kind_set)
312 IF (constr_x_mol(i)%constr(1) == 0) cycle
313 molecule_kind => molecule_kind_set(i)
314 CALL get_molecule_kind(molecule_kind=molecule_kind, &
315 bond_list=bond_list, nbond=nbond, atom_list=atom_list, &
316 molecule_list=molecule_list)
317 molecule => molecule_set(molecule_list(1))
318 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
319 natom = last_atom - first_atom + 1
320 offset = first_atom - 1
321 DO k = 1, nbond
322 ishbond = .false.
323 j = bond_list(k)%a
324 IF (j < 1 .OR. j > natom) cycle
325 atomic_kind => atom_list(j)%atomic_kind
326 CALL get_atomic_kind(atomic_kind=atomic_kind, name=name)
327 is_qm = qmmm_ff_precond_only_qm(id1=name)
328 IF ((name(1:1) == "H") .OR. is_hydrogen(atomic_kind)) ishbond = .true.
329 IF (is_qm .AND. exclude_qm) ishbond = .false.
330 IF (.NOT. (is_qm) .AND. exclude_mm) ishbond = .false.
331 IF (.NOT. ishbond) THEN
332 j = bond_list(k)%b
333 IF (j < 1 .OR. j > natom) cycle
334 atomic_kind => atom_list(j)%atomic_kind
335 CALL get_atomic_kind(atomic_kind=atomic_kind, name=name)
336 is_qm = qmmm_ff_precond_only_qm(id1=name)
337 IF ((name(1:1) == "H") .OR. is_hydrogen(atomic_kind)) ishbond = .true.
338 IF (is_qm .AND. exclude_qm) ishbond = .false.
339 IF (.NOT. (is_qm) .AND. exclude_mm) ishbond = .false.
340 END IF
341 IF (ishbond) THEN
342 nhdist = nhdist + 1
343 rvec = particle_set(offset + bond_list(k)%a)%r - particle_set(offset + bond_list(k)%b)%r
344 rmod = sqrt(dot_product(rvec, rvec))
345 IF (ASSOCIATED(hdist)) THEN
346 IF (SIZE(hdist) > 0) THEN
347 IF (bond_list(k)%a == j) atomic_kind => atom_list(bond_list(k)%b)%atomic_kind
348 IF (bond_list(k)%b == j) atomic_kind => atom_list(bond_list(k)%a)%atomic_kind
349 CALL get_atomic_kind(atomic_kind=atomic_kind, &
350 name=name, element_symbol=element_symbol)
351 ldummy = qmmm_ff_precond_only_qm(id1=name)
352 DO m = 1, SIZE(hdist)
353 IF (trim(name) == trim(atom_typeh(m))) EXIT
354 IF (trim(element_symbol) == trim(atom_typeh(m))) EXIT
355 END DO
356 IF (m <= SIZE(hdist)) THEN
357 rmod = hdist(m)
358 END IF
359 END IF
360 END IF
361 cons_info%const_colv_mol(nhdist + n_start_colv) = i
362 cons_info%const_colv_molname(nhdist + n_start_colv) = "UNDEF"
363 cons_info%const_colv_target(nhdist + n_start_colv) = rmod
364 cons_info%const_colv_target_growth(nhdist + n_start_colv) = 0.0_dp
365 CALL colvar_create(cons_info%colvar_set(nhdist + n_start_colv)%colvar, &
366 dist_colvar_id)
367 cons_info%colvar_set(nhdist + n_start_colv)%colvar%dist_param%i_at = bond_list(k)%a
368 cons_info%colvar_set(nhdist + n_start_colv)%colvar%dist_param%j_at = bond_list(k)%b
369 CALL colvar_setup(cons_info%colvar_set(nhdist + n_start_colv)%colvar)
370 END IF
371 END DO
372 END DO
373 DO j = 1, SIZE(constr_x_mol)
374 DEALLOCATE (constr_x_mol(j)%constr)
375 END DO
376 DEALLOCATE (constr_x_mol)
377 END IF
378
379 CALL timestop(handle2)
380 CALL timeset(routinen//"_3", handle2)
381 !-----------------------------------------------------------------------------
382 !-----------------------------------------------------------------------------
383 ! 3. Set the COLVAR constraint molecule_kind_set(ikind)%colv_list
384 !-----------------------------------------------------------------------------
385 IF (topology%const_colv) THEN
386 ! Post Process of COLVARS..
387 DO ii = 1, SIZE(cons_info%colvar_set)
388 CALL post_process_colvar(cons_info%colvar_set(ii)%colvar, particle_set)
389 END DO
390 ! Real constraint/restraint part..
391 CALL give_constraint_array(cons_info%const_colv_mol, &
392 cons_info%const_colv_molname, &
393 cons_info%colv_intermolecular, &
394 constr_x_mol, &
395 constr_x_glob, &
396 molecule_kind_set, &
397 cons_info%colv_exclude_qm, &
398 cons_info%colv_exclude_mm)
399 ! Intramolecular constraints
400 gind = 0
401 cind = 0
402 DO ii = 1, SIZE(molecule_kind_set)
403 molecule_kind => molecule_kind_set(ii)
404 CALL get_molecule_kind(molecule_kind=molecule_kind, &
405 nmolecule=nmolecule, molecule_list=molecule_list)
406 ncolv_mol = SIZE(constr_x_mol(ii)%constr)
407 ALLOCATE (colv_list(ncolv_mol))
408 ! Starting index of the first molecule of this kind.
409 ! We need the index if no target is provided in the input file
410 ! for the collective variable.. The target will be computed on the
411 ! first molecule of the kind...
412 molecule => molecule_set(molecule_list(1))
413 CALL get_molecule(molecule, first_atom=first_atom)
414 CALL setup_colv_list(colv_list, constr_x_mol(ii)%constr, gind, &
415 cons_info, topology, particle_set, restart_restraint_clv, &
416 colvar_rest, first_atom)
417 CALL setup_colvar_counters(colv_list, ncolv)
418 CALL set_molecule_kind(molecule_kind, colv_list=colv_list, ncolv=ncolv)
419 DO j = 1, nmolecule
420 molecule => molecule_set(molecule_list(j))
421 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
422 ALLOCATE (lcolv(ncolv_mol))
423 CALL setup_lcolv(lcolv, constr_x_mol(ii)%constr, first_atom, last_atom, &
424 cons_info, particle_set, colvar_func_info, use_clv_info, cind)
425 CALL set_molecule(molecule=molecule, lcolv=lcolv)
426 END DO
427 END DO
428 DO j = 1, SIZE(constr_x_mol)
429 DEALLOCATE (constr_x_mol(j)%constr)
430 END DO
431 DEALLOCATE (constr_x_mol)
432 ! Intermolecular constraints
433 ncolv_glob = 0
434 IF (ASSOCIATED(constr_x_glob)) THEN
435 ncolv_glob = SIZE(constr_x_glob)
436 ALLOCATE (colv_list(ncolv_glob))
437 CALL setup_colv_list(colv_list, constr_x_glob, gind, cons_info, &
438 topology, particle_set, restart_restraint_clv, colvar_rest, &
439 first_atom=1)
440 CALL setup_colvar_counters(colv_list, ncolv)
441 ALLOCATE (lcolv(ncolv_glob))
442 CALL setup_lcolv(lcolv, constr_x_glob, 1, SIZE(particle_set), cons_info, &
443 particle_set, colvar_func_info, use_clv_info, cind)
444 gci%colv_list => colv_list
445 gci%lcolv => lcolv
446 gci%ncolv = ncolv
447 ! Total number of Intermolecular constraints
448 gci%ntot = gci%ncolv%ntot + gci%ntot
449 DEALLOCATE (constr_x_glob)
450 END IF
451 END IF
452
453 CALL timestop(handle2)
454 CALL timeset(routinen//"_4", handle2)
455 !-----------------------------------------------------------------------------
456 !-----------------------------------------------------------------------------
457 ! 4. Set the group 3x3 constraint g3x3_list
458 !-----------------------------------------------------------------------------
459 IF (topology%const_33) THEN
460 CALL give_constraint_array(cons_info%const_g33_mol, &
461 cons_info%const_g33_molname, &
462 cons_info%g33_intermolecular, &
463 constr_x_mol, &
464 constr_x_glob, &
465 molecule_kind_set, &
466 cons_info%g33_exclude_qm, &
467 cons_info%g33_exclude_mm)
468 ! Intramolecular constraints
469 DO ii = 1, SIZE(molecule_kind_set)
470 molecule_kind => molecule_kind_set(ii)
471 CALL get_molecule_kind(molecule_kind=molecule_kind, &
472 nmolecule=nmolecule, &
473 molecule_list=molecule_list)
474 ng3x3 = SIZE(constr_x_mol(ii)%constr)
475 ALLOCATE (g3x3_list(ng3x3))
476 CALL setup_g3x3_list(g3x3_list, constr_x_mol(ii)%constr, cons_info, ng3x3_restraint)
477 CALL set_molecule_kind(molecule_kind, ng3x3=ng3x3, ng3x3_restraint=ng3x3_restraint, g3x3_list=g3x3_list)
478 DO j = 1, nmolecule
479 molecule => molecule_set(molecule_list(j))
480 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
481 ALLOCATE (lg3x3(ng3x3))
482 CALL setup_lg3x3(lg3x3, g3x3_list, first_atom, last_atom)
483 CALL set_molecule(molecule=molecule, lg3x3=lg3x3)
484 END DO
485 END DO
486 DO j = 1, SIZE(constr_x_mol)
487 DEALLOCATE (constr_x_mol(j)%constr)
488 END DO
489 DEALLOCATE (constr_x_mol)
490 ! Intermolecular constraints
491 IF (ASSOCIATED(constr_x_glob)) THEN
492 ng3x3 = SIZE(constr_x_glob)
493 ALLOCATE (g3x3_list(ng3x3))
494 CALL setup_g3x3_list(g3x3_list, constr_x_glob, cons_info, ng3x3_restraint)
495 ALLOCATE (lg3x3(ng3x3))
496 CALL setup_lg3x3(lg3x3, g3x3_list, first_atom, last_atom)
497 gci%g3x3_list => g3x3_list
498 gci%lg3x3 => lg3x3
499 gci%ng3x3 = ng3x3
500 gci%ng3x3_restraint = ng3x3_restraint
501 ! Total number of Intermolecular constraints
502 gci%ntot = 3*gci%ng3x3 + gci%ntot
503 DEALLOCATE (constr_x_glob)
504 END IF
505 END IF
506
507 CALL timestop(handle2)
508 CALL timeset(routinen//"_5", handle2)
509 !-----------------------------------------------------------------------------
510 !-----------------------------------------------------------------------------
511 ! 5. Set the group 4x6 constraint g4x6_list
512 !-----------------------------------------------------------------------------
513 IF (topology%const_46) THEN
514 CALL give_constraint_array(cons_info%const_g46_mol, &
515 cons_info%const_g46_molname, &
516 cons_info%g46_intermolecular, &
517 constr_x_mol, &
518 constr_x_glob, &
519 molecule_kind_set, &
520 cons_info%g46_exclude_qm, &
521 cons_info%g46_exclude_mm)
522 ! Intramolecular constraints
523 DO ii = 1, SIZE(molecule_kind_set)
524 molecule_kind => molecule_kind_set(ii)
525 CALL get_molecule_kind(molecule_kind=molecule_kind, &
526 nmolecule=nmolecule, molecule_list=molecule_list)
527 ng4x6 = SIZE(constr_x_mol(ii)%constr)
528 ALLOCATE (g4x6_list(ng4x6))
529 CALL setup_g4x6_list(g4x6_list, constr_x_mol(ii)%constr, cons_info, ng4x6_restraint)
530 CALL set_molecule_kind(molecule_kind, ng4x6=ng4x6, ng4x6_restraint=ng4x6_restraint, g4x6_list=g4x6_list)
531 DO j = 1, nmolecule
532 molecule => molecule_set(molecule_list(j))
533 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
534 ALLOCATE (lg4x6(ng4x6))
535 CALL setup_lg4x6(lg4x6, g4x6_list, first_atom, last_atom)
536 CALL set_molecule(molecule=molecule, lg4x6=lg4x6)
537 END DO
538 END DO
539 DO j = 1, SIZE(constr_x_mol)
540 DEALLOCATE (constr_x_mol(j)%constr)
541 END DO
542 DEALLOCATE (constr_x_mol)
543 ! Intermolecular constraints
544 IF (ASSOCIATED(constr_x_glob)) THEN
545 ng4x6 = SIZE(constr_x_glob)
546 ALLOCATE (g4x6_list(ng4x6))
547 CALL setup_g4x6_list(g4x6_list, constr_x_glob, cons_info, ng4x6_restraint)
548 ALLOCATE (lg4x6(ng4x6))
549 CALL setup_lg4x6(lg4x6, g4x6_list, first_atom, last_atom)
550 gci%g4x6_list => g4x6_list
551 gci%lg4x6 => lg4x6
552 gci%ng4x6 = ng4x6
553 gci%ng4x6_restraint = ng4x6_restraint
554 ! Total number of Intermolecular constraints
555 gci%ntot = 6*gci%ng4x6 + gci%ntot
556 DEALLOCATE (constr_x_glob)
557 END IF
558 END IF
559
560 CALL timestop(handle2)
561 CALL timeset(routinen//"_6", handle2)
562 !-----------------------------------------------------------------------------
563 !-----------------------------------------------------------------------------
564 ! 6. Set the group vsite constraint vsite_list
565 !-----------------------------------------------------------------------------
566 IF (topology%const_vsite) THEN
567 CALL give_constraint_array(cons_info%const_vsite_mol, &
568 cons_info%const_vsite_molname, &
569 cons_info%vsite_intermolecular, &
570 constr_x_mol, &
571 constr_x_glob, &
572 molecule_kind_set, &
573 cons_info%vsite_exclude_qm, &
574 cons_info%vsite_exclude_mm)
575 ! Intramolecular constraints
576 DO ii = 1, SIZE(molecule_kind_set)
577 molecule_kind => molecule_kind_set(ii)
578 CALL get_molecule_kind(molecule_kind=molecule_kind, &
579 nmolecule=nmolecule, molecule_list=molecule_list)
580 nvsite = SIZE(constr_x_mol(ii)%constr)
581 ALLOCATE (vsite_list(nvsite))
582 CALL setup_vsite_list(vsite_list, constr_x_mol(ii)%constr, cons_info, nvsite_restraint)
583 CALL set_molecule_kind(molecule_kind, nvsite=nvsite, nvsite_restraint=nvsite_restraint, &
584 vsite_list=vsite_list)
585 END DO
586 DO j = 1, SIZE(constr_x_mol)
587 DEALLOCATE (constr_x_mol(j)%constr)
588 END DO
589 DEALLOCATE (constr_x_mol)
590 ! Intermolecular constraints
591 IF (ASSOCIATED(constr_x_glob)) THEN
592 nvsite = SIZE(constr_x_glob)
593 ALLOCATE (vsite_list(nvsite))
594 CALL setup_vsite_list(vsite_list, constr_x_glob, cons_info, nvsite_restraint)
595 gci%vsite_list => vsite_list
596 gci%nvsite = nvsite
597 gci%nvsite_restraint = nvsite_restraint
598 ! Total number of Intermolecular constraints
599 gci%ntot = gci%nvsite + gci%ntot
600 DEALLOCATE (constr_x_glob)
601 END IF
602 END IF
603 CALL timestop(handle2)
604 CALL timeset(routinen//"_7", handle2)
605 !-----------------------------------------------------------------------------
606 !-----------------------------------------------------------------------------
607 ! 7. Set the group fixed_atom constraint fixd_list
608 !-----------------------------------------------------------------------------
609 IF (topology%const_atom) THEN
610 ALLOCATE (fixd_list_gci(SIZE(particle_set)))
611 nfixd_list_gci = 0
612 ALLOCATE (missed_molname(SIZE(cons_info%fixed_molnames, 1)))
613 missed_molname = .true.
614 nfixd_restart = 0
615 DO i = 1, SIZE(molecule_kind_set)
616 molecule_kind => molecule_kind_set(i)
617 CALL get_molecule_kind(molecule_kind=molecule_kind, &
618 nmolecule=nmolecule, molecule_list=molecule_list, name=molname)
619 is_qm = qmmm_ff_precond_only_qm(id1=molname)
620 WHERE (molname .EQ. cons_info%fixed_molnames)
621 missed_molname = .false.
622 END WHERE
623 ! Try to figure out how many atoms of the list belong to this molecule_kind
624 nfixed_atoms = 0
625 DO j = 1, nmolecule
626 molecule => molecule_set(molecule_list(j))
627 CALL get_molecule(molecule, first_atom=first, last_atom=last)
628 fix_atom_molname = .false.
629 IF (ASSOCIATED(cons_info%fixed_molnames)) THEN
630 DO k = 1, SIZE(cons_info%fixed_molnames)
631 IF (cons_info%fixed_molnames(k) .EQ. molname) THEN
632 fix_atom_molname = .true.
633 IF (is_qm .AND. cons_info%fixed_exclude_qm(k)) fix_atom_molname = .false.
634 IF ((.NOT. is_qm) .AND. cons_info%fixed_exclude_mm(k)) fix_atom_molname = .false.
635 END IF
636 END DO
637 END IF
638 DO k = first, last
639 fix_atom_qmmm = .false.
640 IF (PRESENT(qmmm_env)) THEN
641 SELECT CASE (cons_info%freeze_qm)
642 CASE (do_constr_atomic)
643 IF (any(qmmm_env%qm_atom_index == k)) fix_atom_qmmm = .true.
644 CASE (do_constr_molec)
645 IF (any(qmmm_env%qm_molecule_index == molecule_list(j))) fix_atom_qmmm = .true.
646 END SELECT
647 SELECT CASE (cons_info%freeze_mm)
648 CASE (do_constr_atomic)
649 IF (all(qmmm_env%qm_atom_index /= k)) fix_atom_qmmm = .true.
650 CASE (do_constr_molec)
651 IF (all(qmmm_env%qm_molecule_index /= molecule_list(j))) fix_atom_qmmm = .true.
652 END SELECT
653 END IF
654 IF (any(cons_info%fixed_atoms == k) .OR. fix_atom_qmmm .OR. fix_atom_molname) THEN
655 nfixed_atoms = nfixed_atoms + 1
656 END IF
657 END DO
658 END DO
659 ALLOCATE (fixd_list(nfixed_atoms))
660 kk = 0
661 nfixd_restraint = 0
662 IF (nfixed_atoms /= 0) THEN
663 DO j = 1, nmolecule
664 molecule => molecule_set(molecule_list(j))
665 CALL get_molecule(molecule, first_atom=first, last_atom=last)
666 fix_atom_molname = .false.
667 IF (ASSOCIATED(cons_info%fixed_molnames)) THEN
668 DO k1loc = 1, SIZE(cons_info%fixed_molnames)
669 IF (cons_info%fixed_molnames(k1loc) .EQ. molname) THEN
670 fix_atom_molname = .true.
671 itype = cons_info%fixed_mol_type(k1loc)
672 EXIT
673 END IF
674 END DO
675 END IF
676 DO k = first, last
677 ! FIXED LIST ATOMS
678 fix_fixed_atom = .false.
679 DO k2loc = 1, SIZE(cons_info%fixed_atoms)
680 IF (cons_info%fixed_atoms(k2loc) == k) THEN
681 fix_fixed_atom = .true.
682 itype = cons_info%fixed_type(k2loc)
683 EXIT
684 END IF
685 END DO
686 ! QMMM FIXED ATOMS (QM OR MM)
687 fix_atom_qmmm = .false.
688 fix_atom_mm = .false.
689 fix_atom_qm = .false.
690 IF (PRESENT(qmmm_env)) THEN
691 SELECT CASE (cons_info%freeze_qm)
692 CASE (do_constr_atomic)
693 IF (any(qmmm_env%qm_atom_index == k)) THEN
694 fix_atom_qmmm = .true.
695 fix_atom_qm = .true.
696 itype = cons_info%freeze_qm_type
697 END IF
698 CASE (do_constr_molec)
699 IF (any(qmmm_env%qm_molecule_index == molecule_list(j))) THEN
700 fix_atom_qmmm = .true.
701 fix_atom_qm = .true.
702 itype = cons_info%freeze_qm_type
703 END IF
704 END SELECT
705 SELECT CASE (cons_info%freeze_mm)
706 CASE (do_constr_atomic)
707 IF (all(qmmm_env%qm_atom_index /= k)) THEN
708 fix_atom_qmmm = .true.
709 fix_atom_mm = .true.
710 itype = cons_info%freeze_mm_type
711 END IF
712 CASE (do_constr_molec)
713 IF (all(qmmm_env%qm_molecule_index /= molecule_list(j))) THEN
714 fix_atom_qmmm = .true.
715 fix_atom_mm = .true.
716 itype = cons_info%freeze_mm_type
717 END IF
718 END SELECT
719 ! We should never reach this point but let's check it anyway
720 IF (fix_atom_qm .AND. fix_atom_mm) THEN
721 CALL cp_abort(__location__, &
722 "Atom number: "//cp_to_string(k)// &
723 " has been defined both QM and MM. General Error!")
724 END IF
725 END IF
726 ! Check that the fixed atom constraint/restraint is unique
727 IF ((fix_fixed_atom .AND. fix_atom_qmmm) .OR. (fix_fixed_atom .AND. fix_atom_molname) &
728 .OR. (fix_atom_qmmm .AND. fix_atom_molname)) THEN
729 CALL cp_abort(__location__, &
730 "Atom number: "//cp_to_string(k)// &
731 " has been constrained/restrained to be fixed in more than one"// &
732 " input section. Check and correct your input file!")
733 END IF
734 ! Let's store the atom index
735 IF (fix_fixed_atom .OR. fix_atom_qmmm .OR. fix_atom_molname) THEN
736 kk = kk + 1
737 fixd_list(kk)%fixd = k
738 fixd_list(kk)%coord = particle_set(k)%r
739 fixd_list(kk)%itype = itype
740 ! Possibly Restraint
741 IF (fix_fixed_atom) THEN
742 fixd_list(kk)%restraint%active = cons_info%fixed_restraint(k2loc)
743 fixd_list(kk)%restraint%k0 = cons_info%fixed_k0(k2loc)
744 ELSEIF (fix_atom_qm) THEN
745 fixd_list(kk)%restraint%active = cons_info%fixed_qm_restraint
746 fixd_list(kk)%restraint%k0 = cons_info%fixed_qm_k0
747 ELSEIF (fix_atom_mm) THEN
748 fixd_list(kk)%restraint%active = cons_info%fixed_mm_restraint
749 fixd_list(kk)%restraint%k0 = cons_info%fixed_mm_k0
750 ELSEIF (fix_atom_molname) THEN
751 fixd_list(kk)%restraint%active = cons_info%fixed_mol_restraint(k1loc)
752 fixd_list(kk)%restraint%k0 = cons_info%fixed_mol_k0(k1loc)
753 ELSE
754 ! Should never reach this point
755 cpabort("")
756 END IF
757 IF (fixd_list(kk)%restraint%active) THEN
758 nfixd_restraint = nfixd_restraint + 1
759 nfixd_restart = nfixd_restart + 1
760 ! Check that we use the components that we really want..
761 SELECT CASE (itype)
762 CASE (use_perd_x)
763 fixd_list(kk)%coord(2) = huge(0.0_dp)
764 fixd_list(kk)%coord(3) = huge(0.0_dp)
765 CASE (use_perd_y)
766 fixd_list(kk)%coord(1) = huge(0.0_dp)
767 fixd_list(kk)%coord(3) = huge(0.0_dp)
768 CASE (use_perd_z)
769 fixd_list(kk)%coord(1) = huge(0.0_dp)
770 fixd_list(kk)%coord(2) = huge(0.0_dp)
771 CASE (use_perd_xy)
772 fixd_list(kk)%coord(3) = huge(0.0_dp)
773 CASE (use_perd_xz)
774 fixd_list(kk)%coord(2) = huge(0.0_dp)
775 CASE (use_perd_yz)
776 fixd_list(kk)%coord(1) = huge(0.0_dp)
777 END SELECT
778 IF (restart_restraint_pos) THEN
779 ! Read coord0 value for restraint
780 CALL section_vals_val_get(fixd_restr_rest, "_DEFAULT_KEYWORD_", &
781 i_rep_val=nfixd_restart, r_vals=r)
782 SELECT CASE (itype)
783 CASE (use_perd_x)
784 cpassert(SIZE(r) == 1)
785 fixd_list(kk)%coord(1) = r(1)
786 CASE (use_perd_y)
787 cpassert(SIZE(r) == 1)
788 fixd_list(kk)%coord(2) = r(1)
789 CASE (use_perd_z)
790 cpassert(SIZE(r) == 1)
791 fixd_list(kk)%coord(3) = r(1)
792 CASE (use_perd_xy)
793 cpassert(SIZE(r) == 2)
794 fixd_list(kk)%coord(1) = r(1)
795 fixd_list(kk)%coord(2) = r(2)
796 CASE (use_perd_xz)
797 cpassert(SIZE(r) == 2)
798 fixd_list(kk)%coord(1) = r(1)
799 fixd_list(kk)%coord(3) = r(2)
800 CASE (use_perd_yz)
801 cpassert(SIZE(r) == 2)
802 fixd_list(kk)%coord(2) = r(1)
803 fixd_list(kk)%coord(3) = r(2)
804 CASE (use_perd_xyz)
805 cpassert(SIZE(r) == 3)
806 fixd_list(kk)%coord(1:3) = r(1:3)
807 END SELECT
808 ELSE
809 ! Write coord0 value for restraint
810 SELECT CASE (itype)
811 CASE (use_perd_x)
812 ALLOCATE (r(1))
813 r(1) = fixd_list(kk)%coord(1)
814 CASE (use_perd_y)
815 ALLOCATE (r(1))
816 r(1) = fixd_list(kk)%coord(2)
817 CASE (use_perd_z)
818 ALLOCATE (r(1))
819 r(1) = fixd_list(kk)%coord(3)
820 CASE (use_perd_xy)
821 ALLOCATE (r(2))
822 r(1) = fixd_list(kk)%coord(1)
823 r(2) = fixd_list(kk)%coord(2)
824 CASE (use_perd_xz)
825 ALLOCATE (r(2))
826 r(1) = fixd_list(kk)%coord(1)
827 r(2) = fixd_list(kk)%coord(3)
828 CASE (use_perd_yz)
829 ALLOCATE (r(2))
830 r(1) = fixd_list(kk)%coord(1)
831 r(2) = fixd_list(kk)%coord(3)
832 CASE (use_perd_xyz)
833 ALLOCATE (r(3))
834 r(1:3) = fixd_list(kk)%coord(1:3)
835 END SELECT
836 CALL section_vals_val_set(fixd_restr_rest, "_DEFAULT_KEYWORD_", &
837 i_rep_val=nfixd_restart, r_vals_ptr=r)
838 END IF
839 END IF
840 END IF
841 END DO
842 END DO
843 END IF
844 IF (iw > 0) THEN
845 WRITE (iw, *) "MOLECULE KIND:", i, " NR. FIXED ATOMS:", SIZE(fixd_list(:)%fixd), " LIST::", fixd_list(:)%fixd
846 END IF
847 CALL set_molecule_kind(molecule_kind, nfixd=nfixed_atoms, nfixd_restraint=nfixd_restraint, &
848 fixd_list=fixd_list)
849 fixd_list_gci(nfixd_list_gci + 1:nfixd_list_gci + nfixed_atoms) = fixd_list
850 nfixd_list_gci = nfixd_list_gci + nfixed_atoms
851 END DO
852 IF (iw > 0) THEN
853 WRITE (iw, *) "TOTAL NUMBER OF FIXED ATOMS:", nfixd_list_gci
854 END IF
855 cpassert(count(missed_molname) == 0)
856 DEALLOCATE (missed_molname)
857 ! Intermolecular constraints
858 IF (gci%ntot /= 0) THEN
859 ALLOCATE (fixd_list(nfixd_list_gci))
860 fixd_list(1:nfixd_list_gci) = fixd_list_gci(1:nfixd_list_gci)
861 gci%fixd_list => fixd_list
862 END IF
863 DEALLOCATE (fixd_list_gci)
864 END IF
865 ! Final setup of the number of possible restraints
866 gci%nrestraint = gci%ng3x3_restraint + &
867 gci%ng4x6_restraint + &
868 gci%nvsite_restraint + &
869 gci%ncolv%nrestraint
870 CALL cp_print_key_finished_output(iw, logger, subsys_section, &
871 "PRINT%TOPOLOGY_INFO/UTIL_INFO")
872 CALL timestop(handle2)
873 CALL timestop(handle)
874 END SUBROUTINE topology_constraint_pack
875
876! **************************************************************************************************
877!> \brief Setup the colv_list for the packing of constraints
878!> \param colv_list ...
879!> \param ilist ...
880!> \param gind ...
881!> \param cons_info ...
882!> \param topology ...
883!> \param particle_set ...
884!> \param restart_restraint_clv ...
885!> \param colvar_rest ...
886!> \param first_atom ...
887!> \par History
888!> Updated 2007 for intermolecular constraints
889!> \author Teodoro Laino [2007]
890! **************************************************************************************************
891 SUBROUTINE setup_colv_list(colv_list, ilist, gind, cons_info, topology, &
892 particle_set, restart_restraint_clv, colvar_rest, first_atom)
893
894 TYPE(colvar_constraint_type), DIMENSION(:), &
895 POINTER :: colv_list
896 INTEGER, DIMENSION(:), POINTER :: ilist
897 INTEGER, INTENT(INOUT) :: gind
898 TYPE(constraint_info_type), POINTER :: cons_info
899 TYPE(topology_parameters_type), INTENT(INOUT) :: topology
900 TYPE(particle_type), DIMENSION(:), OPTIONAL, &
901 POINTER :: particle_set
902 LOGICAL, INTENT(IN) :: restart_restraint_clv
903 TYPE(section_vals_type), POINTER :: colvar_rest
904 INTEGER, INTENT(IN) :: first_atom
905
906 INTEGER :: j, kdim, kk, ncolv_mol
907 REAL(kind=dp) :: rmod
908 TYPE(colvar_type), POINTER :: local_colvar
909
910 ncolv_mol = 0
911 DO kk = 1, SIZE(ilist)
912 j = ilist(kk)
913 ncolv_mol = ncolv_mol + 1
914 kdim = SIZE(cons_info%colvar_set(j)%colvar%i_atom)
915 ALLOCATE (colv_list(ncolv_mol)%i_atoms(kdim))
916 colv_list(ncolv_mol)%inp_seq_num = j
917 colv_list(ncolv_mol)%type_id = cons_info%colvar_set(j)%colvar%type_id
918 colv_list(ncolv_mol)%i_atoms = cons_info%colvar_set(j)%colvar%i_atom
919 colv_list(ncolv_mol)%use_points = cons_info%colvar_set(j)%colvar%use_points
920 ! Restraint
921 colv_list(ncolv_mol)%restraint%active = cons_info%colv_restraint(j)
922 colv_list(ncolv_mol)%restraint%k0 = cons_info%colv_k0(j)
923 IF (cons_info%const_colv_target(j) == -huge(0.0_dp)) THEN
924 ! Let's compute the value..
925 NULLIFY (local_colvar)
926 CALL colvar_clone(local_colvar, cons_info%colvar_set(j)%colvar, &
927 i_atom_offset=first_atom - 1)
928 CALL colvar_eval_mol_f(local_colvar, topology%cell, particle_set)
929 colv_list(ncolv_mol)%expected_value = local_colvar%ss
930 CALL colvar_release(local_colvar)
931 ELSE
932 colv_list(ncolv_mol)%expected_value = cons_info%const_colv_target(j)
933 END IF
934 colv_list(ncolv_mol)%expected_value_growth_speed = cons_info%const_colv_target_growth(j)
935 ! In case of Restraint let's check for possible restart values
936 IF (colv_list(ncolv_mol)%restraint%active .AND. &
937 (colv_list(ncolv_mol)%expected_value_growth_speed == 0.0_dp)) THEN
938 gind = gind + 1
939 IF (restart_restraint_clv) THEN
940 CALL section_vals_val_get(colvar_rest, "_DEFAULT_KEYWORD_", &
941 i_rep_val=gind, r_val=rmod)
942 colv_list(ncolv_mol)%expected_value = rmod
943 ELSE
944 rmod = colv_list(ncolv_mol)%expected_value
945 CALL section_vals_val_set(colvar_rest, "_DEFAULT_KEYWORD_", &
946 i_rep_val=gind, r_val=rmod)
947 END IF
948 END IF
949 ! Only if torsion let's take into account the singularity in the definition
950 ! of the dihedral
951 IF (cons_info%colvar_set(j)%colvar%type_id == torsion_colvar_id) THEN
952 cons_info%colvar_set(j)%colvar%torsion_param%o0 = colv_list(ncolv_mol)%expected_value
953 END IF
954 END DO
955 END SUBROUTINE setup_colv_list
956
957! **************************************************************************************************
958!> \brief Setup the g3x3_list for the packing of constraints
959!> \param g3x3_list ...
960!> \param ilist ...
961!> \param cons_info ...
962!> \param ng3x3_restraint ...
963!> \par History
964!> Updated 2007 for intermolecular constraints
965!> \author Teodoro Laino [2007]
966! **************************************************************************************************
967 SUBROUTINE setup_g3x3_list(g3x3_list, ilist, cons_info, ng3x3_restraint)
968 TYPE(g3x3_constraint_type), DIMENSION(:), POINTER :: g3x3_list
969 INTEGER, DIMENSION(:), POINTER :: ilist
970 TYPE(constraint_info_type), POINTER :: cons_info
971 INTEGER, INTENT(OUT) :: ng3x3_restraint
972
973 INTEGER :: j, ng3x3
974
975 ng3x3_restraint = 0
976 DO ng3x3 = 1, SIZE(ilist)
977 j = ilist(ng3x3)
978 g3x3_list(ng3x3)%a = cons_info%const_g33_a(j)
979 g3x3_list(ng3x3)%b = cons_info%const_g33_b(j)
980 g3x3_list(ng3x3)%c = cons_info%const_g33_c(j)
981 g3x3_list(ng3x3)%dab = cons_info%const_g33_dab(j)
982 g3x3_list(ng3x3)%dac = cons_info%const_g33_dac(j)
983 g3x3_list(ng3x3)%dbc = cons_info%const_g33_dbc(j)
984 ! Restraint
985 g3x3_list(ng3x3)%restraint%active = cons_info%g33_restraint(j)
986 g3x3_list(ng3x3)%restraint%k0 = cons_info%g33_k0(j)
987 IF (g3x3_list(ng3x3)%restraint%active) ng3x3_restraint = ng3x3_restraint + 1
988 END DO
989
990 END SUBROUTINE setup_g3x3_list
991
992! **************************************************************************************************
993!> \brief Setup the g4x6_list for the packing of constraints
994!> \param g4x6_list ...
995!> \param ilist ...
996!> \param cons_info ...
997!> \param ng4x6_restraint ...
998!> \par History
999!> Updated 2007 for intermolecular constraints
1000!> \author Teodoro Laino [2007]
1001! **************************************************************************************************
1002 SUBROUTINE setup_g4x6_list(g4x6_list, ilist, cons_info, ng4x6_restraint)
1003 TYPE(g4x6_constraint_type), DIMENSION(:), POINTER :: g4x6_list
1004 INTEGER, DIMENSION(:), POINTER :: ilist
1005 TYPE(constraint_info_type), POINTER :: cons_info
1006 INTEGER, INTENT(OUT) :: ng4x6_restraint
1007
1008 INTEGER :: j, ng4x6
1009
1010 ng4x6 = 0
1011 ng4x6_restraint = 0
1012 DO ng4x6 = 1, SIZE(ilist)
1013 j = ilist(ng4x6)
1014 g4x6_list(ng4x6)%a = cons_info%const_g46_a(j)
1015 g4x6_list(ng4x6)%b = cons_info%const_g46_b(j)
1016 g4x6_list(ng4x6)%c = cons_info%const_g46_c(j)
1017 g4x6_list(ng4x6)%d = cons_info%const_g46_d(j)
1018 g4x6_list(ng4x6)%dab = cons_info%const_g46_dab(j)
1019 g4x6_list(ng4x6)%dac = cons_info%const_g46_dac(j)
1020 g4x6_list(ng4x6)%dbc = cons_info%const_g46_dbc(j)
1021 g4x6_list(ng4x6)%dad = cons_info%const_g46_dad(j)
1022 g4x6_list(ng4x6)%dbd = cons_info%const_g46_dbd(j)
1023 g4x6_list(ng4x6)%dcd = cons_info%const_g46_dcd(j)
1024 ! Restraint
1025 g4x6_list(ng4x6)%restraint%active = cons_info%g46_restraint(j)
1026 g4x6_list(ng4x6)%restraint%k0 = cons_info%g46_k0(j)
1027 IF (g4x6_list(ng4x6)%restraint%active) ng4x6_restraint = ng4x6_restraint + 1
1028 END DO
1029
1030 END SUBROUTINE setup_g4x6_list
1031
1032! **************************************************************************************************
1033!> \brief Setup the vsite_list for the packing of constraints
1034!> \param vsite_list ...
1035!> \param ilist ...
1036!> \param cons_info ...
1037!> \param nvsite_restraint ...
1038!> \par History
1039!> \author Marcel Baer [2008]
1040! **************************************************************************************************
1041 SUBROUTINE setup_vsite_list(vsite_list, ilist, cons_info, nvsite_restraint)
1042 TYPE(vsite_constraint_type), DIMENSION(:), POINTER :: vsite_list
1043 INTEGER, DIMENSION(:), POINTER :: ilist
1044 TYPE(constraint_info_type), POINTER :: cons_info
1045 INTEGER, INTENT(OUT) :: nvsite_restraint
1046
1047 INTEGER :: j, nvsite
1048
1049 nvsite = 0
1050 nvsite_restraint = 0
1051 DO nvsite = 1, SIZE(ilist)
1052 j = ilist(nvsite)
1053 vsite_list(nvsite)%a = cons_info%const_vsite_a(j)
1054 vsite_list(nvsite)%b = cons_info%const_vsite_b(j)
1055 vsite_list(nvsite)%c = cons_info%const_vsite_c(j)
1056 vsite_list(nvsite)%d = cons_info%const_vsite_d(j)
1057 vsite_list(nvsite)%wbc = cons_info%const_vsite_wbc(j)
1058 vsite_list(nvsite)%wdc = cons_info%const_vsite_wdc(j)
1059 ! Restraint
1060 vsite_list(nvsite)%restraint%active = cons_info%vsite_restraint(j)
1061 vsite_list(nvsite)%restraint%k0 = cons_info%vsite_k0(j)
1062 IF (vsite_list(nvsite)%restraint%active) nvsite_restraint = nvsite_restraint + 1
1063 END DO
1064
1065 END SUBROUTINE setup_vsite_list
1066! **************************************************************************************************
1067!> \brief Setup the lcolv for the packing of constraints
1068!> \param lcolv ...
1069!> \param ilist ...
1070!> \param first_atom ...
1071!> \param last_atom ...
1072!> \param cons_info ...
1073!> \param particle_set ...
1074!> \param colvar_func_info ...
1075!> \param use_clv_info ...
1076!> \param cind ...
1077!> \par History
1078!> Updated 2007 for intermolecular constraints
1079!> \author Teodoro Laino [2007]
1080! **************************************************************************************************
1081 SUBROUTINE setup_lcolv(lcolv, ilist, first_atom, last_atom, cons_info, &
1082 particle_set, colvar_func_info, use_clv_info, &
1083 cind)
1084 TYPE(local_colvar_constraint_type), DIMENSION(:), &
1085 POINTER :: lcolv
1086 INTEGER, DIMENSION(:), POINTER :: ilist
1087 INTEGER, INTENT(IN) :: first_atom, last_atom
1088 TYPE(constraint_info_type), POINTER :: cons_info
1089 TYPE(particle_type), DIMENSION(:), OPTIONAL, &
1090 POINTER :: particle_set
1091 TYPE(section_vals_type), POINTER :: colvar_func_info
1092 LOGICAL, INTENT(IN) :: use_clv_info
1093 INTEGER, INTENT(INOUT) :: cind
1094
1095 INTEGER :: ind, k, kk
1096 REAL(kind=dp), DIMENSION(:), POINTER :: r_vals
1097
1098 DO kk = 1, SIZE(ilist)
1099 k = ilist(kk)
1100 lcolv(kk)%init = .false.
1101 lcolv(kk)%lambda = 0.0_dp
1102 lcolv(kk)%sigma = 0.0_dp
1103
1104 ! Set Up colvar variable
1105 NULLIFY (lcolv(kk)%colvar, lcolv(kk)%colvar_old)
1106 ! Colvar
1107 CALL colvar_clone(lcolv(kk)%colvar, cons_info%colvar_set(k)%colvar, &
1108 i_atom_offset=first_atom - 1)
1109
1110 ! Some COLVARS may need additional information for evaluating the
1111 ! functional form: this is the case for COLVARS which depend on the
1112 ! initial position of the atoms: This information is stored in a proper
1113 ! container in the COLVAR_RESTART section..
1114 IF ((lcolv(kk)%colvar%type_id == xyz_diag_colvar_id) .OR. &
1115 (lcolv(kk)%colvar%type_id == xyz_outerdiag_colvar_id)) THEN
1116 cind = cind + 1
1117 IF (use_clv_info) THEN
1118 CALL section_vals_val_get(colvar_func_info, "_DEFAULT_KEYWORD_", &
1119 i_rep_val=cind, r_vals=r_vals)
1120 SELECT CASE (lcolv(kk)%colvar%type_id)
1121 CASE (xyz_diag_colvar_id)
1122 cpassert(SIZE(r_vals) == 3)
1123 lcolv(kk)%colvar%xyz_diag_param%r0 = r_vals
1124 CASE (xyz_outerdiag_colvar_id)
1125 cpassert(SIZE(r_vals) == 6)
1126 lcolv(kk)%colvar%xyz_outerdiag_param%r0(:, 1) = r_vals(1:3)
1127 lcolv(kk)%colvar%xyz_outerdiag_param%r0(:, 2) = r_vals(4:6)
1128 END SELECT
1129 ELSE
1130 SELECT CASE (lcolv(kk)%colvar%type_id)
1131 CASE (xyz_diag_colvar_id)
1132 ALLOCATE (r_vals(3))
1133 ind = first_atom - 1 + lcolv(kk)%colvar%xyz_diag_param%i_atom
1134 r_vals = particle_set(ind)%r
1135 lcolv(kk)%colvar%xyz_diag_param%r0 = r_vals
1136 CASE (xyz_outerdiag_colvar_id)
1137 ALLOCATE (r_vals(6))
1138 ind = first_atom - 1 + lcolv(kk)%colvar%xyz_outerdiag_param%i_atoms(1)
1139 r_vals(1:3) = particle_set(ind)%r
1140 ind = first_atom - 1 + lcolv(kk)%colvar%xyz_outerdiag_param%i_atoms(2)
1141 r_vals(4:6) = particle_set(ind)%r
1142 lcolv(kk)%colvar%xyz_outerdiag_param%r0(:, 1) = r_vals(1:3)
1143 lcolv(kk)%colvar%xyz_outerdiag_param%r0(:, 2) = r_vals(4:6)
1144 END SELECT
1145 CALL section_vals_val_set(colvar_func_info, "_DEFAULT_KEYWORD_", &
1146 i_rep_val=cind, r_vals_ptr=r_vals)
1147 END IF
1148 END IF
1149
1150 ! Setup Colvar_old
1151 CALL colvar_clone(lcolv(kk)%colvar_old, lcolv(kk)%colvar)
1152
1153 ! Check for consistency in the constraint definition
1154 IF (any(lcolv(kk)%colvar%i_atom > last_atom) .OR. &
1155 any(lcolv(kk)%colvar%i_atom < first_atom)) THEN
1156 WRITE (*, '(T5,"|",T8,A)') "Error in constraints setup!"
1157 WRITE (*, '(T5,"|",T8,A)') "A constraint has been defined for a molecule type", &
1158 " but the atoms specified in the constraint and the atoms defined for", &
1159 " the molecule DO NOT match!", &
1160 "This could be very probable due to a wrong connectivity, or an error", &
1161 " in the constraint specification in the input file.", &
1162 " Please check it carefully!"
1163 cpabort("")
1164 END IF
1165 END DO
1166 END SUBROUTINE setup_lcolv
1167
1168! **************************************************************************************************
1169!> \brief Setup the lg3x3 for the packing of constraints
1170!> \param lg3x3 ...
1171!> \param g3x3_list ...
1172!> \param first_atom ...
1173!> \param last_atom ...
1174!> \par History
1175!> Updated 2007 for intermolecular constraints
1176!> \author Teodoro Laino [2007]
1177! **************************************************************************************************
1178 SUBROUTINE setup_lg3x3(lg3x3, g3x3_list, first_atom, last_atom)
1179 TYPE(local_g3x3_constraint_type), DIMENSION(:), &
1180 POINTER :: lg3x3
1181 TYPE(g3x3_constraint_type), DIMENSION(:), POINTER :: g3x3_list
1182 INTEGER, INTENT(IN) :: first_atom, last_atom
1183
1184 INTEGER :: kk
1185
1186 DO kk = 1, SIZE(lg3x3)
1187 lg3x3(kk)%init = .false.
1188 lg3x3(kk)%scale = 0.0_dp
1189 lg3x3(kk)%scale_old = 0.0_dp
1190 lg3x3(kk)%fa = 0.0_dp
1191 lg3x3(kk)%fb = 0.0_dp
1192 lg3x3(kk)%fc = 0.0_dp
1193 lg3x3(kk)%ra_old = 0.0_dp
1194 lg3x3(kk)%rb_old = 0.0_dp
1195 lg3x3(kk)%rc_old = 0.0_dp
1196 lg3x3(kk)%va = 0.0_dp
1197 lg3x3(kk)%vb = 0.0_dp
1198 lg3x3(kk)%vc = 0.0_dp
1199 lg3x3(kk)%lambda = 0.0_dp
1200 IF ((g3x3_list(kk)%a + first_atom - 1 < first_atom) .OR. &
1201 (g3x3_list(kk)%b + first_atom - 1 < first_atom) .OR. &
1202 (g3x3_list(kk)%c + first_atom - 1 < first_atom) .OR. &
1203 (g3x3_list(kk)%a + first_atom - 1 > last_atom) .OR. &
1204 (g3x3_list(kk)%b + first_atom - 1 > last_atom) .OR. &
1205 (g3x3_list(kk)%c + first_atom - 1 > last_atom)) THEN
1206 WRITE (*, '(T5,"|",T8,A)') "Error in constraints setup!"
1207 WRITE (*, '(T5,"|",T8,A)') "A constraint has been defined for a molecule type", &
1208 " but the atoms specified in the constraint and the atoms defined for", &
1209 " the molecule DO NOT match!", &
1210 "This could be very probable due to a wrong connectivity, or an error", &
1211 " in the constraint specification in the input file.", &
1212 " Please check it carefully!"
1213 cpabort("")
1214 END IF
1215 END DO
1216
1217 END SUBROUTINE setup_lg3x3
1218
1219! **************************************************************************************************
1220!> \brief Setup the lg4x6 for the packing of constraints
1221!> \param lg4x6 ...
1222!> \param g4x6_list ...
1223!> \param first_atom ...
1224!> \param last_atom ...
1225!> \par History
1226!> Updated 2007 for intermolecular constraints
1227!> \author Teodoro Laino [2007]
1228! **************************************************************************************************
1229 SUBROUTINE setup_lg4x6(lg4x6, g4x6_list, first_atom, last_atom)
1230 TYPE(local_g4x6_constraint_type), DIMENSION(:), &
1231 POINTER :: lg4x6
1232 TYPE(g4x6_constraint_type), DIMENSION(:), POINTER :: g4x6_list
1233 INTEGER, INTENT(IN) :: first_atom, last_atom
1234
1235 INTEGER :: kk
1236
1237 DO kk = 1, SIZE(lg4x6)
1238 lg4x6(kk)%init = .false.
1239 lg4x6(kk)%scale = 0.0_dp
1240 lg4x6(kk)%scale_old = 0.0_dp
1241 lg4x6(kk)%fa = 0.0_dp
1242 lg4x6(kk)%fb = 0.0_dp
1243 lg4x6(kk)%fc = 0.0_dp
1244 lg4x6(kk)%fd = 0.0_dp
1245 lg4x6(kk)%fe = 0.0_dp
1246 lg4x6(kk)%ff = 0.0_dp
1247 lg4x6(kk)%ra_old = 0.0_dp
1248 lg4x6(kk)%rb_old = 0.0_dp
1249 lg4x6(kk)%rc_old = 0.0_dp
1250 lg4x6(kk)%rd_old = 0.0_dp
1251 lg4x6(kk)%re_old = 0.0_dp
1252 lg4x6(kk)%rf_old = 0.0_dp
1253 lg4x6(kk)%va = 0.0_dp
1254 lg4x6(kk)%vb = 0.0_dp
1255 lg4x6(kk)%vc = 0.0_dp
1256 lg4x6(kk)%vd = 0.0_dp
1257 lg4x6(kk)%ve = 0.0_dp
1258 lg4x6(kk)%vf = 0.0_dp
1259 lg4x6(kk)%lambda = 0.0_dp
1260 IF ((g4x6_list(kk)%a + first_atom - 1 < first_atom) .OR. &
1261 (g4x6_list(kk)%b + first_atom - 1 < first_atom) .OR. &
1262 (g4x6_list(kk)%c + first_atom - 1 < first_atom) .OR. &
1263 (g4x6_list(kk)%d + first_atom - 1 < first_atom) .OR. &
1264 (g4x6_list(kk)%a + first_atom - 1 > last_atom) .OR. &
1265 (g4x6_list(kk)%b + first_atom - 1 > last_atom) .OR. &
1266 (g4x6_list(kk)%c + first_atom - 1 > last_atom) .OR. &
1267 (g4x6_list(kk)%d + first_atom - 1 > last_atom)) THEN
1268 WRITE (*, '(T5,"|",T8,A)') "Error in constraints setup!"
1269 WRITE (*, '(T5,"|",T8,A)') "A constrained has been defined for a molecule type", &
1270 " but the atoms specified in the constraint and the atoms defined for", &
1271 " the molecule DO NOT match!", &
1272 "This could be very probable due to a wrong connectivity, or an error", &
1273 " in the constraint specification in the input file.", &
1274 " Please check it carefully!"
1275 cpabort("")
1276 END IF
1277 END DO
1278
1279 END SUBROUTINE setup_lg4x6
1280
1281! **************************************************************************************************
1282!> \brief Gives back a list of molecule to which apply the constraint
1283!> \param const_mol ...
1284!> \param const_molname ...
1285!> \param const_intermolecular ...
1286!> \param constr_x_mol ...
1287!> \param constr_x_glob ...
1288!> \param molecule_kind_set ...
1289!> \param exclude_qm ...
1290!> \param exclude_mm ...
1291!> \par History
1292!> Updated 2007 for intermolecular constraints
1293!> \author Teodoro Laino [2006]
1294! **************************************************************************************************
1295 SUBROUTINE give_constraint_array(const_mol, const_molname, const_intermolecular, &
1296 constr_x_mol, constr_x_glob, molecule_kind_set, exclude_qm, exclude_mm)
1297
1298 INTEGER, DIMENSION(:), POINTER :: const_mol
1299 CHARACTER(LEN=default_string_length), &
1300 DIMENSION(:), POINTER :: const_molname
1301 LOGICAL, DIMENSION(:), POINTER :: const_intermolecular
1302 TYPE(constr_list_type), DIMENSION(:), POINTER :: constr_x_mol
1303 INTEGER, DIMENSION(:), POINTER :: constr_x_glob
1304 TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
1305 LOGICAL, DIMENSION(:), POINTER :: exclude_qm, exclude_mm
1306
1307 CHARACTER(len=*), PARAMETER :: routinen = 'give_constraint_array'
1308
1309 CHARACTER(LEN=default_string_length) :: myname, name
1310 INTEGER :: handle, i, iglob, isize, k
1311 LOGICAL :: found_molname, is_qm
1312 TYPE(molecule_kind_type), POINTER :: molecule_kind
1313
1314 CALL timeset(routinen, handle)
1315 NULLIFY (molecule_kind)
1316 ALLOCATE (constr_x_mol(SIZE(molecule_kind_set)))
1317 DO i = 1, SIZE(constr_x_mol)
1318 NULLIFY (constr_x_mol(i)%constr)
1319 ALLOCATE (constr_x_mol(i)%constr(0))
1320 END DO
1321 cpassert(SIZE(const_mol) == SIZE(const_molname))
1322 iglob = 0
1323 DO i = 1, SIZE(const_mol)
1324 IF (const_intermolecular(i)) THEN
1325 ! Intermolecular constraint
1326 iglob = iglob + 1
1327 CALL reallocate(constr_x_glob, 1, iglob)
1328 constr_x_glob(iglob) = i
1329 ELSE
1330 ! Intramolecular constraint
1331 IF (const_mol(i) /= 0) THEN
1332 k = const_mol(i)
1333 IF (k > SIZE(molecule_kind_set)) &
1334 CALL cp_abort(__location__, &
1335 "A constraint has been specified providing the molecule index. But the"// &
1336 " molecule index ("//cp_to_string(k)//") is out of range of the possible"// &
1337 " molecule kinds ("//cp_to_string(SIZE(molecule_kind_set))//").")
1338 isize = SIZE(constr_x_mol(k)%constr)
1339 CALL reallocate(constr_x_mol(k)%constr, 1, isize + 1)
1340 constr_x_mol(k)%constr(isize + 1) = i
1341 ELSE
1342 myname = const_molname(i)
1343 found_molname = .false.
1344 DO k = 1, SIZE(molecule_kind_set)
1345 molecule_kind => molecule_kind_set(k)
1346 name = molecule_kind%name
1347 is_qm = qmmm_ff_precond_only_qm(id1=name)
1348 IF (is_qm .AND. exclude_qm(i)) cycle
1349 IF (.NOT. is_qm .AND. exclude_mm(i)) cycle
1350 IF (name == myname) THEN
1351 isize = SIZE(constr_x_mol(k)%constr)
1352 CALL reallocate(constr_x_mol(k)%constr, 1, isize + 1)
1353 constr_x_mol(k)%constr(isize + 1) = i
1354 found_molname = .true.
1355 END IF
1356 END DO
1357 CALL print_warning_molname(found_molname, myname)
1358 END IF
1359 END IF
1360 END DO
1361 CALL timestop(handle)
1362 END SUBROUTINE give_constraint_array
1363
1364! **************************************************************************************************
1365!> \brief Prints a warning message if undefined molnames are used to define constraints
1366!> \param found ...
1367!> \param name ...
1368!> \author Teodoro Laino [2007] - Zurich University
1369! **************************************************************************************************
1370 SUBROUTINE print_warning_molname(found, name)
1371 LOGICAL, INTENT(IN) :: found
1372 CHARACTER(LEN=*), INTENT(IN) :: name
1373
1374 IF (.NOT. found) &
1375 CALL cp_warn(__location__, &
1376 " MOLNAME ("//trim(name)//") was defined for constraints, but this molecule name "// &
1377 "is not defined. Please check carefully your PDB, PSF (has priority over PDB) or "// &
1378 "input driven CP2K coordinates. In case you may not find the reason for this warning "// &
1379 "it may be a good idea to print all molecule information (including kind name) activating "// &
1380 "the print_key MOLECULES specific of the SUBSYS%PRINT section. ")
1381
1382 END SUBROUTINE print_warning_molname
1383
1384END MODULE topology_constraint_util
cp_log_handling::cp_to_string
Definition cp_log_handling.F:90
memory_utilities::reallocate
Definition memory_utilities.F:27
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind
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.
Definition atomic_kind_types.F:138
atomic_kind_types::is_hydrogen
elemental logical function, public is_hydrogen(atomic_kind)
Determines if the atomic_kind is HYDROGEN.
Definition atomic_kind_types.F:380
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
cell_types::use_perd_xyz
integer, parameter, public use_perd_xyz
Definition cell_types.F:42
cell_types::use_perd_y
integer, parameter, public use_perd_y
Definition cell_types.F:42
cell_types::use_perd_xz
integer, parameter, public use_perd_xz
Definition cell_types.F:42
cell_types::use_perd_x
integer, parameter, public use_perd_x
Definition cell_types.F:42
cell_types::use_perd_z
integer, parameter, public use_perd_z
Definition cell_types.F:42
cell_types::use_perd_yz
integer, parameter, public use_perd_yz
Definition cell_types.F:42
cell_types::use_perd_xy
integer, parameter, public use_perd_xy
Definition cell_types.F:42
colvar_methods
defines collective variables s({R}) and the derivative of this variable wrt R these can then be used ...
Definition colvar_methods.F:16
colvar_methods::colvar_eval_mol_f
subroutine, public colvar_eval_mol_f(colvar, cell, particles, pos, fixd_list)
evaluates the derivatives (dsdr) given and due to the given colvar variables in a molecular environme...
Definition colvar_methods.F:1631
colvar_types
Initialize the collective variables types.
Definition colvar_types.F:15
colvar_types::colvar_release
recursive subroutine, public colvar_release(colvar)
releases the memory that might have been allocated by the colvar
Definition colvar_types.F:1176
colvar_types::xyz_outerdiag_colvar_id
integer, parameter, public xyz_outerdiag_colvar_id
Definition colvar_types.F:56
colvar_types::colvar_create
subroutine, public colvar_create(colvar, colvar_id)
initializes a colvar_param type
Definition colvar_types.F:421
colvar_types::colvar_clone
recursive subroutine, public colvar_clone(colvar_out, colvar_in, i_atom_offset)
Clone a colvar type.
Definition colvar_types.F:1353
colvar_types::colvar_p_reallocate
subroutine, public colvar_p_reallocate(colvar_set, lb1_new, ub1_new)
Change the dimension of a colvar_p_type.
Definition colvar_types.F:1791
colvar_types::dist_colvar_id
integer, parameter, public dist_colvar_id
Definition colvar_types.F:56
colvar_types::colvar_setup
subroutine, public colvar_setup(colvar)
Finalize the setup of the collective variable.
Definition colvar_types.F:502
colvar_types::torsion_colvar_id
integer, parameter, public torsion_colvar_id
Definition colvar_types.F:56
colvar_types::xyz_diag_colvar_id
integer, parameter, public xyz_diag_colvar_id
Definition colvar_types.F:56
colvar_utils
evaluations of colvar for internal coordinates schemes
Definition colvar_utils.F:14
colvar_utils::post_process_colvar
subroutine, public post_process_colvar(colvar, particles)
Complete the description of the COORDINATION colvar when defined using KINDS.
Definition colvar_utils.F:566
cp_log_handling
various routines to log and control the output. The idea is that decisions about where to log should ...
Definition cp_log_handling.F:41
cp_log_handling::cp_get_default_logger
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
Definition cp_log_handling.F:234
cp_output_handling
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
Definition cp_output_handling.F:25
cp_output_handling::cp_print_key_unit_nr
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)
...
Definition cp_output_handling.F:803
cp_output_handling::cp_print_key_finished_output
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,...
Definition cp_output_handling.F:1013
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:17
input_constants::do_constr_atomic
integer, parameter, public do_constr_atomic
Definition input_constants.F:79
input_constants::do_constr_molec
integer, parameter, public do_constr_molec
Definition input_constants.F:79
input_section_types
objects that represent the structure of input sections and the data contained in an input section
Definition input_section_types.F:15
input_section_types::section_vals_val_set
subroutine, public section_vals_val_set(section_vals, keyword_name, i_rep_section, i_rep_val, val, l_val, i_val, r_val, c_val, l_vals_ptr, i_vals_ptr, r_vals_ptr, c_vals_ptr)
sets the requested value
Definition input_section_types.F:1216
input_section_types::section_vals_get_subs_vals
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
Definition input_section_types.F:724
input_section_types::section_vals_get
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
Definition input_section_types.F:697
input_section_types::section_vals_val_get
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Definition input_section_types.F:1040
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_string_length
integer, parameter, public default_string_length
Definition kinds.F:57
memory_utilities
Utility routines for the memory handling.
Definition memory_utilities.F:14
molecule_kind_types
Define the molecule kind structure types and the corresponding functionality.
Definition molecule_kind_types.F:16
molecule_kind_types::get_molecule_kind
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.
Definition molecule_kind_types.F:510
molecule_kind_types::setup_colvar_counters
subroutine, public setup_colvar_counters(colv_list, ncolv)
...
Definition molecule_kind_types.F:239
molecule_kind_types::set_molecule_kind
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.
Definition molecule_kind_types.F:779
molecule_types
Define the data structure for the molecule information.
Definition molecule_types.F:16
molecule_types::get_molecule
subroutine, public get_molecule(molecule, molecule_kind, lmi, lci, lg3x3, lg4x6, lcolv, first_atom, last_atom, first_shell, last_shell)
Get components from a molecule data set.
Definition molecule_types.F:277
molecule_types::set_molecule
subroutine, public set_molecule(molecule, molecule_kind, lmi, lci, lcolv, lg3x3, lg4x6)
Set a molecule data set.
Definition molecule_types.F:337
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
qmmm_ff_fist
Definition qmmm_ff_fist.F:11
qmmm_ff_fist::qmmm_ff_precond_only_qm
logical function, public qmmm_ff_precond_only_qm(id1, id2, id3, id4, is_link)
This function handles the atom names and modifies the "_QM_" prefix, in order to find the parameters ...
Definition qmmm_ff_fist.F:39
qmmm_types_low
Definition qmmm_types_low.F:13
topology_constraint_util
Collection of subroutine needed for topology related things.
Definition topology_constraint_util.F:13
topology_constraint_util::topology_constraint_pack
subroutine, public topology_constraint_pack(molecule_kind_set, molecule_set, topology, qmmm_env, particle_set, input_file, subsys_section, gci)
Pack in all the information needed for the constraints.
Definition topology_constraint_util.F:87
topology_types
Definition topology_types.F:18
topology
Control for reading in different topologies and coordinates.
Definition topology.F:13
atomic_kind_types::atomic_kind_type
Provides all information about an atomic kind.
Definition atomic_kind_types.F:45
colvar_types::colvar_counters
Definition colvar_types.F:388
colvar_types::colvar_type
parameters for a collective variable
Definition colvar_types.F:346
cp_log_handling::cp_logger_type
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Definition cp_log_handling.F:140
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:126
molecule_kind_types::atom_type
Definition molecule_kind_types.F:56
molecule_kind_types::bond_type
Definition molecule_kind_types.F:69
molecule_kind_types::colvar_constraint_type
Definition molecule_kind_types.F:117
molecule_kind_types::fixd_constraint_type
Definition molecule_kind_types.F:149
molecule_kind_types::g3x3_constraint_type
Definition molecule_kind_types.F:128
molecule_kind_types::g4x6_constraint_type
Definition molecule_kind_types.F:135
molecule_kind_types::molecule_kind_type
Definition molecule_kind_types.F:161
molecule_kind_types::vsite_constraint_type
Definition molecule_kind_types.F:142
molecule_types::global_constraint_type
Definition molecule_types.F:93
molecule_types::local_colvar_constraint_type
Definition molecule_types.F:41
molecule_types::local_constraint_type
Definition molecule_types.F:86
molecule_types::local_g3x3_constraint_type
Definition molecule_types.F:49
molecule_types::local_g4x6_constraint_type
Definition molecule_types.F:63
molecule_types::molecule_type
Definition molecule_types.F:110
particle_types::particle_type
Definition particle_types.F:35
qmmm_types_low::qmmm_env_mm_type
...
Definition qmmm_types_low.F:238
topology_types::constr_list_type
Definition topology_types.F:179
topology_types::constraint_info_type
Definition topology_types.F:64
topology_types::topology_parameters_type
Definition topology_types.F:142