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