ROL
ROL_BoundConstraint.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Rapid Optimization Library (ROL) Package
5 // Copyright (2014) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact lead developers:
38 // Drew Kouri (dpkouri@sandia.gov) and
39 // Denis Ridzal (dridzal@sandia.gov)
40 //
41 // ************************************************************************
42 // @HEADER
43 
44 #ifndef ROL_BOUND_CONSTRAINT_H
45 #define ROL_BOUND_CONSTRAINT_H
46 
47 #include "ROL_Vector.hpp"
48 #include "ROL_Types.hpp"
49 #include <iostream>
50 
69 namespace ROL {
70 
71 template <class Real>
73 private:
74  bool activated_;
75 
76 public:
77 
78  virtual ~BoundConstraint() {}
79 
84  BoundConstraint(void) : activated_(true) {}
85 
93  virtual void update( const Vector<Real> &x, bool flag = true, int iter = -1 ) {}
94 
103  virtual void project( Vector<Real> &x ) {}
104 
116  virtual void pruneUpperActive( Vector<Real> &v, const Vector<Real> &x, Real eps = 0.0 ) {}
117 
131  virtual void pruneUpperActive( Vector<Real> &v, const Vector<Real> &g, const Vector<Real> &x, Real eps = 0.0 ) {}
132 
144  virtual void pruneLowerActive( Vector<Real> &v, const Vector<Real> &x, Real eps = 0.0 ) {}
145 
159  virtual void pruneLowerActive( Vector<Real> &v, const Vector<Real> &g, const Vector<Real> &x, Real eps = 0.0 ) {}
160 
166  virtual void setVectorToUpperBound( Vector<Real> &u ) {}
167 
173  virtual void setVectorToLowerBound( Vector<Real> &l ) {}
174 
186  virtual void pruneActive( Vector<Real> &v, const Vector<Real> &x, Real eps = 0.0 ) {
187  this->pruneUpperActive(v,x,eps);
188  this->pruneLowerActive(v,x,eps);
189  }
190 
203  virtual void pruneActive( Vector<Real> &v, const Vector<Real> &g, const Vector<Real> &x, Real eps = 0.0 ) {
204  this->pruneUpperActive(v,g,x,eps);
205  this->pruneLowerActive(v,g,x,eps);
206  }
207 
213  virtual bool isFeasible( const Vector<Real> &v ) {
214  if ( this->activated_ ) {
215  Teuchos::RCP<Vector<Real> > pv = v.clone();
216  pv->set(v);
217  this->project(*pv);
218  pv->axpy(-1.0,v);
219  Real norm = pv->norm();
220  if ( norm < ROL_EPSILON ) {
221  return true;
222  }
223  else {
224  return false;
225  }
226  }
227  else {
228  return true;
229  }
230  }
231 
236  void activate(void) { this->activated_ = true; }
237 
242  void deactivate(void) { this->activated_ = false; }
243 
248  bool isActivated(void) { return this->activated_; }
249 
257  void pruneInactive( Vector<Real> &v, const Vector<Real> &x, Real eps = 0.0 ) {
258  Teuchos::RCP<Vector<Real> > tmp = x.clone();
259  tmp->set(v);
260  this->pruneActive(*tmp,x,eps);
261  v.axpy(-1.0,*tmp);
262  }
263 
272  void pruneInactive( Vector<Real> &v, const Vector<Real> &g, const Vector<Real> &x, Real eps = 0.0 ) {
273  Teuchos::RCP<Vector<Real> > tmp = x.clone();
274  tmp->set(v);
275  this->pruneActive(*tmp,g,x,eps);
276  v.axpy(-1.0,*tmp);
277  }
278 
286  Teuchos::RCP<Vector<Real> > tmp = g.clone();
287  tmp->set(g);
288  this->pruneActive(g,*tmp,x);
289  }
290 
298  v.plus(x);
299  this->project(v);
300  v.axpy(-1.0,x);
301  }
302 
303 }; // class BoundConstraint
304 
305 } // namespace ROL
306 
307 #endif
bool activated_
Flag that determines whether or not the constraints are being used.
virtual void plus(const Vector &x)=0
Compute , where .
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:141
void activate(void)
Turn on bounds.
Contains definitions of custom data types in ROL.
virtual void pruneLowerActive(Vector< Real > &v, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the lower -active set.
virtual void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update bounds.
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual void setVectorToUpperBound(Vector< Real > &u)
Set the input vector to the upper bound.
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:72
virtual void pruneActive(Vector< Real > &v, const Vector< Real > &g, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the -binding set.
bool isActivated(void)
Check if bounds are on.
virtual void pruneUpperActive(Vector< Real > &v, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the upper -active set.
void computeProjectedStep(Vector< Real > &v, const Vector< Real > &x)
Compute projected step.
void pruneInactive(Vector< Real > &v, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the -inactive set.
BoundConstraint(void)
Default constructor.
virtual void setVectorToLowerBound(Vector< Real > &l)
Set the input vector to the lower bound.
Provides the interface to apply upper and lower bound constraints.
void computeProjectedGradient(Vector< Real > &g, const Vector< Real > &x)
Compute projected gradient.
virtual void pruneUpperActive(Vector< Real > &v, const Vector< Real > &g, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the upper -binding set.
virtual void pruneLowerActive(Vector< Real > &v, const Vector< Real > &g, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the lower -binding set.
virtual void pruneActive(Vector< Real > &v, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the -active set.
void deactivate(void)
Turn off bounds.
void pruneInactive(Vector< Real > &v, const Vector< Real > &g, const Vector< Real > &x, Real eps=0.0)
Set variables to zero if they correspond to the -nonbinding set.
virtual bool isFeasible(const Vector< Real > &v)
Check if the vector, v, is feasible.
virtual void project(Vector< Real > &x)
Project optimization variables onto the bounds.
static const double ROL_EPSILON
Platform-dependent machine epsilon.
Definition: ROL_Types.hpp:115