(git:374b731)
Loading...
Searching...
No Matches
topology_gromos.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! **************************************************************************************************
20 USE cp_units, ONLY: cp_unit_to_cp2k
25 USE kinds, ONLY: default_string_length,&
26 dp
29 USE string_table, ONLY: s2s,&
30 str2id
34#include "./base/base_uses.f90"
35
36 IMPLICIT NONE
37
38 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'topology_gromos'
39
40 PRIVATE
42
43CONTAINS
44
45! **************************************************************************************************
46!> \brief Read GROMOS topology file
47!> \param file_name ...
48!> \param topology ...
49!> \param para_env ...
50!> \param subsys_section ...
51! **************************************************************************************************
52 SUBROUTINE read_topology_gromos(file_name, topology, para_env, subsys_section)
53 CHARACTER(LEN=*), INTENT(IN) :: file_name
54 TYPE(topology_parameters_type), INTENT(INOUT) :: topology
55 TYPE(mp_para_env_type), POINTER :: para_env
56 TYPE(section_vals_type), POINTER :: subsys_section
57
58 CHARACTER(len=*), PARAMETER :: routinen = 'read_topology_gromos'
59
60 CHARACTER :: ctemp
61 CHARACTER(LEN=default_string_length) :: label, string
62 CHARACTER(LEN=default_string_length), &
63 DIMENSION(21) :: avail_section
64 CHARACTER(LEN=default_string_length), POINTER :: namearray1(:), namearray2(:)
65 INTEGER :: begin, handle, i, iatom, ibond, icon, ii(50), index_now, iresid, isolvent, itemp, &
66 itype, iw, jatom, natom, natom_prev, nbond, nbond_prev, ncon, nsolute, nsolvent, ntype, &
67 offset, offset2, stat
68 INTEGER, POINTER :: ba(:), bb(:), na(:)
69 LOGICAL :: found
70 REAL(kind=dp) :: ftemp
71 REAL(kind=dp), POINTER :: ac(:), am(:)
72 TYPE(atom_info_type), POINTER :: atom_info
73 TYPE(connectivity_info_type), POINTER :: conn_info
74 TYPE(cp_logger_type), POINTER :: logger
75 TYPE(cp_parser_type) :: parser
76
77 NULLIFY (logger)
78 NULLIFY (namearray1, namearray2)
79 logger => cp_get_default_logger()
80 iw = cp_print_key_unit_nr(logger, subsys_section, "PRINT%TOPOLOGY_INFO/GTOP_INFO", &
81 extension=".subsysLog")
82 CALL timeset(routinen, handle)
83
84 avail_section(1) = "TITLE"
85 avail_section(2) = "TOPPHYSCON"
86 avail_section(3) = "TOPVERSION"
87 avail_section(4) = "ATOMTYPENAME"
88 avail_section(5) = "RESNAME"
89 avail_section(6) = "SOLUTEATOM"
90 avail_section(7) = "BONDTYPE"
91 avail_section(8) = "BONDH"
92 avail_section(9) = "BOND"
93 avail_section(10) = "BONDANGLETYPE"
94 avail_section(11) = "BONDANGLEH"
95 avail_section(12) = "BONDANGLE"
96 avail_section(13) = "IMPDIHEDRALTYPE"
97 avail_section(14) = "IMPDIHEDRALH"
98 avail_section(15) = "IMPDIHEDRAL"
99 avail_section(16) = "DIHEDRALTYPE"
100 avail_section(17) = "DIHEDRALH"
101 avail_section(18) = "DIHEDRAL"
102 avail_section(19) = "LJPARAMETERS"
103 avail_section(20) = "SOLVENTATOM"
104 avail_section(21) = "SOLVENTCONSTR"
105
106 atom_info => topology%atom_info
107 conn_info => topology%conn_info
108
109 natom_prev = 0
110 IF (ASSOCIATED(atom_info%id_molname)) natom_prev = SIZE(atom_info%id_molname)
111 ! TITLE SECTION
112 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the TITLE section'
113 CALL parser_create(parser, file_name, para_env=para_env)
114 label = trim(avail_section(1))
115 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
116 IF (found) THEN
117 DO
118 CALL parser_get_next_line(parser, 1)
119 CALL parser_get_object(parser, string, string_length=80)
120 IF (string == trim("END")) EXIT
121 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| ", trim(string)
122 END DO
123 END IF
124 CALL parser_release(parser)
125 ! TOPPHYSCON SECTION
126 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the TOPPHYSCON section'
127 CALL parser_create(parser, file_name, para_env=para_env)
128 label = trim(avail_section(2))
129 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
130 IF (found) THEN
131 DO
132 CALL parser_get_next_line(parser, 1)
133 CALL parser_get_object(parser, string)
134 IF (string == trim("END")) EXIT
135 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| ", trim(string)
136 END DO
137 END IF
138 CALL parser_release(parser)
139 ! TOPVERSION SECTION
140 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the TOPVERSION section'
141 CALL parser_create(parser, file_name, para_env=para_env)
142 label = trim(avail_section(3))
143 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
144 IF (found) THEN
145 DO
146 CALL parser_get_next_line(parser, 1)
147 CALL parser_get_object(parser, string)
148 IF (string == trim("END")) EXIT
149 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| ", trim(string)
150 END DO
151 END IF
152 CALL parser_release(parser)
153 ! ATOMTYPENAME SECTION
154 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the ATOMTYPENAME section'
155 CALL parser_create(parser, file_name, para_env=para_env)
156 label = trim(avail_section(4))
157 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
158 IF (found) THEN
159 CALL parser_get_next_line(parser, 1)
160 CALL parser_get_object(parser, ntype)
161 CALL reallocate(namearray1, 1, ntype)
162 DO itype = 1, ntype
163 CALL parser_get_next_line(parser, 1)
164 CALL parser_get_object(parser, namearray1(itype))
165 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| ", trim(namearray1(itype))
166 END DO
167 END IF
168 CALL parser_release(parser)
169 ! RESNAME SECTION
170 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the RESNAME section'
171 CALL parser_create(parser, file_name, para_env=para_env)
172 label = trim(avail_section(5))
173 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
174 IF (found) THEN
175 CALL parser_get_next_line(parser, 1)
176 CALL parser_get_object(parser, ntype)
177 CALL reallocate(namearray2, 1, ntype)
178 DO itype = 1, ntype
179 CALL parser_get_next_line(parser, 1)
180 CALL parser_get_object(parser, namearray2(itype))
181 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| ", trim(namearray2(itype))
182 END DO
183 END IF
184 CALL parser_release(parser)
185 ! SOLUTEATOM SECTION
186 iresid = 1
187 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the SOLUTEATOM section'
188 CALL parser_create(parser, file_name, para_env=para_env)
189 label = trim(avail_section(6))
190 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
191 IF (found) THEN
192 CALL parser_get_next_line(parser, 1)
193 CALL parser_get_object(parser, natom)
194 CALL reallocate(atom_info%id_molname, 1, natom_prev + natom)
195 CALL reallocate(atom_info%resid, 1, natom_prev + natom)
196 CALL reallocate(atom_info%id_resname, 1, natom_prev + natom)
197 CALL reallocate(atom_info%id_atmname, 1, natom_prev + natom)
198 CALL reallocate(atom_info%id_element, 1, natom_prev + natom)
199 CALL reallocate(atom_info%atm_charge, 1, natom_prev + natom)
200 CALL reallocate(atom_info%atm_mass, 1, natom_prev + natom)
201 CALL parser_get_next_line(parser, 1)
202 DO iatom = 1, natom
203 index_now = iatom + natom_prev
204 CALL parser_get_object(parser, itemp)
205 CALL parser_get_object(parser, itemp)
206 atom_info%resid(index_now) = itemp
207 atom_info%id_molname(index_now) = str2id(s2s(namearray2(itemp)))
208 atom_info%id_resname(index_now) = str2id(s2s(namearray2(itemp)))
209 CALL parser_get_object(parser, string)
210 CALL parser_get_object(parser, itemp)
211 atom_info%id_atmname(index_now) = str2id(s2s(namearray1(itemp)))
212 atom_info%id_element(index_now) = str2id(s2s(namearray1(itemp)))
213 CALL parser_get_object(parser, atom_info%atm_mass(index_now))
214 CALL parser_get_object(parser, atom_info%atm_charge(index_now))
215 CALL parser_get_object(parser, itemp)
216 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT SOLUTEATOM INFO HERE!!!!"
217 CALL parser_get_object(parser, ntype)
218 DO i = 1, 50
219 ii(i) = -1
220 END DO
221 stat = -1
222 begin = 1
223 IF (ntype /= 0) THEN
224 DO
225 IF (begin .EQ. 1) THEN
226 READ (parser%input_line, iostat=stat, fmt=*) itemp, itemp, ctemp, itemp, ftemp, ftemp, &
227 itemp, itemp, (ii(i), i=begin, ntype)
228 ELSE IF (begin .GT. 1) THEN
229 CALL parser_get_next_line(parser, 1)
230 READ (parser%input_line, iostat=stat, fmt=*) (ii(i), i=begin, ntype)
231 END IF
232 DO i = ntype, 1, -1
233 IF (ii(i) .LT. 0) THEN
234 begin = i
235 END IF
236 END DO
237 IF (stat .EQ. 0) EXIT
238 END DO
239 END IF
240 ! 1-4 list
241 CALL parser_get_next_line(parser, 1)
242 CALL parser_get_object(parser, ntype)
243 IF (ntype /= 0) THEN
244 itemp = (itemp - 1)/6 + 1
245 offset = 0
246 IF (ASSOCIATED(conn_info%onfo_a)) offset = SIZE(conn_info%onfo_a)
247 CALL reallocate(conn_info%onfo_a, 1, offset + ntype)
248 CALL reallocate(conn_info%onfo_b, 1, offset + ntype)
249 conn_info%onfo_a(offset + 1:offset + ntype) = index_now
250 DO i = 1, 50
251 ii(i) = -1
252 END DO
253 stat = -1
254 begin = 1
255 IF (ntype /= 0) THEN
256 DO
257 IF (begin .EQ. 1) THEN
258 READ (parser%input_line, iostat=stat, fmt=*) itemp, (ii(i), i=begin, ntype)
259 ELSE IF (begin .GT. 1) THEN
260 CALL parser_get_next_line(parser, 1)
261 READ (parser%input_line, iostat=stat, fmt=*) (ii(i), i=begin, ntype)
262 END IF
263 DO i = ntype, 1, -1
264 IF (ii(i) .LT. 0) THEN
265 begin = i
266 END IF
267 END DO
268 IF (stat .EQ. 0) EXIT
269 END DO
270 DO i = 1, ntype
271 conn_info%onfo_b(offset + i) = ii(i)
272 END DO
273 END IF
274 CALL parser_get_next_line(parser, 1)
275 ELSE
276 CALL parser_get_next_line(parser, 1)
277 END IF
278 END DO
279 END IF
280 nsolute = natom
281 CALL parser_release(parser)
282
283 CALL parser_create(parser, file_name, para_env=para_env)
284 ! BONDH SECTION
285 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the BONDH section'
286 label = trim(avail_section(8))
287 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
288 IF (found) THEN
289 CALL parser_get_next_line(parser, 1)
290 CALL parser_get_object(parser, ntype)
291 offset = 0
292 IF (ASSOCIATED(conn_info%bond_a)) offset = SIZE(conn_info%bond_a)
293 CALL reallocate(conn_info%bond_a, 1, offset + ntype)
294 CALL reallocate(conn_info%bond_b, 1, offset + ntype)
295 CALL reallocate(conn_info%bond_type, 1, offset + ntype)
296 DO itype = 1, ntype
297 CALL parser_get_next_line(parser, 1)
298 CALL parser_get_object(parser, conn_info%bond_a(offset + itype))
299 CALL parser_get_object(parser, conn_info%bond_b(offset + itype))
300 CALL parser_get_object(parser, itemp)
301 conn_info%bond_type(offset + itype) = itemp
302 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT BONDH INFO HERE!!!!"
303 END DO
304 conn_info%bond_a(offset + 1:offset + ntype) = conn_info%bond_a(offset + 1:offset + ntype) + natom_prev
305 conn_info%bond_b(offset + 1:offset + ntype) = conn_info%bond_b(offset + 1:offset + ntype) + natom_prev
306 END IF
307 ! BOND SECTION
308 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the BOND section'
309 label = trim(avail_section(9))
310 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
311 IF (found) THEN
312 CALL parser_get_next_line(parser, 1)
313 CALL parser_get_object(parser, ntype)
314 offset = 0
315 IF (ASSOCIATED(conn_info%bond_a)) offset = SIZE(conn_info%bond_a)
316 CALL reallocate(conn_info%bond_a, 1, offset + ntype)
317 CALL reallocate(conn_info%bond_b, 1, offset + ntype)
318 CALL reallocate(conn_info%bond_type, 1, offset + ntype)
319 DO itype = 1, ntype
320 CALL parser_get_next_line(parser, 1)
321 CALL parser_get_object(parser, conn_info%bond_a(offset + itype))
322 CALL parser_get_object(parser, conn_info%bond_b(offset + itype))
323 CALL parser_get_object(parser, itemp)
324 conn_info%bond_type(offset + itype) = itemp
325 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT BOND INFO HERE!!!!"
326 END DO
327 conn_info%bond_a(offset + 1:offset + ntype) = conn_info%bond_a(offset + 1:offset + ntype) + natom_prev
328 conn_info%bond_b(offset + 1:offset + ntype) = conn_info%bond_b(offset + 1:offset + ntype) + natom_prev
329 END IF
330 ! BONDANGLEH SECTION
331 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the BONDANGLEH section'
332 label = trim(avail_section(11))
333 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
334 IF (found) THEN
335 CALL parser_get_next_line(parser, 1)
336 CALL parser_get_object(parser, ntype)
337 offset = 0
338 IF (ASSOCIATED(conn_info%theta_a)) offset = SIZE(conn_info%theta_a)
339 CALL reallocate(conn_info%theta_a, 1, offset + ntype)
340 CALL reallocate(conn_info%theta_b, 1, offset + ntype)
341 CALL reallocate(conn_info%theta_c, 1, offset + ntype)
342 CALL reallocate(conn_info%theta_type, 1, offset + ntype)
343 DO itype = 1, ntype
344 CALL parser_get_next_line(parser, 1)
345 CALL parser_get_object(parser, conn_info%theta_a(offset + itype))
346 CALL parser_get_object(parser, conn_info%theta_b(offset + itype))
347 CALL parser_get_object(parser, conn_info%theta_c(offset + itype))
348 CALL parser_get_object(parser, itemp)
349 conn_info%theta_type(offset + itype) = itemp
350 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT BONDANGLEH INFO HERE!!!!"
351 END DO
352 conn_info%theta_a(offset + 1:offset + ntype) = conn_info%theta_a(offset + 1:offset + ntype) + natom_prev
353 conn_info%theta_b(offset + 1:offset + ntype) = conn_info%theta_b(offset + 1:offset + ntype) + natom_prev
354 conn_info%theta_c(offset + 1:offset + ntype) = conn_info%theta_c(offset + 1:offset + ntype) + natom_prev
355 END IF
356 ! BONDANGLE SECTION
357 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the BONDANGLE section'
358 label = trim(avail_section(12))
359 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
360 IF (found) THEN
361 CALL parser_get_next_line(parser, 1)
362 CALL parser_get_object(parser, ntype)
363 offset = 0
364 IF (ASSOCIATED(conn_info%theta_a)) offset = SIZE(conn_info%theta_a)
365 CALL reallocate(conn_info%theta_a, 1, offset + ntype)
366 CALL reallocate(conn_info%theta_b, 1, offset + ntype)
367 CALL reallocate(conn_info%theta_c, 1, offset + ntype)
368 CALL reallocate(conn_info%theta_type, 1, offset + ntype)
369 DO itype = 1, ntype
370 CALL parser_get_next_line(parser, 1)
371 CALL parser_get_object(parser, conn_info%theta_a(offset + itype))
372 CALL parser_get_object(parser, conn_info%theta_b(offset + itype))
373 CALL parser_get_object(parser, conn_info%theta_c(offset + itype))
374 CALL parser_get_object(parser, itemp)
375 conn_info%theta_type(offset + itype) = itemp
376 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT BONDANGLE INFO HERE!!!!"
377 END DO
378 conn_info%theta_a(offset + 1:offset + ntype) = conn_info%theta_a(offset + 1:offset + ntype) + natom_prev
379 conn_info%theta_b(offset + 1:offset + ntype) = conn_info%theta_b(offset + 1:offset + ntype) + natom_prev
380 conn_info%theta_c(offset + 1:offset + ntype) = conn_info%theta_c(offset + 1:offset + ntype) + natom_prev
381 END IF
382 ! IMPDIHEDRALH SECTION
383 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the IMPDIHEDRALH section'
384 label = trim(avail_section(14))
385 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
386 IF (found) THEN
387 CALL parser_get_next_line(parser, 1)
388 CALL parser_get_object(parser, ntype)
389 offset = 0
390 IF (ASSOCIATED(conn_info%impr_a)) offset = SIZE(conn_info%impr_a)
391 CALL reallocate(conn_info%impr_a, 1, offset + ntype)
392 CALL reallocate(conn_info%impr_b, 1, offset + ntype)
393 CALL reallocate(conn_info%impr_c, 1, offset + ntype)
394 CALL reallocate(conn_info%impr_d, 1, offset + ntype)
395 CALL reallocate(conn_info%impr_type, 1, offset + ntype)
396 DO itype = 1, ntype
397 CALL parser_get_next_line(parser, 1)
398 CALL parser_get_object(parser, conn_info%impr_a(offset + itype))
399 CALL parser_get_object(parser, conn_info%impr_b(offset + itype))
400 CALL parser_get_object(parser, conn_info%impr_c(offset + itype))
401 CALL parser_get_object(parser, conn_info%impr_d(offset + itype))
402 CALL parser_get_object(parser, itemp)
403 conn_info%impr_type(offset + itype) = itemp
404 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT IMPDIHEDRALH INFO HERE!!!!"
405 END DO
406 conn_info%impr_a(offset + 1:offset + ntype) = conn_info%impr_a(offset + 1:offset + ntype) + natom_prev
407 conn_info%impr_b(offset + 1:offset + ntype) = conn_info%impr_b(offset + 1:offset + ntype) + natom_prev
408 conn_info%impr_c(offset + 1:offset + ntype) = conn_info%impr_c(offset + 1:offset + ntype) + natom_prev
409 conn_info%impr_d(offset + 1:offset + ntype) = conn_info%impr_d(offset + 1:offset + ntype) + natom_prev
410 END IF
411 ! IMPDIHEDRAL SECTION
412 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the IMPDIHEDRAL section'
413 label = trim(avail_section(15))
414 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
415 IF (found) THEN
416 CALL parser_get_next_line(parser, 1)
417 CALL parser_get_object(parser, ntype)
418 offset = 0
419 IF (ASSOCIATED(conn_info%impr_a)) offset = SIZE(conn_info%impr_a)
420 CALL reallocate(conn_info%impr_a, 1, offset + ntype)
421 CALL reallocate(conn_info%impr_b, 1, offset + ntype)
422 CALL reallocate(conn_info%impr_c, 1, offset + ntype)
423 CALL reallocate(conn_info%impr_d, 1, offset + ntype)
424 CALL reallocate(conn_info%impr_type, 1, offset + ntype)
425 DO itype = 1, ntype
426 CALL parser_get_next_line(parser, 1)
427 CALL parser_get_object(parser, conn_info%impr_a(offset + itype))
428 CALL parser_get_object(parser, conn_info%impr_b(offset + itype))
429 CALL parser_get_object(parser, conn_info%impr_c(offset + itype))
430 CALL parser_get_object(parser, conn_info%impr_d(offset + itype))
431 CALL parser_get_object(parser, itemp)
432 conn_info%impr_type(offset + itype) = itemp
433 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT IMPDIHEDRAL INFO HERE!!!!"
434 END DO
435 conn_info%impr_a(offset + 1:offset + ntype) = conn_info%impr_a(offset + 1:offset + ntype) + natom_prev
436 conn_info%impr_b(offset + 1:offset + ntype) = conn_info%impr_b(offset + 1:offset + ntype) + natom_prev
437 conn_info%impr_c(offset + 1:offset + ntype) = conn_info%impr_c(offset + 1:offset + ntype) + natom_prev
438 conn_info%impr_d(offset + 1:offset + ntype) = conn_info%impr_d(offset + 1:offset + ntype) + natom_prev
439 END IF
440 ! DIHEDRALH SECTION
441 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the DIHEDRALH section'
442 label = trim(avail_section(17))
443 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
444 IF (found) THEN
445 CALL parser_get_next_line(parser, 1)
446 CALL parser_get_object(parser, ntype)
447 offset = 0
448 IF (ASSOCIATED(conn_info%phi_a)) offset = SIZE(conn_info%phi_a)
449 CALL reallocate(conn_info%phi_a, 1, offset + ntype)
450 CALL reallocate(conn_info%phi_b, 1, offset + ntype)
451 CALL reallocate(conn_info%phi_c, 1, offset + ntype)
452 CALL reallocate(conn_info%phi_d, 1, offset + ntype)
453 CALL reallocate(conn_info%phi_type, 1, offset + ntype)
454 DO itype = 1, ntype
455 CALL parser_get_next_line(parser, 1)
456 CALL parser_get_object(parser, conn_info%phi_a(offset + itype))
457 CALL parser_get_object(parser, conn_info%phi_b(offset + itype))
458 CALL parser_get_object(parser, conn_info%phi_c(offset + itype))
459 CALL parser_get_object(parser, conn_info%phi_d(offset + itype))
460 CALL parser_get_object(parser, itemp)
461 conn_info%phi_type(offset + itype) = itemp
462 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT DIHEDRALH INFO HERE!!!!"
463 END DO
464 conn_info%phi_a(offset + 1:offset + ntype) = conn_info%phi_a(offset + 1:offset + ntype) + natom_prev
465 conn_info%phi_b(offset + 1:offset + ntype) = conn_info%phi_b(offset + 1:offset + ntype) + natom_prev
466 conn_info%phi_c(offset + 1:offset + ntype) = conn_info%phi_c(offset + 1:offset + ntype) + natom_prev
467 conn_info%phi_d(offset + 1:offset + ntype) = conn_info%phi_d(offset + 1:offset + ntype) + natom_prev
468 END IF
469 ! DIHEDRAL SECTION
470 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the DIHEDRAL section'
471 label = trim(avail_section(18))
472 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
473 IF (found) THEN
474 CALL parser_get_next_line(parser, 1)
475 CALL parser_get_object(parser, ntype)
476 offset = 0
477 IF (ASSOCIATED(conn_info%phi_a)) offset = SIZE(conn_info%phi_a)
478 CALL reallocate(conn_info%phi_a, 1, offset + ntype)
479 CALL reallocate(conn_info%phi_b, 1, offset + ntype)
480 CALL reallocate(conn_info%phi_c, 1, offset + ntype)
481 CALL reallocate(conn_info%phi_d, 1, offset + ntype)
482 CALL reallocate(conn_info%phi_type, 1, offset + ntype)
483 DO itype = 1, ntype
484 CALL parser_get_next_line(parser, 1)
485 CALL parser_get_object(parser, conn_info%phi_a(offset + itype))
486 CALL parser_get_object(parser, conn_info%phi_b(offset + itype))
487 CALL parser_get_object(parser, conn_info%phi_c(offset + itype))
488 CALL parser_get_object(parser, conn_info%phi_d(offset + itype))
489 CALL parser_get_object(parser, itemp)
490 conn_info%phi_type(offset + itype) = itemp
491 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT DIHEDRAL INFO HERE!!!!"
492 END DO
493 conn_info%phi_a(offset + 1:offset + ntype) = conn_info%phi_a(offset + 1:offset + ntype) + natom_prev
494 conn_info%phi_b(offset + 1:offset + ntype) = conn_info%phi_b(offset + 1:offset + ntype) + natom_prev
495 conn_info%phi_c(offset + 1:offset + ntype) = conn_info%phi_c(offset + 1:offset + ntype) + natom_prev
496 conn_info%phi_d(offset + 1:offset + ntype) = conn_info%phi_d(offset + 1:offset + ntype) + natom_prev
497 END IF
498 CALL parser_release(parser)
499
500 ! SOLVENTATOM and SOLVENTCONSTR SECTION
501 IF (iw > 0) WRITE (iw, '(T2,A)') 'GTOP_INFO| Parsing the SOLVENTATOM section'
502 nsolvent = (SIZE(atom_info%r(1, :)) - nsolute)/3
503
504 NULLIFY (na, am, ac, ba, bb)
505 CALL parser_create(parser, file_name, para_env=para_env)
506 label = trim(avail_section(20))
507 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
508 IF (found) THEN
509 CALL parser_get_next_line(parser, 1)
510 CALL parser_get_object(parser, natom)
511 CALL reallocate(na, 1, natom)
512 CALL reallocate(am, 1, natom)
513 CALL reallocate(ac, 1, natom)
514 DO iatom = 1, natom
515 CALL parser_get_next_line(parser, 1)
516 CALL parser_get_object(parser, itemp)
517 CALL parser_get_object(parser, string)
518 CALL parser_get_object(parser, na(iatom))
519 CALL parser_get_object(parser, am(iatom))
520 CALL parser_get_object(parser, ac(iatom))
521 IF (iw > 0) WRITE (iw, *) "GTOP_INFO| PUT SOLVENTATOM INFO HERE!!!!"
522 END DO
523 END IF
524 label = trim(avail_section(21))
525 ncon = 0
526 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
527 IF (found) THEN
528 CALL parser_get_next_line(parser, 1)
529 CALL parser_get_object(parser, ncon)
530 CALL reallocate(ba, 1, ncon)
531 CALL reallocate(bb, 1, ncon)
532 DO icon = 1, ncon
533 CALL parser_get_next_line(parser, 1)
534 CALL parser_get_object(parser, ba(icon))
535 CALL parser_get_object(parser, bb(icon))
536 END DO
537 END IF
538 CALL parser_release(parser)
539
540 offset = 0
541 IF (ASSOCIATED(atom_info%id_molname)) offset = SIZE(atom_info%id_molname)
542 CALL reallocate(atom_info%id_molname, 1, offset + nsolvent*natom)
543 CALL reallocate(atom_info%resid, 1, offset + nsolvent*natom)
544 CALL reallocate(atom_info%id_resname, 1, offset + nsolvent*natom)
545 CALL reallocate(atom_info%id_atmname, 1, offset + nsolvent*natom)
546 CALL reallocate(atom_info%id_element, 1, offset + nsolvent*natom)
547 CALL reallocate(atom_info%atm_charge, 1, offset + nsolvent*natom)
548 CALL reallocate(atom_info%atm_mass, 1, offset + nsolvent*natom)
549 DO isolvent = 1, nsolvent
550 offset = nsolute + natom*isolvent - natom
551 DO iatom = 1, natom
552 index_now = offset + iatom
553 atom_info%id_atmname(index_now) = str2id(s2s(namearray1(na(iatom))))
554 atom_info%id_element(index_now) = str2id(s2s(namearray1(na(iatom))))
555 atom_info%id_molname(index_now) = str2id(s2s("SOL"))
556 atom_info%id_resname(index_now) = str2id(s2s("SOL"))
557 atom_info%resid(index_now) = isolvent
558 atom_info%atm_mass(index_now) = am(iatom)
559 atom_info%atm_charge(index_now) = ac(iatom)
560 END DO
561 END DO
562
563 offset = 0
564 IF (ASSOCIATED(conn_info%bond_a)) offset = SIZE(conn_info%bond_a)
565 offset2 = maxval(conn_info%bond_type(:))
566 CALL reallocate(conn_info%bond_a, 1, offset + ncon*nsolvent)
567 CALL reallocate(conn_info%bond_b, 1, offset + ncon*nsolvent)
568 CALL reallocate(conn_info%bond_type, 1, offset + ncon*nsolvent)
569 offset = offset - ncon
570 DO isolvent = 1, nsolvent
571 offset = offset + ncon
572 DO icon = 1, ncon
573 conn_info%bond_a(offset + icon) = nsolute + isolvent*ncon - ncon + ba(icon)
574 conn_info%bond_b(offset + icon) = nsolute + isolvent*ncon - ncon + bb(icon)
575 conn_info%bond_type(offset + icon) = offset2 + isolvent*ncon - ncon + icon
576 END DO
577 END DO
578 ! PARA_RES structure
579 i = 0
580 nbond_prev = 0
581 IF (ASSOCIATED(conn_info%c_bond_a)) i = SIZE(conn_info%c_bond_a)
582 nbond = SIZE(conn_info%bond_a)
583 DO ibond = 1 + nbond_prev, nbond + nbond_prev
584 iatom = conn_info%bond_a(ibond)
585 jatom = conn_info%bond_b(ibond)
586 IF (topology%para_res) THEN
587 IF ((atom_info%id_molname(iatom) /= atom_info%id_molname(jatom)) .OR. &
588 (atom_info%resid(iatom) /= atom_info%resid(jatom)) .OR. &
589 (atom_info%id_resname(iatom) /= atom_info%id_resname(jatom))) THEN
590 IF (iw > 0) WRITE (iw, '(T2,A,2I3)') "GTOP_INFO| PARA_RES, bond between molecules atom ", &
591 iatom, jatom
592 i = i + 1
593 CALL reallocate(conn_info%c_bond_a, 1, i)
594 CALL reallocate(conn_info%c_bond_b, 1, i)
595 CALL reallocate(conn_info%c_bond_type, 1, i)
596 conn_info%c_bond_a(i) = iatom
597 conn_info%c_bond_b(i) = jatom
598 conn_info%c_bond_type(i) = conn_info%bond_type(ibond)
599 END IF
600 ELSE
601 IF (atom_info%id_molname(iatom) /= atom_info%id_molname(jatom)) THEN
602 cpabort("")
603 END IF
604 END IF
605 END DO
606
607 DEALLOCATE (namearray1)
608 DEALLOCATE (namearray2)
609 DEALLOCATE (na)
610 DEALLOCATE (am)
611 DEALLOCATE (ac)
612 IF (ASSOCIATED(ba)) &
613 DEALLOCATE (ba)
614 IF (ASSOCIATED(bb)) &
615 DEALLOCATE (bb)
616
617 CALL timestop(handle)
618 CALL cp_print_key_finished_output(iw, logger, subsys_section, &
619 "PRINT%TOPOLOGY_INFO/GTOP_INFO")
620
621 END SUBROUTINE read_topology_gromos
622
623! **************************************************************************************************
624!> \brief ...
625!> \param topology ...
626!> \param para_env ...
627!> \param subsys_section ...
628! **************************************************************************************************
629 SUBROUTINE read_coordinate_g96(topology, para_env, subsys_section)
631 TYPE(mp_para_env_type), POINTER :: para_env
632 TYPE(section_vals_type), POINTER :: subsys_section
633
634 CHARACTER(len=*), PARAMETER :: routinen = 'read_coordinate_g96'
635 INTEGER, PARAMETER :: nblock = 1000
636
637 CHARACTER(LEN=default_string_length) :: label, string, strtmp, strtmp2
638 CHARACTER(LEN=default_string_length), DIMENSION(5) :: avail_section
639 INTEGER :: handle, itemp, iw, natom, newsize
640 LOGICAL :: found
641 REAL(kind=dp) :: pfactor
642 REAL(kind=dp), DIMENSION(:, :), POINTER :: velocity
643 TYPE(atom_info_type), POINTER :: atom_info
644 TYPE(cp_logger_type), POINTER :: logger
645 TYPE(cp_parser_type) :: parser
646 TYPE(section_vals_type), POINTER :: velocity_section
647
648 NULLIFY (logger, velocity)
649 logger => cp_get_default_logger()
650 iw = cp_print_key_unit_nr(logger, subsys_section, "PRINT%TOPOLOGY_INFO/G96_INFO", &
651 extension=".subsysLog")
652 CALL timeset(routinen, handle)
653
654 pfactor = section_get_rval(subsys_section, "TOPOLOGY%MEMORY_PROGRESSION_FACTOR")
655 atom_info => topology%atom_info
656 IF (iw > 0) WRITE (iw, *) " Reading in G96 file ", trim(topology%coord_file_name)
657 avail_section(1) = "TITLE"
658 avail_section(2) = "TIMESTEP"
659 avail_section(3) = "POSITION"
660 avail_section(4) = "VELOCITY"
661 avail_section(5) = "BOX"
662
663 ! Element is assigned on the basis of the atm_name
664 topology%aa_element = .true.
665
666 CALL reallocate(atom_info%id_molname, 1, nblock)
667 CALL reallocate(atom_info%id_resname, 1, nblock)
668 CALL reallocate(atom_info%resid, 1, nblock)
669 CALL reallocate(atom_info%id_atmname, 1, nblock)
670 CALL reallocate(atom_info%id_element, 1, nblock)
671 CALL reallocate(atom_info%r, 1, 3, 1, nblock)
672 CALL reallocate(atom_info%atm_mass, 1, nblock)
673 CALL reallocate(atom_info%atm_charge, 1, nblock)
674 CALL reallocate(atom_info%occup, 1, nblock)
675 CALL reallocate(atom_info%beta, 1, nblock)
676 ! TITLE SECTION
677 IF (iw > 0) WRITE (iw, '(T2,A)') 'G96_INFO| Parsing the TITLE section'
678 CALL parser_create(parser, topology%coord_file_name, para_env=para_env)
679 label = trim(avail_section(1))
680 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
681 IF (found) THEN
682 DO
683 CALL parser_get_next_line(parser, 1)
684 CALL parser_get_object(parser, string, string_length=default_string_length)
685 IF (string == trim("END")) EXIT
686 IF (iw > 0) WRITE (iw, *) "G96_INFO| ", trim(string)
687 END DO
688 END IF
689 CALL parser_release(parser)
690 ! POSITION SECTION
691 IF (iw > 0) WRITE (iw, '(T2,A)') 'G96_INFO| Parsing the POSITION section'
692 CALL parser_create(parser, topology%coord_file_name, para_env=para_env)
693 label = trim(avail_section(3))
694 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
695 IF (found) THEN
696 natom = 0
697 CALL parser_get_next_line(parser, 1)
698 CALL parser_get_object(parser, string, string_length=default_string_length)
699 DO
700 IF (string == trim("END")) EXIT
701 natom = natom + 1
702 IF (natom > SIZE(atom_info%id_molname)) THEN
703 newsize = int(pfactor*natom)
704 CALL reallocate(atom_info%id_molname, 1, newsize)
705 CALL reallocate(atom_info%id_resname, 1, newsize)
706 CALL reallocate(atom_info%resid, 1, newsize)
707 CALL reallocate(atom_info%id_atmname, 1, newsize)
708 CALL reallocate(atom_info%id_element, 1, newsize)
709 CALL reallocate(atom_info%r, 1, 3, 1, newsize)
710 CALL reallocate(atom_info%atm_mass, 1, newsize)
711 CALL reallocate(atom_info%atm_charge, 1, newsize)
712 CALL reallocate(atom_info%occup, 1, newsize)
713 CALL reallocate(atom_info%beta, 1, newsize)
714 END IF
715 READ (string, *) &
716 atom_info%resid(natom), strtmp, strtmp2, &
717 itemp, atom_info%r(1, natom), atom_info%r(2, natom), atom_info%r(3, natom)
718 atom_info%id_resname(natom) = str2id(s2s(strtmp))
719 atom_info%id_atmname(natom) = str2id(s2s(strtmp2))
720 atom_info%id_molname(natom) = atom_info%id_resname(natom)
721 atom_info%id_element(natom) = atom_info%id_atmname(natom)
722 atom_info%r(1, natom) = cp_unit_to_cp2k(atom_info%r(1, natom), "nm")
723 atom_info%r(2, natom) = cp_unit_to_cp2k(atom_info%r(2, natom), "nm")
724 atom_info%r(3, natom) = cp_unit_to_cp2k(atom_info%r(3, natom), "nm")
725 IF (iw > 0) WRITE (iw, *) "G96_INFO| PUT POSITION INFO HERE!!!!"
726 CALL parser_get_next_line(parser, 1)
727 CALL parser_get_object(parser, string, string_length=default_string_length)
728 END DO
729 END IF
730 CALL parser_release(parser)
731 CALL reallocate(velocity, 1, 3, 1, natom)
732 ! VELOCITY SECTION
733 IF (topology%use_g96_velocity) THEN
734 IF (iw > 0) WRITE (iw, '(T2,A)') 'G96_INFO| Parsing the VELOCITY section'
735 CALL parser_create(parser, topology%coord_file_name, para_env=para_env)
736 label = trim(avail_section(4))
737 CALL parser_search_string(parser, label, .true., found, begin_line=.true.)
738 IF (found) THEN
739 natom = 0
740 CALL parser_get_next_line(parser, 1)
741 CALL parser_get_object(parser, string, string_length=default_string_length)
742 DO
743 IF (string == trim("END")) EXIT
744 natom = natom + 1
745 READ (string, *) &
746 atom_info%resid(natom), strtmp, strtmp2, &
747 itemp, velocity(1, natom), velocity(2, natom), velocity(3, natom)
748 atom_info%id_resname(natom) = str2id(strtmp)
749 atom_info%id_atmname(natom) = str2id(strtmp2)
750 atom_info%id_molname(natom) = atom_info%id_resname(natom)
751 atom_info%id_element(natom) = atom_info%id_atmname(natom)
752 velocity(1, natom) = cp_unit_to_cp2k(velocity(1, natom), "nm*ps^-1")
753 velocity(2, natom) = cp_unit_to_cp2k(velocity(2, natom), "nm*ps^-1")
754 velocity(3, natom) = cp_unit_to_cp2k(velocity(3, natom), "nm*ps^-1")
755 IF (iw > 0) WRITE (iw, *) "G96_INFO| PUT VELOCITY INFO HERE!!!!"
756 CALL parser_get_next_line(parser, 1)
757 CALL parser_get_object(parser, string, string_length=default_string_length)
758 END DO
759 CALL parser_release(parser)
760 velocity_section => section_vals_get_subs_vals(subsys_section, "VELOCITY")
761 CALL section_velocity_val_set(velocity_section, velocity=velocity, &
762 conv_factor=1.0_dp)
763 END IF
764 END IF
765 DEALLOCATE (velocity)
766
767 CALL reallocate(atom_info%id_molname, 1, natom)
768 CALL reallocate(atom_info%id_resname, 1, natom)
769 CALL reallocate(atom_info%resid, 1, natom)
770 CALL reallocate(atom_info%id_atmname, 1, natom)
771 CALL reallocate(atom_info%id_element, 1, natom)
772 CALL reallocate(atom_info%r, 1, 3, 1, natom)
773 CALL reallocate(atom_info%atm_mass, 1, natom)
774 CALL reallocate(atom_info%atm_charge, 1, natom)
775 CALL reallocate(atom_info%occup, 1, natom)
776 CALL reallocate(atom_info%beta, 1, natom)
777
778 IF (.NOT. topology%para_res) atom_info%resid(:) = 1
779
780 topology%natoms = natom
781 CALL cp_print_key_finished_output(iw, logger, subsys_section, &
782 "PRINT%TOPOLOGY_INFO/G96_INFO")
783 CALL timestop(handle)
784
785 END SUBROUTINE read_coordinate_g96
786
787END MODULE topology_gromos
788
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,...
Utility routines to read data from files. Kept as close as possible to the old parser because.
subroutine, public parser_get_next_line(parser, nline, at_end)
Read the next input line and broadcast the input information. Skip (nline-1) lines and skip also all ...
subroutine, public parser_search_string(parser, string, ignore_case, found, line, begin_line, search_from_begin_of_file)
Search a string pattern in a file defined by its logical unit number "unit". A case sensitive search ...
Utility routines to read data from files. Kept as close as possible to the old parser because.
subroutine, public parser_release(parser)
releases the parser
subroutine, public parser_create(parser, file_name, unit_nr, para_env, end_section_label, separator_chars, comment_char, continuation_char, quote_char, section_char, parse_white_lines, initial_variables, apply_preprocessing)
Start a parser run. Initial variables allow to @SET stuff before opening the file.
unit conversion facility
Definition cp_units.F:30
real(kind=dp) function, public cp_unit_to_cp2k(value, unit_str, defaults, power)
converts to the internal cp2k units to the given unit
Definition cp_units.F:1150
subroutine, public section_velocity_val_set(velocity_section, particles, velocity, conv_factor)
routine to dump velocities.. fast implementation
objects that represent the structure of input sections and the data contained in an input section
real(kind=dp) function, public section_get_rval(section_vals, keyword_name)
...
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
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
Utility routines for the memory handling.
Interface to the message passing library MPI.
generates a unique id number for a string (str2id) that can be used two compare two strings....
character(len=default_string_length) function, public s2s(str)
converts a string in a string of default_string_length
integer function, public str2id(str)
returns a unique id for a given string, and stores the string for later retrieval using the id.
subroutine, public read_coordinate_g96(topology, para_env, subsys_section)
...
subroutine, public read_topology_gromos(file_name, topology, para_env, subsys_section)
Read GROMOS topology file.
Control for reading in different topologies and coordinates.
Definition topology.F:13
type of a logger, at the moment it contains just a print level starting at which level it should be l...
stores all the informations relevant to an mpi environment