(git:374b731)
Loading...
Searching...
No Matches
eip_environment_types.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 The environment for the empirical interatomic potential methods.
10!> \par History
11!> 03.2006 initial create [tdk]
12!> \author Thomas D. Kuehne (tkuehne@phys.chem.ethz.ch)
13! **************************************************************************************************
19 USE cell_types, ONLY: cell_release,&
29 USE kinds, ONLY: dp
42 USE virial_types, ONLY: virial_type
43#include "./base/base_uses.f90"
44
45 IMPLICIT NONE
46 PRIVATE
47
48 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'eip_environment_types'
49
50 ! *** Public data types ***
51 PUBLIC :: eip_environment_type
52
53 ! *** Public subroutines ***
54 PUBLIC :: eip_env_release, &
58
59! **************************************************************************************************
60!> \brief The empirical interatomic potential environment
61!> \param eip_model Specifies which EIP model is in use.
62!> \param eip_kinetic_energy The EIP kinetic energy
63!> \param eip_potential_energy The EIP potential energy
64!> \param eip_energy The total eip energy
65!> \param eip_energy_var Variance of the energy/atom
66!> \param eip_forces The final eip forces [eV/A]
67!> \param coord_avg The average coordination number
68!> \param coord_var The variance of the coordination number
69!> \param count Counts how often the empirical interatomic potential function
70!> is called. Don't ask why this is a real!
71!> \param subsystem the particles, molecules,... of this environment
72!> \param eip_input Pointer to the EIP input section
73!> \param force_env_input Pointer to the force_env input section
74!> \param cell The simulation cell
75!> \param cell_ref The reference simulation cell
76!> \param use_ref_cell Logical which indicates if reference
77!> simulation cell is used
78!> \param virial Dummy virial pointer
79!> \par History
80!> 03.2006 initial create [tdk]
81!> \author Thomas D. Kuehne (tkuehne@phys.chem.ethz.ch)
82! **************************************************************************************************
84 INTEGER :: eip_model = 0
85 REAL(kind=dp) :: eip_energy = 0.0_dp, &
86 eip_kinetic_energy = 0.0_dp, &
87 eip_potential_energy = 0.0_dp, &
88 eip_energy_var = 0.0_dp
89 REAL(kind=dp), DIMENSION(:, :), POINTER :: eip_forces => null()
90 REAL(kind=dp) :: coord_avg = 0.0_dp, &
91 coord_var = 0.0_dp, &
92 count = 0.0_dp
93 TYPE(cp_subsys_type), POINTER :: subsys => null()
94 TYPE(section_vals_type), POINTER :: eip_input => null(), &
95 force_env_input => null()
96 TYPE(cell_type), POINTER :: cell_ref => null()
97 LOGICAL :: use_ref_cell = .false.
99
100CONTAINS
101
102! **************************************************************************************************
103!> \brief Releases the given eip environment (see doc/ReferenceCounting.html)
104!> \param eip_env The eip environment to release
105!> \par History
106!> 03.2006 initial create [tdk]
107!> \author Thomas D. Kuehne (tkuehne@phys.chem.ethz.ch)
108! **************************************************************************************************
109 SUBROUTINE eip_env_release(eip_env)
110
111 TYPE(eip_environment_type), INTENT(INOUT) :: eip_env
112
113 IF (ASSOCIATED(eip_env%eip_forces)) THEN
114 DEALLOCATE (eip_env%eip_forces)
115 END IF
116 IF (ASSOCIATED(eip_env%subsys)) THEN
117 CALL cp_subsys_release(eip_env%subsys)
118 END IF
119 IF (ASSOCIATED(eip_env%subsys)) THEN
120 CALL cp_subsys_release(eip_env%subsys)
121 END IF
122 !IF (ASSOCIATED(eip_env%eip_input)) THEN
123 ! CALL section_vals_release(eip_env%eip_input)
124 !END IF
125 !IF (ASSOCIATED(eip_env%force_env_input)) THEN
126 ! CALL section_vals_release(eip_env%force_env_input)
127 !END IF
128 IF (ASSOCIATED(eip_env%cell_ref)) THEN
129 CALL cell_release(eip_env%cell_ref)
130 END IF
131 END SUBROUTINE eip_env_release
132
133! **************************************************************************************************
134!> \brief Returns various attributes of the eip environment
135!> \param eip_env The enquired eip environment
136!> \param eip_model Specifies which EIP model is in use.
137!> \param eip_energy The total eip energy
138!> \param eip_energy_var Variance of the energy/atom
139!> \param eip_forces The final eip forces [eV/A]
140!> \param coord_avg The average coordination number
141!> \param coord_var The variance of the coordination number
142!> \param count Counts how often the empirical interatomic potential function
143!> is called. Don't ask why this is a real!
144!> \param subsys the particles, molecules,... of this environment
145!> \param atomic_kind_set The set of all atomic kinds involved
146!> \param particle_set The set of all particles
147!> \param local_particles All particles on this particular node
148!> \param molecule_kind_set The set of all different molecule kinds involved
149!> \param molecule_set The set of all molecules
150!> \param local_molecules All molecules on this particular node
151!> \param eip_input the pointer to the EIP input section
152!> \param force_env_input Pointer to the force_env input section
153!> \param cell The simulation cell
154!> \param cell_ref The reference simulation cell
155!> \param use_ref_cell Logical which indicates if reference
156!> simulation cell is used
157!> \param eip_kinetic_energy The EIP kinetic energy
158!> \param eip_potential_energy The EIP potential energy
159!> \param virial Dummy virial pointer
160!>
161!> For possible missing arguments see the attributes of
162!> eip_environment_type
163!> \par History
164!> 03.2006 initial create [tdk]
165!> \author Thomas D. Kuehne (tkuehne@phys.chem.ethz.ch)
166! **************************************************************************************************
167 SUBROUTINE eip_env_get(eip_env, eip_model, eip_energy, eip_energy_var, &
168 eip_forces, coord_avg, coord_var, count, subsys, &
169 atomic_kind_set, particle_set, local_particles, &
170 molecule_kind_set, molecule_set, local_molecules, &
171 eip_input, force_env_input, cell, cell_ref, &
172 use_ref_cell, eip_kinetic_energy, eip_potential_energy, &
173 virial)
174
175 TYPE(eip_environment_type), INTENT(IN) :: eip_env
176 INTEGER, INTENT(OUT), OPTIONAL :: eip_model
177 REAL(kind=dp), INTENT(OUT), OPTIONAL :: eip_energy, eip_energy_var
178 REAL(kind=dp), DIMENSION(:, :), OPTIONAL, POINTER :: eip_forces
179 REAL(kind=dp), INTENT(OUT), OPTIONAL :: coord_avg, coord_var, count
180 TYPE(cp_subsys_type), OPTIONAL, POINTER :: subsys
181 TYPE(atomic_kind_type), DIMENSION(:), OPTIONAL, &
182 POINTER :: atomic_kind_set
183 TYPE(particle_type), DIMENSION(:), OPTIONAL, &
184 POINTER :: particle_set
185 TYPE(distribution_1d_type), OPTIONAL, POINTER :: local_particles
186 TYPE(molecule_kind_type), DIMENSION(:), OPTIONAL, &
187 POINTER :: molecule_kind_set
188 TYPE(molecule_type), DIMENSION(:), OPTIONAL, &
189 POINTER :: molecule_set
190 TYPE(distribution_1d_type), OPTIONAL, POINTER :: local_molecules
191 TYPE(section_vals_type), OPTIONAL, POINTER :: eip_input, force_env_input
192 TYPE(cell_type), OPTIONAL, POINTER :: cell, cell_ref
193 LOGICAL, INTENT(OUT), OPTIONAL :: use_ref_cell
194 REAL(kind=dp), INTENT(OUT), OPTIONAL :: eip_kinetic_energy, eip_potential_energy
195 TYPE(virial_type), OPTIONAL, POINTER :: virial
196
197 TYPE(atomic_kind_list_type), POINTER :: atomic_kinds
198 TYPE(molecule_kind_list_type), POINTER :: molecule_kinds
199 TYPE(molecule_list_type), POINTER :: molecules
200 TYPE(particle_list_type), POINTER :: particles
201
202! ------------------------------------------------------------------------
203
204 NULLIFY (atomic_kinds, particles, molecules, molecule_kinds)
205
206 IF (PRESENT(eip_model)) eip_model = eip_env%eip_model
207 IF (PRESENT(eip_kinetic_energy)) THEN
208 eip_kinetic_energy = eip_env%eip_kinetic_energy
209 END IF
210 IF (PRESENT(eip_potential_energy)) THEN
211 eip_potential_energy = eip_env%eip_potential_energy
212 END IF
213 IF (PRESENT(eip_energy)) eip_energy = eip_env%eip_energy
214 IF (PRESENT(eip_energy_var)) eip_energy_var = eip_env%eip_energy_var
215 IF (PRESENT(eip_forces)) eip_forces = eip_env%eip_forces
216 IF (PRESENT(coord_avg)) coord_avg = eip_env%coord_avg
217 IF (PRESENT(coord_var)) coord_var = eip_env%coord_var
218 IF (PRESENT(count)) count = eip_env%count
219 IF (PRESENT(subsys)) subsys => eip_env%subsys
220 CALL cp_subsys_get(eip_env%subsys, &
221 atomic_kinds=atomic_kinds, &
222 particles=particles, &
223 molecule_kinds=molecule_kinds, &
224 molecules=molecules, &
225 local_molecules=local_molecules, &
226 local_particles=local_particles, &
227 virial=virial, &
228 cell=cell)
229 IF (PRESENT(atomic_kind_set)) atomic_kind_set => atomic_kinds%els
230 IF (PRESENT(particle_set)) particle_set => particles%els
231 IF (PRESENT(molecule_kind_set)) molecule_kind_set => molecule_kinds%els
232 IF (PRESENT(molecule_set)) molecule_set => molecules%els
233
234 IF (PRESENT(eip_input)) eip_input => eip_env%eip_input
235 IF (PRESENT(force_env_input)) force_env_input => eip_env%force_env_input
236 IF (PRESENT(cell_ref)) cell_ref => eip_env%cell_ref
237 IF (PRESENT(use_ref_cell)) use_ref_cell = eip_env%use_ref_cell
238
239 END SUBROUTINE eip_env_get
240
241! **************************************************************************************************
242!> \brief Sets various attributes of the eip environment
243!> \param eip_env The enquired eip environment
244!> \param eip_model Specifies which EIP model is in use
245!> \param eip_energy The total eip energy
246!> \param eip_energy_var Variance of the energy/atom
247!> \param eip_forces The final eip forces [eV/A]
248!> \param coord_avg The average coordination number
249!> \param coord_var The variance of the coordination number
250!> \param count Counts how often the empirical interatomic potential function
251!> is called. Don't ask why this is a real!
252!> \param subsys the particles, molecules,... of this environment
253!> \param atomic_kind_set The set of all atomic kinds involved
254!> \param particle_set The set of all particles
255!> \param local_particles All particles on this particular node
256!> \param molecule_kind_set The set of all different molecule kinds involved
257!> \param molecule_set The set of all molecules
258!> \param local_molecules All molecules on this particular node
259!> \param eip_input the pointer to the EIP input section
260!> \param force_env_input Pointer to the force_env input section
261!> \param cell_ref The reference simulation cell
262!> \param use_ref_cell Logical which indicates if reference
263!> simulation cell is used
264!> \param eip_kinetic_energy The EIP kinetic energy
265!> \param eip_potential_energy The EIP potential energy
266!> \par History
267!> 03.2006 initial create [tdk]
268!> \author Thomas D. Kuehne (tkuehne@phys.chem.ethz.ch)
269!> \note
270!> For possible missing arguments see the attributes of eip_environment_type
271! **************************************************************************************************
272 SUBROUTINE eip_env_set(eip_env, eip_model, eip_energy, eip_energy_var, &
273 eip_forces, coord_avg, coord_var, count, subsys, &
274 atomic_kind_set, particle_set, local_particles, &
275 molecule_kind_set, molecule_set, local_molecules, &
276 eip_input, force_env_input, cell_ref, &
277 use_ref_cell, eip_kinetic_energy, eip_potential_energy)
278
279 TYPE(eip_environment_type), INTENT(INOUT) :: eip_env
280 INTEGER, INTENT(IN), OPTIONAL :: eip_model
281 REAL(kind=dp), INTENT(IN), OPTIONAL :: eip_energy, eip_energy_var
282 REAL(kind=dp), DIMENSION(:, :), OPTIONAL, POINTER :: eip_forces
283 REAL(kind=dp), INTENT(IN), OPTIONAL :: coord_avg, coord_var, count
284 TYPE(cp_subsys_type), OPTIONAL, POINTER :: subsys
285 TYPE(atomic_kind_type), DIMENSION(:), OPTIONAL, &
286 POINTER :: atomic_kind_set
287 TYPE(particle_type), DIMENSION(:), OPTIONAL, &
288 POINTER :: particle_set
289 TYPE(distribution_1d_type), OPTIONAL, POINTER :: local_particles
290 TYPE(molecule_kind_type), DIMENSION(:), OPTIONAL, &
291 POINTER :: molecule_kind_set
292 TYPE(molecule_type), DIMENSION(:), OPTIONAL, &
293 POINTER :: molecule_set
294 TYPE(distribution_1d_type), OPTIONAL, POINTER :: local_molecules
295 TYPE(section_vals_type), OPTIONAL, POINTER :: eip_input, force_env_input
296 TYPE(cell_type), OPTIONAL, POINTER :: cell_ref
297 LOGICAL, INTENT(IN), OPTIONAL :: use_ref_cell
298 REAL(kind=dp), INTENT(IN), OPTIONAL :: eip_kinetic_energy, eip_potential_energy
299
300 TYPE(atomic_kind_list_type), POINTER :: atomic_kinds
301 TYPE(molecule_kind_list_type), POINTER :: molecule_kinds
302 TYPE(molecule_list_type), POINTER :: molecules
303 TYPE(particle_list_type), POINTER :: particles
304
305 IF (PRESENT(eip_model)) eip_env%eip_model = eip_model
306 IF (PRESENT(eip_kinetic_energy)) THEN
307 eip_env%eip_kinetic_energy = eip_kinetic_energy
308 END IF
309 IF (PRESENT(eip_potential_energy)) THEN
310 eip_env%eip_potential_energy = eip_potential_energy
311 END IF
312 IF (PRESENT(eip_energy)) eip_env%eip_energy = eip_energy
313 IF (PRESENT(eip_energy_var)) eip_env%eip_energy_var = eip_energy_var
314 IF (PRESENT(eip_forces)) eip_env%eip_forces = eip_forces
315 IF (PRESENT(coord_avg)) eip_env%coord_avg = coord_avg
316 IF (PRESENT(coord_var)) eip_env%coord_var = coord_var
317 IF (PRESENT(count)) eip_env%count = count
318 IF (PRESENT(subsys)) THEN
319 IF (ASSOCIATED(eip_env%subsys)) THEN
320 IF (.NOT. ASSOCIATED(eip_env%subsys, subsys)) THEN
321 CALL cp_subsys_release(eip_env%subsys)
322 END IF
323 END IF
324 eip_env%subsys => subsys
325 END IF
326 IF (PRESENT(atomic_kind_set)) THEN
327 CALL atomic_kind_list_create(atomic_kinds, els_ptr=atomic_kind_set)
328 CALL cp_subsys_set(eip_env%subsys, atomic_kinds=atomic_kinds)
329 CALL atomic_kind_list_release(atomic_kinds)
330 END IF
331 IF (PRESENT(particle_set)) THEN
332 CALL particle_list_create(particles, els_ptr=particle_set)
333 CALL cp_subsys_set(eip_env%subsys, particles=particles)
334 CALL particle_list_release(particles)
335 END IF
336 IF (PRESENT(molecule_kind_set)) THEN
337 CALL molecule_kind_list_create(molecule_kinds, els_ptr=molecule_kind_set)
338 CALL cp_subsys_set(eip_env%subsys, molecule_kinds=molecule_kinds)
339 CALL molecule_kind_list_release(molecule_kinds)
340 END IF
341 IF (PRESENT(molecule_set)) THEN
342 CALL molecule_list_create(molecules, els_ptr=molecule_set)
343 CALL cp_subsys_set(eip_env%subsys, molecules=molecules)
344 CALL molecule_list_release(molecules)
345 END IF
346 IF (PRESENT(local_particles)) THEN
347 CALL cp_subsys_set(eip_env%subsys, local_particles=local_particles)
348 END IF
349 IF (PRESENT(local_molecules)) THEN
350 CALL cp_subsys_set(eip_env%subsys, local_molecules=local_molecules)
351 END IF
352
353 IF (PRESENT(eip_input)) eip_env%eip_input => eip_input
354 IF (PRESENT(force_env_input)) THEN
355 eip_env%force_env_input => force_env_input
356 END IF
357 IF (PRESENT(cell_ref)) THEN
358 CALL cell_retain(cell_ref)
359 CALL cell_release(eip_env%cell_ref)
360 eip_env%cell_ref => cell_ref
361 END IF
362 IF (PRESENT(use_ref_cell)) eip_env%use_ref_cell = use_ref_cell
363 END SUBROUTINE eip_env_set
364
365! **************************************************************************************************
366!> \brief Reinitializes the eip environment
367!> \param eip_env The eip environment to be reinitialized
368!> \par History
369!> 03.2006 initial create [tdk]
370!> \author Thomas D. Kuehne (tkuehne@phys.chem.ethz.ch)
371! **************************************************************************************************
372 SUBROUTINE eip_env_clear(eip_env)
373
374 TYPE(eip_environment_type), INTENT(INOUT) :: eip_env
375
376 eip_env%eip_model = 0
377 eip_env%eip_kinetic_energy = 0.0_dp
378 eip_env%eip_potential_energy = 0.0_dp
379 eip_env%eip_energy = 0.0_dp
380 eip_env%eip_energy_var = 0.0_dp
381 eip_env%coord_avg = 0.0_dp
382 eip_env%coord_var = 0.0_dp
383 eip_env%count = 0.0_dp
384 IF (ASSOCIATED(eip_env%eip_forces)) THEN
385 eip_env%eip_forces(:, :) = 0.0_dp
386 END IF
387 IF (ASSOCIATED(eip_env%subsys)) THEN
388 CALL cp_subsys_release(eip_env%subsys)
389 END IF
390 IF (ASSOCIATED(eip_env%eip_input)) THEN
391 CALL section_vals_release(eip_env%eip_input)
392 END IF
393 IF (ASSOCIATED(eip_env%force_env_input)) THEN
394 CALL section_vals_release(eip_env%force_env_input)
395 END IF
396 IF (ASSOCIATED(eip_env%cell_ref)) THEN
397 CALL cell_release(eip_env%cell_ref)
398 END IF
399 END SUBROUTINE eip_env_clear
400
401! **************************************************************************************************
402!> \brief Creates the eip environment
403!> \param eip_env The eip environment to be created
404!> \par History
405!> 03.2006 initial create [tdk]
406!> \author Thomas D. Kuehne (tkuehne@phys.chem.ethz.ch)
407! **************************************************************************************************
408 SUBROUTINE eip_env_create(eip_env)
409
410 TYPE(eip_environment_type), INTENT(OUT) :: eip_env
411
412 NULLIFY (eip_env%eip_forces)
413 NULLIFY (eip_env%subsys)
414 NULLIFY (eip_env%eip_input)
415 NULLIFY (eip_env%force_env_input)
416 NULLIFY (eip_env%cell_ref)
417
418 eip_env%use_ref_cell = .false.
419 CALL eip_env_clear(eip_env)
420 END SUBROUTINE eip_env_create
421
422END MODULE eip_environment_types
represent a simple array based list of the given type
subroutine, public atomic_kind_list_release(list)
releases a list (see doc/ReferenceCounting.html)
subroutine, public atomic_kind_list_create(list, els_ptr, owns_els, n_els)
creates a list
Define the atomic kind types and their sub types.
Handles all functions related to the CELL.
Definition cell_types.F:15
subroutine, public cell_release(cell)
releases the given cell (see doc/ReferenceCounting.html)
Definition cell_types.F:559
subroutine, public cell_retain(cell)
retains the given cell (see doc/ReferenceCounting.html)
Definition cell_types.F:542
types that represent a subsys, i.e. a part of the system
subroutine, public cp_subsys_release(subsys)
releases a subsys (see doc/ReferenceCounting.html)
subroutine, public cp_subsys_set(subsys, atomic_kinds, particles, local_particles, molecules, molecule_kinds, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, results, cell)
sets various propreties of the subsys
subroutine, public cp_subsys_get(subsys, ref_count, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell)
returns information about various attributes of the given subsys
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
The environment for the empirical interatomic potential methods.
subroutine, public eip_env_release(eip_env)
Releases the given eip environment (see doc/ReferenceCounting.html)
subroutine, public eip_env_get(eip_env, eip_model, eip_energy, eip_energy_var, eip_forces, coord_avg, coord_var, count, subsys, atomic_kind_set, particle_set, local_particles, molecule_kind_set, molecule_set, local_molecules, eip_input, force_env_input, cell, cell_ref, use_ref_cell, eip_kinetic_energy, eip_potential_energy, virial)
Returns various attributes of the eip environment.
subroutine, public eip_env_create(eip_env)
Creates the eip environment.
subroutine, public eip_env_set(eip_env, eip_model, eip_energy, eip_energy_var, eip_forces, coord_avg, coord_var, count, subsys, atomic_kind_set, particle_set, local_particles, molecule_kind_set, molecule_set, local_molecules, eip_input, force_env_input, cell_ref, use_ref_cell, eip_kinetic_energy, eip_potential_energy)
Sets various attributes of the eip environment.
objects that represent the structure of input sections and the data contained in an input section
recursive subroutine, public section_vals_release(section_vals)
releases the given object
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
represent a simple array based list of the given type
subroutine, public molecule_kind_list_create(list, els_ptr, owns_els, n_els)
creates a list
subroutine, public molecule_kind_list_release(list)
releases a list (see doc/ReferenceCounting.html)
Define the molecule kind structure types and the corresponding functionality.
represent a simple array based list of the given type
subroutine, public molecule_list_create(list, els_ptr, owns_els, n_els)
creates a list
subroutine, public molecule_list_release(list)
releases a list (see doc/ReferenceCounting.html)
Define the data structure for the molecule information.
represent a simple array based list of the given type
subroutine, public particle_list_create(list, els_ptr, owns_els, n_els)
creates a list
subroutine, public particle_list_release(list)
releases a list (see doc/ReferenceCounting.html)
Define the data structure for the particle information.
Provides all information about an atomic kind.
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
represents a system: atoms, molecules, their pos,vel,...
structure to store local (to a processor) ordered lists of integers.
The empirical interatomic potential environment.