(git:ccc2433)
dimer_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 Contains types used for a Dimer Method calculations
10 !> \par History
11 !> Luca Bellucci 11.2017 added kdimer and beta
12 !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
13 ! **************************************************************************************************
15 
16  USE cell_types, ONLY: use_perd_x,&
17  use_perd_xy,&
18  use_perd_xyz,&
19  use_perd_xz,&
20  use_perd_y,&
21  use_perd_yz,&
24  USE cp_subsys_types, ONLY: cp_subsys_get,&
25  cp_subsys_type
26  USE global_types, ONLY: global_environment_type
30  section_vals_type,&
32  USE kinds, ONLY: dp
33  USE molecule_kind_list_types, ONLY: molecule_kind_list_type
34  USE molecule_kind_types, ONLY: fixd_constraint_type,&
36  molecule_kind_type
37 #include "../base/base_uses.f90"
38 
39  IMPLICIT NONE
40  PRIVATE
41 
42  LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .true.
43  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'dimer_types'
44 
45  PUBLIC :: dimer_env_type, &
50 
51 ! **************************************************************************************************
52 !> \brief Type containing all informations abour the rotation of the Dimer
53 !> \par History
54 !> none
55 !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
56 ! **************************************************************************************************
57  TYPE dimer_rotational_type
58  ! Rotational parameters
59  INTEGER :: rotation_step = 0
60  LOGICAL :: interpolate_gradient = .false.
61  REAL(KIND=dp) :: angle_tol = 0.0_dp, angle1 = 0.0_dp, angle2 = 0.0_dp, &
62  dcdp = 0.0_dp, curvature = 0.0_dp
63  REAL(KIND=dp), POINTER, DIMENSION(:) :: g0 => null(), g1 => null(), g1p => null()
64  END TYPE dimer_rotational_type
65 
66 ! **************************************************************************************************
67 !> \brief Type containing all informations abour the translation of the Dimer
68 !> \par History
69 !> none
70 !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
71 ! **************************************************************************************************
72  TYPE dimer_translational_type
73  ! Translational parameters
74  REAL(KIND=dp), POINTER, DIMENSION(:) :: tls_vec => null()
75  END TYPE dimer_translational_type
76 
77 ! **************************************************************************************************
78 !> \brief Conjugate Directions type
79 !> \par History
80 !> none
81 !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
82 ! **************************************************************************************************
83  TYPE dimer_cg_rot_type
84  REAL(KIND=dp) :: norm_theta = 0.0_dp, norm_theta_old = 0.0_dp, norm_h = 0.0_dp
85  REAL(KIND=dp), POINTER, DIMENSION(:) :: nvec_old => null()
86  END TYPE dimer_cg_rot_type
87 
88 ! **************************************************************************************************
89 !> \brief Defines the environment for a Dimer Method calculation
90 !> \par History
91 !> none
92 !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
93 ! **************************************************************************************************
94  TYPE dimer_env_type
95  INTEGER :: ref_count = 0
96  REAL(KIND=dp) :: dr = 0.0_dp
97  REAL(KIND=dp), POINTER, DIMENSION(:) :: nvec => null()
98  REAL(KIND=dp) :: beta = 0.0_dp
99  TYPE(dimer_rotational_type) :: rot = dimer_rotational_type()
100  TYPE(dimer_translational_type) :: tsl = dimer_translational_type()
101  TYPE(dimer_cg_rot_type) :: cg_rot = dimer_cg_rot_type()
102  LOGICAL :: kdimer = .false.
103  END TYPE dimer_env_type
104 
105 CONTAINS
106 
107 ! **************************************************************************************************
108 !> \brief ...
109 !> \param dimer_env ...
110 !> \param subsys ...
111 !> \param globenv ...
112 !> \param dimer_section ...
113 !> \par History
114 !> Luca Bellucci 11.2017 added K-DIMER and BETA
115 !> 2016/03/03 [LTong] changed input natom to subsys
116 !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
117 ! **************************************************************************************************
118  SUBROUTINE dimer_env_create(dimer_env, subsys, globenv, dimer_section)
119  TYPE(dimer_env_type), POINTER :: dimer_env
120  TYPE(cp_subsys_type), POINTER :: subsys
121  TYPE(global_environment_type), POINTER :: globenv
122  TYPE(section_vals_type), POINTER :: dimer_section
123 
124  INTEGER :: i, isize, j, k, n_rep_val, natom, unit_nr
125  LOGICAL :: explicit
126  REAL(kind=dp) :: norm, xval(3)
127  REAL(kind=dp), DIMENSION(:), POINTER :: array
128  TYPE(section_vals_type), POINTER :: nvec_section
129 
131  cpassert(.NOT. ASSOCIATED(dimer_env))
132  ALLOCATE (dimer_env)
133  dimer_env%ref_count = 1
134  ! Setup NVEC
135  ! get natom
136  CALL cp_subsys_get(subsys=subsys, natom=natom)
137  ! Allocate the working arrays
138  ALLOCATE (dimer_env%nvec(natom*3))
139  ALLOCATE (dimer_env%rot%g0(natom*3))
140  ALLOCATE (dimer_env%rot%g1(natom*3))
141  ALLOCATE (dimer_env%rot%g1p(natom*3))
142  ! Check if the dimer vector is available in the input or not..
143  nvec_section => section_vals_get_subs_vals(dimer_section, "DIMER_VECTOR")
144  CALL section_vals_get(nvec_section, explicit=explicit)
145  IF (explicit) THEN
146  IF (unit_nr > 0) WRITE (unit_nr, *) "Reading Dimer Vector from file!"
147  NULLIFY (array)
148  CALL section_vals_val_get(nvec_section, "_DEFAULT_KEYWORD_", n_rep_val=n_rep_val)
149  isize = 0
150  DO i = 1, n_rep_val
151  CALL section_vals_val_get(nvec_section, "_DEFAULT_KEYWORD_", r_vals=array, i_rep_val=i)
152  DO j = 1, SIZE(array)
153  isize = isize + 1
154  dimer_env%nvec(isize) = array(j)
155  END DO
156  END DO
157  cpassert(isize == SIZE(dimer_env%nvec))
158  ELSE
159  CALL globenv%gaussian_rng_stream%fill(dimer_env%nvec)
160  END IF
161  ! Check for translation in the dimer vector and remove them
162  IF (natom > 1) THEN
163  xval = 0.0_dp
164  DO j = 1, natom
165  DO k = 1, 3
166  i = (j - 1)*3 + k
167  xval(k) = xval(k) + dimer_env%nvec(i)
168  END DO
169  END DO
170  ! Subtract net translations
171  xval = xval/real(natom*3, kind=dp)
172  DO j = 1, natom
173  DO k = 1, 3
174  i = (j - 1)*3 + k
175  dimer_env%nvec(i) = dimer_env%nvec(i) - xval(k)
176  END DO
177  END DO
178  END IF
179  ! set nvec components to zero for the corresponding constraints
180  CALL dimer_fixed_atom_control(dimer_env%nvec, subsys)
181  norm = sqrt(sum(dimer_env%nvec**2))
182  IF (norm <= epsilon(0.0_dp)) &
183  cpabort("The norm of the dimer vector is 0! Calculation cannot proceed further.")
184  dimer_env%nvec = dimer_env%nvec/norm
185  dimer_env%rot%rotation_step = do_first_rotation_step
186  CALL section_vals_val_get(dimer_section, "DR", r_val=dimer_env%dr)
187  CALL section_vals_val_get(dimer_section, "INTERPOLATE_GRADIENT", &
188  l_val=dimer_env%rot%interpolate_gradient)
189  CALL section_vals_val_get(dimer_section, "ANGLE_TOLERANCE", &
190  r_val=dimer_env%rot%angle_tol)
191  CALL section_vals_val_get(dimer_section, "K-DIMER", &
192  l_val=dimer_env%kdimer)
193  CALL section_vals_val_get(dimer_section, "BETA", &
194  r_val=dimer_env%beta)
195  ! initialise values
196  dimer_env%cg_rot%norm_h = 1.0_dp
197  dimer_env%rot%g0 = 0.0_dp
198  dimer_env%rot%g1 = 0.0_dp
199  dimer_env%rot%g1p = 0.0_dp
200  ALLOCATE (dimer_env%cg_rot%nvec_old(natom*3))
201  END SUBROUTINE dimer_env_create
202 
203 ! **************************************************************************************************
204 !> \brief ...
205 !> \param dimer_env ...
206 !> \par History
207 !> none
208 !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
209 ! **************************************************************************************************
210  SUBROUTINE dimer_env_retain(dimer_env)
211  TYPE(dimer_env_type), POINTER :: dimer_env
212 
213  cpassert(ASSOCIATED(dimer_env))
214  cpassert(dimer_env%ref_count > 0)
215  dimer_env%ref_count = dimer_env%ref_count + 1
216  END SUBROUTINE dimer_env_retain
217 
218 ! **************************************************************************************************
219 !> \brief ...
220 !> \param dimer_env ...
221 !> \par History
222 !> none
223 !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
224 ! **************************************************************************************************
225  SUBROUTINE dimer_env_release(dimer_env)
226  TYPE(dimer_env_type), POINTER :: dimer_env
227 
228  IF (ASSOCIATED(dimer_env)) THEN
229  cpassert(dimer_env%ref_count > 0)
230  dimer_env%ref_count = dimer_env%ref_count - 1
231  IF (dimer_env%ref_count == 0) THEN
232  IF (ASSOCIATED(dimer_env%nvec)) THEN
233  DEALLOCATE (dimer_env%nvec)
234  END IF
235  IF (ASSOCIATED(dimer_env%rot%g0)) THEN
236  DEALLOCATE (dimer_env%rot%g0)
237  END IF
238  IF (ASSOCIATED(dimer_env%rot%g1)) THEN
239  DEALLOCATE (dimer_env%rot%g1)
240  END IF
241  IF (ASSOCIATED(dimer_env%rot%g1p)) THEN
242  DEALLOCATE (dimer_env%rot%g1p)
243  END IF
244  IF (ASSOCIATED(dimer_env%cg_rot%nvec_old)) THEN
245  DEALLOCATE (dimer_env%cg_rot%nvec_old)
246  END IF
247  ! No need to deallocate tls_vec (just a pointer to aother local array)
248  NULLIFY (dimer_env%tsl%tls_vec)
249  DEALLOCATE (dimer_env)
250  END IF
251  END IF
252  END SUBROUTINE dimer_env_release
253 
254 ! **************************************************************************************************
255 !> \brief Set parts of a given array vec to zero according to fixed atom constraints.
256 !> When atoms are (partially) fixed then the relevant components of
257 !> nvec should be set to zero. Furthermore, the relevant components
258 !> of the gradient in CG should also be set to zero.
259 !> \param vec : vector to be modified
260 !> \param subsys : subsys type object used by CP2k
261 !> \par History
262 !> 2016/03/03 [LTong] created
263 !> \author Lianheng Tong [LTong]
264 ! **************************************************************************************************
265  SUBROUTINE dimer_fixed_atom_control(vec, subsys)
266  REAL(kind=dp), DIMENSION(:), INTENT(INOUT) :: vec
267  TYPE(cp_subsys_type), POINTER :: subsys
268 
269  INTEGER :: ii, ikind, ind, iparticle, nfixed_atoms, &
270  nkinds
271  TYPE(fixd_constraint_type), DIMENSION(:), POINTER :: fixd_list
272  TYPE(molecule_kind_list_type), POINTER :: molecule_kinds
273  TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
274  TYPE(molecule_kind_type), POINTER :: molecule_kind
275 
276  NULLIFY (molecule_kinds, molecule_kind, fixd_list)
277 
278  ! need to get constraint information from molecule information
279  CALL cp_subsys_get(subsys=subsys, &
280  molecule_kinds=molecule_kinds)
281  molecule_kind_set => molecule_kinds%els
282 
283  ! get total number of fixed atoms
284  ! nkinds is the kinds of molecules, not atoms
285  nkinds = molecule_kinds%n_els
286  DO ikind = 1, nkinds
287  molecule_kind => molecule_kind_set(ikind)
288  CALL get_molecule_kind(molecule_kind, &
289  nfixd=nfixed_atoms, &
290  fixd_list=fixd_list)
291  IF (ASSOCIATED(fixd_list)) THEN
292  DO ii = 1, nfixed_atoms
293  IF (.NOT. fixd_list(ii)%restraint%active) THEN
294  iparticle = fixd_list(ii)%fixd
295  ind = (iparticle - 1)*3
296  ! apply constraint to nvec
297  SELECT CASE (fixd_list(ii)%itype)
298  CASE (use_perd_x)
299  vec(ind + 1) = 0.0_dp
300  CASE (use_perd_y)
301  vec(ind + 2) = 0.0_dp
302  CASE (use_perd_z)
303  vec(ind + 3) = 0.0_dp
304  CASE (use_perd_xy)
305  vec(ind + 1) = 0.0_dp
306  vec(ind + 2) = 0.0_dp
307  CASE (use_perd_xz)
308  vec(ind + 1) = 0.0_dp
309  vec(ind + 3) = 0.0_dp
310  CASE (use_perd_yz)
311  vec(ind + 2) = 0.0_dp
312  vec(ind + 3) = 0.0_dp
313  CASE (use_perd_xyz)
314  vec(ind + 1) = 0.0_dp
315  vec(ind + 2) = 0.0_dp
316  vec(ind + 3) = 0.0_dp
317  END SELECT
318  END IF ! .NOT.fixd_list(ii)%restraint%active
319  END DO ! ii
320  END IF ! ASSOCIATED(fixd_list)
321  END DO ! ikind
322  END SUBROUTINE dimer_fixed_atom_control
323 
324 END MODULE dimer_types
Handles all functions related to the CELL.
Definition: cell_types.F:15
integer, parameter, public use_perd_xyz
Definition: cell_types.F:42
integer, parameter, public use_perd_y
Definition: cell_types.F:42
integer, parameter, public use_perd_xz
Definition: cell_types.F:42
integer, parameter, public use_perd_x
Definition: cell_types.F:42
integer, parameter, public use_perd_z
Definition: cell_types.F:42
integer, parameter, public use_perd_yz
Definition: cell_types.F:42
integer, parameter, public use_perd_xy
Definition: cell_types.F:42
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...
types that represent a subsys, i.e. a part of the system
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
Contains types used for a Dimer Method calculations.
Definition: dimer_types.F:14
subroutine, public dimer_env_retain(dimer_env)
...
Definition: dimer_types.F:211
subroutine, public dimer_env_release(dimer_env)
...
Definition: dimer_types.F:226
subroutine, public dimer_fixed_atom_control(vec, subsys)
Set parts of a given array vec to zero according to fixed atom constraints. When atoms are (partially...
Definition: dimer_types.F:266
subroutine, public dimer_env_create(dimer_env, subsys, globenv, dimer_section)
...
Definition: dimer_types.F:119
Define type storing the global information of a run. Keep the amount of stored data small....
Definition: global_types.F:21
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public do_first_rotation_step
objects that represent the structure of input sections and the data contained in an input section
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 dp
Definition: kinds.F:34
represent a simple array based list of the given type
Define the molecule kind structure types and the corresponding functionality.
subroutine, public get_molecule_kind(molecule_kind, atom_list, bond_list, bend_list, ub_list, impr_list, opbend_list, colv_list, fixd_list, g3x3_list, g4x6_list, vsite_list, torsion_list, shell_list, name, mass, charge, kind_number, natom, nbend, nbond, nub, nimpr, nopbend, nconstraint, nconstraint_fixd, nfixd, ncolv, ng3x3, ng4x6, nvsite, nfixd_restraint, ng3x3_restraint, ng4x6_restraint, nvsite_restraint, nrestraints, nmolecule, nsgf, nshell, ntorsion, molecule_list, nelectron, nelectron_alpha, nelectron_beta, bond_kind_set, bend_kind_set, ub_kind_set, impr_kind_set, opbend_kind_set, torsion_kind_set, molname_generated)
Get informations about a molecule kind.