ROL
ROL_PoissonInversion.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 
49 #ifndef USE_HESSVEC
50 #define USE_HESSVEC 1
51 #endif
52 
53 #ifndef ROL_POISSONINVERSION_HPP
54 #define ROL_POISSONINVERSION_HPP
55 
56 #include "ROL_StdVector.hpp"
57 #include "ROL_Objective.hpp"
58 #include "ROL_HelperFunctions.hpp"
59 
60 #include "Teuchos_LAPACK.hpp"
61 
62 namespace ROL {
63 namespace ZOO {
64 
67 template<class Real>
68 class Objective_PoissonInversion : public Objective<Real> {
69 
70  typedef std::vector<Real> vector;
71  typedef Vector<Real> V;
73 
74  typedef typename vector::size_type uint;
75 
76 private:
79 
80  Real hu_;
81  Real hz_;
82 
83  Real alpha_;
84 
85  Real eps_;
86  int reg_type_;
87 
88  Teuchos::RCP<const vector> getVector( const V& x ) {
89  using Teuchos::dyn_cast;
90  return dyn_cast<const SV>(x).getVector();
91  }
92 
93  Teuchos::RCP<vector> getVector( V& x ) {
94  using Teuchos::dyn_cast;
95  return dyn_cast<SV>(x).getVector();
96  }
97 
98 public:
99 
100  /* CONSTRUCTOR */
101  Objective_PoissonInversion(uint nz = 32, Real alpha = 1.e-4)
102  : nu_(nz-1), nz_(nz), hu_(1./((Real)nu_+1.)), hz_(hu_),
103  alpha_(alpha), eps_(1.e-4), reg_type_(2) {}
104 
105  /* REGULARIZATION DEFINITIONS */
106  Real reg_value(const Vector<Real> &z) {
107  using Teuchos::RCP;
108 
109  RCP<const vector> zp = getVector(z);
110 
111  Real val = 0.0;
112  for (uint i = 0; i < nz_; i++) {
113  if ( reg_type_ == 2 ) {
114  val += alpha_/2.0 * hz_ * (*zp)[i]*(*zp)[i];
115  }
116  else if ( reg_type_ == 1 ) {
117  val += alpha_ * hz_ * std::sqrt((*zp)[i]*(*zp)[i] + eps_);
118  }
119  else if ( reg_type_ == 0 ) {
120  if ( i < nz_-1 ) {
121  val += alpha_ * std::sqrt(std::pow((*zp)[i]-(*zp)[i+1],2.0)+eps_);
122  }
123  }
124  }
125  return val;
126  }
127 
129  using Teuchos::RCP;
130 
131  if ( reg_type_ == 2 ) {
132  g.set(z);
133  g.scale(alpha_*hz_);
134  }
135  else if ( reg_type_ == 1 ) {
136  RCP<const vector> zp = getVector(z);
137  RCP<vector > gp = getVector(g);
138 
139  for (uint i = 0; i < nz_; i++) {
140  (*gp)[i] = alpha_ * hz_ * (*zp)[i]/std::sqrt(std::pow((*zp)[i],2.0)+eps_);
141  }
142  }
143  else if ( reg_type_ == 0 ) {
144  RCP<const vector> zp = getVector(z);
145  RCP<vector> gp = getVector(g);
146 
147  Real diff = 0.0;
148  for (uint i = 0; i < nz_; i++) {
149  if ( i == 0 ) {
150  diff = (*zp)[i]-(*zp)[i+1];
151  (*gp)[i] = alpha_ * diff/std::sqrt(std::pow(diff,2.0)+eps_);
152  }
153  else if ( i == nz_-1 ) {
154  diff = (*zp)[i-1]-(*zp)[i];
155  (*gp)[i] = -alpha_ * diff/std::sqrt(std::pow(diff,2.0)+eps_);
156  }
157  else {
158  diff = (*zp)[i]-(*zp)[i+1];
159  (*gp)[i] = alpha_ * diff/std::sqrt(std::pow(diff,2.0)+eps_);
160  diff = (*zp)[i-1]-(*zp)[i];
161  (*gp)[i] -= alpha_ * diff/std::sqrt(std::pow(diff,2.0)+eps_);
162  }
163  }
164  }
165  }
166 
167  void reg_hessVec(Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &z) {
168 
169  using Teuchos::RCP;
170 
171  if ( reg_type_ == 2 ) {
172  hv.set(v);
173  hv.scale(alpha_*hz_);
174  }
175  else if ( reg_type_ == 1 ) {
176  RCP<const vector> zp = getVector(z);
177  RCP<const vector> vp = getVector(v);
178  RCP<vector> hvp = getVector(hv);
179 
180  for (uint i = 0; i < nz_; i++) {
181  (*hvp)[i] = alpha_*hz_*(*vp)[i]*eps_/std::pow(std::pow((*zp)[i],2.0)+eps_,3.0/2.0);
182  }
183  }
184  else if ( reg_type_ == 0 ) {
185  RCP<const vector> zp = getVector(z);
186  RCP<const vector> vp = getVector(v);
187  RCP<vector> hvp = getVector(hv);
188 
189  Real diff1 = 0.0;
190  Real diff2 = 0.0;
191  for (uint i = 0; i < nz_; i++) {
192  if ( i == 0 ) {
193  diff1 = (*zp)[i]-(*zp)[i+1];
194  diff1 = eps_/std::pow(std::pow(diff1,2.0)+eps_,3.0/2.0);
195  (*hvp)[i] = alpha_* ((*vp)[i]*diff1 - (*vp)[i+1]*diff1);
196  }
197  else if ( i == nz_-1 ) {
198  diff2 = (*zp)[i-1]-(*zp)[i];
199  diff2 = eps_/std::pow(std::pow(diff2,2.0)+eps_,3.0/2.0);
200  (*hvp)[i] = alpha_* (-(*vp)[i-1]*diff2 + (*vp)[i]*diff2);
201  }
202  else {
203  diff1 = (*zp)[i]-(*zp)[i+1];
204  diff1 = eps_/std::pow(std::pow(diff1,2.0)+eps_,3.0/2.0);
205  diff2 = (*zp)[i-1]-(*zp)[i];
206  diff2 = eps_/std::pow(std::pow(diff2,2.0)+eps_,3.0/2.0);
207  (*hvp)[i] = alpha_* (-(*vp)[i-1]*diff2 + (*vp)[i]*(diff1 + diff2) - (*vp)[i+1]*diff1);
208  }
209  }
210  }
211  }
212 
213  /* FINITE ELEMENT DEFINTIONS */
214  void apply_mass(Vector<Real> &Mf, const Vector<Real> &f ) {
215 
216  using Teuchos::RCP;
217  RCP<const vector> fp = getVector(f);
218  RCP<vector> Mfp = getVector(Mf);
219 
220  for (uint i = 0; i < nu_; i++) {
221  if ( i == 0 ) {
222  (*Mfp)[i] = hu_/6.0*(4.0*(*fp)[i] + (*fp)[i+1]);
223  }
224  else if ( i == nu_-1 ) {
225  (*Mfp)[i] = hu_/6.0*((*fp)[i-1] + 4.0*(*fp)[i]);
226  }
227  else {
228  (*Mfp)[i] = hu_/6.0*((*fp)[i-1] + 4.0*(*fp)[i] + (*fp)[i+1]);
229  }
230  }
231  }
232 
234 
235  using Teuchos::RCP;
236  RCP<const vector> zp = getVector(z);
237  RCP<vector> up = getVector(u);
238  RCP<vector> bp = getVector(b);
239 
240  // Get Diagonal and Off-Diagonal Entries of PDE Jacobian
241  vector d(nu_,1.0);
242  vector o(nu_-1,1.0);
243  for ( uint i = 0; i < nu_; i++ ) {
244  d[i] = (std::exp((*zp)[i]) + std::exp((*zp)[i+1]))/hu_;
245  if ( i < nu_-1 ) {
246  o[i] *= -std::exp((*zp)[i+1])/hu_;
247  }
248  }
249 
250  // Solve Tridiagonal System Using LAPACK's SPD Tridiagonal Solver
251  Teuchos::LAPACK<int,Real> lp;
252  int info;
253  int ldb = nu_;
254  int nhrs = 1;
255  lp.PTTRF(nu_,&d[0],&o[0],&info);
256  lp.PTTRS(nu_,nhrs,&d[0],&o[0],&(*bp)[0],ldb,&info);
257  u.set(b);
258  }
259 
260  Real evaluate_target(Real x) {
261  return x*(1.0-x);
262  }
263 
265  const Vector<Real> &d, const Vector<Real> &u ) {
266 
267  using Teuchos::RCP;
268 
269  RCP<const vector> zp = getVector(z);
270  RCP<const vector> up = getVector(u);
271  RCP<const vector> dp = getVector(d);
272  RCP<vector> Bdp = getVector(Bd);
273 
274  for (uint i = 0; i < nu_; i++) {
275  if ( i == 0 ) {
276  (*Bdp)[i] = 1.0/hu_*( std::exp((*zp)[i])*(*up)[i]*(*dp)[i]
277  + std::exp((*zp)[i+1])*((*up)[i]-(*up)[i+1])*(*dp)[i+1] );
278  }
279  else if ( i == nu_-1 ) {
280  (*Bdp)[i] = 1.0/hu_*( std::exp((*zp)[i])*((*up)[i]-(*up)[i-1])*(*dp)[i]
281  + std::exp((*zp)[i+1])*(*up)[i]*(*dp)[i+1] );
282  }
283  else {
284  (*Bdp)[i] = 1.0/hu_*( std::exp((*zp)[i])*((*up)[i]-(*up)[i-1])*(*dp)[i]
285  + std::exp((*zp)[i+1])*((*up)[i]-(*up)[i+1])*(*dp)[i+1] );
286  }
287  }
288  }
289 
291  const Vector<Real> &d, const Vector<Real> &u ) {
292  using Teuchos::RCP;
293 
294  RCP<const vector> zp = getVector(z);
295  RCP<const vector> up = getVector(u);
296  RCP<const vector> dp = getVector(d);
297  RCP<vector> Bdp = getVector(Bd);
298 
299  for (uint i = 0; i < nz_; i++) {
300  if ( i == 0 ) {
301  (*Bdp)[i] = std::exp((*zp)[i])/hu_*(*up)[i]*(*dp)[i];
302  }
303  else if ( i == nz_-1 ) {
304  (*Bdp)[i] = std::exp((*zp)[i])/hu_*(*up)[i-1]*(*dp)[i-1];
305  }
306  else {
307  (*Bdp)[i] = std::exp((*zp)[i])/hu_*( ((*up)[i]-(*up)[i-1])*((*dp)[i]-(*dp)[i-1]) );
308  }
309  }
310  }
311 
313  const Vector<Real> &d, const Vector<Real> &u ) {
314  using Teuchos::RCP;
315  RCP<const vector> zp = getVector(z);
316  RCP<const vector> vp = getVector(v);
317  RCP<const vector> up = getVector(u);
318  RCP<const vector> dp = getVector(d);
319  RCP<vector> Bdp = getVector(Bd);
320 
321  for (uint i = 0; i < nz_; i++) {
322  if ( i == 0 ) {
323  (*Bdp)[i] = (*vp)[i]*std::exp((*zp)[i])/hu_*(*up)[i]*(*dp)[i];
324  }
325  else if ( i == nz_-1 ) {
326  (*Bdp)[i] = (*vp)[i]*std::exp((*zp)[i])/hu_*(*up)[i-1]*(*dp)[i-1];
327  }
328  else {
329  (*Bdp)[i] = (*vp)[i]*std::exp((*zp)[i])/hu_*( ((*up)[i]-(*up)[i-1])*((*dp)[i]-(*dp)[i-1]) );
330  }
331  }
332  }
333 
334  /* STATE AND ADJOINT EQUATION DEFINTIONS */
336 
337  using Teuchos::RCP;
338  using Teuchos::rcp;
339 
340  Real k1 = 1.0;
341  Real k2 = 2.0;
342  // Right Hand Side
343  RCP<vector> bp = rcp( new vector(nu_, 0.0) );
344  for ( uint i = 0; i < nu_; i++ ) {
345  if ( (Real)(i+1)*hu_ < 0.5 ) {
346  (*bp)[i] = 2.0*k1*hu_;
347  }
348  else if ( std::abs((Real)(i+1)*hu_ - 0.5) < ROL_EPSILON<Real>() ) {
349  (*bp)[i] = (k1+k2)*hu_;
350  }
351  else if ( (Real)(i+1)*hu_ > 0.5 ) {
352  (*bp)[i] = 2.0*k2*hu_;
353  }
354  }
355 
356  SV b(bp);
357  // Solve Equation
358  solve_poisson(u,z,b);
359  }
360 
362 
363  using Teuchos::RCP;
364  using Teuchos::rcp;
365 
366  RCP<const vector> up = getVector(u);
367  RCP<vector> rp = rcp( new vector(nu_,0.0) );
368  SV res(rp);
369 
370  for ( uint i = 0; i < nu_; i++) {
371  (*rp)[i] = -((*up)[i]-evaluate_target((Real)(i+1)*hu_));
372  }
373  StdVector<Real> Mres( Teuchos::rcp( new std::vector<Real>(nu_,0.0) ) );
374  apply_mass(Mres,res);
375  solve_poisson(p,z,Mres);
376  }
377 
379  const Vector<Real> &u, const Vector<Real> &z) {
380 
381  using Teuchos::rcp;
382  SV b( rcp( new vector(nu_,0.0) ) );
384  solve_poisson(w,z,b);
385  }
386 
388  const Vector<Real> &p, const Vector<Real> &u, const Vector<Real> &z) {
389 
390  using Teuchos::rcp;
391 
392  SV res( rcp( new vector(nu_,0.0) ) );
393  apply_mass(res,w);
394  SV res1( rcp( new vector(nu_,0.0) ) );
396  res.axpy(-1.0,res1);
397  solve_poisson(q,z,res);
398  }
399 
400  /* OBJECTIVE FUNCTION DEFINITIONS */
401  Real value( const Vector<Real> &z, Real &tol ) {
402 
403  using Teuchos::RCP;
404  using Teuchos::rcp;
405 
406  // SOLVE STATE EQUATION
407  RCP<vector> up = rcp( new vector(nu_,0.0) );
408  SV u( up );
409 
411 
412  // COMPUTE MISFIT
413  RCP<vector> rp = rcp( new vector(nu_,0.0) );
414  SV res( rp );
415 
416  for ( uint i = 0; i < nu_; i++) {
417  (*rp)[i] = ((*up)[i]-evaluate_target((Real)(i+1)*hu_));
418  }
419 
420  RCP<V> Mres = res.clone();
421  apply_mass(*Mres,res);
422  return 0.5*Mres->dot(res) + reg_value(z);
423  }
424 
425  void gradient( Vector<Real> &g, const Vector<Real> &z, Real &tol ) {
426 
427  using Teuchos::rcp;
428 
429  // SOLVE STATE EQUATION
430  SV u( rcp( new vector(nu_,0.0) ) );
432 
433  // SOLVE ADJOINT EQUATION
434  SV p( Teuchos::rcp( new std::vector<Real>(nu_,0.0) ) );
435  solve_adjoint_equation(p,u,z);
436 
437  // Apply Transpose of Linearized Control Operator
439 
440  // Regularization gradient
441  SV g_reg( rcp( new vector(nz_,0.0) ) );
442  reg_gradient(g_reg,z);
443 
444  // Build Gradient
445  g.plus(g_reg);
446  }
447 #if USE_HESSVEC
448  void hessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &z, Real &tol ) {
449 
450  using Teuchos::rcp;
451 
452  // SOLVE STATE EQUATION
453  SV u( rcp( new vector(nu_,0.0) ) );
455 
456  // SOLVE ADJOINT EQUATION
457  SV p( rcp( new vector(nu_,0.0) ) );
458  solve_adjoint_equation(p,u,z);
459 
460  // SOLVE STATE SENSITIVITY EQUATION
461  SV w( rcp( new vector(nu_,0.0) ) );
463 
464  // SOLVE ADJOINT SENSITIVITY EQUATION
465  SV q( rcp( new vector(nu_,0.0) ) );
467 
468  // Apply Transpose of Linearized Control Operator
470 
471  // Apply Transpose of Linearized Control Operator
472  SV tmp( rcp( new vector(nz_,0.0) ) );
474  hv.axpy(-1.0,tmp);
475 
476  // Apply Transpose of 2nd Derivative of Control Operator
477  tmp.zero();
479  hv.plus(tmp);
480 
481  // Regularization hessVec
482  SV hv_reg( rcp( new vector(nz_,0.0) ) );
483  reg_hessVec(hv_reg,v,z);
484 
485  // Build hessVec
486  hv.plus(hv_reg);
487  }
488 #endif
489 
490  void invHessVec( Vector<Real> &hv, const Vector<Real> &v, const Vector<Real> &x, Real &tol ) {
491 
492  using Teuchos::RCP;
493  using Teuchos::rcp;
494 
495  // Cast hv and v vectors to std::vector.
496  RCP<vector> hvp = getVector(hv);
497 
498  std::vector<Real> vp(*getVector(v));
499 
500  int dim = static_cast<int>(vp.size());
501 
502 
503  // Compute dense Hessian.
504  Teuchos::SerialDenseMatrix<int, Real> H(dim, dim);
505  Objective_PoissonInversion<Real> & obj = *this;
506  H = computeDenseHessian<Real>(obj, x);
507 
508  // Compute eigenvalues, sort real part.
509  std::vector<vector> eigenvals = computeEigenvalues<Real>(H);
510  std::sort((eigenvals[0]).begin(), (eigenvals[0]).end());
511 
512  // Perform 'inertia' correction.
513  Real inertia = (eigenvals[0])[0];
514  Real correction = 0.0;
515  if ( inertia <= 0.0 ) {
516  correction = 2.0*std::abs(inertia);
517  if ( inertia == 0.0 ) {
518  int cnt = 0;
519  while ( eigenvals[0][cnt] == 0.0 ) {
520  cnt++;
521  }
522  correction = 0.5*eigenvals[0][cnt];
523  if ( cnt == dim-1 ) {
524  correction = 1.0;
525  }
526  }
527  for (int i=0; i<dim; i++) {
528  H(i,i) += correction;
529  }
530  }
531 
532  // Compute dense inverse Hessian.
533  Teuchos::SerialDenseMatrix<int, Real> invH = computeInverse<Real>(H);
534 
535  // Apply dense inverse Hessian.
536  Teuchos::SerialDenseVector<int, Real> hv_teuchos(Teuchos::View, &((*hvp)[0]), dim);
537  const Teuchos::SerialDenseVector<int, Real> v_teuchos(Teuchos::View, &(vp[0]), dim);
538  hv_teuchos.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1.0, invH, v_teuchos, 0.0);
539  }
540 
541 };
542 
543 template<class Real>
544 void getPoissonInversion( Teuchos::RCP<Objective<Real> > &obj,
545  Teuchos::RCP<Vector<Real> > &x0,
546  Teuchos::RCP<Vector<Real> > &x ) {
547  // Problem dimension
548  int n = 128;
549 
550  // Get Initial Guess
551  Teuchos::RCP<std::vector<Real> > x0p = Teuchos::rcp(new std::vector<Real>(n,0.0));
552  for ( int i = 0; i < n; i++ ) {
553  (*x0p)[i] = 1.5;
554  }
555  x0 = Teuchos::rcp(new StdVector<Real>(x0p));
556 
557  // Get Solution
558  Teuchos::RCP<std::vector<Real> > xp = Teuchos::rcp(new std::vector<Real>(n,0.0));
559  Real h = 1.0/((Real)n+1), pt = 0.0, k1 = 1.0, k2 = 2.0;
560  for( int i = 0; i < n; i++ ) {
561  pt = (Real)(i+1)*h;
562  if ( pt >= 0.0 && pt < 0.5 ) {
563  (*xp)[i] = std::log(k1);
564  }
565  else if ( pt >= 0.5 && pt < 1.0 ) {
566  (*xp)[i] = std::log(k2);
567  }
568  }
569  x = Teuchos::rcp(new StdVector<Real>(xp));
570 
571  // Instantiate Objective Function
572  obj = Teuchos::rcp(new Objective_PoissonInversion<Real>(n,1.e-6));
573 }
574 
575 } // End ZOO Namespace
576 } // End ROL Namespace
577 
578 #endif
Provides the interface to evaluate objective functions.
void apply_transposed_linearized_control_operator(Vector< Real > &Bd, const Vector< Real > &z, const Vector< Real > &d, const Vector< Real > &u)
virtual void scale(const Real alpha)=0
Compute where .
void axpy(const Real alpha, const Vector< Real > &x)
Compute where .
void solve_state_sensitivity_equation(Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &u, const Vector< Real > &z)
virtual void plus(const Vector &x)=0
Compute , where .
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:143
void gradient(Vector< Real > &g, const Vector< Real > &z, Real &tol)
Compute gradient.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
Contains definitions for helper functions in ROL.
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:74
void solve_poisson(Vector< Real > &u, const Vector< Real > &z, Vector< Real > &b)
void solve_adjoint_sensitivity_equation(Vector< Real > &q, const Vector< Real > &w, const Vector< Real > &v, const Vector< Real > &p, const Vector< Real > &u, const Vector< Real > &z)
void apply_transposed_linearized_control_operator_2(Vector< Real > &Bd, const Vector< Real > &z, const Vector< Real > &v, const Vector< Real > &d, const Vector< Real > &u)
void invHessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply inverse Hessian approximation to vector.
void solve_state_equation(Vector< Real > &u, const Vector< Real > &z)
void reg_hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &z)
void getPoissonInversion(Teuchos::RCP< Objective< Real > > &obj, Teuchos::RCP< Vector< Real > > &x0, Teuchos::RCP< Vector< Real > > &x)
void apply_mass(Vector< Real > &Mf, const Vector< Real > &f)
void solve_adjoint_equation(Vector< Real > &p, const Vector< Real > &u, const Vector< Real > &z)
Teuchos::RCP< const vector > getVector(const V &x)
void apply_linearized_control_operator(Vector< Real > &Bd, const Vector< Real > &z, const Vector< Real > &d, const Vector< Real > &u)
virtual Teuchos::RCP< Vector< Real > > clone() const
Clone to make a new (uninitialized) vector.
Real value(const Vector< Real > &z, Real &tol)
Compute value.
Objective_PoissonInversion(uint nz=32, Real alpha=1.e-4)
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:196
Real reg_value(const Vector< Real > &z)
void reg_gradient(Vector< Real > &g, const Vector< Real > &z)
Teuchos::RCP< vector > getVector(V &x)