44 #ifndef ROL_AUGMENTEDLAGRANGIAN_H 45 #define ROL_AUGMENTEDLAGRANGIAN_H 51 #include "Teuchos_RCP.hpp" 67 Teuchos::RCP<Objective<Real> >
obj_;
68 Teuchos::RCP<EqualityConstraint<Real> >
con_;
69 Teuchos::RCP<Vector<Real> >
lam_;
70 Teuchos::RCP<Vector<Real> >
dlam_;
71 Teuchos::RCP<Vector<Real> >
x_;
72 Teuchos::RCP<Vector<Real> >
g_;
73 Teuchos::RCP<Vector<Real> >
c_;
74 Teuchos::RCP<Vector<Real> >
dc1_;
75 Teuchos::RCP<Vector<Real> >
dc2_;
96 Teuchos::ParameterList &parlist)
97 : fval_(0.0), ncval_(0), nfval_(0), ngval_(0), penaltyParameter_(mu),
98 isConstraintComputed_(false), isValueComputed_(false), isGradientComputed_(false) {
99 obj_ = Teuchos::rcp(&obj,
false);
100 con_ = Teuchos::rcp(&con,
false);
103 g_ = x.
dual().clone();
104 dc1_ = x.
dual().clone();
111 Teuchos::ParameterList& sublist = parlist.sublist(
"Step").sublist(
"Augmented Lagrangian");
112 scaleLagrangian_ = sublist.get(
"Use Scaled Augmented Lagrangian",
false);
113 HessianLevel_ = sublist.get(
"Level of Hessian Approximation", 0);
118 if ( !isValueComputed_ ) {
120 fval_ = obj_->value(x,tol);
122 isValueComputed_ =
true;
129 if ( !isGradientComputed_ ) {
131 obj_->gradient(*g_,x,tol);
133 isGradientComputed_ =
true;
140 if ( !isConstraintComputed_ ) {
142 con_->value(*c_,x,tol);
144 isConstraintComputed_ =
true;
162 ncval_ = 0.; nfval_ = 0.; ngval_ = 0.;
164 penaltyParameter_ = penaltyParameter;
175 obj_->update(x,flag,iter);
176 con_->update(x,flag,iter);
178 isConstraintComputed_ =
false;
179 isValueComputed_ =
false;
180 isGradientComputed_ =
false;
191 if ( !isConstraintComputed_ ) {
193 con_->value(*c_,x,tol);
195 isConstraintComputed_ =
true;
197 if ( !isValueComputed_ ) {
199 fval_ = obj_->value(x,tol);
201 isValueComputed_ =
true;
204 Real cval = lam_->dot(c_->dual());
206 Real pval = c_->dot(*c_);
209 if (scaleLagrangian_) {
210 val = (fval_ + cval)/penaltyParameter_ + 0.5*pval;
213 val = fval_ + cval + 0.5*penaltyParameter_*pval;
226 if ( !isConstraintComputed_ ) {
228 con_->value(*c_,x,tol);
230 isConstraintComputed_ =
true;
232 if ( !isGradientComputed_ ) {
234 obj_->gradient(*g_,x,tol);
236 isGradientComputed_ =
true;
240 dlam_->set(c_->dual());
241 if ( scaleLagrangian_ ) {
242 g.
scale(1./penaltyParameter_);
243 dlam_->axpy(1./penaltyParameter_,*lam_);
246 dlam_->scale(penaltyParameter_);
249 con_->applyAdjointJacobian(*dc1_,*dlam_,x,tol);
263 obj_->hessVec(hv,v,x,tol);
264 if (HessianLevel_ < 2) {
265 con_->applyJacobian(*dc2_,v,x,tol);
266 con_->applyAdjointJacobian(*dc1_,dc2_->dual(),x,tol);
267 if (scaleLagrangian_) {
268 hv.
scale(1./penaltyParameter_);
272 hv.
axpy(penaltyParameter_,*dc1_);
275 if (HessianLevel_ == 0) {
276 if ( !isConstraintComputed_ ) {
278 con_->value(*c_,x,tol);
280 isConstraintComputed_ =
true;
283 dlam_->set(c_->dual());
284 if ( scaleLagrangian_ ) {
285 dlam_->axpy(1./penaltyParameter_,*lam_);
288 dlam_->scale(penaltyParameter_);
291 con_->applyAdjointHessian(*dc1_,*dlam_,v,x,tol);
Provides the interface to evaluate objective functions.
Provides the interface to evaluate the augmented Lagrangian.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
virtual void scale(const Real alpha)=0
Compute where .
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update augmented Lagrangian function.
virtual void plus(const Vector &x)=0
Compute , where .
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Teuchos::RCP< Vector< Real > > g_
Contains definitions of custom data types in ROL.
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
Teuchos::RCP< Objective< Real > > obj_
bool isConstraintComputed_
Defines the linear algebra or vector space interface.
Teuchos::RCP< Vector< Real > > c_
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Real getObjectiveValue(const Vector< Real > &x)
Teuchos::RCP< Vector< Real > > dc2_
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Defines the equality constraint operator interface.
Teuchos::RCP< EqualityConstraint< Real > > con_
int getNumberFunctionEvaluations(void)
void reset(const Vector< Real > &lam, const Real penaltyParameter)
void getConstraintVec(Vector< Real > &c, const Vector< Real > &x)
Real value(const Vector< Real > &x, Real &tol)
Compute value.
int getNumberConstraintEvaluations(void)
Teuchos::RCP< Vector< Real > > dc1_
void getObjectiveGradient(Vector< Real > &g, const Vector< Real > &x)
virtual void set(const Vector &x)
Set where .
int getNumberGradientEvaluations(void)
Teuchos::RCP< Vector< Real > > lam_
Teuchos::RCP< Vector< Real > > x_
AugmentedLagrangian(Objective< Real > &obj, EqualityConstraint< Real > &con, const Vector< Real > &x, const Vector< Real > &c, const Vector< Real > &l, const Real mu, Teuchos::ParameterList &parlist)
Teuchos::RCP< Vector< Real > > dlam_
static const double ROL_EPSILON
Platform-dependent machine epsilon.