ROL
ROL_LinMore.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 
44 #ifndef ROL_LINMORE_H
45 #define ROL_LINMORE_H
46 
51 #include "ROL_TrustRegion.hpp"
52 #include "ROL_LinMoreModel.hpp"
53 #include "ROL_Elementwise_Function.hpp"
54 #include "ROL_Elementwise_Reduce.hpp"
55 
56 namespace ROL {
57 
58 template<class Real>
59 class LinMore : public TrustRegion<Real> {
60 private:
61 
62  Ptr<Vector<Real>> x_, s_, g_;
63  Ptr<Vector<Real>> pwa1_, pwa2_, dwa1_, dwa2_;
64 
65  Real tol1_, tol2_, alpha_;
66  int maxit_;
67 
68  unsigned verbosity_;
69 
70  class LowerBreakPoint : public Elementwise::BinaryFunction<Real> {
71  public:
72  Real apply( const Real &x, const Real &y) const {
73  const Real zero(0), one(1);
74  return (x > zero && y < zero) ? -x/y : -one;
75  }
76  } lbp_;
77 
78  class UpperBreakPoint : public Elementwise::BinaryFunction<Real> {
79  public:
80  Real apply( const Real &x, const Real &y) const {
81  const Real zero(0), one(1);
82  return (x > zero && y > zero) ? x/y : -one;
83  }
84  } ubp_;
85 
86  class PositiveMin : public Elementwise::ReductionOp<Real> {
87  public:
88  void reduce(const Real &input, Real &output) const {
89  const Real zero(0);
90  output = (input<output && input>zero) ? input : output;
91  }
92  void reduce( const volatile Real &input, Real volatile &output ) const {
93  const Real zero(0);
94  output = (input<output && input>zero) ? input : output;
95  }
96  Real initialValue() const {
97  return ROL_INF<Real>();
98  }
99  Elementwise::EReductionType reductionType() const {
100  return Elementwise::REDUCE_MIN;
101  }
102  } pmin_;
103 
104  class PositiveMax : public Elementwise::ReductionOp<Real> {
105  public:
106  void reduce(const Real &input, Real &output) const {
107  const Real zero(0);
108  output = (input>output && input>zero) ? input : output;
109  }
110  void reduce( const volatile Real &input, Real volatile &output ) const {
111  const Real zero(0);
112  output = (input>output && input>zero) ? input : output;
113  }
114  Real initialValue() const {
115  return static_cast<Real>(0);
116  }
117  Elementwise::EReductionType reductionType() const {
118  return Elementwise::REDUCE_MAX;
119  }
120  } pmax_;
121 
122 public:
123 
124  // Constructor
125  LinMore( ROL::ParameterList &parlist ) : TrustRegion<Real>(parlist), alpha_(1) {
126  // Unravel Parameter List
127  Real em4(1e-4), em2(1e-2);
128  maxit_ = parlist.sublist("General").sublist("Krylov").get("Iteration Limit",20);
129  tol1_ = parlist.sublist("General").sublist("Krylov").get("Absolute Tolerance",em4);
130  tol2_ = parlist.sublist("General").sublist("Krylov").get("Relative Tolerance",em2);
131  // Get verbosity level
132  verbosity_ = parlist.sublist("General").get("Print Verbosity", 0);
133  }
134 
135  void initialize( const Vector<Real> &x, const Vector<Real> &s, const Vector<Real> &g) {
137  x_ = x.clone(); s_ = x.clone(); g_ = g.clone();
138  pwa1_ = x.clone(); pwa2_ = x.clone();
139  dwa1_ = g.clone(); dwa2_ = g.clone();
140  }
141 
142  void run( Vector<Real> &s,
143  Real &snorm,
144  int &iflag,
145  int &iter,
146  const Real del,
147  TrustRegionModel<Real> &model ) {
148  const Real zero(0), half(0.5), one(1);
149  Real tol0 = std::sqrt(ROL_EPSILON<Real>());
150  Real gfnorm(0), gfnormf(0), tol(0), stol(0);
151  int dim = s.dimension();
152  // Compute Cauchy point (TRON notation: x_ = x[1])
153  snorm = dcauchy(*s_,alpha_,*model.getIterate(),model.getGradient()->dual(),
154  del,model,*pwa1_,*pwa2_,*dwa1_); // Solve 1D optimization problem for alpha_
155  x_->set(*model.getIterate()); // TRON notation: model.getIterate() = x[0]
156  x_->plus(*s_); // Set x_ = x[0] + alpha_*g
157  model.getBoundConstraint()->project(*x_); // Project x_ onto bounds
158 
159  // Model gradient at s = x[1] - x[0]
160  s.set(*x_); s.axpy(-one,*model.getIterate()); // s_ = x[i+1]-x[0]
161  dynamic_cast<LinMoreModel<Real>&>(model).applyFullHessian(*g_,s,tol0);
162  g_->plus(*model.getGradient());
163  model.getBoundConstraint()->pruneActive(*g_,*x_,zero);
164  gfnorm = g_->norm();
165  if (verbosity_ > 0) {
166  std::cout << std::endl;
167  std::cout << " Computation of Cauchy point" << std::endl;
168  std::cout << " Step length (alpha): " << alpha_ << std::endl;
169  std::cout << " Step length (alpha*g): " << snorm << std::endl;
170  std::cout << " Norm of Cauchy point: " << x_->norm() << std::endl;
171  std::cout << " Norm of free gradient components: " << gfnorm << std::endl;
172  }
173 
174  // Main step computation loop
175  // There are at most dim iterations since at least one
176  // face becomes active at each iteration
177  iter = 0;
178  for (int i = 0; i < dim; ++i) {
179  // Run Truncated CG
180  int flagCG = 0, iterCG = 0;
181  tol = tol2_*gfnorm;
182  stol = zero;
183  snorm = dtrpcg(*s_,flagCG,iterCG,*g_,*x_,del,model,
184  tol,stol,maxit_,
185  *pwa1_,*dwa1_,*pwa2_,*dwa2_);
186  iter += iterCG;
187  if (verbosity_ > 0) {
188  std::cout << std::endl;
189  std::cout << " Computation of CG step" << std::endl;
190  std::cout << " Number of faces: " << dim << std::endl;
191  std::cout << " Current face (i): " << i << std::endl;
192  std::cout << " CG step length: " << snorm << std::endl;
193  std::cout << " Number of CG iterations: " << iterCG << std::endl;
194  std::cout << " CG flag: " << flagCG << std::endl;
195  std::cout << " Total number of iterations: " << iter << std::endl;
196  }
197 
198  // Projected search
199  snorm = dprsrch(*x_,*s_,g_->dual(),model,*pwa1_,*dwa1_);
200  if (verbosity_ > 0) {
201  std::cout << " Step length (beta*s): " << snorm << std::endl;
202  std::cout << " Iterate length: " << x_->norm() << std::endl;
203  }
204 
205  // Model gradient at s = x[i+1] - x[0]
206  s.set(*x_); s.axpy(-one,*model.getIterate()); // s_ = x[i+1]-x[0]
207  dynamic_cast<LinMoreModel<Real>&>(model).applyFullHessian(*g_,s,tol0);
208  g_->plus(*model.getGradient());
209  model.getBoundConstraint()->pruneActive(*g_,*x_,zero);
210  gfnormf = g_->norm();
211  if (verbosity_ > 0) {
212  std::cout << std::endl;
213  std::cout << " Update model gradient" << std::endl;
214  std::cout << " Step length: " << s.norm() << std::endl;
215  std::cout << " Norm of free gradient components: " << gfnormf << std::endl;
216  std::cout << std::endl;
217  }
218 
219  // Termination check
220  if (gfnormf <= tol2_*gfnorm) {
221  iflag = 0;
222  break;
223  }
224  else if (iter >= maxit_) {
225  iflag = 1;
226  break;
227  }
228  else if (flagCG == 2) {
229  iflag = 2;
230  break;
231  }
232  else if (flagCG == 3) {
233  iflag = 3;
234  break;
235  }
236 
237  // Update free gradient norm
238  gfnorm = gfnormf;
239  }
240  // Update norm of step and update model predicted reduction
241  snorm = s.norm();
242  Real gs(0), q(0);
243  dynamic_cast<LinMoreModel<Real>&>(model).applyFullHessian(*dwa1_,s,tol0);
244  gs = s.dot(model.getGradient()->dual());
245  q = half * s.dot(dwa1_->dual()) + gs;
247  }
248 
249 private:
250 
251  // Compute the projected step s = P(x + alpha*w) - x
252  // Returns the norm of the projected step s
253  // s -- The projected step upon return
254  // w -- The direction vector w (unchanged)
255  // x -- The anchor vector x (unchanged)
256  // alpha -- The step size (unchanged)
257  // model -- Contains the bound constraint information
258  Real dgpstep(Vector<Real> &s, const Vector<Real> &w,
259  const Vector<Real> &x, const Real alpha,
260  TrustRegionModel<Real> &model) const {
261  s.set(x); s.axpy(alpha,w);
262  model.getBoundConstraint()->project(s);
263  s.axpy(static_cast<Real>(-1),x);
264  return s.norm();
265  }
266 
267  // Compute minimal and maximal break points of x+alpha*s
268  // with in the interval [xl,xu] specified by the bound constraint
269  // x -- The anchor vector x (unchanged)
270  // s -- The descent vector s (unchanged)
271  // model -- Contains the bound constraint information
272  // bpmin -- The minimum break point
273  // bpmax -- The maximum break point
274  // pwa -- A primal working vector
275  void dbreakpt(const Vector<Real> &x, const Vector<Real> &s,
276  TrustRegionModel<Real> &model,
277  Real &bpmin, Real &bpmax,
278  Vector<Real> &pwa) {
279  const Real zero(0), one(1);
280  bpmin = one; bpmax = zero;
281  Real lbpmin = one, lbpmax = zero, ubpmin = one, ubpmax = zero;
282  // Compute lower break points
283  if (model.getBoundConstraint()->isLowerActivated()) {
284  pwa.set(x);
285  pwa.axpy(-one,*model.getBoundConstraint()->getLowerBound());
286  pwa.applyBinary(lbp_,s);
287  if (pwa.norm() != zero) {
288  lbpmin = pwa.reduce(pmin_);
289  lbpmax = pwa.reduce(pmax_);
290  }
291  }
292  // Compute upper break points
293  if (model.getBoundConstraint()->isUpperActivated()) {
294  pwa.set(*model.getBoundConstraint()->getUpperBound());
295  pwa.axpy(-one,x);
296  pwa.applyBinary(ubp_,s);
297  if (pwa.norm() != zero) {
298  ubpmin = pwa.reduce(pmin_);
299  ubpmax = pwa.reduce(pmax_);
300  }
301  }
302  bpmin = std::min(lbpmin,ubpmin);
303  bpmax = std::max(lbpmax,ubpmax);
304  if (bpmin > bpmax) {
305  bpmin = zero;
306  bpmax = zero;
307  }
308  if (verbosity_ > 0) {
309  std::cout << std::endl;
310  std::cout << " Computation of break points" << std::endl;
311  std::cout << " Minimum break point: " << bpmin << std::endl;
312  std::cout << " Maximum break point: " << bpmax << std::endl;
313  }
314  }
315 
316  // Compute Cauchy point, i.e., the minimizer of q(P(x - alpha*g)-x)
317  // subject to the trust region constraint ||P(x - alpha*g)-x|| <= del
318  // s -- The Cauchy point upon return
319  // alpha -- The step length for the Cauchy point upon return
320  // x -- The anchor vector x (unchanged)
321  // g -- The (dual) gradient vector g (unchanged)
322  // del -- The trust region radius (unchanged)
323  // model -- Contains the objective and bound constraint information
324  // pwa1 -- Primal working array
325  // pwa2 -- Primal working array
326  // dwa -- Dual working array
327  Real dcauchy(Vector<Real> &s, Real &alpha,
328  const Vector<Real> &x, const Vector<Real> &g,
329  const Real del, TrustRegionModel<Real> &model,
330  Vector<Real> &pwa1, Vector<Real> &pwa2, Vector<Real> &dwa) {
331  const Real half(0.5), one(1), mu0(0.01), interpf(0.1), extrapf(10);
332  // const Real zero(0); // Unused
333  Real tol = std::sqrt(ROL_EPSILON<Real>());
334  bool interp = false;
335  Real q(0), gs(0), bpmin(0), bpmax(0), snorm(0);
336  // Compute minimal and maximal break points of x[0] - alpha g[0]
337  pwa1.set(g); pwa1.scale(-one);
338  dbreakpt(x,pwa1,model,bpmin,bpmax,pwa2);
339  // Compute s = P(x[0] - alpha g[0].dual())
340  snorm = dgpstep(s,g,x,-alpha,model);
341  if (snorm > del) {
342  interp = true;
343  }
344  else {
345  dynamic_cast<LinMoreModel<Real>&>(model).applyFullHessian(dwa,s,tol);
346  gs = s.dot(g);
347  q = half * s.dot(dwa.dual()) + gs;
348  interp = (q > mu0*gs);
349  }
350  // Either increase or decrease alpha to find approximate Cauchy point
351  if (interp) {
352  bool search = true;
353  while (search) {
354  alpha *= interpf;
355  snorm = dgpstep(s,g,x,-alpha,model);
356  if (snorm <= del) {
357  dynamic_cast<LinMoreModel<Real>&>(model).applyFullHessian(dwa,s,tol);
358  gs = s.dot(g);
359  q = half * s.dot(dwa.dual()) + gs;
360  search = (q > mu0*gs);
361  }
362  }
363  }
364  else {
365  bool search = true;
366  Real alphas = alpha;
367  while (search && alpha <= bpmax) {
368  alpha *= extrapf;
369  snorm = dgpstep(s,g,x,-alpha,model);
370  if (snorm <= del) {
371  dynamic_cast<LinMoreModel<Real>&>(model).applyFullHessian(dwa,s,tol);
372  gs = s.dot(g);
373  q = half * s.dot(dwa.dual()) + gs;
374  if (q < mu0*gs) {
375  search = true;
376  alphas = alpha;
377  }
378  }
379  else {
380  search = false;
381  }
382  }
383  alpha = alphas;
384  snorm = dgpstep(s,g,x,-alpha,model);
385  }
386  return snorm;
387  }
388 
389  // Perform projected search to determine beta such that
390  // q(P(x + beta*s)-x) <= mu0*g'(P(x + beta*s)-x) for mu0 in (0,1)
391  // x -- The anchor vector x, upon return x = P(x + beta*s)
392  // s -- The direction vector s, upon return s = P(x + beta*s) - x
393  // g -- The free components of the gradient vector g (unchanged)
394  // model -- Contains objective and bound constraint information
395  // pwa -- Primal working array
396  // dwa -- Dual working array
398  const Vector<Real> &g, TrustRegionModel<Real> &model,
399  Vector<Real> &pwa, Vector<Real> &dwa) {
400  const Real half(0.5), one(1), mu0(0.01), interpf(0.5);
401  Real tol = std::sqrt(ROL_EPSILON<Real>());
402  Real beta(1), snorm(0), q(0), gs(0), bpmin(0), bpmax(0);
403  int nsteps = 0;
404  // Compute break points of x+beta*s;
405  dbreakpt(x,s,model,bpmin,bpmax,pwa);
406  // Reduce beta until sufficient decrease is satisfied
407  bool search = true;
408  while (search && beta > bpmin) {
409  nsteps++;
410  snorm = dgpstep(pwa,s,x,beta,model);
411  dynamic_cast<LinMoreModel<Real>&>(model).applyFreeHessian(dwa,pwa,x,tol);
412  gs = pwa.dot(g);
413  q = half * s.dot(dwa.dual()) + gs;
414  if (q <= mu0*gs) {
415  search = false;
416  }
417  else {
418  beta *= interpf;
419  }
420  }
421  if (beta < one && beta < bpmin) {
422  beta = bpmin;
423  }
424  snorm = dgpstep(pwa,s,x,beta,model);
425  s.set(pwa);
426  x.plus(s);
427  if (verbosity_ > 0) {
428  std::cout << std::endl;
429  std::cout << " Projected search" << std::endl;
430  std::cout << " Step length (beta): " << beta << std::endl;
431  }
432  return snorm;
433  }
434 
435  // Compute sigma such that ||x+sigma*p||_inv(M) = del. This is called
436  // if dtrpcg detects negative curvature or if the step violates
437  // the trust region bound
438  // xtx -- The dot product <x, inv(M)x> (unchanged)
439  // ptp -- The dot product <p, inv(M)p> (unchanged)
440  // ptx -- The dot product <p, inv(M)x> (unchanged)
441  // del -- Trust region radius (unchanged)
442  Real dtrqsol(const Real xtx, const Real ptp, const Real ptx, const Real del) {
443  const Real zero(0);
444  Real dsq = del*del;
445  Real rad = ptx*ptx + ptp*(dsq-xtx);
446  rad = std::sqrt(std::max(rad,zero));
447  Real sigma(0);
448  if (ptx > zero) {
449  sigma = (dsq-xtx)/(ptx+rad);
450  }
451  else if (rad > zero) {
452  sigma = (rad-ptx)/ptp;
453  }
454  else {
455  sigma = zero;
456  }
457  return sigma;
458  }
459 
460  // Solve the trust region subproblem: minimize q(w) subject to the
461  // trust region constraint ||w||_inv(M) <= del using the Steihaug-Toint
462  // Conjugate Gradients algorithm
463  // w -- The step vector to be computed
464  // iflag -- Termination flag
465  // iter -- Number of CG iterations
466  // del -- Trust region radius (unchanged)
467  // model -- Contains the objective and bound constraint information
468  // tol -- Residual stopping tolerance (unchanged)
469  // stol -- Preconditioned residual stopping tolerance (unchanged)
470  // itermax -- Maximum number of iterations
471  // p -- Primal working array that stores the CG step
472  // q -- Dual working array that stores the Hessian applied to p
473  // r -- Primal working array that stores the preconditioned residual
474  // t -- Dual working array that stores the residual
475  Real dtrpcg(Vector<Real> &w, int &iflag, int &iter,
476  const Vector<Real> &g, const Vector<Real> &x,
477  const Real del, TrustRegionModel<Real> &model,
478  const Real tol, const Real stol, const int itermax,
480  Vector<Real> &t) {
481  Real tol0 = std::sqrt(ROL_EPSILON<Real>());
482  const Real zero(0), one(1), two(2);
483  Real rho(0), tnorm(0), rnorm(0), rnorm0(0), kappa(0), beta(0), sigma(0), alpha(0), rtr(0);
484  Real sMs(0), pMp(0), sMp(0);
485  iter = 0; iflag = 0;
486  // Initialize step
487  w.zero();
488  // Compute residual
489  t.set(g); t.scale(-one);
490  // Preconditioned residual
491  dynamic_cast<LinMoreModel<Real>&>(model).applyFreePrecond(r,t,x,tol0);
492  rho = r.dot(t.dual());
493  rnorm0 = std::sqrt(rho);
494  if ( rnorm0 == zero ) {
495  return zero;
496  }
497  // Initialize direction
498  p.set(r);
499  pMp = rho;
500  // Iterate CG
501  for (iter = 0; iter < itermax; ++iter) {
502  // Apply Hessian to direction dir
503  dynamic_cast<LinMoreModel<Real>&>(model).applyFreeHessian(q,p,x,tol0);
504  // Compute sigma such that ||s+sigma*dir|| = del
505  kappa = p.dot(q.dual());
506  alpha = (kappa>zero) ? rho/kappa : zero;
507  sigma = dtrqsol(sMs,pMp,sMp,del);
508  // Check for negative curvature or if iterate exceeds trust region
509  if (kappa <= zero || alpha >= sigma) {
510  w.axpy(sigma,p);
511  iflag = (kappa<=zero) ? 2 : 3;
512  break;
513  }
514  // Update iterate and residuals
515  w.axpy(alpha,p);
516  t.axpy(-alpha,q);
517  dynamic_cast<LinMoreModel<Real>&>(model).applyFreePrecond(r,t,x,tol0);
518  // Exit if residual tolerance is met
519  rtr = r.dot(t.dual());
520  rnorm = std::sqrt(rtr);
521  tnorm = t.norm();
522  if (rnorm <= stol || tnorm <= tol) {
523  iflag = 0;
524  break;
525  }
526  // Compute p = r + beta * p
527  beta = rtr/rho;
528  p.scale(beta); p.plus(r);
529  rho = rtr;
530  // Update dot products
531  // sMs = <s, inv(M)s>
532  // sMp = <s, inv(M)p>
533  // pMp = <p, inv(M)p>
534  sMs = sMs + two*alpha*sMp + alpha*alpha*pMp;
535  sMp = beta*(sMp + alpha*pMp);
536  pMp = rho + beta*beta*pMp;
537  }
538  // Check iteration count
539  if (iter == itermax) {
540  iflag = 1;
541  }
542  if (iflag != 1) {
543  iter++;
544  }
545  return w.norm();
546  }
547 
548 };
549 
550 }
551 
552 #endif
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Definition: ROL_Vector.hpp:226
void dbreakpt(const Vector< Real > &x, const Vector< Real > &s, TrustRegionModel< Real > &model, Real &bpmin, Real &bpmax, Vector< Real > &pwa)
Provides interface for truncated CG trust-region subproblem solver.
Definition: ROL_LinMore.hpp:59
virtual void scale(const Real alpha)=0
Compute where .
Real apply(const Real &x, const Real &y) const
Definition: ROL_LinMore.hpp:80
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual int dimension() const
Return dimension of the vector space.
Definition: ROL_Vector.hpp:196
virtual void plus(const Vector &x)=0
Compute , where .
void run(Vector< Real > &s, Real &snorm, int &iflag, int &iter, const Real del, TrustRegionModel< Real > &model)
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Definition: ROL_Vector.hpp:153
Ptr< Vector< Real > > dwa2_
Definition: ROL_LinMore.hpp:63
virtual Real reduce(const Elementwise::ReductionOp< Real > &r) const
Definition: ROL_Vector.hpp:243
virtual void applyBinary(const Elementwise::BinaryFunction< Real > &f, const Vector &x)
Definition: ROL_Vector.hpp:236
virtual void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g)
LinMore(ROL::ParameterList &parlist)
Provides the interface to evaluate projected trust-region model functions from the Kelley-Sachs bound...
Ptr< Vector< Real > > g_
Definition: ROL_LinMore.hpp:62
void reduce(const Real &input, Real &output) const
Provides interface for and implements trust-region subproblem solvers.
Provides the interface to evaluate trust-region model functions.
virtual void zero()
Set to zero vector.
Definition: ROL_Vector.hpp:167
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:80
virtual const Ptr< const Vector< Real > > getGradient(void) const
Real dcauchy(Vector< Real > &s, Real &alpha, const Vector< Real > &x, const Vector< Real > &g, const Real del, TrustRegionModel< Real > &model, Vector< Real > &pwa1, Vector< Real > &pwa2, Vector< Real > &dwa)
ROL::LinMore::PositiveMax pmax_
ROL::LinMore::UpperBreakPoint ubp_
Ptr< Vector< Real > > pwa2_
Definition: ROL_LinMore.hpp:63
virtual Real dot(const Vector &x) const =0
Compute where .
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Real dgpstep(Vector< Real > &s, const Vector< Real > &w, const Vector< Real > &x, const Real alpha, TrustRegionModel< Real > &model) const
void reduce(const Real &input, Real &output) const
Definition: ROL_LinMore.hpp:88
void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g)
void setPredictedReduction(const Real pRed)
void reduce(const volatile Real &input, Real volatile &output) const
Real dtrpcg(Vector< Real > &w, int &iflag, int &iter, const Vector< Real > &g, const Vector< Real > &x, const Real del, TrustRegionModel< Real > &model, const Real tol, const Real stol, const int itermax, Vector< Real > &p, Vector< Real > &q, Vector< Real > &r, Vector< Real > &t)
ROL::LinMore::LowerBreakPoint lbp_
Real dprsrch(Vector< Real > &x, Vector< Real > &s, const Vector< Real > &g, TrustRegionModel< Real > &model, Vector< Real > &pwa, Vector< Real > &dwa)
Ptr< Vector< Real > > s_
Definition: ROL_LinMore.hpp:62
Ptr< Vector< Real > > pwa1_
Definition: ROL_LinMore.hpp:63
Real initialValue() const
Definition: ROL_LinMore.hpp:96
Elementwise::EReductionType reductionType() const
Ptr< Vector< Real > > x_
Definition: ROL_LinMore.hpp:62
Ptr< Vector< Real > > dwa1_
Definition: ROL_LinMore.hpp:63
Elementwise::EReductionType reductionType() const
Definition: ROL_LinMore.hpp:99
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:209
virtual Real norm() const =0
Returns where .
unsigned verbosity_
Definition: ROL_LinMore.hpp:68
Real dtrqsol(const Real xtx, const Real ptp, const Real ptx, const Real del)
ROL::LinMore::PositiveMin pmin_
virtual const Ptr< BoundConstraint< Real > > getBoundConstraint(void) const
void reduce(const volatile Real &input, Real volatile &output) const
Definition: ROL_LinMore.hpp:92
virtual const Ptr< const Vector< Real > > getIterate(void) const
Real apply(const Real &x, const Real &y) const
Definition: ROL_LinMore.hpp:72