ROL
ROL_TypeB_PrimalDualActiveSetAlgorithm.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_TYPEB_PRIMALDUALACTIVESETALGORITHM_HPP
45#define ROL_TYPEB_PRIMALDUALACTIVESETALGORITHM_HPP
46
48#include "ROL_KrylovFactory.hpp"
49#include "ROL_SecantFactory.hpp"
50
55namespace ROL {
56namespace TypeB {
57
58template<typename Real>
60private:
61 Ptr<Secant<Real>> secant_;
63 std::string secantName_;
64
65 Ptr<Krylov<Real>> krylov_;
67 std::string krylovName_;
68
72
75
76 int maxit_;
77 int iter_;
78 int flag_;
79 Real stol_;
80 Real gtol_;
81 Real scale_;
82 Real neps_;
83 Real itol_;
87 bool feasible_;
88
92
93 class HessianPDAS : public LinearOperator<Real> {
94 private:
95 const Ptr<Objective<Real>> obj_;
96 const Ptr<BoundConstraint<Real>> bnd_;
97 const Ptr<const Vector<Real>> x_;
98 const Ptr<const Vector<Real>> xlam_;
99 const Real eps_;
100 const Ptr<Secant<Real>> secant_;
101 const bool useSecant_;
102 const Ptr<Vector<Real>> pwa_;
103 public:
105 const Ptr<BoundConstraint<Real>> &bnd,
106 const Ptr<const Vector<Real>> &x,
107 const Ptr<const Vector<Real>> &xlam,
108 Real eps,
109 const Ptr<Secant<Real>> &secant,
110 bool useSecant,
111 const Ptr<Vector<Real>> &pwa)
112 : obj_(obj), bnd_(bnd), x_(x), xlam_(xlam), eps_(eps),
113 secant_(secant), useSecant_(useSecant), pwa_(pwa) {}
114 void apply(Vector<Real> &Hv, const Vector<Real> &v, Real &tol) const {
115 pwa_->set(v);
116 bnd_->pruneActive(*pwa_,*xlam_,eps_);
117 if (!useSecant_) obj_->hessVec(Hv,*pwa_,*x_,tol);
118 else secant_->applyB(Hv,*pwa_);
119 bnd_->pruneActive(Hv,*xlam_,eps_);
120 }
121 };
122
123 class PrecondPDAS : public LinearOperator<Real> {
124 private:
125 const Ptr<Objective<Real>> obj_;
126 const Ptr<BoundConstraint<Real>> bnd_;
127 const Ptr<const Vector<Real>> x_;
128 const Ptr<const Vector<Real>> xlam_;
129 const Real eps_;
130 const Ptr<Secant<Real>> secant_;
131 const bool useSecant_;
132 const Ptr<Vector<Real>> dwa_;
133 public:
135 const Ptr<BoundConstraint<Real>> &bnd,
136 const Ptr<const Vector<Real>> &x,
137 const Ptr<const Vector<Real>> &xlam,
138 Real eps,
139 const Ptr<Secant<Real>> &secant,
140 bool useSecant,
141 const Ptr<Vector<Real>> &dwa)
142 : obj_(obj), bnd_(bnd), x_(x), xlam_(xlam), eps_(eps),
143 secant_(secant), useSecant_(useSecant), dwa_(dwa) {}
144 void apply(Vector<Real> &Hv, const Vector<Real> &v, Real &tol) const {
145 Hv.set(v.dual());
146 }
147 void applyInverse(Vector<Real> &Hv, const Vector<Real> &v, Real &tol) const {
148 dwa_->set(v);
149 bnd_->pruneActive(*dwa_,*xlam_,eps_);
150 if ( useSecant_ ) secant_->applyH(Hv,*dwa_);
151 else obj_->precond(Hv,*dwa_,*x_,tol);
152 bnd_->pruneActive(Hv,*xlam_,eps_);
153 dwa_->set(v);
154 bnd_->pruneInactive(*dwa_,*xlam_,eps_);
155 Hv.plus(dwa_->dual());
156 }
157 };
158
159 class HessianPDAS_Poly : public LinearOperator<Real> {
160 private:
161 const Ptr<Objective<Real>> obj_;
162 const Ptr<BoundConstraint<Real>> bnd_;
163 const Ptr<Constraint<Real>> con_;
164 const Ptr<const Vector<Real>> x_;
165 const Ptr<const Vector<Real>> xlam_;
166 const Real eps_;
167 const Ptr<Secant<Real>> secant_;
168 const bool useSecant_;
169 const Ptr<Vector<Real>> pwa_, dwa_;
170 public:
172 const Ptr<BoundConstraint<Real>> &bnd,
173 const Ptr<Constraint<Real>> &con,
174 const Ptr<const Vector<Real>> &x,
175 const Ptr<const Vector<Real>> &xlam,
176 Real eps,
177 const Ptr<Secant<Real>> &secant,
178 bool useSecant,
179 const Ptr<Vector<Real>> &pwa,
180 const Ptr<Vector<Real>> &dwa)
181 : obj_(obj), bnd_(bnd), con_(con), x_(x), xlam_(xlam), eps_(eps),
182 secant_(secant), useSecant_(useSecant), pwa_(pwa), dwa_(dwa) {}
183 void apply(Vector<Real> &Hv, const Vector<Real> &v, Real &tol) const {
184 PartitionedVector<Real> &Hvp = dynamic_cast<PartitionedVector<Real>&>(Hv);
185 const PartitionedVector<Real> &vp = dynamic_cast<const PartitionedVector<Real>&>(v);
186 pwa_->set(*vp.get(0));
187 bnd_->pruneActive(*pwa_,*xlam_,eps_);
188 if (!useSecant_) obj_->hessVec(*Hvp.get(0),*pwa_,*x_,tol);
189 else secant_->applyB(*Hvp.get(0),*pwa_);
190 con_->applyAdjointJacobian(*dwa_,*vp.get(1),*x_,tol);
191 Hvp.get(0)->plus(*dwa_);
192 bnd_->pruneActive(*Hvp.get(0),*xlam_,eps_);
193 con_->applyJacobian(*Hvp.get(1),*pwa_,*x_,tol);
194 }
195 };
196
197 class PrecondPDAS_Poly : public LinearOperator<Real> {
198 private:
199 const Ptr<Objective<Real>> obj_;
200 const Ptr<BoundConstraint<Real>> bnd_;
201 const Ptr<const Vector<Real>> x_;
202 const Ptr<const Vector<Real>> xlam_;
203 const Real eps_;
204 const Ptr<Secant<Real>> secant_;
205 const bool useSecant_;
206 const Ptr<Vector<Real>> dwa_;
207 public:
209 const Ptr<BoundConstraint<Real>> &bnd,
210 const Ptr<const Vector<Real>> &x,
211 const Ptr<const Vector<Real>> &xlam,
212 Real eps,
213 const Ptr<Secant<Real>> &secant,
214 bool useSecant,
215 const Ptr<Vector<Real>> &dwa)
216 : obj_(obj), bnd_(bnd), x_(x), xlam_(xlam), eps_(eps),
217 secant_(secant), useSecant_(useSecant), dwa_(dwa) {}
218 void apply(Vector<Real> &Hv, const Vector<Real> &v, Real &tol) const {
219 Hv.set(v.dual());
220 }
221 void applyInverse(Vector<Real> &Hv, const Vector<Real> &v, Real &tol) const {
222 PartitionedVector<Real> &Hvp = dynamic_cast<PartitionedVector<Real>&>(Hv);
223 const PartitionedVector<Real> &vp = dynamic_cast<const PartitionedVector<Real>&>(v);
224 dwa_->set(*vp.get(0));
225 bnd_->pruneActive(*dwa_,*xlam_,eps_);
226 if ( useSecant_ ) secant_->applyH(*Hvp.get(0),*dwa_);
227 else obj_->precond(*Hvp.get(0),*dwa_,*x_,tol);
228 bnd_->pruneActive(*Hvp.get(0),*xlam_,eps_);
229 dwa_->set(*vp.get(0));
230 bnd_->pruneInactive(*dwa_,*xlam_,eps_);
231 Hvp.get(0)->plus(dwa_->dual());
232 Hvp.get(1)->set(vp.get(1)->dual());
233 }
234 };
235
236 using TypeB::Algorithm<Real>::status_;
237 using TypeB::Algorithm<Real>::state_;
238 using TypeB::Algorithm<Real>::proj_;
239
240 void initialize(Vector<Real> &x,
241 const Vector<Real> &g,
242 Objective<Real> &obj,
244 std::ostream &outStream = std::cout);
245
246public:
247
248 PrimalDualActiveSetAlgorithm(ParameterList &list, const Ptr<Secant<Real>> &secant = nullPtr);
249
250 using TypeB::Algorithm<Real>::run;
251 void run( Vector<Real> &x,
252 const Vector<Real> &g,
253 Objective<Real> &obj,
255 std::ostream &outStream = std::cout) override;
256
257 void writeHeader( std::ostream& os ) const override;
258
259 void writeName( std::ostream& os ) const override;
260
261 void writeOutput( std::ostream& os, bool write_header = false ) const override;
262
263}; // class ROL::TypeB::PrimalDualActiveSetAlgorithm
264
265} // namespace TypeB
266} // namespace ROL
267
269
270#endif
Provides the interface to apply upper and lower bound constraints.
Defines the general constraint operator interface.
Provides the interface to apply a linear operator.
Provides the interface to evaluate objective functions.
Defines the linear algebra of vector space on a generic partitioned vector.
ROL::Ptr< const Vector< Real > > get(size_type i) const
Provides interface for and implements limited-memory secant operators.
Definition: ROL_Secant.hpp:79
Provides an interface to run bound constrained optimization algorithms.
Ptr< PolyhedralProjection< Real > > proj_
const Ptr< AlgorithmState< Real > > state_
const Ptr< CombinedStatusTest< Real > > status_
HessianPDAS_Poly(const Ptr< Objective< Real > > &obj, const Ptr< BoundConstraint< Real > > &bnd, const Ptr< Constraint< Real > > &con, const Ptr< const Vector< Real > > &x, const Ptr< const Vector< Real > > &xlam, Real eps, const Ptr< Secant< Real > > &secant, bool useSecant, const Ptr< Vector< Real > > &pwa, const Ptr< Vector< Real > > &dwa)
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
HessianPDAS(const Ptr< Objective< Real > > &obj, const Ptr< BoundConstraint< Real > > &bnd, const Ptr< const Vector< Real > > &x, const Ptr< const Vector< Real > > &xlam, Real eps, const Ptr< Secant< Real > > &secant, bool useSecant, const Ptr< Vector< Real > > &pwa)
PrecondPDAS_Poly(const Ptr< Objective< Real > > &obj, const Ptr< BoundConstraint< Real > > &bnd, const Ptr< const Vector< Real > > &x, const Ptr< const Vector< Real > > &xlam, Real eps, const Ptr< Secant< Real > > &secant, bool useSecant, const Ptr< Vector< Real > > &dwa)
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
void applyInverse(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply inverse of linear operator.
PrecondPDAS(const Ptr< Objective< Real > > &obj, const Ptr< BoundConstraint< Real > > &bnd, const Ptr< const Vector< Real > > &x, const Ptr< const Vector< Real > > &xlam, Real eps, const Ptr< Secant< Real > > &secant, bool useSecant, const Ptr< Vector< Real > > &dwa)
void apply(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply linear operator.
void applyInverse(Vector< Real > &Hv, const Vector< Real > &v, Real &tol) const
Apply inverse of linear operator.
Provides an interface to run the projected secant algorithm.
int flagKrylov_
Termination flag for Krylov method (used for inexact Newton)
int iterKrylov_
Number of Krylov iterations (used for inexact Newton)
Real stol_
PDAS minimum step size stopping tolerance (default: 1e-8)
void writeName(std::ostream &os) const override
Print step name.
Real scale_
Scale for dual variables in the active set, (default: 1)
bool feasible_
Flag whether the current iterate is feasible or not.
Ptr< Secant< Real > > secant_
Secant object (used for quasi-Newton)
Real atolKrylov_
Absolute tolerance for Krylov solve (default: 1e-4)
int totalKrylov_
Total number of Krylov iterations per PDAS iteration.
Real gtol_
PDAS gradient stopping tolerance (default: 1e-6)
void initialize(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &bnd, std::ostream &outStream=std::cout)
void writeOutput(std::ostream &os, bool write_header=false) const override
Print iterate status.
void run(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &bnd, std::ostream &outStream=std::cout) override
Run algorithm on bound constrained problems (Type-B). This general interface supports the use of dual...
int maxit_
Maximum number of PDAS steps (default: 10)
Ptr< Krylov< Real > > krylov_
Krylov solver object (used for inexact Newton)
void writeHeader(std::ostream &os) const override
Print iterate header.
int maxitKrylov_
Maximum number of Krylov iterations (default: 100)
bool useSecantHessVec_
Whether or not to use to a secant approximation as the Hessian.
Real rtolKrylov_
Relative tolerance for Krylov solve (default: 1e-2)
bool useSecantPrecond_
Whether or not to use a secant approximation to precondition inexact Newton.
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 const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis,...
Definition: ROL_Vector.hpp:226
virtual void plus(const Vector &x)=0
Compute , where .