Newton.h

Go to the documentation of this file.

 
#ifndef PBAT_MATH_OPTIMIZATION_NEWTON_H
#define PBAT_MATH_OPTIMIZATION_NEWTON_H
 
#include "LineSearch.h"
#include "pbat/Aliases.h"
 
#include <Eigen/Core>
#include <optional>
#include <type_traits>
 
namespace pbat::math::optimization {

template <class TScalar = Scalar>
struct Newton
{
    int nMaxIters; 
    TScalar gtol2; 
    Eigen::Vector<TScalar, Eigen::Dynamic> dxk; 
    Eigen::Vector<TScalar, Eigen::Dynamic> gk;  

    Newton(int nMaxIters = 10, TScalar gtol = TScalar(1e-4), Index n = 0);
    template <
        class FPrepareDerivatives,
        class FObjective,
        class FGradient,
        class FHessianInverseProduct,
        class TDerivedX>
    TScalar Solve(
        FPrepareDerivatives prepareDerivatives,
        FObjective f,
        FGradient g,
        FHessianInverseProduct Hinv,
        Eigen::MatrixBase<TDerivedX>& xk,
        std::optional<BackTrackingLineSearch<TScalar>> lineSearch = std::nullopt);
};
 
template <class TScalar>
inline Newton<TScalar>::Newton(int nMaxItersIn, TScalar gtol, Index n)
    : nMaxIters(nMaxItersIn), gtol2(gtol * gtol), dxk(n), gk(n)
{
}
 
template <class TScalar>
template <
    class FPrepareDerivatives,
    class FObjective,
    class FGradient,
    class FHessianInverseProduct,
    class TDerivedX>
inline TScalar Newton<TScalar>::Solve(
    FPrepareDerivatives prepareDerivatives,
    FObjective f,
    FGradient g,
    FHessianInverseProduct Hinv,
    Eigen::MatrixBase<TDerivedX>& xk,
    std::optional<BackTrackingLineSearch<TScalar>> lineSearch)
{
    TScalar gnorm2{0};
    prepareDerivatives(xk);
    gk = g(xk);
    for (auto k = 0; k < nMaxIters; ++k)
    {
        gnorm2 = gk.squaredNorm();
        if (gnorm2 < gtol2)
            break;
        dxk = -Hinv(xk, gk);
        if (lineSearch)
            lineSearch->Solve(f, gk, dxk, xk);
        else
            xk += dxk;
        prepareDerivatives(xk);
        gk = g(xk);
    }
    return gnorm2;
}
 
} // namespace pbat::math::optimization
 
#endif // PBAT_MATH_OPTIMIZATION_NEWTON_H