Line data    Source code

       1             : #pragma once
       2             : 
       3             : #include "Functional.h"
       4             : #include "RegularizationTerm.h"
       5             : 
       6             : #include <vector>
       7             : 
       8             : namespace elsa
       9             : {
      10             :     /**
      11             :      * @brief Class representing a generic optimization problem consisting of data term and
      12             :      * regularization term(s).
      13             :      *
      14             :      * @author Matthias Wieczorek - initial code
      15             :      * @author Maximilian Hornung - modularization
      16             :      * @author Tobias Lasser - rewrite, modernization
      17             :      *
      18             :      * @tparam data_t data type for the domain and range of the problem, defaulting to real_t
      19             :      *
      20             :      * This class represents a generic optimization problem, which consists of a data term and
      21             :      * (optionally) of one (or many) regularization terms, \f$ \argmin_x D(x) + \sum_{i=1}^n
      22             :      * \lambda_i R(x) \f$. Here, the data term \f$ D(x) \f$ is represented through a Functional (or
      23             :      * it derivatives), the regularization terms are represented by RegularizationTerms, which
      24             :      * encapsulate regularization parameters \f$ \lambda_i \f$ (scalar values) and the actual
      25             :      * regularization terms \f$ R(x) \f$ (Functionals or its derivatives).
      26             :      */
      27             :     template <typename data_t = real_t>
      28             :     class Problem : public Cloneable<Problem<data_t>>
      29             :     {
      30             :     public:
      31             :         /**
      32             :          * @brief Constructor for optimization problem, accepting a data and multiple regularization
      33             :          * terms, and an initial guess x0.
      34             :          *
      35             :          * @param[in] dataTerm functional expressing the data term
      36             :          * @param[in] regTerms vector of RegularizationTerms (weight and functional)
      37             :          * @param[in] x0 initial value for the current estimated solution
      38             :          * @param[in] lipschitzConstant if non-null the known lipschitz constant of the
      39             :          * problem. If null the lipschitz constant will be computed using power-iteration. Useful in
      40             :          * cases where the numerical approximation is not accurate and the constant is known.
      41             :          */
      42             :         Problem(const Functional<data_t>& dataTerm,
      43             :                 const std::vector<RegularizationTerm<data_t>>& regTerms,
      44             :                 const DataContainer<data_t>& x0, std::optional<data_t> lipschitzConstant = {});
      45             : 
      46             :         /**
      47             :          * @brief Constructor for optimization problem, accepting a data and multiple regularization
      48             :          * terms.
      49             :          *
      50             :          * @param[in] dataTerm functional expressing the data term
      51             :          * @param[in] regTerms vector of RegularizationTerms (weight and functional)
      52             :          * @param[in] lipschitzConstant if non-null the known lipschitz constant of the
      53             :          * problem. If null the lipschitz constant will be computed using power-iteration. Useful in
      54             :          * cases where the numerical approximation is not accurate and the constant is known.
      55             :          */
      56             :         Problem(const Functional<data_t>& dataTerm,
      57             :                 const std::vector<RegularizationTerm<data_t>>& regTerms,
      58             :                 std::optional<data_t> lipschitzConstant = {});
      59             : 
      60             :         /**
      61             :          * @brief Constructor for optimization problem, accepting a data and one regularization
      62             :          * term, and an initial guess x0.
      63             :          *
      64             :          * @param[in] dataTerm functional expressing the data term
      65             :          * @param[in] regTerm RegularizationTerm (weight and functional)
      66             :          * @param[in] x0 initial value for the current estimated solution
      67             :          * @param[in] lipschitzConstant if non-null the known lipschitz constant of the
      68             :          * problem. If null the lipschitz constant will be computed using power-iteration. Useful in
      69             :          * cases where the numerical approximation is not accurate and the constant is known.
      70             :          */
      71             :         Problem(const Functional<data_t>& dataTerm, const RegularizationTerm<data_t>& regTerm,
      72             :                 const DataContainer<data_t>& x0, std::optional<data_t> lipschitzConstant = {});
      73             : 
      74             :         /**
      75             :          * @brief Constructor for optimization problem, accepting a data and one regularization
      76             :          * term.
      77             :          *
      78             :          * @param[in] dataTerm functional expressing the data term
      79             :          * @param[in] regTerm RegularizationTerm (weight and functional)
      80             :          * @param[in] lipschitzConstant if non-null the known lipschitz constant of the
      81             :          * problem. If null the lipschitz constant will be computed using power-iteration. Useful in
      82             :          * cases where the numerical approximation is not accurate and the constant is known.
      83             :          */
      84             :         Problem(const Functional<data_t>& dataTerm, const RegularizationTerm<data_t>& regTerm,
      85             :                 std::optional<data_t> lipschitzConstant = {});
      86             : 
      87             :         /**
      88             :          * @brief Constructor for optimization problem, accepting a data term and an initial guess
      89             :          * x0.
      90             :          *
      91             :          * @param[in] dataTerm functional expressing the data term
      92             :          * @param[in] x0 initial value for the current estimated solution
      93             :          * @param[in] lipschitzConstant if non-null the known lipschitz constant of the
      94             :          * problem. If null the lipschitz constant will be computed using power-iteration. Useful in
      95             :          * cases where the numerical approximation is not accurate and the constant is known.
      96             :          */
      97             :         Problem(const Functional<data_t>& dataTerm, const DataContainer<data_t>& x0,
      98             :                 std::optional<data_t> lipschitzConstant = {});
      99             : 
     100             :         /**
     101             :          * @brief Constructor for optimization problem, accepting a data term.
     102             :          *
     103             :          * @param[in] dataTerm functional expressing the data term
     104             :          * @param[in] lipschitzConstant if non-null the known lipschitz constant of the
     105             :          * problem. If null the lipschitz constant will be computed using power-iteration. Useful in
     106             :          * cases where the numerical approximation is not accurate and the constant is known.
     107             :          */
     108             :         explicit Problem(const Functional<data_t>& dataTerm,
     109             :                          std::optional<data_t> lipschitzConstant = {});
     110             : 
     111             :         /// default destructor
     112           0 :         ~Problem() override = default;
     113             : 
     114             :         /// return the data term
     115             :         const Functional<data_t>& getDataTerm() const;
     116             : 
     117             :         /// return the vector of regularization terms
     118             :         const std::vector<RegularizationTerm<data_t>>& getRegularizationTerms() const;
     119             : 
     120             :         /// return the current estimated solution (const version)
     121             :         const DataContainer<data_t>& getCurrentSolution() const;
     122             : 
     123             :         /// return the current estimated solution
     124             :         DataContainer<data_t>& getCurrentSolution();
     125             : 
     126             :         /**
     127             :          * @brief evaluate the problem at the current estimated solution
     128             :          *
     129             :          * @returns the value of the problem evaluated at the current estimated solution
     130             :          *
     131             :          * Please note: this method calls the method evaluateImpl that has to be overridden in
     132             :          * derived classes.
     133             :          */
     134             :         data_t evaluate();
     135             : 
     136             :         /**
     137             :          * @brief return the gradient of the problem at the current estimated solution
     138             :          *
     139             :          * @returns DataContainer (in the domain of the problem) containing the result of the
     140             :          * gradient at the current solution
     141             :          *
     142             :          * Please note: this method used getGradient(result) to perform the actual operation.
     143             :          */
     144             :         DataContainer<data_t> getGradient();
     145             : 
     146             :         /**
     147             :          * @brief compute the gradient of the problem at the current estimated solution
     148             :          *
     149             :          * @param[out] result output DataContainer containing the gradient (in the domain of the
     150             :          * problem)
     151             :          *
     152             :          * Please note: this method calls the method getGradientImpl that has to be overridden in
     153             :          * derived classes.
     154             :          */
     155             :         void getGradient(DataContainer<data_t>& result);
     156             : 
     157             :         /**
     158             :          * @brief return the Hessian of the problem at the current estimated solution
     159             :          *
     160             :          * @returns a LinearOperator (the Hessian)
     161             :          *
     162             :          * Please note: this method calls the method getHessianImpl that has to be overridden in
     163             :          * derived classes.
     164             :          */
     165             :         LinearOperator<data_t> getHessian() const;
     166             : 
     167             :         /**
     168             :          * @brief return the Lipschitz Constant of the problem at the current estimated solution. If
     169             :          * an explicit lipschitz constant has been passed to the problem it will be returned here.
     170             :          *
     171             :          * @param[in] nIterations number of iterations to compute the lipschitz constant using
     172             :          * power iteration.
     173             :          *
     174             :          * @returns data_t (the Lipschitz Constant)
     175             :          *
     176             :          * Please note: this method calls the method getLipschitzConstantImpl that has to be
     177             :          * overridden in derived classes which want to provide a more specific way of computing
     178             :          * the Lipschitz constant, e.g. by not using power iteration or where the hessian is already
     179             :          * approximated as a diagonal matrix.
     180             :          */
     181             :         data_t getLipschitzConstant(index_t nIterations = 5) const;
     182             : 
     183             :     protected:
     184             :         /// the data term
     185             :         std::unique_ptr<Functional<data_t>> _dataTerm{};
     186             : 
     187             :         /// the regularization terms
     188             :         std::vector<RegularizationTerm<data_t>> _regTerms{};
     189             : 
     190             :         /// the current estimated solution
     191             :         DataContainer<data_t> _currentSolution;
     192             : 
     193             :         /// the known lipschitz constant for a problem, if not given will be computed on demand
     194             :         std::optional<data_t> _lipschitzConstant = {};
     195             : 
     196             :         /// protected copy constructor, simplifies cloning (of the subclasses primarily)
     197             :         Problem(const Problem<data_t>& problem);
     198             : 
     199             :         /// the evaluation of the optimization problem
     200             :         virtual data_t evaluateImpl();
     201             : 
     202             :         /// the getGradient method for the optimization problem
     203             :         virtual void getGradientImpl(DataContainer<data_t>& result);
     204             : 
     205             :         /// the getHessian method for the optimization problem
     206             :         virtual LinearOperator<data_t> getHessianImpl() const;
     207             : 
     208             :         /// the getLipschitzConstant method for the optimization problem
     209             :         virtual data_t getLipschitzConstantImpl(index_t nIterations) const;
     210             : 
     211             :         /// implement the polymorphic clone operation
     212             :         Problem<data_t>* cloneImpl() const override;
     213             : 
     214             :         /// implement the polymorphic comparison operation
     215             :         bool isEqual(const Problem<data_t>& other) const override;
     216             :     };
     217             : } // namespace elsa