44 #ifndef ROL_TYPEB_LINMOREALGORITHM_DEF_HPP
45 #define ROL_TYPEB_LINMOREALGORITHM_DEF_HPP
50 template<
typename Real>
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 ParameterList &glist = list.sublist(
"General");
95 useInexact_.push_back(glist.get(
"Inexact Objective Function",
false));
96 useInexact_.push_back(glist.get(
"Inexact Gradient",
false));
97 useInexact_.push_back(glist.get(
"Inexact Hessian-Times-A-Vector",
false));
99 ParameterList &ilist = trlist.sublist(
"Inexact").sublist(
"Gradient");
100 scale0_ = ilist.get(
"Tolerance Scaling", static_cast<Real>(0.1));
101 scale1_ = ilist.get(
"Relative Tolerance", static_cast<Real>(2));
103 ParameterList &vlist = trlist.sublist(
"Inexact").sublist(
"Value");
104 scale_ = vlist.get(
"Tolerance Scaling", static_cast<Real>(1.e-1));
105 omega_ = vlist.get(
"Exponent", static_cast<Real>(0.9));
106 force_ = vlist.get(
"Forcing Sequence Initial Value", static_cast<Real>(1.0));
107 updateIter_ = vlist.get(
"Forcing Sequence Update Frequency", static_cast<int>(10));
108 forceFactor_ = vlist.get(
"Forcing Sequence Reduction Factor", static_cast<Real>(0.1));
110 verbosity_ = list.sublist(
"General").get(
"Output Level",0);
111 writeHeader_ = verbosity_ > 2;
113 useSecantPrecond_ = list.sublist(
"General").sublist(
"Secant").get(
"Use as Preconditioner",
false);
114 useSecantHessVec_ = list.sublist(
"General").sublist(
"Secant").get(
"Use as Hessian",
false);
119 model_ = makePtr<TrustRegionModel_U<Real>>(list,secant,mode);
120 if (secant == nullPtr) {
121 esec_ =
StringToESecant(list.sublist(
"General").sublist(
"Secant").get(
"Type",
"Limited-Memory BFGS"));
125 template<
typename Real>
130 std::ostream &outStream) {
133 if (proj_ == nullPtr) {
134 proj_ = makePtr<PolyhedralProjection<Real>>(makePtrFromRef(bnd));
141 Real ftol =
static_cast<Real
>(0.1)*ROL_OVERFLOW<Real>();
142 proj_->project(x,outStream); state_->nproj++;
143 state_->iterateVec->set(x);
145 state_->value = obj.
value(x,ftol);
148 computeGradient(x,*state_->gradientVec,*state_->stepVec,state_->searchSize,obj,
true,gtol_,state_->gnorm,outStream);
155 state_->snorm = ROL_INF<Real>();
157 if (normAlpha_) alpha_ /= state_->gradientVec->norm();
159 if ( state_->searchSize <= static_cast<Real>(0) )
160 state_->searchSize = state_->gradientVec->norm();
163 rcon_ = makePtr<ReducedLinearConstraint<Real>>(proj_->getLinearConstraint(),
166 ns_ = makePtr<NullSpaceOperator<Real>>(rcon_,x,
167 *proj_->getResidual());
171 template<
typename Real>
180 outTol = std::sqrt(ROL_EPSILON<Real>());
181 if ( useInexact_[0] ) {
182 if (!(iter%updateIter_) && (iter!=0)) force_ *= forceFactor_;
184 Real eta =
static_cast<Real
>(0.999)*std::min(eta1_,one-eta2_);
185 outTol = scale_*std::pow(eta*std::min(pRed,force_),one/omega_);
186 if (inTol > outTol) fold = obj.
value(xold,outTol);
190 Real fval = obj.
value(x,outTol);
194 template<
typename Real>
203 std::ostream &outStream)
const {
204 if ( useInexact_[1] ) {
206 Real gtol0 = scale0_*del;
207 if (accept) gtol = gtol0 + one;
208 else gtol0 = scale0_*std::min(gnorm,del);
209 while ( gtol > gtol0 ) {
213 gtol0 = scale0_*std::min(gnorm,del);
218 gtol = std::sqrt(ROL_EPSILON<Real>());
279 template<
typename Real>
284 std::ostream &outStream ) {
286 Real tol0 = std::sqrt(ROL_EPSILON<Real>());
287 Real inTol =
static_cast<Real
>(0.1)*ROL_OVERFLOW<Real>(), outTol(inTol);
288 Real gfnorm(0), gfnormf(0), tol(0), stol(0), snorm(0);
289 Real ftrial(0), pRed(0), rho(1), q(0), delta(0);
290 int flagCG(0), iterCG(0), maxit(0);
292 initialize(x,g,obj,bnd,outStream);
293 Ptr<Vector<Real>> s = x.
clone();
294 Ptr<Vector<Real>> gmod = g.
clone(), gfree = g.
clone();
295 Ptr<Vector<Real>> pwa1 = x.
clone(), pwa2 = x.
clone(), pwa3 = x.
clone();
296 Ptr<Vector<Real>> dwa1 = g.
clone(), dwa2 = g.
clone(), dwa3 = g.
clone();
298 Real rhoNM(0), sigmac(0), sigmar(0);
299 Real fr(state_->value), fc(state_->value), fmin(state_->value);
304 if (verbosity_ > 0) writeOutput(outStream,
true);
306 while (status_->check(*state_)) {
308 model_->setData(obj,*state_->iterateVec,*state_->gradientVec,gtol_);
312 snorm = dcauchy(*state_->stepVec,alpha_,q,*state_->iterateVec,
313 state_->gradientVec->dual(),state_->searchSize,
314 *model_,*dwa1,*dwa2,outStream);
315 x.
plus(*state_->stepVec);
316 state_->snorm = snorm;
317 delta = state_->searchSize - snorm;
322 gmod->plus(*state_->gradientVec);
325 pwa1->set(gfree->dual());
327 gfree->set(pwa1->dual());
329 applyFreePrecond(*pwa1,*gfree,x,*model_,bnd,tol0,*dwa1,*pwa2);
330 gfnorm = pwa1->norm();
333 gfnorm = gfree->norm();
335 SPiter_ = 0; SPflag_ = 0;
336 if (verbosity_ > 1) {
337 outStream <<
" Norm of free gradient components: " << gfnorm << std::endl;
338 outStream << std::endl;
343 for (
int i = 0; i < minit_; ++i) {
345 flagCG = 0; iterCG = 0;
347 tol = std::min(tol1_,tol2_*std::pow(gfnorm,spexp_));
349 if (gfnorm > zero && delta > zero) {
350 snorm = dtrpcg(*s,flagCG,iterCG,*gfree,x,
351 delta,*model_,bnd,tol,stol,maxit,
352 *pwa1,*dwa1,*pwa2,*dwa2,*pwa3,*dwa3,outStream);
355 if (verbosity_ > 1) {
356 outStream <<
" Computation of CG step" << std::endl;
357 outStream <<
" Current face (i): " << i << std::endl;
358 outStream <<
" CG step length: " << snorm << std::endl;
359 outStream <<
" Number of CG iterations: " << iterCG << std::endl;
360 outStream <<
" CG flag: " << flagCG << std::endl;
361 outStream <<
" Total number of iterations: " << SPiter_ << std::endl;
362 outStream << std::endl;
366 snorm = dprsrch(x,*s,q,gmod->
dual(),*model_,bnd,*pwa1,*dwa1,outStream);
370 state_->stepVec->plus(*s);
371 state_->snorm = state_->stepVec->norm();
372 delta = state_->searchSize - state_->snorm;
376 pwa1->set(gfree->dual());
377 bnd.pruneActive(*pwa1,x,zero);
378 gfree->set(pwa1->dual());
380 applyFreePrecond(*pwa1,*gfree,x,*model_,bnd,tol0,*dwa1,*pwa2);
381 gfnormf = pwa1->norm();
384 gfnormf = gfree->norm();
386 if (verbosity_ > 1) {
387 outStream <<
" Norm of free gradient components: " << gfnormf << std::endl;
388 outStream << std::endl;
393 if (gfnormf <= tol) {
397 else if (SPiter_ >= maxit_) {
401 else if (flagCG == 2) {
405 else if (delta <= zero) {
418 ftrial = computeValue(inTol,outTol,pRed,state_->value,state_->iter,x,*state_->iterateVec,obj);
423 TRUtils::analyzeRatio<Real>(rho,TRflag_,state_->value,ftrial,pRed,eps_,outStream,verbosity_>1);
425 TRUtils::analyzeRatio<Real>(rhoNM,TRflagNM,fr,ftrial,pRed+sigmar,eps_,outStream,verbosity_>1);
426 TRflag_ = (rho < rhoNM ? TRflagNM : TRflag_);
427 rho = (rho < rhoNM ? rhoNM : rho );
434 x.
set(*state_->iterateVec);
438 state_->searchSize = TRUtils::interpolateRadius<Real>(*state_->gradientVec,*state_->stepVec,
439 state_->snorm,pRed,state_->value,ftrial,state_->searchSize,gamma0_,gamma1_,eta2_,
440 outStream,verbosity_>1);
443 state_->searchSize = gamma1_*std::min(state_->snorm,state_->searchSize);
445 computeGradient(x,*state_->gradientVec,*pwa1,state_->searchSize,obj,
false,gtol_,state_->gnorm,outStream);
449 state_->value = ftrial;
452 sigmac += pRed; sigmar += pRed;
453 if (ftrial < fmin) { fmin = ftrial; fc = fmin; sigmac =
zero; L = 0; }
456 if (ftrial > fc) { fc = ftrial; sigmac =
zero; }
457 if (L == storageNM_) { fr = fc; sigmar = sigmac; }
461 if (rho >= eta2_) state_->searchSize = std::min(gamma2_*state_->searchSize, delMax_);
463 dwa1->set(*state_->gradientVec);
465 computeGradient(x,*state_->gradientVec,*pwa1,state_->searchSize,obj,
true,gtol_,state_->gnorm,outStream);
468 state_->iterateVec->set(x);
470 model_->update(x,*state_->stepVec,*dwa1,*state_->gradientVec,
471 state_->snorm,state_->iter);
475 if (verbosity_ > 0) writeOutput(outStream,writeHeader_);
480 template<
typename Real>
483 std::ostream &outStream)
const {
485 proj_->project(s,outStream); state_->nproj++;
486 s.
axpy(static_cast<Real>(-1),x);
490 template<
typename Real>
499 std::ostream &outStream) {
500 const Real half(0.5);
501 Real tol = std::sqrt(ROL_EPSILON<Real>());
503 Real gs(0), snorm(0);
505 snorm = dgpstep(s,g,x,-alpha,outStream);
510 model.
hessVec(dwa,s,x,tol); nhess_++;
512 q = half * s.
apply(dwa) + gs;
513 interp = (q > mu0_*gs);
519 while (search && cnt < redlim_) {
521 snorm = dgpstep(s,g,x,-alpha,outStream);
523 model.
hessVec(dwa,s,x,tol); nhess_++;
525 q = half * s.
apply(dwa) + gs;
526 search = (q > mu0_*gs);
530 if (cnt >= redlim_ && q > mu0_*gs) {
531 outStream <<
"Cauchy point: The interpolation limit was met without producing sufficient decrease." << std::endl;
532 outStream <<
" Lin-More trust-region algorithm may not converge!" << std::endl;
542 snorm = dgpstep(s,g,x,-alpha,outStream);
543 if (snorm <= del && cnt < explim_) {
544 model.
hessVec(dwa,s,x,tol); nhess_++;
546 q = half * s.
apply(dwa) + gs;
547 if (q <= mu0_*gs && std::abs(q-qs) > qtol_*std::abs(qs)) {
565 snorm = dgpstep(s,g,x,-alpha,outStream);
567 if (verbosity_ > 1) {
568 outStream <<
" Cauchy point" << std::endl;
569 outStream <<
" Step length (alpha): " << alpha << std::endl;
570 outStream <<
" Step length (alpha*g): " << snorm << std::endl;
571 outStream <<
" Model decrease (pRed): " << -q << std::endl;
573 outStream <<
" Number of extrapolation steps: " << cnt << std::endl;
579 template<
typename Real>
585 std::ostream &outStream) {
586 const Real
zero(0.0), half(0.5);
587 Real tol = std::sqrt(ROL_EPSILON<Real>());
588 Real beta(1), snorm(0), gs(0);
594 snorm = dgpstep(pwa,s,x,beta,outStream);
595 model.
hessVec(dwa,pwa,x,tol); nhess_++;
598 q = half * pwa.
apply(dwa) + gs;
599 if (q <= mu0_*std::min(gs,
zero) || nsteps > pslim_) {
608 if (verbosity_ > 1) {
609 outStream << std::endl;
610 outStream <<
" Projected search" << std::endl;
611 outStream <<
" Step length (beta): " << beta << std::endl;
612 outStream <<
" Step length (beta*s): " << snorm << std::endl;
613 outStream <<
" Model decrease (pRed): " << -q << std::endl;
614 outStream <<
" Number of steps: " << nsteps << std::endl;
619 template<
typename Real>
623 const Real del)
const {
626 Real rad = ptx*ptx + ptp*(dsq-xtx);
627 rad = std::sqrt(std::max(rad,zero));
630 sigma = (dsq-xtx)/(ptx+rad);
632 else if (rad > zero) {
633 sigma = (rad-ptx)/ptp;
641 template<
typename Real>
646 const Real tol,
const Real stol,
const int itermax,
649 std::ostream &outStream)
const {
654 Real tol0 = std::sqrt(ROL_EPSILON<Real>());
655 const Real
zero(0), one(1), two(2);
656 Real rho(0), kappa(0), beta(0), sigma(0), alpha(0);
657 Real rtr(0), tnorm(0), sMs(0), pMp(0), sMp(0);
664 applyFreePrecond(r,t,x,model,bnd,tol0,dwa,pwa);
669 pMp = (!hasEcon_ ? rho : p.
dot(p));
671 for (iter = 0; iter < itermax; ++iter) {
673 applyFreeHessian(q,p,x,model,bnd,tol0,pwa);
677 alpha = (kappa>
zero) ? rho/kappa :
zero;
678 sigma = dtrqsol(sMs,pMp,sMp,del);
680 if (kappa <= zero || alpha >= sigma) {
683 iflag = (kappa<=
zero) ? 2 : 3;
689 applyFreePrecond(r,t,x,model,bnd,tol0,dwa,pwa);
694 if (rtr <= stol*stol || tnorm <= tol) {
695 sMs = sMs + two*alpha*sMp + alpha*alpha*pMp;
707 sMs = sMs + two*alpha*sMp + alpha*alpha*pMp;
708 sMp = beta*(sMp + alpha*pMp);
709 pMp = (!hasEcon_ ? rho : p.
dot(p)) + beta*beta*pMp;
712 if (iter == itermax) {
718 return std::sqrt(sMs);
721 template<
typename Real>
732 model.
hessVec(hv,pwa,x,tol); nhess_++;
742 template<
typename Real>
765 rcon_->setX(makePtrFromRef(x));
768 ns_->apply(hv,pwa,tol);
772 template<
typename Real>
774 std::ios_base::fmtflags osFlags(os.flags());
775 if (verbosity_ > 1) {
776 os << std::string(114,
'-') << std::endl;
777 os <<
" Lin-More trust-region method status output definitions" << std::endl << std::endl;
778 os <<
" iter - Number of iterates (steps taken)" << std::endl;
779 os <<
" value - Objective function value" << std::endl;
780 os <<
" gnorm - Norm of the gradient" << std::endl;
781 os <<
" snorm - Norm of the step (update to optimization vector)" << std::endl;
782 os <<
" delta - Trust-Region radius" << std::endl;
783 os <<
" #fval - Number of times the objective function was evaluated" << std::endl;
784 os <<
" #grad - Number of times the gradient was computed" << std::endl;
785 os <<
" #hess - Number of times the Hessian was applied" << std::endl;
786 os <<
" #proj - Number of times the projection was applied" << std::endl;
788 os <<
" tr_flag - Trust-Region flag" << std::endl;
795 os <<
" iterCG - Number of Truncated CG iterations" << std::endl << std::endl;
796 os <<
" flagGC - Trust-Region Truncated CG flag" << std::endl;
802 os << std::string(114,
'-') << std::endl;
805 os << std::setw(6) << std::left <<
"iter";
806 os << std::setw(15) << std::left <<
"value";
807 os << std::setw(15) << std::left <<
"gnorm";
808 os << std::setw(15) << std::left <<
"snorm";
809 os << std::setw(15) << std::left <<
"delta";
810 os << std::setw(10) << std::left <<
"#fval";
811 os << std::setw(10) << std::left <<
"#grad";
812 os << std::setw(10) << std::left <<
"#hess";
813 os << std::setw(10) << std::left <<
"#proj";
814 os << std::setw(10) << std::left <<
"tr_flag";
816 os << std::setw(10) << std::left <<
"iterCG";
817 os << std::setw(10) << std::left <<
"flagCG";
823 template<
typename Real>
825 std::ios_base::fmtflags osFlags(os.flags());
826 os << std::endl <<
"Lin-More Trust-Region Method (Type B, Bound Constraints)" << std::endl;
830 template<
typename Real>
832 std::ios_base::fmtflags osFlags(os.flags());
833 os << std::scientific << std::setprecision(6);
834 if ( state_->iter == 0 ) writeName(os);
835 if ( write_header ) writeHeader(os);
836 if ( state_->iter == 0 ) {
838 os << std::setw(6) << std::left << state_->iter;
839 os << std::setw(15) << std::left << state_->value;
840 os << std::setw(15) << std::left << state_->gnorm;
841 os << std::setw(15) << std::left <<
"---";
842 os << std::setw(15) << std::left << state_->searchSize;
843 os << std::setw(10) << std::left << state_->nfval;
844 os << std::setw(10) << std::left << state_->ngrad;
845 os << std::setw(10) << std::left << nhess_;
846 os << std::setw(10) << std::left << state_->nproj;
847 os << std::setw(10) << std::left <<
"---";
849 os << std::setw(10) << std::left <<
"---";
850 os << std::setw(10) << std::left <<
"---";
856 os << std::setw(6) << std::left << state_->iter;
857 os << std::setw(15) << std::left << state_->value;
858 os << std::setw(15) << std::left << state_->gnorm;
859 os << std::setw(15) << std::left << state_->snorm;
860 os << std::setw(15) << std::left << state_->searchSize;
861 os << std::setw(10) << std::left << state_->nfval;
862 os << std::setw(10) << std::left << state_->ngrad;
863 os << std::setw(10) << std::left << nhess_;
864 os << std::setw(10) << std::left << state_->nproj;
865 os << std::setw(10) << std::left << TRflag_;
867 os << std::setw(10) << std::left << SPiter_;
868 os << std::setw(10) << std::left << SPflag_;
std::string ECGFlagToString(ECGFlag cgf)
Provides the interface to evaluate objective functions.
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Real computeValue(Real inTol, Real &outTol, Real pRed, Real &fold, int iter, const Vector< Real > &x, const Vector< Real > &xold, Objective< Real > &obj)
virtual void scale(const Real alpha)=0
Compute where .
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual Real apply(const Vector< Real > &x) const
Apply to a dual vector. This is equivalent to the call .
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)
virtual void plus(const Vector &x)=0
Compute , where .
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
LinMoreAlgorithm(ParameterList &list, const Ptr< Secant< Real >> &secant=nullPtr)
virtual void writeExitStatus(std::ostream &os) const
void applyFreePrecond(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, TrustRegionModel_U< Real > &model, BoundConstraint< Real > &bnd, Real &tol, Vector< Real > &dwa, Vector< Real > &pwa) const
ESecant StringToESecant(std::string s)
Real dtrpcg(Vector< Real > &w, int &iflag, int &iter, const Vector< Real > &g, const Vector< Real > &x, const Real del, TrustRegionModel_U< Real > &model, BoundConstraint< Real > &bnd, const Real tol, const Real stol, const int itermax, Vector< Real > &p, Vector< Real > &q, Vector< Real > &r, Vector< Real > &t, Vector< Real > &pwa, Vector< Real > &dwa, std::ostream &outStream=std::cout) const
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
virtual void update(const Vector< Real > &x, UpdateType type, int iter=-1)
Update objective function.
virtual Real dot(const Vector &x) const =0
Compute where .
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Real dgpstep(Vector< Real > &s, const Vector< Real > &w, const Vector< Real > &x, const Real alpha, std::ostream &outStream=std::cout) const
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol) override
Apply Hessian approximation to vector.
ETRFlag
Enumation of flags used by trust-region solvers.
virtual void gradient(Vector< Real > &g, const Vector< Real > &x, Real &tol)
Compute gradient.
Real dprsrch(Vector< Real > &x, Vector< Real > &s, Real &q, const Vector< Real > &g, TrustRegionModel_U< Real > &model, BoundConstraint< Real > &bnd, Vector< Real > &pwa, Vector< Real > &dwa, std::ostream &outStream=std::cout)
std::string NumberToString(T Number)
Provides the interface to evaluate trust-region model functions.
void writeHeader(std::ostream &os) const override
Print iterate header.
Real dtrqsol(const Real xtx, const Real ptp, const Real ptx, const Real del) const
void applyFreeHessian(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, TrustRegionModel_U< Real > &model, BoundConstraint< Real > &bnd, Real &tol, Vector< Real > &pwa) const
Provides interface for and implements limited-memory secant operators.
void initialize(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &bnd, std::ostream &outStream=std::cout)
Provides an interface to check status of optimization algorithms.
std::string ETRFlagToString(ETRFlag trf)
void writeName(std::ostream &os) const override
Print step name.
void pruneActive(Vector< Real > &v, const Vector< Real > &x, Real eps=Real(0))
Set variables to zero if they correspond to the -active set.
Provides the interface to apply upper and lower bound constraints.
Real optimalityCriterion(const Vector< Real > &x, const Vector< Real > &g, Vector< Real > &primal, std::ostream &outStream=std::cout) const
void initialize(const Vector< Real > &x, const Vector< Real > &g)
void writeOutput(std::ostream &os, const bool write_header=false) const override
Print iterate status.
virtual void set(const Vector &x)
Set where .
virtual Real norm() const =0
Returns where .
void computeGradient(const Vector< Real > &x, Vector< Real > &g, Vector< Real > &pwa, Real del, Objective< Real > &obj, bool accept, Real >ol, Real &gnorm, 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...
virtual void precond(Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &s, Real &tol) override
Apply preconditioner to vector.