ROL
dual-spaces/simple-eq-constr-1/example_01.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 
50 #include "ROL_CompositeStepSQP.hpp"
51 #include "ROL_Algorithm.hpp"
52 #include "Teuchos_oblackholestream.hpp"
53 #include "Teuchos_GlobalMPISession.hpp"
54 
55 #include <iostream>
56 
57 typedef double RealT;
58 
59 /*** Declare four vector spaces. ***/
60 
61 // Forward declarations:
62 
63 template <class Real, class Element=Real>
64 class OptStdVector; // Optimization space.
65 
66 template <class Real, class Element=Real>
67 class OptDualStdVector; // Dual optimization space.
68 
69 template <class Real, class Element=Real>
70 class ConStdVector; // Constraint space.
71 
72 template <class Real, class Element=Real>
73 class ConDualStdVector; // Dual constraint space.
74 
75 // Vector space definitions:
76 
77 // Optimization space.
78 template <class Real, class Element>
79 class OptStdVector : public ROL::Vector<Real> {
80 
81 private:
82 Teuchos::RCP<std::vector<Element> > std_vec_;
83 mutable Teuchos::RCP<OptDualStdVector<Real> > dual_vec_;
84 
85 public:
86 
87 OptStdVector(const Teuchos::RCP<std::vector<Element> > & std_vec) : std_vec_(std_vec), dual_vec_(Teuchos::null) {}
88 
89 void plus( const ROL::Vector<Real> &x ) {
90  OptStdVector &ex = Teuchos::dyn_cast<OptStdVector>(const_cast <ROL::Vector<Real> &>(x));
91  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
92  unsigned dimension = std_vec_->size();
93  for (unsigned i=0; i<dimension; i++) {
94  (*std_vec_)[i] += (*xvalptr)[i];
95  }
96 }
97 
98 void scale( const Real alpha ) {
99  unsigned dimension = std_vec_->size();
100  for (unsigned i=0; i<dimension; i++) {
101  (*std_vec_)[i] *= alpha;
102  }
103 }
104 
105 Real dot( const ROL::Vector<Real> &x ) const {
106  Real val = 0;
107  OptStdVector<Real, Element> & ex = Teuchos::dyn_cast<OptStdVector<Real, Element> >(const_cast <ROL::Vector<Real> &>(x));
108  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
109  unsigned dimension = std_vec_->size();
110  for (unsigned i=0; i<dimension; i++) {
111  val += (*std_vec_)[i]*(*xvalptr)[i];
112  }
113  return val;
114 }
115 
116 Real norm() const {
117  Real val = 0;
118  val = std::sqrt( dot(*this) );
119  return val;
120 }
121 
122 Teuchos::RCP<ROL::Vector<Real> > clone() const {
123  return Teuchos::rcp( new OptStdVector( Teuchos::rcp( new std::vector<Element>(std_vec_->size()) ) ) );
124 }
125 
126 Teuchos::RCP<const std::vector<Element> > getVector() const {
127  return std_vec_;
128 }
129 
130 Teuchos::RCP<ROL::Vector<Real> > basis( const int i ) const {
131  Teuchos::RCP<OptStdVector> e = Teuchos::rcp( new OptStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size(), 0.0)) ) );
132  (const_cast <std::vector<Element> &> (*e->getVector()))[i]= 1.0;
133  return e;
134 }
135 
136 int dimension() const {return std_vec_->size();}
137 
138 const ROL::Vector<Real> & dual() const {
139  dual_vec_ = Teuchos::rcp( new OptDualStdVector<Real>( Teuchos::rcp( new std::vector<Element>(*std_vec_) ) ) );
140  return *dual_vec_;
141 }
142 
143 }; // class OptStdVector
144 
145 
146 // Dual optimization space.
147 template <class Real, class Element>
148 class OptDualStdVector : public ROL::Vector<Real> {
149 
150 private:
151 Teuchos::RCP<std::vector<Element> > std_vec_;
152 mutable Teuchos::RCP<OptStdVector<Real> > dual_vec_;
153 
154 public:
155 
156 OptDualStdVector(const Teuchos::RCP<std::vector<Element> > & std_vec) : std_vec_(std_vec), dual_vec_(Teuchos::null) {}
157 
158 void plus( const ROL::Vector<Real> &x ) {
159  OptDualStdVector &ex = Teuchos::dyn_cast<OptDualStdVector>(const_cast <ROL::Vector<Real> &>(x));
160  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
161  unsigned dimension = std_vec_->size();
162  for (unsigned i=0; i<dimension; i++) {
163  (*std_vec_)[i] += (*xvalptr)[i];
164  }
165 }
166 
167 void scale( const Real alpha ) {
168  unsigned dimension = std_vec_->size();
169  for (unsigned i=0; i<dimension; i++) {
170  (*std_vec_)[i] *= alpha;
171  }
172 }
173 
174 Real dot( const ROL::Vector<Real> &x ) const {
175  Real val = 0;
176  OptDualStdVector<Real, Element> & ex = Teuchos::dyn_cast<OptDualStdVector<Real, Element> >(const_cast <ROL::Vector<Real> &>(x));
177  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
178  unsigned dimension = std_vec_->size();
179  for (unsigned i=0; i<dimension; i++) {
180  val += (*std_vec_)[i]*(*xvalptr)[i];
181  }
182  return val;
183 }
184 
185 Real norm() const {
186  Real val = 0;
187  val = std::sqrt( dot(*this) );
188  return val;
189 }
190 
191 Teuchos::RCP<ROL::Vector<Real> > clone() const {
192  return Teuchos::rcp( new OptDualStdVector( Teuchos::rcp( new std::vector<Element>(std_vec_->size()) ) ) );
193 }
194 
195 Teuchos::RCP<const std::vector<Element> > getVector() const {
196  return std_vec_;
197 }
198 
199 Teuchos::RCP<ROL::Vector<Real> > basis( const int i ) const {
200  Teuchos::RCP<OptDualStdVector> e = Teuchos::rcp( new OptDualStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size(), 0.0)) ) );
201  (const_cast <std::vector<Element> &> (*e->getVector()))[i]= 1.0;
202  return e;
203 }
204 
205 int dimension() const {return std_vec_->size();}
206 
207 const ROL::Vector<Real> & dual() const {
208  dual_vec_ = Teuchos::rcp( new OptStdVector<Real>( Teuchos::rcp( new std::vector<Element>(*std_vec_) ) ) );
209  return *dual_vec_;
210 }
211 
212 }; // class OptDualStdVector
213 
214 
215 // Constraint space.
216 template <class Real, class Element>
217 class ConStdVector : public ROL::Vector<Real> {
218 
219 private:
220 Teuchos::RCP<std::vector<Element> > std_vec_;
221 mutable Teuchos::RCP<ConDualStdVector<Real> > dual_vec_;
222 
223 public:
224 
225 ConStdVector(const Teuchos::RCP<std::vector<Element> > & std_vec) : std_vec_(std_vec), dual_vec_(Teuchos::null) {}
226 
227 void plus( const ROL::Vector<Real> &x ) {
228  ConStdVector &ex = Teuchos::dyn_cast<ConStdVector>(const_cast <ROL::Vector<Real> &>(x));
229  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
230  unsigned dimension = std_vec_->size();
231  for (unsigned i=0; i<dimension; i++) {
232  (*std_vec_)[i] += (*xvalptr)[i];
233  }
234 }
235 
236 void scale( const Real alpha ) {
237  unsigned dimension = std_vec_->size();
238  for (unsigned i=0; i<dimension; i++) {
239  (*std_vec_)[i] *= alpha;
240  }
241 }
242 
243 Real dot( const ROL::Vector<Real> &x ) const {
244  Real val = 0;
245  ConStdVector<Real, Element> & ex = Teuchos::dyn_cast<ConStdVector<Real, Element> >(const_cast <ROL::Vector<Real> &>(x));
246  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
247  unsigned dimension = std_vec_->size();
248  for (unsigned i=0; i<dimension; i++) {
249  val += (*std_vec_)[i]*(*xvalptr)[i];
250  }
251  return val;
252 }
253 
254 Real norm() const {
255  Real val = 0;
256  val = std::sqrt( dot(*this) );
257  return val;
258 }
259 
260 Teuchos::RCP<ROL::Vector<Real> > clone() const {
261  return Teuchos::rcp( new ConStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size())) ) );
262 }
263 
264 Teuchos::RCP<const std::vector<Element> > getVector() const {
265  return std_vec_;
266 }
267 
268 Teuchos::RCP<ROL::Vector<Real> > basis( const int i ) const {
269  Teuchos::RCP<ConStdVector> e = Teuchos::rcp( new ConStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size(), 0.0)) ) );
270  (const_cast <std::vector<Element> &> (*e->getVector()))[i]= 1.0;
271  return e;
272 }
273 
274 int dimension() const {return std_vec_->size();}
275 
276 const ROL::Vector<Real> & dual() const {
277  dual_vec_ = Teuchos::rcp( new ConDualStdVector<Real>( Teuchos::rcp( new std::vector<Element>(*std_vec_) ) ) );
278  return *dual_vec_;
279 }
280 
281 }; // class ConStdVector
282 
283 
284 // Dual constraint space.
285 template <class Real, class Element>
286 class ConDualStdVector : public ROL::Vector<Real> {
287 private:
288 
289 Teuchos::RCP<std::vector<Element> > std_vec_;
290 mutable Teuchos::RCP<ConStdVector<Real> > dual_vec_;
291 
292 public:
293 
294 ConDualStdVector(const Teuchos::RCP<std::vector<Element> > & std_vec) : std_vec_(std_vec), dual_vec_(Teuchos::null) {}
295 
296 void plus( const ROL::Vector<Real> &x ) {
297  ConDualStdVector &ex = Teuchos::dyn_cast<ConDualStdVector>(const_cast <ROL::Vector<Real> &>(x));
298  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
299  unsigned dimension = std_vec_->size();
300  for (unsigned i=0; i<dimension; i++) {
301  (*std_vec_)[i] += (*xvalptr)[i];
302  }
303 }
304 
305 void scale( const Real alpha ) {
306  unsigned dimension = std_vec_->size();
307  for (unsigned i=0; i<dimension; i++) {
308  (*std_vec_)[i] *= alpha;
309  }
310 }
311 
312 Real dot( const ROL::Vector<Real> &x ) const {
313  Real val = 0;
314  ConDualStdVector<Real, Element> & ex = Teuchos::dyn_cast<ConDualStdVector<Real, Element> >(const_cast <ROL::Vector<Real> &>(x));
315  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
316  unsigned dimension = std_vec_->size();
317  for (unsigned i=0; i<dimension; i++) {
318  val += (*std_vec_)[i]*(*xvalptr)[i];
319  }
320  return val;
321 }
322 
323 Real norm() const {
324  Real val = 0;
325  val = std::sqrt( dot(*this) );
326  return val;
327 }
328 
329 Teuchos::RCP<ROL::Vector<Real> > clone() const {
330  return Teuchos::rcp( new ConDualStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size())) ) );
331 }
332 
333 Teuchos::RCP<const std::vector<Element> > getVector() const {
334  return std_vec_;
335 }
336 
337 Teuchos::RCP<ROL::Vector<Real> > basis( const int i ) const {
338  Teuchos::RCP<ConDualStdVector> e = Teuchos::rcp( new ConDualStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size(), 0.0)) ) );
339  (const_cast <std::vector<Element> &> (*e->getVector()))[i]= 1.0;
340  return e;
341 }
342 
343 int dimension() const {return std_vec_->size();}
344 
345 const ROL::Vector<Real> & dual() const {
346  dual_vec_ = Teuchos::rcp( new ConStdVector<Real>( Teuchos::rcp( new std::vector<Element>(*std_vec_) ) ) );
347  return *dual_vec_;
348 }
349 
350 }; // class ConDualStdVector
351 
352 /*** End of declaration of four vector spaces. ***/
353 
354 
355 int main(int argc, char *argv[]) {
356 
357  Teuchos::GlobalMPISession mpiSession(&argc, &argv);
358 
359  // This little trick lets us print to std::cout only if a (dummy) command-line argument is provided.
360  int iprint = argc - 1;
361  Teuchos::RCP<std::ostream> outStream;
362  Teuchos::oblackholestream bhs; // outputs nothing
363  if (iprint > 0)
364  outStream = Teuchos::rcp(&std::cout, false);
365  else
366  outStream = Teuchos::rcp(&bhs, false);
367 
368  int errorFlag = 0;
369 
370  // *** Example body.
371 
372  try {
373 
374  Teuchos::RCP<ROL::Objective<RealT> > obj;
375  Teuchos::RCP<ROL::EqualityConstraint<RealT> > constr;
376  Teuchos::RCP<std::vector<RealT> > x_rcp = Teuchos::rcp( new std::vector<RealT> (0, 0.0) );
377  Teuchos::RCP<std::vector<RealT> > sol_rcp = Teuchos::rcp( new std::vector<RealT> (0, 0.0) );
378  OptStdVector<RealT> x(x_rcp); // Iteration vector.
379  OptStdVector<RealT> sol(sol_rcp); // Reference solution vector.
380 
381  // Retrieve objective, constraint, iteration vector, solution vector.
382  ROL::ZOO::getSimpleEqConstrained <RealT, OptStdVector<RealT>, OptDualStdVector<RealT>, ConStdVector<RealT>, ConDualStdVector<RealT> > (obj, constr, x, sol);
383 
384  Teuchos::ParameterList parlist;
385  // Define Step
386  parlist.set("Nominal SQP Optimality Solver Tolerance", 1.e-2);
387  ROL::CompositeStepSQP<RealT> step(parlist);
388 
389  // Run derivative checks, etc.
390  int dim = 5;
391  int nc = 3;
392  RealT left = -1e0, right = 1e0;
393  Teuchos::RCP<std::vector<RealT> > xtest_rcp = Teuchos::rcp( new std::vector<RealT> (dim, 0.0) );
394  Teuchos::RCP<std::vector<RealT> > g_rcp = Teuchos::rcp( new std::vector<RealT> (dim, 0.0) );
395  Teuchos::RCP<std::vector<RealT> > d_rcp = Teuchos::rcp( new std::vector<RealT> (dim, 0.0) );
396  Teuchos::RCP<std::vector<RealT> > gd_rcp = Teuchos::rcp( new std::vector<RealT> (dim, 0.0) );
397  Teuchos::RCP<std::vector<RealT> > v_rcp = Teuchos::rcp( new std::vector<RealT> (dim, 0.0) );
398  Teuchos::RCP<std::vector<RealT> > vc_rcp = Teuchos::rcp( new std::vector<RealT> (nc, 0.0) );
399  Teuchos::RCP<std::vector<RealT> > vl_rcp = Teuchos::rcp( new std::vector<RealT> (nc, 0.0) );
400  OptStdVector<RealT> xtest(xtest_rcp);
401  OptDualStdVector<RealT> g(g_rcp);
402  OptStdVector<RealT> d(d_rcp);
403  OptDualStdVector<RealT> gd(gd_rcp);
404  OptStdVector<RealT> v(v_rcp);
405  ConStdVector<RealT> vc(vc_rcp);
406  ConDualStdVector<RealT> vl(vl_rcp);
407  // set xtest, d, v
408  for (int i=0; i<dim; i++) {
409  (*xtest_rcp)[i] = ( (RealT)rand() / (RealT)RAND_MAX ) * (right - left) + left;
410  (*d_rcp)[i] = ( (RealT)rand() / (RealT)RAND_MAX ) * (right - left) + left;
411  (*gd_rcp)[i] = ( (RealT)rand() / (RealT)RAND_MAX ) * (right - left) + left;
412  (*v_rcp)[i] = ( (RealT)rand() / (RealT)RAND_MAX ) * (right - left) + left;
413  }
414  // set vc, vl
415  for (int i=0; i<nc; i++) {
416  (*vc_rcp)[i] = ( (RealT)rand() / (RealT)RAND_MAX ) * (right - left) + left;
417  (*vl_rcp)[i] = ( (RealT)rand() / (RealT)RAND_MAX ) * (right - left) + left;
418  }
419  obj->checkGradient(xtest, g, d, true, *outStream); *outStream << "\n";
420  obj->checkHessVec(xtest, g, v, true, *outStream); *outStream << "\n";
421  obj->checkHessSym(xtest, g, d, v, true, *outStream); *outStream << "\n";
422  constr->checkApplyJacobian(xtest, v, vc, true, *outStream); *outStream << "\n";
423  constr->checkApplyAdjointJacobian(xtest, vl, vc, g, true, *outStream); *outStream << "\n";
424  constr->checkApplyAdjointHessian(xtest, vl, d, g, true, *outStream); *outStream << "\n";
425 
426  Teuchos::RCP<std::vector<RealT> > v1_rcp = Teuchos::rcp( new std::vector<RealT> (dim, 0.0) );
427  Teuchos::RCP<std::vector<RealT> > v2_rcp = Teuchos::rcp( new std::vector<RealT> (nc, 0.0) );
428  OptStdVector<RealT> v1(v1_rcp);
429  ConDualStdVector<RealT> v2(v2_rcp);
430  RealT augtol = 1e-8;
431  constr->solveAugmentedSystem(v1, v2, gd, vc, xtest, augtol);
432 
433  // Define Status Test
434  RealT gtol = 1e-12; // norm of gradient tolerance
435  RealT ctol = 1e-12; // norm of constraint tolerance
436  RealT stol = 1e-18; // norm of step tolerance
437  int maxit = 1000; // maximum number of iterations
438  ROL::StatusTestSQP<RealT> status(gtol, ctol, stol, maxit);
439 
440  // Define Algorithm
441  ROL::DefaultAlgorithm<RealT> algo(step, status, false);
442 
443  // Run Algorithm
444  vl.zero();
445  //(*x_rcp)[0] = 3.0; (*x_rcp)[1] = 2.0; (*x_rcp)[2] = 2.0; (*x_rcp)[3] = 1.0; (*x_rcp)[4] = 1.0;
446  //(*x_rcp)[0] = -5.0; (*x_rcp)[1] = -5.0; (*x_rcp)[2] = -5.0; (*x_rcp)[3] = -6.0; (*x_rcp)[4] = -6.0;
447 
448  std::vector<std::string> output = algo.run(x, g, vl, vc, *obj, *constr, false);
449  for ( unsigned i = 0; i < output.size(); i++ ) {
450  *outStream << output[i];
451  }
452 
453  // Compute Error
454  *outStream << "\nReference solution x_r =\n";
455  *outStream << std::scientific << " " << (*sol_rcp)[0] << "\n";
456  *outStream << std::scientific << " " << (*sol_rcp)[1] << "\n";
457  *outStream << std::scientific << " " << (*sol_rcp)[2] << "\n";
458  *outStream << std::scientific << " " << (*sol_rcp)[3] << "\n";
459  *outStream << std::scientific << " " << (*sol_rcp)[4] << "\n";
460  *outStream << "\nOptimal solution x =\n";
461  *outStream << std::scientific << " " << (*x_rcp)[0] << "\n";
462  *outStream << std::scientific << " " << (*x_rcp)[1] << "\n";
463  *outStream << std::scientific << " " << (*x_rcp)[2] << "\n";
464  *outStream << std::scientific << " " << (*x_rcp)[3] << "\n";
465  *outStream << std::scientific << " " << (*x_rcp)[4] << "\n";
466  x.axpy(-1.0, sol);
467  RealT abserr = x.norm();
468  RealT relerr = abserr/sol.norm();
469  *outStream << std::scientific << "\n Absolute Error: " << abserr;
470  *outStream << std::scientific << "\n Relative Error: " << relerr << "\n";
471  if ( relerr > sqrt(ROL::ROL_EPSILON) ) {
472  errorFlag += 1;
473  }
474  }
475  catch (std::logic_error err) {
476  *outStream << err.what() << "\n";
477  errorFlag = -1000;
478  }; // end try
479 
480  if (errorFlag != 0)
481  std::cout << "End Result: TEST FAILED\n";
482  else
483  std::cout << "End Result: TEST PASSED\n";
484 
485  return 0;
486 
487 }
488 
Teuchos::RCP< std::vector< Element > > std_vec_
Teuchos::RCP< OptDualStdVector< Real > > dual_vec_
void scale(const Real alpha)
Compute where .
OptStdVector(const Teuchos::RCP< std::vector< Element > > &std_vec)
Real dot(const ROL::Vector< Real > &x) const
Compute where .
Teuchos::RCP< ROL::Vector< Real > > basis(const int i) const
Return i-th basis vector.
ConStdVector(const Teuchos::RCP< std::vector< Element > > &std_vec)
const ROL::Vector< Real > & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
const ROL::Vector< Real > & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Teuchos::RCP< const std::vector< Element > > getVector() const
Teuchos::RCP< ROL::Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
int dimension() const
Return dimension of the vector space.
Teuchos::RCP< OptStdVector< Real > > dual_vec_
Teuchos::RCP< ROL::Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
Teuchos::RCP< std::vector< Element > > std_vec_
virtual void zero()
Set to zero vector.
Definition: ROL_Vector.hpp:155
Teuchos::RCP< const std::vector< Element > > getVector() const
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:72
Teuchos::RCP< ROL::Vector< Real > > basis(const int i) const
Return i-th basis vector.
void plus(const ROL::Vector< Real > &x)
Compute , where .
Teuchos::RCP< ConDualStdVector< Real > > dual_vec_
const ROL::Vector< Real > & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
const ROL::Vector< Real > & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
OptDualStdVector(const Teuchos::RCP< std::vector< Element > > &std_vec)
Teuchos::RCP< ROL::Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
void scale(const Real alpha)
Compute where .
void scale(const Real alpha)
Compute where .
virtual std::vector< std::string > run(Vector< Real > &x, Objective< Real > &obj, bool print=false, std::ostream &outStream=std::cout)
Run algorithm on unconstrained problems (Type-U). This is the primary Type-U interface.
Teuchos::RCP< const std::vector< Element > > getVector() const
Real dot(const ROL::Vector< Real > &x) const
Compute where .
Implements the computation of optimization steps with composite-step trust-region SQP methods...
int dimension() const
Return dimension of the vector space.
int dimension() const
Return dimension of the vector space.
Teuchos::RCP< ROL::Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
Teuchos::RCP< ROL::Vector< Real > > basis(const int i) const
Return i-th basis vector.
Teuchos::RCP< ROL::Vector< Real > > basis(const int i) const
Return i-th basis vector.
Contains definitions for the equality constrained NLP from Nocedal/Wright, 2nd edition, page 574, example 18.2; note the typo in reversing the initial guess and the solution.
int dimension() const
Return dimension of the vector space.
Teuchos::RCP< std::vector< Element > > std_vec_
Real norm() const
Returns where .
ConDualStdVector(const Teuchos::RCP< std::vector< Element > > &std_vec)
Teuchos::RCP< ConStdVector< Real > > dual_vec_
Real dot(const ROL::Vector< Real > &x) const
Compute where .
Real dot(const ROL::Vector< Real > &x) const
Compute where .
double RealT
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:194
Teuchos::RCP< const std::vector< Element > > getVector() const
void plus(const ROL::Vector< Real > &x)
Compute , where .
int main(int argc, char *argv[])
Teuchos::RCP< std::vector< Element > > std_vec_
void plus(const ROL::Vector< Real > &x)
Compute , where .
void plus(const ROL::Vector< Real > &x)
Compute , where .
Real norm() const
Returns where .
void scale(const Real alpha)
Compute where .
static const double ROL_EPSILON
Platform-dependent machine epsilon.
Definition: ROL_Types.hpp:115