45 #ifndef ROL_BOUNDFLETCHER_H
46 #define ROL_BOUNDFLETCHER_H
53 #include "ROL_Ptr.hpp"
66 Ptr<const Vector<Real> >
low_;
67 Ptr<const Vector<Real> >
upp_;
96 Ptr<Vector<Real> >
Q_;
98 Ptr<Vector<Real> >
DQ_;
119 Ptr<Vector<Real> >
w_;
120 Ptr<Vector<Real> >
v_;
146 class DiffLower :
public Elementwise::BinaryFunction<Real> {
149 Real
apply(
const Real& x,
const Real& y)
const {
150 const Real NINF(ROL_NINF<Real>());
151 return (y <= NINF ? static_cast<Real>(-1.) : x - y);
155 class DiffUpper :
public Elementwise::BinaryFunction<Real> {
158 Real
apply(
const Real& x,
const Real& y)
const {
159 const Real INF(ROL_INF<Real>());
160 return (y >= INF ? static_cast<Real>(-1.) : y - x);
164 class FormQ :
public Elementwise::BinaryFunction<Real> {
167 Real
apply(
const Real& x,
const Real& y)
const {
169 if( x < zero && y < zero) {
170 return static_cast<Real
>(1);
172 if( x < zero && y >= zero ) {
175 if( x >= zero && y < zero ) {
178 return std::min(x, y);
182 class FormDQ :
public Elementwise::BinaryFunction<Real> {
185 Real
apply(
const Real& x,
const Real& y)
const {
186 Real
zero(0.), one(1.), mone(-1.);
190 if( x < zero && y >=
zero ) {
208 const Ptr<Constraint<Real> >
con_;
209 const Ptr<const Vector<Real> >
x_;
211 const Ptr<Vector<Real> >
Qv_;
224 con_->applyAdjointJacobian(*(Hvp.
get(0)), *(vp.
get(1)), *
x_, tol);
225 Hvp.
get(0)->applyBinary(Elementwise::Multiply<Real>(), *
Qsqrt_);
226 Hvp.
get(0)->plus(*(vp.
get(0)));
229 Qv_->applyBinary(Elementwise::Multiply<Real>(), *
Qsqrt_);
237 const Ptr<Constraint<Real> >
con_;
238 const Ptr<const Vector<Real> >
x_;
239 const Ptr<Vector<Real> >
Q_;
240 const Ptr<Vector<Real> >
Qv_;
253 con_->applyAdjointJacobian(*(Hvp.
get(0)), *(vp.
get(1)), *
x_, tol);
254 Hvp.
get(0)->plus(*(vp.
get(0)));
257 Qv_->applyBinary(Elementwise::Multiply<Real>(), *
Q_);
265 const Ptr<Constraint<Real> >
con_;
266 const Ptr<const Vector<Real> >
x_;
281 con_->applyPreconditioner(*(Hvp.
get(1)),*(vp.
get(1)),*
x_,*
x_, zero);
291 ROL::ParameterList &parlist)
294 low_ = bnd->getLowerBound();
295 upp_ = bnd->getUpperBound();
298 y_ = conVec.
dual().clone();
299 g_ = optVec.
dual().clone();
315 w_ = optVec.
dual().clone();
316 v_ = conVec.
dual().clone();
327 vv_ = makePtr<PartitionedVector<Real>>(std::vector<Ptr<Vector<Real>> >({
v1_,
v2_}));
331 bb_ = makePtr<PartitionedVector<Real>>(std::vector<Ptr<Vector<Real>> >({
b1_,
b2_}));
333 ROL::ParameterList& sublist = parlist.sublist(
"Step").sublist(
"Fletcher");
334 HessianApprox_ = sublist.get(
"Level of Hessian Approximation", 0);
336 AugSolve_ = sublist.get(
"Type of Augmented System Solve", 0);
341 delta_ = sublist.get(
"Regularization Parameter", 0.0);
343 useInexact_ = sublist.get(
"Inexact Solves",
false);
345 ROL::ParameterList krylovList;
346 Real atol =
static_cast<Real
>(1e-12);
347 Real rtol =
static_cast<Real
>(1e-2);
348 krylovList.sublist(
"General").sublist(
"Krylov").set(
"Type",
"GMRES");
349 krylovList.sublist(
"General").sublist(
"Krylov").set(
"Absolute Tolerance", atol);
350 krylovList.sublist(
"General").sublist(
"Krylov").set(
"Relative Tolerance", rtol);
351 krylovList.sublist(
"General").sublist(
"Krylov").set(
"Iteration Limit", 200);
352 krylov_ = KrylovFactory<Real>(krylovList);
356 obj_->update(x,flag,iter);
357 con_->update(x,flag,iter);
409 w_->applyBinary(Elementwise::Multiply<Real>(), *
Qsqrt_);
410 con_->applyAdjointHessian( *
gPhi_, *
y_, *
w_, x, tol2 ); tol2 = origTol;
411 obj_->hessVec( *
Tv_, *
w_, x, tol2 ); tol2 = origTol;
412 gPhi_->axpy( static_cast<Real>(-1), *
Tv_ );
414 con_->applyAdjointJacobian( *
Tv_, *
v_, x, tol2); tol2 = origTol;
417 Tv_->applyBinary(Elementwise::Multiply<Real>(), *
DQgL_);
420 con_->applyAdjointHessian( *
Tv_, *
v_, *
QgL_, x, tol2 ); tol2 = origTol;
434 Tv_->applyBinary( Elementwise::Multiply<Real>(), *
DQgL_ );
435 gPhi_->axpy(static_cast<Real>(-1), *
Tv_);
437 w_->applyBinary( Elementwise::Multiply<Real>(), *
Q_ );
438 obj_->hessVec( *
Tv_, *
w_, x, tol2); tol2 = origTol;
439 gPhi_->axpy( static_cast<Real>(-1), *
Tv_ );
440 con_->applyAdjointHessian( *
Tv_, *
y_, *
w_, x, tol2 ); tol2 = origTol;
443 con_->applyAdjointHessian( *
Tv_, *
v_, *
QgL_, x, tol2 ); tol2 = origTol;
462 value(x, tol2); tol2 = origTol;
469 obj_->hessVec( hv, v, x, tol2 ); tol2 = origTol;
470 con_->applyAdjointHessian( *
Tv_, *
y_, v, x, tol2 ); tol2 = origTol;
471 hv.
axpy(static_cast<Real>(-1), *
Tv_ );
475 htmp1_->applyBinary(Elementwise::Multiply<Real>(), *
Qsqrt_);
476 htmp1_->scale(static_cast<Real>(-1));
479 Tv_->applyBinary( Elementwise::Multiply<Real>(), v );
481 con_->applyJacobian( *
htmp2_, *
Tv_, x, tol2); tol2 = origTol;
484 con_->applyAdjointJacobian( *
Tv_, *
v_, x, tol2 ); tol2 = origTol;
487 con_->applyJacobian( *
htmp2_, v, x, tol2 ); tol2 = origTol;
489 con_->applyAdjointJacobian( *
Tv_, *
v_, x, tol2 ); tol2 = origTol;
491 Tv_->applyBinary( Elementwise::Multiply<Real>(), *
DQgL_ );
494 w_->applyBinary( Elementwise::Multiply<Real>(), *
Qsqrt_ );
495 obj_->hessVec( *
Tv_, *
w_, x, tol2 ); tol2 = origTol;
496 hv.
axpy( static_cast<Real>(-1), *
Tv_);
497 con_->applyAdjointHessian( *
Tv_, *
y_, *
w_, x, tol2 ); tol2 = origTol;
503 obj_->hessVec( hv, v, x, tol2 ); tol2 = origTol;
504 con_->applyAdjointHessian( *
Tv_, *
y_, v, x, tol2 ); tol2 = origTol;
505 hv.
axpy(static_cast<Real>(-1), *
Tv_ );
509 Tv_->applyBinary( Elementwise::Multiply<Real>(), *
DQgL_ );
511 Tv_->scale( static_cast<Real>(-1) );
512 con_->applyJacobian( *
htmp2_, *
Tv_, x, tol2 ); tol2 = origTol;
516 con_->applyJacobian( *
htmp2_, v, x, tol2 ); tol2 = origTol;
520 Tv_->applyBinary( Elementwise::Multiply<Real>(), *
DQgL_ );
521 hv.
axpy( static_cast<Real>(-1), *
Tv_ );
523 w_->applyBinary( Elementwise::Multiply<Real>(), *
Q_ );
524 obj_->hessVec( *
Tv_, *
w_, x, tol2 ); tol2 = tol;
525 hv.
axpy( static_cast<Real>(-1), *
Tv_);
526 con_->applyAdjointHessian( *
Tv_, *
y_, *
w_, x, tol2 ); tol2 = origTol;
540 ROL::Ptr<LinearOperator<Real> > K;
547 K = ROL::makePtr<AugSystemNonSym>(
con_, makePtrFromRef(x),
Q_,
Qv_,
delta_);
551 ROL::Ptr<LinearOperator<Real> > P
552 = ROL::makePtr<AugSystemPrecond>(
con_, makePtrFromRef(x));
562 krylov_->resetAbsoluteTolerance(tol);
584 Qsg_->applyBinary(Elementwise::Multiply<Real>(), *
Qsqrt_);
590 gL_->applyBinary(Elementwise::Divide<Real>(), *
Qsqrt_);
592 QgL_->applyBinary(Elementwise::Multiply<Real>(), *
Qsqrt_);
599 QgL_->applyBinary(Elementwise::Multiply<Real>(), *
Q_);
605 DQgL_->applyBinary(Elementwise::Multiply<Real>(), *
DQ_);
621 Qsqrt_->applyUnary(Elementwise::SquareRoot<Real>());
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...
Ptr< Vector< Real > > gPhi_
Ptr< Vector< Real > > Qsg_
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
const Ptr< Constraint< Real > > con_
ROL::Ptr< const Vector< Real > > get(size_type i) const
virtual void plus(const Vector &x)=0
Compute , where .
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Defines the linear algebra of vector space on a generic partitioned vector.
AugSystemSym(const Ptr< Constraint< Real > > &con, const Ptr< const Vector< Real > > &x, const Ptr< Vector< Real > > &Qsqrt, const Ptr< Vector< Real > > &Qv, const Real delta)
void objValue(const Vector< Real > &x, Real &tol)
const Ptr< Objective< Real > > obj_
Contains definitions of custom data types in ROL.
const Ptr< Constraint< Real > > con_
Ptr< Vector< Real > > czeros_
Real apply(const Real &x, const Real &y) const
Real apply(const Real &x, const Real &y) const
Defines the linear algebra or vector space interface.
Ptr< Krylov< Real > > krylov_
const Ptr< Vector< Real > > Qv_
const Ptr< Vector< Real > > Q_
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
void applyInverse(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply inverse of linear operator.
const Ptr< const Vector< Real > > x_
void solveAugmentedSystem(Vector< Real > &v1, Vector< Real > &v2, const Vector< Real > &b1, const Vector< Real > &b2, const Vector< Real > &x, Real &tol)
BoundFletcher(const ROL::Ptr< Objective< Real > > &obj, const ROL::Ptr< Constraint< Real > > &con, const ROL::Ptr< BoundConstraint< Real > > &bnd, const Vector< Real > &optVec, const Vector< Real > &conVec, ROL::ParameterList &parlist)
Ptr< Vector< Real > > QgL_
Ptr< Vector< Real > > gL_
Ptr< Vector< Real > > DQgL_
Ptr< const Vector< Real > > low_
Ptr< Vector< Real > > Tv_
Ptr< Vector< Real > > v1_
Ptr< Vector< Real > > b1_
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
Ptr< Vector< Real > > Qv_
Real value(const Vector< Real > &x, Real &tol)
Compute value.
bool isMultiplierComputed_
void computeMultipliers(const Vector< Real > &x, const Real tol)
const Ptr< const Vector< Real > > x_
Ptr< Vector< Real > > scaledc_
const Ptr< Constraint< Real > > con_
void conValue(const Vector< Real > &x, Real &tol)
Ptr< Vector< Real > > Qsqrt_
const Ptr< Constraint< Real > > con_
void computeDQ(const Vector< Real > &x)
void set(const V &x)
Set where .
const Ptr< Vector< Real > > Qsqrt_
Ptr< PartitionedVector< Real > > vv_
AugSystemPrecond(const Ptr< Constraint< Real > > con, const Ptr< const Vector< Real > > x)
Provides the interface to apply a linear operator.
Provides the interface to apply upper and lower bound constraints.
Ptr< const Vector< Real > > upp_
AugSystemNonSym(const Ptr< Constraint< Real > > &con, const Ptr< const Vector< Real > > &x, const Ptr< Vector< Real > > &Q, const Ptr< Vector< Real > > &Qv, const Real delta)
Ptr< Vector< Real > > umx_
Ptr< Vector< Real > > DQ_
const Ptr< Vector< Real > > Qv_
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
virtual void set(const Vector &x)
Set where .
Ptr< Vector< Real > > v2_
Ptr< Vector< Real > > QsgL_
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Ptr< Vector< Real > > b2_
Ptr< PartitionedVector< Real > > bb_
void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
Ptr< Vector< Real > > htmp1_
const Ptr< const Vector< Real > > x_
void computeQ(const Vector< Real > &x)
Ptr< Vector< Real > > htmp2_
void objGrad(const Vector< Real > &x, Real &tol)
Defines the general constraint operator interface.
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Ptr< Vector< Real > > xzeros_