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