ROL
ROL_ObjectiveFromBoundConstraint.hpp
Go to the documentation of this file.
1// Redistribution and use in source and binary forms, with or without
2// modification, are permitted provided that the following conditions are
3// met:
4//
5// 1. Redistributions of source code must retain the above copyright
6// notice, this list of conditions and the following disclaimer.
7//
8// 2. Redistributions in binary form must reproduce the above copyright
9// notice, this list of conditions and the following disclaimer in the
10// documentation and/or other materials provided with the distribution.
11//
12// 3. Neither the name of the Corporation nor the names of the
13// contributors may be used to endorse or promote products derived from
14// this software without specific prior written permission.
15//
16// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
17// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
18// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
19// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
20// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
21// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
22// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
23// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
24// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
25// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
26// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
27//
28// Questions? Contact lead developers:
29// Drew Kouri (dpkouri@sandia.gov) and
30// Denis Ridzal (dridzal@sandia.gov)
31//
32// ************************************************************************
33// @HEADER
34
35#ifndef ROL_OBJECTIVE_FROM_BOUND_CONSTRAINT_H
36#define ROL_OBJECTIVE_FROM_BOUND_CONSTRAINT_H
37
38#include "ROL_Objective.hpp"
40
41#include "ROL_ParameterList.hpp"
42
43namespace ROL {
44
45
51template <class Real>
53
54 typedef Vector<Real> V;
55
56 typedef Elementwise::Fill<Real> Fill;
57 typedef Elementwise::Reciprocal<Real> Reciprocal;
58 typedef Elementwise::Power<Real> Power;
59 typedef Elementwise::Logarithm<Real> Logarithm;
60 typedef Elementwise::Multiply<Real> Multiply;
61 typedef Elementwise::Heaviside<Real> Heaviside;
62 typedef Elementwise::ThresholdUpper<Real> ThresholdUpper;
63 typedef Elementwise::ThresholdLower<Real> ThresholdLower;
64 typedef Elementwise::ReductionSum<Real> Sum;
65 typedef Elementwise::UnaryFunction<Real> UnaryFunction;
66
67
68
75
76 inline std::string EBarrierToString( EBarrierType type ) {
77 std::string retString;
78 switch(type) {
80 retString = "Logarithmic";
81 break;
83 retString = "Quadratic";
84 break;
86 retString = "Double Well";
87 break;
88 case BARRIER_LAST:
89 retString = "Last Type (Dummy)";
90 break;
91 default:
92 retString = "Invalid EBarrierType";
93 break;
94 }
95 return retString;
96 }
97
98 inline EBarrierType StringToEBarrierType( std::string s ) {
99 s = removeStringFormat(s);
101 for( int to = BARRIER_LOGARITHM; to != BARRIER_LAST; ++to ) {
102 type = static_cast<EBarrierType>(to);
103 if( !s.compare(removeStringFormat(EBarrierToString(type))) ) {
104 return type;
105 }
106 }
107 return type;
108 }
109
110private:
111 const ROL::Ptr<const V> lo_;
112 const ROL::Ptr<const V> up_;
113 ROL::Ptr<V> a_; // scratch vector
114 ROL::Ptr<V> b_; // scratch vector
118
119public:
120
122 ROL::ParameterList &parlist ) :
123 lo_( bc.getLowerBound() ),
124 up_( bc.getUpperBound() ) {
125
128
129 a_ = lo_->clone();
130 b_ = up_->clone();
131
132 std::string bfstring = parlist.sublist("Barrier Function").get("Type","Logarithmic");
133 btype_ = StringToEBarrierType(bfstring);
134 }
135
137 lo_( bc.getLowerBound() ),
138 up_( bc.getUpperBound() ),
140
143
144 a_ = lo_->clone();
145 b_ = up_->clone();
146 }
147
148
149 Real value( const Vector<Real> &x, Real &tol ) {
150 const Real zero(0), one(1), two(2);
151
152 ROL::Ptr<UnaryFunction> func;
153
154 a_->zero(); b_->zero();
155 switch(btype_) {
157
158 if ( isLowerActivated_ ) {
159 a_->set(x); // a = x
160 a_->axpy(-one,*lo_); // a = x-l
161 a_->applyUnary(Logarithm()); // a = log(x-l)
162 }
163
164 if ( isUpperActivated_ ) {
165 b_->set(*up_); // b = u
166 b_->axpy(-one,x); // b = u-x
167 b_->applyUnary(Logarithm()); // b = log(u-x)
168 }
169
170 b_->plus(*a_); // b = log(x-l)+log(u-x)
171 b_->scale(-one); // b = -log(x-l)-log(u-x)
172
173 break;
174
176
177 if ( isLowerActivated_ ) {
178 a_->set(x); // a = x
179 a_->axpy(-one,*lo_); // a = x-l
180 a_->applyUnary(ThresholdLower(zero)); // a = min(x-l,0)
181 a_->applyUnary(Power(two)); // a = min(x-l,0)^2
182 }
183
184 if ( isUpperActivated_ ) {
185 b_->set(*up_); // b = u
186 b_->axpy(-one,x); // b = u-x
187 b_->applyUnary(ThresholdUpper(zero)); // b = max(x-u,0)
188 b_->applyUnary(Power(two)); // b = max(x-u,0)^2
189 }
190
191 b_->plus(*a_); // b = min(x-l,0)^2 + max(x-u,0)^2
192
193 break;
194
196
197 if ( isLowerActivated_ ) {
198 a_->set(x); // a = x
199 a_->axpy(-one,*lo_); // a = x-l
200 a_->applyUnary(Power(two)); // a = (x-l)^2
201 }
202 else {
203 a_->applyUnary(Fill(one));
204 }
205
206 if ( isUpperActivated_ ) {
207 b_->set(*up_); // b = u
208 b_->axpy(-one,x); // b = u-x
209 b_->applyUnary(Power(two)); // b = (u-x)^2
210 }
211 else {
212 b_->applyUnary(Fill(one));
213 }
214
215 b_->applyBinary(Multiply(),*a_); // b = (x-l)^2*(u-x)^2
216
217 break;
218
219 default:
220
221 ROL_TEST_FOR_EXCEPTION(true,std::invalid_argument,
222 ">>>(ObjectiveFromBoundConstraint::value): Undefined barrier function type!");
223
224 }
225
226 Real result = b_->reduce(Sum());
227 return result;
228
229 }
230
231 void gradient( Vector<Real> &g, const Vector<Real> &x, Real &tol ) {
232 const Real zero(0), one(1), two(2);
233
234 a_->zero(); b_->zero();
235 switch(btype_) {
237
238 if ( isLowerActivated_ ) {
239 a_->set(*lo_); // a = l
240 a_->axpy(-one,x); // a = l-x
241 a_->applyUnary(Reciprocal()); // a = -1/(x-l)
242 }
243
244 if ( isUpperActivated_ ) {
245 b_->set(*up_); // b = u
246 b_->axpy(-one,x); // b = u-x
247 b_->applyUnary(Reciprocal()); // b = 1/(u-x)
248 }
249
250 b_->plus(*a_); // b = -1/(x-l)+1/(u-x)
251
252 break;
253
255
256 if ( isLowerActivated_ ) {
257 a_->set(x); // a = x
258 a_->axpy(-one,*lo_); // a = x-l
259 a_->applyUnary(ThresholdLower(zero)); // a = min(x-l,0)
260 }
261
262 if ( isUpperActivated_ ) {
263 b_->set(*up_); // b = u
264 b_->axpy(-one,x); // b = u-x
265 b_->applyUnary(ThresholdUpper(zero)); // b = max(x-u,0)
266 }
267
268 b_->plus(*a_); // b = max(x-u,0) + min(x-l,0)
269 b_->scale(two); // b = 2*max(x-u,0) + 2*min(x-l,0)
270 break;
271
273
274 if ( isLowerActivated_ ) {
275 a_->set(x); // a = l
276 a_->axpy(-one,*lo_); // a = x-l
277 }
278 else {
279 a_->applyUnary(Fill(one));
280 }
281
282 if ( isUpperActivated_ ) {
283 b_->set(*up_); // b = x
284 b_->axpy(-one,x); // b = u-x
285 }
286 else {
287 b_->applyUnary(Fill(one));
288 }
289
290 b_->applyBinary(Multiply(),*a_); // b = (x-l)*(u-x)
291 b_->scale(two); // b = 2*(x-l)*(u-x)
292
294 a_->set(*up_); // a = u
295 a_->axpy(-two,x); // a = (u-x)-x
296 a_->plus(*lo_); // a = (u-x)-(x-l)
297 b_->applyBinary(Multiply(),*b_); // b = 2*(x-l)*(u-x)*[(u-x)-(x-l)]
298 }
299
300 break;
301
302 default:
303
304 ROL_TEST_FOR_EXCEPTION(true,std::invalid_argument,
305 ">>>(ObjectiveFromBoundConstraint::gradient): Undefined barrier function type!");
306
307 }
308
309 g.set(*b_);
310
311 }
312
313 void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
314 const Real one(1), two(2), eight(8);
315
316 switch(btype_) {
318
319 if ( isLowerActivated_ ) {
320 a_->set(x); // a = x
321 a_->axpy(-one,*lo_); // a = x-l
322 a_->applyUnary(Reciprocal()); // a = 1/(x-l)
323 a_->applyUnary(Power(two)); // a = 1/(x-l)^2
324 }
325
326 if ( isUpperActivated_ ) {
327 b_->set(*up_); // b = u
328 b_->axpy(-one,x); // b = u-x
329 b_->applyUnary(Reciprocal()); // b = 1/(u-x)
330 b_->applyUnary(Power(two)); // b = 1/(u-x)^2
331 }
332
333 b_->plus(*a_); // b = 1/(x-l)^2 + 1/(u-x)^2
334
335 break;
336
338
339 if ( isLowerActivated_ ) {
340 a_->set(*lo_); // a = l
341 a_->axpy(-one,x); // a = l-x
342 a_->applyUnary(Heaviside()); // a = theta(l-x)
343 }
344
345 if ( isUpperActivated_ ) {
346 b_->set(x); // b = x
347 b_->axpy(-one,*up_); // b = x-u
348 b_->applyUnary(Heaviside()); // b = theta(x-u)
349 }
350
351 b_->plus(*a_); // b = theta(l-x) + theta(x-u)
352 b_->scale(two); // b = 2*theta(l-x) + 2*theta(x-u)
353
354 break;
355
357
359 a_->set(x); // a = x
360 a_->axpy(-one,*lo_); // a = x-l
361
362 b_->set(*up_); // b = u
363 b_->axpy(-one,x); // b = u-x
364
365 b_->applyBinary(Multiply(),*a_); // b = (u-x)*(x-l)
366 b_->scale(-eight); // b = -8*(u-x)*(x-l)
367
368 a_->applyUnary(Power(two)); // a = (x-l)^2
369 a_->scale(two); // a = 2*(x-l)^2
370
371 b_->plus(*a_); // b = 2*(x-l)^2-8*(u-x)*(x-l)
372
373 a_->set(*up_); // a = u
374 a_->axpy(-one,x); // a = u-x
375 a_->applyUnary(Power(two)); // a = (u-x)^2
376 a_->scale(two); // a = 2*(u-x)^2
377
378 b_->plus(*a_); // b = 2*(u-x)^2-8*(u-x)*(x-l)+2*(x-l)^2
379 }
380 else {
381 b_->applyUnary(Fill(two));
382 }
383
384 break;
385
386 default:
387
388 ROL_TEST_FOR_EXCEPTION(true,std::invalid_argument,
389 ">>>(ObjectiveFromBoundConstraint::hessVec): Undefined barrier function type!");
390
391 }
392
393 hv.set(v);
394 hv.applyBinary(Multiply(),*b_);
395
396 }
397
398 // For testing purposes
399 ROL::Ptr<Vector<Real> > getBarrierVector(void) {
400 return b_;
401 }
402
403
404}; // class ObjectiveFromBoundConstraint
405
406}
407
408
409#endif // ROL_OBJECTIVE_FROM_BOUND_CONSTRAINT_H
410
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Provides the interface to apply upper and lower bound constraints.
bool isLowerActivated(void) const
Check if lower bound are on.
bool isUpperActivated(void) const
Check if upper bound are on.
Create a penalty objective from upper and lower bound vectors.
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Elementwise::ThresholdLower< Real > ThresholdLower
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
enum ROL::ObjectiveFromBoundConstraint::EBarrierType eBarrierType_
ObjectiveFromBoundConstraint(const BoundConstraint< Real > &bc)
ROL::Ptr< Vector< Real > > getBarrierVector(void)
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Elementwise::ThresholdUpper< Real > ThresholdUpper
Elementwise::UnaryFunction< Real > UnaryFunction
ObjectiveFromBoundConstraint(const BoundConstraint< Real > &bc, ROL::ParameterList &parlist)
Provides the interface to evaluate objective functions.
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:84
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:209
virtual void applyBinary(const Elementwise::BinaryFunction< Real > &f, const Vector &x)
Definition: ROL_Vector.hpp:248
std::string removeStringFormat(std::string s)
Definition: ROL_Types.hpp:249