57 ParameterList &trlist = list.sublist(
"Step").sublist(
"Trust Region");
59 state_->searchSize = trlist.get(
"Initial Radius", -1.0);
60 delMax_ = trlist.get(
"Maximum Radius", ROL_INF<Real>());
61 eta0_ = trlist.get(
"Step Acceptance Threshold", 0.05);
62 eta1_ = trlist.get(
"Radius Shrinking Threshold", 0.05);
63 eta2_ = trlist.get(
"Radius Growing Threshold", 0.9);
64 gamma0_ = trlist.get(
"Radius Shrinking Rate (Negative rho)", 0.0625);
65 gamma1_ = trlist.get(
"Radius Shrinking Rate (Positive rho)", 0.25);
66 gamma2_ = trlist.get(
"Radius Growing Rate", 2.5);
67 TRsafe_ = trlist.get(
"Safeguard Size", 100.0);
68 eps_ = TRsafe_*ROL_EPSILON<Real>();
69 interpRad_ = trlist.get(
"Use Radius Interpolation",
false);
71 storageNM_ = trlist.get(
"Nonmonotone Storage Size", 0);
72 useNM_ = (storageNM_ <= 0 ? false :
true);
74 maxit_ = list.sublist(
"General").sublist(
"Krylov").get(
"Iteration Limit", 20);
75 tol1_ = list.sublist(
"General").sublist(
"Krylov").get(
"Absolute Tolerance", 1e-4);
76 tol2_ = list.sublist(
"General").sublist(
"Krylov").get(
"Relative Tolerance", 1e-2);
78 ROL::ParameterList &lmlist = trlist.sublist(
"Lin-More");
79 minit_ = lmlist.get(
"Maximum Number of Minor Iterations", 10);
80 mu0_ = lmlist.get(
"Sufficient Decrease Parameter", 1e-2);
81 spexp_ = lmlist.get(
"Relative Tolerance Exponent", 1.0);
82 spexp_ = std::max(
static_cast<Real
>(1),std::min(spexp_,
static_cast<Real
>(2)));
83 redlim_ = lmlist.sublist(
"Cauchy Point").get(
"Maximum Number of Reduction Steps", 10);
84 explim_ = lmlist.sublist(
"Cauchy Point").get(
"Maximum Number of Expansion Steps", 10);
85 alpha_ = lmlist.sublist(
"Cauchy Point").get(
"Initial Step Size", 1.0);
86 normAlpha_ = lmlist.sublist(
"Cauchy Point").get(
"Normalize Initial Step Size",
false);
87 interpf_ = lmlist.sublist(
"Cauchy Point").get(
"Reduction Rate", 0.1);
88 extrapf_ = lmlist.sublist(
"Cauchy Point").get(
"Expansion Rate", 10.0);
89 qtol_ = lmlist.sublist(
"Cauchy Point").get(
"Decrease Tolerance", 1e-8);
90 interpfPS_ = lmlist.sublist(
"Projected Search").get(
"Backtracking Rate", 0.5);
91 pslim_ = lmlist.sublist(
"Projected Search").get(
"Maximum Number of Steps", 20);
93 verbosity_ = list.sublist(
"General").get(
"Output Level",0);
94 writeHeader_ = verbosity_ > 2;
96 useSecantPrecond_ = list.sublist(
"General").sublist(
"Secant").get(
"Use as Preconditioner",
false);
97 useSecantHessVec_ = list.sublist(
"General").sublist(
"Secant").get(
"Use as Hessian",
false);
102 model_ = makePtr<TrustRegionModel_U<Real>>(list,secant,mode);
103 if (secant == nullPtr) {
104 esec_ =
StringToESecant(list.sublist(
"General").sublist(
"Secant").get(
"Type",
"Limited-Memory BFGS"));
113 std::ostream &outStream) {
116 if (proj_ == nullPtr) {
117 proj_ = makePtr<PolyhedralProjection<Real>>(makePtrFromRef(bnd));
124 Real ftol =
static_cast<Real
>(0.1)*ROL_OVERFLOW<Real>();
125 proj_->project(x,outStream); state_->nproj++;
126 state_->iterateVec->set(x);
128 state_->value = obj.
value(x,ftol);
130 obj.
gradient(*state_->gradientVec,x,ftol);
132 state_->stepVec->set(x);
133 state_->stepVec->axpy(-one,state_->gradientVec->dual());
134 proj_->project(*state_->stepVec,outStream); state_->nproj++;
135 state_->stepVec->axpy(-one,x);
136 state_->gnorm = state_->stepVec->norm();
137 state_->snorm = ROL_INF<Real>();
140 alpha_ /= state_->gradientVec->norm();
143 if ( state_->searchSize <=
static_cast<Real
>(0) ) {
144 state_->searchSize = state_->gradientVec->norm();
148 rcon_ = makePtr<ReducedLinearConstraint<Real>>(proj_->getLinearConstraint(),
151 ns_ = makePtr<NullSpaceOperator<Real>>(rcon_,x,
152 *proj_->getResidual());
161 std::ostream &outStream ) {
163 Real tol0 = std::sqrt(ROL_EPSILON<Real>());
164 Real gfnorm(0), gfnormf(0), tol(0), stol(0), snorm(0);
165 Real ftrial(0), pRed(0), rho(1), q(0), delta(0);
166 int flagCG(0), iterCG(0), maxit(0);
168 initialize(x,g,obj,bnd,outStream);
169 Ptr<Vector<Real>> s = x.
clone();
170 Ptr<Vector<Real>> gmod = g.
clone(), gfree = g.
clone();
171 Ptr<Vector<Real>> pwa1 = x.
clone(), pwa2 = x.
clone(), pwa3 = x.
clone();
172 Ptr<Vector<Real>> dwa1 = g.
clone(), dwa2 = g.
clone(), dwa3 = g.
clone();
174 Real rhoNM(0), sigmac(0), sigmar(0);
175 Real fr(state_->value), fc(state_->value), fmin(state_->value);
180 if (verbosity_ > 0) writeOutput(outStream,
true);
182 while (status_->check(*state_)) {
184 model_->setData(obj,*state_->iterateVec,*state_->gradientVec);
188 snorm = dcauchy(*state_->stepVec,alpha_,q,*state_->iterateVec,
189 state_->gradientVec->dual(),state_->searchSize,
190 *model_,*dwa1,*dwa2,outStream);
191 x.
plus(*state_->stepVec);
192 state_->snorm = snorm;
193 delta = state_->searchSize - snorm;
198 gmod->plus(*state_->gradientVec);
201 pwa1->set(gfree->dual());
203 gfree->set(pwa1->dual());
205 applyFreePrecond(*pwa1,*gfree,x,*model_,bnd,tol0,*dwa1,*pwa2);
206 gfnorm = pwa1->norm();
209 gfnorm = gfree->norm();
211 SPiter_ = 0; SPflag_ = 0;
212 if (verbosity_ > 1) {
213 outStream <<
" Norm of free gradient components: " << gfnorm << std::endl;
214 outStream << std::endl;
219 for (
int i = 0; i < minit_; ++i) {
221 flagCG = 0; iterCG = 0;
223 tol = std::min(tol1_,tol2_*std::pow(gfnorm,spexp_));
225 if (gfnorm >
zero && delta >
zero) {
226 snorm = dtrpcg(*s,flagCG,iterCG,*gfree,x,
227 delta,*model_,bnd,tol,stol,maxit,
228 *pwa1,*dwa1,*pwa2,*dwa2,*pwa3,*dwa3,outStream);
231 if (verbosity_ > 1) {
232 outStream <<
" Computation of CG step" << std::endl;
233 outStream <<
" Current face (i): " << i << std::endl;
234 outStream <<
" CG step length: " << snorm << std::endl;
235 outStream <<
" Number of CG iterations: " << iterCG << std::endl;
236 outStream <<
" CG flag: " << flagCG << std::endl;
237 outStream <<
" Total number of iterations: " << SPiter_ << std::endl;
238 outStream << std::endl;
242 snorm = dprsrch(x,*s,q,gmod->
dual(),*model_,bnd,*pwa1,*dwa1,outStream);
246 state_->stepVec->plus(*s);
247 state_->snorm = state_->stepVec->norm();
248 delta = state_->searchSize - state_->snorm;
252 pwa1->set(gfree->dual());
254 gfree->set(pwa1->dual());
256 applyFreePrecond(*pwa1,*gfree,x,*model_,bnd,tol0,*dwa1,*pwa2);
257 gfnormf = pwa1->norm();
260 gfnormf = gfree->norm();
262 if (verbosity_ > 1) {
263 outStream <<
" Norm of free gradient components: " << gfnormf << std::endl;
264 outStream << std::endl;
269 if (gfnormf <= tol) {
273 else if (SPiter_ >= maxit_) {
277 else if (flagCG == 2) {
281 else if (delta <=
zero) {
293 ftrial = obj.
value(x,tol0);
298 TRUtils::analyzeRatio<Real>(rho,TRflag_,state_->value,ftrial,pRed,eps_,outStream,verbosity_>1);
300 TRUtils::analyzeRatio<Real>(rhoNM,TRflagNM,fr,ftrial,pRed+sigmar,eps_,outStream,verbosity_>1);
301 TRflag_ = (rho < rhoNM ? TRflagNM : TRflag_);
302 rho = (rho < rhoNM ? rhoNM : rho );
309 x.
set(*state_->iterateVec);
313 state_->searchSize = TRUtils::interpolateRadius<Real>(*state_->gradientVec,*state_->stepVec,
314 state_->snorm,pRed,state_->value,ftrial,state_->searchSize,gamma0_,gamma1_,eta2_,
315 outStream,verbosity_>1);
318 state_->searchSize = gamma1_*std::min(state_->snorm,state_->searchSize);
323 state_->value = ftrial;
326 sigmac += pRed; sigmar += pRed;
327 if (ftrial < fmin) { fmin = ftrial; fc = fmin; sigmac =
zero; L = 0; }
330 if (ftrial > fc) { fc = ftrial; sigmac =
zero; }
331 if (L == storageNM_) { fr = fc; sigmar = sigmac; }
335 if (rho >= eta2_) state_->searchSize = std::min(gamma2_*state_->searchSize, delMax_);
337 dwa1->set(*state_->gradientVec);
338 obj.
gradient(*state_->gradientVec,x,tol0);
341 state_->iterateVec->set(x);
343 model_->update(x,*state_->stepVec,*dwa1,*state_->gradientVec,
344 state_->snorm,state_->iter);
348 if (verbosity_ > 0) writeOutput(outStream,writeHeader_);
519 const Real tol,
const Real stol,
const int itermax,
522 std::ostream &outStream)
const {
527 Real tol0 = std::sqrt(ROL_EPSILON<Real>());
528 const Real
zero(0), one(1), two(2);
529 Real rho(0), kappa(0), beta(0), sigma(0), alpha(0);
530 Real rtr(0), tnorm(0), sMs(0), pMp(0), sMp(0);
537 applyFreePrecond(r,t,x,model,bnd,tol0,dwa,pwa);
542 pMp = (!hasEcon_ ? rho : p.
dot(p));
544 for (iter = 0; iter < itermax; ++iter) {
546 applyFreeHessian(q,p,x,model,bnd,tol0,pwa);
550 alpha = (kappa>
zero) ? rho/kappa :
zero;
551 sigma = dtrqsol(sMs,pMp,sMp,del);
553 if (kappa <= zero || alpha >= sigma) {
556 iflag = (kappa<=
zero) ? 2 : 3;
562 applyFreePrecond(r,t,x,model,bnd,tol0,dwa,pwa);
567 if (rtr <= stol*stol || tnorm <= tol) {
568 sMs = sMs + two*alpha*sMp + alpha*alpha*pMp;
580 sMs = sMs + two*alpha*sMp + alpha*alpha*pMp;
581 sMp = beta*(sMp + alpha*pMp);
582 pMp = (!hasEcon_ ? rho : p.
dot(p)) + beta*beta*pMp;
585 if (iter == itermax) {
591 return std::sqrt(sMs);
705 std::stringstream hist;
706 hist << std::scientific << std::setprecision(6);
707 if ( state_->iter == 0 ) writeName(os);
708 if ( write_header ) writeHeader(os);
709 if ( state_->iter == 0 ) {
711 hist << std::setw(6) << std::left << state_->iter;
712 hist << std::setw(15) << std::left << state_->value;
713 hist << std::setw(15) << std::left << state_->gnorm;
714 hist << std::setw(15) << std::left <<
"---";
715 hist << std::setw(15) << std::left << state_->searchSize;
716 hist << std::setw(10) << std::left << state_->nfval;
717 hist << std::setw(10) << std::left << state_->ngrad;
718 hist << std::setw(10) << std::left << nhess_;
719 hist << std::setw(10) << std::left << state_->nproj;
720 hist << std::setw(10) << std::left <<
"---";
722 hist << std::setw(10) << std::left <<
"---";
723 hist << std::setw(10) << std::left <<
"---";
729 hist << std::setw(6) << std::left << state_->iter;
730 hist << std::setw(15) << std::left << state_->value;
731 hist << std::setw(15) << std::left << state_->gnorm;
732 hist << std::setw(15) << std::left << state_->snorm;
733 hist << std::setw(15) << std::left << state_->searchSize;
734 hist << std::setw(10) << std::left << state_->nfval;
735 hist << std::setw(10) << std::left << state_->ngrad;
736 hist << std::setw(10) << std::left << nhess_;
737 hist << std::setw(10) << std::left << state_->nproj;
738 hist << std::setw(10) << std::left << TRflag_;
740 hist << std::setw(10) << std::left << SPiter_;
741 hist << std::setw(10) << std::left << SPflag_;