(git:0de0cc2)
mixed_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 !> \author fschiff SEPT-11-06
10 ! **************************************************************************************************
14  atomic_kind_list_type
15  USE atomic_kind_types, ONLY: atomic_kind_type
16  USE cell_types, ONLY: cell_release,&
17  cell_retain,&
18  cell_type
19  USE cp_log_handling, ONLY: cp_logger_p_type,&
21  USE cp_result_types, ONLY: cp_result_type
22  USE cp_subsys_types, ONLY: cp_subsys_get,&
25  cp_subsys_type
26  USE distribution_1d_types, ONLY: distribution_1d_type
29  section_vals_type
30  USE kinds, ONLY: default_path_length,&
32  dp
33  USE message_passing, ONLY: mp_para_env_p_type,&
35  mp_para_env_type
36  USE mixed_cdft_types, ONLY: mixed_cdft_type,&
39  mixed_energy_type
42  molecule_kind_list_type
43  USE molecule_kind_types, ONLY: molecule_kind_type
46  molecule_list_type
47  USE molecule_types, ONLY: molecule_type
50  particle_list_type
51  USE particle_types, ONLY: particle_type
52  USE qs_rho_types, ONLY: qs_rho_p_type
53 #include "./base/base_uses.f90"
54 
55  IMPLICIT NONE
56  PRIVATE
57 
58 ! **************************************************************************************************
59 !> \param mixed_env the pointer to the mixed_env
60 !> \par History
61 !> 11/06 Created [fschiff]
62 !> 12/15-12/16 Mixed CDFT [Nico Holmberg]
63 ! **************************************************************************************************
64  TYPE mixed_environment_type
65  TYPE(cell_type), POINTER :: cell_ref => null()
66  TYPE(mixed_energy_type), POINTER :: mixed_energy => null()
67  TYPE(mp_para_env_type), POINTER :: para_env => null()
68  TYPE(cp_subsys_type), POINTER :: subsys => null()
69  TYPE(section_vals_type), POINTER :: input => null()
70  REAL(KIND=dp), DIMENSION(:), POINTER :: energies => null()
71  ! Parallelization of multiple force_eval
72  INTEGER :: ngroups = -1
73  INTEGER, DIMENSION(:), POINTER :: group_distribution => null()
74  TYPE(mp_para_env_p_type), DIMENSION(:), POINTER :: sub_para_env => null()
75  TYPE(cp_logger_p_type), DIMENSION(:), POINTER :: sub_logger => null()
76  REAL(KIND=dp), POINTER, DIMENSION(:) :: val => null()
77  CHARACTER(LEN=default_string_length), &
78  DIMENSION(:), POINTER :: par => null()
79  REAL(KIND=dp) :: dx = 0.0_dp, lerr = 0.0_dp
80  CHARACTER(default_path_length) :: coupling_function = ""
81  ! Mixed CDFT control parameters
82  LOGICAL :: do_mixed_cdft = .false., do_mixed_et = .false., &
83  do_mixed_qmmm_cdft = .false.
84  INTEGER :: et_freq = -1
85  REAL(KIND=dp), DIMENSION(:, :), POINTER :: strength => null()
86  TYPE(mixed_cdft_type), POINTER :: cdft_control => null()
87  ! Densities from sunbsystem
88  TYPE(qs_rho_p_type), DIMENSION(:), ALLOCATABLE :: subsys_dens
89  END TYPE mixed_environment_type
90 
91 ! *** Public data types ***
92 
93  PUBLIC :: mixed_environment_type
94 
95 ! *** Public subroutines ***
96 
97  PUBLIC :: get_mixed_env, &
98  set_mixed_env, &
101 
102  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'mixed_environment_types'
103 
104 CONTAINS
105 
106 ! **************************************************************************************************
107 !> \brief Get the MIXED environment.
108 !> \param mixed_env the pointer to the mixed_env
109 !> \param atomic_kind_set ...
110 !> \param particle_set ...
111 !> \param local_particles ...
112 !> \param local_molecules ...
113 !> \param molecule_kind_set ...
114 !> \param molecule_set ...
115 !> \param cell ...
116 !> \param cell_ref ...
117 !> \param mixed_energy ...
118 !> \param para_env ...
119 !> \param sub_para_env ...
120 !> \param subsys ...
121 !> \param input ...
122 !> \param results ...
123 !> \param cdft_control ...
124 ! **************************************************************************************************
125  SUBROUTINE get_mixed_env(mixed_env, atomic_kind_set, particle_set, &
126  local_particles, local_molecules, molecule_kind_set, &
127  molecule_set, cell, cell_ref, &
128  mixed_energy, para_env, sub_para_env, subsys, &
129  input, results, cdft_control)
130 
131  TYPE(mixed_environment_type), INTENT(IN) :: mixed_env
132  TYPE(atomic_kind_type), OPTIONAL, POINTER :: atomic_kind_set(:)
133  TYPE(particle_type), OPTIONAL, POINTER :: particle_set(:)
134  TYPE(distribution_1d_type), OPTIONAL, POINTER :: local_particles, local_molecules
135  TYPE(molecule_kind_type), OPTIONAL, POINTER :: molecule_kind_set(:)
136  TYPE(molecule_type), OPTIONAL, POINTER :: molecule_set(:)
137  TYPE(cell_type), OPTIONAL, POINTER :: cell, cell_ref
138  TYPE(mixed_energy_type), OPTIONAL, POINTER :: mixed_energy
139  TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env
140  TYPE(mp_para_env_p_type), DIMENSION(:), OPTIONAL, &
141  POINTER :: sub_para_env
142  TYPE(cp_subsys_type), OPTIONAL, POINTER :: subsys
143  TYPE(section_vals_type), OPTIONAL, POINTER :: input
144  TYPE(cp_result_type), OPTIONAL, POINTER :: results
145  TYPE(mixed_cdft_type), OPTIONAL, POINTER :: cdft_control
146 
147  TYPE(atomic_kind_list_type), POINTER :: atomic_kinds
148  TYPE(molecule_kind_list_type), POINTER :: molecule_kinds
149  TYPE(molecule_list_type), POINTER :: molecules
150  TYPE(particle_list_type), POINTER :: particles
151 
152  NULLIFY (atomic_kinds, particles, molecules, molecule_kinds)
153  cpassert(ASSOCIATED(mixed_env%subsys))
154 
155  IF (PRESENT(input)) input => mixed_env%input
156  IF (PRESENT(cell_ref)) cell_ref => mixed_env%cell_ref
157  IF (PRESENT(mixed_energy)) mixed_energy => mixed_env%mixed_energy
158  IF (PRESENT(para_env)) para_env => mixed_env%para_env
159  IF (PRESENT(sub_para_env)) sub_para_env => mixed_env%sub_para_env
160  IF (PRESENT(cdft_control)) cdft_control => mixed_env%cdft_control
161  IF (PRESENT(subsys)) subsys => mixed_env%subsys
162  CALL cp_subsys_get(mixed_env%subsys, &
163  atomic_kinds=atomic_kinds, &
164  local_molecules=local_molecules, &
165  local_particles=local_particles, &
166  particles=particles, &
167  molecule_kinds=molecule_kinds, &
168  molecules=molecules, &
169  results=results, &
170  cell=cell)
171  IF (PRESENT(atomic_kind_set)) atomic_kind_set => atomic_kinds%els
172  IF (PRESENT(particle_set)) particle_set => particles%els
173  IF (PRESENT(molecule_kind_set)) molecule_kind_set => molecule_kinds%els
174  IF (PRESENT(molecule_set)) molecule_set => molecules%els
175 
176  END SUBROUTINE get_mixed_env
177 
178 ! **************************************************************************************************
179 !> \brief Set the MIXED environment.
180 !> \param mixed_env the pointer to the mixed_env
181 !> \param atomic_kind_set ...
182 !> \param particle_set ...
183 !> \param local_particles ...
184 !> \param local_molecules ...
185 !> \param molecule_kind_set ...
186 !> \param molecule_set ...
187 !> \param cell_ref ...
188 !> \param mixed_energy ...
189 !> \param subsys ...
190 !> \param input ...
191 !> \param sub_para_env ...
192 !> \param cdft_control ...
193 ! **************************************************************************************************
194  SUBROUTINE set_mixed_env(mixed_env, atomic_kind_set, particle_set, &
195  local_particles, local_molecules, molecule_kind_set, &
196  molecule_set, cell_ref, mixed_energy, subsys, &
197  input, sub_para_env, cdft_control)
198 
199  TYPE(mixed_environment_type), INTENT(INOUT) :: mixed_env
200  TYPE(atomic_kind_type), OPTIONAL, POINTER :: atomic_kind_set(:)
201  TYPE(particle_type), OPTIONAL, POINTER :: particle_set(:)
202  TYPE(distribution_1d_type), OPTIONAL, POINTER :: local_particles, local_molecules
203  TYPE(molecule_kind_type), OPTIONAL, POINTER :: molecule_kind_set(:)
204  TYPE(molecule_type), OPTIONAL, POINTER :: molecule_set(:)
205  TYPE(cell_type), OPTIONAL, POINTER :: cell_ref
206  TYPE(mixed_energy_type), OPTIONAL, POINTER :: mixed_energy
207  TYPE(cp_subsys_type), OPTIONAL, POINTER :: subsys
208  TYPE(section_vals_type), OPTIONAL, POINTER :: input
209  TYPE(mp_para_env_p_type), DIMENSION(:), OPTIONAL, &
210  POINTER :: sub_para_env
211  TYPE(mixed_cdft_type), OPTIONAL, POINTER :: cdft_control
212 
213  TYPE(atomic_kind_list_type), POINTER :: atomic_kinds
214  TYPE(molecule_kind_list_type), POINTER :: molecule_kinds
215  TYPE(molecule_list_type), POINTER :: molecules
216  TYPE(particle_list_type), POINTER :: particles
217 
218  IF (PRESENT(cell_ref)) THEN
219  CALL cell_retain(cell_ref)
220  CALL cell_release(mixed_env%cell_ref)
221  mixed_env%cell_ref => cell_ref
222  END IF
223  IF (PRESENT(input)) THEN
224  CALL section_vals_retain(input)
225  CALL section_vals_release(mixed_env%input)
226  mixed_env%input => input
227  END IF
228  IF (PRESENT(mixed_energy)) mixed_env%mixed_energy => mixed_energy
229  IF (PRESENT(subsys)) THEN
230  IF (ASSOCIATED(mixed_env%subsys)) THEN
231  IF (.NOT. ASSOCIATED(mixed_env%subsys, subsys)) THEN
232  CALL cp_subsys_release(mixed_env%subsys)
233  END IF
234  END IF
235  mixed_env%subsys => subsys
236  END IF
237  IF (PRESENT(sub_para_env)) THEN
238  mixed_env%sub_para_env => sub_para_env
239  END IF
240  IF (PRESENT(cdft_control)) mixed_env%cdft_control => cdft_control
241  IF (PRESENT(atomic_kind_set)) THEN
242  CALL atomic_kind_list_create(atomic_kinds, &
243  els_ptr=atomic_kind_set)
244  CALL cp_subsys_set(mixed_env%subsys, &
245  atomic_kinds=atomic_kinds)
246  CALL atomic_kind_list_release(atomic_kinds)
247  END IF
248  IF (PRESENT(particle_set)) THEN
249  CALL particle_list_create(particles, &
250  els_ptr=particle_set)
251  CALL cp_subsys_set(mixed_env%subsys, &
252  particles=particles)
253  CALL particle_list_release(particles)
254  END IF
255  IF (PRESENT(local_particles)) THEN
256  CALL cp_subsys_set(mixed_env%subsys, &
257  local_particles=local_particles)
258  END IF
259  IF (PRESENT(local_molecules)) THEN
260  CALL cp_subsys_set(mixed_env%subsys, &
261  local_molecules=local_molecules)
262  END IF
263  IF (PRESENT(molecule_kind_set)) THEN
264  CALL molecule_kind_list_create(molecule_kinds, els_ptr=molecule_kind_set)
265  CALL cp_subsys_set(mixed_env%subsys, molecule_kinds=molecule_kinds)
266  CALL molecule_kind_list_release(molecule_kinds)
267  END IF
268  IF (PRESENT(molecule_set)) THEN
269  CALL molecule_list_create(molecules, els_ptr=molecule_set)
270  CALL cp_subsys_set(mixed_env%subsys, molecules=molecules)
271  CALL molecule_list_release(molecules)
272  END IF
273 
274  END SUBROUTINE set_mixed_env
275 
276 ! **************************************************************************************************
277 !> \brief allocates and intitializes a mixed_env
278 !> \param mixed_env the object to create
279 !> \param para_env the parallel environment for the qs_env
280 !> \author fschiff 11.06
281 ! **************************************************************************************************
282  SUBROUTINE mixed_env_create(mixed_env, para_env)
283  TYPE(mixed_environment_type), INTENT(OUT) :: mixed_env
284  TYPE(mp_para_env_type), INTENT(IN), TARGET :: para_env
285 
286  mixed_env%para_env => para_env
287  CALL mixed_env%para_env%retain()
288  END SUBROUTINE mixed_env_create
289 
290 ! **************************************************************************************************
291 !> \brief releases the given mixed_env (see doc/ReferenceCounting.html)
292 !> \param mixed_env the object to release
293 !> \author fschiff 11.06
294 ! **************************************************************************************************
295  SUBROUTINE mixed_env_release(mixed_env)
296  TYPE(mixed_environment_type), INTENT(INOUT) :: mixed_env
297 
298  INTEGER :: i, ngroups
299 
300  ngroups = SIZE(mixed_env%sub_para_env)
301  DO i = 1, ngroups
302  IF (ASSOCIATED(mixed_env%sub_para_env(i)%para_env)) THEN
303  CALL cp_logger_release(mixed_env%sub_logger(i)%p)
304  CALL mp_para_env_release(mixed_env%sub_para_env(i)%para_env)
305  END IF
306  END DO
307  DEALLOCATE (mixed_env%sub_para_env)
308  DEALLOCATE (mixed_env%sub_logger)
309  DEALLOCATE (mixed_env%energies)
310  IF (ASSOCIATED(mixed_env%par)) THEN
311  DEALLOCATE (mixed_env%par)
312  END IF
313  IF (ASSOCIATED(mixed_env%val)) THEN
314  DEALLOCATE (mixed_env%val)
315  END IF
316  CALL cell_release(mixed_env%cell_ref)
317  CALL mp_para_env_release(mixed_env%para_env)
318  CALL deallocate_mixed_energy(mixed_env%mixed_energy)
319  CALL cp_subsys_release(mixed_env%subsys)
320  CALL section_vals_release(mixed_env%input)
321  IF (ASSOCIATED(mixed_env%group_distribution)) THEN
322  DEALLOCATE (mixed_env%group_distribution)
323  END IF
324  IF (ASSOCIATED(mixed_env%cdft_control)) &
325  CALL mixed_cdft_type_release(mixed_env%cdft_control)
326  IF (ASSOCIATED(mixed_env%strength)) &
327  DEALLOCATE (mixed_env%strength)
328 
329  END SUBROUTINE mixed_env_release
330 
331 END MODULE mixed_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
various routines to log and control the output. The idea is that decisions about where to log should ...
subroutine, public cp_logger_release(logger)
releases this logger
set of type/routines to handle the storage of results in force_envs
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...
objects that represent the structure of input sections and the data contained in an input section
subroutine, public section_vals_retain(section_vals)
retains the given section values (see doc/ReferenceCounting.html)
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
integer, parameter, public default_string_length
Definition: kinds.F:57
integer, parameter, public default_path_length
Definition: kinds.F:58
Interface to the message passing library MPI.
subroutine, public mp_para_env_release(para_env)
releases the para object (to be called when you don't want anymore the shared copy of this object)
Types for mixed CDFT calculations.
subroutine, public mixed_cdft_type_release(cdft_control)
releases the given mixed_cdft_type
subroutine, public deallocate_mixed_energy(mixed_energy)
Deallocate a mixed energy data structure.
subroutine, public mixed_env_release(mixed_env)
releases the given mixed_env (see doc/ReferenceCounting.html)
subroutine, public get_mixed_env(mixed_env, atomic_kind_set, particle_set, local_particles, local_molecules, molecule_kind_set, molecule_set, cell, cell_ref, mixed_energy, para_env, sub_para_env, subsys, input, results, cdft_control)
Get the MIXED environment.
subroutine, public set_mixed_env(mixed_env, atomic_kind_set, particle_set, local_particles, local_molecules, molecule_kind_set, molecule_set, cell_ref, mixed_energy, subsys, input, sub_para_env, cdft_control)
Set the MIXED environment.
subroutine, public mixed_env_create(mixed_env, para_env)
allocates and intitializes a mixed_env
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.
superstucture that hold various representations of the density and keeps track of which ones are vali...
Definition: qs_rho_types.F:18