37 use,
intrinsic :: iso_fortran_env
77 type(filter_t),
save :: filter
78 type(oct_t),
save :: oct
79 type(oct_iterator_t),
save :: iterator
80 type(target_t),
save :: oct_target
81 type(opt_control_state_t),
save :: initial_st
85 type(controlfunction_t),
save :: par_
86 type(electrons_t),
pointer :: sys_
87 type(hamiltonian_elec_t),
pointer :: hm_
88 type(td_t),
pointer :: td_
89 real(real64),
allocatable :: x_(:)
96 class(*),
intent(inout) :: system
102 message(1) =
"CalculationMode = opt_control not implemented for multi-system calculations"
114 type(electrons_t),
target,
intent(inout) :: sys
116 type(td_t),
target :: td
117 type(controlfunction_t) :: par, par_new, par_prev
120 type(oct_prop_t) :: prop_chi, prop_psi
121 type(states_elec_t) :: psi
123 type(lasers_t),
pointer :: lasers
127 if (sys%hm%pcm%run_pcm)
then
131 if (sys%kpoints%use_symmetries)
then
132 call messages_experimental(
"KPoints symmetries with CalculationMode = opt_control", namespace=sys%namespace)
141 call td_init(td, sys%namespace, sys%space, sys%gr, sys%ions, sys%st, sys%ks, sys%hm, sys%ext_partners, sys%outp)
159 if(
associated(lasers))
call laser_write_info(lasers%lasers, namespace=sys%namespace)
170 (oct%algorithm == option__octscheme__oct_zbr98), &
171 (oct%algorithm == option__octscheme__oct_cg) .or. &
172 (oct%algorithm == option__octscheme__oct_bfgs) .or. &
173 (oct%algorithm == option__octscheme__oct_nlopt_lbfgs))
177 call filter_init(td%max_iter, sys%namespace, td%dt, filter)
183 call target_init(sys%gr, sys%kpoints, sys%namespace, sys%space, sys%ions, initial_st, td, &
192 call output_states(sys%outp, sys%namespace, sys%space,
oct_dir//
'initial', psi, sys%gr, sys%ions, sys%hm, -1)
193 call target_output(oct_target, sys%namespace, sys%space, sys%gr,
oct_dir//
'target', sys%ions, sys%hm, sys%outp)
198 select case (oct%algorithm)
199 case (option__octscheme__oct_zbr98)
200 message(1) =
"Info: Starting OCT iteration using scheme: ZBR98"
203 case (option__octscheme__oct_wg05)
204 message(1) =
"Info: Starting OCT iteration using scheme: WG05"
207 case (option__octscheme__oct_zr98)
208 message(1) =
"Info: Starting OCT iteration using scheme: ZR98"
211 case (option__octscheme__oct_mt03)
212 message(1) =
"Info: Starting OCT iteration using scheme: MT03"
215 case (option__octscheme__oct_krotov)
216 message(1) =
"Info: Starting OCT iteration using scheme: KROTOV"
219 case (option__octscheme__oct_straight_iteration)
220 message(1) =
"Info: Starting OCT iterations using scheme: STRAIGHT ITERATION"
223 case (option__octscheme__oct_cg)
224 message(1) =
"Info: Starting OCT iterations using scheme: CONJUGATE GRADIENTS"
227 case (option__octscheme__oct_bfgs)
228 message(1) =
"Info: Starting OCT iterations using scheme: BFGS"
231 case (option__octscheme__oct_direct)
232 message(1) =
"Info: Starting OCT iterations using scheme: DIRECT OPTIMIZATION (NELDER-MEAD)"
235 case (option__octscheme__oct_nlopt_bobyqa)
236 message(1) =
"Info: Starting OCT iterations using scheme: DIRECT OPTIMIZATION (NLOPT - BOBYQA)"
239 case (option__octscheme__oct_nlopt_lbfgs)
240 message(1) =
"Info: Starting OCT iterations using scheme: DIRECT OPTIMIZATION (NLOPT - LBFGS)"
274 if (
clean_stop(sys%mc%master_comm) .or. stop_loop)
exit ctr_loop
291 call oct_prop_init(prop_chi, sys%namespace,
"chi", sys%gr, sys%mc)
292 call oct_prop_init(prop_psi, sys%namespace,
"psi", sys%gr, sys%mc)
297 call f_iter(sys, td, psi, par, prop_psi, prop_chi, j1)
299 if (
clean_stop(sys%mc%master_comm) .or. stop_loop)
exit ctr_loop
317 call oct_prop_init(prop_chi, sys%namespace,
"chi", sys%gr, sys%mc)
318 call oct_prop_init(prop_psi, sys%namespace,
"psi", sys%gr, sys%mc)
320 if (oct%mode_fixed_fluence)
then
329 call f_wg05(sys, td, psi, par, prop_psi, prop_chi, j1)
331 if (
clean_stop(sys%mc%master_comm) .or. stop_loop)
exit ctr_loop
351 call oct_prop_init(prop_chi, sys%namespace,
"chi", sys%gr, sys%mc)
352 call oct_prop_init(prop_psi, sys%namespace,
"psi", sys%gr, sys%mc)
356 j1 =
target_j1(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi)
369 call f_zbr98(sys, td, qcpsi, prop_psi, prop_chi, par)
370 j1 =
target_j1(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi)
372 if (
clean_stop(sys%mc%master_comm) .or. stop_loop)
exit ctr_loop
387 integer :: dof, ierr, maxiter
388 real(real64):: step, minvalue
389 real(real64),
allocatable :: theta(:)
390 real(real64),
allocatable :: x(:)
419 safe_allocate(x(1:dof))
420 safe_allocate(theta(1:dof))
424 step = oct%direct_step *
m_pi
427 select case (oct%algorithm)
428 case (option__octscheme__oct_bfgs)
430 real(oct_iterator_tolerance(iterator), real64),
real(oct_iterator_tolerance(iterator), real64), &
431 maxiter, opt_control_cg_calc, opt_control_cg_write_info, minvalue, ierr)
432 case (option__octscheme__oct_cg)
434 real(oct_iterator_tolerance(iterator), real64),
real(oct_iterator_tolerance(iterator), real64), &
435 maxiter, opt_control_cg_calc, opt_control_cg_write_info, minvalue, ierr)
439 if (ierr <= 1024)
then
440 message(1) =
"Error occurred during the GSL minimization procedure:"
444 message(1) =
"The optimization did not meet the convergence criterion."
451 safe_deallocate_a(theta)
459 integer :: ierr, maxiter
460 real(real64):: minvalue, step
461 real(real64),
allocatable :: theta(:)
462 real(real64),
allocatable :: x(:)
484 safe_allocate(x(1:dim))
485 safe_allocate(theta(1:dim))
500 step = oct%direct_step *
m_pi
504 real(oct_iterator_tolerance(iterator), real64) , maxiter, &
505 opt_control_direct_calc, opt_control_direct_message_info, minvalue, ierr)
508 if (ierr <= 1024)
then
509 message(1) =
"Error occurred during the GSL minimization procedure:"
513 message(1) =
"The OCT direct optimization did not meet the convergence criterion."
520 safe_deallocate_a(theta)
528#if defined(HAVE_NLOPT)
529 integer :: method, dim, maxiter, ierr
530 real(real64),
allocatable :: x(:), xl(:), xu(:)
531 real(real64) :: step, toldr, minimum, f
550 safe_allocate(x(dim))
551 safe_allocate(xl(1:dim))
552 safe_allocate(xu(1:dim))
567 step = oct%direct_step
568 select case (oct%algorithm)
569 case (option__octscheme__oct_nlopt_bobyqa)
571 case (option__octscheme__oct_nlopt_lbfgs)
576 call minimize_multidim_nlopt(ierr, method, dim, x, step, toldr, maxiter,
opt_control_nlopt_func, minimum, &
578 if (ierr < 1 .or. ierr > 4)
then
579 message(1) =
"The nlopt minimization procedure did not find convergence, or found an error"
580 write(
message(2),
'(a,i5)')
"Error code =", ierr
585 safe_deallocate_a(xl)
586 safe_deallocate_a(xu)
598 subroutine f_zbr98(sys, td, qcpsi, prop_psi, prop_chi, par)
600 type(
td_t),
intent(inout) :: td
616 call bwd_step(sys, td, oct_target, par, par_chi, qcchi, prop_chi, prop_psi)
618 call fwd_step(sys, td, oct_target, par, par_chi, qcpsi, prop_chi, prop_psi)
628 subroutine f_wg05(sys, td, qcpsi, par, prop_psi, prop_chi, j1)
630 type(
td_t),
intent(inout) :: td
635 real(real64),
intent(out) :: j1
637 real(real64) :: new_penalty
646 j1 =
target_j1(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi)
655 call target_chi(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi, qcchi, sys%ions)
656 call bwd_step(sys, td, oct_target, par, parp, qcchi, prop_chi, prop_psi)
661 if (oct%mode_fixed_fluence)
then
673 j1 =
target_j1(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi)
685 type(
td_t),
intent(inout) :: td
687 real(real64),
intent(out) :: j1
696 call oct_prop_init(prop_chi, sys%namespace,
"chi", sys%gr, sys%mc)
697 call oct_prop_init(prop_psi, sys%namespace,
"psi", sys%gr, sys%mc)
709 j1 =
target_j1(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi, sys%ions)
714 call target_chi(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi, qcchi, sys%ions)
718 call bwd_step_2(sys, td, oct_target, par, par_chi, qcchi, prop_chi, prop_psi)
736 subroutine f_iter(sys, td, qcpsi, par, prop_psi, prop_chi, j1)
738 type(
td_t),
intent(inout) :: td
743 real(real64),
intent(out) :: j1
753 j1 =
target_j1(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi)
762 call target_chi(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi, qcchi, sys%ions)
763 call bwd_step(sys, td, oct_target, par, par_chi, qcchi, prop_chi, prop_psi)
766 call fwd_step(sys, td, oct_target, par, par_chi, qcpsi, prop_chi, prop_psi)
768 j1 =
target_j1(oct_target, sys%namespace, sys%gr, sys%kpoints, qcpsi)
776#include "opt_control_c_inc.F90"
777#include "check_input_inc.F90"
778#include "finalcheck_inc.F90"
Define which routines can be seen from the outside.
double sqrt(double __x) __attribute__((__nothrow__
subroutine minimize_multidim(method, dim, x, step, line_tol, tolgrad, toldr, maxiter, f, write_iter_info, minimum, ierr)
This module contains the definition of the data type that holds a "control function" used for OCT run...
subroutine, public controlfunction_filter(par, filter)
subroutine, public controlfunction_apply_envelope(cp)
subroutine, public controlfunction_set_fluence(par)
subroutine, public controlfunction_bounds(par, lower_bounds, upper_bounds)
subroutine, public controlfunction_to_realtime(par)
integer pure function, public controlfunction_dof(par)
real(real64) function, public controlfunction_fluence(par)
subroutine, public controlfunction_set_alpha(par, alpha)
subroutine, public controlfunction_write(filename, cp, namespace)
subroutine, public controlfunction_prepare_initial(par)
"Prepares" the initial guess control field: maybe it has to be normalized to a certain fluence,...
real(real64) pure function, public controlfunction_w0(par)
subroutine, public controlfunction_end(cp)
subroutine, public controlfunction_get_theta(par, theta)
real(real64) function, public controlfunction_j2(par)
subroutine, public controlfunction_randomize(par)
real(real64) pure function, public controlfunction_alpha(par, ipar)
real(real64) pure function, public controlfunction_targetfluence()
subroutine, public controlfunction_copy(cp_out, cp_in)
subroutine, public controlfunction_set_rep(par)
Transforms the control function to frequency space, if this is the space in which the functions are d...
subroutine, public controlfunction_mod_close()
subroutine, public controlfunction_to_h(cp, ext_partners)
subroutine, public controlfunction_init(cp, dt, ntiter)
Before using an controlfunction_t variable, it needs to be initialized, either by calling controlfunc...
subroutine, public controlfunction_mod_init(ext_partners, namespace, dt, max_iter, mode_fixed_fluence)
Initializes the module, should be the first subroutine to be called (the last one should be controlfu...
subroutine, public controlfunction_set(cp, ext_partners)
The external fields defined in epot_t "ep" are transferred to the control functions described in "cp"...
type(lasers_t) function, pointer, public list_get_lasers(partners)
subroutine, public filter_init(steps, namespace, dt, filter)
subroutine, public filter_write(filter, namespace)
subroutine, public filter_end(filter)
real(real64), parameter, public m_pi
some mathematical constants
character(len= *), parameter, public oct_dir
This module implements the underlying real-space grid.
subroutine, public initial_state_init(sys, qcstate)
subroutine, public io_mkdir(fname, namespace, parents)
A module to handle KS potential, without the external potential.
integer, parameter, public independent_particles
subroutine, public laser_write_info(lasers, namespace, dt, max_iter, iunit)
This module is intended to contain "only mathematical" functions and procedures.
subroutine, public messages_print_with_emphasis(msg, iunit, namespace)
subroutine, public messages_not_implemented(feature, namespace)
character(len=512), private msg
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
subroutine, public messages_input_error(namespace, var, details, row, column)
subroutine, public messages_experimental(name, namespace)
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
integer minmethod_nmsimplex
integer minmethod_nlopt_lbfgs
integer minmethod_nlopt_bobyqa
This module implements the basic mulsisystem class, a container system for other systems.
This module contains the definition of the oct_t data type, which contains some of the basic informat...
subroutine, public oct_read_inp(oct, namespace)
Reads, from the inp file, some global information about how the QOCT run should be.
logical function, public iteration_manager(namespace, j1, par, par_prev, iterator)
subroutine, public oct_iterator_init(iterator, namespace, par)
integer pure function, public oct_iterator_maxiter(iterator)
subroutine, public oct_iterator_end(iterator, namespace)
integer pure function, public oct_iterator_current(iterator)
subroutine, public iteration_manager_direct(j, par, iterator, sys, dx)
real(real64) pure function, public oct_iterator_tolerance(iterator)
This module contains the main procedure ("opt_control_run") that is used when optimal control runs ar...
subroutine, public opt_control_cg_calc(n, x, f, getgrad, df)
subroutine check_faulty_runmodes(sys, tr)
subroutine, public opt_control_run(system)
subroutine oct_finalcheck(sys, td)
subroutine, public opt_control_cg_write_info(iter, n, val, maxdx, maxdf, x)
interface is required by its being passed as dummy routine to minimize_multidim
subroutine f_striter(sys, td, par, j1)
subroutine, public opt_control_direct_message_info(iter, n, val, maxdx, x)
subroutine f_zbr98(sys, td, qcpsi, prop_psi, prop_chi, par)
subroutine f_iter(sys, td, qcpsi, par, prop_psi, prop_chi, j1)
subroutine opt_control_run_legacy(sys)
This is the main procedure for all types of optimal control runs. It is called from the "run" procedu...
subroutine f_wg05(sys, td, qcpsi, par, prop_psi, prop_chi, j1)
subroutine, public opt_control_function_forward(x, f)
subroutine opt_control_nlopt_func(val, n, x, grad, need_gradient, f_data)
subroutine, public opt_control_direct_calc(n, x, f)
No intents here is unfortunately required because this will be passed to newuoa routines as a dummy f...
This module holds the "opt_control_state_t" datatype, which contains a quantum-classical state.
subroutine, public opt_control_state_end(ocs)
subroutine, public opt_control_state_null(ocs)
subroutine, public opt_control_state_copy(ocsout, ocsin)
subroutine, public opt_control_state_init(ocs, qstate, ions)
subroutine, public opt_control_get_qs(qstate, ocs)
this module contains the output system
subroutine, public output_states(outp, namespace, space, dir, st, gr, ions, hm, iter)
subroutine, public propagation_mod_init(niter, eta, delta, number_checkpoints, zbr98, gradients)
This subroutine must be called before any QOCT propagations are done. It simply stores in the module ...
subroutine, public oct_prop_end(prop)
subroutine, public bwd_step(sys, td, tg, par, par_chi, qcchi, prop_chi, prop_psi)
subroutine, public propagate_forward(sys, td, par, tg, qcpsi, prop, write_iter)
subroutine, public fwd_step(sys, td, tg, par, par_chi, qcpsi, prop_chi, prop_psi)
subroutine, public oct_prop_init(prop, namespace, dirname, mesh, mc)
subroutine, public bwd_step_2(sys, td, tg, par, par_chi, qcchi, prop_chi, prop_psi)
subroutine, public propagator_elec_set_scf_prop(tr, threshold)
logical function, public clean_stop(comm)
returns true if a file named stop exists
This module handles spin dimensions of the states and the k-point distribution.
subroutine, public states_elec_end(st)
finalize the states_elec_t object
real(real64) function, public target_j1(tg, namespace, gr, kpoints, qcpsi, ions)
Calculates the J1 functional, i.e.: in the time-independent case, or else in the time-dependent cas...
subroutine, public target_init(gr, kpoints, namespace, space, ions, qcs, td, w0, tg, oct, ep, mc)
The target is initialized, mainly by reading from the inp file.
subroutine, public target_get_state(tg, st)
This just copies the states_elec_t variable present in target, into st.
subroutine, public target_output(tg, namespace, space, gr, dir, ions, hm, outp)
subroutine, public target_chi(tg, namespace, gr, kpoints, qcpsi_in, qcchi_out, ions)
Calculate .
subroutine, public target_end(tg, oct)
subroutine, public td_end(td)
subroutine, public td_init(td, namespace, space, gr, ions, st, ks, hm, ext_partners, outp)
subroutine scheme_zbr98()
subroutine scheme_direct()
subroutine scheme_straight_iteration()
subroutine scheme_nlopt()
This is the data type used to hold a control function.
Class describing the electron system.
Container class for lists of system_oct_m::system_t.
This is the datatype that contains the objects that are propagated: in principle this could be both t...
The states_elec_t class contains all electronic wave functions.