ROL
gross-pitaevskii/example_03.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 
68 #include <iostream>
69 
70 #include "Teuchos_oblackholestream.hpp"
71 #include "Teuchos_GlobalMPISession.hpp"
72 #include "Teuchos_XMLParameterListHelpers.hpp"
73 
74 #include "ROL_StdVector.hpp"
75 #include "ROL_Objective.hpp"
77 #include "ROL_CompositeStepSQP.hpp"
78 #include "ROL_Algorithm.hpp"
79 
80 #include "numerics/NodalBasis.hpp"
82 
83 
84 
85 using namespace ROL;
86 
87 template <class Real, class Element=Real>
88 class OptStdVector; // Optimization space.
89 
90 template <class Real, class Element=Real>
91 class OptDualStdVector; // Dual optimization space.
92 
93 template <class Real, class Element=Real>
94 class ConStdVector; // Constraint space.
95 
96 template <class Real, class Element=Real>
97 class ConDualStdVector; // Dual constraint space.
98 
99 // Vector space definitions:
100 
101 
102 // Optimization space.
103 template <class Real, class Element>
104 class OptStdVector : public Vector<Real> {
105 
106 private:
107 Teuchos::RCP<std::vector<Element> > std_vec_;
108 mutable Teuchos::RCP<OptDualStdVector<Real> > dual_vec_;
109 bool useRiesz_;
110 Teuchos::RCP<InnerProductMatrix<Real> > ipmat_;
111 
112 public:
113 
114 OptStdVector(const Teuchos::RCP<std::vector<Element> > & std_vec, bool useRiesz, Teuchos::RCP<InnerProductMatrix<Real> > ipmat) :
115  std_vec_(std_vec), dual_vec_(Teuchos::null), useRiesz_(useRiesz), ipmat_(ipmat) {}
116 
117 void plus( const Vector<Real> &x ) {
118  OptStdVector &ex = Teuchos::dyn_cast<OptStdVector>(const_cast <Vector<Real> &>(x));
119  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
120  unsigned dimension = std_vec_->size();
121  for (unsigned i=0; i<dimension; i++) {
122  (*std_vec_)[i] += (*xvalptr)[i];
123  }
124 }
125 
126 void scale( const Real alpha ) {
127  unsigned dimension = std_vec_->size();
128  for (unsigned i=0; i<dimension; i++) {
129  (*std_vec_)[i] *= alpha;
130  }
131 }
132 
133 
135 Real dot( const Vector<Real> &x ) const {
136  Real val = 0;
137  OptStdVector<Real, Element> & ex = Teuchos::dyn_cast<OptStdVector<Real, Element> >(const_cast <Vector<Real> &>(x));
138  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
139  unsigned dimension = std_vec_->size();
140  if(useRiesz_){
141  val = this->ipmat_->inner(std_vec_,xvalptr);
142  }
143  else{
144  for (unsigned i=0; i<dimension; i++) {
145  val += (*std_vec_)[i]*(*xvalptr)[i];
146  }
147  }
148  return val;
149 }
150 
151 Real norm() const {
152  Real val = 0;
153  val = std::sqrt( dot(*this) );
154  return val;
155 }
156 
157 Teuchos::RCP<Vector<Real> > clone() const {
158  return Teuchos::rcp( new OptStdVector( Teuchos::rcp( new std::vector<Element>(std_vec_->size()) ),useRiesz_,ipmat_ ) );
159 }
160 
161 Teuchos::RCP<const std::vector<Element> > getVector() const {
162  return std_vec_;
163 }
164 
165 Teuchos::RCP<Vector<Real> > basis( const int i ) const {
166  Teuchos::RCP<OptStdVector> e = Teuchos::rcp( new OptStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size(), 0.0)),useRiesz_, ipmat_ ) );
167  (const_cast <std::vector<Element> &> (*e->getVector()))[i]= 1.0;
168  return e;
169 }
170 
171 int dimension() const {return std_vec_->size();}
172 
173 
175 const Vector<Real> & dual() const {
176 
177  dual_vec_ = Teuchos::rcp( new OptDualStdVector<Real>( Teuchos::rcp( new std::vector<Element>(*std_vec_) ), useRiesz_, ipmat_ ) );
178 
179  if(useRiesz_){
180  this->ipmat_->apply(std_vec_, Teuchos::rcp_const_cast<std::vector<Real> >(dual_vec_->getVector()));
181  }
182  return *dual_vec_;
183 }
184 
185 }; // class OptStdVector
186 
187 
188 
189 
190 // Dual optimization space.
191 template <class Real, class Element>
192 class OptDualStdVector : public Vector<Real> {
193 
194 private:
195 Teuchos::RCP<std::vector<Element> > std_vec_;
196 mutable Teuchos::RCP<OptStdVector<Real> > dual_vec_;
198 Teuchos::RCP<InnerProductMatrix<Real> >ipmat_;
199 
200 public:
201 
202 OptDualStdVector(const Teuchos::RCP<std::vector<Element> > & std_vec, bool useRiesz, Teuchos::RCP<InnerProductMatrix<Real> > ipmat) :
203  std_vec_(std_vec), dual_vec_(Teuchos::null), useRiesz_(useRiesz), ipmat_(ipmat) {}
204 
205 void plus( const Vector<Real> &x ) {
206  OptDualStdVector &ex = Teuchos::dyn_cast<OptDualStdVector>(const_cast <Vector<Real> &>(x));
207  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
208  unsigned dimension = std_vec_->size();
209  for (unsigned i=0; i<dimension; i++) {
210  (*std_vec_)[i] += (*xvalptr)[i];
211  }
212 }
213 
214 void scale( const Real alpha ) {
215  unsigned dimension = std_vec_->size();
216  for (unsigned i=0; i<dimension; i++) {
217  (*std_vec_)[i] *= alpha;
218  }
219 }
220 
221 Real dot( const Vector<Real> &x ) const {
222  Real val = 0;
223  OptDualStdVector<Real, Element> & ex = Teuchos::dyn_cast<OptDualStdVector<Real, Element> >(const_cast <Vector<Real> &>(x));
224  Teuchos::RCP<const std::vector<Element> > kxvalptr = ex.getVector();
225  unsigned dimension = std_vec_->size();
226  if(useRiesz_)
227  {
228  val = ipmat_->inv_inner(std_vec_,kxvalptr);
229  }
230  else {
231  for (unsigned i=0; i<dimension; i++) {
232  val += (*std_vec_)[i]*(*kxvalptr)[i];
233  }
234 
235  }
236  return val;
237 }
238 
239 Real norm() const {
240  Real val = 0;
241  val = std::sqrt( dot(*this) );
242  return val;
243 }
244 
245 Teuchos::RCP<Vector<Real> > clone() const {
246  return Teuchos::rcp( new OptDualStdVector( Teuchos::rcp( new std::vector<Element>(std_vec_->size()) ), useRiesz_, ipmat_ ) );
247 }
248 
249 Teuchos::RCP<const std::vector<Element> > getVector() const {
250  return std_vec_;
251 }
252 
253 Teuchos::RCP<Vector<Real> > basis( const int i ) const {
254  Teuchos::RCP<OptDualStdVector> e = Teuchos::rcp( new OptDualStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size(), 0.0)),useRiesz_, ipmat_ ) );
255  (const_cast <std::vector<Element> &> (*e->getVector()))[i]= 1.0;
256  return e;
257 }
258 
259 int dimension() const {return std_vec_->size();}
260 
261 const Vector<Real> & dual() const {
262  dual_vec_ = Teuchos::rcp( new OptStdVector<Real>( Teuchos::rcp( new std::vector<Element>(*std_vec_) ), useRiesz_, ipmat_ ) );
263  if(useRiesz_){
264  this->ipmat_->solve(std_vec_, Teuchos::rcp_const_cast<std::vector<Real> >(dual_vec_->getVector()));
265  }
266  return *dual_vec_;
267 }
268 
269 }; // class OptDualStdVector
270 
271 
272 
273 
274 
275 // Constraint space.
276 template <class Real, class Element>
277 class ConStdVector : public Vector<Real> {
278 
279 private:
280 Teuchos::RCP<std::vector<Element> > std_vec_;
281 mutable Teuchos::RCP<ConDualStdVector<Real> > dual_vec_;
282 
283 public:
284 
285 ConStdVector(const Teuchos::RCP<std::vector<Element> > & std_vec) : std_vec_(std_vec), dual_vec_(Teuchos::null) {}
286 
287 void plus( const Vector<Real> &x ) {
288  ConStdVector &ex = Teuchos::dyn_cast<ConStdVector>(const_cast <Vector<Real> &>(x));
289  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
290  unsigned dimension = std_vec_->size();
291  for (unsigned i=0; i<dimension; i++) {
292  (*std_vec_)[i] += (*xvalptr)[i];
293  }
294 }
295 
296 void scale( const Real alpha ) {
297  unsigned dimension = std_vec_->size();
298  for (unsigned i=0; i<dimension; i++) {
299  (*std_vec_)[i] *= alpha;
300  }
301 }
302 
303 Real dot( const Vector<Real> &x ) const {
304  Real val = 0;
305  ConStdVector<Real, Element> & ex = Teuchos::dyn_cast<ConStdVector<Real, Element> >(const_cast <Vector<Real> &>(x));
306  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
307 
308 
309 
310  unsigned dimension = std_vec_->size();
311  for (unsigned i=0; i<dimension; i++) {
312  val += (*std_vec_)[i]*(*xvalptr)[i];
313  }
314  return val;
315 }
316 
317 Real norm() const {
318  Real val = 0;
319  val = std::sqrt( dot(*this) );
320  return val;
321 }
322 
323 Teuchos::RCP<Vector<Real> > clone() const {
324  return Teuchos::rcp( new ConStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size())) ) );
325 }
326 
327 Teuchos::RCP<const std::vector<Element> > getVector() const {
328  return std_vec_;
329 }
330 
331 Teuchos::RCP<Vector<Real> > basis( const int i ) const {
332  Teuchos::RCP<ConStdVector> e = Teuchos::rcp( new ConStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size(), 0.0)) ) );
333  (const_cast <std::vector<Element> &> (*e->getVector()))[i]= 1.0;
334  return e;
335 }
336 
337 int dimension() const {return std_vec_->size();}
338 
339 const Vector<Real> & dual() const {
340  dual_vec_ = Teuchos::rcp( new ConDualStdVector<Real>( Teuchos::rcp( new std::vector<Element>(*std_vec_) ) ) );
341  return *dual_vec_;
342 }
343 
344 }; // class ConStdVector
345 
346 
347 // Dual constraint space.
348 template <class Real, class Element>
349 class ConDualStdVector : public Vector<Real> {
350 private:
351 
352 Teuchos::RCP<std::vector<Element> > std_vec_;
353 mutable Teuchos::RCP<ConStdVector<Real> > dual_vec_;
354 
355 public:
356 
357 ConDualStdVector(const Teuchos::RCP<std::vector<Element> > & std_vec) : std_vec_(std_vec), dual_vec_(Teuchos::null) {}
358 
359 void plus( const Vector<Real> &x ) {
360  ConDualStdVector &ex = Teuchos::dyn_cast<ConDualStdVector>(const_cast <Vector<Real> &>(x));
361  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
362  unsigned dimension = std_vec_->size();
363  for (unsigned i=0; i<dimension; i++) {
364  (*std_vec_)[i] += (*xvalptr)[i];
365  }
366 }
367 
368 void scale( const Real alpha ) {
369  unsigned dimension = std_vec_->size();
370  for (unsigned i=0; i<dimension; i++) {
371  (*std_vec_)[i] *= alpha;
372  }
373 }
374 
375 Real dot( const Vector<Real> &x ) const {
376  Real val = 0;
377  ConDualStdVector<Real, Element> & ex = Teuchos::dyn_cast<ConDualStdVector<Real, Element> >(const_cast <Vector<Real> &>(x));
378  Teuchos::RCP<const std::vector<Element> > xvalptr = ex.getVector();
379  unsigned dimension = std_vec_->size();
380  for (unsigned i=0; i<dimension; i++) {
381  val += (*std_vec_)[i]*(*xvalptr)[i];
382  }
383  return val;
384 }
385 
386 Real norm() const {
387  Real val = 0;
388  val = std::sqrt( dot(*this) );
389  return val;
390 }
391 
392 Teuchos::RCP<Vector<Real> > clone() const {
393  return Teuchos::rcp( new ConDualStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size())) ) );
394 }
395 
396 Teuchos::RCP<const std::vector<Element> > getVector() const {
397  return std_vec_;
398 }
399 
400 Teuchos::RCP<Vector<Real> > basis( const int i ) const {
401  Teuchos::RCP<ConDualStdVector> e = Teuchos::rcp( new ConDualStdVector( Teuchos::rcp(new std::vector<Element>(std_vec_->size(), 0.0)) ) );
402  (const_cast <std::vector<Element> &> (*e->getVector()))[i]= 1.0;
403  return e;
404 }
405 
406 int dimension() const {return std_vec_->size();}
407 
408 const Vector<Real> & dual() const {
409  dual_vec_ = Teuchos::rcp( new ConStdVector<Real>( Teuchos::rcp( new std::vector<Element>(*std_vec_) ) ) );
410  return *dual_vec_;
411 }
412 
413 }; // class ConDualStdVector
414 
415 
416 /*** End of declaration of four vector spaces. ***/
417 
418 
419 
420 
421 
422 
424 template<class Real, class XPrim=StdVector<Real>, class XDual=StdVector<Real> >
425 class Objective_GrossPitaevskii : public Objective<Real> {
426 
427  private:
428 
431  const int ni_;
432  const Real gnl_;
433  Teuchos::RCP<NodalBasis<Real> > nb_;
434 
435  Teuchos::RCP<InnerProductMatrix<Real> > kinetic_;
436  Teuchos::RCP<InnerProductMatrix<Real> > potential_;
437  Teuchos::RCP<InnerProductMatrix<Real> > nonlinear_;
438 
439  void updateNonlinear(Teuchos::RCP<const std::vector<Real> > psip) {
440 
441  // This should be a member variable instead of reallocating each time
442  Teuchos::RCP<std::vector<Real> > psi2q =
443  Teuchos::rcp( new std::vector<Real> (nb_->nq_, 0.0) );
444 
445  nb_->lagrange_->dinterp(*psip,*psi2q);
446 
447  for(int i=0;i<nb_->nq_;++i) {
448  (*psi2q)[i] *= (*psi2q)[i];
449  }
450  nonlinear_->update(*psi2q);
451  }
452 
453  public:
454 
455  Objective_GrossPitaevskii(const int ni, const Real gnl,
456  Teuchos::RCP<NodalBasis<Real> > nb,
457  Teuchos::RCP<InnerProductMatrix<Real> > kinetic,
458  Teuchos::RCP<InnerProductMatrix<Real> > potential,
459  Teuchos::RCP<InnerProductMatrix<Real> > nonlinear):
460  ni_(ni), gnl_(gnl), nb_(nb), kinetic_(kinetic),
461  potential_(potential), nonlinear_(nonlinear) {}
462 
463 
464  Real value( const Vector<Real> &psi, Real &tol ) {
465 
466  // Pointer to opt vector
467  Teuchos::RCP<const std::vector<Real> > psip =
468  (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(psi))).getVector();
469 
470  this->updateNonlinear(psip);
471 
472  // Compute objective function value
473  Real J = 0.5*( kinetic_->inner(psip,psip) +
474  potential_->inner(psip,psip) +
475  gnl_*nonlinear_->inner(psip,psip) );
476 
477  return J;
478  }
479 
480 
481  void gradient( Vector<Real> &g, const Vector<Real> &psi, Real &tol ) {
482 
483  // Pointer to opt vector
484  Teuchos::RCP<const std::vector<Real> > psip =
485  (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(psi))).getVector();
486 
487  // Pointer to gradient vector
488  Teuchos::RCP<std::vector<Real> > gp =
489  Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(g)).getVector());
490 
491  this->updateNonlinear(psip);
492 
493  kinetic_->apply(psip,gp);
494  potential_->applyadd(psip,gp);
495  nonlinear_->applyaddtimes(psip,gp,2*gnl_);
496  }
497 
498  void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &psi, Real &tol ) {
499  // Pointer to opt vector
500  Teuchos::RCP<const std::vector<Real> > psip =
501  (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(psi))).getVector();
502 
503  // Pointer to direction vector
504  Teuchos::RCP<const std::vector<Real> > vp = (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(v))).getVector();
505 
506  // Pointer to action of Hessian on direction vector
507  Teuchos::RCP<std::vector<Real> > hvp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(hv)).getVector());
508 
509  kinetic_->apply(vp,hvp);
510  potential_->applyadd(vp,hvp);
511  nonlinear_->applyaddtimes(vp,hvp,6*gnl_);
512 
513  }
514 
515 };
516 
517 
518 
520 template<class Real, class XPrim=StdVector<Real>, class XDual=StdVector<Real>, class CPrim=StdVector<Real>, class CDual=StdVector<Real> >
521 class Normalization_Constraint : public EqualityConstraint<Real> {
522 
523  private:
524  Teuchos::RCP<InnerProductMatrix<Real> > mass_;
525 
526  public:
528  mass_(mass) {}
529 
530  void value(Vector<Real> &c, const Vector<Real> &psi, Real &tol){
531 
532  // Pointer to constraint vector (only one element)
533  Teuchos::RCP<std::vector<Real> > cp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<CPrim>(c)).getVector());
534 
535  // Pointer to opt vector
536  Teuchos::RCP<const std::vector<Real> > psip =
537  (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(psi))).getVector();
538 
539  (*cp)[0] = mass_->inner(psip,psip)-1.0;
540 
541  }
542 
543  void applyJacobian(Vector<Real> &jv, const Vector<Real> &v, const Vector<Real> &psi, Real &tol){
544 
545  // Pointer to action of Jacobian of constraint on direction vector (yields scalar)
546  Teuchos::RCP<std::vector<Real> > jvp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<CPrim>(jv)).getVector());
547 
548  // Pointer to direction vector
549  Teuchos::RCP<const std::vector<Real> > vp = (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(v))).getVector();
550 
551  // Pointer to opt vector
552  Teuchos::RCP<const std::vector<Real> > psip =
553  (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(psi))).getVector();
554 
555  (*jvp)[0] = 2.0*mass_->inner(psip,vp);
556 
557  }
558 
559  void applyAdjointJacobian(Vector<Real> &ajv, const Vector<Real> &v, const Vector<Real> &psi, Real &tol){
560 
561  // Pointer to action of adjoint of Jacobian of constraint on direction vector (yields vector)
562  Teuchos::RCP<std::vector<Real> > ajvp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(ajv)).getVector());
563 
564  // Pointer to direction vector
565  Teuchos::RCP<const std::vector<Real> > vp = (Teuchos::dyn_cast<CDual>(const_cast<Vector<Real> &>(v))).getVector();
566 
567  // Pointer to opt vector
568  Teuchos::RCP<const std::vector<Real> > psip =
569  (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(psi))).getVector();
570 
571  mass_->apply(psip,ajvp);
572  for(unsigned i=0;i<psip->size();++i) {
573  (*ajvp)[i] *= 2.0*(*vp)[0];
574  }
575 
576  }
577 
579  const Vector<Real> &psi, Real &tol){
580 
581  // The pointer to action of constraint Hessian in u,v inner product
582  Teuchos::RCP<std::vector<Real> > ahuvp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<XDual>(ahuv)).getVector());
583 
584  // Pointer to direction vector u
585  Teuchos::RCP<const std::vector<Real> > up = (Teuchos::dyn_cast<CDual>(const_cast<Vector<Real> &>(u))).getVector();
586 
587  // Pointer to direction vector v
588  Teuchos::RCP<const std::vector<Real> > vp = (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(v))).getVector();
589 
590  // Pointer to opt vector
591  Teuchos::RCP<const std::vector<Real> > psip =
592  (Teuchos::dyn_cast<XPrim>(const_cast<Vector<Real> &>(psi))).getVector();
593 
594  mass_->apply(vp,ahuvp);
595  for(unsigned i=0;i<psip->size();++i) {
596  (*ahuvp)[i] *= 2.0*(*up)[0];
597  }
598  }
599 
600 };
601 
602 
Normalization_Constraint(Teuchos::RCP< InnerProductMatrix< Real > > mass)
Provides the interface to evaluate objective functions.
Real value(const Vector< Real > &psi, Real &tol)
Compute value.
void scale(const Real alpha)
Compute where .
void applyAdjointJacobian(Vector< Real > &ajv, const Vector< Real > &v, const Vector< Real > &psi, Real &tol)
Apply the adjoint of the the constraint Jacobian at , , to vector .
Teuchos::RCP< Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
ConStdVector(const Teuchos::RCP< std::vector< Element > > &std_vec)
Real norm() const
Returns where .
OptDualStdVector(const Teuchos::RCP< std::vector< Element > > &std_vec, bool useRiesz, Teuchos::RCP< InnerProductMatrix< Real > > ipmat)
void value(Vector< Real > &c, const Vector< Real > &psi, Real &tol)
Evaluate the constraint operator at .
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &psi, Real &tol)
Apply Hessian approximation to vector.
Real dot(const Vector< Real > &x) const
Compute where .
Teuchos::RCP< Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
Teuchos::RCP< const std::vector< Element > > getVector() const
Teuchos::RCP< InnerProductMatrix< Real > > kinetic_
Teuchos::RCP< Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
int dimension() const
Return dimension of the vector space.
Teuchos::RCP< const std::vector< Element > > getVector() const
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:72
const Vector< Real > & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Teuchos::RCP< InnerProductMatrix< Real > > nonlinear_
void applyJacobian(Vector< Real > &jv, const Vector< Real > &v, const Vector< Real > &psi, Real &tol)
Apply the constraint Jacobian at , , to vector .
const Vector< Real > & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Teuchos::RCP< Vector< Real > > basis(const int i) const
Return i-th basis vector.
Teuchos::RCP< Vector< Real > > basis(const int i) const
Return i-th basis vector.
Defines the equality constraint operator interface.
void scale(const Real alpha)
Compute where .
Teuchos::RCP< Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
Real dot(const Vector< Real > &x) const
Compute where .
void scale(const Real alpha)
Compute where .
Teuchos::RCP< InnerProductMatrix< Real > > ipmat_
Teuchos::RCP< const std::vector< Element > > getVector() const
Teuchos::RCP< NodalBasis< Real > > nb_
const Vector< Real > & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Objective_GrossPitaevskii(const int ni, const Real gnl, Teuchos::RCP< NodalBasis< Real > > nb, Teuchos::RCP< InnerProductMatrix< Real > > kinetic, Teuchos::RCP< InnerProductMatrix< Real > > potential, Teuchos::RCP< InnerProductMatrix< Real > > nonlinear)
int dimension() const
Return dimension of the vector space.
int dimension() const
Return dimension of the vector space.
Teuchos::RCP< InnerProductMatrix< Real > > ipmat_
Teuchos::RCP< Vector< Real > > basis(const int i) const
Return i-th basis vector.
Real norm() const
Returns where .
void plus(const Vector< Real > &x)
Compute , where .
int dimension() const
Return dimension of the vector space.
Real dot(const Vector< Real > &x) const
Compute where .
OptStdVector(const Teuchos::RCP< std::vector< Element > > &std_vec, bool useRiesz, Teuchos::RCP< InnerProductMatrix< Real > > ipmat)
void plus(const Vector< Real > &x)
Compute , where .
const Vector< Real > & dual() const
Modify the dual of vector u to be .
Real norm() const
Returns where .
Real dot(const Vector< Real > &x) const
Modify the dot product between primal variables to be .
ConDualStdVector(const Teuchos::RCP< std::vector< Element > > &std_vec)
Teuchos::RCP< InnerProductMatrix< Real > > mass_
Teuchos::RCP< InnerProductMatrix< Real > > potential_
void gradient(Vector< Real > &g, const Vector< Real > &psi, Real &tol)
Compute gradient.
Teuchos::RCP< const std::vector< Element > > getVector() const
void plus(const Vector< Real > &x)
Compute , where .
void updateNonlinear(Teuchos::RCP< const std::vector< Real > > psip)
Teuchos::RCP< Vector< Real > > basis(const int i) const
Return i-th basis vector.
void plus(const Vector< Real > &x)
Compute , where .
Real norm() const
Returns where .
void applyAdjointHessian(Vector< Real > &ahuv, const Vector< Real > &u, const Vector< Real > &v, const Vector< Real > &psi, Real &tol)
Apply the derivative of the adjoint of the constraint Jacobian at to vector in direction ...
void scale(const Real alpha)
Compute where .