gsGradientDescent.h

Go to the documentation of this file.

 
#pragma once
 
#include <gsCore/gsLinearAlgebra.h>
#include <gsIO/gsOptionList.h>
#include <gsOptimizer/gsOptimizer.h>
#include <gsOptimizer/gsOptProblem.h>
#define Eigen gsEigen
#include "gdcpp.h"
#undef Eigen
//#include "lsqcpp.h"
 
namespace gismo
{
 
template<typename T>
struct gsGradientDescentObjective
{
    typedef gsEigen::Matrix<T, gsEigen::Dynamic, 1> Vector;
    // typedef gsEigen::Matrix<T, gsEigen::Dynamic, gsEigen::Dynamic> Matrix;
 
    gsGradientDescentObjective(gsOptProblem<T>* objective)
    :
    obj(objective)
    { }
 
    gsGradientDescentObjective()
    :
    obj(nullptr)
    { }
 
 
    T operator()(const Vector & vx, Vector & vfgrad) const
    {
        vfgrad.resize(obj->numDesignVars());
        gsAsConstVector<T> xvec(vx.data(),vx.size());
        gsAsVector<T> gvec(vfgrad.data(),vfgrad.size());
        obj->gradObj_into(xvec,gvec);
        return obj->evalObj(xvec);
    }
 
    gsOptProblem<T> * obj;
};
 
template<typename T = real_t,
         typename StepSize=gdc::BarzilaiBorwein<T>,
         typename Callback=gdc::NoCallback<T>,
         typename FiniteDifferences=gdc::CentralDifferences<T> >
class gsGradientDescent : public gsOptimizer<T>
{
    using Base = gsOptimizer<T>;
 
    typedef typename gdc::GradientDescent<T, gsGradientDescentObjective<T>, StepSize, Callback, FiniteDifferences>::Result Result;
 
public:
    gsGradientDescent(gsOptProblem<T> * problem)
    :
    Base(problem),
    m_solver()
    {
        this->defaultOptions();
        gsGradientDescentObjective<T> obj(m_op);
        m_solver.setObjective(obj);
    }
 
 
public:
    // const gsMatrix<T> & lambda() const { return m_lambda; }
 
    void minimize(const gsMatrix<T> &initialGuess)
    {
        m_result = m_solver.minimize(initialGuess);
        m_numIterations = m_result.iterations;
        m_finalObjective = m_result.fval;
        m_curDesign = m_result.xval;
    }
 
    Result result() { return m_result; };
 
protected:
    void defaultOptions()
    {
        Base::defaultOptions();
        m_options.addReal("MinGradientLength","Minimal gradient length",1e-9);
        m_options.addReal("MinStepLength","Minimal step length",1e-9);
    }
 
    void getOptions()
    {
        Base::getOptions();
        m_minGradientLength = m_options.getReal("MinGradientLength");
        m_minStepLength = m_options.getReal("MinStepLength");
 
        m_solver.setMaxIterations(m_maxIterations);
        m_solver.setMinGradientLength(m_minGradientLength);
        m_solver.setMinStepLength(m_minStepLength);
        m_solver.setVerbosity(m_verbose);
    }
 
public:
 
    using Base::solve;
    void solve(const gsMatrix<T> &initialGuess)
    {
        this->getOptions();
        this->minimize(initialGuess);
    }
 
protected:
    using Base::m_op;
    using Base::m_numIterations;
    using Base::m_finalObjective;
    using Base::m_curDesign;
    using Base::m_options;
    using Base::m_verbose;
    using Base::m_maxIterations;
 
    using Base::defaultOptions;
    using Base::getOptions;
 
    Result m_result;
 
protected:
    T m_minGradientLength;
    T m_minStepLength;
 
    gdc::GradientDescent<T, gsGradientDescentObjective<T>, StepSize, Callback, FiniteDifferences> m_solver;
 
};
 
 
// using gsGradientDescent = gdc::GradientDescent<T, Objective, StepSize, Callback, FiniteDifferences>;
 
} //namespace gismo