10 #ifndef ROL_TRUSTREGIONMODEL_H
11 #define ROL_TRUSTREGIONMODEL_H
36 Ptr<BoundConstraint<Real>>
bnd_;
37 Ptr<const Vector<Real>>
x_,
g_;
62 obj_->hessVec(hv,v,*
x_,tol);
71 obj_->invHessVec(hv,v,*
x_,tol);
80 obj_->precond(Pv,v,*
x_,tol);
94 const bool useSecantPrecond =
false,
const bool useSecantHessVec =
false)
95 :
obj_(makePtrFromRef(obj)),
bnd_(makePtrFromRef(bnd)),
96 x_(makePtrFromRef(x)),
g_(makePtrFromRef(g)),
107 obj_ = makePtrFromRef(obj);
108 bnd_ = makePtrFromRef(bnd);
109 x_ = makePtrFromRef(x);
110 g_ = makePtrFromRef(g);
120 dual_->scale(static_cast<Real>(0.5));
152 virtual const Ptr<const Vector<Real>>
getIterate(
void)
const {
161 if (!
bnd_->isActivated()) {
Provides the interface to evaluate objective functions.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
const bool useSecantHessVec_
Ptr< Secant< Real > > secant_
virtual void plus(const Vector &x)=0
Compute , where .
Ptr< Objective< Real > > obj_
virtual void update(Objective< Real > &obj, BoundConstraint< Real > &bnd, const Vector< Real > &x, const Vector< Real > &g, const Ptr< Secant< Real >> &secant=nullPtr)
virtual ~TrustRegionModel()
Ptr< const Vector< Real > > g_
Provides the interface to evaluate trust-region model functions.
Defines the linear algebra or vector space interface.
virtual const Ptr< const Vector< Real > > getGradient(void) const
virtual void updatePredictedReduction(Real &pred, const Vector< Real > &s)
virtual void dualTransform(Vector< Real > &tv, const Vector< Real > &v)
void applyPrecond(Vector< Real > &Pv, const Vector< Real > &v, Real &tol)
void applyHessian(Vector< Real > &hv, const Vector< Real > &v, Real &tol)
TrustRegionModel(Objective< Real > &obj, BoundConstraint< Real > &bnd, const Vector< Real > &x, const Vector< Real > &g, const Ptr< Secant< Real >> &secant=nullPtr, const bool useSecantPrecond=false, const bool useSecantHessVec=false)
virtual const Ptr< Objective< Real > > getObjective(void) const
void applyInvHessian(Vector< Real > &hv, const Vector< Real > &v, Real &tol)
virtual Real value(const Vector< Real > &s, Real &tol)
Compute value.
Provides interface for and implements limited-memory secant operators.
const bool useSecantPrecond_
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol)
Apply Hessian approximation to vector.
void initialize(const Vector< Real > &s)
virtual void invHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol)
Apply inverse Hessian approximation to vector.
virtual void gradient(Vector< Real > &g, const Vector< Real > &s, Real &tol)
Compute gradient.
Provides the interface to apply upper and lower bound constraints.
Ptr< BoundConstraint< Real > > bnd_
virtual void set(const Vector &x)
Set where .
Ptr< Vector< Real > > dual_
virtual void updateActualReduction(Real &ared, const Vector< Real > &s)
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &s, Real &tol)
Apply preconditioner to vector.
Ptr< const Vector< Real > > x_
virtual const Ptr< BoundConstraint< Real > > getBoundConstraint(void) const
virtual void primalTransform(Vector< Real > &tv, const Vector< Real > &v)
virtual const Ptr< const Vector< Real > > getIterate(void) const