opt_control_global.F90

Go to the documentation of this file.

!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
!!
!! This program is free software; you can redistribute it and/or modify
!! it under the terms of the GNU General Public License as published by
!! the Free Software Foundation; either version 2, or (at your option)
!! any later version.
!!
!! This program is distributed in the hope that it will be useful,
!! but WITHOUT ANY WARRANTY; without even the implied warranty of
!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!! GNU General Public License for more details.
!!
!! You should have received a copy of the GNU General Public License
!! along with this program; if not, write to the Free Software
!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
!! 02110-1301, USA.
!!
 
#include "global.h"
 
module opt_control_global_oct_m
  use debug_oct_m
  use global_oct_m
  use, intrinsic :: iso_fortran_env
  use messages_oct_m
  use namespace_oct_m
  use parser_oct_m
  use varinfo_oct_m
 
  implicit none
 
  private
  public ::                   &
    oct_t,                    &
    oct_read_inp,             &
    oct_algorithm_is_direct
 
  type oct_t
    ! Components are public by default
    integer :: algorithm
    logical :: mode_fixed_fluence
    real(real64)   :: eta, delta
    real(real64)   :: direct_step
    logical :: oct_double_check
    real(real64)   :: check_gradient
    integer :: number_checkpoints
    logical :: random_initial_guess
  end type oct_t
 
contains
 
  subroutine oct_read_inp(oct, namespace)
    type(oct_t),       intent(inout) :: oct
    type(namespace_t), intent(in)    :: namespace
 
    push_sub(oct_read_inp)
 
    call messages_print_with_emphasis(msg="OCT run mode", namespace=namespace)
    call messages_obsolete_variable(namespace, 'OCTControlRepresentation')
 
    !%Variable OCTScheme
    !%Type integer
    !%Section Calculation Modes::Optimal Control
    !%Default oct_zbr98
    !%Description
    !% Optimal Control Theory can be performed with <tt>Octopus</tt> with a variety of different
    !% algorithms. Not all of them can be used with any choice of target or control function
    !% representation. For example, some algorithms cannot be used if
    !% <tt>OCTControlRepresentation = control_function_real_time</tt>
    !% (<tt>OCTScheme</tt> > <tt>oct_straight_iteration</tt>), and others cannot be used
    !% if <tt>OCTControlRepresentation = control_function_parametrized</tt>
    !% (<tt>OCTScheme</tt>  <  <tt>oct_straight_iteration</tt>).
    !%Option oct_zbr98 1
    !% Backward-Forward-Backward scheme described in <i>JCP</i> <b>108</b>, 1953 (1998).
    !% Only possible if target operator is a projection operator.
    !% Provides fast, stable and monotonic convergence.
    !%Option oct_zr98  2
    !% Forward-Backward-Forward scheme described in <i>JCP</i> <b>109</b>, 385 (1998).
    !% Works for projection and more general target operators also. The convergence is
    !% stable but slower than ZBR98.
    !% Note that local operators show an extremely slow convergence. It ensures monotonic
    !% convergence.
    !%Option oct_wg05  3
    !% Forward-Backward scheme described in <i>J. Opt. B.</i> <b>7</b>, 300 (2005).
    !% Works for all kinds of target operators, can be used with all kinds of filters, and
    !% allows a fixed fluence.
    !% The price is a rather unstable convergence.
    !% If the restrictions set by the filter and fluence are reasonable, a good overlap can be
    !% expected within 20 iterations.
    !% No monotonic convergence.
    !%Option oct_mt03 4
    !% Basically an improved and generalized scheme.
    !% Comparable to ZBR98/ZR98. See [Y. Maday and G. Turinici, <i>J. Chem. Phys.</i> <b>118</b>, 8191 (2003)].
    !%Option oct_krotov 5
    !% The procedure reported in [D. Tannor, V. Kazakov and V.
    !% Orlov, in <i>Time-Dependent Quantum Molecular Dynamics</i>, edited by J. Broeckhove
    !% and L. Lathouweres (Plenum, New York, 1992), pp. 347-360].
    !%Option oct_straight_iteration 6
    !% Straight iteration: one forward and one backward propagation is performed at each
    !% iteration, both with the same control field. An output field is calculated with the
    !% resulting wavefunctions.
    !%Option oct_cg 7
    !% Conjugate-gradients, as implemented in the GNU GSL library. In particular, the
    !% Fletcher-Reeves version.
    !% The seed for the random number generator can be modified by setting
    !% <tt>GSL_RNG_SEED</tt> environment variable.
    !%Option oct_bfgs 8
    !% The methods use the vector Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm.
    !% Also, it calls the GNU GSL library version of the algorithm. It is a quasi-Newton
    !% method which builds up an approximation to the second derivatives of the function using
    !% the difference between successive gradient vectors.  By combining the first and second
    !% derivatives the algorithm is able to take Newton-type steps towards the function minimum,
    !% assuming quadratic behavior in that region. We have chosen to implement the "bfgs2" version,
    !% as GSL calls it, which is supposed to be the most efficient version available, and a faithful
    !% implementation of the line minimization scheme described in "Practical Methods of Optimization",
    !% (Fletcher), Algorithms 2.6.2 and 2.6.4.
    !%Option oct_direct 9
    !% This is a "direct" optimization scheme. This means that we do not make use of the
    !% "usual" QOCT equations (backward-forward propagations, etc), but we use some gradient-free
    !% maximization algorithm for the function that we want to optimize. In this case, the
    !% maximization algorithm is the Nelder-Mead algorithm as implemented in the GSL. The function
    !% values are obtained by successive forward propagations.
    !% The seed for the random number generator can be modified by setting
    !% <tt>GSL_RNG_SEED</tt> environment variable.
    !%Option oct_nlopt_bobyqa 11
    !% The BOBYQA algorithm, as implemented in the NLOPT library -- therefore, octopus has to
    !% be compiled with it in order to be able to use this option.
    !% The seed for the random number generator can be modified by setting
    !% <tt>GSL_RNG_SEED</tt> environment variable.
    !%Option oct_nlopt_lbfgs 12
    !% The local BFGS, as implemented in the NLOPT library -- therefore, octopus has to
    !% be compiled with it in order to be able to use this option.
    !% The seed for the random number generator can be modified by setting
    !% <tt>GSL_RNG_SEED</tt> environment variable.
    !%End
    call parse_variable(namespace, 'OCTScheme', option__octscheme__oct_zr98, oct%algorithm)
    if (.not. varinfo_valid_option('OCTScheme', oct%algorithm)) call messages_input_error(namespace, 'OCTScheme')
    ! We must check that the algorithm is consistent with OCTControlRepresentation, i.e.
    ! some algorithms only make sense if the control functions are handled directly in real
    ! time, some others make only sense if the control functions are parameterized.
    ! This check cannot be here any more, and it should be placed somewhere else.
    call messages_print_var_option("OCTScheme", oct%algorithm, namespace=namespace)
    select case (oct%algorithm)
    case (option__octscheme__oct_mt03)
      oct%delta = m_two
      oct%eta = m_zero
      !%Variable OCTEta
      !%Type float
      !%Section Calculation Modes::Optimal Control
      !%Default 1.0
      !%Description
      !% If <tt>OCTScheme = oct_mt03</tt>, then you can supply the "eta" and "delta" parameters
      !% described in [Y. Maday and G. Turinici, <i>J. Chem. Phys.</i> <b>118</b>, 8191 (2003)], using the
      !% <tt>OCTEta</tt> and <tt>OCTDelta</tt> variables.
      !%End
      call parse_variable(namespace, 'OCTEta', m_one, oct%eta)
      !%Variable OCTDelta
      !%Type float
      !%Section Calculation Modes::Optimal Control
      !%Default 0.0
      !%Description
      !% If <tt>OCTScheme = oct_mt03</tt>, then you can supply the "eta" and "delta" parameters
      !% described in [Y. Maday and G. Turinici, <i>J. Chem. Phys.</i> <b>118</b>, 8191 (2003)], using the
      !% <tt>OCTEta</tt> and <tt>OCTDelta</tt> variables.
      !%End
      call parse_variable(namespace, 'OCTDelta', m_zero, oct%delta)
 
    case (option__octscheme__oct_zr98)
      oct%delta = m_one
      oct%eta = m_one
    case (option__octscheme__oct_krotov)
      oct%delta = m_one
      oct%eta = m_zero
    case (option__octscheme__oct_straight_iteration)
      oct%delta = m_zero
      oct%eta = m_one
    case (option__octscheme__oct_cg)
      oct%delta = m_zero
      oct%eta = m_one
    case (option__octscheme__oct_bfgs)
      oct%delta = m_zero
      oct%eta = m_one
    case (option__octscheme__oct_direct)
      ! The use of these variables for the direct and bobyqa schemes remain undocumented for the moment.
      call parse_variable(namespace, 'OCTEta', m_one, oct%eta)
      call parse_variable(namespace, 'OCTDelta', m_zero, oct%delta)
    case (option__octscheme__oct_nlopt_bobyqa, option__octscheme__oct_nlopt_lbfgs)
#if defined(HAVE_NLOPT)
      !WARNING: not clear if this is needed, probably not.
      call parse_variable(namespace, 'OCTEta', m_one, oct%eta)
      call parse_variable(namespace, 'OCTDelta', m_zero, oct%delta)
#else
      write(message(1), '(a)') '"OCTScheme = oct_nlopt_bobyqa" or "OCTScheme = oct_nlopt_lbfgs" are'
      write(message(2), '(a)') ' only possible if the nlopt library has been compiled.'
      call messages_fatal(2, namespace=namespace)
#endif
    case default
      oct%delta = m_one
      oct%eta = m_one
    end select
 
    !%Variable OCTDoubleCheck
    !%Type logical
    !%Section Calculation Modes::Optimal Control
    !%Default true
    !%Description
    !% In order to make sure that the optimized field indeed does its job, the code
    !% may run a normal propagation after the optimization using the optimized field.
    !%End
    call parse_variable(namespace, 'OCTDoubleCheck', .true., oct%oct_double_check)
    call messages_print_var_value("OCTDoubleCheck", oct%oct_double_check, namespace=namespace)
 
 
    !%Variable OCTCheckGradient
    !%Type float
    !%Section Calculation Modes::Optimal Control
    !%Default 0.0
    !%Description
    !% When doing QOCT with the conjugate-gradient optimization scheme, the gradient is
    !% computed thanks to a forward-backwards propagation. For debugging purposes, this
    !% gradient can be compared with the value obtained "numerically" (<i>i.e.</i> by doing
    !% successive forward propagations with control fields separated by small finite
    !% differences).
    !%
    !% In order to activate this feature, set <tt>OCTCheckGradient</tt> to some non-zero value,
    !% which will be the finite difference used to numerically compute the gradient.
    !%End
    call parse_variable(namespace, 'OCTCheckGradient', 0.0_real64, oct%check_gradient)
    call messages_print_var_value("OCTCheckGradient", oct%check_gradient, namespace=namespace)
 
 
    !%Variable OCTDirectStep
    !%Type float
    !%Section Calculation Modes::Optimal Control
    !%Default 0.25
    !%Description
    !% If you choose <tt>OCTScheme = oct_direct</tt> or <tt>OCTScheme = oct_nlopt_bobyqa</tt>,
    !% the algorithms necessitate an initial "step" to perform the direct search for the
    !% optimal value. The precise meaning of this "step" differs.
    !%End
    call parse_variable(namespace, 'OCTDirectStep', 0.25_real64, oct%direct_step)
    call messages_print_var_value("OCTDirectStep", oct%direct_step, namespace=namespace)
 
    !%Variable OCTNumberCheckPoints
    !%Type integer
    !%Section Calculation Modes::Optimal Control
    !%Default 0
    !%Description
    !% During an OCT propagation, the code may write the wavefunctions at some time steps (the
    !% "check points"). When the inverse backward or forward propagation
    !% is performed in a following step, the wavefunction should reverse its path
    !% (almost) exactly. This can be checked to make sure that it is the case -- otherwise
    !% one should try reducing the time-step, or altering in some other way the
    !% variables that control the propagation.
    !%
    !% If the backward (or forward) propagation is not retracing the steps of the previous
    !% forward (or backward) propagation, the code will write a warning.
    !%End
    call parse_variable(namespace, 'OCTNumberCheckPoints', 0, oct%number_checkpoints)
    call messages_print_var_value("OCTNumberCheckPoints", oct%number_checkpoints, namespace=namespace)
 
    !%Variable OCTRandomInitialGuess
    !%Type logical
    !%Section Calculation Modes::Optimal Control
    !%Default false
    !%Description
    !% The initial field to start the optimization search is usually given in the <tt>inp</tt> file,
    !% through a <tt>TDExternalFields</tt> block. However, you can start from a random guess if you
    !% set this variable to true.
    !%
    !% Note, however, that this is only valid for the "direct" optimization schemes; moreover
    !% you still need to provide a <tt>TDExternalFields</tt> block.
    !%End
    call parse_variable(namespace, 'OCTRandomInitialGuess', .false., oct%random_initial_guess)
    call messages_print_var_value("OCTRandomInitialGuess", oct%random_initial_guess, namespace=namespace)
 
    call messages_print_with_emphasis(namespace=namespace)
 
    pop_sub(oct_read_inp)
  end subroutine oct_read_inp
  ! ---------------------------------------------------------
 
 
 
  logical pure function oct_algorithm_is_direct(oct)
    type(oct_t), intent(in) :: oct
    oct_algorithm_is_direct = (oct%algorithm == option__octscheme__oct_direct) .or. &
      (oct%algorithm == option__octscheme__oct_nlopt_bobyqa)
  end function oct_algorithm_is_direct
  ! ---------------------------------------------------------
 
 
end module opt_control_global_oct_m
 
!! Local Variables:
!! mode: f90
!! coding: utf-8
!! End: