ROL
ROL_SimpleEqConstrained.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 
51 #ifndef ROL_SIMPLEEQCONSTRAINED_HPP
52 #define ROL_SIMPLEEQCONSTRAINED_HPP
53 
54 #include "ROL_Objective.hpp"
55 #include "ROL_StdVector.hpp"
57 #include "Teuchos_SerialDenseVector.hpp"
58 #include "Teuchos_SerialDenseSolver.hpp"
59 
60 namespace ROL {
61 namespace ZOO {
62 
66  template< class Real, class XPrim=StdVector<Real>, class XDual=StdVector<Real> >
68 
69  public:
71 
72  Real value( const Vector<Real> &x, Real &tol ) {
73  Teuchos::RCP<const std::vector<Real> > xp =
74  (Teuchos::dyn_cast<const XPrim>(x)).getVector();
75 
76  int n = xp->size();
77  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, objective value): "
78  "Primal vector x must be of length 5.");
79 
80  Real x1 = (*xp)[0];
81  Real x2 = (*xp)[1];
82  Real x3 = (*xp)[2];
83  Real x4 = (*xp)[3];
84  Real x5 = (*xp)[4];
85 
86  Real val = exp(x1*x2*x3*x4*x5) - 0.5 * pow( (pow(x1,3)+pow(x2,3)+1.0), 2);
87 
88  return val;
89  }
90 
91  void gradient( Vector<Real> &g, const Vector<Real> &x, Real &tol ) {
92  Teuchos::RCP<const std::vector<Real> > xp =
93  (Teuchos::dyn_cast<const XPrim>(x)).getVector();
94  Teuchos::RCP<std::vector<Real> > gp =
95  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(g)).getVector());
96 
97  int n = xp->size();
98  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, objective gradient): "
99  " Primal vector x must be of length 5.");
100 
101  n = gp->size();
102  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, objective gradient): "
103  "Gradient vector g must be of length 5.");
104 
105  Real x1 = (*xp)[0];
106  Real x2 = (*xp)[1];
107  Real x3 = (*xp)[2];
108  Real x4 = (*xp)[3];
109  Real x5 = (*xp)[4];
110 
111  Real expxi = exp(x1*x2*x3*x4*x5);
112 
113  (*gp)[0] = x2*x3*x4*x5 * expxi - 3*pow(x1,2) * (pow(x1,3) + pow(x2,3) + 1);
114  (*gp)[1] = x1*x3*x4*x5 * expxi - 3*pow(x2,2) * (pow(x1,3) + pow(x2,3) + 1);
115  (*gp)[2] = x1*x2*x4*x5 * expxi;
116  (*gp)[3] = x1*x2*x3*x5 * expxi;
117  (*gp)[4] = x1*x2*x3*x4 * expxi;
118  }
119 
120  void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
121  Teuchos::RCP<const std::vector<Real> > xp =
122  (Teuchos::dyn_cast<const XPrim>(x)).getVector();
123  Teuchos::RCP<const std::vector<Real> > vp =
124  (Teuchos::dyn_cast<const XPrim>(v)).getVector();
125  Teuchos::RCP<std::vector<Real> > hvp =
126  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(hv)).getVector());
127 
128  int n = xp->size();
129  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, objective hessVec): "
130  "Primal vector x must be of length 5.");
131 
132  n = vp->size();
133  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, objective hessVec): "
134  "Input vector v must be of length 5.");
135 
136  n = hvp->size();
137  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, objective hessVec): "
138  "Output vector hv must be of length 5.");
139 
140  Real x1 = (*xp)[0];
141  Real x2 = (*xp)[1];
142  Real x3 = (*xp)[2];
143  Real x4 = (*xp)[3];
144  Real x5 = (*xp)[4];
145 
146  Real v1 = (*vp)[0];
147  Real v2 = (*vp)[1];
148  Real v3 = (*vp)[2];
149  Real v4 = (*vp)[3];
150  Real v5 = (*vp)[4];
151 
152  Real expxi = exp(x1*x2*x3*x4*x5);
153 
154  (*hvp)[0] = ( pow(x2,2)*pow(x3,2)*pow(x4,2)*pow(x5,2)*expxi-9.0*pow(x1,4)-6.0*(pow(x1,3)+pow(x2,3)+1.0)*x1 ) * v1 +
155  ( x3*x4*x5*expxi+x2*pow(x3,2)*pow(x4,2)*pow(x5,2)*x1*expxi-9.0*pow(x2,2)*pow(x1,2) ) * v2 +
156  ( x2*x4*x5*expxi+pow(x2,2)*x3*pow(x4,2)*pow(x5,2)*x1*expxi ) * v3 +
157  ( x2*x3*x5*expxi+pow(x2,2)*pow(x3,2)*x4*pow(x5,2)*x1*expxi ) * v4 +
158  ( x2*x3*x4*expxi+pow(x2,2)*pow(x3,2)*pow(x4,2)*x5*x1*expxi ) * v5;
159 
160  (*hvp)[1] = ( x3*x4*x5*expxi+x2*pow(x3,2)*pow(x4,2)*pow(x5,2)*x1*expxi-9.0*pow(x2,2)*pow(x1,2) ) * v1 +
161  ( pow(x1,2)*pow(x3,2)*pow(x4,2)*pow(x5,2)*expxi-9.0*pow(x2,4)-6.0*(pow(x1,3)+pow(x2,3)+1.0)*x2 ) * v2 +
162  ( x1*x4*x5*expxi+pow(x1,2)*x3*pow(x4,2)*pow(x5,2)*x2*expxi ) * v3 +
163  ( x1*x3*x5*expxi+pow(x1,2)*pow(x3,2)*x4*pow(x5,2)*x2*expxi ) * v4 +
164  ( x1*x3*x4*expxi+pow(x1,2)*pow(x3,2)*pow(x4,2)*x5*x2*expxi ) * v5;
165 
166  (*hvp)[2] = ( x2*x4*x5*expxi+pow(x2,2)*x3*pow(x4,2)*pow(x5,2)*x1*expxi ) * v1 +
167  ( x1*x4*x5*expxi+pow(x1,2)*x3*pow(x4,2)*pow(x5,2)*x2*expxi ) * v2 +
168  ( pow(x1,2)*pow(x2,2)*pow(x4,2)*pow(x5,2)*expxi ) * v3 +
169  ( x1*x2*x5*expxi+pow(x1,2)*pow(x2,2)*x4*pow(x5,2)*x3*expxi ) * v4 +
170  ( x1*x2*x4*expxi+pow(x1,2)*pow(x2,2)*pow(x4,2)*x5*x3*expxi ) * v5;
171 
172  (*hvp)[3] = ( x2*x3*x5*expxi+pow(x2,2)*pow(x3,2)*x4*pow(x5,2)*x1*expxi ) * v1 +
173  ( x1*x3*x5*expxi+pow(x1,2)*pow(x3,2)*x4*pow(x5,2)*x2*expxi ) * v2 +
174  ( x1*x2*x5*expxi+pow(x1,2)*pow(x2,2)*x4*pow(x5,2)*x3*expxi ) * v3 +
175  ( pow(x1,2)*pow(x2,2)*pow(x3,2)*pow(x5,2)*expxi ) * v4 +
176  ( x1*x2*x3*expxi+pow(x1,2)*pow(x2,2)*pow(x3,2)*x5*x4*expxi ) * v5;
177 
178  (*hvp)[4] = ( x2*x3*x4*expxi+pow(x2,2)*pow(x3,2)*pow(x4,2)*x5*x1*expxi ) * v1 +
179  ( x1*x3*x4*expxi+pow(x1,2)*pow(x3,2)*pow(x4,2)*x5*x2*expxi ) * v2 +
180  ( x1*x2*x4*expxi+pow(x1,2)*pow(x2,2)*pow(x4,2)*x5*x3*expxi ) * v3 +
181  ( x1*x2*x3*expxi+pow(x1,2)*pow(x2,2)*pow(x3,2)*x5*x4*expxi ) * v4 +
182  ( pow(x1,2)*pow(x2,2)*pow(x3,2)*pow(x4,2)*expxi ) * v5;
183  }
184 
185  };
186 
187 
193  template<class Real, class XPrim=StdVector<Real>, class XDual=StdVector<Real>, class CPrim=StdVector<Real>, class CDual=StdVector<Real> >
195 
196  public:
198 
199  void value( Vector<Real> &c, const Vector<Real> &x, Real &tol ) {
200  Teuchos::RCP<const std::vector<Real> > xp =
201  (Teuchos::dyn_cast<const XPrim>(x)).getVector();
202  Teuchos::RCP<std::vector<Real> > cp =
203  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<CPrim>(c)).getVector());
204 
205  int n = xp->size();
206  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint value): "
207  "Primal vector x must be of length 5.");
208 
209  int m = cp->size();
210  TEUCHOS_TEST_FOR_EXCEPTION( (m != 3), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint value): "
211  "Constraint vector c must be of length 3.");
212 
213  Real x1 = (*xp)[0];
214  Real x2 = (*xp)[1];
215  Real x3 = (*xp)[2];
216  Real x4 = (*xp)[3];
217  Real x5 = (*xp)[4];
218 
219  (*cp)[0] = x1*x1+x2*x2+x3*x3+x4*x4+x5*x5 - 10.0;
220  (*cp)[1] = x2*x3 - 5.0*x4*x5;
221  (*cp)[2] = x1*x1*x1 + x2*x2*x2 + 1.0;
222  }
223 
224  void applyJacobian( Vector<Real> &jv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
225  Teuchos::RCP<const std::vector<Real> > xp =
226  (Teuchos::dyn_cast<const XPrim>(x)).getVector();
227  Teuchos::RCP<const std::vector<Real> > vp =
228  (Teuchos::dyn_cast<const XPrim>(v)).getVector();
229  Teuchos::RCP<std::vector<Real> > jvp =
230  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<CPrim>(jv)).getVector());
231 
232  int n = xp->size();
233  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyJacobian): "
234  "Primal vector x must be of length 5.");
235 
236  int d = vp->size();
237  TEUCHOS_TEST_FOR_EXCEPTION( (d != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyJacobian): "
238  "Input vector v must be of length 5.");
239  d = jvp->size();
240  TEUCHOS_TEST_FOR_EXCEPTION( (d != 3), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyJacobian): "
241  "Output vector jv must be of length 3.");
242 
243  Real x1 = (*xp)[0];
244  Real x2 = (*xp)[1];
245  Real x3 = (*xp)[2];
246  Real x4 = (*xp)[3];
247  Real x5 = (*xp)[4];
248 
249  Real v1 = (*vp)[0];
250  Real v2 = (*vp)[1];
251  Real v3 = (*vp)[2];
252  Real v4 = (*vp)[3];
253  Real v5 = (*vp)[4];
254 
255  (*jvp)[0] = 2.0*(x1*v1+x2*v2+x3*v3+x4*v4+x5*v5);
256  (*jvp)[1] = x3*v2+x2*v3-5.0*x5*v4-5.0*x4*v5;
257  (*jvp)[2] = 3.0*x1*x1*v1+3.0*x2*x2*v2;
258 
259  } //applyJacobian
260 
261  void applyAdjointJacobian( Vector<Real> &ajv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
262  Teuchos::RCP<const std::vector<Real> > xp =
263  (Teuchos::dyn_cast<const XPrim>(x)).getVector();
264  Teuchos::RCP<const std::vector<Real> > vp =
265  (Teuchos::dyn_cast<const CDual>(v)).getVector();
266  Teuchos::RCP<std::vector<Real> > ajvp =
267  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(ajv)).getVector());
268 
269  int n = xp->size();
270  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyAdjointJacobian): "
271  "Primal vector x must be of length 5.");
272 
273  int d = vp->size();
274  TEUCHOS_TEST_FOR_EXCEPTION( (d != 3), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyAdjointJacobian): "
275  "Input vector v must be of length 3.");
276 
277  d = ajvp->size();
278  TEUCHOS_TEST_FOR_EXCEPTION( (d != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyAdjointJacobian): "
279  "Output vector ajv must be of length 5.");
280 
281  Real x1 = (*xp)[0];
282  Real x2 = (*xp)[1];
283  Real x3 = (*xp)[2];
284  Real x4 = (*xp)[3];
285  Real x5 = (*xp)[4];
286 
287  Real v1 = (*vp)[0];
288  Real v2 = (*vp)[1];
289  Real v3 = (*vp)[2];
290 
291  (*ajvp)[0] = 2.0*x1*v1+3.0*x1*x1*v3;
292  (*ajvp)[1] = 2.0*x2*v1+x3*v2+3.0*x2*x2*v3;
293  (*ajvp)[2] = 2.0*x3*v1+x2*v2;
294  (*ajvp)[3] = 2.0*x4*v1-5.0*x5*v2;
295  (*ajvp)[4] = 2.0*x5*v1-5.0*x4*v2;
296 
297  } //applyAdjointJacobian
298 
299  void applyAdjointHessian( Vector<Real> &ahuv, const Vector<Real> &u, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
300  Teuchos::RCP<const std::vector<Real> > xp =
301  (Teuchos::dyn_cast<const XPrim>(x)).getVector();
302  Teuchos::RCP<const std::vector<Real> > up =
303  (Teuchos::dyn_cast<const CDual>(u)).getVector();
304  Teuchos::RCP<const std::vector<Real> > vp =
305  (Teuchos::dyn_cast<const XPrim>(v)).getVector();
306  Teuchos::RCP<std::vector<Real> > ahuvp =
307  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(ahuv)).getVector());
308 
309  int n = xp->size();
310  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyAdjointHessian): "
311  "Primal vector x must be of length 5.");
312 
313  n = vp->size();
314  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyAdjointHessian): "
315  "Direction vector v must be of length 5.");
316 
317  n = ahuvp->size();
318  TEUCHOS_TEST_FOR_EXCEPTION( (n != 5), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyAdjointHessian): "
319  "Output vector ahuv must be of length 5.");
320  int d = up->size();
321  TEUCHOS_TEST_FOR_EXCEPTION( (d != 3), std::invalid_argument, ">>> ERROR (ROL_SimpleEqConstrained, constraint applyAdjointHessian): "
322  "Dual constraint vector u must be of length 3.");
323 
324  Real x1 = (*xp)[0];
325  Real x2 = (*xp)[1];
326 
327  Real v1 = (*vp)[0];
328  Real v2 = (*vp)[1];
329  Real v3 = (*vp)[2];
330  Real v4 = (*vp)[3];
331  Real v5 = (*vp)[4];
332 
333  Real u1 = (*up)[0];
334  Real u2 = (*up)[1];
335  Real u3 = (*up)[2];
336 
337  (*ahuvp)[0] = 2.0*u1*v1 + 6.0*u3*x1*v1;
338  (*ahuvp)[1] = 2.0*u1*v2 + u2*v3 + 6.0*u3*x2*v2;
339  (*ahuvp)[2] = 2.0*u1*v3 + u2*v2;
340  (*ahuvp)[3] = 2.0*u1*v4 - 5.0*u2*v5;
341  (*ahuvp)[4] = 2.0*u1*v5 - 5.0*u2*v4;
342 
343  } //applyAdjointHessian
344 
345  /*std::vector<Real> solveAugmentedSystem(Vector<Real> &v1, Vector<Real> &v2, const Vector<Real> &b1, const Vector<Real> &b2, const Vector<Real> &x, Real &tol) {
346  Teuchos::RCP<std::vector<Real> > v1p =
347  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XPrim>(v1)).getVector());
348  Teuchos::RCP<std::vector<Real> > v2p =
349  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<CDual>(v2)).getVector());
350  Teuchos::RCP<const std::vector<Real> > b1p =
351  (Teuchos::dyn_cast<XDual>(const_cast<Vector<Real> &>(b1))).getVector();
352  Teuchos::RCP<const std::vector<Real> > b2p =
353  (Teuchos::dyn_cast<CPrim>(const_cast<Vector<Real> &>(b2))).getVector();
354  Teuchos::RCP<const std::vector<Real> > xp =
355  (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(x))).getVector();
356 
357  Real x1 = (*xp)[0];
358  Real x2 = (*xp)[1];
359  Real x3 = (*xp)[2];
360  Real x4 = (*xp)[3];
361  Real x5 = (*xp)[4];
362 
363  int i = 0;
364 
365  // Assemble augmented system.
366  Teuchos::SerialDenseMatrix<int, Real> augmat(8,8);
367  for (i=0; i<5; i++) {
368  augmat(i,i) = 1.0;
369  }
370  augmat(5,0) = 2.0*x1; augmat(5,1) = 2.0*x2; augmat(5,2) = 2.0*x3; augmat(5,3) = 2.0*x4; augmat(5,4) = 2.0*x5;
371  augmat(6,1) = x3; augmat(6,2) = x2; augmat(6,3) = -5.0*x5; augmat(6,4) = -5.0*x4;
372  augmat(7,0) = 3.0*x1*x1; augmat(7,1) = 3.0*x2*x2;
373  augmat(0,5) = 2.0*x1; augmat(1,5) = 2.0*x2; augmat(2,5) = 2.0*x3; augmat(3,5) = 2.0*x4; augmat(4,5) = 2.0*x5;
374  augmat(1,6) = x3; augmat(2,6) = x2; augmat(3,6) = -5.0*x5; augmat(4,6) = -5.0*x4;
375  augmat(0,7) = 3.0*x1*x1; augmat(1,7) = 3.0*x2*x2;
376  Teuchos::SerialDenseVector<int, Real> lhs(8);
377  Teuchos::SerialDenseVector<int, Real> rhs(8);
378  for (i=0; i<5; i++) {
379  rhs(i) = (*b1p)[i];
380  }
381  for (i=5; i<8; i++) {
382  rhs(i) = (*b2p)[i-5];
383  }
384 
385  // Solve augmented system.
386  Teuchos::SerialDenseSolver<int, Real> augsolver;
387  augsolver.setMatrix(Teuchos::rcp(&augmat, false));
388  augsolver.setVectors(Teuchos::rcp(&lhs, false), Teuchos::rcp(&rhs, false));
389  augsolver.solve();
390 
391  // Retrieve solution.
392  for (i=0; i<5; i++) {
393  (*v1p)[i] = lhs(i);
394  }
395  for (i=0; i<3; i++) {
396  (*v2p)[i] = lhs(i+5);
397  }
398 
399  return std::vector<Real>(0);
400 
401  }*/ //solveAugmentedSystem
402 
403  };
404 
405 
406  template<class Real, class XPrim, class XDual, class CPrim, class CDual>
407  void getSimpleEqConstrained( Teuchos::RCP<Objective<Real> > &obj,
408  Teuchos::RCP<EqualityConstraint<Real> > &constr,
409  Vector<Real> &x0,
410  Vector<Real> &sol ) {
411  // Cast initial guess and solution vectors.
412  Teuchos::RCP<std::vector<Real> > x0p =
413  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XPrim>(x0)).getVector());
414  Teuchos::RCP<std::vector<Real> > solp =
415  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XPrim>(sol)).getVector());
416  int n = 5;
417  // Resize vectors.
418  x0p->resize(n);
419  solp->resize(n);
420  // Instantiate objective function.
421  obj = Teuchos::rcp( new Objective_SimpleEqConstrained<Real, XPrim, XDual> );
422  // Instantiate constraints.
424  // later we will bundle equality constraints into constraints ...
425  //std::vector<Teuchos::RCP<EqualityConstraint<Real> > > eqc( 1, Teuchos::rcp( new EqualityConstraint_SimpleEqConstrained<Real> ) );
426  //constr = Teuchos::rcp( new Constraints<Real>(eqc) );
427 
428  // Get initial guess.
429  (*x0p)[0] = -1.8;
430  (*x0p)[1] = 1.7;
431  (*x0p)[2] = 1.9;
432  (*x0p)[3] = -0.8;
433  (*x0p)[4] = -0.8;
434  // Get solution.
435  (*solp)[0] = -1.717143570394391e+00;
436  (*solp)[1] = 1.595709690183565e+00;
437  (*solp)[2] = 1.827245752927178e+00;
438  (*solp)[3] = -7.636430781841294e-01;
439  (*solp)[4] = -7.636430781841294e-01;
440  }
441 
442 } // End ZOO Namespace
443 } // End ROL Namespace
444 
445 #endif
Provides the interface to evaluate objective functions.
Equality constraints c_i(x) = 0, where: c1(x) = x1^2+x2^2+x3^2+x4^2+x5^2 - 10 c2(x) = x2*x3-5*x4*x5 c...
void value(Vector< Real > &c, const Vector< Real > &x, Real &tol)
Evaluate the constraint operator at .
void applyJacobian(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply the constraint Jacobian at , , to vector .
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:72
void applyAdjointJacobian(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply the adjoint of the the constraint Jacobian at , , to vector .
Defines the equality constraint operator interface.
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
void getSimpleEqConstrained(Teuchos::RCP< Objective< Real > > &obj, Teuchos::RCP< EqualityConstraint< Real > > &constr, Vector< Real > &x0, Vector< Real > &sol)
void applyAdjointHessian(Vector< Real > &ahuv, const Vector< Real > &u, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply the derivative of the adjoint of the constraint Jacobian at to vector in direction ...
Real value(const Vector< Real > &x, Real &tol)
Compute value.
Objective function: f(x) = exp(x1*x2*x3*x4*x5) + 0.5*(x1^3+x2^3+1)^2.