Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 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 : ! **************************************************************************************************
14 : MODULE dimer_types
15 :
16 : USE cell_types, ONLY: &
17 : cell_transform_input_cartesian, cell_type, use_perd_x, use_perd_xy, use_perd_xyz, &
18 : use_perd_xz, use_perd_y, use_perd_yz, use_perd_z
19 : USE cp_log_handling, ONLY: cp_logger_get_default_io_unit,&
20 : cp_to_string
21 : USE cp_parser_methods, ONLY: parser_get_next_line,&
22 : parser_search_string
23 : USE cp_parser_types, ONLY: cp_parser_type,&
24 : parser_create,&
25 : parser_release
26 : USE cp_subsys_types, ONLY: cp_subsys_get,&
27 : cp_subsys_type
28 : USE force_env_types, ONLY: force_env_type
29 : USE global_types, ONLY: global_environment_type
30 : USE input_constants, ONLY: dimer_init_molden,&
31 : dimer_init_random,&
32 : do_first_rotation_step
33 : USE input_section_types, ONLY: section_vals_get,&
34 : section_vals_get_subs_vals,&
35 : section_vals_type,&
36 : section_vals_val_get
37 : USE kinds, ONLY: default_path_length,&
38 : dp
39 : USE message_passing, ONLY: mp_para_env_type
40 : USE molecule_kind_list_types, ONLY: molecule_kind_list_type
41 : USE molecule_kind_types, ONLY: fixd_constraint_type,&
42 : get_molecule_kind,&
43 : molecule_kind_type
44 : #include "../base/base_uses.f90"
45 :
46 : IMPLICIT NONE
47 : PRIVATE
48 :
49 : LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .TRUE.
50 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'dimer_types'
51 :
52 : PUBLIC :: dimer_env_type, &
53 : dimer_env_create, &
54 : dimer_env_retain, &
55 : dimer_env_release, &
56 : dimer_fixed_atom_control, &
57 : dimer_init_vector, &
58 : vib_get_index_weight
59 :
60 : ! **************************************************************************************************
61 : !> \brief Type containing all informations abour the rotation of the Dimer
62 : !> \par History
63 : !> none
64 : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
65 : ! **************************************************************************************************
66 : TYPE dimer_rotational_type
67 : ! Rotational parameters
68 : INTEGER :: rotation_step = 0
69 : LOGICAL :: interpolate_gradient = .FALSE.
70 : REAL(KIND=dp) :: angle_tol = 0.0_dp, angle1 = 0.0_dp, angle2 = 0.0_dp, &
71 : dCdp = 0.0_dp, curvature = 0.0_dp
72 : REAL(KIND=dp), POINTER, DIMENSION(:) :: g0 => NULL(), g1 => NULL(), g1p => NULL()
73 : END TYPE dimer_rotational_type
74 :
75 : ! **************************************************************************************************
76 : !> \brief Type containing all informations abour the translation of the Dimer
77 : !> \par History
78 : !> none
79 : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
80 : ! **************************************************************************************************
81 : TYPE dimer_translational_type
82 : ! Translational parameters
83 : REAL(KIND=dp), POINTER, DIMENSION(:) :: tls_vec => NULL()
84 : END TYPE dimer_translational_type
85 :
86 : ! **************************************************************************************************
87 : !> \brief Conjugate Directions type
88 : !> \par History
89 : !> none
90 : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
91 : ! **************************************************************************************************
92 : TYPE dimer_cg_rot_type
93 : REAL(KIND=dp) :: norm_theta = 0.0_dp, norm_theta_old = 0.0_dp, norm_h = 0.0_dp
94 : REAL(KIND=dp), POINTER, DIMENSION(:) :: nvec_old => NULL()
95 : END TYPE dimer_cg_rot_type
96 :
97 : ! **************************************************************************************************
98 : !> \brief Defines the environment for a Dimer Method calculation
99 : !> \par History
100 : !> none
101 : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
102 : ! **************************************************************************************************
103 : TYPE dimer_env_type
104 : INTEGER :: ref_count = 0
105 : REAL(KIND=dp) :: dr = 0.0_dp
106 : REAL(KIND=dp), POINTER, DIMENSION(:) :: nvec => NULL()
107 : REAL(KIND=dp) :: beta = 0.0_dp
108 : TYPE(dimer_rotational_type) :: rot = dimer_rotational_type()
109 : TYPE(dimer_translational_type) :: tsl = dimer_translational_type()
110 : TYPE(dimer_cg_rot_type) :: cg_rot = dimer_cg_rot_type()
111 : LOGICAL :: kdimer = .FALSE.
112 : END TYPE dimer_env_type
113 :
114 : CONTAINS
115 :
116 : ! **************************************************************************************************
117 : !> \brief ...
118 : !> \param dimer_env ...
119 : !> \param subsys ...
120 : !> \param globenv ...
121 : !> \param dimer_section ...
122 : !> \param force_env ...
123 : !> \par History
124 : !> Luca Bellucci 11.2017 added K-DIMER and BETA
125 : !> 2016/03/03 [LTong] changed input natom to subsys
126 : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
127 : ! **************************************************************************************************
128 24 : SUBROUTINE dimer_env_create(dimer_env, subsys, globenv, dimer_section, force_env)
129 : TYPE(dimer_env_type), POINTER :: dimer_env
130 : TYPE(cp_subsys_type), POINTER :: subsys
131 : TYPE(global_environment_type), POINTER :: globenv
132 : TYPE(section_vals_type), POINTER :: dimer_section
133 : TYPE(force_env_type), POINTER :: force_env
134 :
135 : INTEGER :: i, j, k, natom, unit_nr
136 : REAL(KIND=dp) :: norm, xval(3)
137 : TYPE(cell_type), POINTER :: cell
138 :
139 24 : NULLIFY (cell)
140 24 : unit_nr = cp_logger_get_default_io_unit()
141 24 : CPASSERT(.NOT. ASSOCIATED(dimer_env))
142 24 : ALLOCATE (dimer_env)
143 24 : dimer_env%ref_count = 1
144 : ! Setup NVEC
145 : ! get natom
146 24 : CALL cp_subsys_get(subsys=subsys, cell=cell, natom=natom)
147 : ! Allocate the working arrays
148 72 : ALLOCATE (dimer_env%rot%g0(natom*3))
149 48 : ALLOCATE (dimer_env%rot%g1(natom*3))
150 48 : ALLOCATE (dimer_env%rot%g1p(natom*3))
151 : ! Read dimer vector from input
152 : CALL dimer_init_vector(dimer_env, dimer_section, natom, &
153 24 : unit_nr, globenv, force_env)
154 24 : IF (unit_nr > 0) THEN
155 : WRITE (unit_nr, "(/,T2,A,T9,A,T71,I10)") &
156 12 : "DIMER|", "Dimension of dimer vector (natom * 3)", natom*3
157 53 : DO j = 1, natom
158 : WRITE (unit_nr, "(T2,A,T9,3(F12.6,3X))") &
159 176 : "DIMER|", dimer_env%nvec(3*j - 2:3*j)
160 : END DO
161 : END IF
162 : ! Transform input cell
163 106 : DO i = 1, natom
164 328 : xval(:) = dimer_env%nvec(3*i - 2:3*i)
165 82 : CALL cell_transform_input_cartesian(cell, xval)
166 352 : dimer_env%nvec(3*i - 2:3*i) = xval(:)
167 : END DO
168 : ! Check for translation in the dimer vector and remove them
169 24 : IF (natom > 1) THEN
170 22 : xval = 0.0_dp
171 102 : DO j = 1, natom
172 342 : DO k = 1, 3
173 240 : i = (j - 1)*3 + k
174 320 : xval(k) = xval(k) + dimer_env%nvec(i)
175 : END DO
176 : END DO
177 22 : IF (unit_nr > 0) THEN
178 : WRITE (unit_nr, "(/,T2,A,T9,A)") &
179 11 : "DIMER|", "Overall translation to be removed from the initial dimer vector"
180 : WRITE (unit_nr, "(T2,A,T9,3(F12.6,3X))") &
181 11 : "DIMER|", xval(1:3)
182 : END IF
183 : ! Subtract net translations
184 88 : xval = xval/REAL(natom*3, KIND=dp)
185 102 : DO j = 1, natom
186 342 : DO k = 1, 3
187 240 : i = (j - 1)*3 + k
188 320 : dimer_env%nvec(i) = dimer_env%nvec(i) - xval(k)
189 : END DO
190 : END DO
191 : END IF
192 : ! set nvec components to zero for the corresponding constraints
193 24 : CALL dimer_fixed_atom_control(dimer_env%nvec, subsys, unit_nr)
194 : ! Normalize dimer vector
195 270 : norm = SQRT(SUM(dimer_env%nvec**2))
196 24 : IF (norm <= EPSILON(0.0_dp)) &
197 0 : CPABORT("The norm of the dimer vector is 0! Calculation cannot proceed further.")
198 24 : IF (unit_nr > 0) THEN
199 : WRITE (unit_nr, "(T2,A,T9,A,T69,F12.6)") &
200 12 : "DIMER|", "Norm of dimer vector to be normalized by rescaling", norm
201 : END IF
202 270 : dimer_env%nvec = dimer_env%nvec/norm
203 24 : dimer_env%rot%rotation_step = do_first_rotation_step
204 24 : CALL section_vals_val_get(dimer_section, "DR", r_val=dimer_env%dr)
205 : CALL section_vals_val_get(dimer_section, "INTERPOLATE_GRADIENT", &
206 24 : l_val=dimer_env%rot%interpolate_gradient)
207 : CALL section_vals_val_get(dimer_section, "ANGLE_TOLERANCE", &
208 24 : r_val=dimer_env%rot%angle_tol)
209 : CALL section_vals_val_get(dimer_section, "K-DIMER", &
210 24 : l_val=dimer_env%kdimer)
211 : CALL section_vals_val_get(dimer_section, "BETA", &
212 24 : r_val=dimer_env%beta)
213 : ! initialise values
214 24 : dimer_env%cg_rot%norm_h = 1.0_dp
215 270 : dimer_env%rot%g0 = 0.0_dp
216 270 : dimer_env%rot%g1 = 0.0_dp
217 270 : dimer_env%rot%g1p = 0.0_dp
218 48 : ALLOCATE (dimer_env%cg_rot%nvec_old(natom*3))
219 24 : END SUBROUTINE dimer_env_create
220 :
221 : ! **************************************************************************************************
222 : !> \brief ...
223 : !> \param dimer_env ...
224 : !> \par History
225 : !> none
226 : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
227 : ! **************************************************************************************************
228 24 : SUBROUTINE dimer_env_retain(dimer_env)
229 : TYPE(dimer_env_type), POINTER :: dimer_env
230 :
231 24 : CPASSERT(ASSOCIATED(dimer_env))
232 24 : CPASSERT(dimer_env%ref_count > 0)
233 24 : dimer_env%ref_count = dimer_env%ref_count + 1
234 24 : END SUBROUTINE dimer_env_retain
235 :
236 : ! **************************************************************************************************
237 : !> \brief ...
238 : !> \param dimer_env ...
239 : !> \par History
240 : !> none
241 : !> \author Luca Bellucci and Teodoro Laino - created [tlaino] - 01.2008
242 : ! **************************************************************************************************
243 1121 : SUBROUTINE dimer_env_release(dimer_env)
244 : TYPE(dimer_env_type), POINTER :: dimer_env
245 :
246 1121 : IF (ASSOCIATED(dimer_env)) THEN
247 48 : CPASSERT(dimer_env%ref_count > 0)
248 48 : dimer_env%ref_count = dimer_env%ref_count - 1
249 48 : IF (dimer_env%ref_count == 0) THEN
250 24 : IF (ASSOCIATED(dimer_env%nvec)) THEN
251 24 : DEALLOCATE (dimer_env%nvec)
252 : END IF
253 24 : IF (ASSOCIATED(dimer_env%rot%g0)) THEN
254 24 : DEALLOCATE (dimer_env%rot%g0)
255 : END IF
256 24 : IF (ASSOCIATED(dimer_env%rot%g1)) THEN
257 24 : DEALLOCATE (dimer_env%rot%g1)
258 : END IF
259 24 : IF (ASSOCIATED(dimer_env%rot%g1p)) THEN
260 24 : DEALLOCATE (dimer_env%rot%g1p)
261 : END IF
262 24 : IF (ASSOCIATED(dimer_env%cg_rot%nvec_old)) THEN
263 24 : DEALLOCATE (dimer_env%cg_rot%nvec_old)
264 : END IF
265 : ! No need to deallocate tls_vec (just a pointer to aother local array)
266 24 : NULLIFY (dimer_env%tsl%tls_vec)
267 24 : DEALLOCATE (dimer_env)
268 : END IF
269 : END IF
270 1121 : END SUBROUTINE dimer_env_release
271 :
272 : ! **************************************************************************************************
273 : !> \brief Set parts of a given array vec to zero according to fixed atom constraints.
274 : !> When atoms are (partially) fixed then the relevant components of
275 : !> nvec should be set to zero. Furthermore, the relevant components
276 : !> of the gradient in CG should also be set to zero.
277 : !> \param vec : vector to be modified
278 : !> \param subsys : subsys type object used by CP2k
279 : !> \param unit : unit to write message to
280 : !> \par History
281 : !> 2016/03/03 [LTong] created
282 : !> \author Lianheng Tong [LTong]
283 : ! **************************************************************************************************
284 432 : SUBROUTINE dimer_fixed_atom_control(vec, subsys, unit)
285 : REAL(KIND=dp), DIMENSION(:), INTENT(INOUT) :: vec
286 : TYPE(cp_subsys_type), POINTER :: subsys
287 : INTEGER, INTENT(IN), OPTIONAL :: unit
288 :
289 : INTEGER :: ii, ikind, ind, iparticle, nfixed_atoms, &
290 : nkinds
291 432 : TYPE(fixd_constraint_type), DIMENSION(:), POINTER :: fixd_list
292 : TYPE(molecule_kind_list_type), POINTER :: molecule_kinds
293 432 : TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
294 : TYPE(molecule_kind_type), POINTER :: molecule_kind
295 :
296 432 : NULLIFY (molecule_kinds, molecule_kind, fixd_list)
297 :
298 : ! need to get constraint information from molecule information
299 : CALL cp_subsys_get(subsys=subsys, &
300 432 : molecule_kinds=molecule_kinds)
301 432 : molecule_kind_set => molecule_kinds%els
302 :
303 : ! get total number of fixed atoms
304 : ! nkinds is the kinds of molecules, not atoms
305 432 : nkinds = molecule_kinds%n_els
306 1158 : DO ikind = 1, nkinds
307 726 : molecule_kind => molecule_kind_set(ikind)
308 : CALL get_molecule_kind(molecule_kind, &
309 : nfixd=nfixed_atoms, &
310 726 : fixd_list=fixd_list)
311 1158 : IF (ASSOCIATED(fixd_list)) THEN
312 168 : IF (PRESENT(unit) .AND. unit > 0) THEN
313 : WRITE (unit, "(/,T2,A,T9,A,T71,I10)") &
314 12 : "DIMER|", "Number of fixed atoms to adjust dimer vector", nfixed_atoms
315 : END IF
316 224 : DO ii = 1, nfixed_atoms
317 224 : IF (.NOT. fixd_list(ii)%restraint%active) THEN
318 56 : iparticle = fixd_list(ii)%fixd
319 56 : ind = (iparticle - 1)*3
320 : ! apply constraint to nvec
321 56 : SELECT CASE (fixd_list(ii)%itype)
322 : CASE (use_perd_x)
323 0 : vec(ind + 1) = 0.0_dp
324 : CASE (use_perd_y)
325 0 : vec(ind + 2) = 0.0_dp
326 : CASE (use_perd_z)
327 0 : vec(ind + 3) = 0.0_dp
328 : CASE (use_perd_xy)
329 0 : vec(ind + 1) = 0.0_dp
330 0 : vec(ind + 2) = 0.0_dp
331 : CASE (use_perd_xz)
332 0 : vec(ind + 1) = 0.0_dp
333 0 : vec(ind + 3) = 0.0_dp
334 : CASE (use_perd_yz)
335 0 : vec(ind + 2) = 0.0_dp
336 0 : vec(ind + 3) = 0.0_dp
337 : CASE (use_perd_xyz)
338 56 : vec(ind + 1) = 0.0_dp
339 56 : vec(ind + 2) = 0.0_dp
340 56 : vec(ind + 3) = 0.0_dp
341 : END SELECT
342 : END IF ! .NOT.fixd_list(ii)%restraint%active
343 : END DO ! ii
344 : END IF ! ASSOCIATED(fixd_list)
345 : END DO ! ikind
346 432 : END SUBROUTINE dimer_fixed_atom_control
347 :
348 : ! **************************************************************************************************
349 : !> \brief Read dimer vector from input and store into dimer_env%nvec
350 : !> \param dimer_env ...
351 : !> \param dimer_section ...
352 : !> \param natom ...
353 : !> \param unit ...
354 : !> \param globenv ...
355 : !> \param force_env ...
356 : !> \par History
357 : !> 2026/05 created
358 : !> \author HE Zilong
359 : ! **************************************************************************************************
360 24 : SUBROUTINE dimer_init_vector(dimer_env, dimer_section, natom, unit, globenv, force_env)
361 : TYPE(dimer_env_type), INTENT(INOUT), POINTER :: dimer_env
362 : TYPE(section_vals_type), INTENT(INOUT), POINTER :: dimer_section
363 : INTEGER, INTENT(IN) :: natom, unit
364 : TYPE(global_environment_type), INTENT(IN), POINTER :: globenv
365 : TYPE(force_env_type), POINTER :: force_env
366 :
367 : CHARACTER(LEN=17) :: vib_name
368 : CHARACTER(LEN=default_path_length) :: molden_name
369 : INTEGER :: dimer_init_method, i, ierr, isize, j, &
370 : n_rep_val, nvec_size, vib_list_size
371 24 : INTEGER, ALLOCATABLE, DIMENSION(:) :: vib_id_list
372 : LOGICAL :: explicit, found
373 48 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: array_r, vib_wt_list
374 24 : REAL(KIND=dp), DIMENSION(:), POINTER :: array
375 : TYPE(cp_parser_type) :: parser
376 : TYPE(mp_para_env_type), POINTER :: para_env
377 : TYPE(section_vals_type), POINTER :: nvec_section
378 :
379 24 : para_env => force_env%para_env
380 48 : nvec_section => section_vals_get_subs_vals(dimer_section, "DIMER_VECTOR")
381 24 : NULLIFY (array)
382 24 : nvec_size = natom*3
383 :
384 72 : ALLOCATE (dimer_env%nvec(nvec_size))
385 270 : dimer_env%nvec(:) = 0.0_dp
386 24 : CALL section_vals_get(nvec_section, explicit=explicit)
387 24 : IF (explicit) THEN
388 4 : IF (unit > 0) WRITE (unit, "(/,T2,A,T9,A)") &
389 2 : "DIMER|", "Initial dimer vector read from input DIMER_VECTOR section"
390 4 : isize = 0
391 4 : CALL section_vals_val_get(nvec_section, "_DEFAULT_KEYWORD_", n_rep_val=n_rep_val)
392 10 : DO i = 1, n_rep_val
393 6 : CALL section_vals_val_get(nvec_section, "_DEFAULT_KEYWORD_", r_vals=array, i_rep_val=i)
394 34 : DO j = 1, SIZE(array)
395 24 : isize = isize + 1
396 30 : IF (isize <= nvec_size) THEN
397 24 : dimer_env%nvec(isize) = array(j)
398 : ELSE
399 0 : CPABORT("Size of input DIMER_VECTOR more than natom * 3")
400 : END IF
401 : END DO
402 : END DO
403 4 : IF (isize /= nvec_size) THEN
404 0 : CPABORT("Size of input DIMER_VECTOR inconsistent with natom * 3")
405 : END IF
406 : ELSE
407 20 : CALL section_vals_val_get(dimer_section, "INITIALIZATION_METHOD", i_val=dimer_init_method)
408 18 : SELECT CASE (dimer_init_method)
409 : CASE (dimer_init_random)
410 18 : IF (unit > 0) WRITE (unit, "(/,T2,A,T9,A)") &
411 9 : "DIMER|", "Initial dimer vector generated randomly"
412 18 : CALL globenv%gaussian_rng_stream%fill(dimer_env%nvec)
413 : CASE (dimer_init_molden)
414 2 : CALL section_vals_val_get(dimer_section, "VIB_MOLDEN_NAME", explicit=explicit)
415 2 : IF (.NOT. explicit) &
416 0 : CPABORT("INITIALIZATION_METHOD MOLDEN requires VIB_MOLDEN_NAME")
417 2 : CALL section_vals_val_get(dimer_section, "VIB_MOLDEN_NAME", c_val=molden_name)
418 2 : IF (unit > 0) THEN
419 : WRITE (unit, "(/,T2,A,T9,A)") &
420 1 : "DIMER|", "Initial dimer vector by linear combination of normal modes from:"
421 : WRITE (unit, "(T2,A,T9,A)") &
422 1 : "DIMER|", TRIM(ADJUSTL(molden_name))
423 : END IF
424 : ! Build lists of indices and weights for linear combination
425 : CALL vib_get_index_weight(dimer_section, vib_id_list, vib_wt_list, &
426 2 : vib_list_size, unit)
427 : ! Find vibrational modes in molden and do linear combination
428 4 : ALLOCATE (array_r(nvec_size))
429 2 : CALL parser_create(parser, molden_name, para_env=para_env, apply_preprocessing=.FALSE.)
430 4 : Vib_modes: DO i = 1, vib_list_size
431 : ! Format below is from molden_utils.F
432 2 : WRITE (UNIT=vib_name, FMT='(T2,A,1X,I6)') "vibration", vib_id_list(i)
433 : CALL parser_search_string(parser, TRIM(vib_name), &
434 : ignore_case=.FALSE., found=found, &
435 2 : begin_line=.TRUE., search_from_begin_of_file=.TRUE.)
436 2 : IF (.NOT. found) CALL cp_abort(__LOCATION__, &
437 0 : "Could not found <"//vib_name//"> from molden file")
438 12 : DO j = 1, natom
439 10 : CALL parser_get_next_line(parser, 1)
440 : READ (UNIT=parser%input_line, FMT=*, IOSTAT=ierr) &
441 10 : array_r(3*j - 2), array_r(3*j - 1), array_r(3*j)
442 10 : IF (ierr /= 0) &
443 : CALL cp_abort(__LOCATION__, &
444 : "Error while reading MOLDEN file: cannot parse the line "// &
445 : TRIM(ADJUSTL(cp_to_string(j)))//" of <"//vib_name//"> "// &
446 2 : "for components of the normal mode")
447 : END DO
448 34 : dimer_env%nvec(:) = dimer_env%nvec(:) + array_r(:)*vib_wt_list(i)
449 : END DO Vib_modes
450 2 : CALL parser_release(parser)
451 2 : IF (ALLOCATED(array_r)) DEALLOCATE (array_r)
452 2 : IF (ALLOCATED(vib_id_list)) DEALLOCATE (vib_id_list)
453 4 : IF (ALLOCATED(vib_wt_list)) DEALLOCATE (vib_wt_list)
454 : CASE DEFAULT
455 20 : CPABORT("Invalid or not yet implemented dimer initialization method")
456 : END SELECT
457 : END IF
458 :
459 120 : END SUBROUTINE dimer_init_vector
460 :
461 : ! **************************************************************************************************
462 : !> \brief Read VIB_INDEX and VIB_WEIGHT under DIMER section
463 : !> \param dimer_section ...
464 : !> \param vib_id_list ...
465 : !> \param vib_wt_list ...
466 : !> \param vib_list_size ...
467 : !> \param unit ...
468 : !> \par History
469 : !> 2026/05 created
470 : !> \author HE Zilong
471 : ! **************************************************************************************************
472 2 : SUBROUTINE vib_get_index_weight(dimer_section, vib_id_list, vib_wt_list, vib_list_size, unit)
473 : TYPE(section_vals_type), INTENT(IN), POINTER :: dimer_section
474 : INTEGER, ALLOCATABLE, DIMENSION(:), INTENT(OUT) :: vib_id_list
475 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:), &
476 : INTENT(OUT) :: vib_wt_list
477 : INTEGER, INTENT(OUT) :: vib_list_size
478 : INTEGER, INTENT(IN) :: unit
479 :
480 : INTEGER :: i, isize, n_rep_val, vib_id_size, &
481 : vib_wt_size
482 2 : INTEGER, DIMENSION(:), POINTER :: tmpilist
483 2 : REAL(KIND=dp), DIMENSION(:), POINTER :: tmprlist
484 :
485 2 : vib_id_size = 0
486 2 : CALL section_vals_val_get(dimer_section, "VIB_INDEX", n_rep_val=n_rep_val)
487 4 : DO i = 1, n_rep_val
488 : CALL section_vals_val_get(dimer_section, "VIB_INDEX", &
489 2 : i_vals=tmpilist, i_rep_val=i)
490 4 : vib_id_size = vib_id_size + SIZE(tmpilist)
491 : END DO
492 6 : ALLOCATE (vib_id_list(vib_id_size))
493 2 : isize = 0
494 4 : DO i = 1, n_rep_val
495 : CALL section_vals_val_get(dimer_section, "VIB_INDEX", &
496 2 : i_vals=tmpilist, i_rep_val=i)
497 4 : vib_id_list(isize + 1:isize + SIZE(tmpilist)) = tmpilist
498 4 : isize = isize + SIZE(tmpilist)
499 : END DO
500 :
501 2 : vib_wt_size = 0
502 2 : CALL section_vals_val_get(dimer_section, "VIB_WEIGHT", n_rep_val=n_rep_val)
503 4 : DO i = 1, n_rep_val
504 : CALL section_vals_val_get(dimer_section, "VIB_WEIGHT", &
505 2 : r_vals=tmprlist, i_rep_val=i)
506 4 : vib_wt_size = vib_wt_size + SIZE(tmprlist)
507 : END DO
508 6 : ALLOCATE (vib_wt_list(vib_wt_size))
509 2 : isize = 0
510 4 : DO i = 1, n_rep_val
511 : CALL section_vals_val_get(dimer_section, "VIB_WEIGHT", &
512 2 : r_vals=tmprlist, i_rep_val=i)
513 4 : vib_wt_list(isize + 1:isize + SIZE(tmprlist)) = tmprlist
514 4 : isize = isize + SIZE(tmprlist)
515 : END DO
516 :
517 2 : IF (vib_id_size /= vib_wt_size) THEN
518 : CALL cp_abort(__LOCATION__, &
519 : "Inconsistent count of values speficied in input between "// &
520 : "VIB_INDEX ("//TRIM(ADJUSTL(cp_to_string(vib_id_size)))//") "// &
521 0 : "and VIB_WEIGHT ("//TRIM(ADJUSTL(cp_to_string(vib_wt_size)))//")")
522 : ELSE
523 2 : vib_list_size = vib_id_size
524 : END IF
525 2 : IF (unit > 0) THEN
526 : WRITE (unit, "(/,T2,A,T9,A)") &
527 1 : "DIMER|", "Indices and weights of vibrational modes in linear combination"
528 2 : DO i = 1, vib_list_size
529 : WRITE (unit, "(T2,A,T10,I6,T18,F12.6)") &
530 2 : "DIMER|", vib_id_list(i), vib_wt_list(i)
531 : END DO
532 : END IF
533 :
534 4 : END SUBROUTINE vib_get_index_weight
535 :
536 0 : END MODULE dimer_types
|
