Crocoddyl
 
Loading...
Searching...
No Matches
SolverIntro Class Reference
Inheritance diagram for SolverIntro:
SolverFDDP

Public Member Functions

EIGEN_MAKE_ALIGNED_OPERATOR_NEW SolverIntro (std::shared_ptr< ShootingProblem > problem)
 Initialize the INTRO solver.
 
virtual double calcDiff ()
 
virtual void computeGains (const std::size_t t)
 
virtual void computeValueFunction (const std::size_t t, const std::shared_ptr< ActionModelAbstract > &model)
 
EqualitySolverType get_equality_solver () const
 Return the type of solver used for handling the equality constraints.
 
const std::vector< Eigen::MatrixXd > & get_Hy () const
 Return span-projected Jacobian of the equality-constraint with respect to the control.
 
const std::vector< Eigen::VectorXd > & get_ks () const
 Return feedforward term related to the equality constraints.
 
const std::vector< Eigen::MatrixXd > & get_Ks () const
 Return feedback gain related to the equality constraints.
 
const std::vector< Eigen::VectorXd > & get_kz () const
 Return feedforward term related to the nullspace of \(\mathbf{H_u}\).
 
const std::vector< Eigen::MatrixXd > & get_Kz () const
 Return feedback gain related to the nullspace of \(\mathbf{H_u}\).
 
const std::vector< Eigen::MatrixXd > & get_Quz () const
 Return Hessian of the reduced Hamiltonian \(\mathbf{Q_{uz}}\).
 
const std::vector< Eigen::MatrixXd > & get_Qxz () const
 Return Hessian of the reduced Hamiltonian \(\mathbf{Q_{xz}}\).
 
const std::vector< Eigen::VectorXd > & get_Qz () const
 Return Jacobian of the reduced Hamiltonian \(\mathbf{Q_{z}}\).
 
const std::vector< Eigen::MatrixXd > & get_Qzz () const
 Return Hessian of the reduced Hamiltonian \(\mathbf{Q_{zz}}\).
 
double get_rho () const
 Return the rho parameter used in the merit function.
 
double get_th_feas () const
 Return the threshold for switching to feasibility.
 
double get_upsilon () const
 Return the estimated penalty parameter that balances relative contribution of the cost function and equality constraints.
 
const std::vector< Eigen::MatrixXd > & get_YZ () const
 Return the rank of control-equality constraints \(\mathbf{H_u}\f */ const std::vector<std::size_t>& get_Hu_rank() const; /** @brief Return the span and kernel of control-equality constraints \)\mathbf{H_u}\f.
 
bool get_zero_upsilon () const
 Return the zero-upsilon label.
 
virtual void resizeData ()
 
void set_equality_solver (const EqualitySolverType type)
 Modify the type of solver used for handling the equality constraints.
 
void set_rho (const double rho)
 Modify the rho parameter used in the merit function.
 
void set_th_feas (const double th_feas)
 Modify the threshold for switching to feasibility.
 
void set_zero_upsilon (const bool zero_upsilon)
 Modify the zero-upsilon label.
 
virtual bool solve (const std::vector< Eigen::VectorXd > &init_xs=DEFAULT_VECTOR, const std::vector< Eigen::VectorXd > &init_us=DEFAULT_VECTOR, const std::size_t maxiter=100, const bool is_feasible=false, const double init_reg=NAN)
 
virtual double stoppingCriteria ()
 
virtual double tryStep (const double step_length=1)
 
- Public Member Functions inherited from SolverFDDP
EIGEN_MAKE_ALIGNED_OPERATOR_NEW SolverFDDP (std::shared_ptr< ShootingProblem > problem)
 Initialize the FDDP solver.
 
virtual const Eigen::Vector2d & expectedImprovement ()
 Return the expected improvement \(dV_{exp}\) from a given current search direction \((\delta\mathbf{x}^k,\delta\mathbf{u}^k)\).
 
virtual void forwardPass (const double stepLength)
 
double get_th_acceptnegstep () const
 Return the threshold used for accepting step along ascent direction.
 
void set_th_acceptnegstep (const double th_acceptnegstep)
 Modify the threshold used for accepting step along ascent direction.
 
void updateExpectedImprovement ()
 Update internal values for computing the expected improvement.
 

Protected Attributes

enum EqualitySolverType eq_solver_
 Strategy used for handling the equality constraints.
 
std::vector< Eigen::FullPivLU< Eigen::MatrixXd > > Hu_lu_
 
std::vector< Eigen::ColPivHouseholderQR< Eigen::MatrixXd > > Hu_qr_
 
std::vector< std::size_t > Hu_rank_
 Rank of the control Jacobian of the equality constraints.
 
std::vector< Eigen::MatrixXd > Hy_
 
std::vector< Eigen::PartialPivLU< Eigen::MatrixXd > > Hy_lu_
 
std::vector< Eigen::MatrixXd > KQuu_tmp_
 
std::vector< Eigen::VectorXd > ks_
 Feedforward term related to the equality constraints.
 
std::vector< Eigen::MatrixXd > Ks_
 Feedback gain related to the equality constraints.
 
std::vector< Eigen::VectorXd > kz_
 Feedforward term in the nullspace of \(\mathbf{H_u}\).
 
std::vector< Eigen::MatrixXd > Kz_
 Feedback gain in the nullspace of \(\mathbf{H_u}\).
 
std::vector< Eigen::MatrixXd > QuuinvHuT_
 
std::vector< Eigen::MatrixXd > Quz_
 Hessian of the reduced Hamiltonian \(\mathbf{Q_{uz}}\).
 
std::vector< Eigen::MatrixXd > Qxz_
 Hessian of the reduced Hamiltonian \(\mathbf{Q_{xz}}\).
 
std::vector< Eigen::VectorXd > Qz_
 Jacobian of the reduced Hamiltonian \(\mathbf{Q_{z}}\).
 
std::vector< Eigen::MatrixXd > Qzz_
 Hessian of the reduced Hamiltonian \(\mathbf{Q_{zz}}\).
 
std::vector< Eigen::LLT< Eigen::MatrixXd > > Qzz_llt_
 Cholesky LLT solver.
 
double rho_
 
double th_feas_
 Threshold for switching to feasibility.
 
double upsilon_
 
std::vector< Eigen::MatrixXd > YZ_
 
bool zero_upsilon_
 
- Protected Attributes inherited from SolverFDDP
double dg_
 Internal data for computing the expected improvement.
 
double dq_
 Internal data for computing the expected improvement.
 
double dv_
 Internal data for computing the expected improvement.
 
double th_acceptnegstep_
 

Detailed Description

Definition at line 18 of file intro.hpp.

Constructor & Destructor Documentation

◆ SolverIntro()

SolverIntro ( std::shared_ptr< ShootingProblem problem)
explicit

Initialize the INTRO solver.

Parameters
[in]problemShooting problem
[in]reducedUse the reduced Schur-complement approach (default true)

Definition at line 15 of file intro.cpp.

◆ ~SolverIntro()

~SolverIntro ( )
virtual

Definition at line 68 of file intro.cpp.

Member Function Documentation

◆ solve()

bool solve ( const std::vector< Eigen::VectorXd > &  init_xs = DEFAULT_VECTOR,
const std::vector< Eigen::VectorXd > &  init_us = DEFAULT_VECTOR,
const std::size_t  maxiter = 100,
const bool  is_feasible = false,
const double  init_reg = NAN 
)
virtual

Reimplemented from SolverFDDP.

Definition at line 70 of file intro.cpp.

◆ tryStep()

double tryStep ( const double  step_length = 1)
virtual

Definition at line 187 of file intro.cpp.

◆ stoppingCriteria()

double stoppingCriteria ( )
virtual

Definition at line 193 of file intro.cpp.

◆ resizeData()

void resizeData ( )
virtual

Definition at line 198 of file intro.cpp.

◆ calcDiff()

double calcDiff ( )
virtual

Definition at line 226 of file intro.cpp.

◆ computeValueFunction()

void computeValueFunction ( const std::size_t  t,
const std::shared_ptr< ActionModelAbstract > &  model 
)
virtual

Definition at line 297 of file intro.cpp.

◆ computeGains()

void computeGains ( const std::size_t  t)
virtual

Definition at line 329 of file intro.cpp.

◆ get_equality_solver()

EqualitySolverType get_equality_solver ( ) const

Return the type of solver used for handling the equality constraints.

Definition at line 401 of file intro.cpp.

◆ get_th_feas()

double get_th_feas ( ) const

Return the threshold for switching to feasibility.

Definition at line 405 of file intro.cpp.

◆ get_rho()

double get_rho ( ) const

Return the rho parameter used in the merit function.

Definition at line 407 of file intro.cpp.

◆ get_upsilon()

double get_upsilon ( ) const

Return the estimated penalty parameter that balances relative contribution of the cost function and equality constraints.

Definition at line 409 of file intro.cpp.

◆ get_YZ()

const std::vector< Eigen::MatrixXd > & get_YZ ( ) const

Return the rank of control-equality constraints \(\mathbf{H_u}\f */ const std::vector<std::size_t>& get_Hu_rank() const; /** @brief Return the span and kernel of control-equality constraints \)\mathbf{H_u}\f.

Definition at line 417 of file intro.cpp.

◆ get_Qzz()

const std::vector< Eigen::MatrixXd > & get_Qzz ( ) const

Return Hessian of the reduced Hamiltonian \(\mathbf{Q_{zz}}\).

Definition at line 419 of file intro.cpp.

◆ get_Qxz()

const std::vector< Eigen::MatrixXd > & get_Qxz ( ) const

Return Hessian of the reduced Hamiltonian \(\mathbf{Q_{xz}}\).

Definition at line 423 of file intro.cpp.

◆ get_Quz()

const std::vector< Eigen::MatrixXd > & get_Quz ( ) const

Return Hessian of the reduced Hamiltonian \(\mathbf{Q_{uz}}\).

Definition at line 427 of file intro.cpp.

◆ get_Qz()

const std::vector< Eigen::VectorXd > & get_Qz ( ) const

Return Jacobian of the reduced Hamiltonian \(\mathbf{Q_{z}}\).

Definition at line 431 of file intro.cpp.

◆ get_Hy()

const std::vector< Eigen::MatrixXd > & get_Hy ( ) const

Return span-projected Jacobian of the equality-constraint with respect to the control.

Definition at line 433 of file intro.cpp.

◆ get_kz()

const std::vector< Eigen::VectorXd > & get_kz ( ) const

Return feedforward term related to the nullspace of \(\mathbf{H_u}\).

Definition at line 435 of file intro.cpp.

◆ get_Kz()

const std::vector< Eigen::MatrixXd > & get_Kz ( ) const

Return feedback gain related to the nullspace of \(\mathbf{H_u}\).

Definition at line 437 of file intro.cpp.

◆ get_ks()

const std::vector< Eigen::VectorXd > & get_ks ( ) const

Return feedforward term related to the equality constraints.

Definition at line 439 of file intro.cpp.

◆ get_Ks()

const std::vector< Eigen::MatrixXd > & get_Ks ( ) const

Return feedback gain related to the equality constraints.

Definition at line 441 of file intro.cpp.

◆ get_zero_upsilon()

bool get_zero_upsilon ( ) const

Return the zero-upsilon label.

True if we set the estimated penalty parameter (upsilon) to zero when solve is called.

Definition at line 411 of file intro.cpp.

◆ set_equality_solver()

void set_equality_solver ( const EqualitySolverType  type)

Modify the type of solver used for handling the equality constraints.

Note that the default solver is nullspace LU. When we enable parallelization, this strategy is generally faster than others for medium to large systems.

Definition at line 443 of file intro.cpp.

◆ set_th_feas()

void set_th_feas ( const double  th_feas)

Modify the threshold for switching to feasibility.

Definition at line 447 of file intro.cpp.

◆ set_rho()

void set_rho ( const double  rho)

Modify the rho parameter used in the merit function.

Definition at line 449 of file intro.cpp.

◆ set_zero_upsilon()

void set_zero_upsilon ( const bool  zero_upsilon)

Modify the zero-upsilon label.

Parameters
zero_upsilonTrue if we set estimated penalty parameter (upsilon) to zero when solve is called.

Definition at line 456 of file intro.cpp.

Member Data Documentation

◆ eq_solver_

enum EqualitySolverType eq_solver_
protected

Strategy used for handling the equality constraints.

Definition at line 161 of file intro.hpp.

◆ th_feas_

double th_feas_
protected

Threshold for switching to feasibility.

Definition at line 162 of file intro.hpp.

◆ rho_

double rho_
protected

Parameter used in the merit function to predict the expected reduction

Definition at line 163 of file intro.hpp.

◆ upsilon_

double upsilon_
protected

Estimated penalty parameter that balances relative contribution of the cost function and equality constraints

Definition at line 166 of file intro.hpp.

◆ zero_upsilon_

bool zero_upsilon_
protected

True if we wish to set estimated penalty parameter (upsilon) to zero when solve is called.

Definition at line 168 of file intro.hpp.

◆ Hu_rank_

std::vector<std::size_t> Hu_rank_
protected

Rank of the control Jacobian of the equality constraints.

Definition at line 172 of file intro.hpp.

◆ KQuu_tmp_

std::vector<Eigen::MatrixXd> KQuu_tmp_
protected

Definition at line 173 of file intro.hpp.

◆ YZ_

std::vector<Eigen::MatrixXd> YZ_
protected

Span \(\mathbf{Y}\in\mathbb{R}^{rank}\) and kernel \(\mathbf{Z}\in\mathbb{R}^{nullity}\) of the control-equality constraints \(\mathbf{H_u}\)

Definition at line 175 of file intro.hpp.

◆ Hy_

std::vector<Eigen::MatrixXd> Hy_
protected

Span-projected Jacobian of the equality-constraint with respect to the control

Definition at line 179 of file intro.hpp.

◆ Qz_

std::vector<Eigen::VectorXd> Qz_
protected

Jacobian of the reduced Hamiltonian \(\mathbf{Q_{z}}\).

Definition at line 182 of file intro.hpp.

◆ Qzz_

std::vector<Eigen::MatrixXd> Qzz_
protected

Hessian of the reduced Hamiltonian \(\mathbf{Q_{zz}}\).

Definition at line 184 of file intro.hpp.

◆ Qxz_

std::vector<Eigen::MatrixXd> Qxz_
protected

Hessian of the reduced Hamiltonian \(\mathbf{Q_{xz}}\).

Definition at line 186 of file intro.hpp.

◆ Quz_

std::vector<Eigen::MatrixXd> Quz_
protected

Hessian of the reduced Hamiltonian \(\mathbf{Q_{uz}}\).

Definition at line 188 of file intro.hpp.

◆ kz_

std::vector<Eigen::VectorXd> kz_
protected

Feedforward term in the nullspace of \(\mathbf{H_u}\).

Definition at line 190 of file intro.hpp.

◆ Kz_

std::vector<Eigen::MatrixXd> Kz_
protected

Feedback gain in the nullspace of \(\mathbf{H_u}\).

Definition at line 192 of file intro.hpp.

◆ ks_

std::vector<Eigen::VectorXd> ks_
protected

Feedforward term related to the equality constraints.

Definition at line 194 of file intro.hpp.

◆ Ks_

std::vector<Eigen::MatrixXd> Ks_
protected

Feedback gain related to the equality constraints.

Definition at line 196 of file intro.hpp.

◆ QuuinvHuT_

std::vector<Eigen::MatrixXd> QuuinvHuT_
protected

Definition at line 197 of file intro.hpp.

◆ Qzz_llt_

std::vector<Eigen::LLT<Eigen::MatrixXd> > Qzz_llt_
protected

Cholesky LLT solver.

Definition at line 198 of file intro.hpp.

◆ Hu_lu_

std::vector<Eigen::FullPivLU<Eigen::MatrixXd> > Hu_lu_
protected

Full-pivot LU solvers used for computing the span and nullspace matrices

Definition at line 200 of file intro.hpp.

◆ Hu_qr_

std::vector<Eigen::ColPivHouseholderQR<Eigen::MatrixXd> > Hu_qr_
protected

Column-pivot QR solvers used for computing the span and nullspace matrices

Definition at line 203 of file intro.hpp.

◆ Hy_lu_

std::vector<Eigen::PartialPivLU<Eigen::MatrixXd> > Hy_lu_
protected

Partial-pivot LU solvers used for computing the feedforward and feedback gain related to the equality constraint

Definition at line 206 of file intro.hpp.


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