10 #ifndef ROL_INTERIORPOINTOBJECTIVE_H
11 #define ROL_INTERIORPOINTOBJECTIVE_H
15 #include "ROL_ParameterList.hpp"
36 const Ptr<Objective<Real>>
obj_;
37 const Ptr<BoundConstraint<Real>>
bnd_;
38 const Ptr<const Vector<Real>>
lo_;
39 const Ptr<const Vector<Real>>
up_;
53 Ptr<ScalarController<Real,int>>
fval_;
63 Real
apply(
const Real &x )
const {
64 const Real
zero(0), NINF(ROL_NINF<Real>());
65 return (x>
zero) ? std::log(x) : NINF;
74 Real
apply(
const Real &x )
const {
75 const Real
zero(0), one(1);
86 Real
apply(
const Real &x,
const Real &y )
const {
88 return (x>zero) ? y/x :
zero;
97 class Mask :
public Elementwise::BinaryFunction<Real> {
102 Real
apply(
const Real &x,
const Real &y )
const {
109 const Real
zero(0), one(1);
111 fval_ = makePtr<ScalarController<Real,int>>();
112 gradient_ = makePtr<VectorController<Real,int>>();
115 ValueSet isBoundedBelow( ROL_NINF<Real>(), ValueSet::GREATER_THAN, one,
zero );
116 ValueSet isBoundedAbove( ROL_INF<Real>(), ValueSet::LESS_THAN, one,
zero );
143 const bool useLinearDamping,
146 :
obj_(obj),
bnd_(bnd),
lo_(bnd->getLowerBound()),
up_(bnd->getUpperBound()),
156 ParameterList &parlist )
157 :
obj_(obj),
bnd_(bnd),
lo_(bnd->getLowerBound()),
up_(bnd->getUpperBound()),
159 ParameterList &iplist = parlist.sublist(
"Step").sublist(
"Primal Dual Interior Point");
160 ParameterList &lblist = iplist.sublist(
"Barrier Objective");
163 kappaD_ = lblist.get(
"Linear Damping Coefficient", 1.e-4);
164 mu_ = lblist.get(
"Initial Barrier Parameter", 0.1);
172 bool isComputed =
fval_->get(val,key);
202 obj_->update(x,type,iter);
203 fval_->objectiveUpdate(type);
208 const Real
zero(0), one(1);
209 Real linearTerm =
zero;
214 Elementwise::ReductionSum<Real> sum;
215 Elementwise::Multiply<Real> mult;
223 pwa_->applyUnary(mlog);
232 pwa_->applyUnary(mlog);
235 fval -=
mu_*(aval+bval);
250 pwa_->applyUnary(mrec);
259 pwa_->applyUnary(mrec);
268 const Real one(1), two(2);
270 obj_->hessVec(hv,v,x,tol);
274 Elementwise::Multiply<Real> mult;
275 Elementwise::Power<Real> square(two);
279 pwa_->applyUnary(mrec);
281 pwa_->applyUnary(square);
283 pwa_->applyBinary(mult,v);
288 pwa_->applyUnary(mrec);
290 pwa_->applyUnary(square);
292 pwa_->applyBinary(mult,v);
300 #endif // ROL_INTERIORPOINTOBJECTIVE_H
Provides the interface to evaluate objective functions.
InteriorPointObjective(const Ptr< Objective< Real >> &obj, const Ptr< BoundConstraint< Real >> &bnd, const Vector< Real > &x, const Vector< Real > &g, ParameterList &parlist)
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
const Ptr< const Vector< Real > > lo_
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Ptr< Vector< Real > > maskL0_
Defines the linear algebra or vector space interface.
void initialize(const Vector< Real > &x, const Vector< Real > &g)
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
void updatePenalty(const Real mu)
Ptr< ScalarController< Real, int > > fval_
Real apply(const Real &x) const
Ptr< Vector< Real > > maskL_
Real apply(const Real &x, const Real &y) const
int getNumberGradientEvaluations(void) const
Ptr< Vector< Real > > maskU0_
void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
const Ptr< Objective< Real > > obj_
Ptr< Vector< Real > > pwa_
Provides the interface to apply upper and lower bound constraints.
Real apply(const Real &x) const
Ptr< Vector< Real > > maskU_
Elementwise::ValueSet< Real > ValueSet
const Ptr< BoundConstraint< Real > > bnd_
Real value(const Vector< Real > &x, Real &tol)
Compute value.
const Ptr< const Vector< Real > > getObjectiveGradient(const Vector< Real > &x, Real &tol)
const Ptr< const Vector< Real > > up_
int getNumberFunctionEvaluations(void) const
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
virtual void set(const Vector &x)
Set where .
Ptr< VectorController< Real, int > > gradient_
InteriorPointObjective(const Ptr< Objective< Real >> &obj, const Ptr< BoundConstraint< Real >> &bnd, const Vector< Real > &x, const Vector< Real > &g, const bool useLinearDamping, const Real kappaD, const Real mu)
Real apply(const Real &x, const Real &y) const
Real getObjectiveValue(const Vector< Real > &x, Real &tol)