nmpc_cgmres
Public Types | Public Member Functions | Public Attributes | List of all members
nmpc_cgmres::Gmres Class Reference

GMRES method to solve a linear equation. More...

#include <Gmres.h>

Public Types

using AmulFunc = std::function< Eigen::VectorXd(const Eigen::Ref< const Eigen::VectorXd > &)>
 Type of function that returns x multiplied by A. More...
 

Public Member Functions

EIGEN_MAKE_ALIGNED_OPERATOR_NEW Gmres ()
 Constructor. More...
 
void solve (const Eigen::Ref< const Eigen::MatrixXd > &A, const Eigen::Ref< const Eigen::VectorXd > &b, Eigen::Ref< Eigen::VectorXd > x, int k_max=100, double eps=1e-10)
 Solve. More...
 
void solve (const AmulFunc &Amul_func, const Eigen::Ref< const Eigen::VectorXd > &b, Eigen::Ref< Eigen::VectorXd > x, int k_max=100, double eps=1e-10)
 Solve. More...
 

Public Attributes

bool make_triangular_ = true
 
bool apply_reorth_ = true
 
Eigen::MatrixXd H_
 
Eigen::VectorXd g_
 
std::vector< double > err_list_
 
std::vector< Eigen::VectorXd > basis_
 

Detailed Description

GMRES method to solve a linear equation.

See the following articles about the GMRES method:

Definition at line 21 of file Gmres.h.

Member Typedef Documentation

◆ AmulFunc

using nmpc_cgmres::Gmres::AmulFunc = std::function<Eigen::VectorXd(const Eigen::Ref<const Eigen::VectorXd> &)>

Type of function that returns x multiplied by A.

Definition at line 25 of file Gmres.h.

Constructor & Destructor Documentation

◆ Gmres()

EIGEN_MAKE_ALIGNED_OPERATOR_NEW nmpc_cgmres::Gmres::Gmres ( )
inline

Constructor.

Definition at line 31 of file Gmres.h.

Member Function Documentation

◆ solve() [1/2]

void nmpc_cgmres::Gmres::solve ( const AmulFunc Amul_func,
const Eigen::Ref< const Eigen::VectorXd > &  b,
Eigen::Ref< Eigen::VectorXd >  x,
int  k_max = 100,
double  eps = 1e-10 
)
inline

Solve.

Parameters
Amul_functhe function to return $ A * v $ where $ v $ is given
bthe vector of linear equation
xthe initial guess of solution, which is overwritten by the final solution
k_maxthe maximum number of GMRES iteration
epsthe required solution tolerance

Solve the linear equation: $ A x = b $.

Refer to the Algorithm 3.5.1. in [1] for the case that make_triangular_ is true. Refer to the Algorithm 3.4.2. in [1] for the case that make_triangular_ is false. [1] Kelley, Carl T. Iterative methods for linear and nonlinear equations. Society for Industrial and Applied Mathematics, 1995.

Definition at line 67 of file Gmres.h.

◆ solve() [2/2]

void nmpc_cgmres::Gmres::solve ( const Eigen::Ref< const Eigen::MatrixXd > &  A,
const Eigen::Ref< const Eigen::VectorXd > &  b,
Eigen::Ref< Eigen::VectorXd >  x,
int  k_max = 100,
double  eps = 1e-10 
)
inline

Solve.

Parameters
Athe matrix of linear equation
bthe vector of linear equation
xthe initial guess of solution, which is overwritten by the final solution
k_maxthe maximum number of GMRES iteration
epsthe required solution tolerance

Solve the linear equation: $ A x = b $.

Definition at line 42 of file Gmres.h.

Member Data Documentation

◆ apply_reorth_

bool nmpc_cgmres::Gmres::apply_reorth_ = true

Definition at line 196 of file Gmres.h.

◆ basis_

std::vector<Eigen::VectorXd> nmpc_cgmres::Gmres::basis_

Definition at line 203 of file Gmres.h.

◆ err_list_

std::vector<double> nmpc_cgmres::Gmres::err_list_

Definition at line 201 of file Gmres.h.

◆ g_

Eigen::VectorXd nmpc_cgmres::Gmres::g_

Definition at line 199 of file Gmres.h.

◆ H_

Eigen::MatrixXd nmpc_cgmres::Gmres::H_

Definition at line 198 of file Gmres.h.

◆ make_triangular_

bool nmpc_cgmres::Gmres::make_triangular_ = true

Definition at line 195 of file Gmres.h.


The documentation for this class was generated from the following file: