(git:3919fb7)
Loading...
Searching...
No Matches
cp_lbfgs_optimizer_gopt.F
Go to the documentation of this file.
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 routines that optimize a functional using the limited memory bfgs
10!> quasi-newton method.
11!> The process set up so that a master runs the real optimizer and the
12!> others help then to calculate the objective function.
13!> The arguments for the objective function are physically present in
14!> every processor (nedeed in the actual implementation of pao).
15!> In the future tha arguments themselves could be distributed.
16!> \par History
17!> 09.2003 globenv->para_env, retain/release, better parallel behaviour
18!> 01.2020 Space Group Symmetry introduced by Pierre-André Cazade [pcazade]
19!> \author Fawzi Mohamed
20!> @version 2.2002
21! **************************************************************************************************
22MODULE cp_lbfgs_optimizer_gopt
23 USE cp_lbfgs, ONLY: setulb
24 USE cp_log_handling, ONLY: cp_get_default_logger,&
25 cp_logger_type,&
26 cp_to_string
27 USE cp_output_handling, ONLY: cp_print_key_finished_output,&
28 cp_print_key_unit_nr
29 USE cp_subsys_types, ONLY: cp_subsys_type
30 USE force_env_types, ONLY: force_env_get,&
31 force_env_type
32 USE gopt_f_methods, ONLY: cp_eval_at,&
33 gopt_f_io
34 USE gopt_f_types, ONLY: gopt_f_release,&
35 gopt_f_retain,&
36 gopt_f_type
37 USE gopt_param_types, ONLY: gopt_param_type
38 USE input_section_types, ONLY: section_vals_type
39 USE kinds, ONLY: dp
40 USE machine, ONLY: m_walltime
41 USE message_passing, ONLY: mp_para_env_release,&
42 mp_para_env_type
43 USE space_groups, ONLY: spgr_apply_rotations_coord,&
44 spgr_apply_rotations_force
45 USE space_groups_types, ONLY: spgr_type
46#include "../base/base_uses.f90"
47
48 IMPLICIT NONE
49 PRIVATE
50
51 LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .true.
52 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cp_lbfgs_optimizer_gopt'
53
54 ! types
55 PUBLIC :: cp_lbfgs_opt_gopt_type
56
57 ! core methods
58
59 ! special methos
60
61 ! underlying functions
62 PUBLIC :: cp_opt_gopt_create, cp_opt_gopt_release, &
63 cp_opt_gopt_next, &
64 cp_opt_gopt_stop
65
66! **************************************************************************************************
67!> \brief info for the optimizer (see the description of this module)
68!> \param task the actual task of the optimizer (in the master it is up to
69!> date, in case of error also the minions one get updated.
70!> \param csave internal character string used by the lbfgs optimizer,
71!> meaningful only in the master
72!> \param lsave logical array used by the lbfgs optimizer, updated only
73!> in the master
74!> On exit with task = 'NEW_X', the following information is
75!> available:
76!> lsave(1) = .true. the initial x did not satisfy the bounds;
77!> lsave(2) = .true. the problem contains bounds;
78!> lsave(3) = .true. each variable has upper and lower bounds.
79!> \param ref_count reference count (see doc/ReferenceCounting.html)
80!> \param m the dimension of the subspace used to approximate the second
81!> derivative
82!> \param print_every every how many iterations output should be written.
83!> if 0 only at end, if print_every<0 never
84!> \param master the pid of the master processor
85!> \param max_f_per_iter the maximum number of function evaluations per
86!> iteration
87!> \param status 0: just initialized, 1: f g calculation,
88!> 2: begin new iteration, 3: ended iteration,
89!> 4: normal (converged) exit, 5: abnormal (error) exit,
90!> 6: daellocated
91!> \param n_iter the actual iteration number
92!> \param kind_of_bound an array with 0 (no bound), 1 (lower bound),
93!> 2 (both bounds), 3 (upper bound), to describe the bounds
94!> of every variable
95!> \param i_work_array an integer workarray of dimension 3*n, present only
96!> in the master
97!> \param isave is an INTEGER working array of dimension 44.
98!> On exit with task = 'NEW_X', it contains information that
99!> the user may want to access:
100!> \param isave (30) = the current iteration number;
101!> \param isave (34) = the total number of function and gradient
102!> evaluations;
103!> \param isave (36) = the number of function value or gradient
104!> evaluations in the current iteration;
105!> \param isave (38) = the number of free variables in the current
106!> iteration;
107!> \param isave (39) = the number of active constraints at the current
108!> iteration;
109!> \param f the actual best value of the object function
110!> \param wanted_relative_f_delta the wanted relative error on f
111!> (to be multiplied by epsilon), 0.0 -> no check
112!> \param wanted_projected_gradient the wanted error on the projected
113!> gradient (hessian times the gradient), 0.0 -> no check
114!> \param last_f the value of f in the last iteration
115!> \param projected_gradient the value of the sup norm of the projected
116!> gradient
117!> \param x the actual evaluation point (best one if converged or stopped)
118!> \param lower_bound the lower bounds
119!> \param upper_bound the upper bounds
120!> \param gradient the actual gradient
121!> \param dsave info date for lbfgs (master only)
122!> \param work_array a work array for lbfgs (master only)
123!> \param para_env the parallel environment for this optimizer
124!> \param obj_funct the objective function to be optimized
125!> \par History
126!> none
127!> \author Fawzi Mohamed
128!> @version 2.2002
129! **************************************************************************************************
130 TYPE cp_lbfgs_opt_gopt_type
131 CHARACTER(len=60) :: task = ""
132 CHARACTER(len=60) :: csave = ""
133 LOGICAL :: lsave(4) = .false.
134 INTEGER :: m = 0, print_every = 0, master = 0, max_f_per_iter = 0, status = 0, n_iter = 0
135 INTEGER, DIMENSION(:), POINTER :: kind_of_bound => null(), i_work_array => null(), isave => null()
136 REAL(kind=dp) :: f = 0.0_dp, wanted_relative_f_delta = 0.0_dp, wanted_projected_gradient = 0.0_dp, &
137 last_f = 0.0_dp, projected_gradient = 0.0_dp, eold = 0.0_dp, emin = 0.0_dp, trust_radius = 0.0_dp
138 REAL(kind=dp), DIMENSION(:), POINTER :: x => null(), lower_bound => null(), upper_bound => null(), &
139 gradient => null(), dsave => null(), work_array => null()
140 TYPE(mp_para_env_type), POINTER :: para_env => null()
141 TYPE(gopt_f_type), POINTER :: obj_funct => null()
142 END TYPE cp_lbfgs_opt_gopt_type
143
144CONTAINS
145
146! **************************************************************************************************
147!> \brief initializes the optimizer
148!> \param optimizer ...
149!> \param para_env ...
150!> \param obj_funct ...
151!> \param x0 ...
152!> \param m ...
153!> \param print_every ...
154!> \param wanted_relative_f_delta ...
155!> \param wanted_projected_gradient ...
156!> \param lower_bound ...
157!> \param upper_bound ...
158!> \param kind_of_bound ...
159!> \param master ...
160!> \param max_f_per_iter ...
161!> \param trust_radius ...
162!> \par History
163!> 02.2002 created [fawzi]
164!> 09.2003 refactored (retain/release,para_env) [fawzi]
165!> \author Fawzi Mohamed
166!> \note
167!> redirects the lbfgs output the the default unit
168! **************************************************************************************************
169 SUBROUTINE cp_opt_gopt_create(optimizer, para_env, obj_funct, x0, m, print_every, &
170 wanted_relative_f_delta, wanted_projected_gradient, lower_bound, upper_bound, &
171 kind_of_bound, master, max_f_per_iter, trust_radius)
172 TYPE(cp_lbfgs_opt_gopt_type), INTENT(OUT) :: optimizer
173 TYPE(mp_para_env_type), POINTER :: para_env
174 TYPE(gopt_f_type), POINTER :: obj_funct
175 REAL(kind=dp), DIMENSION(:), INTENT(in) :: x0
176 INTEGER, INTENT(in), OPTIONAL :: m, print_every
177 REAL(kind=dp), INTENT(in), OPTIONAL :: wanted_relative_f_delta, &
178 wanted_projected_gradient
179 REAL(kind=dp), DIMENSION(SIZE(x0)), INTENT(in), &
180 OPTIONAL :: lower_bound, upper_bound
181 INTEGER, DIMENSION(SIZE(x0)), INTENT(in), OPTIONAL :: kind_of_bound
182 INTEGER, INTENT(in), OPTIONAL :: master, max_f_per_iter
183 REAL(kind=dp), INTENT(in), OPTIONAL :: trust_radius
184
185 CHARACTER(len=*), PARAMETER :: routinen = 'cp_opt_gopt_create'
186
187 INTEGER :: handle, lenwa, n
188
189 CALL timeset(routinen, handle)
190
191 NULLIFY (optimizer%kind_of_bound, &
192 optimizer%i_work_array, &
193 optimizer%isave, &
194 optimizer%x, &
195 optimizer%lower_bound, &
196 optimizer%upper_bound, &
197 optimizer%gradient, &
198 optimizer%dsave, &
199 optimizer%work_array, &
200 optimizer%para_env, &
201 optimizer%obj_funct)
202 n = SIZE(x0)
203 optimizer%m = 4
204 IF (PRESENT(m)) optimizer%m = m
205 optimizer%master = para_env%source
206 optimizer%para_env => para_env
207 CALL para_env%retain()
208 optimizer%obj_funct => obj_funct
209 CALL gopt_f_retain(obj_funct)
210 optimizer%max_f_per_iter = 20
211 IF (PRESENT(max_f_per_iter)) optimizer%max_f_per_iter = max_f_per_iter
212 optimizer%print_every = -1
213 optimizer%n_iter = 0
214 optimizer%f = -1.0_dp
215 optimizer%last_f = -1.0_dp
216 optimizer%projected_gradient = -1.0_dp
217 IF (PRESENT(print_every)) optimizer%print_every = print_every
218 IF (PRESENT(master)) optimizer%master = master
219 IF (optimizer%master == optimizer%para_env%mepos) THEN
220 !MK This has to be adapted for a new L-BFGS version possibly
221 lenwa = 2*optimizer%m*n + 5*n + 11*optimizer%m*optimizer%m + 8*optimizer%m
222 ALLOCATE (optimizer%kind_of_bound(n), optimizer%i_work_array(3*n), &
223 optimizer%isave(44))
224 ALLOCATE (optimizer%x(n), optimizer%lower_bound(n), &
225 optimizer%upper_bound(n), optimizer%gradient(n), &
226 optimizer%dsave(29), optimizer%work_array(lenwa))
227 optimizer%x = x0
228 optimizer%task = 'START'
229 optimizer%i_work_array = 0
230 optimizer%isave = 0
231 optimizer%lower_bound = 0.0_dp
232 optimizer%upper_bound = 0.0_dp
233 optimizer%gradient = 0.0_dp
234 optimizer%dsave = 0.0_dp
235 optimizer%work_array = 0.0_dp
236 IF (PRESENT(wanted_relative_f_delta)) THEN
237 optimizer%wanted_relative_f_delta = wanted_relative_f_delta
238 END IF
239 IF (PRESENT(wanted_projected_gradient)) THEN
240 optimizer%wanted_projected_gradient = wanted_projected_gradient
241 END IF
242 optimizer%kind_of_bound = 0
243 IF (PRESENT(kind_of_bound)) optimizer%kind_of_bound = kind_of_bound
244 IF (PRESENT(lower_bound)) optimizer%lower_bound = lower_bound
245 IF (PRESENT(upper_bound)) optimizer%upper_bound = upper_bound
246 IF (PRESENT(trust_radius)) optimizer%trust_radius = trust_radius
247
248 CALL setulb(SIZE(optimizer%x), optimizer%m, optimizer%x, &
249 optimizer%lower_bound, optimizer%upper_bound, &
250 optimizer%kind_of_bound, optimizer%f, optimizer%gradient, &
251 optimizer%wanted_relative_f_delta, &
252 optimizer%wanted_projected_gradient, optimizer%work_array, &
253 optimizer%i_work_array, optimizer%task, optimizer%print_every, &
254 optimizer%csave, optimizer%lsave, optimizer%isave, &
255 optimizer%dsave, optimizer%trust_radius)
256 ELSE
257 NULLIFY ( &
258 optimizer%kind_of_bound, optimizer%i_work_array, optimizer%isave, &
259 optimizer%lower_bound, optimizer%upper_bound, optimizer%gradient, &
260 optimizer%dsave, optimizer%work_array)
261 ALLOCATE (optimizer%x(n))
262 optimizer%x(:) = 0.0_dp
263 ALLOCATE (optimizer%gradient(n))
264 optimizer%gradient(:) = 0.0_dp
265 END IF
266 CALL optimizer%para_env%bcast(optimizer%x, optimizer%master)
267 optimizer%status = 0
268
269 CALL timestop(handle)
270
271 END SUBROUTINE cp_opt_gopt_create
272
273! **************************************************************************************************
274!> \brief releases the optimizer (see doc/ReferenceCounting.html)
275!> \param optimizer the object that should be released
276!> \par History
277!> 02.2002 created [fawzi]
278!> 09.2003 dealloc_ref->release [fawzi]
279!> \author Fawzi Mohamed
280! **************************************************************************************************
281 SUBROUTINE cp_opt_gopt_release(optimizer)
282 TYPE(cp_lbfgs_opt_gopt_type), INTENT(INOUT) :: optimizer
283
284 CHARACTER(len=*), PARAMETER :: routinen = 'cp_opt_gopt_release'
285
286 INTEGER :: handle
287
288 CALL timeset(routinen, handle)
289
290 IF (ASSOCIATED(optimizer%kind_of_bound)) THEN
291 DEALLOCATE (optimizer%kind_of_bound)
292 END IF
293 IF (ASSOCIATED(optimizer%i_work_array)) THEN
294 DEALLOCATE (optimizer%i_work_array)
295 END IF
296 IF (ASSOCIATED(optimizer%isave)) THEN
297 DEALLOCATE (optimizer%isave)
298 END IF
299 IF (ASSOCIATED(optimizer%x)) THEN
300 DEALLOCATE (optimizer%x)
301 END IF
302 IF (ASSOCIATED(optimizer%lower_bound)) THEN
303 DEALLOCATE (optimizer%lower_bound)
304 END IF
305 IF (ASSOCIATED(optimizer%upper_bound)) THEN
306 DEALLOCATE (optimizer%upper_bound)
307 END IF
308 IF (ASSOCIATED(optimizer%gradient)) THEN
309 DEALLOCATE (optimizer%gradient)
310 END IF
311 IF (ASSOCIATED(optimizer%dsave)) THEN
312 DEALLOCATE (optimizer%dsave)
313 END IF
314 IF (ASSOCIATED(optimizer%work_array)) THEN
315 DEALLOCATE (optimizer%work_array)
316 END IF
317 CALL mp_para_env_release(optimizer%para_env)
318 CALL gopt_f_release(optimizer%obj_funct)
319
320 CALL timestop(handle)
321 END SUBROUTINE cp_opt_gopt_release
322
323! **************************************************************************************************
324!> \brief takes different valuse from the optimizer
325!> \param optimizer ...
326!> \param para_env ...
327!> \param obj_funct ...
328!> \param m ...
329!> \param print_every ...
330!> \param wanted_relative_f_delta ...
331!> \param wanted_projected_gradient ...
332!> \param x ...
333!> \param lower_bound ...
334!> \param upper_bound ...
335!> \param kind_of_bound ...
336!> \param master ...
337!> \param actual_projected_gradient ...
338!> \param n_var ...
339!> \param n_iter ...
340!> \param status ...
341!> \param max_f_per_iter ...
342!> \param at_end ...
343!> \param is_master ...
344!> \param last_f ...
345!> \param f ...
346!> \par History
347!> none
348!> \author Fawzi Mohamed
349!> @version 2.2002
350! **************************************************************************************************
351 SUBROUTINE cp_opt_gopt_get(optimizer, para_env, &
352 obj_funct, m, print_every, &
353 wanted_relative_f_delta, wanted_projected_gradient, &
354 x, lower_bound, upper_bound, kind_of_bound, master, &
355 actual_projected_gradient, &
356 n_var, n_iter, status, max_f_per_iter, at_end, &
357 is_master, last_f, f)
358 TYPE(cp_lbfgs_opt_gopt_type), INTENT(IN) :: optimizer
359 TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env
360 TYPE(gopt_f_type), OPTIONAL, POINTER :: obj_funct
361 INTEGER, INTENT(out), OPTIONAL :: m, print_every
362 REAL(kind=dp), INTENT(out), OPTIONAL :: wanted_relative_f_delta, &
363 wanted_projected_gradient
364 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: x, lower_bound, upper_bound
365 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: kind_of_bound
366 INTEGER, INTENT(out), OPTIONAL :: master
367 REAL(kind=dp), INTENT(out), OPTIONAL :: actual_projected_gradient
368 INTEGER, INTENT(out), OPTIONAL :: n_var, n_iter, status, max_f_per_iter
369 LOGICAL, INTENT(out), OPTIONAL :: at_end, is_master
370 REAL(kind=dp), INTENT(out), OPTIONAL :: last_f, f
371
372 IF (PRESENT(is_master)) is_master = optimizer%master == optimizer%para_env%mepos
373 IF (PRESENT(master)) master = optimizer%master
374 IF (PRESENT(status)) status = optimizer%status
375 IF (PRESENT(para_env)) para_env => optimizer%para_env
376 IF (PRESENT(obj_funct)) obj_funct = optimizer%obj_funct
377 IF (PRESENT(m)) m = optimizer%m
378 IF (PRESENT(max_f_per_iter)) max_f_per_iter = optimizer%max_f_per_iter
379 IF (PRESENT(wanted_projected_gradient)) THEN
380 wanted_projected_gradient = optimizer%wanted_projected_gradient
381 END IF
382 IF (PRESENT(wanted_relative_f_delta)) THEN
383 wanted_relative_f_delta = optimizer%wanted_relative_f_delta
384 END IF
385 IF (PRESENT(print_every)) print_every = optimizer%print_every
386 IF (PRESENT(x)) x => optimizer%x
387 IF (PRESENT(n_var)) n_var = SIZE(x)
388 IF (PRESENT(lower_bound)) lower_bound => optimizer%lower_bound
389 IF (PRESENT(upper_bound)) upper_bound => optimizer%upper_bound
390 IF (PRESENT(kind_of_bound)) kind_of_bound => optimizer%kind_of_bound
391 IF (PRESENT(n_iter)) n_iter = optimizer%n_iter
392 IF (PRESENT(last_f)) last_f = optimizer%last_f
393 IF (PRESENT(f)) f = optimizer%f
394 IF (PRESENT(at_end)) at_end = optimizer%status > 3
395 IF (PRESENT(actual_projected_gradient)) THEN
396 actual_projected_gradient = optimizer%projected_gradient
397 END IF
398 IF (optimizer%master == optimizer%para_env%mepos) THEN
399 IF (optimizer%isave(30) > 1 .AND. (optimizer%task(1:5) == "NEW_X" .OR. &
400 optimizer%task(1:4) == "STOP" .AND. optimizer%task(7:9) == "CPU")) THEN
401 ! nr iterations >1 .and. dsave contains the wanted data
402 IF (PRESENT(last_f)) last_f = optimizer%dsave(2)
403 IF (PRESENT(actual_projected_gradient)) THEN
404 actual_projected_gradient = optimizer%dsave(13)
405 END IF
406 ELSE
407 cpassert(.NOT. PRESENT(last_f))
408 cpassert(.NOT. PRESENT(actual_projected_gradient))
409 END IF
410 ELSE IF (PRESENT(lower_bound) .OR. PRESENT(upper_bound) .OR. PRESENT(kind_of_bound)) THEN
411 cpwarn("asked undefined types")
412 END IF
413
414 END SUBROUTINE cp_opt_gopt_get
415
416! **************************************************************************************************
417!> \brief does one optimization step
418!> \param optimizer ...
419!> \param n_iter ...
420!> \param f ...
421!> \param last_f ...
422!> \param projected_gradient ...
423!> \param converged ...
424!> \param geo_section ...
425!> \param force_env ...
426!> \param gopt_param ...
427!> \param spgr ...
428!> \par History
429!> 01.2020 modified [pcazade]
430!> \author Fawzi Mohamed
431!> @version 2.2002
432!> \note
433!> use directly mainlb in place of setulb ??
434! **************************************************************************************************
435 SUBROUTINE cp_opt_gopt_step(optimizer, n_iter, f, last_f, &
436 projected_gradient, converged, geo_section, force_env, &
437 gopt_param, spgr)
438 TYPE(cp_lbfgs_opt_gopt_type), INTENT(INOUT) :: optimizer
439 INTEGER, INTENT(out), OPTIONAL :: n_iter
440 REAL(kind=dp), INTENT(out), OPTIONAL :: f, last_f, projected_gradient
441 LOGICAL, INTENT(out), OPTIONAL :: converged
442 TYPE(section_vals_type), POINTER :: geo_section
443 TYPE(force_env_type), POINTER :: force_env
444 TYPE(gopt_param_type), POINTER :: gopt_param
445 TYPE(spgr_type), OPTIONAL, POINTER :: spgr
446
447 CHARACTER(len=*), PARAMETER :: routinen = 'cp_opt_gopt_step'
448
449 CHARACTER(LEN=5) :: wildcard
450 INTEGER :: dataunit, handle, its
451 LOGICAL :: conv, is_master, justentred, &
452 keep_space_group
453 REAL(kind=dp) :: t_diff, t_now, t_old
454 REAL(kind=dp), DIMENSION(:), POINTER :: xold
455 TYPE(cp_logger_type), POINTER :: logger
456 TYPE(cp_subsys_type), POINTER :: subsys
457
458 NULLIFY (logger, xold)
459 logger => cp_get_default_logger()
460 CALL timeset(routinen, handle)
461 justentred = .true.
462 is_master = optimizer%master == optimizer%para_env%mepos
463 IF (PRESENT(converged)) converged = optimizer%status == 4
464 ALLOCATE (xold(SIZE(optimizer%x)))
465
466 ! collecting subsys
467 CALL force_env_get(force_env, subsys=subsys)
468
469 keep_space_group = .false.
470 IF (PRESENT(spgr)) THEN
471 IF (ASSOCIATED(spgr)) keep_space_group = spgr%keep_space_group
472 END IF
473
474 ! applies rotation matrices to coordinates
475 IF (keep_space_group) THEN
476 CALL spgr_apply_rotations_coord(spgr, optimizer%x)
477 END IF
478
479 xold = optimizer%x
480 t_old = m_walltime()
481
482 IF (optimizer%status >= 4) THEN
483 cpwarn("status>=4, trying to restart")
484 optimizer%status = 0
485 dataunit = cp_print_key_unit_nr(logger, geo_section, &
486 "PRINT%PROGRAM_RUN_INFO", extension=".geoLog")
487 IF (is_master) THEN
488 optimizer%task = 'START'
489 CALL setulb(SIZE(optimizer%x), optimizer%m, optimizer%x, &
490 optimizer%lower_bound, optimizer%upper_bound, &
491 optimizer%kind_of_bound, optimizer%f, optimizer%gradient, &
492 optimizer%wanted_relative_f_delta, &
493 optimizer%wanted_projected_gradient, optimizer%work_array, &
494 optimizer%i_work_array, optimizer%task, optimizer%print_every, &
495 optimizer%csave, optimizer%lsave, optimizer%isave, &
496 optimizer%dsave, optimizer%trust_radius, spgr=spgr, iwunit=dataunit)
497 END IF
498 CALL cp_print_key_finished_output(dataunit, logger, geo_section, &
499 "PRINT%PROGRAM_RUN_INFO")
500 END IF
501
502 DO
503 dataunit = cp_print_key_unit_nr(logger, geo_section, &
504 "PRINT%PROGRAM_RUN_INFO", extension=".geoLog")
505 ifmaster: IF (is_master) THEN
506 IF (optimizer%task(1:7) == 'RESTART') THEN
507 ! restart the optimizer
508 optimizer%status = 0
509 optimizer%task = 'START'
510 ! applies rotation matrices to coordinates and forces
511 IF (keep_space_group) THEN
512 CALL spgr_apply_rotations_coord(spgr, optimizer%x)
513 CALL spgr_apply_rotations_force(spgr, optimizer%gradient)
514 END IF
515 CALL setulb(SIZE(optimizer%x), optimizer%m, optimizer%x, &
516 optimizer%lower_bound, optimizer%upper_bound, &
517 optimizer%kind_of_bound, optimizer%f, optimizer%gradient, &
518 optimizer%wanted_relative_f_delta, &
519 optimizer%wanted_projected_gradient, optimizer%work_array, &
520 optimizer%i_work_array, optimizer%task, optimizer%print_every, &
521 optimizer%csave, optimizer%lsave, optimizer%isave, &
522 optimizer%dsave, optimizer%trust_radius, spgr=spgr, iwunit=dataunit)
523 IF (keep_space_group) THEN
524 CALL spgr_apply_rotations_coord(spgr, optimizer%x)
525 CALL spgr_apply_rotations_force(spgr, optimizer%gradient)
526 END IF
527 END IF
528 IF (optimizer%task(1:2) == 'FG') THEN
529 IF (optimizer%isave(36) > optimizer%max_f_per_iter) THEN
530 optimizer%task = 'STOP: CPU, hit max f eval in iter'
531 optimizer%status = 5 ! anormal exit
532 CALL setulb(SIZE(optimizer%x), optimizer%m, optimizer%x, &
533 optimizer%lower_bound, optimizer%upper_bound, &
534 optimizer%kind_of_bound, optimizer%f, optimizer%gradient, &
535 optimizer%wanted_relative_f_delta, &
536 optimizer%wanted_projected_gradient, optimizer%work_array, &
537 optimizer%i_work_array, optimizer%task, optimizer%print_every, &
538 optimizer%csave, optimizer%lsave, optimizer%isave, &
539 optimizer%dsave, optimizer%trust_radius, spgr=spgr, iwunit=dataunit)
540 ELSE
541 optimizer%status = 1
542 END IF
543 ELSE IF (optimizer%task(1:5) == 'NEW_X') THEN
544 IF (justentred) THEN
545 optimizer%status = 2
546 ! applies rotation matrices to coordinates and forces
547 IF (keep_space_group) THEN
548 CALL spgr_apply_rotations_coord(spgr, optimizer%x)
549 CALL spgr_apply_rotations_force(spgr, optimizer%gradient)
550 END IF
551 CALL setulb(SIZE(optimizer%x), optimizer%m, optimizer%x, &
552 optimizer%lower_bound, optimizer%upper_bound, &
553 optimizer%kind_of_bound, optimizer%f, optimizer%gradient, &
554 optimizer%wanted_relative_f_delta, &
555 optimizer%wanted_projected_gradient, optimizer%work_array, &
556 optimizer%i_work_array, optimizer%task, optimizer%print_every, &
557 optimizer%csave, optimizer%lsave, optimizer%isave, &
558 optimizer%dsave, optimizer%trust_radius, spgr=spgr, iwunit=dataunit)
559 IF (keep_space_group) THEN
560 CALL spgr_apply_rotations_coord(spgr, optimizer%x)
561 CALL spgr_apply_rotations_force(spgr, optimizer%gradient)
562 END IF
563 ELSE
564 ! applies rotation matrices to coordinates and forces
565 IF (keep_space_group) THEN
566 CALL spgr_apply_rotations_coord(spgr, optimizer%x)
567 CALL spgr_apply_rotations_force(spgr, optimizer%gradient)
568 END IF
569 optimizer%status = 3
570 END IF
571 ELSE IF (optimizer%task(1:4) == 'CONV') THEN
572 optimizer%status = 4
573 ELSE IF (optimizer%task(1:4) == 'STOP') THEN
574 optimizer%status = 5
575 cpwarn("task became stop in an unknown way")
576 ELSE IF (optimizer%task(1:5) == 'ERROR') THEN
577 optimizer%status = 5
578 ELSE
579 cpwarn("unknown task '"//optimizer%task//"'")
580 END IF
581 END IF ifmaster
582 CALL cp_print_key_finished_output(dataunit, logger, geo_section, &
583 "PRINT%PROGRAM_RUN_INFO")
584 CALL optimizer%para_env%bcast(optimizer%status, optimizer%master)
585 ! Dump info
586 IF (optimizer%status == 3) THEN
587 its = 0
588 IF (is_master) THEN
589 ! Iteration level is taken into account in the optimizer external loop
590 its = optimizer%isave(30)
591 END IF
592 END IF
593 !
594 SELECT CASE (optimizer%status)
595 CASE (1)
596 !op=1 evaluate f and g
597 CALL cp_eval_at(optimizer%obj_funct, x=optimizer%x, &
598 f=optimizer%f, &
599 gradient=optimizer%gradient, &
600 master=optimizer%master, para_env=optimizer%para_env)
601 ! do not use keywords?
602 dataunit = cp_print_key_unit_nr(logger, geo_section, &
603 "PRINT%PROGRAM_RUN_INFO", extension=".geoLog")
604 IF (is_master) THEN
605 ! applies rotation matrices to coordinates and forces
606 IF (keep_space_group) THEN
607 CALL spgr_apply_rotations_coord(spgr, optimizer%x)
608 CALL spgr_apply_rotations_force(spgr, optimizer%gradient)
609 END IF
610 CALL setulb(SIZE(optimizer%x), optimizer%m, optimizer%x, &
611 optimizer%lower_bound, optimizer%upper_bound, &
612 optimizer%kind_of_bound, optimizer%f, optimizer%gradient, &
613 optimizer%wanted_relative_f_delta, &
614 optimizer%wanted_projected_gradient, optimizer%work_array, &
615 optimizer%i_work_array, optimizer%task, optimizer%print_every, &
616 optimizer%csave, optimizer%lsave, optimizer%isave, &
617 optimizer%dsave, optimizer%trust_radius, spgr=spgr, iwunit=dataunit)
618 IF (keep_space_group) THEN
619 CALL spgr_apply_rotations_coord(spgr, optimizer%x)
620 CALL spgr_apply_rotations_force(spgr, optimizer%gradient)
621 END IF
622 END IF
623 CALL cp_print_key_finished_output(dataunit, logger, geo_section, &
624 "PRINT%PROGRAM_RUN_INFO")
625 CALL optimizer%para_env%bcast(optimizer%x, optimizer%master)
626 CASE (2)
627 !op=2 begin new iter
628 CALL optimizer%para_env%bcast(optimizer%x, optimizer%master)
629 t_old = m_walltime()
630 CASE (3)
631 !op=3 ended iter
632 wildcard = "LBFGS"
633 dataunit = cp_print_key_unit_nr(logger, geo_section, &
634 "PRINT%PROGRAM_RUN_INFO", extension=".geoLog")
635 IF (is_master) its = optimizer%isave(30)
636 CALL optimizer%para_env%bcast(its, optimizer%master)
637
638 ! Some IO and Convergence check
639 t_now = m_walltime()
640 t_diff = t_now - t_old
641 t_old = t_now
642 CALL gopt_f_io(optimizer%obj_funct, force_env, force_env%root_section, &
643 its, optimizer%f, dataunit, optimizer%eold, optimizer%emin, wildcard, gopt_param, &
644 SIZE(optimizer%x), optimizer%x - xold, optimizer%gradient, conv, used_time=t_diff)
645 CALL optimizer%para_env%bcast(conv, optimizer%master)
646 CALL cp_print_key_finished_output(dataunit, logger, geo_section, &
647 "PRINT%PROGRAM_RUN_INFO")
648 optimizer%eold = optimizer%f
649 optimizer%emin = min(optimizer%emin, optimizer%eold)
650 xold = optimizer%x
651 IF (PRESENT(converged)) converged = conv
652 EXIT
653 CASE (4)
654 !op=4 (convergence - normal exit)
655 ! Specific L-BFGS convergence criteria.. overrides the convergence criteria on
656 ! stepsize and gradients
657 dataunit = cp_print_key_unit_nr(logger, geo_section, &
658 "PRINT%PROGRAM_RUN_INFO", extension=".geoLog")
659 IF (dataunit > 0) THEN
660 WRITE (dataunit, '(T2,A)') ""
661 WRITE (dataunit, '(T2,A)') "************************************************"
662 WRITE (dataunit, '(T2,A)') "* Specific L-BFGS convergence criteria *"
663 WRITE (dataunit, '(T2,A)') "* WANTED_PROJ_GRADIENT and WANTED_REL_F_ERROR *"
664 WRITE (dataunit, '(T2,A)') "* satisfied .... run CONVERGED! *"
665 WRITE (dataunit, '(T2,A)') "* * * * *"
666 WRITE (dataunit, '(T2,A)') "* General convergence criteria on stepsize and *"
667 WRITE (dataunit, '(T2,A)') "* gradients may or may not have been satisfied *"
668 WRITE (dataunit, '(T2,A)') "* yet; if unsatisfactory, try tightening the *"
669 WRITE (dataunit, '(T2,A)') "* L-BFGS convergence criteria and restart run. *"
670 WRITE (dataunit, '(T2,A)') "************************************************"
671 WRITE (dataunit, '(T2,A)') ""
672 END IF
673 CALL cp_print_key_finished_output(dataunit, logger, geo_section, &
674 "PRINT%PROGRAM_RUN_INFO")
675 IF (PRESENT(converged)) converged = .true.
676 EXIT
677 CASE (5)
678 ! op=5 abnormal exit ()
679 CALL optimizer%para_env%bcast(optimizer%task, optimizer%master)
680 CASE (6)
681 ! deallocated
682 cpabort("step on a deallocated opt structure ")
683 CASE default
684 CALL cp_abort(__location__, &
685 "unknown status "//cp_to_string(optimizer%status))
686 optimizer%status = 5
687 EXIT
688 END SELECT
689 IF (optimizer%status == 1 .AND. justentred) THEN
690 optimizer%eold = optimizer%f
691 optimizer%emin = optimizer%eold
692 END IF
693 justentred = .false.
694 END DO
695
696 CALL optimizer%para_env%bcast(optimizer%x, optimizer%master)
697 CALL cp_opt_gopt_bcast_res(optimizer, &
698 n_iter=optimizer%n_iter, &
699 f=optimizer%f, last_f=optimizer%last_f, &
700 projected_gradient=optimizer%projected_gradient)
701
702 DEALLOCATE (xold)
703 IF (PRESENT(f)) f = optimizer%f
704 IF (PRESENT(last_f)) last_f = optimizer%last_f
705 IF (PRESENT(projected_gradient)) projected_gradient = optimizer%projected_gradient
706 IF (PRESENT(n_iter)) n_iter = optimizer%n_iter
707 CALL timestop(handle)
708
709 END SUBROUTINE cp_opt_gopt_step
710
711! **************************************************************************************************
712!> \brief returns the results (and broadcasts them)
713!> \param optimizer the optimizer object the info is taken from
714!> \param n_iter the number of iterations
715!> \param f the actual value of the objective function (f)
716!> \param last_f the last value of f
717!> \param projected_gradient the infinity norm of the projected gradient
718!> \par History
719!> none
720!> \author Fawzi Mohamed
721!> @version 2.2002
722!> \note
723!> private routine
724! **************************************************************************************************
725 SUBROUTINE cp_opt_gopt_bcast_res(optimizer, n_iter, f, last_f, &
726 projected_gradient)
727 TYPE(cp_lbfgs_opt_gopt_type), INTENT(IN) :: optimizer
728 INTEGER, INTENT(out), OPTIONAL :: n_iter
729 REAL(kind=dp), INTENT(inout), OPTIONAL :: f, last_f, projected_gradient
730
731 REAL(kind=dp), DIMENSION(4) :: results
732
733 IF (optimizer%master == optimizer%para_env%mepos) THEN
734 results = [real(optimizer%isave(30), kind=dp), &
735 optimizer%f, optimizer%dsave(2), optimizer%dsave(13)]
736 END IF
737 CALL optimizer%para_env%bcast(results, optimizer%master)
738 IF (PRESENT(n_iter)) n_iter = nint(results(1))
739 IF (PRESENT(f)) f = results(2)
740 IF (PRESENT(last_f)) last_f = results(3)
741 IF (PRESENT(projected_gradient)) projected_gradient = results(4)
742
743 END SUBROUTINE cp_opt_gopt_bcast_res
744
745! **************************************************************************************************
746!> \brief goes to the next optimal point (after an optimizer iteration)
747!> returns true if converged
748!> \param optimizer the optimizer that goes to the next point
749!> \param n_iter ...
750!> \param f ...
751!> \param last_f ...
752!> \param projected_gradient ...
753!> \param converged ...
754!> \param geo_section ...
755!> \param force_env ...
756!> \param gopt_param ...
757!> \param spgr ...
758!> \return ...
759!> \par History
760!> 01.2020 modified [pcazade]
761!> \author Fawzi Mohamed
762!> @version 2.2002
763!> \note
764!> if you deactivate convergence control it returns never false
765! **************************************************************************************************
766 FUNCTION cp_opt_gopt_next(optimizer, n_iter, f, last_f, &
767 projected_gradient, converged, geo_section, force_env, &
768 gopt_param, spgr) RESULT(res)
769 TYPE(cp_lbfgs_opt_gopt_type), INTENT(INOUT) :: optimizer
770 INTEGER, INTENT(out), OPTIONAL :: n_iter
771 REAL(kind=dp), INTENT(out), OPTIONAL :: f, last_f, projected_gradient
772 LOGICAL, INTENT(out) :: converged
773 TYPE(section_vals_type), POINTER :: geo_section
774 TYPE(force_env_type), POINTER :: force_env
775 TYPE(gopt_param_type), POINTER :: gopt_param
776 TYPE(spgr_type), OPTIONAL, POINTER :: spgr
777 LOGICAL :: res
778
779 ! passes spgr structure if present
780 CALL cp_opt_gopt_step(optimizer, n_iter=n_iter, f=f, &
781 last_f=last_f, projected_gradient=projected_gradient, &
782 converged=converged, geo_section=geo_section, &
783 force_env=force_env, gopt_param=gopt_param, spgr=spgr)
784 res = (optimizer%status < 40) .AND. .NOT. converged
785
786 END FUNCTION cp_opt_gopt_next
787
788! **************************************************************************************************
789!> \brief stops the optimization
790!> \param optimizer ...
791!> \par History
792!> none
793!> \author Fawzi Mohamed
794!> @version 2.2002
795! **************************************************************************************************
796 SUBROUTINE cp_opt_gopt_stop(optimizer)
797 TYPE(cp_lbfgs_opt_gopt_type), INTENT(INOUT) :: optimizer
798
799 optimizer%task = 'STOPPED on user request'
800 optimizer%status = 4 ! normal exit
801 IF (optimizer%master == optimizer%para_env%mepos) THEN
802 CALL setulb(SIZE(optimizer%x), optimizer%m, optimizer%x, &
803 optimizer%lower_bound, optimizer%upper_bound, &
804 optimizer%kind_of_bound, optimizer%f, optimizer%gradient, &
805 optimizer%wanted_relative_f_delta, &
806 optimizer%wanted_projected_gradient, optimizer%work_array, &
807 optimizer%i_work_array, optimizer%task, optimizer%print_every, &
808 optimizer%csave, optimizer%lsave, optimizer%isave, &
809 optimizer%dsave, optimizer%trust_radius)
810 END IF
811
812 END SUBROUTINE cp_opt_gopt_stop
813
814END MODULE cp_lbfgs_optimizer_gopt
cp_log_handling::cp_to_string
Definition cp_log_handling.F:90
cp_lbfgs_optimizer_gopt
routines that optimize a functional using the limited memory bfgs quasi-newton method....
Definition cp_lbfgs_optimizer_gopt.F:22
cp_lbfgs_optimizer_gopt::cp_opt_gopt_create
subroutine, public cp_opt_gopt_create(optimizer, para_env, obj_funct, x0, m, print_every, wanted_relative_f_delta, wanted_projected_gradient, lower_bound, upper_bound, kind_of_bound, master, max_f_per_iter, trust_radius)
initializes the optimizer
Definition cp_lbfgs_optimizer_gopt.F:172
cp_lbfgs_optimizer_gopt::cp_opt_gopt_next
logical function, public cp_opt_gopt_next(optimizer, n_iter, f, last_f, projected_gradient, converged, geo_section, force_env, gopt_param, spgr)
goes to the next optimal point (after an optimizer iteration) returns true if converged
Definition cp_lbfgs_optimizer_gopt.F:769
cp_lbfgs_optimizer_gopt::cp_opt_gopt_stop
subroutine, public cp_opt_gopt_stop(optimizer)
stops the optimization
Definition cp_lbfgs_optimizer_gopt.F:797
cp_lbfgs_optimizer_gopt::cp_opt_gopt_release
subroutine, public cp_opt_gopt_release(optimizer)
releases the optimizer (see doc/ReferenceCounting.html)
Definition cp_lbfgs_optimizer_gopt.F:282
cp_lbfgs
LBFGS-B routine (version 3.0, April 25, 2011).
Definition cp_lbfgs.F:19
cp_lbfgs::setulb
subroutine, public setulb(n, m, x, lower_bound, upper_bound, nbd, f, g, factr, pgtol, wa, iwa, task, iprint, csave, lsave, isave, dsave, trust_radius, spgr, iwunit)
This subroutine partitions the working arrays wa and iwa, and then uses the limited memory BFGS metho...
Definition cp_lbfgs.F:188
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:1073
cp_subsys_types
types that represent a subsys, i.e. a part of the system
Definition cp_subsys_types.F:16
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:344
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::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:505
gopt_f_methods::cp_eval_at
subroutine, public 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_f_methods.F:173
gopt_f_types
contains a functional that calculates the energy and its derivatives for the geometry optimizer
Definition gopt_f_types.F:15
gopt_f_types::gopt_f_retain
subroutine, public gopt_f_retain(gopt_env)
...
Definition gopt_f_types.F:203
gopt_f_types::gopt_f_release
recursive subroutine, public gopt_f_release(gopt_env)
...
Definition gopt_f_types.F:217
gopt_param_types
contains typo and related routines to handle parameters controlling the GEO_OPT module
Definition gopt_param_types.F:14
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
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
machine
Machine interface based on Fortran 2003 and POSIX.
Definition machine.F:17
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:141
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:27
message_passing::mp_para_env_release
subroutine, public mp_para_env_release(para_env)
releases the para object (to be called when you don't want anymore the shared copy of this object)
Definition message_passing.F:2163
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::spgr_apply_rotations_coord
subroutine, public spgr_apply_rotations_coord(spgr, coord)
routine applies the rotation matrices to the coordinates.
Definition space_groups.F:610
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:671
cp_lbfgs_optimizer_gopt::cp_lbfgs_opt_gopt_type
info for the optimizer (see the description of this module)
Definition cp_lbfgs_optimizer_gopt.F:130
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
gopt_f_types::gopt_f_type
calculates the potential energy of a system, and its derivatives
Definition gopt_f_types.F:59
gopt_param_types::gopt_param_type
Definition gopt_param_types.F:55
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:768
space_groups_types::spgr_type
Definition space_groups_types.F:27