10 #ifndef ROL_OBJECTIVE_FROM_BOUND_CONSTRAINT_H
11 #define ROL_OBJECTIVE_FROM_BOUND_CONSTRAINT_H
16 #include "ROL_ParameterList.hpp"
31 typedef Elementwise::Fill<Real>
Fill;
33 typedef Elementwise::Power<Real>
Power;
39 typedef Elementwise::ReductionSum<Real>
Sum;
52 std::string retString;
55 retString =
"Logarithmic";
58 retString =
"Quadratic";
61 retString =
"Double Well";
64 retString =
"Last Type (Dummy)";
67 retString =
"Invalid EBarrierType";
86 const ROL::Ptr<const V>
lo_;
87 const ROL::Ptr<const V>
up_;
97 ROL::ParameterList &parlist ) :
98 lo_( bc.getLowerBound() ),
99 up_( bc.getUpperBound() ) {
107 std::string bfstring = parlist.sublist(
"Barrier Function").get(
"Type",
"Logarithmic");
112 lo_( bc.getLowerBound() ),
113 up_( bc.getUpperBound() ),
125 const Real
zero(0), one(1), two(2);
127 ROL::Ptr<UnaryFunction> func;
129 a_->zero();
b_->zero();
178 a_->applyUnary(
Fill(one));
187 b_->applyUnary(
Fill(one));
196 ROL_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
197 ">>>(ObjectiveFromBoundConstraint::value): Undefined barrier function type!");
201 Real result =
b_->reduce(
Sum());
207 const Real
zero(0), one(1), two(2);
209 a_->zero();
b_->zero();
254 a_->applyUnary(
Fill(one));
262 b_->applyUnary(
Fill(one));
279 ROL_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
280 ">>>(ObjectiveFromBoundConstraint::gradient): Undefined barrier function type!");
289 const Real one(1), two(2), eight(8);
356 b_->applyUnary(
Fill(two));
363 ROL_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
364 ">>>(ObjectiveFromBoundConstraint::hessVec): Undefined barrier function type!");
384 #endif // ROL_OBJECTIVE_FROM_BOUND_CONSTRAINT_H
Provides the interface to evaluate objective functions.
const ROL::Ptr< const V > up_
virtual void applyBinary(const Elementwise::BinaryFunction< Real > &f, const Vector &x)
Elementwise::Multiply< Real > Multiply
std::string removeStringFormat(std::string s)
Elementwise::ThresholdUpper< Real > ThresholdUpper
ROL::Ptr< Vector< Real > > getBarrierVector(void)
Elementwise::Logarithm< Real > Logarithm
ObjectiveFromBoundConstraint(const BoundConstraint< Real > &bc)
Elementwise::Power< Real > Power
Defines the linear algebra or vector space interface.
ObjectiveFromBoundConstraint(const BoundConstraint< Real > &bc, ROL::ParameterList &parlist)
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Elementwise::ReductionSum< Real > Sum
enum ROL::ObjectiveFromBoundConstraint::EBarrierType eBarrierType_
std::string EBarrierToString(EBarrierType type)
Elementwise::UnaryFunction< Real > UnaryFunction
EBarrierType StringToEBarrierType(std::string s)
Elementwise::Reciprocal< Real > Reciprocal
Provides the interface to apply upper and lower bound constraints.
Elementwise::Heaviside< Real > Heaviside
Real value(const Vector< Real > &x, Real &tol)
Compute value.
const ROL::Ptr< const V > lo_
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
virtual void set(const Vector &x)
Set where .
bool isLowerActivated(void) const
Check if lower bound are on.
Elementwise::Fill< Real > Fill
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Elementwise::ThresholdLower< Real > ThresholdLower
Create a penalty objective from upper and lower bound vectors.
bool isUpperActivated(void) const
Check if upper bound are on.