(git:e7e05ae)
topology_cp2k.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 ! **************************************************************************************************
9 
10  USE cell_types, ONLY: cell_type,&
13  cp_logger_type
16  USE cp_parser_methods, ONLY: parser_get_object,&
19  USE cp_parser_types, ONLY: cp_parser_type,&
22  USE cp_units, ONLY: cp_unit_to_cp2k
27  section_vals_type,&
30  USE kinds, ONLY: default_string_length,&
31  dp,&
33  USE memory_utilities, ONLY: reallocate
34  USE message_passing, ONLY: mp_para_env_type
35  USE periodic_table, ONLY: nelem,&
36  ptable
37  USE string_table, ONLY: id2str,&
38  s2s,&
39  str2id
40  USE topology_types, ONLY: atom_info_type,&
42 #include "./base/base_uses.f90"
43 
44  IMPLICIT NONE
45 
46  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'topology_cp2k'
47 
48  PRIVATE
49  PUBLIC :: read_coordinate_cp2k
50 
51 CONTAINS
52 
53 ! **************************************************************************************************
54 !> \brief Read the CP2K &COORD section from an external file, i.e. read
55 !> atomic coordinates and molecule/residue information in CP2K format.
56 !> \param topology ...
57 !> \param para_env ...
58 !> \param subsys_section ...
59 !> \date 17.01.2011 (Creation, MK)
60 !> \author Matthias Krack (MK)
61 !> \version 1.0
62 ! **************************************************************************************************
63  SUBROUTINE read_coordinate_cp2k(topology, para_env, subsys_section)
64 
66  TYPE(mp_para_env_type), POINTER :: para_env
67  TYPE(section_vals_type), POINTER :: subsys_section
68 
69  CHARACTER(LEN=*), PARAMETER :: routinen = 'read_coordinate_cp2k'
70 
71  CHARACTER(LEN=default_string_length) :: string
72  CHARACTER(LEN=max_line_length) :: error_message
73  INTEGER :: handle, i, ian, iw, natom, newsize, &
74  number_of_atoms
75  LOGICAL :: eof, explicit, scaled_coordinates
76  REAL(kind=dp) :: pfactor, unit_conv
77  REAL(kind=dp), DIMENSION(3) :: r
78  TYPE(atom_info_type), POINTER :: atom_info
79  TYPE(cell_type), POINTER :: cell
80  TYPE(cp_logger_type), POINTER :: logger
81  TYPE(cp_parser_type), POINTER :: parser
82  TYPE(section_vals_type), POINTER :: coord_section
83 
84  CALL timeset(routinen, handle)
85 
86  NULLIFY (coord_section)
87  NULLIFY (logger)
88  NULLIFY (parser)
89  logger => cp_get_default_logger()
90  iw = cp_print_key_unit_nr(logger, subsys_section, "PRINT%TOPOLOGY_INFO/XYZ_INFO", &
91  extension=".subsysLog")
92 
93  ! Check if there is a &COORD section
94  coord_section => section_vals_get_subs_vals(subsys_section, "COORD")
95  CALL section_vals_get(coord_section, explicit=explicit)
96  IF (explicit) THEN
97  CALL section_vals_val_get(coord_section, "UNIT", c_val=string)
98  CALL section_vals_val_get(coord_section, "SCALED", l_val=scaled_coordinates)
99  ELSE
100  ! The default is Cartesian coordinates in Angstrom
101  scaled_coordinates = .false.
102  string = "angstrom"
103  END IF
104  unit_conv = cp_unit_to_cp2k(1.0_dp, trim(string))
105 
106  atom_info => topology%atom_info
107  cell => topology%cell_muc
108 
109  IF (iw > 0) THEN
110  WRITE (unit=iw, fmt="(T2,A)") &
111  "BEGIN of COORD section data read from file "//trim(topology%coord_file_name)
112  END IF
113 
114  pfactor = section_get_rval(subsys_section, "TOPOLOGY%MEMORY_PROGRESSION_FACTOR")
115  number_of_atoms = section_get_ival(subsys_section, "TOPOLOGY%NUMBER_OF_ATOMS")
116  IF (number_of_atoms < 1) THEN
117  newsize = 1000
118  ELSE
119  newsize = number_of_atoms
120  END IF
121 
122  CALL reallocate(atom_info%id_molname, 1, newsize)
123  CALL reallocate(atom_info%id_resname, 1, newsize)
124  CALL reallocate(atom_info%resid, 1, newsize)
125  CALL reallocate(atom_info%id_atmname, 1, newsize)
126  CALL reallocate(atom_info%r, 1, 3, 1, newsize)
127  CALL reallocate(atom_info%atm_mass, 1, newsize)
128  CALL reallocate(atom_info%atm_charge, 1, newsize)
129  CALL reallocate(atom_info%occup, 1, newsize)
130  CALL reallocate(atom_info%beta, 1, newsize)
131  CALL reallocate(atom_info%id_element, 1, newsize)
132 
133  topology%molname_generated = .false.
134  ! Element is assigned on the basis of the atm_name
135  topology%aa_element = .true.
136 
137  ALLOCATE (parser)
138  CALL parser_create(parser, topology%coord_file_name, para_env=para_env)
139 
140  natom = 0
141  DO
142  CALL parser_get_object(parser, object=string, newline=.true., at_end=eof)
143  IF (eof) EXIT
144  natom = natom + 1
145  IF (natom > SIZE(atom_info%id_atmname)) THEN
146  newsize = int(pfactor*natom)
147  CALL reallocate(atom_info%id_molname, 1, newsize)
148  CALL reallocate(atom_info%id_resname, 1, newsize)
149  CALL reallocate(atom_info%resid, 1, newsize)
150  CALL reallocate(atom_info%id_atmname, 1, newsize)
151  CALL reallocate(atom_info%r, 1, 3, 1, newsize)
152  CALL reallocate(atom_info%atm_mass, 1, newsize)
153  CALL reallocate(atom_info%atm_charge, 1, newsize)
154  CALL reallocate(atom_info%occup, 1, newsize)
155  CALL reallocate(atom_info%beta, 1, newsize)
156  CALL reallocate(atom_info%id_element, 1, newsize)
157  END IF
158  error_message = ""
159  CALL read_integer_object(string, ian, error_message)
160  IF (len_trim(error_message) == 0) THEN
161  ! Integer value found: assume atomic number, check it, and load
162  ! the corresponding element symbol if valid
163  IF ((ian < 0) .OR. (ian > nelem)) THEN
164  error_message = "Invalid atomic number <"//trim(string)// &
165  "> found in the xyz file <"//trim(topology%coord_file_name)//">!"
166  cpabort(trim(error_message))
167  ELSE
168  atom_info%id_atmname(natom) = str2id(s2s(ptable(ian)%symbol))
169  END IF
170  ELSE
171  atom_info%id_atmname(natom) = str2id(s2s(string))
172  END IF
173  ! Read x, y, and z coordinate of the current atom
174  DO i = 1, 3
175  CALL parser_get_object(parser, object=r(i))
176  END DO
177  IF (scaled_coordinates) THEN
178  CALL scaled_to_real(atom_info%r(1:3, natom), r, cell)
179  ELSE
180  atom_info%r(1:3, natom) = r(1:3)*unit_conv
181  END IF
182  IF (parser_test_next_token(parser) /= "EOL") THEN
183  CALL parser_get_object(parser, object=string)
184  atom_info%id_molname(natom) = str2id(s2s(string))
185  IF (parser_test_next_token(parser) /= "EOL") THEN
186  CALL parser_get_object(parser, object=string)
187  atom_info%id_resname(natom) = str2id(s2s(string))
188  ELSE
189  atom_info%id_resname(natom) = atom_info%id_molname(natom)
190  END IF
191  ELSE
192  string = ""
193  WRITE (unit=string, fmt="(I0)") natom
194  atom_info%id_molname(natom) = str2id(s2s(trim(id2str(atom_info%id_atmname(natom)))//trim(string)))
195  atom_info%id_resname(natom) = atom_info%id_molname(natom)
196  topology%molname_generated = .true.
197  END IF
198  atom_info%resid(natom) = 1
199  atom_info%id_element(natom) = atom_info%id_atmname(natom)
200  atom_info%atm_mass(natom) = huge(0.0_dp)
201  atom_info%atm_charge(natom) = -huge(0.0_dp)
202  IF (iw > 0) THEN
203  WRITE (unit=iw, fmt="(T2,A,3F20.8,2(2X,A))") &
204  trim(id2str(atom_info%id_atmname(natom))), atom_info%r(1:3, natom), &
205  adjustl(trim(id2str(atom_info%id_molname(natom)))), &
206  adjustl(trim(id2str(atom_info%id_resname(natom))))
207  END IF
208  IF (natom == number_of_atoms) EXIT
209  END DO
210 
211  CALL parser_release(parser)
212  DEALLOCATE (parser)
213 
214  topology%natoms = natom
215  CALL reallocate(atom_info%id_molname, 1, natom)
216  CALL reallocate(atom_info%id_resname, 1, natom)
217  CALL reallocate(atom_info%resid, 1, natom)
218  CALL reallocate(atom_info%id_atmname, 1, natom)
219  CALL reallocate(atom_info%r, 1, 3, 1, natom)
220  CALL reallocate(atom_info%atm_mass, 1, natom)
221  CALL reallocate(atom_info%atm_charge, 1, natom)
222  CALL reallocate(atom_info%occup, 1, natom)
223  CALL reallocate(atom_info%beta, 1, natom)
224  CALL reallocate(atom_info%id_element, 1, natom)
225 
226  CALL section_vals_val_set(subsys_section, "TOPOLOGY%NUMBER_OF_ATOMS", i_val=natom)
227 
228  IF (iw > 0) THEN
229  WRITE (unit=iw, fmt="(T2,A)") &
230  "END of COORD section data read from file "//trim(topology%coord_file_name)
231  END IF
232 
233  CALL cp_print_key_finished_output(iw, logger, subsys_section, &
234  "PRINT%TOPOLOGY_INFO/XYZ_INFO")
235 
236  CALL timestop(handle)
237 
238  END SUBROUTINE read_coordinate_cp2k
239 
240 END MODULE topology_cp2k
Handles all functions related to the CELL.
Definition: cell_types.F:15
subroutine, public scaled_to_real(r, s, cell)
Transform scaled cell coordinates real coordinates. r=h*s.
Definition: cell_types.F:516
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.
elemental subroutine, public read_integer_object(string, object, error_message)
Returns an integer number read from a string including products of integer numbers like iz1*iz2*iz3.
character(len=3) function, public parser_test_next_token(parser, string_length)
Test next input object.
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
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)
...
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
integer function, public section_get_ival(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
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
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
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public max_line_length
Definition: kinds.F:59
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.
Periodic Table related data definitions.
type(atom), dimension(0:nelem), public ptable
integer, parameter, public nelem
generates a unique id number for a string (str2id) that can be used two compare two strings....
Definition: string_table.F:22
character(len=default_string_length) function, public s2s(str)
converts a string in a string of default_string_length
Definition: string_table.F:141
integer function, public str2id(str)
returns a unique id for a given string, and stores the string for later retrieval using the id.
Definition: string_table.F:72
character(len=default_string_length) function, public id2str(id)
returns the string associated with a given id
Definition: string_table.F:115
subroutine, public read_coordinate_cp2k(topology, para_env, subsys_section)
Read the CP2K &COORD section from an external file, i.e. read atomic coordinates and molecule/residue...
Definition: topology_cp2k.F:64
Control for reading in different topologies and coordinates.
Definition: topology.F:13