Broyden.h

Go to the documentation of this file.

/*
 * Thermal-FIST package
 * 
 * Copyright (c) 2018-2019 Volodymyr Vovchenko
 *
 * GNU General Public License (GPLv3 or later)
 */

#ifndef BROYDEN_H
#define BROYDEN_H
 
#include <cstddef>
#include <vector>
 
namespace thermalfist {

  class BroydenEquations
  {
  public:
    BroydenEquations() = default;

    virtual ~BroydenEquations(void) = default;

    virtual int Dimension() const { return m_N; }

    void SetDimension(int dim);

    virtual std::vector<double> Equations(const std::vector<double> &x) = 0;
 
  protected:
    int m_N;           
  };
 

  class BroydenJacobian
  {
  public:
    static const double EPS;
    BroydenJacobian(BroydenEquations *eqs = NULL) : m_Equations(eqs), m_dx(EPS) { }

    virtual ~BroydenJacobian(void) = default;

    virtual std::vector<double> Jacobian(const std::vector<double> &x);

    void SetDx(double dx) { m_dx = dx; }

    double CurrentDx() const { return m_dx; }
 
  private:
    BroydenEquations *m_Equations;  
    double           m_dx;
  };

  class Broyden
  {
  public:
    class BroydenSolutionCriterium {
    public:
      BroydenSolutionCriterium(double maximum_error = TOL) : m_MaximumError(maximum_error) { }
      
      virtual ~BroydenSolutionCriterium() = default;

      virtual bool IsSolved(const std::vector<double>& x, const std::vector<double>& f, const std::vector<double>& xdelta = std::vector<double>()) const;
    protected:
      double m_MaximumError;
    };
  public:
    static const double TOL;       
    static const int    MAX_ITERS; 

    Broyden(BroydenEquations *eqs = NULL, BroydenJacobian *jaco = NULL) : m_Equations(eqs), m_Jacobian(jaco), m_Iterations(0), m_MaxDifference(0.), m_UseNewton(false) { }
    
    virtual ~Broyden(void) = default;

    virtual std::vector<double> Solve(const std::vector<double> &x0, BroydenSolutionCriterium *solcrit = NULL, int max_iterations = MAX_ITERS);

    int Iterations() const { return m_Iterations; }

    int MaxIterations() const { return m_MaxIterations; }

    double MaxDifference() const { return m_MaxDifference; }

    void UseNewton(bool flag) { m_UseNewton = flag; }

    bool UseNewton() const { return m_UseNewton; }
 
  protected:
    BroydenEquations *m_Equations; 
    BroydenJacobian  *m_Jacobian;   
    int m_Iterations;              
    int m_MaxIterations;      
    double m_MaxDifference;         
    bool m_UseNewton;
  };
 
} // namespace thermalfist
 
#endif