51 #include "./base/base_uses.f90"
57 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'atoms_input'
72 LOGICAL,
INTENT(IN),
OPTIONAL :: overwrite
73 TYPE(section_vals_type),
POINTER :: subsys_section
74 LOGICAL,
INTENT(IN),
OPTIONAL :: save_mem
76 CHARACTER(len=*),
PARAMETER :: routinen =
'read_atoms_input'
78 CHARACTER(len=2*default_string_length) :: line_att
79 CHARACTER(len=default_string_length) :: error_message, my_default_index, strtmp, &
81 INTEGER :: default_id, end_c, handle, iatom, j, &
82 natom, output_unit, start_c, wrd
83 LOGICAL :: explicit, is_ok, my_overwrite, &
84 my_save_mem, scaled_coordinates
85 REAL(kind=
dp) :: r0(3), unit_conv
87 TYPE(cell_type),
POINTER :: cell
88 TYPE(cp_sll_val_type),
POINTER ::
list
89 TYPE(section_vals_type),
POINTER :: coord_section
90 TYPE(val_type),
POINTER :: val
92 my_overwrite = .false.
96 IF (
PRESENT(overwrite)) my_overwrite = overwrite
97 IF (
PRESENT(save_mem)) my_save_mem = save_mem
98 NULLIFY (coord_section)
101 IF (.NOT. explicit)
RETURN
103 CALL timeset(routinen, handle)
121 IF (my_overwrite)
THEN
122 cpassert(
SIZE(atom_info%r, 2) == natom)
123 CALL cp_warn(__location__, &
124 "Overwriting coordinates. Active coordinates read from &COORD section."// &
125 " Active coordinates READ from &COORD section ")
129 CALL val_get(val, c_val=line_att)
133 DO j = start_c, len(line_att)
134 IF (line_att(j:j) /=
' ')
THEN
139 end_c = len(line_att) + 1
140 DO j = start_c, len(line_att)
141 IF (line_att(j:j) ==
' ')
THEN
146 IF (len_trim(line_att(start_c:end_c - 1)) == 0) &
147 cpabort(
"incorrectly formatted line in coord section'"//line_att//
"'")
149 atom_info%id_atmname(iatom) =
str2id(
s2s(line_att(start_c:end_c - 1)))
151 READ (line_att(start_c:end_c - 1), *) atom_info%r(wrd - 1, iatom)
160 CALL reallocate(atom_info%id_molname, 1, natom)
161 CALL reallocate(atom_info%id_resname, 1, natom)
162 CALL reallocate(atom_info%resid, 1, natom)
163 CALL reallocate(atom_info%id_atmname, 1, natom)
164 CALL reallocate(atom_info%id_element, 1, natom)
165 CALL reallocate(atom_info%r, 1, 3, 1, natom)
166 CALL reallocate(atom_info%atm_mass, 1, natom)
167 CALL reallocate(atom_info%atm_charge, 1, natom)
173 CALL val_get(val, c_val=line_att)
175 atom_info%id_molname(iatom) = default_id
176 atom_info%id_resname(iatom) = default_id
177 atom_info%resid(iatom) = 1
178 atom_info%id_atmname(iatom) = default_id
179 atom_info%id_element(iatom) = default_id
184 DO j = start_c, len(line_att)
185 IF (line_att(j:j) /=
' ')
THEN
190 end_c = len(line_att) + 1
191 DO j = start_c, len(line_att)
192 IF (line_att(j:j) ==
' ')
THEN
197 IF (len_trim(line_att(start_c:end_c - 1)) == 0) &
198 CALL cp_abort(__location__, &
199 "Incorrectly formatted input line for atom "// &
200 trim(adjustl(cp_to_string(iatom)))// &
201 " found in COORD section. Input line: <"// &
202 trim(line_att)//
"> ")
205 atom_info%id_atmname(iatom) =
str2id(
s2s(line_att(start_c:end_c - 1)))
208 atom_info%r(wrd - 1, iatom), error_message)
209 IF (len_trim(error_message) /= 0) &
210 CALL cp_abort(__location__, &
211 "Incorrectly formatted input line for atom "// &
212 trim(adjustl(cp_to_string(iatom)))// &
213 " found in COORD section. "//trim(error_message)// &
214 " Input line: <"//trim(line_att)//
"> ")
216 READ (line_att(start_c:end_c - 1), *) strtmp
217 atom_info%id_molname(iatom) =
str2id(strtmp)
218 atom_info%id_resname(iatom) = atom_info%id_molname(iatom)
219 topology%molname_generated = .false.
221 READ (line_att(start_c:end_c - 1), *) strtmp
222 atom_info%id_resname(iatom) =
str2id(strtmp)
225 IF (start_c > len_trim(line_att))
EXIT
227 IF (
topology%molname_generated)
THEN
229 WRITE (my_default_index,
'(I0)') iatom
230 atom_info%id_molname(iatom) =
str2id(
s2s(trim(
id2str(atom_info%id_atmname(iatom)))//trim(my_default_index)))
231 atom_info%id_resname(iatom) = atom_info%id_molname(iatom)
233 atom_info%id_element(iatom) = atom_info%id_atmname(iatom)
234 atom_info%atm_mass(iatom) = 0.0_dp
235 atom_info%atm_charge(iatom) = -huge(0.0_dp)
243 IF (scaled_coordinates)
THEN
244 r0 = atom_info%r(:, iatom)
247 atom_info%r(:, iatom) = atom_info%r(:, iatom)*unit_conv
252 CALL timestop(handle)
266 subsys_section, core_particle_set, save_mem)
268 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set, shell_particle_set
269 TYPE(cell_type),
POINTER :: cell
270 TYPE(section_vals_type),
POINTER :: subsys_section
271 TYPE(particle_type),
DIMENSION(:),
OPTIONAL, &
272 POINTER :: core_particle_set
273 LOGICAL,
INTENT(IN),
OPTIONAL :: save_mem
275 CHARACTER(len=*),
PARAMETER :: routinen =
'read_shell_coord_input'
277 CHARACTER(len=2*default_string_length) :: line_att
278 CHARACTER(len=default_string_length) :: name_kind, unit_str
279 CHARACTER(len=default_string_length), &
280 ALLOCATABLE,
DIMENSION(:) :: at_name, at_name_c
281 INTEGER :: end_c, handle, ishell, j, nshell, &
282 output_unit, sh_index, start_c, wrd
283 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: at_index, at_index_c
284 LOGICAL :: core_scaled_coordinates, explicit, &
285 is_ok, is_shell, my_save_mem, &
286 shell_scaled_coordinates
287 REAL(kind=
dp) :: dab, mass_com, rab(3), unit_conv_core, &
289 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: r, rc
290 TYPE(atomic_kind_type),
POINTER :: atomic_kind
291 TYPE(cp_sll_val_type),
POINTER ::
list
292 TYPE(section_vals_type),
POINTER :: core_coord_section, shell_coord_section
293 TYPE(shell_kind_type),
POINTER :: shell
294 TYPE(val_type),
POINTER :: val
296 my_save_mem = .false.
297 NULLIFY (atomic_kind,
list, shell_coord_section, shell, val)
300 IF (
PRESENT(save_mem)) my_save_mem = save_mem
301 NULLIFY (shell_coord_section, core_coord_section)
304 IF (.NOT. explicit)
RETURN
306 CALL timeset(routinen, handle)
307 cpassert(
ASSOCIATED(particle_set))
318 IF (
ASSOCIATED(shell_particle_set))
THEN
319 cpassert((
SIZE(shell_particle_set, 1) == nshell))
320 ALLOCATE (r(3, nshell), at_name(nshell), at_index(nshell))
321 CALL cp_warn(__location__, &
322 "Overwriting shell coordinates. "// &
323 "Active coordinates READ from &SHELL_COORD section. ")
325 DO ishell = 1, nshell
328 CALL val_get(val, c_val=line_att)
331 DO j = start_c, len(line_att)
332 IF (line_att(j:j) /=
' ')
THEN
337 end_c = len(line_att) + 1
338 DO j = start_c, len(line_att)
339 IF (line_att(j:j) ==
' ')
THEN
344 IF (wrd /= 5 .AND. end_c >= len(line_att) + 1) &
345 cpabort(
"incorrectly formatted line in coord section'"//line_att//
"'")
347 at_name(ishell) = line_att(start_c:end_c - 1)
349 ELSE IF (wrd == 5)
THEN
350 READ (line_att(start_c:end_c - 1), *) at_index(ishell)
352 READ (line_att(start_c:end_c - 1), *) r(wrd - 1, ishell)
358 IF (
PRESENT(core_particle_set))
THEN
359 cpassert(
ASSOCIATED(core_particle_set))
369 cpassert((
SIZE(core_particle_set, 1) == nshell))
370 ALLOCATE (rc(3, nshell), at_name_c(nshell), at_index_c(nshell))
371 CALL cp_warn(__location__, &
372 "Overwriting cores coordinates. "// &
373 "Active coordinates READ from &CORE_COORD section. ")
375 DO ishell = 1, nshell
378 CALL val_get(val, c_val=line_att)
381 DO j = start_c, len(line_att)
382 IF (line_att(j:j) /=
' ')
THEN
387 end_c = len(line_att) + 1
388 DO j = start_c, len(line_att)
389 IF (line_att(j:j) ==
' ')
THEN
394 IF (wrd /= 5 .AND. end_c >= len(line_att) + 1) &
395 cpabort(
"incorrectly formatted line in coord section'"//line_att//
"'")
397 at_name_c(ishell) = line_att(start_c:end_c - 1)
399 ELSE IF (wrd == 5)
THEN
400 READ (line_att(start_c:end_c - 1), *) at_index_c(ishell)
402 READ (line_att(start_c:end_c - 1), *) rc(wrd - 1, ishell)
414 DO ishell = 1, nshell
415 atomic_kind => particle_set(at_index(ishell))%atomic_kind
417 name=name_kind, shell_active=is_shell, mass=mass_com, shell=shell)
419 IF ((trim(at_name(ishell)) == trim(name_kind)) .AND. is_shell)
THEN
420 sh_index = particle_set(at_index(ishell))%shell_index
421 IF (shell_scaled_coordinates)
THEN
422 CALL scaled_to_real(r(:, ishell), shell_particle_set(sh_index)%r(:), cell)
424 shell_particle_set(sh_index)%r(:) = r(:, ishell)*unit_conv_shell
426 shell_particle_set(sh_index)%atom_index = at_index(ishell)
428 IF (
PRESENT(core_particle_set) .AND. .NOT. explicit)
THEN
429 core_particle_set(sh_index)%r(1) = (mass_com*particle_set(at_index(ishell))%r(1) - &
430 shell%mass_shell*shell_particle_set(sh_index)%r(1))/shell%mass_core
431 core_particle_set(sh_index)%r(2) = (mass_com*particle_set(at_index(ishell))%r(2) - &
432 shell%mass_shell*shell_particle_set(sh_index)%r(2))/shell%mass_core
433 core_particle_set(sh_index)%r(3) = (mass_com*particle_set(at_index(ishell))%r(3) - &
434 shell%mass_shell*shell_particle_set(sh_index)%r(3))/shell%mass_core
435 core_particle_set(sh_index)%atom_index = at_index(ishell)
436 rab =
pbc(shell_particle_set(sh_index)%r, core_particle_set(sh_index)%r, cell)
437 ELSE IF (explicit)
THEN
438 IF (core_scaled_coordinates)
THEN
439 CALL scaled_to_real(rc(:, ishell), core_particle_set(sh_index)%r(:), cell)
441 core_particle_set(sh_index)%r(:) = rc(:, ishell)*unit_conv_core
443 core_particle_set(sh_index)%atom_index = at_index_c(ishell)
444 rab =
pbc(shell_particle_set(sh_index)%r, core_particle_set(sh_index)%r, cell)
445 cpassert(trim(at_name(ishell)) == trim(at_name_c(ishell)))
446 cpassert(at_index(ishell) == at_index_c(ishell))
448 rab =
pbc(shell_particle_set(sh_index)%r, particle_set(at_index(ishell))%r, cell)
451 dab = sqrt(rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3))
452 IF (shell%max_dist > 0.0_dp .AND. shell%max_dist < dab)
THEN
453 IF (output_unit > 0)
THEN
454 WRITE (output_unit, *)
"WARNING : shell and core for atom ", at_index(ishell),
" seem to be too distant. "
459 cpabort(
"shell coordinate assigned to the wrong atom. check the shell indexes in the input")
462 DEALLOCATE (r, at_index, at_name)
463 DEALLOCATE (rc, at_index_c, at_name_c)
469 CALL timestop(handle)
subroutine pbc(r, r_pbc, s, s_pbc, a, b, c, alpha, beta, gamma, debug, info, pbc0, h, hinv)
...
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Handles all functions related to the CELL.
subroutine, public scaled_to_real(r, s, cell)
Transform scaled cell coordinates real coordinates. r=h*s.
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Utility routines to read data from files. Kept as close as possible to the old parser because.
elemental subroutine, public read_float_object(string, object, error_message)
Returns a floating point number read from a string including fraction like z1/z2.
real(kind=dp) function, public cp_unit_to_cp2k(value, unit_str, defaults, power)
converts to the internal cp2k units to the given unit
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
An array-based list which grows on demand. When the internal array is full, a new array of twice the ...
Utility routines for the memory handling.
Define the data structure for the particle information.
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.
character(len=default_string_length) function, public id2str(id)
returns the string associated with a given id
Utilities for string manipulations.
elemental subroutine, public uppercase(string)
Convert all lower case characters in a string to upper case.
Control for reading in different topologies and coordinates.