NLPInterfacePack: C++ Interfaces and Implementation for Non-Linear Programs
Version of the Day
|
NLP interface class that adds gradient information for the objective function {abstract}. More...
#include <NLPInterfacePack_NLPObjGrad.hpp>
Classes | |
struct | ObjGradInfo |
Struct for gradient (objective), objective and constriants (pointers) More... | |
Protected Member Functions | |
const ObjGradInfo | obj_grad_info () const |
Return objective gradient and zero order information. More... | |
Protected Member Functions inherited from NLPInterfacePack::NLP | |
template<class T > | |
void | assert_ref_set (T *p, std::string info) const |
Assert referece has been set for a quanity. More... | |
virtual void | imp_calc_f (const Vector &x, bool newx, const ZeroOrderInfo &zero_order_info) const =0 |
Overridden to compute f(x) (and perhaps other quantities if set). More... | |
virtual void | imp_calc_c (const Vector &x, bool newx, const ZeroOrderInfo &zero_order_info) const =0 |
Overridden to compute c(x) and perhaps f(x) and/or h(x) (if multiple calculaiton = true). More... | |
virtual void | imp_calc_c_breve (const Vector &x, bool newx, const ZeroOrderInfo &zero_order_info_breve) const |
Overridden to compute c_breve(x_breve) and perhaps f(x) and/or h_breve(x_breve) More... | |
virtual void | imp_calc_h_breve (const Vector &x, bool newx, const ZeroOrderInfo &zero_order_info_breve) const |
Overridden to compute h_breve(x_breve) and perhaps f(x) and/or c_breve(x_breve). More... | |
Constructors | |
NLPObjGrad () | |
Initialize to no reference set to calculation quanities. More... | |
NLP initialization | |
void | initialize (bool test_setup) |
Initialize the NLP for its first use. More... | |
Information | |
virtual bool | supports_Gf () const |
Determine if the objective gradient is supported or not. More... | |
virtual bool | supports_Gf_prod () const |
Determine if the objective gradient product is supported or not. More... | |
<<std aggr>> members for the gradient of the objective function Gf(x) | |
virtual void | set_Gf (VectorMutable *Gf) |
Set a pointer to a vector to be updated when this->calc_Gf() is called. More... | |
virtual VectorMutable * | get_Gf () |
Return pointer passed to this->set_Gf() . More... | |
virtual VectorMutable & | Gf () |
Returns non-const *this->get_Gf() . More... | |
virtual const Vector & | Gf () const |
Returns const *this->get_Gf() . More... | |
Unset calculation quantities | |
void | unset_quantities () |
Call to unset all storage quantities (both in this class and all subclasses). More... | |
Calculation Members | |
virtual void | calc_Gf (const Vector &x, bool newx=true) const |
Update the vector for Gf at the point x and put it in the stored reference. More... | |
virtual value_type | calc_Gf_prod (const Vector &x, const Vector &d, bool newx=true) const |
Calculate the inner product Gf(x)'*d at the point x and put it in the stored reference. More... | |
Function evaluation counts. | |
virtual size_type | num_Gf_evals () const |
Objective gradient evaluations count. More... | |
Protected methods to be overridden by subclasses | |
virtual void | imp_calc_Gf (const Vector &x, bool newx, const ObjGradInfo &obj_grad_info) const =0 |
Overridden to compute f(x) and perhaps c(x) (if multiple calculaiton = true). More... | |
Additional Inherited Members | |
Public Types inherited from NLPInterfacePack::NLP | |
typedef Teuchos::RCP< const VectorSpace > | vec_space_ptr_t |
typedef Teuchos::RCP< const OptionsFromStreamPack::OptionsFromStream > | options_ptr_t |
Public Member Functions inherited from NLPInterfacePack::NLP | |
const ZeroOrderInfo | zero_order_info () const |
Return pointer to set quantities. More... | |
const ZeroOrderInfo | zero_order_info_breve () const |
Return pointer to set hat quantities. More... | |
NLP () | |
Initialize to no reference set to calculation quanities. More... | |
virtual | ~NLP () |
Destructor that cleans all the memory it owns. More... | |
virtual void | force_xinit_in_bounds (bool force_xinit_in_bounds)=0 |
Set if the initial point must be within the bounds. More... | |
virtual bool | force_xinit_in_bounds () const =0 |
Returns if the initial point must be within the bounds. More... | |
virtual void | set_options (const options_ptr_t &options) |
Set the options that this NLP may be interested in. More... | |
virtual const options_ptr_t & | get_options () const |
Get the OptionsFromStream object being used to extract the options from. More... | |
virtual bool | is_initialized () const =0 |
Return if this is initialized. More... | |
virtual size_type | n () const |
Return the number of variables. More... | |
virtual size_type | m () const |
Return the number of general equality constraints. More... | |
virtual vec_space_ptr_t | space_x () const =0 |
Vector space object for unknown variables x (dimension n). More... | |
virtual vec_space_ptr_t | space_c () const =0 |
Vector space object for general equality constraints c(x) (dimension m). More... | |
virtual size_type | num_bounded_x () const =0 |
Returns the number of variables in x(i) for which xl(i)> -infinite_bound() or xu(i) < +infinite_bound() . More... | |
virtual const Vector & | xl () const =0 |
Returns the lower bounds on the variables x . More... | |
virtual const Vector & | xu () const =0 |
Returns a reference to the vector of upper bounds on the variables x . More... | |
virtual value_type | max_var_bounds_viol () const =0 |
Set the maximum absolute value for which the variable bounds may be violated by when computing function and gradient values. More... | |
virtual const Vector & | xinit () const =0 |
Returns a reference to the vector of the initial guess for the solution x . More... | |
virtual void | get_init_lagrange_mult (VectorMutable *lambda, VectorMutable *nu) const |
Get the initial value of the Lagrange multipliers lambda. More... | |
virtual void | set_f (value_type *f) |
Set a pointer to an value to be updated when this->calc_f() is called. More... | |
virtual value_type * | get_f () |
Return pointer passed to this->set_f() . More... | |
virtual value_type & | f () |
Returns non-const *this->get_f() . More... | |
virtual const value_type & | f () const |
Returns const *this->get_f() . More... | |
virtual void | set_c (VectorMutable *c) |
Set a pointer to a vector to be updated when this->calc_c() is called. More... | |
virtual VectorMutable * | get_c () |
Return pointer passed to this->set_c() . More... | |
virtual VectorMutable & | c () |
Returns non-const *this->get_c() . More... | |
virtual const Vector & | c () const |
Returns const *this->get_c() . More... | |
virtual void | scale_f (value_type scale_f)=0 |
Set the scaling of the objective function. More... | |
virtual value_type | scale_f () const =0 |
Return the scaling being used for the objective function. More... | |
virtual void | calc_f (const Vector &x, bool newx=true) const |
Update the value for the objective f at the point x and put it in the stored reference. More... | |
virtual void | calc_c (const Vector &x, bool newx=true) const |
Update the constraint residual vector for c at the point x and put it in the stored reference. More... | |
virtual void | report_final_solution (const Vector &x, const Vector *lambda, const Vector *nu, bool is_optimal) |
Used by the solver to report the final solution and multipliers. More... | |
virtual size_type | num_f_evals () const |
Gives the number of object function f(x) evaluations called by the solver since initialize() was called. More... | |
virtual size_type | num_c_evals () const |
Gives the number of constraint function c(x) evaluations called by the solver since initialize() was called. Throws exception if this->m() == 0 . More... | |
virtual size_type | ns () const |
Return the number of slack variables (i.e. number of general inequalities). More... | |
virtual vec_space_ptr_t | space_c_breve () const |
Vector space object for the original equalities c_breve(x_breve) More... | |
virtual vec_space_ptr_t | space_h_breve () const |
Vector space object for the original inequalities h_breve(x_breve) More... | |
virtual const Vector & | hl_breve () const |
Returns a reference to the vector of lower bounds on the general inequality constraints h_breve(x_breve) . More... | |
virtual const Vector & | hu_breve () const |
Returns a reference to the vector of upper bounds on the general inequality constraints h_breve(x_breve) . More... | |
virtual void | set_c_breve (VectorMutable *c_breve) |
Set a pointer to a vector to be updated when this->calc_c_breve() is called. More... | |
virtual VectorMutable * | get_c_breve () |
Return pointer passed to this->set_c_breve() . More... | |
virtual VectorMutable & | c_breve () |
Returns non-const *this->get_c_breve() . More... | |
virtual const Vector & | c_breve () const |
Returns const *this->get_c_breve() . More... | |
virtual void | set_h_breve (VectorMutable *h_breve) |
Set a pointer to a vector to be updated when this->calc_h_breve() is called. More... | |
virtual VectorMutable * | get_h_breve () |
Return pointer passed to this->set_h_breve() . More... | |
virtual VectorMutable & | h_breve () |
Returns non-const *this->get_h_breve() . More... | |
virtual const Vector & | h_breve () const |
Returns const *this->get_h_breve() . More... | |
virtual const Permutation & | P_var () const |
Return the permutation object for the variables. More... | |
virtual const Permutation & | P_equ () const |
Return the permutation object for the constraints. More... | |
virtual void | calc_c_breve (const Vector &x, bool newx=true) const |
Update the constraint residual vector for c_breve at the point x and put it in the stored reference. More... | |
virtual void | calc_h_breve (const Vector &x, bool newx=true) const |
Update the constraint residual vector for h_breve at the point x and put it in the stored reference. More... | |
Static Public Member Functions inherited from NLPInterfacePack::NLP | |
static value_type | infinite_bound () |
Value for an infinite bound. More... | |
NLP interface class that adds gradient information for the objective function {abstract}.
Overview:
This class adds the ability to compute the gradient of the objective function Gf(x)
to the basic information given in the NLP
interface class. Note that Gf
is in the vector space space_x()
.
Client Usage:
As with the NLP
base interface, the initialize()
method must be called before the NLP object can be used. The method set_Gf()
is used to set a pointer to a vector to update when the gradient of the objective Gf
is computed when calc_Gf()
is called.
The number of evaluations of Gf
using calc_Gf()
is returned by num_Gf_evals()
.
Subclass developer's notes:
In addition to the methods that must be overridden by the NLP
interface (see) the following methods must also be overridden: imp_calc_Gf()
.
In addition to the methods that should be overridden from NLP
by most subclasses (see), the following additional methods should be overridden: initialize()
.
The following methods should never have to be overridden by most subclasses except in some very strange situations: set_Gf()
, get_Gf()
, Gf()
, num_Gf_evals()
.
Definition at line 79 of file NLPInterfacePack_NLPObjGrad.hpp.
NLPInterfacePack::NLPObjGrad::NLPObjGrad | ( | ) |
Initialize to no reference set to calculation quanities.
Definition at line 53 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Initialize the NLP for its first use.
This function implementation should be called by subclass implementations in order to reset counts for f(x)
, c(x)
, h(x)
and Gf(x)
evaluations. This implementation calls this->NLP::initialize()
Postconditions:
NLP::initialize()
this->num_Gf_evals() == 0
Reimplemented from NLPInterfacePack::NLP.
Reimplemented in NLPInterfacePack::NLPSerialPreprocess, NLPInterfacePack::NLPSerialPreprocessExplJac, NLPInterfacePack::NLPSecondOrder, and NLPInterfacePack::NLPBarrier.
Definition at line 57 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Determine if the objective gradient is supported or not.
The default implementation returns true
.
Definition at line 64 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Determine if the objective gradient product is supported or not.
The default implementation returns true
.
Definition at line 69 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Set a pointer to a vector to be updated when this->calc_Gf()
is called.
Gf | [in] Pointer to gradient vector. May be NULL . |
Preconditions:
this->is_initialized() == true
(throw NotInitialized
) this->supports_Gf()
Gf != NULL
] Gf->space().is_compatible(*this->space_x()) == true
(throw VectorSpace::IncompatibleVectorSpaces
) Postconditions:
this->get_Gf() == Gf
Reimplemented in NLPInterfacePack::NLPBarrier.
Definition at line 76 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Return pointer passed to this->set_Gf()
.
Preconditions:
this->is_initialized() == true
(throw NotInitialized
) this->supports_Gf()
Reimplemented in NLPInterfacePack::NLPBarrier.
Definition at line 81 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Returns non-const
*this->get_Gf()
.
Preconditions:
this->is_initialized() == true
(throw NotInitialized
) this->supports_Gf()
this->get_Gf() != NULL
(throw NoRefSet
) Reimplemented in NLPInterfacePack::NLPBarrier.
Definition at line 86 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Returns const
*this->get_Gf()
.
Preconditions:
this->is_initialized() == true
(throw NotInitialized
) this->supports_Gf()
this->get_Gf() != NULL
(throw NoRefSet
) Reimplemented in NLPInterfacePack::NLPBarrier.
Definition at line 91 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Call to unset all storage quantities (both in this class and all subclasses).
Preconditions:
this->is_initialized() == true
(throw NotInitialized
) Postconditions:
NLP::unset_quantities()
this->get_Gf() == NULL
This method must be called by all subclasses that override it.
Reimplemented from NLPInterfacePack::NLP.
Reimplemented in NLPInterfacePack::NLPSecondOrder.
Definition at line 96 of file NLPInterfacePack_NLPObjGrad.cpp.
Update the vector for Gf
at the point x
and put it in the stored reference.
x | [in] Point at which to calculate the gradient of the objective Gf(x) . |
newx | [in] (default true ) If false , the values in x are assumed to be the same as the last call to a this->calc_*(x,newx) member. If true , the values in x are assumed to not be the same as the last call to a this->calc_*(x,newx) member. |
Preconditions:
this->is_initialized() == true
(throw NotInitialized
) this->supports_Gf()
x.space().is_compatible(*this->space_x()) == true
(throw VectorSpace::IncompatibleVectorSpaces
) this->get_Gf() != NULL
(throw NoRefSet
) Postconditions:
this->Gf()
is updated to Gf(x)
If set_multi_calc(true)
was called then referenced storage for f
and/or c
may also be updated but are not guaranteed to be. But no other quanities from possible subclasses are allowed to be updated as a side effect (i.e. no higher order derivatives).
Reimplemented in NLPInterfacePack::NLPBarrier.
Definition at line 104 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Calculate the inner product Gf(x)'*d
at the point x
and put it in the stored reference.
x | [in] Base point |
d | [in] Direction to compute the product along. |
newx | [in] (default true ) If false , the values in x are assumed to be the same as the last call to a this->calc_*(x,newx) member. If true , the values in x are assumed to not be the same as the last call to a this->calc_*(x,newx) member. |
Preconditions:
this->is_initialized() == true
(throw NotInitialized
) this->supports_Gf()
x.space().is_compatible(*this->space_x()) == true
(throw VectorSpace::IncompatibleVectorSpaces
) Postconditions:
return
gives the desired product. If set_multi_calc(true)
was called then referenced storage for f
and/or c
may also be updated but are not guaranteed to be. But no other quanities from possible subclasses are allowed to be updated as a side effect (i.e. no higher order derivatives).
Definition at line 111 of file NLPInterfacePack_NLPObjGrad.cpp.
|
virtual |
Objective gradient evaluations count.
This function can be called to find out how many evaluations this->calc_Gf()
the client requested since this->initialize()
was called.
Reimplemented in NLPInterfacePack::NLPBarrier.
Definition at line 128 of file NLPInterfacePack_NLPObjGrad.cpp.
|
inlineprotected |
Return objective gradient and zero order information.
Definition at line 332 of file NLPInterfacePack_NLPObjGrad.hpp.
|
protectedpure virtual |
Overridden to compute f(x) and perhaps c(x) (if multiple calculaiton = true).
Preconditions:
x.space().is_compatible(*this->space_x())
(throw IncompatibleType
) obj_grad_info.Gf != NULL
(throw std::invalid_argument
) Postconditions:
*obj_grad_info.Gf
is updated to Gf(x). x | [in] Unknown vector (size n). |
newx | [in] (default true ) If false , the values in x are assumed to be the same as the last call to a this->imp_calc_*(x,newx) member. If true , the values in x are assumed to not be the same as the last call to a this->imp_calc_*(x,newx) member. |
obj_grad_info | [out] Pointers to f , c and Gf . On output *obj_grad_info.Gf is updated to Gf(x). Any of the other objects pointed to in obj_grad_info may be set if this->multi_calc() == true but are now guaranteed to be. |
Implemented in NLPInterfacePack::NLPSerialPreprocess, and NLPInterfacePack::NLPBarrier.