45 #ifndef ROL_FLETCHEROBJECTIVEEDEF_H
46 #define ROL_FLETCHEROBJECTIVEEDEF_H
50 template<
typename Real>
57 ROL::ParameterList &parlist)
70 template<
typename Real>
74 bool isComputed = fPhi_->get(val,key);
75 if( isComputed && multSolverError_*cnorm_ <= tol) {
76 tol = multSolverError_*cnorm_;
84 multSolverError_ = origTol / (
static_cast<Real
>(2) * std::max(static_cast<Real>(1), cnorm_));
86 tol = multSolverError_*cnorm_;
88 val = fval - cprim_->apply(*cdual_);
89 if( quadPenaltyParameter_ > static_cast<Real>(0) )
90 val +=
static_cast<Real
>(0.5)*quadPenaltyParameter_*(cprim_->dot(*cprim_));
97 template<
typename Real>
100 bool isComputed = gPhi_->get(g,key);
101 if( isComputed && gradSolveError_ <= tol) {
102 tol = gradSolveError_;
109 gradSolveError_ = origTol /
static_cast<Real
>(2);
111 gL_->set(gLdual_->dual());
112 bool refine = isComputed;
114 solveAugmentedSystem( *wdual_, *vg_, *xzeros_, *cprim_, x, gradSolveError_, refine );
115 gradSolveError_ += multSolverError_;
116 tol = gradSolveError_;
117 wg_->set(wdual_->dual());
118 con_->applyAdjointHessian( g, *cdual_, *wdual_, x, tol2 ); tol2 = origTol;
119 g.
axpy( sigma_, *wg_ );
120 obj_->hessVec( *Tv_, *wdual_, x, tol2 ); tol2 = origTol;
121 g.
axpy( static_cast<Real>(-1), *Tv_ );
122 con_->applyAdjointHessian( *Tv_, *vg_, *gLdual_, x, tol2 ); tol2 = origTol;
125 if( quadPenaltyParameter_ > static_cast<Real>(0) ) {
126 con_->applyAdjointJacobian( *Tv_, *cprim_, x, tol2 ); tol2 = origTol;
127 g.
axpy( quadPenaltyParameter_, *Tv_ );
133 template<
typename Real>
139 bool isComputed = y_->get(*cdual_,key);
140 if( !isComputed || !useInexact_) {
149 obj_->hessVec( hv, v, x, tol2 ); tol2 = origTol;
150 con_->applyAdjointHessian( *Tv_, *cdual_, v, x, tol2 ); tol2 = origTol;
151 hv.
axpy(static_cast<Real>(-1), *Tv_ );
154 solveAugmentedSystem( *w_, *v_, hv, *czeros_, x, tol2 ); tol2 = origTol;
155 hv.
scale( static_cast<Real>(-1) );
160 solveAugmentedSystem( *w_, *v_, *Tv_, *czeros_, x, tol2 ); tol2 = origTol;
161 hv.
axpy(static_cast<Real>(-2)*sigma_, *w_);
163 wdual_->set(w_->dual());
165 obj_->hessVec( *Tv_, *wdual_, x, tol2 ); tol2 = origTol;
167 con_->applyAdjointHessian( *Tv_, *cdual_, *wdual_, x, tol2 ); tol2 = origTol;
168 hv.
axpy( static_cast<Real>(-1), *Tv_ );
170 hv.
axpy( static_cast<Real>(2)*sigma_, v );
172 if( quadPenaltyParameter_ > static_cast<Real>(0) ) {
173 con_->applyJacobian( *b2_, v, x, tol2 ); tol2 = origTol;
174 con_->applyAdjointJacobian( *Tv_, *b2_, x, tol2 ); tol2 = origTol;
175 hv.
axpy( quadPenaltyParameter_, *Tv_ );
176 con_->applyAdjointHessian( *Tv_, *cprim_, v, x, tol2); tol2 = origTol;
177 hv.
axpy( -quadPenaltyParameter_, *Tv_ );
181 template<
typename Real>
190 ROL::Ptr<LinearOperator<Real>>
191 K = ROL::makePtr<AugSystem>(con_, makePtrFromRef(x), delta_);
192 ROL::Ptr<LinearOperator<Real>>
193 P = ROL::makePtr<AugSystemPrecond>(con_, makePtrFromRef(x), makePtrFromRef(b1));
202 K->apply(*bb_, *ww_, tol); tol = origTol;
203 bb_->scale(static_cast<Real>(-1));
209 if( useInexact_ ) krylov_->resetAbsoluteTolerance(tol);
213 tol = krylov_->run(*vv_,*K,*bb_,*P,iterKrylov_,flagKrylov_);
Provides the interface to evaluate objective functions.
void solveAugmentedSystem(Vector< Real > &v1, Vector< Real > &v2, const Vector< Real > &b1, const Vector< Real > &b2, const Vector< Real > &x, Real &tol, bool refine=false) override
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 computeMultipliers(Vector< Real > &y, Vector< Real > &gL, const Vector< Real > &x, Vector< Real > &g, Vector< Real > &c, Real tol)
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
Real objValue(const Vector< Real > &x, Real &tol)
virtual void plus(const Vector &x)=0
Compute , where .
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Ptr< Vector< Real > > wg_
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol) override
Compute gradient.
Defines the linear algebra or vector space interface.
Ptr< Vector< Real > > czeros_
FletcherObjectiveE(const ROL::Ptr< Objective< Real >> &obj, const ROL::Ptr< Constraint< Real >> &con, const Vector< Real > &xprim, const Vector< Real > &xdual, const Vector< Real > &cprim, const Vector< Real > &cdual, ROL::ParameterList &parlist)
Ptr< Vector< Real > > vg_
Ptr< Vector< Real > > Tv_
Ptr< Vector< Real > > xzeros_
Real value(const Vector< Real > &x, Real &tol) override
Compute value.
Ptr< Vector< Real > > wdual_
virtual void set(const Vector &x)
Set where .
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol) override
Apply Hessian approximation to vector.
Defines the general constraint operator interface.