(git:d657d38)
Loading...
Searching...
No Matches
bfgs_optimizer.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Routines for Geometry optimization using BFGS algorithm
10!> \par History
11!> Module modified by Pierre-André Cazade [pcazade] 01.2020 - University of Limerick.
12!> Modifications enable Space Group Symmetry.
13! **************************************************************************************************
14MODULE bfgs_optimizer
15
16 USE atomic_kind_list_types, ONLY: atomic_kind_list_type
17 USE atomic_kind_types, ONLY: get_atomic_kind, &
18 get_atomic_kind_set
19 USE cell_types, ONLY: cell_type, &
20 pbc
21 USE constraint_fxd, ONLY: fix_atom_control
22 USE cp_blacs_env, ONLY: cp_blacs_env_create, &
23 cp_blacs_env_release, &
24 cp_blacs_env_type
25 USE cp_external_control, ONLY: external_control
26 USE cp_files, ONLY: close_file, &
27 open_file
28 USE cp_fm_basic_linalg, ONLY: cp_fm_column_scale, &
29 cp_fm_transpose
30 USE cp_fm_diag, ONLY: choose_eigv_solver
31 USE cp_fm_struct, ONLY: cp_fm_struct_create, &
32 cp_fm_struct_release, &
33 cp_fm_struct_type
34 USE cp_fm_types, ONLY: &
35 cp_fm_get_info, &
36 cp_fm_read_unformatted, &
37 cp_fm_set_all, &
38 cp_fm_to_fm, &
39 cp_fm_type, &
40 cp_fm_write_unformatted, cp_fm_create, cp_fm_release
41 USE parallel_gemm_api, ONLY: parallel_gemm
42 USE cp_log_handling, ONLY: cp_get_default_logger, &
43 cp_logger_type, &
44 cp_to_string
45 USE cp_output_handling, ONLY: cp_iterate, &
46 cp_p_file, &
47 cp_print_key_finished_output, &
48 cp_print_key_should_output, &
49 cp_print_key_unit_nr
50 USE message_passing, ONLY: mp_para_env_type
51 USE cp_subsys_types, ONLY: cp_subsys_get, &
52 cp_subsys_type
53 USE force_env_types, ONLY: force_env_get, &
54 force_env_type
55 USE global_types, ONLY: global_environment_type
56 USE gopt_f_methods, ONLY: gopt_f_ii, &
57 gopt_f_io, &
58 gopt_f_io_finalize, &
59 gopt_f_io_init, &
60 print_geo_opt_header, &
61 print_geo_opt_nc
62 USE gopt_f_types, ONLY: gopt_f_type
63 USE gopt_param_types, ONLY: gopt_param_type
64 USE input_constants, ONLY: default_cell_method_id, &
65 default_ts_method_id
66 USE input_section_types, ONLY: section_vals_get_subs_vals, &
67 section_vals_type, &
68 section_vals_val_get, &
69 section_vals_val_set
70 USE kinds, ONLY: default_path_length, &
71 dp
72 USE machine, ONLY: m_flush, &
73 m_walltime
74 USE particle_list_types, ONLY: particle_list_type
75 USE space_groups, ONLY: identify_space_group, &
76 print_spgr, &
77 spgr_apply_rotations_coord, &
78 spgr_apply_rotations_force
79 USE space_groups_types, ONLY: spgr_type
80#include "../base/base_uses.f90"
81
82 IMPLICIT NONE
83 PRIVATE
84
85
86! **************************************************************************************************
87!> \brief evaluate the potential energy and its gradients using an array
88!> with same dimension as the particle_set
89!> \param gopt_env the geometry optimization environment
90!> \param x the position where the function should be evaluated
91!> \param f the function value
92!> \param gradient the value of its gradient
93!> \par History
94!> none
95!> \author Teodoro Laino [tlaino] - University of Zurich - 01.2008
96! **************************************************************************************************
97INTERFACE
98
99 SUBROUTINE cp_eval_at(gopt_env, x, f, gradient, master, &
100 final_evaluation, para_env)
101
102 USE message_passing, ONLY: mp_para_env_type
103 USE gopt_f_types, ONLY: gopt_f_type
104 USE kinds, ONLY: dp
105
106 TYPE(gopt_f_type), POINTER :: gopt_env
107 REAL(KIND=dp), DIMENSION(:), POINTER :: x
108 REAL(KIND=dp), INTENT(out), OPTIONAL :: f
109 REAL(KIND=dp), DIMENSION(:), OPTIONAL, &
110 POINTER :: gradient
111 INTEGER, INTENT(IN) :: master
112 LOGICAL, INTENT(IN), OPTIONAL :: final_evaluation
113 TYPE(mp_para_env_type), POINTER :: para_env
114
115 END SUBROUTINE cp_eval_at
116
117END INTERFACE
118
119 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'bfgs_optimizer'
120 LOGICAL, PARAMETER :: debug_this_module = .true.
121
122 PUBLIC :: geoopt_bfgs
123
124CONTAINS
125
126! **************************************************************************************************
127!> \brief Main driver for BFGS geometry optimizations
128!> \param force_env ...
129!> \param gopt_param ...
130!> \param globenv ...
131!> \param geo_section ...
132!> \param gopt_env ...
133!> \param x0 ...
134!> \par History
135!> 01.2020 modified to perform Space Group Symmetry [pcazade]
136! **************************************************************************************************
137 RECURSIVE SUBROUTINE geoopt_bfgs(force_env, gopt_param, globenv, geo_section, gopt_env, x0)
138
139 TYPE(force_env_type), POINTER :: force_env
140 TYPE(gopt_param_type), POINTER :: gopt_param
141 TYPE(global_environment_type), POINTER :: globenv
142 TYPE(section_vals_type), POINTER :: geo_section
143 TYPE(gopt_f_type), POINTER :: gopt_env
144 REAL(kind=dp), DIMENSION(:), POINTER :: x0
145
146 CHARACTER(len=*), PARAMETER :: routinen = 'geoopt_bfgs'
147 REAL(kind=dp), PARAMETER :: one = 1.0_dp, zero = 0.0_dp
148
149 CHARACTER(LEN=5) :: wildcard
150 CHARACTER(LEN=default_path_length) :: hes_filename
151 INTEGER :: handle, hesunit_read, indf, info, &
152 iter_nr, its, maxiter, ndf, nfree, &
153 output_unit
154 LOGICAL :: conv, hesrest, ionode, shell_present, &
155 should_stop, use_mod_hes, use_rfo
156 REAL(kind=dp) :: ediff, emin, eold, etot, pred, rad, rat, &
157 step, t_diff, t_now, t_old
158 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: dg, dr, dx, eigval, gold, work, xold
159 REAL(kind=dp), DIMENSION(:), POINTER :: g
160 TYPE(atomic_kind_list_type), POINTER :: atomic_kinds
161 TYPE(cp_blacs_env_type), POINTER :: blacs_env
162 TYPE(cp_fm_struct_type), POINTER :: fm_struct_hes
163 TYPE(cp_fm_type) :: eigvec_mat, hess_mat, hess_tmp
164 TYPE(cp_logger_type), POINTER :: logger
165 TYPE(mp_para_env_type), POINTER :: para_env
166 TYPE(cp_subsys_type), POINTER :: subsys
167 TYPE(section_vals_type), POINTER :: print_key, root_section
168 TYPE(spgr_type), POINTER :: spgr
169
170 NULLIFY (logger, g, blacs_env, spgr)
171 logger => cp_get_default_logger()
172 para_env => force_env%para_env
173 root_section => force_env%root_section
174 spgr => gopt_env%spgr
175 t_old = m_walltime()
176
177 CALL timeset(routinen, handle)
178 CALL section_vals_val_get(geo_section, "BFGS%TRUST_RADIUS", r_val=rad)
179 print_key => section_vals_get_subs_vals(geo_section, "BFGS%RESTART")
180 ionode = para_env%is_source()
181 maxiter = gopt_param%max_iter
182 conv = .false.
183 rat = 0.0_dp
184 wildcard = " BFGS"
185 hes_filename = ""
186
187 ! Stop if not yet implemented
188 SELECT CASE (gopt_env%type_id)
189 CASE (default_ts_method_id)
190 cpabort("BFGS method not yet working with DIMER")
191 END SELECT
192
193 CALL section_vals_val_get(geo_section, "BFGS%USE_RAT_FUN_OPT", l_val=use_rfo)
194 CALL section_vals_val_get(geo_section, "BFGS%USE_MODEL_HESSIAN", l_val=use_mod_hes)
195 CALL section_vals_val_get(geo_section, "BFGS%RESTART_HESSIAN", l_val=hesrest)
196 output_unit = cp_print_key_unit_nr(logger, geo_section, "PRINT%PROGRAM_RUN_INFO", &
197 extension=".geoLog")
198 IF (output_unit > 0) THEN
199 IF (use_rfo) THEN
200 WRITE (unit=output_unit, fmt="(/,T2,A,T78,A3)") &
201 "BFGS| Use rational function optimization for step estimation: ", "YES"
202 ELSE
203 WRITE (unit=output_unit, fmt="(/,T2,A,T78,A3)") &
204 "BFGS| Use rational function optimization for step estimation: ", " NO"
205 END IF
206 IF (use_mod_hes) THEN
207 WRITE (unit=output_unit, fmt="(T2,A,T78,A3)") &
208 "BFGS| Use model Hessian for initial guess: ", "YES"
209 ELSE
210 WRITE (unit=output_unit, fmt="(T2,A,T78,A3)") &
211 "BFGS| Use model Hessian for initial guess: ", " NO"
212 END IF
213 IF (hesrest) THEN
214 WRITE (unit=output_unit, fmt="(T2,A,T78,A3)") &
215 "BFGS| Restart Hessian: ", "YES"
216 ELSE
217 WRITE (unit=output_unit, fmt="(T2,A,T78,A3)") &
218 "BFGS| Restart Hessian: ", " NO"
219 END IF
220 WRITE (unit=output_unit, fmt="(T2,A,T61,F20.3)") &
221 "BFGS| Trust radius: ", rad
222 END IF
223
224 ndf = SIZE(x0)
225 nfree = gopt_env%nfree
226 IF (ndf > 3000) &
227 CALL cp_warn(__location__, &
228 "The dimension of the Hessian matrix ("// &
229 trim(adjustl(cp_to_string(ndf)))//") is greater than 3000. "// &
230 "The diagonalisation of the full Hessian matrix needed for BFGS "// &
231 "is computationally expensive. You should consider to use the linear "// &
232 "scaling variant L-BFGS instead.")
233
234 ! Initialize hessian (hes = unitary matrix or model hessian )
235 CALL cp_blacs_env_create(blacs_env, para_env, globenv%blacs_grid_layout, &
236 globenv%blacs_repeatable)
237 CALL cp_fm_struct_create(fm_struct_hes, para_env=para_env, context=blacs_env, &
238 nrow_global=ndf, ncol_global=ndf)
239 CALL cp_fm_create(hess_mat, fm_struct_hes, name="hess_mat")
240 CALL cp_fm_create(hess_tmp, fm_struct_hes, name="hess_tmp")
241 CALL cp_fm_create(eigvec_mat, fm_struct_hes, name="eigvec_mat")
242 ALLOCATE (eigval(ndf))
243 eigval(:) = 0.0_dp
244
245 CALL force_env_get(force_env=force_env, subsys=subsys)
246 CALL cp_subsys_get(subsys, atomic_kinds=atomic_kinds)
247 CALL get_atomic_kind_set(atomic_kind_set=atomic_kinds%els, shell_present=shell_present)
248 IF (use_mod_hes) THEN
249 IF (shell_present) THEN
250 CALL cp_warn(__location__, &
251 "No model Hessian is available for core-shell models. "// &
252 "A unit matrix is used as the initial Hessian.")
253 use_mod_hes = .false.
254 END IF
255 IF (gopt_env%type_id == default_cell_method_id) THEN
256 CALL cp_warn(__location__, &
257 "No model Hessian is available for cell optimizations. "// &
258 "A unit matrix is used as the initial Hessian.")
259 use_mod_hes = .false.
260 END IF
261 END IF
262
263 IF (use_mod_hes) THEN
264 CALL cp_fm_set_all(hess_mat, alpha=zero, beta=0.00_dp)
265 CALL construct_initial_hess(gopt_env%force_env, hess_mat)
266 CALL cp_fm_to_fm(hess_mat, hess_tmp)
267 CALL choose_eigv_solver(hess_tmp, eigvec_mat, eigval, info=info)
268 ! In rare cases the diagonalization of hess_mat fails (bug in scalapack?)
269 IF (info /= 0) THEN
270 CALL cp_fm_set_all(hess_mat, alpha=zero, beta=one)
271 IF (output_unit > 0) WRITE (output_unit, *) &
272 "BFGS: Matrix diagonalization failed, using unity as model Hessian."
273 ELSE
274 DO its = 1, SIZE(eigval)
275 IF (eigval(its) .LT. 0.1_dp) eigval(its) = 0.1_dp
276 END DO
277 CALL cp_fm_to_fm(eigvec_mat, hess_tmp)
278 CALL cp_fm_column_scale(eigvec_mat, eigval)
279 CALL parallel_gemm("N", "T", ndf, ndf, ndf, one, hess_tmp, eigvec_mat, zero, hess_mat)
280 END IF
281 ELSE
282 CALL cp_fm_set_all(hess_mat, alpha=zero, beta=one)
283 END IF
284
285 ALLOCATE (xold(ndf))
286 xold(:) = x0(:)
287
288 ALLOCATE (g(ndf))
289 g(:) = 0.0_dp
290
291 ALLOCATE (gold(ndf))
292 gold(:) = 0.0_dp
293
294 ALLOCATE (dx(ndf))
295 dx(:) = 0.0_dp
296
297 ALLOCATE (dg(ndf))
298 dg(:) = 0.0_dp
299
300 ALLOCATE (work(ndf))
301 work(:) = 0.0_dp
302
303 ALLOCATE (dr(ndf))
304 dr(:) = 0.0_dp
305
306 ! find space_group
307 CALL section_vals_val_get(geo_section, "KEEP_SPACE_GROUP", l_val=spgr%keep_space_group)
308 IF (spgr%keep_space_group) THEN
309 CALL identify_space_group(subsys, geo_section, gopt_env, output_unit)
310 CALL spgr_apply_rotations_coord(spgr, x0)
311 CALL print_spgr(spgr)
312 END IF
313
314 ! Geometry optimization starts now
315 CALL cp_iterate(logger%iter_info, increment=0, iter_nr_out=iter_nr)
316 CALL print_geo_opt_header(gopt_env, output_unit, wildcard)
317
318 ! Calculate Energy & Gradients
319 CALL cp_eval_at(gopt_env, x0, etot, g, gopt_env%force_env%para_env%mepos, &
320 .false., gopt_env%force_env%para_env)
321
322 ! Symmetrize coordinates and forces
323 IF (spgr%keep_space_group) THEN
324 CALL spgr_apply_rotations_coord(spgr, x0)
325 CALL spgr_apply_rotations_force(spgr, g)
326 END IF
327
328 ! Print info at time 0
329 emin = etot
330 t_now = m_walltime()
331 t_diff = t_now - t_old
332 t_old = t_now
333 CALL gopt_f_io_init(gopt_env, output_unit, etot, wildcard=wildcard, its=iter_nr, used_time=t_diff)
334 DO its = iter_nr + 1, maxiter
335 CALL cp_iterate(logger%iter_info, last=(its == maxiter))
336 CALL section_vals_val_set(geo_section, "STEP_START_VAL", i_val=its)
337 CALL gopt_f_ii(its, output_unit)
338
339 ! Hessian update/restarting
340 IF (((its - iter_nr) == 1) .AND. hesrest) THEN
341 IF (ionode) THEN
342 CALL section_vals_val_get(geo_section, "BFGS%RESTART_FILE_NAME", c_val=hes_filename)
343 IF (len_trim(hes_filename) == 0) THEN
344 ! Set default Hessian restart file name if no file name is defined
345 hes_filename = trim(logger%iter_info%project_name)//"-BFGS.Hessian"
346 END IF
347 IF (output_unit > 0) THEN
348 WRITE (unit=output_unit, fmt="(/,T2,A)") &
349 "BFGS| Checking for Hessian restart file <"//trim(adjustl(hes_filename))//">"
350 END IF
351 CALL open_file(file_name=trim(hes_filename), file_status="OLD", &
352 file_form="UNFORMATTED", file_action="READ", unit_number=hesunit_read)
353 IF (output_unit > 0) THEN
354 WRITE (unit=output_unit, fmt="(T2,A)") &
355 "BFGS| Hessian restart file read"
356 END IF
357 END IF
358 CALL cp_fm_read_unformatted(hess_mat, hesunit_read)
359 IF (ionode) CALL close_file(unit_number=hesunit_read)
360 ELSE
361 IF ((its - iter_nr) > 1) THEN
362 ! Symmetrize old coordinates and old forces
363 IF (spgr%keep_space_group) THEN
364 CALL spgr_apply_rotations_coord(spgr, xold)
365 CALL spgr_apply_rotations_force(spgr, gold)
366 END IF
367
368 DO indf = 1, ndf
369 dx(indf) = x0(indf) - xold(indf)
370 dg(indf) = g(indf) - gold(indf)
371 END DO
372
373 CALL bfgs(ndf, dx, dg, hess_mat, work, para_env)
374
375 ! Symmetrize coordinates and forces change
376 IF (spgr%keep_space_group) THEN
377 CALL spgr_apply_rotations_force(spgr, dx)
378 CALL spgr_apply_rotations_force(spgr, dg)
379 END IF
380
381 !Possibly dump the Hessian file
382 IF (btest(cp_print_key_should_output(logger%iter_info, print_key), cp_p_file)) THEN
383 CALL write_bfgs_hessian(geo_section, hess_mat, logger)
384 END IF
385 END IF
386 END IF
387
388 ! Symmetrize coordinates and forces
389 IF (spgr%keep_space_group) THEN
390 CALL spgr_apply_rotations_coord(spgr, x0)
391 CALL spgr_apply_rotations_force(spgr, g)
392 END IF
393
394 ! Setting the present positions & gradients as old
395 xold(:) = x0
396 gold(:) = g
397
398 ! Copying hessian hes to (ndf x ndf) matrix hes_mat for diagonalization
399 CALL cp_fm_to_fm(hess_mat, hess_tmp)
400
401 CALL choose_eigv_solver(hess_tmp, eigvec_mat, eigval, info=info)
402
403 ! In rare cases the diagonalization of hess_mat fails (bug in scalapack?)
404 IF (info /= 0) THEN
405 IF (output_unit > 0) WRITE (output_unit, *) &
406 "BFGS: Matrix diagonalization failed, resetting Hessian to unity."
407 CALL cp_fm_set_all(hess_mat, alpha=zero, beta=one)
408 CALL cp_fm_to_fm(hess_mat, hess_tmp)
409 CALL choose_eigv_solver(hess_tmp, eigvec_mat, eigval)
410 END IF
411
412 IF (use_rfo) THEN
413 CALL set_hes_eig(ndf, eigval, work)
414 dx(:) = eigval
415 CALL rat_fun_opt(ndf, dg, eigval, work, eigvec_mat, g, para_env)
416 END IF
417 CALL geoopt_get_step(ndf, eigval, eigvec_mat, hess_tmp, dr, g, para_env, use_rfo)
418
419 ! Symmetrize dr
420 IF (spgr%keep_space_group) THEN
421 CALL spgr_apply_rotations_force(spgr, dr)
422 END IF
423
424 CALL trust_radius(ndf, step, rad, rat, dr, output_unit)
425
426 ! Update the atomic positions
427 x0 = x0 + dr
428
429 ! Symmetrize coordinates
430 IF (spgr%keep_space_group) THEN
431 CALL spgr_apply_rotations_coord(spgr, x0)
432 END IF
433
434 CALL energy_predict(ndf, work, hess_mat, dr, g, conv, pred, para_env)
435 eold = etot
436
437 ! Energy & Gradients at new step
438 CALL cp_eval_at(gopt_env, x0, etot, g, gopt_env%force_env%para_env%mepos, &
439 .false., gopt_env%force_env%para_env)
440
441 ediff = etot - eold
442
443 ! Symmetrize forces
444 IF (spgr%keep_space_group) THEN
445 CALL spgr_apply_rotations_force(spgr, g)
446 END IF
447
448 ! check for an external exit command
449 CALL external_control(should_stop, "GEO", globenv=globenv)
450 IF (should_stop) EXIT
451
452 ! Some IO and Convergence check
453 t_now = m_walltime()
454 t_diff = t_now - t_old
455 t_old = t_now
456 CALL gopt_f_io(gopt_env, force_env, root_section, its, etot, output_unit, &
457 eold, emin, wildcard, gopt_param, ndf, dr, g, conv, pred, rat, &
458 step, rad, used_time=t_diff)
459
460 IF (conv .OR. (its == maxiter)) EXIT
461 IF (etot < emin) emin = etot
462 IF (use_rfo) CALL update_trust_rad(rat, rad, step, ediff)
463 END DO
464
465 IF (its == maxiter .AND. (.NOT. conv)) THEN
466 CALL print_geo_opt_nc(gopt_env, output_unit)
467 END IF
468
469 ! show space_group
470 CALL section_vals_val_get(geo_section, "SHOW_SPACE_GROUP", l_val=spgr%show_space_group)
471 IF (spgr%show_space_group) THEN
472 CALL identify_space_group(subsys, geo_section, gopt_env, output_unit)
473 CALL print_spgr(spgr)
474 END IF
475
476 ! Write final information, if converged
477 CALL cp_iterate(logger%iter_info, last=.true., increment=0)
478 CALL write_bfgs_hessian(geo_section, hess_mat, logger)
479 CALL gopt_f_io_finalize(gopt_env, force_env, x0, conv, its, root_section, &
480 gopt_env%force_env%para_env, gopt_env%force_env%para_env%mepos, output_unit)
481
482 CALL cp_fm_struct_release(fm_struct_hes)
483 CALL cp_fm_release(hess_mat)
484 CALL cp_fm_release(eigvec_mat)
485 CALL cp_fm_release(hess_tmp)
486
487 CALL cp_blacs_env_release(blacs_env)
488 DEALLOCATE (xold)
489 DEALLOCATE (g)
490 DEALLOCATE (gold)
491 DEALLOCATE (dx)
492 DEALLOCATE (dg)
493 DEALLOCATE (eigval)
494 DEALLOCATE (work)
495 DEALLOCATE (dr)
496
497 CALL cp_print_key_finished_output(output_unit, logger, geo_section, &
498 "PRINT%PROGRAM_RUN_INFO")
499 CALL timestop(handle)
500
501 END SUBROUTINE geoopt_bfgs
502
503! **************************************************************************************************
504!> \brief ...
505!> \param ndf ...
506!> \param dg ...
507!> \param eigval ...
508!> \param work ...
509!> \param eigvec_mat ...
510!> \param g ...
511!> \param para_env ...
512! **************************************************************************************************
513 SUBROUTINE rat_fun_opt(ndf, dg, eigval, work, eigvec_mat, g, para_env)
514
515 INTEGER, INTENT(IN) :: ndf
516 REAL(kind=dp), INTENT(INOUT) :: dg(ndf), eigval(ndf), work(ndf)
517 TYPE(cp_fm_type), INTENT(IN) :: eigvec_mat
518 REAL(kind=dp), INTENT(INOUT) :: g(ndf)
519 TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env
520
521 CHARACTER(LEN=*), PARAMETER :: routinen = 'rat_fun_opt'
522 REAL(kind=dp), PARAMETER :: one = 1.0_dp
523
524 INTEGER :: handle, i, indf, iref, iter, j, k, l, &
525 maxit, ncol_local, nrow_local
526 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
527 LOGICAL :: bisec, fail, set
528 REAL(kind=dp) :: fun, fun1, fun2, fun3, fung, lam1, lam2, &
529 ln, lp, ssize, step, stol
530 REAL(kind=dp), CONTIGUOUS, DIMENSION(:, :), POINTER :: local_data
531
532 CALL timeset(routinen, handle)
533
534 stol = 1.0e-8_dp
535 ssize = 0.2_dp
536 maxit = 999
537 fail = .false.
538 bisec = .false.
539
540 dg = 0._dp
541
542 CALL cp_fm_get_info(eigvec_mat, row_indices=row_indices, col_indices=col_indices, &
543 local_data=local_data, nrow_local=nrow_local, ncol_local=ncol_local)
544
545 DO i = 1, nrow_local
546 j = row_indices(i)
547 DO k = 1, ncol_local
548 l = col_indices(k)
549 dg(l) = dg(l) + local_data(i, k)*g(j)
550 END DO
551 END DO
552 CALL para_env%sum(dg)
553
554 set = .false.
555
556 DO
557
558! calculating Lambda
559
560 lp = 0.0_dp
561 iref = 1
562 ln = 0.0_dp
563 IF (eigval(iref) < 0.0_dp) ln = eigval(iref) - 0.01_dp
564
565 iter = 0
566 DO
567 iter = iter + 1
568 fun = 0.0_dp
569 fung = 0.0_dp
570 DO indf = 1, ndf
571 fun = fun + dg(indf)**2/(ln - eigval(indf))
572 fung = fung - dg(indf)**2/(ln - eigval(indf)**2)
573 END DO
574 fun = fun - ln
575 fung = fung - one
576 step = fun/fung
577 ln = ln - step
578 IF (abs(step) < stol) GOTO 200
579 IF (iter >= maxit) EXIT
580 END DO
581100 CONTINUE
582 bisec = .true.
583 iter = 0
584 maxit = 9999
585 lam1 = 0.0_dp
586 IF (eigval(iref) < 0.0_dp) lam1 = eigval(iref) - 0.01_dp
587 fun1 = 0.0_dp
588 DO indf = 1, ndf
589 fun1 = fun1 + dg(indf)**2/(lam1 - eigval(indf))
590 END DO
591 fun1 = fun1 - lam1
592 step = abs(lam1)/1000.0_dp
593 IF (step < ssize) step = ssize
594 DO
595 iter = iter + 1
596 IF (iter > maxit) THEN
597 ln = 0.0_dp
598 lp = 0.0_dp
599 fail = .true.
600 GOTO 300
601 END IF
602 fun2 = 0.0_dp
603 lam2 = lam1 - iter*step
604 DO indf = 1, ndf
605 fun2 = fun2 + eigval(indf)**2/(lam2 - eigval(indf))
606 END DO
607 fun2 = fun2 - lam2
608 IF (fun2*fun1 < 0.0_dp) THEN
609 iter = 0
610 DO
611 iter = iter + 1
612 IF (iter > maxit) THEN
613 ln = 0.0_dp
614 lp = 0.0_dp
615 fail = .true.
616 GOTO 300
617 END IF
618 step = (lam1 + lam2)/2
619 fun3 = 0.0_dp
620 DO indf = 1, ndf
621 fun3 = fun3 + dg(indf)**2/(step - eigval(indf))
622 END DO
623 fun3 = fun3 - step
624
625 IF (abs(step - lam2) < stol) THEN
626 ln = step
627 GOTO 200
628 END IF
629
630 IF (fun3*fun1 < stol) THEN
631 lam2 = step
632 ELSE
633 lam1 = step
634 END IF
635 END DO
636 END IF
637 END DO
638
639200 CONTINUE
640 IF ((ln > eigval(iref)) .OR. ((ln > 0.0_dp) .AND. &
641 (eigval(iref) > 0.0_dp))) THEN
642
643 IF (.NOT. bisec) GOTO 100
644 ln = 0.0_dp
645 lp = 0.0_dp
646 fail = .true.
647 END IF
648
649300 CONTINUE
650
651 IF (fail .AND. .NOT. set) THEN
652 set = .true.
653 DO indf = 1, ndf
654 eigval(indf) = eigval(indf)*work(indf)
655 END DO
656 cycle
657 END IF
658
659 IF (.NOT. set) THEN
660 work(1:ndf) = one
661 END IF
662
663 DO indf = 1, ndf
664 eigval(indf) = eigval(indf) - ln
665 END DO
666 EXIT
667 END DO
668
669 CALL timestop(handle)
670
671 END SUBROUTINE rat_fun_opt
672
673! **************************************************************************************************
674!> \brief ...
675!> \param ndf ...
676!> \param dx ...
677!> \param dg ...
678!> \param hess_mat ...
679!> \param work ...
680!> \param para_env ...
681! **************************************************************************************************
682 SUBROUTINE bfgs(ndf, dx, dg, hess_mat, work, para_env)
683 INTEGER, INTENT(IN) :: ndf
684 REAL(kind=dp), INTENT(INOUT) :: dx(ndf), dg(ndf)
685 TYPE(cp_fm_type), INTENT(IN) :: hess_mat
686 REAL(kind=dp), INTENT(INOUT) :: work(ndf)
687 TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env
688
689 CHARACTER(LEN=*), PARAMETER :: routinen = 'bfgs'
690 REAL(kind=dp), PARAMETER :: one = 1.0_dp, zero = 0.0_dp
691
692 INTEGER :: handle, i, j, k, l, ncol_local, &
693 nrow_local
694 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
695 REAL(kind=dp) :: ddot, dxw, gdx
696 REAL(kind=dp), CONTIGUOUS, DIMENSION(:, :), POINTER :: local_hes
697
698 CALL timeset(routinen, handle)
699
700 CALL cp_fm_get_info(hess_mat, row_indices=row_indices, col_indices=col_indices, &
701 local_data=local_hes, nrow_local=nrow_local, ncol_local=ncol_local)
702
703 work = zero
704 DO i = 1, nrow_local
705 j = row_indices(i)
706 DO k = 1, ncol_local
707 l = col_indices(k)
708 work(j) = work(j) + local_hes(i, k)*dx(l)
709 END DO
710 END DO
711
712 CALL para_env%sum(work)
713
714 gdx = ddot(ndf, dg, 1, dx, 1)
715 gdx = one/gdx
716 dxw = ddot(ndf, dx, 1, work, 1)
717 dxw = one/dxw
718
719 DO i = 1, nrow_local
720 j = row_indices(i)
721 DO k = 1, ncol_local
722 l = col_indices(k)
723 local_hes(i, k) = local_hes(i, k) + gdx*dg(j)*dg(l) - &
724 dxw*work(j)*work(l)
725 END DO
726 END DO
727
728 CALL timestop(handle)
729
730 END SUBROUTINE bfgs
731
732! **************************************************************************************************
733!> \brief ...
734!> \param ndf ...
735!> \param eigval ...
736!> \param work ...
737! **************************************************************************************************
738 SUBROUTINE set_hes_eig(ndf, eigval, work)
739 INTEGER, INTENT(IN) :: ndf
740 REAL(kind=dp), INTENT(INOUT) :: eigval(ndf), work(ndf)
741
742 CHARACTER(LEN=*), PARAMETER :: routinen = 'set_hes_eig'
743 REAL(kind=dp), PARAMETER :: max_neg = -0.5_dp, max_pos = 5.0_dp, &
744 min_eig = 0.005_dp, one = 1.0_dp
745
746 INTEGER :: handle, indf
747 LOGICAL :: neg
748
749 CALL timeset(routinen, handle)
750
751 DO indf = 1, ndf
752 IF (eigval(indf) < 0.0_dp) neg = .true.
753 IF (eigval(indf) > 1000.0_dp) eigval(indf) = 1000.0_dp
754 END DO
755 DO indf = 1, ndf
756 IF (eigval(indf) < 0.0_dp) THEN
757 IF (eigval(indf) < max_neg) THEN
758 eigval(indf) = max_neg
759 ELSE IF (eigval(indf) > -min_eig) THEN
760 eigval(indf) = -min_eig
761 END IF
762 ELSE IF (eigval(indf) < 1000.0_dp) THEN
763 IF (eigval(indf) < min_eig) THEN
764 eigval(indf) = min_eig
765 ELSE IF (eigval(indf) > max_pos) THEN
766 eigval(indf) = max_pos
767 END IF
768 END IF
769 END DO
770
771 DO indf = 1, ndf
772 IF (eigval(indf) < 0.0_dp) THEN
773 work(indf) = -one
774 ELSE
775 work(indf) = one
776 END IF
777 END DO
778
779 CALL timestop(handle)
780
781 END SUBROUTINE set_hes_eig
782
783! **************************************************************************************************
784!> \brief ...
785!> \param ndf ...
786!> \param eigval ...
787!> \param eigvec_mat ...
788!> \param hess_tmp ...
789!> \param dr ...
790!> \param g ...
791!> \param para_env ...
792!> \param use_rfo ...
793! **************************************************************************************************
794 SUBROUTINE geoopt_get_step(ndf, eigval, eigvec_mat, hess_tmp, dr, g, para_env, use_rfo)
795
796 INTEGER, INTENT(IN) :: ndf
797 REAL(kind=dp), INTENT(INOUT) :: eigval(ndf)
798 TYPE(cp_fm_type), INTENT(IN) :: eigvec_mat, hess_tmp
799 REAL(kind=dp), INTENT(INOUT) :: dr(ndf), g(ndf)
800 TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env
801 LOGICAL :: use_rfo
802
803 REAL(kind=dp), PARAMETER :: one = 1.0_dp, zero = 0.0_dp
804
805 INTEGER :: i, indf, j, k, l, ncol_local, nrow_local
806 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
807 REAL(kind=dp), CONTIGUOUS, DIMENSION(:, :), POINTER :: local_data
808 TYPE(cp_fm_struct_type), POINTER :: matrix_struct
809 TYPE(cp_fm_type) :: tmp
810
811 CALL cp_fm_to_fm(eigvec_mat, hess_tmp)
812 IF (use_rfo) THEN
813 DO indf = 1, ndf
814 eigval(indf) = one/eigval(indf)
815 END DO
816 ELSE
817 DO indf = 1, ndf
818 eigval(indf) = one/max(0.0001_dp, eigval(indf))
819 END DO
820 END IF
821
822 CALL cp_fm_column_scale(hess_tmp, eigval)
823 CALL cp_fm_get_info(eigvec_mat, matrix_struct=matrix_struct)
824 CALL cp_fm_create(tmp, matrix_struct, name="tmp")
825
826 CALL parallel_gemm("N", "T", ndf, ndf, ndf, one, hess_tmp, eigvec_mat, zero, tmp)
827
828 CALL cp_fm_transpose(tmp, hess_tmp)
829 CALL cp_fm_release(tmp)
830
831! ** New step **
832
833 CALL cp_fm_get_info(hess_tmp, row_indices=row_indices, col_indices=col_indices, &
834 local_data=local_data, nrow_local=nrow_local, ncol_local=ncol_local)
835
836 dr = 0.0_dp
837 DO i = 1, nrow_local
838 j = row_indices(i)
839 DO k = 1, ncol_local
840 l = col_indices(k)
841 dr(j) = dr(j) - local_data(i, k)*g(l)
842 END DO
843 END DO
844
845 CALL para_env%sum(dr)
846
847 END SUBROUTINE geoopt_get_step
848
849! **************************************************************************************************
850!> \brief ...
851!> \param ndf ...
852!> \param step ...
853!> \param rad ...
854!> \param rat ...
855!> \param dr ...
856!> \param output_unit ...
857! **************************************************************************************************
858 SUBROUTINE trust_radius(ndf, step, rad, rat, dr, output_unit)
859 INTEGER, INTENT(IN) :: ndf
860 REAL(kind=dp), INTENT(INOUT) :: step, rad, rat, dr(ndf)
861 INTEGER, INTENT(IN) :: output_unit
862
863 CHARACTER(LEN=*), PARAMETER :: routinen = 'trust_radius'
864 REAL(kind=dp), PARAMETER :: one = 1.0_dp
865
866 INTEGER :: handle
867 REAL(kind=dp) :: scal
868
869 CALL timeset(routinen, handle)
870
871 step = maxval(abs(dr))
872 scal = max(one, rad/step)
873
874 IF (step > rad) THEN
875 rat = rad/step
876 CALL dscal(ndf, rat, dr, 1)
877 step = rad
878 IF (output_unit > 0) THEN
879 WRITE (unit=output_unit, fmt="(/,T2,A,F8.5)") &
880 " Step is scaled; Scaling factor = ", rat
881 CALL m_flush(output_unit)
882 END IF
883 END IF
884 CALL timestop(handle)
885
886 END SUBROUTINE trust_radius
887
888! **************************************************************************************************
889!> \brief ...
890!> \param ndf ...
891!> \param work ...
892!> \param hess_mat ...
893!> \param dr ...
894!> \param g ...
895!> \param conv ...
896!> \param pred ...
897!> \param para_env ...
898! **************************************************************************************************
899 SUBROUTINE energy_predict(ndf, work, hess_mat, dr, g, conv, pred, para_env)
900
901 INTEGER, INTENT(IN) :: ndf
902 REAL(kind=dp), INTENT(INOUT) :: work(ndf)
903 TYPE(cp_fm_type), INTENT(IN) :: hess_mat
904 REAL(kind=dp), INTENT(INOUT) :: dr(ndf), g(ndf)
905 LOGICAL, INTENT(INOUT) :: conv
906 REAL(kind=dp), INTENT(INOUT) :: pred
907 TYPE(mp_para_env_type), POINTER :: para_env
908
909 CHARACTER(LEN=*), PARAMETER :: routinen = 'energy_predict'
910 REAL(kind=dp), PARAMETER :: zero = 0.0_dp
911
912 INTEGER :: handle, i, j, k, l, ncol_local, &
913 nrow_local
914 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
915 REAL(kind=dp) :: ddot, ener1, ener2
916 REAL(kind=dp), CONTIGUOUS, DIMENSION(:, :), POINTER :: local_data
917
918 CALL timeset(routinen, handle)
919
920 ener1 = ddot(ndf, g, 1, dr, 1)
921
922 CALL cp_fm_get_info(hess_mat, row_indices=row_indices, col_indices=col_indices, &
923 local_data=local_data, nrow_local=nrow_local, ncol_local=ncol_local)
924
925 work = zero
926 DO i = 1, nrow_local
927 j = row_indices(i)
928 DO k = 1, ncol_local
929 l = col_indices(k)
930 work(j) = work(j) + local_data(i, k)*dr(l)
931 END DO
932 END DO
933
934 CALL para_env%sum(work)
935 ener2 = ddot(ndf, dr, 1, work, 1)
936 pred = ener1 + 0.5_dp*ener2
937 conv = .false.
938 CALL timestop(handle)
939
940 END SUBROUTINE energy_predict
941
942! **************************************************************************************************
943!> \brief ...
944!> \param rat ...
945!> \param rad ...
946!> \param step ...
947!> \param ediff ...
948! **************************************************************************************************
949 SUBROUTINE update_trust_rad(rat, rad, step, ediff)
950
951 REAL(kind=dp), INTENT(INOUT) :: rat, rad, step, ediff
952
953 CHARACTER(LEN=*), PARAMETER :: routinen = 'update_trust_rad'
954 REAL(kind=dp), PARAMETER :: max_trust = 1.0_dp, min_trust = 0.1_dp
955
956 INTEGER :: handle
957
958 CALL timeset(routinen, handle)
959
960 IF (rat > 4.0_dp) THEN
961 IF (ediff < 0.0_dp) THEN
962 rad = step*0.5_dp
963 ELSE
964 rad = step*0.25_dp
965 END IF
966 ELSE IF (rat > 2.0_dp) THEN
967 IF (ediff < 0.0_dp) THEN
968 rad = step*0.75_dp
969 ELSE
970 rad = step*0.5_dp
971 END IF
972 ELSE IF (rat > 4.0_dp/3.0_dp) THEN
973 IF (ediff < 0.0_dp) THEN
974 rad = step
975 ELSE
976 rad = step*0.75_dp
977 END IF
978 ELSE IF (rat > 10.0_dp/9.0_dp) THEN
979 IF (ediff < 0.0_dp) THEN
980 rad = step*1.25_dp
981 ELSE
982 rad = step
983 END IF
984 ELSE IF (rat > 0.9_dp) THEN
985 IF (ediff < 0.0_dp) THEN
986 rad = step*1.5_dp
987 ELSE
988 rad = step*1.25_dp
989 END IF
990 ELSE IF (rat > 0.75_dp) THEN
991 IF (ediff < 0.0_dp) THEN
992 rad = step*1.25_dp
993 ELSE
994 rad = step
995 END IF
996 ELSE IF (rat > 0.5_dp) THEN
997 IF (ediff < 0.0_dp) THEN
998 rad = step
999 ELSE
1000 rad = step*0.75_dp
1001 END IF
1002 ELSE IF (rat > 0.25_dp) THEN
1003 IF (ediff < 0.0_dp) THEN
1004 rad = step*0.75_dp
1005 ELSE
1006 rad = step*0.5_dp
1007 END IF
1008 ELSE IF (ediff < 0.0_dp) THEN
1009 rad = step*0.5_dp
1010 ELSE
1011 rad = step*0.25_dp
1012 END IF
1013
1014 rad = max(rad, min_trust)
1015 rad = min(rad, max_trust)
1016 CALL timestop(handle)
1017
1018 END SUBROUTINE update_trust_rad
1019
1020! **************************************************************************************************
1021
1022! **************************************************************************************************
1023!> \brief ...
1024!> \param geo_section ...
1025!> \param hess_mat ...
1026!> \param logger ...
1027! **************************************************************************************************
1028 SUBROUTINE write_bfgs_hessian(geo_section, hess_mat, logger)
1029
1030 TYPE(section_vals_type), POINTER :: geo_section
1031 TYPE(cp_fm_type), INTENT(IN) :: hess_mat
1032 TYPE(cp_logger_type), POINTER :: logger
1033
1034 CHARACTER(LEN=*), PARAMETER :: routinen = 'write_bfgs_hessian'
1035
1036 INTEGER :: handle, hesunit
1037
1038 CALL timeset(routinen, handle)
1039
1040 hesunit = cp_print_key_unit_nr(logger, geo_section, "BFGS%RESTART", &
1041 extension=".Hessian", file_form="UNFORMATTED", file_action="WRITE", &
1042 file_position="REWIND")
1043
1044 CALL cp_fm_write_unformatted(hess_mat, hesunit)
1045
1046 CALL cp_print_key_finished_output(hesunit, logger, geo_section, "BFGS%RESTART")
1047
1048 CALL timestop(handle)
1049
1050 END SUBROUTINE write_bfgs_hessian
1051! **************************************************************************************************
1052
1053! **************************************************************************************************
1054!> \brief ...
1055!> \param force_env ...
1056!> \param hess_mat ...
1057! **************************************************************************************************
1058 SUBROUTINE construct_initial_hess(force_env, hess_mat)
1059
1060 TYPE(force_env_type), POINTER :: force_env
1061 TYPE(cp_fm_type), INTENT(IN) :: hess_mat
1062
1063 INTEGER :: i, iat_col, iat_row, iglobal, iind, j, &
1064 jat_row, jglobal, jind, k, natom, &
1065 ncol_local, nrow_local, z
1066 INTEGER, ALLOCATABLE, DIMENSION(:) :: at_row
1067 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
1068 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: d_ij, rho_ij
1069 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :) :: r_ij
1070 REAL(kind=dp), DIMENSION(3, 3) :: alpha, r0
1071 REAL(kind=dp), CONTIGUOUS, DIMENSION(:, :), POINTER :: fixed, local_data
1072 TYPE(cell_type), POINTER :: cell
1073 TYPE(cp_subsys_type), POINTER :: subsys
1074 TYPE(particle_list_type), POINTER :: particles
1075
1076 CALL force_env_get(force_env=force_env, subsys=subsys, cell=cell)
1077 CALL cp_subsys_get(subsys, &
1078 particles=particles)
1079
1080 alpha(1, :) = (/1._dp, 0.3949_dp, 0.3949_dp/)
1081 alpha(2, :) = (/0.3494_dp, 0.2800_dp, 0.2800_dp/)
1082 alpha(3, :) = (/0.3494_dp, 0.2800_dp, 0.1800_dp/)
1083
1084 r0(1, :) = (/1.35_dp, 2.10_dp, 2.53_dp/)
1085 r0(2, :) = (/2.10_dp, 2.87_dp, 3.40_dp/)
1086 r0(3, :) = (/2.53_dp, 3.40_dp, 3.40_dp/)
1087
1088 CALL cp_fm_get_info(hess_mat, row_indices=row_indices, col_indices=col_indices, &
1089 local_data=local_data, nrow_local=nrow_local, ncol_local=ncol_local)
1090 natom = particles%n_els
1091 ALLOCATE (at_row(natom))
1092 ALLOCATE (rho_ij(natom, natom))
1093 ALLOCATE (d_ij(natom, natom))
1094 ALLOCATE (r_ij(natom, natom, 3))
1095 ALLOCATE (fixed(3, natom))
1096 fixed = 1.0_dp
1097 CALL fix_atom_control(force_env, fixed)
1098 DO i = 1, 3
1099 CALL hess_mat%matrix_struct%para_env%min(fixed(i, :))
1100 END DO
1101 rho_ij = 0
1102 !XXXX insert proper rows !XXX
1103 at_row = 3
1104 DO i = 1, natom
1105 CALL get_atomic_kind(atomic_kind=particles%els(i)%atomic_kind, z=z)
1106 IF (z .LE. 10) at_row(i) = 2
1107 IF (z .LE. 2) at_row(i) = 1
1108 END DO
1109 DO i = 2, natom
1110 iat_row = at_row(i)
1111 DO j = 1, i - 1
1112 jat_row = at_row(j)
1113 !pbc for a distance vector
1114 r_ij(j, i, :) = pbc(particles%els(i)%r, particles%els(j)%r, cell)
1115 r_ij(i, j, :) = -r_ij(j, i, :)
1116 d_ij(j, i) = sqrt(dot_product(r_ij(j, i, :), r_ij(j, i, :)))
1117 d_ij(i, j) = d_ij(j, i)
1118 rho_ij(j, i) = exp(alpha(jat_row, iat_row)*(r0(jat_row, iat_row)**2 - d_ij(j, i)**2))
1119 rho_ij(i, j) = rho_ij(j, i)
1120 END DO
1121 END DO
1122 DO i = 1, ncol_local
1123 iglobal = col_indices(i)
1124 iind = mod(iglobal - 1, 3) + 1
1125 iat_col = (iglobal + 2)/3
1126 IF (iat_col .GT. natom) cycle
1127 DO j = 1, nrow_local
1128 jglobal = row_indices(j)
1129 jind = mod(jglobal - 1, 3) + 1
1130 iat_row = (jglobal + 2)/3
1131 IF (iat_row .GT. natom) cycle
1132 IF (iat_row .NE. iat_col) THEN
1133 IF (d_ij(iat_row, iat_col) .LT. 6.0_dp) &
1134 local_data(j, i) = local_data(j, i) + &
1135 angle_second_deriv(r_ij, d_ij, rho_ij, iind, jind, iat_col, iat_row, natom)
1136 ELSE
1137 local_data(j, i) = local_data(j, i) + &
1138 angle_second_deriv(r_ij, d_ij, rho_ij, iind, jind, iat_col, iat_row, natom)
1139 END IF
1140 IF (iat_col .NE. iat_row) THEN
1141 IF (d_ij(iat_row, iat_col) .LT. 6.0_dp) &
1142 local_data(j, i) = local_data(j, i) - &
1143 dist_second_deriv(r_ij(iat_col, iat_row, :), &
1144 iind, jind, d_ij(iat_row, iat_col), rho_ij(iat_row, iat_col))
1145 ELSE
1146 DO k = 1, natom
1147 IF (k == iat_col) cycle
1148 IF (d_ij(iat_row, k) .LT. 6.0_dp) &
1149 local_data(j, i) = local_data(j, i) + &
1150 dist_second_deriv(r_ij(iat_col, k, :), &
1151 iind, jind, d_ij(iat_row, k), rho_ij(iat_row, k))
1152 END DO
1153 END IF
1154 IF (fixed(jind, iat_row) .LT. 0.5_dp .OR. fixed(iind, iat_col) .LT. 0.5_dp) THEN
1155 local_data(j, i) = 0.0_dp
1156 IF (jind == iind .AND. iat_row == iat_col) local_data(j, i) = 1.0_dp
1157 END IF
1158 END DO
1159 END DO
1160 DEALLOCATE (fixed)
1161 DEALLOCATE (rho_ij)
1162 DEALLOCATE (d_ij)
1163 DEALLOCATE (r_ij)
1164 DEALLOCATE (at_row)
1165
1166 END SUBROUTINE construct_initial_hess
1167
1168! **************************************************************************************************
1169!> \brief ...
1170!> \param r1 ...
1171!> \param i ...
1172!> \param j ...
1173!> \param d ...
1174!> \param rho ...
1175!> \return ...
1176! **************************************************************************************************
1177 FUNCTION dist_second_deriv(r1, i, j, d, rho) RESULT(deriv)
1178 REAL(kind=dp), DIMENSION(3) :: r1
1179 INTEGER :: i, j
1180 REAL(kind=dp) :: d, rho, deriv
1181
1182 deriv = 0.45_dp*rho*(r1(i)*r1(j))/d**2
1183 END FUNCTION
1184
1185! **************************************************************************************************
1186!> \brief ...
1187!> \param r_ij ...
1188!> \param d_ij ...
1189!> \param rho_ij ...
1190!> \param idir ...
1191!> \param jdir ...
1192!> \param iat_der ...
1193!> \param jat_der ...
1194!> \param natom ...
1195!> \return ...
1196! **************************************************************************************************
1197 FUNCTION angle_second_deriv(r_ij, d_ij, rho_ij, idir, jdir, iat_der, jat_der, natom) RESULT(deriv)
1198 REAL(kind=dp), DIMENSION(:, :, :) :: r_ij
1199 REAL(kind=dp), DIMENSION(:, :) :: d_ij, rho_ij
1200 INTEGER :: idir, jdir, iat_der, jat_der, natom
1201 REAL(kind=dp) :: deriv
1202
1203 INTEGER :: i, iat, idr, j, jat, jdr
1204 REAL(kind=dp) :: d12, d23, d31, d_mat(3, 2), denom1, &
1205 denom2, denom3, ka1, ka2, ka3, rho12, &
1206 rho23, rho31, rsst1, rsst2, rsst3
1207 REAL(kind=dp), DIMENSION(3) :: r12, r23, r31
1208
1209 deriv = 0._dp
1210 IF (iat_der == jat_der) THEN
1211 DO i = 1, natom - 1
1212 IF (rho_ij(iat_der, i) .LT. 0.00001) cycle
1213 DO j = i + 1, natom
1214 IF (rho_ij(iat_der, j) .LT. 0.00001) cycle
1215 IF (i == iat_der .OR. j == iat_der) cycle
1216 IF (iat_der .LT. i .OR. iat_der .GT. j) THEN
1217 r12 = r_ij(iat_der, i, :); r23 = r_ij(i, j, :); r31 = r_ij(j, iat_der, :)
1218 d12 = d_ij(iat_der, i); d23 = d_ij(i, j); d31 = d_ij(j, iat_der)
1219 rho12 = rho_ij(iat_der, i); rho23 = rho_ij(i, j); rho31 = rho_ij(j, iat_der)
1220 ELSE
1221 r12 = r_ij(iat_der, j, :); r23 = r_ij(j, i, :); r31 = r_ij(i, iat_der, :)
1222 d12 = d_ij(iat_der, j); d23 = d_ij(j, i); d31 = d_ij(i, iat_der)
1223 rho12 = rho_ij(iat_der, j); rho23 = rho_ij(j, i); rho31 = rho_ij(i, iat_der)
1224 END IF
1225 ka1 = 0.15_dp*rho12*rho23; ka2 = 0.15_dp*rho23*rho31; ka3 = 0.15_dp*rho31*rho12
1226 rsst1 = dot_product(r12, r23); rsst2 = dot_product(r23, r31); rsst3 = dot_product(r31, r12)
1227 denom1 = 1.0_dp - rsst1**2/(d12**2*d23**2); denom2 = 1.0_dp - rsst2**2/(d23**2*d31**2)
1228 denom3 = 1.0_dp - rsst3**2/(d31**2*d12**2)
1229 denom1 = sign(1.0_dp, denom1)*max(abs(denom1), 0.01_dp)
1230 denom2 = sign(1.0_dp, denom2)*max(abs(denom2), 0.01_dp)
1231 denom3 = sign(1.0_dp, denom3)*max(abs(denom3), 0.01_dp)
1232 d_mat(1, 1) = r23(idir)/(d12*d23) - rsst1*r12(idir)/(d12**3*d23)
1233 d_mat(1, 2) = r23(jdir)/(d12*d23) - rsst1*r12(jdir)/(d12**3*d23)
1234 d_mat(2, 1) = -r23(idir)/(d23*d31) + rsst2*r31(idir)/(d23*d31**3)
1235 d_mat(2, 2) = -r23(jdir)/(d23*d31) + rsst2*r31(jdir)/(d23*d31**3)
1236 d_mat(3, 1) = (r31(idir) - r12(idir))/(d31*d12) + rsst3*r31(idir)/(d31**3*d12) - &
1237 rsst3*r12(idir)/(d31*d12**3)
1238 d_mat(3, 2) = (r31(jdir) - r12(jdir))/(d31*d12) + rsst3*r31(jdir)/(d31**3*d12) - &
1239 rsst3*r12(jdir)/(d31*d12**3)
1240 IF (abs(denom1) .LE. 0.011_dp) d_mat(1, 1) = 0.0_dp
1241 IF (abs(denom2) .LE. 0.011_dp) d_mat(2, 1) = 0.0_dp
1242 IF (abs(denom3) .LE. 0.011_dp) d_mat(3, 1) = 0.0_dp
1243 deriv = deriv + ka1*d_mat(1, 1)*d_mat(1, 2)/denom1 + &
1244 ka2*d_mat(2, 1)*d_mat(2, 2)/denom2 + &
1245 ka3*d_mat(3, 1)*d_mat(3, 2)/denom3
1246
1247 END DO
1248 END DO
1249 ELSE
1250 DO i = 1, natom
1251 IF (i == iat_der .OR. i == jat_der) cycle
1252 IF (jat_der .LT. iat_der) THEN
1253 iat = jat_der; jat = iat_der; idr = jdir; jdr = idir
1254 ELSE
1255 iat = iat_der; jat = jat_der; idr = idir; jdr = jdir
1256 END IF
1257 IF (jat .LT. i .OR. iat .GT. i) THEN
1258 r12 = r_ij(iat, jat, :); r23 = r_ij(jat, i, :); r31 = r_ij(i, iat, :)
1259 d12 = d_ij(iat, jat); d23 = d_ij(jat, i); d31 = d_ij(i, iat)
1260 rho12 = rho_ij(iat, jat); rho23 = rho_ij(jat, i); rho31 = rho_ij(i, iat)
1261 ELSE
1262 r12 = r_ij(iat, i, :); r23 = r_ij(i, jat, :); r31 = r_ij(jat, iat, :)
1263 d12 = d_ij(iat, i); d23 = d_ij(i, jat); d31 = d_ij(jat, iat)
1264 rho12 = rho_ij(iat, i); rho23 = rho_ij(i, jat); rho31 = rho_ij(jat, iat)
1265 END IF
1266 ka1 = 0.15_dp*rho12*rho23; ka2 = 0.15_dp*rho23*rho31; ka3 = 0.15_dp*rho31*rho12
1267 rsst1 = dot_product(r12, r23); rsst2 = dot_product(r23, r31); rsst3 = dot_product(r31, r12)
1268 denom1 = 1.0_dp - rsst1**2/(d12**2*d23**2); denom2 = 1.0_dp - rsst2**2/(d23**2*d31**2)
1269 denom3 = 1.0_dp - rsst3**2/(d31**2*d12**2)
1270 denom1 = sign(1.0_dp, denom1)*max(abs(denom1), 0.01_dp)
1271 denom2 = sign(1.0_dp, denom2)*max(abs(denom2), 0.01_dp)
1272 denom3 = sign(1.0_dp, denom3)*max(abs(denom3), 0.01_dp)
1273 d_mat(1, 1) = r23(idr)/(d12*d23) - rsst1*r12(idr)/(d12**3*d23)
1274 d_mat(2, 1) = -r23(idr)/(d23*d31) + rsst2*r31(idr)/(d23*d31**3)
1275 d_mat(3, 1) = (r31(idr) - r12(idr))/(d31*d12) + rsst3*r31(idr)/(d31**3*d12) - &
1276 rsst3*r12(idr)/(d31*d12**3)
1277 IF (jat .LT. i .OR. iat .GT. i) THEN
1278 d_mat(1, 2) = (r12(jdr) - r23(jdr))/(d12*d23) + rsst1*r12(jdr)/(d12**3*d23) - &
1279 rsst1*r23(jdr)/(d12*d23**3)
1280 d_mat(2, 2) = r31(jdr)/(d23*d31) - rsst2*r23(jdr)/(d23**3*d31)
1281 d_mat(3, 2) = -r31(jdr)/(d31*d12) + rsst3*r12(jdr)/(d31*d12**3)
1282 ELSE
1283 d_mat(1, 2) = -r12(jdr)/(d12*d23) + rsst1*r23(jdr)/(d12*d23**3)
1284 d_mat(2, 2) = (r23(jdr) - r31(jdr))/(d23*d31) + rsst2*r23(jdr)/(d23**3*d31) - &
1285 rsst2*r31(jdr)/(d23*d31**3)
1286 d_mat(3, 2) = r12(jdr)/(d31*d12) - rsst3*r31(jdr)/(d31**3*d12)
1287 END IF
1288 IF (abs(denom1) .LE. 0.011_dp) d_mat(1, 1) = 0.0_dp
1289 IF (abs(denom2) .LE. 0.011_dp) d_mat(2, 1) = 0.0_dp
1290 IF (abs(denom3) .LE. 0.011_dp) d_mat(3, 1) = 0.0_dp
1291
1292 deriv = deriv + ka1*d_mat(1, 1)*d_mat(1, 2)/denom1 + &
1293 ka2*d_mat(2, 1)*d_mat(2, 2)/denom2 + &
1294 ka3*d_mat(3, 1)*d_mat(3, 2)/denom3
1295 END DO
1296 END IF
1297 deriv = 0.25_dp*deriv
1298
1299 END FUNCTION angle_second_deriv
1300
1301END MODULE bfgs_optimizer
cp_eval_at
subroutine cp_eval_at(gopt_env, x, f, gradient, master, final_evaluation, para_env)
evaluete the potential energy and its gradients using an array with same dimension as the particle_se...
Definition gopt_f77_methods.F:25
cell_types::pbc
Definition cell_types.F:92
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:87
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:82
cp_log_handling::cp_to_string
Definition cp_log_handling.F:90
parallel_gemm_api::parallel_gemm
Definition parallel_gemm_api.F:38
atomic_kind_list_types
represent a simple array based list of the given type
Definition atomic_kind_list_types.F:15
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind_set
subroutine, public get_atomic_kind_set(atomic_kind_set, atom_of_kind, kind_of, natom_of_kind, maxatom, natom, nshell, fist_potential_present, shell_present, shell_adiabatic, shell_check_distance, damping_present)
Get attributes of an atomic kind set.
Definition atomic_kind_types.F:223
atomic_kind_types::get_atomic_kind
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.
Definition atomic_kind_types.F:138
bfgs_optimizer
Routines for Geometry optimization using BFGS algorithm.
Definition bfgs_optimizer.F:14
bfgs_optimizer::geoopt_bfgs
recursive subroutine, public geoopt_bfgs(force_env, gopt_param, globenv, geo_section, gopt_env, x0)
Main driver for BFGS geometry optimizations.
Definition bfgs_optimizer.F:138
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
constraint_fxd
Definition constraint_fxd.F:12
constraint_fxd::fix_atom_control
subroutine, public fix_atom_control(force_env, w)
allows for fix atom constraints
Definition constraint_fxd.F:65
cp_blacs_env
methods related to the blacs parallel environment
Definition cp_blacs_env.F:15
cp_blacs_env::cp_blacs_env_release
subroutine, public cp_blacs_env_release(blacs_env)
releases the given blacs_env
Definition cp_blacs_env.F:274
cp_blacs_env::cp_blacs_env_create
subroutine, public cp_blacs_env_create(blacs_env, para_env, blacs_grid_layout, blacs_repeatable, row_major, grid_2d)
allocates and initializes a type that represent a blacs context
Definition cp_blacs_env.F:123
cp_external_control
Routines to handle the external control of CP2K.
Definition cp_external_control.F:15
cp_external_control::external_control
subroutine, public external_control(should_stop, flag, globenv, target_time, start_time, force_check)
External manipulations during a run : when the <PROJECT_NAME>.EXIT_$runtype command is sent the progr...
Definition cp_external_control.F:90
cp_files
Utility routines to open and close files. Tracking of preconnections.
Definition cp_files.F:16
cp_files::open_file
subroutine, public open_file(file_name, file_status, file_form, file_action, file_position, file_pad, unit_number, debug, skip_get_unit_number, file_access)
Opens the requested file using a free unit number.
Definition cp_files.F:308
cp_files::close_file
subroutine, public close_file(unit_number, file_status, keep_preconnection)
Close an open file given by its logical unit number. Optionally, keep the file and unit preconnected.
Definition cp_files.F:119
cp_fm_basic_linalg
Basic linear algebra operations for full matrices.
Definition cp_fm_basic_linalg.F:14
cp_fm_basic_linalg::cp_fm_column_scale
subroutine, public cp_fm_column_scale(matrixa, scaling)
scales column i of matrix a with scaling(i)
Definition cp_fm_basic_linalg.F:1885
cp_fm_basic_linalg::cp_fm_transpose
subroutine, public cp_fm_transpose(matrix, matrixt)
transposes a matrix matrixt = matrix ^ T
Definition cp_fm_basic_linalg.F:1695
cp_fm_diag
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition cp_fm_diag.F:17
cp_fm_diag::choose_eigv_solver
subroutine, public choose_eigv_solver(matrix, eigenvectors, eigenvalues, info)
Choose the Eigensolver depending on which library is available ELPA seems to be unstable for small sy...
Definition cp_fm_diag.F:238
cp_fm_struct
represent the structure of a full matrix
Definition cp_fm_struct.F:14
cp_fm_struct::cp_fm_struct_create
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
Definition cp_fm_struct.F:132
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:351
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_info
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
Definition cp_fm_types.F:1009
cp_fm_types::cp_fm_write_unformatted
subroutine, public cp_fm_write_unformatted(fm, unit)
...
Definition cp_fm_types.F:2146
cp_fm_types::cp_fm_read_unformatted
subroutine, public cp_fm_read_unformatted(fm, unit)
...
Definition cp_fm_types.F:2377
cp_fm_types::cp_fm_set_all
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
Definition cp_fm_types.F:528
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
Definition cp_fm_types.F:164
cp_log_handling
various routines to log and control the output. The idea is that decisions about where to log should ...
Definition cp_log_handling.F:41
cp_log_handling::cp_get_default_logger
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
Definition cp_log_handling.F:234
cp_output_handling
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
Definition cp_output_handling.F:25
cp_output_handling::cp_print_key_unit_nr
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
Definition cp_output_handling.F:853
cp_output_handling::cp_print_key_finished_output
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
Definition cp_output_handling.F:1063
cp_output_handling::cp_p_file
integer, parameter, public cp_p_file
Definition cp_output_handling.F:103
cp_output_handling::cp_iterate
subroutine, public cp_iterate(iteration_info, last, iter_nr, increment, iter_nr_out)
adds one to the actual iteration
Definition cp_output_handling.F:574
cp_output_handling::cp_print_key_should_output
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
Definition cp_output_handling.F:345
cp_subsys_types
types that represent a subsys, i.e. a part of the system
Definition cp_subsys_types.F:16
cp_subsys_types::cp_subsys_get
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
Definition cp_subsys_types.F:345
force_env_types
Interface for the force calculations.
Definition force_env_types.F:18
force_env_types::force_env_get
recursive subroutine, public force_env_get(force_env, in_use, fist_env, qs_env, meta_env, fp_env, subsys, para_env, potential_energy, additional_potential, kinetic_energy, harmonic_shell, kinetic_shell, cell, sub_force_env, qmmm_env, qmmmx_env, eip_env, pwdft_env, globenv, input, force_env_section, method_name_id, root_section, mixed_env, nnp_env, embed_env, ipi_env)
returns various attributes about the force environment
Definition force_env_types.F:342
global_types
Define type storing the global information of a run. Keep the amount of stored data small....
Definition global_types.F:21
gopt_f_methods
contains a functional that calculates the energy and its derivatives for the geometry optimizer
Definition gopt_f_methods.F:14
gopt_f_methods::print_geo_opt_header
subroutine, public print_geo_opt_header(gopt_env, output_unit, label)
...
Definition gopt_f_methods.F:986
gopt_f_methods::gopt_f_io_init
subroutine, public gopt_f_io_init(gopt_env, output_unit, opt_energy, wildcard, its, used_time)
Handles the Output during an optimization run.
Definition gopt_f_methods.F:220
gopt_f_methods::gopt_f_io_finalize
recursive subroutine, public gopt_f_io_finalize(gopt_env, force_env, x0, conv, its, root_section, para_env, master, output_unit)
Handles the Output at the end of an optimization run.
Definition gopt_f_methods.F:475
gopt_f_methods::gopt_f_io
subroutine, public gopt_f_io(gopt_env, force_env, root_section, its, opt_energy, output_unit, eold, emin, wildcard, gopt_param, ndf, dx, xi, conv, pred, rat, step, rad, used_time)
Handles the Output during an optimization run.
Definition gopt_f_methods.F:321
gopt_f_methods::print_geo_opt_nc
subroutine, public print_geo_opt_nc(gopt_env, output_unit)
...
Definition gopt_f_methods.F:1025
gopt_f_methods::gopt_f_ii
subroutine, public gopt_f_ii(its, output_unit)
Prints iteration step of the optimization procedure on screen.
Definition gopt_f_methods.F:197
gopt_f_types
contains a functional that calculates the energy and its derivatives for the geometry optimizer
Definition gopt_f_types.F:15
gopt_param_types
contains typo and related routines to handle parameters controlling the GEO_OPT module
Definition gopt_param_types.F:14
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:17
input_constants::default_cell_method_id
integer, parameter, public default_cell_method_id
Definition input_constants.F:537
input_constants::default_ts_method_id
integer, parameter, public default_ts_method_id
Definition input_constants.F:537
input_section_types
objects that represent the structure of input sections and the data contained in an input section
Definition input_section_types.F:15
input_section_types::section_vals_val_set
subroutine, public section_vals_val_set(section_vals, keyword_name, i_rep_section, i_rep_val, val, l_val, i_val, r_val, c_val, l_vals_ptr, i_vals_ptr, r_vals_ptr, c_vals_ptr)
sets the requested value
Definition input_section_types.F:1223
input_section_types::section_vals_get_subs_vals
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
Definition input_section_types.F:731
input_section_types::section_vals_val_get
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
Definition input_section_types.F:1047
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_path_length
integer, parameter, public default_path_length
Definition kinds.F:58
machine
Machine interface based on Fortran 2003 and POSIX.
Definition machine.F:17
machine::m_flush
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition machine.F:130
machine::m_walltime
real(kind=dp) function, public m_walltime()
returns time from a real-time clock, protected against rolling early/easily
Definition machine.F:147
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:23
parallel_gemm_api
basic linear algebra operations for full matrixes
Definition parallel_gemm_api.F:14
particle_list_types
represent a simple array based list of the given type
Definition particle_list_types.F:15
space_groups_types
Space Group Symmetry Type Module (version 1.0, Ferbruary 12, 2021)
Definition space_groups_types.F:14
space_groups
Space Group Symmetry Module (version 1.0, January 16, 2020)
Definition space_groups.F:14
space_groups::print_spgr
subroutine, public print_spgr(spgr)
routine prints Space Group Information.
Definition space_groups.F:830
space_groups::spgr_apply_rotations_coord
subroutine, public spgr_apply_rotations_coord(spgr, coord)
routine applies the rotation matrices to the coordinates.
Definition space_groups.F:618
space_groups::identify_space_group
subroutine, public identify_space_group(subsys, geo_section, gopt_env, iunit)
routine indentifies the space group and finds rotation matrices.
Definition space_groups.F:286
space_groups::spgr_apply_rotations_force
subroutine, public spgr_apply_rotations_force(spgr, force)
routine applies the rotation matrices to the forces.
Definition space_groups.F:679
atomic_kind_list_types::atomic_kind_list_type
represent a list of objects
Definition atomic_kind_list_types.F:45
cell_types::cell_type
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
cp_blacs_env::cp_blacs_env_type
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
Definition cp_blacs_env.F:53
cp_fm_struct::cp_fm_struct_type
keeps the information about the structure of a full matrix
Definition cp_fm_struct.F:82
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:113
cp_log_handling::cp_logger_type
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Definition cp_log_handling.F:140
cp_subsys_types::cp_subsys_type
represents a system: atoms, molecules, their pos,vel,...
Definition cp_subsys_types.F:89
force_env_types::force_env_type
wrapper to abstract the force evaluation of the various methods
Definition force_env_types.F:145
global_types::global_environment_type
contains the initially parsed file and the initial parallel environment
Definition global_types.F:66
gopt_f_types::gopt_f_type
calculates the potential energy of a system, and its derivatives
Definition gopt_f_types.F:60
gopt_param_types::gopt_param_type
Definition gopt_param_types.F:56
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:763
particle_list_types::particle_list_type
represent a list of objects
Definition particle_list_types.F:46
space_groups_types::spgr_type
Definition space_groups_types.F:27