ROL
ROL_TypeB_TrustRegionSPGAlgorithm_Def.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_TYPEB_TRUSTREGIONSPGALGORITHM_DEF_HPP
45#define ROL_TYPEB_TRUSTREGIONSPGALGORITHM_DEF_HPP
46
47#include <deque>
48
49namespace ROL {
50namespace TypeB {
51
52template<typename Real>
54 const Ptr<Secant<Real>> &secant) {
55 // Set status test
56 status_->reset();
57 status_->add(makePtr<StatusTest<Real>>(list));
58
59 ParameterList &trlist = list.sublist("Step").sublist("Trust Region");
60 // Trust-Region Parameters
61 state_->searchSize = trlist.get("Initial Radius", -1.0);
62 delMax_ = trlist.get("Maximum Radius", ROL_INF<Real>());
63 eta0_ = trlist.get("Step Acceptance Threshold", 0.05);
64 eta1_ = trlist.get("Radius Shrinking Threshold", 0.05);
65 eta2_ = trlist.get("Radius Growing Threshold", 0.9);
66 gamma0_ = trlist.get("Radius Shrinking Rate (Negative rho)", 0.0625);
67 gamma1_ = trlist.get("Radius Shrinking Rate (Positive rho)", 0.25);
68 gamma2_ = trlist.get("Radius Growing Rate", 2.5);
69 TRsafe_ = trlist.get("Safeguard Size", 100.0);
70 eps_ = TRsafe_*ROL_EPSILON<Real>();
71 interpRad_ = trlist.get("Use Radius Interpolation", false);
72 verbosity_ = trlist.sublist("General").get("Output Level", 0);
73 // Nonmonotone Parameters
74 storageNM_ = trlist.get("Nonmonotone Storage Size", 0);
75 useNM_ = (storageNM_ <= 0 ? false : true);
76 // Algorithm-Specific Parameters
77 ROL::ParameterList &lmlist = trlist.sublist("SPG");
78 mu0_ = lmlist.get("Sufficient Decrease Parameter", 1e-2);
79 spexp_ = lmlist.get("Relative Tolerance Exponent", 1.0);
80 spexp_ = std::max(static_cast<Real>(1),std::min(spexp_,static_cast<Real>(2)));
81 redlim_ = lmlist.sublist("Cauchy Point").get("Maximum Number of Reduction Steps", 10);
82 explim_ = lmlist.sublist("Cauchy Point").get("Maximum Number of Expansion Steps", 10);
83 alpha_ = lmlist.sublist("Cauchy Point").get("Initial Step Size", 1.0);
84 normAlpha_ = lmlist.sublist("Cauchy Point").get("Normalize Initial Step Size", false);
85 interpf_ = lmlist.sublist("Cauchy Point").get("Reduction Rate", 0.1);
86 extrapf_ = lmlist.sublist("Cauchy Point").get("Expansion Rate", 10.0);
87 qtol_ = lmlist.sublist("Cauchy Point").get("Decrease Tolerance", 1e-8);
88 // Spectral projected gradient parameters
89 lambdaMin_ = lmlist.sublist("Solver").get("Minimum Spectral Step Size", 1e-8);
90 lambdaMax_ = lmlist.sublist("Solver").get("Maximum Spectral Step Size", 1e8);
91 gamma_ = lmlist.sublist("Solver").get("Sufficient Decrease Tolerance", 1e-4);
92 maxSize_ = lmlist.sublist("Solver").get("Maximum Storage Size", 10);
93 maxit_ = lmlist.sublist("Solver").get("Iteration Limit", 25);
94 tol1_ = lmlist.sublist("Solver").get("Absolute Tolerance", 1e-4);
95 tol2_ = lmlist.sublist("Solver").get("Relative Tolerance", 1e-2);
96 useMin_ = lmlist.sublist("Solver").get("Use Smallest Model Iterate", true);
97 useNMSP_ = lmlist.sublist("Solver").get("Use Nonmonotone Search", false);
98 useSimpleSPG_ = !lmlist.sublist("Solver").get("Compute Cauchy Point", true);
99 // Inexactness Information
100 ParameterList &glist = list.sublist("General");
101 useInexact_.clear();
102 useInexact_.push_back(glist.get("Inexact Objective Function", false));
103 useInexact_.push_back(glist.get("Inexact Gradient", false));
104 useInexact_.push_back(glist.get("Inexact Hessian-Times-A-Vector", false));
105 // Trust-Region Inexactness Parameters
106 ParameterList &ilist = trlist.sublist("Inexact").sublist("Gradient");
107 scale0_ = ilist.get("Tolerance Scaling", static_cast<Real>(0.1));
108 scale1_ = ilist.get("Relative Tolerance", static_cast<Real>(2));
109 // Inexact Function Evaluation Information
110 ParameterList &vlist = trlist.sublist("Inexact").sublist("Value");
111 scale_ = vlist.get("Tolerance Scaling", static_cast<Real>(1.e-1));
112 omega_ = vlist.get("Exponent", static_cast<Real>(0.9));
113 force_ = vlist.get("Forcing Sequence Initial Value", static_cast<Real>(1.0));
114 updateIter_ = vlist.get("Forcing Sequence Update Frequency", static_cast<int>(10));
115 forceFactor_ = vlist.get("Forcing Sequence Reduction Factor", static_cast<Real>(0.1));
116 // Output Parameters
117 verbosity_ = list.sublist("General").get("Output Level",0);
118 writeHeader_ = verbosity_ > 2;
119 // Secant Information
120 useSecantPrecond_ = list.sublist("General").sublist("Secant").get("Use as Preconditioner", false);
121 useSecantHessVec_ = list.sublist("General").sublist("Secant").get("Use as Hessian", false);
123 if (useSecantPrecond_ && !useSecantHessVec_) mode = SECANTMODE_INVERSE;
124 else if (useSecantHessVec_ && !useSecantPrecond_) mode = SECANTMODE_FORWARD;
125 // Initialize trust region model
126 model_ = makePtr<TrustRegionModel_U<Real>>(list,secant,mode);
127 if (secant == nullPtr) {
128 esec_ = StringToESecant(list.sublist("General").sublist("Secant").get("Type","Limited-Memory BFGS"));
129 }
130}
131
132template<typename Real>
134 const Vector<Real> &g,
135 Real ftol,
136 Objective<Real> &obj,
138 std::ostream &outStream) {
139 //const Real one(1);
140 if (proj_ == nullPtr)
141 proj_ = makePtr<PolyhedralProjection<Real>>(makePtrFromRef(bnd));
142 // Initialize data
144 nhess_ = 0;
145 // Update approximate gradient and approximate objective function.
146 proj_->project(x,outStream); state_->nproj++;
147 state_->iterateVec->set(x);
148 obj.update(x,UpdateType::Initial,state_->iter);
149 state_->value = obj.value(x,ftol);
150 state_->nfval++;
151 //obj.gradient(*state_->gradientVec,x,ftol);
152 state_->gnorm = computeGradient(x,*state_->gradientVec,*state_->stepVec,state_->searchSize,obj,outStream);
153 state_->ngrad++;
154 //state_->stepVec->set(x);
155 //state_->stepVec->axpy(-one,state_->gradientVec->dual());
156 //proj_->project(*state_->stepVec,outStream); state_->nproj++;
157 //state_->stepVec->axpy(-one,x);
158 //state_->gnorm = state_->stepVec->norm();
159 state_->snorm = ROL_INF<Real>();
160 // Normalize initial CP step length
161 if (normAlpha_) alpha_ /= state_->gradientVec->norm();
162 // Compute initial trust region radius if desired.
163 if ( state_->searchSize <= static_cast<Real>(0) )
164 state_->searchSize = state_->gradientVec->norm();
165}
166
167template<typename Real>
169 Real &outTol,
170 Real pRed,
171 Real &fold,
172 int iter,
173 const Vector<Real> &x,
174 const Vector<Real> &xold,
175 Objective<Real> &obj) {
176 outTol = std::sqrt(ROL_EPSILON<Real>());
177 if ( useInexact_[0] ) {
178 if (!(iter%updateIter_) && (iter!=0)) force_ *= forceFactor_;
179 const Real one(1);
180 Real eta = static_cast<Real>(0.999)*std::min(eta1_,one-eta2_);
181 outTol = scale_*std::pow(eta*std::min(pRed,force_),one/omega_);
182 if (inTol > outTol) fold = obj.value(xold,outTol);
183 }
184 // Evaluate objective function at new iterate
186 Real fval = obj.value(x,outTol);
187 return fval;
188}
189
190template<typename Real>
192 Vector<Real> &g,
193 Vector<Real> &pwa,
194 Real del,
195 Objective<Real> &obj,
196 std::ostream &outStream) const {
197 Real gnorm(0);
198 if ( useInexact_[1] ) {
199 const Real one(1);
200 Real gtol1 = scale0_*del;
201 Real gtol0 = gtol1 + one;
202 while ( gtol0 > gtol1 ) {
203 obj.gradient(g,x,gtol1);
204 gnorm = TypeB::Algorithm<Real>::optimalityCriterion(x,g,pwa,outStream);
205 gtol0 = gtol1;
206 gtol1 = scale0_*std::min(gnorm,del);
207 }
208 }
209 else {
210 Real gtol = std::sqrt(ROL_EPSILON<Real>());
211 obj.gradient(g,x,gtol);
212 gnorm = TypeB::Algorithm<Real>::optimalityCriterion(x,g,pwa,outStream);
213 }
214 return gnorm;
215}
216
217template<typename Real>
219 const Vector<Real> &g,
220 Objective<Real> &obj,
222 std::ostream &outStream ) {
223 const Real zero(0), one(1);
224 //Real tol0 = std::sqrt(ROL_EPSILON<Real>());
225 Real inTol = static_cast<Real>(0.1)*ROL_OVERFLOW<Real>(), outTol(inTol);
226 Real ftrial(0), pRed(0), rho(1), q(0);
227 // Initialize trust-region data
228 std::vector<std::string> output;
229 initialize(x,g,inTol,obj,bnd,outStream);
230 Ptr<Vector<Real>> gmod = g.clone();
231 Ptr<Vector<Real>> pwa1 = x.clone(), pwa2 = x.clone();
232 Ptr<Vector<Real>> pwa3 = x.clone(), pwa4 = x.clone();
233 Ptr<Vector<Real>> pwa5 = x.clone(), pwa6 = x.clone();
234 Ptr<Vector<Real>> pwa7 = x.clone();
235 Ptr<Vector<Real>> dwa1 = g.clone(), dwa2 = g.clone();
236 // Initialize nonmonotone data
237 Real rhoNM(0), sigmac(0), sigmar(0);
238 Real fr(state_->value), fc(state_->value), fmin(state_->value);
239 TRUtils::ETRFlag TRflagNM;
240 int L(0);
241
242 // Output
243 if (verbosity_ > 0) writeOutput(outStream,true);
244
245 while (status_->check(*state_)) {
246 // Build trust-region model
247 model_->setData(obj,*state_->iterateVec,*state_->gradientVec);
248
249 /**** SOLVE TRUST-REGION SUBPROBLEM ****/
250 q = zero;
251 gmod->set(*state_->gradientVec);
252 if (useSimpleSPG_)
253 dpsg_simple(x,q,*gmod,*state_->iterateVec,state_->searchSize,*model_,
254 *pwa1,*pwa2,*dwa1,outStream);
255 else {
256 // Compute Cauchy point (TRON notation: x = x[1])
257 dcauchy(*state_->stepVec,alpha_,q,*state_->iterateVec,
258 state_->gradientVec->dual(),state_->searchSize,
259 *model_,*dwa1,*dwa2,outStream); // Solve 1D optimization problem for alpha
260 x.plus(*state_->stepVec); // Set x = x[0] + alpha*g
261
262 // Model gradient at s = x[1] - x[0]
263 gmod->plus(*dwa1); // hessVec from Cauchy point computation
264
265 // Apply SPG starting from the Cauchy point
266 dpsg(x,q,*gmod,*state_->iterateVec,state_->searchSize,*model_,
267 *pwa1,*pwa2,*pwa3,*pwa4,*pwa5,*pwa6,*pwa7,*dwa1,outStream);
268 }
269
270 // Update storage and compute predicted reduction
271 pRed = -q;
272 state_->stepVec->set(x); state_->stepVec->axpy(-one,*state_->iterateVec);
273 state_->snorm = state_->stepVec->norm();
274
275 // Compute trial objective value
276 ftrial = computeValue(inTol,outTol,pRed,state_->value,state_->iter,x,*state_->iterateVec,obj);
277 state_->nfval++;
278
279 // Compute ratio of acutal and predicted reduction
280 TRflag_ = TRUtils::SUCCESS;
281 TRUtils::analyzeRatio<Real>(rho,TRflag_,state_->value,ftrial,pRed,eps_,outStream,verbosity_>1);
282 if (useNM_) {
283 TRUtils::analyzeRatio<Real>(rhoNM,TRflagNM,fr,ftrial,pRed+sigmar,eps_,outStream,verbosity_>1);
284 TRflag_ = (rho < rhoNM ? TRflagNM : TRflag_);
285 rho = (rho < rhoNM ? rhoNM : rho );
286 }
287
288 // Update algorithm state
289 state_->iter++;
290 // Accept/reject step and update trust region radius
291 if ((rho < eta0_ && TRflag_ == TRUtils::SUCCESS) || (TRflag_ >= 2)) { // Step Rejected
292 x.set(*state_->iterateVec);
293 obj.update(x,UpdateType::Revert,state_->iter);
294 if (interpRad_ && (rho < zero && TRflag_ != TRUtils::TRNAN)) {
295 // Negative reduction, interpolate to find new trust-region radius
296 state_->searchSize = TRUtils::interpolateRadius<Real>(*state_->gradientVec,*state_->stepVec,
297 state_->snorm,pRed,state_->value,ftrial,state_->searchSize,gamma0_,gamma1_,eta2_,
298 outStream,verbosity_>1);
299 }
300 else { // Shrink trust-region radius
301 state_->searchSize = gamma1_*std::min(state_->snorm,state_->searchSize);
302 }
303 }
304 else if ((rho >= eta0_ && TRflag_ != TRUtils::NPOSPREDNEG)
305 || (TRflag_ == TRUtils::POSPREDNEG)) { // Step Accepted
306 state_->value = ftrial;
307 obj.update(x,UpdateType::Accept,state_->iter);
308 inTol = outTol;
309 if (useNM_) {
310 sigmac += pRed; sigmar += pRed;
311 if (ftrial < fmin) { fmin = ftrial; fc = fmin; sigmac = zero; L = 0; }
312 else {
313 L++;
314 if (ftrial > fc) { fc = ftrial; sigmac = zero; }
315 if (L == storageNM_) { fr = fc; sigmar = sigmac; }
316 }
317 }
318 // Increase trust-region radius
319 if (rho >= eta2_) state_->searchSize = std::min(gamma2_*state_->searchSize, delMax_);
320 // Compute gradient at new iterate
321 dwa1->set(*state_->gradientVec);
322 state_->gnorm = computeGradient(x,*state_->gradientVec,*pwa1,state_->searchSize,obj,outStream);
323 state_->ngrad++;
324 state_->iterateVec->set(x);
325 // Update secant information in trust-region model
326 model_->update(x,*state_->stepVec,*dwa1,*state_->gradientVec,
327 state_->snorm,state_->iter);
328 }
329
330 // Update Output
331 if (verbosity_ > 0) writeOutput(outStream,writeHeader_);
332 }
333 if (verbosity_ > 0) TypeB::Algorithm<Real>::writeExitStatus(outStream);
334}
335
336template<typename Real>
338 const Vector<Real> &x, const Real alpha,
339 std::ostream &outStream) const {
340 s.set(x); s.axpy(alpha,w);
341 proj_->project(s,outStream); state_->nproj++;
342 s.axpy(static_cast<Real>(-1),x);
343 return s.norm();
344}
345
346template<typename Real>
348 Real &alpha,
349 Real &q,
350 const Vector<Real> &x,
351 const Vector<Real> &g,
352 const Real del,
354 Vector<Real> &dwa, Vector<Real> &dwa1,
355 std::ostream &outStream) {
356 const Real half(0.5);
357 // const Real zero(0); // Unused
358 Real tol = std::sqrt(ROL_EPSILON<Real>());
359 bool interp = false;
360 Real gs(0), snorm(0);
361 // Compute s = P(x[0] - alpha g[0])
362 snorm = dgpstep(s,g,x,-alpha,outStream);
363 if (snorm > del) {
364 interp = true;
365 }
366 else {
367 model.hessVec(dwa,s,x,tol); nhess_++;
368 gs = s.dot(g);
369 //q = half * s.dot(dwa.dual()) + gs;
370 q = half * s.apply(dwa) + gs;
371 interp = (q > mu0_*gs);
372 }
373 // Either increase or decrease alpha to find approximate Cauchy point
374 int cnt = 0;
375 if (interp) {
376 bool search = true;
377 while (search) {
378 alpha *= interpf_;
379 snorm = dgpstep(s,g,x,-alpha,outStream);
380 if (snorm <= del) {
381 model.hessVec(dwa,s,x,tol); nhess_++;
382 gs = s.dot(g);
383 //q = half * s.dot(dwa.dual()) + gs;
384 q = half * s.apply(dwa) + gs;
385 search = (q > mu0_*gs) && (cnt < redlim_);
386 }
387 cnt++;
388 }
389 }
390 else {
391 bool search = true;
392 Real alphas = alpha;
393 Real qs = q;
394 dwa1.set(dwa);
395 while (search) {
396 alpha *= extrapf_;
397 snorm = dgpstep(s,g,x,-alpha,outStream);
398 if (snorm <= del && cnt < explim_) {
399 model.hessVec(dwa,s,x,tol); nhess_++;
400 gs = s.dot(g);
401 //q = half * s.dot(dwa.dual()) + gs;
402 q = half * s.apply(dwa) + gs;
403 if (q <= mu0_*gs && std::abs(q-qs) > qtol_*std::abs(qs)) {
404 dwa1.set(dwa);
405 search = true;
406 alphas = alpha;
407 qs = q;
408 }
409 else {
410 q = qs;
411 dwa.set(dwa1);
412 search = false;
413 }
414 }
415 else {
416 search = false;
417 }
418 cnt++;
419 }
420 alpha = alphas;
421 snorm = dgpstep(s,g,x,-alpha,outStream);
422 }
423 if (verbosity_ > 1) {
424 outStream << " Cauchy point" << std::endl;
425 outStream << " Step length (alpha): " << alpha << std::endl;
426 outStream << " Step length (alpha*g): " << snorm << std::endl;
427 outStream << " Model decrease (pRed): " << -q << std::endl;
428 if (!interp) {
429 outStream << " Number of extrapolation steps: " << cnt << std::endl;
430 }
431 }
432 return snorm;
433}
434
435template<typename Real>
437 Real &q,
438 Vector<Real> &gmod,
439 const Vector<Real> &x,
440 Real del,
442 Vector<Real> &pwa,
443 Vector<Real> &pwa1,
444 Vector<Real> &dwa,
445 std::ostream &outStream) {
446 // Use SPG to approximately solve TR subproblem:
447 // min 1/2 <H(y-x), (y-x)> + <g, (y-x)> subject to y\in C, ||y|| \le del
448 //
449 // Inpute:
450 // y = Primal vector
451 // x = Current iterate
452 // g = Current gradient
453 const Real half(0.5), one(1), safeguard(1e2*ROL_EPSILON<Real>());
454 Real tol(std::sqrt(ROL_EPSILON<Real>()));
455 Real alpha(1), alphaMax(1), s0s0(0), ss0(0), sHs(0), lambdaTmp(1), snorm(0);
456 pwa1.zero();
457
458 // Set y = x
459 y.set(x);
460
461 // Compute initial step
462 Real coeff = one/gmod.norm();
463 Real lambda = std::max(lambdaMin_,std::min(coeff,lambdaMax_));
464 pwa.set(y); pwa.axpy(-lambda,gmod.dual()); // pwa = x - lambda gmod.dual()
465 proj_->project(pwa,outStream); state_->nproj++; // pwa = P(x - lambda gmod.dual())
466 pwa.axpy(-one,y); // pwa = P(x - lambda gmod.dual()) - x = step
467 Real gs = gmod.apply(pwa); // gs = <step, model gradient>
468 Real ss = pwa.dot(pwa); // Norm squared of step
469 Real gnorm = std::sqrt(ss);
470
471 // Compute initial projected gradient norm
472 const Real gtol = std::min(tol1_,tol2_*gnorm);
473
474 if (verbosity_ > 1)
475 outStream << " Spectral Projected Gradient" << std::endl;
476
477 SPiter_ = 0;
478 while (SPiter_ < maxit_) {
479 SPiter_++;
480
481 // Evaluate model Hessian
482 model.hessVec(dwa,pwa,x,tol); nhess_++; // dwa = H step
483 sHs = dwa.apply(pwa); // sHs = <step, H step>
484
485 // Perform line search
486 alphaMax = 1;
487 if (gnorm >= del-safeguard) { // Trust-region constraint is violated
488 ss0 = pwa1.dot(pwa);
489 alphaMax = std::min(one, (-ss0 + std::sqrt(ss0*ss0 - ss*(s0s0-del*del)))/ss);
490 }
491 if (sHs <= safeguard)
492 alpha = alphaMax;
493 else
494 alpha = std::min(alphaMax, -gs/sHs);
495
496 // Update model quantities
497 q += alpha * (gs + half * alpha * sHs); // Update model value
498 gmod.axpy(alpha,dwa); // Update model gradient
499 y.axpy(alpha,pwa); // New iterate
500
501 // Check trust-region constraint violation
502 pwa1.set(y); pwa1.axpy(-one,x);
503 s0s0 = pwa1.dot(pwa1);
504 snorm = std::sqrt(s0s0);
505
506 if (verbosity_ > 1) {
507 outStream << std::endl;
508 outStream << " Iterate: " << SPiter_ << std::endl;
509 outStream << " Spectral step length (lambda): " << lambda << std::endl;
510 outStream << " Step length (alpha): " << alpha << std::endl;
511 outStream << " Model decrease (pRed): " << -q << std::endl;
512 outStream << " Optimality criterion: " << gnorm << std::endl;
513 outStream << " Step norm: " << snorm << std::endl;
514 outStream << std::endl;
515 }
516
517 if (snorm >= del - safeguard) { SPflag_ = 2; break; }
518
519 // Compute new spectral step
520 lambdaTmp = (sHs <= safeguard ? one/gmod.norm() : ss/sHs);
521 lambda = std::max(lambdaMin_,std::min(lambdaTmp,lambdaMax_));
522 pwa.set(y); pwa.axpy(-lambda,gmod.dual());
523 proj_->project(pwa,outStream); state_->nproj++;
524 pwa.axpy(-one,y);
525 gs = gmod.apply(pwa);
526 ss = pwa.dot(pwa);
527 gnorm = std::sqrt(ss);
528
529 if (gnorm <= gtol) { SPflag_ = 0; break; }
530 }
531 SPflag_ = (SPiter_==maxit_) ? 1 : SPflag_;
532}
533
534template<typename Real>
536 Real &q,
537 Vector<Real> &gmod,
538 const Vector<Real> &x,
539 Real del,
541 Vector<Real> &ymin,
542 Vector<Real> &pwa,
543 Vector<Real> &pwa1,
544 Vector<Real> &pwa2,
545 Vector<Real> &pwa3,
546 Vector<Real> &pwa4,
547 Vector<Real> &pwa5,
548 Vector<Real> &dwa,
549 std::ostream &outStream) {
550 // Use SPG to approximately solve TR subproblem:
551 // min 1/2 <H(y-x), (y-x)> + <g, (y-x)> subject to y\in C, ||y|| \le del
552 //
553 // Inpute:
554 // y = Cauchy step
555 // x = Current iterate
556 // g = Current gradient
557 const Real zero(0), half(0.5), one(1), two(2); //, eps(std::sqrt(ROL_EPSILON<Real>()));
558 Real tol(std::sqrt(ROL_EPSILON<Real>()));
559 Real alpha(1), sHs(0), alphaTmp(1), mmax(0), qmin(0), lambdaTmp(1);
560 std::deque<Real> mqueue; mqueue.push_back(q);
561
562 if (useNMSP_ && useMin_) { qmin = q; ymin.set(y); }
563
564 // Compute initial projected gradient norm
565 pwa1.set(gmod.dual());
566 pwa.set(y); pwa.axpy(-one,pwa1);
567 dproj(pwa,x,del,pwa2,pwa3,pwa4,pwa5,outStream);
568 pwa.axpy(-one,y);
569 Real gnorm = pwa.norm();
570 const Real gtol = std::min(tol1_,tol2_*gnorm);
571
572 // Compute initial step
573 Real coeff = one/gmod.norm();
574 Real lambda = std::max(lambdaMin_,std::min(coeff,lambdaMax_));
575 pwa.set(y); pwa.axpy(-lambda,pwa1); // pwa = y - lambda gmod.dual()
576 dproj(pwa,x,del,pwa2,pwa3,pwa4,pwa5,outStream); // pwa = P(y - lambda gmod.dual())
577 pwa.axpy(-one,y); // pwa = P(y - lambda gmod.dual()) - y = step
578 Real gs = gmod.apply(pwa); // gs = <step, model gradient>
579 Real ss = pwa.dot(pwa); // Norm squared of step
580
581 if (verbosity_ > 1)
582 outStream << " Spectral Projected Gradient" << std::endl;
583
584 SPiter_ = 0;
585 while (SPiter_ < maxit_) {
586 SPiter_++;
587
588 // Evaluate model Hessian
589 model.hessVec(dwa,pwa,x,tol); nhess_++; // dwa = H step
590 sHs = dwa.apply(pwa); // sHs = <step, H step>
591
592 // Perform line search
593 if (useNMSP_) { // Nonmonotone
594 mmax = *std::max_element(mqueue.begin(),mqueue.end());
595 alphaTmp = (-(one-gamma_)*gs + std::sqrt(std::pow((one-gamma_)*gs,two)-two*sHs*(q-mmax)))/sHs;
596 }
597 else { // Exact
598 alphaTmp = -gs/sHs;
599 }
600 alpha = (sHs > zero ? std::min(one,std::max(zero,alphaTmp)) : one);
601
602 // Update model quantities
603 q += alpha * (gs + half * alpha * sHs); // Update model value
604 gmod.axpy(alpha,dwa); // Update model gradient
605 y.axpy(alpha,pwa); // New iterate
606
607 // Update nonmonotone line search information
608 if (useNMSP_) {
609 if (static_cast<int>(mqueue.size())==maxSize_) mqueue.pop_front();
610 mqueue.push_back(q);
611 if (useMin_ && q <= qmin) { qmin = q; ymin.set(y); }
612 }
613
614 // Compute projected gradient norm
615 pwa1.set(gmod.dual());
616 pwa.set(y); pwa.axpy(-one,pwa1);
617 dproj(pwa,x,del,pwa2,pwa3,pwa4,pwa5,outStream);
618 pwa.axpy(-one,y);
619 gnorm = pwa.norm();
620
621 if (verbosity_ > 1) {
622 outStream << std::endl;
623 outStream << " Iterate: " << SPiter_ << std::endl;
624 outStream << " Spectral step length (lambda): " << lambda << std::endl;
625 outStream << " Step length (alpha): " << alpha << std::endl;
626 outStream << " Model decrease (pRed): " << -q << std::endl;
627 outStream << " Optimality criterion: " << gnorm << std::endl;
628 outStream << std::endl;
629 }
630 if (gnorm < gtol) break;
631
632 // Compute new spectral step
633 //lambda = (sHs<=eps ? lambdaMax_ : std::max(lambdaMin_,std::min(ss/sHs,lambdaMax_)));
634 lambdaTmp = (sHs == 0 ? coeff : ss/sHs);
635 lambda = std::max(lambdaMin_,std::min(lambdaTmp,lambdaMax_));
636 pwa.set(y); pwa.axpy(-lambda,pwa1);
637 dproj(pwa,x,del,pwa2,pwa3,pwa4,pwa5,outStream);
638 pwa.axpy(-one,y);
639 gs = gmod.apply(pwa);
640 ss = pwa.dot(pwa);
641 }
642 if (useNMSP_ && useMin_) { q = qmin; y.set(ymin); }
643 SPflag_ = (SPiter_==maxit_) ? 1 : 0;
644}
645
646template<typename Real>
648 const Vector<Real> &x0,
649 Real del,
650 Vector<Real> &y0,
651 Vector<Real> &y1,
652 Vector<Real> &yc,
653 Vector<Real> &pwa,
654 std::ostream &outStream) const {
655 // Solve ||P(t*x0 + (1-t)*(x-x0))-x0|| = del using Brent's method
656 const Real zero(0), half(0.5), one(1), two(2), three(3);
657 const Real eps(ROL_EPSILON<Real>()), tol0(1e1*eps), fudge(1.0-1e-2*sqrt(eps));
658 Real f0(0), f1(0), fc(0), t0(0), t1(1), tc(0), d1(1), d2(1), tol(1);
659 Real p(0), q(0), r(0), s(0), m(0);
660 int cnt(state_->nproj);
661 y1.set(x);
662 proj_->project(y1,outStream); state_->nproj++;
663 pwa.set(y1); pwa.axpy(-one,x0);
664 f1 = pwa.norm();
665 if (f1 <= del) {
666 x.set(y1);
667 return;
668 }
669 y0.set(x0);
670 tc = t0; fc = f0; yc.set(y0);
671 d1 = t1-t0; d2 = d1;
672 int code = 0;
673 while (true) {
674 if (std::abs(fc-del) < std::abs(f1-del)) {
675 t0 = t1; t1 = tc; tc = t0;
676 f0 = f1; f1 = fc; fc = f0;
677 y0.set(y1); y1.set(yc); yc.set(y0);
678 }
679 tol = two*eps*std::abs(t1) + half*tol0;
680 m = half*(tc - t1);
681 if (std::abs(m) <= tol) { code = 1; break; }
682 if ((f1 >= fudge*del && f1 <= del)) break;
683 if (std::abs(d1) < tol || std::abs(f0-del) <= std::abs(f1-del)) {
684 d1 = m; d2 = d1;
685 }
686 else {
687 s = (f1-del)/(f0-del);
688 if (t0 == tc) {
689 p = two*m*s;
690 q = one-s;
691 }
692 else {
693 q = (f0-del)/(fc-del);
694 r = (f1-del)/(fc-del);
695 p = s*(two*m*q*(q-r)-(t1-t0)*(r-one));
696 q = (q-one)*(r-one)*(s-one);
697 }
698 if (p > zero) q = -q;
699 else p = -p;
700 s = d1;
701 d1 = d2;
702 if (two*p < three*m*q-std::abs(tol*q) && p < std::abs(half*s*q)) {
703 d2 = p/q;
704 }
705 else {
706 d1 = m; d2 = d1;
707 }
708 }
709 t0 = t1; f0 = f1; y0.set(y1);
710 if (std::abs(d2) > tol) t1 += d2;
711 else if (m > zero) t1 += tol;
712 else t1 -= tol;
713 y1.set(x); y1.scale(t1); y1.axpy(one-t1,x0);
714 proj_->project(y1,outStream); state_->nproj++;
715 pwa.set(y1); pwa.axpy(-one,x0);
716 f1 = pwa.norm();
717 if ((f1 > del && fc > del) || (f1 <= del && fc <= del)) {
718 tc = t0; fc = f0; yc.set(y0);
719 d1 = t1-t0; d2 = d1;
720 }
721 }
722 if (code==1 && f1>del) x.set(yc);
723 else x.set(y1);
724 if (verbosity_ > 1) {
725 outStream << std::endl;
726 outStream << " Trust-Region Subproblem Projection" << std::endl;
727 outStream << " Number of polyhedral projections: " << state_->nproj-cnt << std::endl;
728 if (code == 1 && f1 > del) {
729 outStream << " Transformed Multiplier: " << tc << std::endl;
730 outStream << " Dual Residual: " << fc-del << std::endl;
731 }
732 else {
733 outStream << " Transformed Multiplier: " << t1 << std::endl;
734 outStream << " Dual Residual: " << f1-del << std::endl;
735 }
736 outStream << " Exit Code: " << code << std::endl;
737 outStream << std::endl;
738 }
739}
740
741// BRACKETING AND BRENTS FOR UNTRANSFORMED MULTIPLIER
742//template<typename Real>
743//void TrustRegionSPGAlgorithm<Real>::dproj(Vector<Real> &x,
744// const Vector<Real> &x0,
745// Real del,
746// Vector<Real> &y0,
747// Vector<Real> &y1,
748// Vector<Real> &yc,
749// Vector<Real> &pwa,
750// std::ostream &outStream) const {
751// // Solve ||P(t*x0 + (1-t)*(x-x0))-x0|| = del using Brent's method
752// const Real zero(0), half(0.5), one(1), two(2), three(3);
753// const Real eps(ROL_EPSILON<Real>()), tol0(1e1*eps), fudge(1.0-1e-2*sqrt(eps));
754// Real f0(0), f1(0), fc(0), u0(0), u1(0), uc(0), t0(1), t1(0), tc(0), d1(1), d2(1), tol(1);
755// Real p(0), q(0), r(0), s(0), m(0);
756// int cnt(state_->nproj);
757// y0.set(x);
758// proj_->project(y0,outStream); state_->nproj++;
759// pwa.set(y0); pwa.axpy(-one,x0);
760// f0 = pwa.norm();
761// if (f0 <= del) {
762// x.set(y0);
763// return;
764// }
765//
766// // Bracketing
767// t1 = static_cast<Real>(1e-1);
768// f1 = one+del;
769// while (f1 >= del) {
770// t1 *= static_cast<Real>(5e-2);
771// y1.set(x); y1.scale(t1); y1.axpy(one-t1,x0);
772// proj_->project(y1,outStream); state_->nproj++;
773// pwa.set(y1); pwa.axpy(-one,x0);
774// f1 = pwa.norm();
775// }
776// u1 = (one-t1)/t1;
777//
778// // Brents
779// uc = u0; tc = t0; fc = f0; yc.set(y0);
780// d1 = u1-u0; d2 = d1;
781// int code = 0;
782// while (true) {
783// if (std::abs(fc-del) < std::abs(f1-del)) {
784// u0 = u1; u1 = uc; uc = u0;
785// t0 = t1; t1 = tc; tc = t0;
786// f0 = f1; f1 = fc; fc = f0;
787// y0.set(y1); y1.set(yc); yc.set(y0);
788// }
789// tol = two*eps*abs(u1) + half*tol0;
790// m = half*(uc - u1);
791// if (std::abs(m) <= tol) { code = 1; break; }
792// if ((f1 >= fudge*del && f1 <= del)) break;
793// if (std::abs(d1) < tol || std::abs(f0-del) <= std::abs(f1-del)) {
794// d1 = m; d2 = d1;
795// }
796// else {
797// s = (f1-del)/(f0-del);
798// if (u0 == uc) {
799// p = two*m*s;
800// q = one-s;
801// }
802// else {
803// q = (f0-del)/(fc-del);
804// r = (f1-del)/(fc-del);
805// p = s*(two*m*q*(q-r)-(u1-u0)*(r-one));
806// q = (q-one)*(r-one)*(s-one);
807// }
808// if (p > zero) q = -q;
809// else p = -p;
810// s = d1;
811// d1 = d2;
812// if (two*p < three*m*q-std::abs(tol*q) && p < std::abs(half*s*q)) {
813// d2 = p/q;
814// }
815// else {
816// d1 = m; d2 = d1;
817// }
818// }
819// u0 = u1; t0 = t1; f0 = f1; y0.set(y1);
820// if (std::abs(d2) > tol) u1 += d2;
821// else if (m > zero) u1 += tol;
822// else u1 -= tol;
823// t1 = one/(one+u1);
824// y1.set(x); y1.scale(t1); y1.axpy(one-t1,x0);
825// proj_->project(y1,outStream); state_->nproj++;
826// pwa.set(y1); pwa.axpy(-one,x0);
827// f1 = pwa.norm();
828// if ((f1 > del && fc > del) || (f1 <= del && fc <= del)) {
829// uc = u0; tc = t0; fc = f0; yc.set(y0);
830// d1 = u1-u0; d2 = d1;
831// }
832// }
833// if (code==1 && f1>del) x.set(yc);
834// else x.set(y1);
835// if (verbosity_ > 1) {
836// outStream << std::endl;
837// outStream << " Trust-Region Subproblem Projection" << std::endl;
838// outStream << " Number of polyhedral projections: " << state_->nproj-cnt << std::endl;
839// if (code == 1 && f1 > del) {
840// outStream << " Multiplier: " << uc << std::endl;
841// outStream << " Transformed Multiplier: " << tc << std::endl;
842// outStream << " Dual Residual: " << fc-del << std::endl;
843// }
844// else {
845// outStream << " Multiplier: " << u1 << std::endl;
846// outStream << " Transformed Multiplier: " << t1 << std::endl;
847// outStream << " Dual Residual: " << f1-del << std::endl;
848// }
849// outStream << " Exit Code: " << code << std::endl;
850// outStream << std::endl;
851// }
852//}
853
854// RIDDERS' METHOD FOR TRUST-REGION PROJECTION
855//template<typename Real>
856//void TrustRegionSPGAlgorithm<Real>::dproj(Vector<Real> &x,
857// const Vector<Real> &x0,
858// Real del,
859// Vector<Real> &y,
860// Vector<Real> &y1,
861// Vector<Real> &yc,
862// Vector<Real> &p,
863// std::ostream &outStream) const {
864// // Solve ||P(t*x0 + (1-t)*(x-x0))-x0|| = del using Ridder's method
865// const Real half(0.5), one(1), tol(1e1*ROL_EPSILON<Real>());
866// const Real fudge(1.0-1e-2*std::sqrt(ROL_EPSILON<Real>()));
867// Real e0(0), e1(0), e2(0), e(0), a0(0), a1(0.5), a2(1), a(0);
868// int cnt(state_->nproj);
869// y.set(x);
870// proj_->project(y,outStream); state_->nproj++;
871// p.set(y); p.axpy(-one,x0);
872// e2 = p.norm();
873// if (e2 <= del) {
874// x.set(y);
875// return;
876// }
877// bool code = 1;
878// while (a2-a0 > tol) {
879// a1 = half*(a0+a2);
880// y.set(x); y.scale(a1); y.axpy(one-a1,x0);
881// proj_->project(y,outStream); state_->nproj++;
882// p.set(y); p.axpy(-one,x0);
883// e1 = p.norm();
884// if (e1 >= fudge*del && e1 <= del) break;
885// a = a1-(a1-a0)*(e1-del)/std::sqrt((e1-del)*(e1-del)-(e0-del)*(e2-del));
886// y.set(x); y.scale(a); y.axpy(one-a,x0);
887// proj_->project(y,outStream); state_->nproj++;
888// p.set(y); p.axpy(-one,x0);
889// e = p.norm();
890// if (e < fudge*del) {
891// if (e1 < fudge*del) { e0 = (a < a1 ? e1 : e); a0 = (a < a1 ? a1 : a); }
892// else { e0 = e; a0 = a; e2 = e1; a2 = a1; };
893// }
894// else if (e > del) {
895// if (e1 < fudge*del) { e0 = e1; a0 = a1; e2 = e; a2 = a; }
896// else { e2 = (a < a1 ? e : e1); a2 = (a < a1 ? a : a1); }
897// }
898// else {
899// code = 0;
900// break; // Exit if fudge*del <= snorm <= del
901// }
902// }
903// x.set(y);
904// if (verbosity_ > 1) {
905// outStream << std::endl;
906// outStream << " Trust-Region Subproblem Projection" << std::endl;
907// outStream << " Number of polyhedral projections: " << state_->nproj-cnt << std::endl;
908// outStream << " Transformed Multiplier: " << a1 << std::endl;
909// outStream << " Dual Residual: " << e1-del << std::endl;
910// outStream << " Exit Code: " << code << std::endl;
911// outStream << std::endl;
912// }
913//}
914
915template<typename Real>
916void TrustRegionSPGAlgorithm<Real>::writeHeader( std::ostream& os ) const {
917 std::stringstream hist;
918 if (verbosity_ > 1) {
919 hist << std::string(114,'-') << std::endl;
920 hist << " SPG trust-region method status output definitions" << std::endl << std::endl;
921 hist << " iter - Number of iterates (steps taken)" << std::endl;
922 hist << " value - Objective function value" << std::endl;
923 hist << " gnorm - Norm of the gradient" << std::endl;
924 hist << " snorm - Norm of the step (update to optimization vector)" << std::endl;
925 hist << " delta - Trust-Region radius" << std::endl;
926 hist << " #fval - Number of times the objective function was evaluated" << std::endl;
927 hist << " #grad - Number of times the gradient was computed" << std::endl;
928 hist << " #hess - Number of times the Hessian was applied" << std::endl;
929 hist << " #proj - Number of times the projection was computed" << std::endl;
930 hist << std::endl;
931 hist << " tr_flag - Trust-Region flag" << std::endl;
932 for( int flag = TRUtils::SUCCESS; flag != TRUtils::UNDEFINED; ++flag ) {
933 hist << " " << NumberToString(flag) << " - "
934 << TRUtils::ETRFlagToString(static_cast<TRUtils::ETRFlag>(flag)) << std::endl;
935 }
936 hist << std::endl;
937 hist << " iterSPG - Number of Spectral Projected Gradient iterations" << std::endl << std::endl;
938 hist << " flagSPG - Trust-Region Truncated CG flag" << std::endl;
939 hist << " 0 - Converged" << std::endl;
940 hist << " 1 - Iteration Limit Exceeded" << std::endl;
941 hist << std::string(114,'-') << std::endl;
942 }
943 hist << " ";
944 hist << std::setw(6) << std::left << "iter";
945 hist << std::setw(15) << std::left << "value";
946 hist << std::setw(15) << std::left << "gnorm";
947 hist << std::setw(15) << std::left << "snorm";
948 hist << std::setw(15) << std::left << "delta";
949 hist << std::setw(10) << std::left << "#fval";
950 hist << std::setw(10) << std::left << "#grad";
951 hist << std::setw(10) << std::left << "#hess";
952 hist << std::setw(10) << std::left << "#proj";
953 hist << std::setw(10) << std::left << "tr_flag";
954 hist << std::setw(10) << std::left << "iterSPG";
955 hist << std::setw(10) << std::left << "flagSPG";
956 hist << std::endl;
957 os << hist.str();
958}
959
960template<typename Real>
961void TrustRegionSPGAlgorithm<Real>::writeName( std::ostream& os ) const {
962 std::stringstream hist;
963 hist << std::endl << "SPG Trust-Region Method (Type B, Bound Constraints)" << std::endl;
964 os << hist.str();
965}
966
967template<typename Real>
968void TrustRegionSPGAlgorithm<Real>::writeOutput( std::ostream& os, bool write_header ) const {
969 std::stringstream hist;
970 hist << std::scientific << std::setprecision(6);
971 if ( state_->iter == 0 ) writeName(os);
972 if ( write_header ) writeHeader(os);
973 if ( state_->iter == 0 ) {
974 hist << " ";
975 hist << std::setw(6) << std::left << state_->iter;
976 hist << std::setw(15) << std::left << state_->value;
977 hist << std::setw(15) << std::left << state_->gnorm;
978 hist << std::setw(15) << std::left << "---";
979 hist << std::setw(15) << std::left << state_->searchSize;
980 hist << std::setw(10) << std::left << state_->nfval;
981 hist << std::setw(10) << std::left << state_->ngrad;
982 hist << std::setw(10) << std::left << nhess_;
983 hist << std::setw(10) << std::left << state_->nproj;
984 hist << std::setw(10) << std::left << "---";
985 hist << std::setw(10) << std::left << "---";
986 hist << std::setw(10) << std::left << "---";
987 hist << std::endl;
988 }
989 else {
990 hist << " ";
991 hist << std::setw(6) << std::left << state_->iter;
992 hist << std::setw(15) << std::left << state_->value;
993 hist << std::setw(15) << std::left << state_->gnorm;
994 hist << std::setw(15) << std::left << state_->snorm;
995 hist << std::setw(15) << std::left << state_->searchSize;
996 hist << std::setw(10) << std::left << state_->nfval;
997 hist << std::setw(10) << std::left << state_->ngrad;
998 hist << std::setw(10) << std::left << nhess_;
999 hist << std::setw(10) << std::left << state_->nproj;
1000 hist << std::setw(10) << std::left << TRflag_;
1001 hist << std::setw(10) << std::left << SPiter_;
1002 hist << std::setw(10) << std::left << SPflag_;
1003 hist << std::endl;
1004 }
1005 os << hist.str();
1006}
1007
1008} // namespace TypeB
1009} // namespace ROL
1010
1011#endif
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0 zero)()
Provides the interface to apply upper and lower bound constraints.
Provides the interface to evaluate objective functions.
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
virtual void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
Provides interface for and implements limited-memory secant operators.
Provides an interface to check status of optimization algorithms.
Provides the interface to evaluate trust-region model functions.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol) override
Apply Hessian approximation to vector.
void initialize(const Vector< Real > &x, const Vector< Real > &g)
Real optimalityCriterion(const Vector< Real > &x, const Vector< Real > &g, Vector< Real > &primal, std::ostream &outStream=std::cout) const
virtual void writeExitStatus(std::ostream &os) const
Real computeValue(Real inTol, Real &outTol, Real pRed, Real &fold, int iter, const Vector< Real > &x, const Vector< Real > &xold, Objective< Real > &obj)
void initialize(Vector< Real > &x, const Vector< Real > &g, Real ftol, Objective< Real > &obj, BoundConstraint< Real > &bnd, std::ostream &outStream=std::cout)
Real dcauchy(Vector< Real > &s, Real &alpha, Real &q, const Vector< Real > &x, const Vector< Real > &g, const Real del, TrustRegionModel_U< Real > &model, Vector< Real > &dwa, Vector< Real > &dwa1, std::ostream &outStream=std::cout)
void dpsg_simple(Vector< Real > &y, Real &q, Vector< Real > &gmod, const Vector< Real > &x, Real del, TrustRegionModel_U< Real > &model, Vector< Real > &pwa, Vector< Real > &pwa1, Vector< Real > &dwa, std::ostream &outStream=std::cout)
void dpsg(Vector< Real > &y, Real &q, Vector< Real > &gmod, const Vector< Real > &x, Real del, TrustRegionModel_U< Real > &model, Vector< Real > &ymin, Vector< Real > &pwa, Vector< Real > &pwa1, Vector< Real > &pwa2, Vector< Real > &pwa3, Vector< Real > &pwa4, Vector< Real > &pwa5, Vector< Real > &dwa, std::ostream &outStream=std::cout)
Real computeGradient(const Vector< Real > &x, Vector< Real > &g, Vector< Real > &pwa, Real del, Objective< Real > &obj, std::ostream &outStream=std::cout) const
void run(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &bnd, std::ostream &outStream=std::cout) override
Run algorithm on bound constrained problems (Type-B). This general interface supports the use of dual...
void writeOutput(std::ostream &os, bool write_header=false) const override
Print iterate status.
TrustRegionSPGAlgorithm(ParameterList &list, const Ptr< Secant< Real > > &secant=nullPtr)
void writeName(std::ostream &os) const override
Print step name.
void writeHeader(std::ostream &os) const override
Print iterate header.
void dproj(Vector< Real > &x, const Vector< Real > &x0, Real del, Vector< Real > &y0, Vector< Real > &y1, Vector< Real > &yc, Vector< Real > &pwa, std::ostream &outStream=std::cout) const
Real dgpstep(Vector< Real > &s, const Vector< Real > &w, const Vector< Real > &x, const Real alpha, std::ostream &outStream=std::cout) const
Defines the linear algebra or vector space interface.
virtual Real apply(const Vector< Real > &x) const
Apply to a dual vector. This is equivalent to the call .
virtual Real norm() const =0
Returns where .
virtual void set(const Vector &x)
Set where .
virtual void scale(const Real alpha)=0
Compute where .
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis,...
virtual void plus(const Vector &x)=0
Compute , where .
virtual void zero()
Set to zero vector.
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
virtual Real dot(const Vector &x) const =0
Compute where .
std::string ETRFlagToString(ETRFlag trf)
std::string NumberToString(T Number)
Definition ROL_Types.hpp:81
ESecant StringToESecant(std::string s)
ESecantMode
@ SECANTMODE_FORWARD
@ SECANTMODE_INVERSE
@ SECANTMODE_BOTH