ROL
function/test_14.cpp
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
49#include "ROL_HS41.hpp"
50#include "ROL_HS53.hpp"
51#include "ROL_HS55.hpp"
52#include "ROL_Bounds.hpp"
54
55#include "ROL_Stream.hpp"
56#include "Teuchos_GlobalMPISession.hpp"
57
58typedef double RealT;
59
60int main(int argc, char *argv[]) {
61
62 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
63
64 // This little trick lets us print to std::cout only if a
65 // (dummy) command-line argument is provided.
66 int iprint = argc - 1;
67 ROL::Ptr<std::ostream> outStream;
68 ROL::nullstream bhs; // outputs nothing
69 if (iprint > 0)
70 outStream = ROL::makePtrFromRef(std::cout);
71 else
72 outStream = ROL::makePtrFromRef(bhs);
73
74 int errorFlag = 0;
75
76 try {
77 RealT tol = std::sqrt(ROL::ROL_EPSILON<RealT>());
78 RealT cnorm(0);
79 ROL::Ptr<ROL::Vector<RealT>> sol, mul, x, lam, l, u, c;
80 ROL::Ptr<ROL::Objective<RealT>> obj;
81 ROL::Ptr<ROL::Constraint<RealT>> con;
82 ROL::Ptr<ROL::BoundConstraint<RealT>> bnd;
83 ROL::Ptr<ROL::PolyhedralProjection<RealT>> proj;
84 ROL::ParameterList list;
85 list.sublist("General").set("Output Level",2);
86 std::vector<RealT> data;
87
88 *outStream << std::endl << "Hock and Schittkowski Problem #41" << std::endl << std::endl;
90 obj = HS41.getObjective();
91 sol = HS41.getInitialGuess();
92 con = HS41.getEqualityConstraint();
93 mul = HS41.getEqualityMultiplier();
94 bnd = HS41.getBoundConstraint();
95
96 lam = mul->clone(); lam->set(*mul);
97 x = sol->clone(); x->set(*sol);
98 l = sol->clone(); l->zero();
99 u = sol->clone(); u->setScalar(static_cast<RealT>(1));
100 c = mul->dual().clone();
101
102 list.sublist("General").sublist("Polyhedral Projection").set("Type","Dai-Fletcher");
103 proj = ROL::PolyhedralProjectionFactory<RealT>(*sol,sol->dual(),bnd,con,*lam,*c,list);
104 proj->project(*x,*outStream);
105
106 con->value(*c,*x,tol);
107 cnorm = c->norm();
108
109 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(sol)->getVector();
110 *outStream << " Initial: x1 = " << data[0] << " x2 = " << data[1]
111 << " x3 = " << data[2] << std::endl;
112 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(x)->getVector();
113 *outStream << " Result: x1 = " << data[0] << " x2 = " << data[1]
114 << " x3 = " << data[2] << std::endl;
115 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(lam)->getVector();
116 *outStream << " Multiplier: l1 = " << data[0] << std::endl;
117
118 *outStream << std::endl;
119 *outStream << " is equality feasible = " << (cnorm<=tol) << std::endl
120 << " are bounds feasible = " << bnd->isFeasible(*x) << std::endl;
121
122 errorFlag += !bnd->isFeasible(*x);
123 errorFlag += (cnorm > tol);
124
125 *outStream << std::endl << "Hock and Schittkowski Problem #41" << std::endl << std::endl;
127 obj = HS41a.getObjective();
128 sol = HS41a.getInitialGuess();
129 con = HS41a.getEqualityConstraint();
130 mul = HS41a.getEqualityMultiplier();
131 bnd = HS41a.getBoundConstraint();
132
133 lam = mul->clone(); lam->set(*mul);
134 x = sol->clone(); x->set(*sol);
135 l = sol->clone(); l->zero();
136 u = sol->clone(); u->setScalar(static_cast<RealT>(1));
137 c = mul->dual().clone();
138
139 list.sublist("General").sublist("Polyhedral Projection").set("Type","Ridders");
140 proj = ROL::PolyhedralProjectionFactory<RealT>(*sol,sol->dual(),bnd,con,*lam,*c,list);
141 proj->project(*x,*outStream);
142
143 con->value(*c,*x,tol);
144 cnorm = c->norm();
145
146 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(sol)->getVector();
147 *outStream << " Initial: x1 = " << data[0] << " x2 = " << data[1]
148 << " x3 = " << data[2] << std::endl;
149 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(x)->getVector();
150 *outStream << " Result: x1 = " << data[0] << " x2 = " << data[1]
151 << " x3 = " << data[2] << std::endl;
152 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(lam)->getVector();
153 *outStream << " Multiplier: l1 = " << data[0] << std::endl;
154
155 *outStream << std::endl;
156 *outStream << " is equality feasible = " << (cnorm<=tol) << std::endl
157 << " are bounds feasible = " << bnd->isFeasible(*x) << std::endl;
158
159 errorFlag += !bnd->isFeasible(*x);
160 errorFlag += (cnorm > tol);
161
162 *outStream << std::endl << "Hock and Schittkowski Problem #41" << std::endl << std::endl;
164 obj = HS41b.getObjective();
165 sol = HS41b.getInitialGuess();
166 con = HS41b.getEqualityConstraint();
167 mul = HS41b.getEqualityMultiplier();
168 bnd = HS41b.getBoundConstraint();
169
170 lam = mul->clone(); lam->set(*mul);
171 x = sol->clone(); x->set(*sol);
172 l = sol->clone(); l->zero();
173 u = sol->clone(); u->setScalar(static_cast<RealT>(1));
174 c = mul->dual().clone();
175
176 list.sublist("General").sublist("Polyhedral Projection").set("Type","Brents");
177 proj = ROL::PolyhedralProjectionFactory<RealT>(*sol,sol->dual(),bnd,con,*lam,*c,list);
178 proj->project(*x,*outStream);
179
180 con->value(*c,*x,tol);
181 cnorm = c->norm();
182
183 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(sol)->getVector();
184 *outStream << " Initial: x1 = " << data[0] << " x2 = " << data[1]
185 << " x3 = " << data[2] << std::endl;
186 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(x)->getVector();
187 *outStream << " Result: x1 = " << data[0] << " x2 = " << data[1]
188 << " x3 = " << data[2] << std::endl;
189 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(lam)->getVector();
190 *outStream << " Multiplier: l1 = " << data[0] << std::endl;
191
192 *outStream << std::endl;
193 *outStream << " is equality feasible = " << (cnorm<=tol) << std::endl
194 << " are bounds feasible = " << bnd->isFeasible(*x) << std::endl;
195
196 errorFlag += !bnd->isFeasible(*x);
197 errorFlag += (cnorm > tol);
198
199 *outStream << std::endl << "Hock and Schittkowski Problem #53" << std::endl << std::endl;
201 obj = HS53.getObjective();
202 sol = HS53.getInitialGuess();
203 con = HS53.getEqualityConstraint();
204 mul = HS53.getEqualityMultiplier();
205 bnd = HS53.getBoundConstraint();
206
207 lam = mul->clone(); lam->set(*mul);
208 x = sol->clone(); x->set(*sol);
209 l = sol->clone(); l->zero();
210 u = sol->clone(); u->setScalar(static_cast<RealT>(1));
211 c = mul->dual().clone();
212
213 list.sublist("General").sublist("Polyhedral Projection").set("Type","Dykstra");
214 proj = ROL::PolyhedralProjectionFactory<RealT>(*sol,sol->dual(),bnd,con,*lam,*c,list);
215 proj->project(*x,*outStream);
216
217 con->value(*c,*x,tol);
218 cnorm = c->norm();
219
220 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(sol)->getVector();
221 *outStream << " Initial: x1 = " << data[0] << " x2 = " << data[1]
222 << " x3 = " << data[2] << " x4 = " << data[3]
223 << " x5 = " << data[4] << std::endl;
224 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(x)->getVector();
225 *outStream << " Result: x1 = " << data[0] << " x2 = " << data[1]
226 << " x3 = " << data[2] << " x4 = " << data[3]
227 << " x5 = " << data[4] << std::endl;
228 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(lam)->getVector();
229 *outStream << " Multiplier: l1 = " << data[0] << " l2 = " << data[1]
230 << " l3 = " << data[2] << std::endl;
231
232 *outStream << std::endl;
233 *outStream << " is equality feasible = " << (cnorm<=tol) << std::endl
234 << " are bounds feasible = " << bnd->isFeasible(*x) << std::endl;
235
236 errorFlag += !bnd->isFeasible(*x);
237 errorFlag += (cnorm > tol);
238
239 *outStream << std::endl << "Hock and Schittkowski Problem #53" << std::endl << std::endl;
241 obj = HS53a.getObjective();
242 sol = HS53a.getInitialGuess();
243 con = HS53a.getEqualityConstraint();
244 mul = HS53a.getEqualityMultiplier();
245 bnd = HS53a.getBoundConstraint();
246
247 lam = mul->clone(); lam->set(*mul);
248 x = sol->clone(); x->set(*sol);
249 l = sol->clone(); l->zero();
250 u = sol->clone(); u->setScalar(static_cast<RealT>(1));
251 c = mul->dual().clone();
252
253 list.sublist("General").sublist("Polyhedral Projection").set("Type","Douglas-Rachford");
254 proj = ROL::PolyhedralProjectionFactory<RealT>(*sol,sol->dual(),bnd,con,*lam,*c,list);
255 proj->project(*x,*outStream);
256
257 con->value(*c,*x,tol);
258 cnorm = c->norm();
259
260 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(sol)->getVector();
261 *outStream << " Initial: x1 = " << data[0] << " x2 = " << data[1]
262 << " x3 = " << data[2] << " x4 = " << data[3]
263 << " x5 = " << data[4] << std::endl;
264 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(x)->getVector();
265 *outStream << " Result: x1 = " << data[0] << " x2 = " << data[1]
266 << " x3 = " << data[2] << " x4 = " << data[3]
267 << " x5 = " << data[4] << std::endl;
268 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(lam)->getVector();
269 *outStream << " Multiplier: l1 = " << data[0] << " l2 = " << data[1]
270 << " l3 = " << data[2] << std::endl;
271
272 *outStream << std::endl;
273 *outStream << " is equality feasible = " << (cnorm<=tol) << std::endl
274 << " are bounds feasible = " << bnd->isFeasible(*x) << std::endl;
275
276 errorFlag += !bnd->isFeasible(*x);
277 errorFlag += (cnorm > tol);
278
279 *outStream << std::endl << "Hock and Schittkowski Problem #55" << std::endl << std::endl;
281 obj = HS55.getObjective();
282 sol = HS55.getInitialGuess();
283 con = HS55.getEqualityConstraint();
284 mul = HS55.getEqualityMultiplier();
285 bnd = HS55.getBoundConstraint();
286
287 //ROL::Ptr<ROL::OptimizationProblem<RealT>> problem;
288 //ROL::Ptr<ROL::Vector<RealT>> xt;
289 //std::vector<ROL::Ptr<ROL::Vector<RealT>>> xv;
290 //HS55.get(problem,xt,xv);
291 //problem->check(*outStream);
292
293 lam = mul->clone(); lam->set(*mul);
294 x = sol->clone(); x->set(*sol);
295 l = sol->clone(); l->zero();
296 u = sol->clone(); u->setScalar(static_cast<RealT>(1));
297 c = mul->dual().clone();
298
299 list.sublist("General").sublist("Polyhedral Projection").set("Type","Semismooth Newton");
300 proj = ROL::PolyhedralProjectionFactory<RealT>(*sol,sol->dual(),bnd,con,*lam,*c,list);
301 proj->project(*x,*outStream);
302
303 con->value(*c,*x,tol);
304 cnorm = c->norm();
305
306 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(sol)->getVector();
307 *outStream << " Initial: x1 = " << data[0] << " x2 = " << data[1]
308 << " x3 = " << data[2] << " x4 = " << data[3]
309 << " x5 = " << data[4] << " x6 = " << data[5] << std::endl;
310 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(x)->getVector();
311 *outStream << " Result: x1 = " << data[0] << " x2 = " << data[1]
312 << " x3 = " << data[2] << " x4 = " << data[3]
313 << " x5 = " << data[4] << " x6 = " << data[5] << std::endl;
314 data = *ROL::staticPtrCast<ROL::StdVector<RealT>>(lam)->getVector();
315 *outStream << " Multiplier: l1 = " << data[0] << " l2 = " << data[1]
316 << " l3 = " << data[2] << " l4 = " << data[3]
317 << " l5 = " << data[4] << " l6 = " << data[5] << std::endl;
318
319 *outStream << std::endl;
320 *outStream << " is equality feasible = " << (cnorm<=tol) << std::endl
321 << " are bounds feasible = " << bnd->isFeasible(*x) << std::endl;
322 *outStream << std::endl;
323
324 errorFlag += !bnd->isFeasible(*x);
325 errorFlag += (cnorm > tol);
326 }
327
328 catch (std::logic_error& err) {
329 *outStream << err.what() << "\n";
330 errorFlag = -1000;
331 }; // end try
332
333 if (errorFlag != 0)
334 std::cout << "End Result: TEST FAILED\n";
335 else
336 std::cout << "End Result: TEST PASSED\n";
337
338 return 0;
339}
340
Contains definitions for W. Hock and K. Schittkowski 41th test function.
Contains definitions for W. Hock and K. Schittkowski 53th test function.
Contains definitions for W. Hock and K. Schittkowski 55th test function.
Defines a no-output stream class ROL::NullStream and a function makeStreamPtr which either wraps a re...
Ptr< BoundConstraint< Real > > getBoundConstraint(void) const
Definition ROL_HS41.hpp:153
Ptr< Constraint< Real > > getEqualityConstraint(void) const
Definition ROL_HS41.hpp:144
Ptr< Vector< Real > > getInitialGuess(void) const
Definition ROL_HS41.hpp:131
Ptr< Objective< Real > > getObjective(void) const
Definition ROL_HS41.hpp:127
Ptr< Vector< Real > > getEqualityMultiplier(void) const
Definition ROL_HS41.hpp:148
Ptr< BoundConstraint< Real > > getBoundConstraint(void) const
Definition ROL_HS53.hpp:166
Ptr< Objective< Real > > getObjective(void) const
Definition ROL_HS53.hpp:137
Ptr< Constraint< Real > > getEqualityConstraint(void) const
Definition ROL_HS53.hpp:157
Ptr< Vector< Real > > getInitialGuess(void) const
Definition ROL_HS53.hpp:141
Ptr< Vector< Real > > getEqualityMultiplier(void) const
Definition ROL_HS53.hpp:161
Ptr< Constraint< Real > > getEqualityConstraint(void) const
Definition ROL_HS55.hpp:175
Ptr< BoundConstraint< Real > > getBoundConstraint(void) const
Definition ROL_HS55.hpp:184
Ptr< Vector< Real > > getEqualityMultiplier(void) const
Definition ROL_HS55.hpp:179
Ptr< Objective< Real > > getObjective(void) const
Definition ROL_HS55.hpp:147
Ptr< Vector< Real > > getInitialGuess(void) const
Definition ROL_HS55.hpp:151
int main(int argc, char *argv[])
double RealT