44#ifndef ROL_TYPEB_NEWTONKRYLOVALGORITHM_DEF_HPP
45#define ROL_TYPEB_NEWTONKRYLOVALGORITHM_DEF_HPP
50template<
typename Real>
57 secantName_ = list.sublist(
"General").sublist(
"Secant").get(
"Type",
"Limited-Memory BFGS");
59 secant_ = SecantFactory<Real>(list);
62 secantName_ = list.sublist(
"General").sublist(
"Secant").get(
"User Defined Secant Name",
63 "Unspecified User Defined Secant Method");
66 krylovName_ = list.sublist(
"General").sublist(
"Krylov").get(
"Type",
"Conjugate Gradients");
68 krylov_ = KrylovFactory<Real>(list);
71template<
typename Real>
79 secantName_ = list.sublist(
"General").sublist(
"Secant").get(
"Type",
"Limited-Memory BFGS");
81 secant_ = SecantFactory<Real>(list);
84 secantName_ = list.sublist(
"General").sublist(
"Secant").get(
"User Defined Secant Name",
85 "Unspecified User Defined Secant Method");
88 krylovName_ = list.sublist(
"General").sublist(
"Krylov").get(
"User Defined Krylov Name",
89 "Unspecified User Defined Krylov Method");
92template<
typename Real>
99 ParameterList &lslist = list.sublist(
"Step").sublist(
"Line Search");
100 maxit_ = lslist.get(
"Function Evaluation Limit", 20);
101 alpha0_ = lslist.get(
"Initial Step Size", 1.0);
102 useralpha_ = lslist.get(
"User Defined Initial Step Size",
false);
103 usePrevAlpha_ = lslist.get(
"Use Previous Step Length as Initial Guess",
false);
104 c1_ = lslist.get(
"Sufficient Decrease Tolerance", 1e-4);
105 rhodec_ = lslist.sublist(
"Line-Search Method").get(
"Backtracking Rate", 0.5);
107 useSecantHessVec_ = list.sublist(
"General").sublist(
"Secant").get(
"Use as Hessian",
false);
108 useSecantPrecond_ = list.sublist(
"General").sublist(
"Secant").get(
"Use as Preconditioner",
false);
110 verbosity_ = list.sublist(
"General").get(
"Output Level", 0);
111 writeHeader_ = verbosity_ > 2;
114template<
typename Real>
119 std::ostream &outStream) {
121 if (proj_ == nullPtr) {
122 proj_ = makePtr<PolyhedralProjection<Real>>(makePtrFromRef(bnd));
127 Real ftol = std::sqrt(ROL_EPSILON<Real>());
128 proj_->project(x,outStream);
129 state_->iterateVec->set(x);
131 state_->value = obj.
value(x,ftol); state_->nfval++;
132 obj.
gradient(*state_->gradientVec,x,ftol); state_->ngrad++;
133 state_->stepVec->set(x);
134 state_->stepVec->axpy(-one,state_->gradientVec->dual());
135 proj_->project(*state_->stepVec,outStream);
136 state_->stepVec->axpy(-one,x);
137 state_->gnorm = state_->stepVec->norm();
138 state_->snorm = ROL_INF<Real>();
139 if (!useralpha_) alpha0_ = one;
140 state_->searchSize = alpha0_;
143template<
typename Real>
148 std::ostream &outStream ) {
151 initialize(x,g,obj,bnd,outStream);
153 Ptr<Vector<Real>> pwa = x.
clone(), pwa1 = x.
clone();
154 Real ftrial(0), gs(0), tol(std::sqrt(ROL_EPSILON<Real>()));
156 Ptr<LinearOperator<Real>> hessian, precond;
159 if (verbosity_ > 0) writeOutput(outStream,
true);
162 gp->set(state_->gradientVec->dual());
163 while (status_->check(*state_)) {
165 hessian = makePtr<HessianPNK>(makePtrFromRef(obj),makePtrFromRef(bnd),
166 state_->iterateVec,state_->gradientVec,state_->gnorm,
167 secant_,useSecantHessVec_,pwa);
168 precond = makePtr<PrecondPNK>(makePtrFromRef(obj),makePtrFromRef(bnd),
169 state_->iterateVec,state_->gradientVec,state_->gnorm,
170 secant_,useSecantPrecond_,pwa1);
172 krylov_->run(*s,*hessian,*state_->gradientVec,*precond,iterKrylov_,flagKrylov_);
173 if (flagKrylov_ == 2 && iterKrylov_ <= 1) {
177 if (!usePrevAlpha_) state_->searchSize = alpha0_;
178 x.
set(*state_->iterateVec);
179 x.
axpy(-state_->searchSize,*s);
180 proj_->project(x,outStream);
182 ftrial = obj.
value(x,tol); ls_nfval_ = 1;
183 state_->stepVec->set(x);
184 state_->stepVec->axpy(-one,*state_->iterateVec);
185 gs = state_->stepVec->dot(*gp);
186 if (verbosity_ > 1) {
187 outStream <<
" In TypeB::NewtonKrylovAlgorithm: Line Search" << std::endl;
188 outStream <<
" Step size: " << state_->searchSize << std::endl;
189 outStream <<
" Trial objective value: " << ftrial << std::endl;
190 outStream <<
" Computed reduction: " << state_->value-ftrial << std::endl;
191 outStream <<
" Dot product of gradient and step: " << gs << std::endl;
192 outStream <<
" Sufficient decrease bound: " << -gs*c1_ << std::endl;
193 outStream <<
" Number of function evaluations: " << ls_nfval_ << std::endl;
195 while ( state_->value - ftrial < -c1_*gs && ls_nfval_ < maxit_ ) {
196 state_->searchSize *= rhodec_;
197 x.
set(*state_->iterateVec);
198 x.
axpy(-state_->searchSize,*s);
199 proj_->project(x,outStream);
201 ftrial = obj.
value(x,tol); ls_nfval_++;
202 state_->stepVec->set(x);
203 state_->stepVec->axpy(-one,*state_->iterateVec);
204 gs = state_->stepVec->dot(*gp);
205 if (verbosity_ > 1) {
206 outStream << std::endl;
207 outStream <<
" Step size: " << state_->searchSize << std::endl;
208 outStream <<
" Trial objective value: " << ftrial << std::endl;
209 outStream <<
" Computed reduction: " << state_->value-ftrial << std::endl;
210 outStream <<
" Dot product of gradient and step: " << gs << std::endl;
211 outStream <<
" Sufficient decrease bound: " << -gs*c1_ << std::endl;
212 outStream <<
" Number of function evaluations: " << ls_nfval_ << std::endl;
215 state_->nfval += ls_nfval_;
218 state_->snorm = state_->stepVec->norm();
221 state_->iterateVec->set(x);
225 state_->value = ftrial;
227 gold->set(*state_->gradientVec);
228 obj.
gradient(*state_->gradientVec,x,tol); state_->ngrad++;
229 gp->set(state_->gradientVec->dual());
232 s->set(x); s->axpy(-one,*gp);
233 proj_->project(*s,outStream);
235 state_->gnorm = s->norm();
238 secant_->updateStorage(x,*state_->gradientVec,*gold,*state_->stepVec,state_->snorm,state_->iter);
241 if (verbosity_ > 0) writeOutput(outStream,writeHeader_);
246template<
typename Real>
248 std::ostream &outStream ) {
253 throw Exception::NotImplemented(
">>> TypeB::NewtonKrylovAlgorithm::run : This algorithm cannot solve problems with linear equality constraints!");
257template<
typename Real>
265 std::ostream &outStream ) {
266 throw Exception::NotImplemented(
">>> TypeB::NewtonKrylovAlgorithm::run : This algorithm cannot solve problems with linear equality constraints!");
269template<
typename Real>
271 std::stringstream hist;
272 if (verbosity_ > 1) {
273 hist << std::string(114,
'-') << std::endl;
274 if (!useSecantHessVec_) {
275 hist <<
"Line-Search Projected Newton";
278 hist <<
"Line-Search Projected Quasi-Newton with " << secantName_ <<
" Hessian approximation";
280 hist <<
" status output definitions" << std::endl << std::endl;
281 hist <<
" iter - Number of iterates (steps taken)" << std::endl;
282 hist <<
" value - Objective function value" << std::endl;
283 hist <<
" gnorm - Norm of the gradient" << std::endl;
284 hist <<
" snorm - Norm of the step (update to optimization vector)" << std::endl;
285 hist <<
" alpha - Line search step length" << std::endl;
286 hist <<
" #fval - Cumulative number of times the objective function was evaluated" << std::endl;
287 hist <<
" #grad - Cumulative number of times the gradient was computed" << std::endl;
288 hist <<
" ls_#fval - Number of the times the objective function was evaluated during the line search" << std::endl;
289 hist <<
" iterCG - Number of Krylov iterations" << std::endl << std::endl;
290 hist <<
" flagGC - Krylov flag" << std::endl;
295 hist << std::string(114,
'-') << std::endl;
299 hist << std::setw(6) << std::left <<
"iter";
300 hist << std::setw(15) << std::left <<
"value";
301 hist << std::setw(15) << std::left <<
"gnorm";
302 hist << std::setw(15) << std::left <<
"snorm";
303 hist << std::setw(15) << std::left <<
"alpha";
304 hist << std::setw(10) << std::left <<
"#fval";
305 hist << std::setw(10) << std::left <<
"#grad";
306 hist << std::setw(10) << std::left <<
"#ls_fval";
307 hist << std::setw(10) << std::left <<
"iterCG";
308 hist << std::setw(10) << std::left <<
"flagCG";
313template<
typename Real>
315 std::stringstream hist;
316 if (!useSecantHessVec_) {
317 hist << std::endl <<
"Line-Search Projected Newton (Type B, Bound Constraints)" << std::endl;
320 hist << std::endl <<
"Line-Search Projected Quasi-Newton with "
321 << secantName_ <<
" Hessian approximation" << std::endl;
326template<
typename Real>
328 std::stringstream hist;
329 hist << std::scientific << std::setprecision(6);
330 if ( state_->iter == 0 ) writeName(os);
331 if ( write_header ) writeHeader(os);
332 if ( state_->iter == 0 ) {
334 hist << std::setw(6) << std::left << state_->iter;
335 hist << std::setw(15) << std::left << state_->value;
336 hist << std::setw(15) << std::left << state_->gnorm;
337 hist << std::setw(15) << std::left <<
"---";
338 hist << std::setw(15) << std::left <<
"---";
339 hist << std::setw(10) << std::left << state_->nfval;
340 hist << std::setw(10) << std::left << state_->ngrad;
341 hist << std::setw(10) << std::left <<
"---";
342 hist << std::setw(10) << std::left <<
"---";
343 hist << std::setw(10) << std::left <<
"---";
348 hist << std::setw(6) << std::left << state_->iter;
349 hist << std::setw(15) << std::left << state_->value;
350 hist << std::setw(15) << std::left << state_->gnorm;
351 hist << std::setw(15) << std::left << state_->snorm;
352 hist << std::setw(15) << std::left << state_->searchSize;
353 hist << std::setw(10) << std::left << state_->nfval;
354 hist << std::setw(10) << std::left << state_->ngrad;
355 hist << std::setw(10) << std::left << ls_nfval_;
356 hist << std::setw(10) << std::left << iterKrylov_;
357 hist << std::setw(10) << std::left << flagKrylov_;
Provides the interface to apply upper and lower bound constraints.
Defines the general constraint operator interface.
Provides definitions for Krylov solvers.
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.
const Ptr< PolyhedralProjection< Real > > & getPolyhedralProjection()
Get the polyhedral projection object. This is a null pointer if no linear constraints and/or bounds a...
Provides interface for and implements limited-memory secant operators.
Provides an interface to check status of optimization algorithms.
virtual void run(Problem< Real > &problem, std::ostream &outStream=std::cout)
Run algorithm on bound constrained problems (Type-B). This is the primary Type-B interface.
void initialize(const Vector< Real > &x, const Vector< Real > &g)
virtual void writeExitStatus(std::ostream &os) const
Ptr< Secant< Real > > secant_
Secant object (used for quasi-Newton)
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...
ESecant esec_
Secant type.
void writeOutput(std::ostream &os, bool write_header=false) const override
Print iterate status.
void parseParameterList(ParameterList &list)
std::string secantName_
Secant name.
NewtonKrylovAlgorithm(ParameterList &list, const Ptr< Secant< Real > > &secant=nullPtr)
std::string krylovName_
Krylov name.
void writeName(std::ostream &os) const override
Print step name.
Ptr< Krylov< Real > > krylov_
Krylov solver object (used for inexact Newton)
void initialize(Vector< Real > &x, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &bnd, std::ostream &outStream=std::cout)
void writeHeader(std::ostream &os) const override
Print iterate header.
Defines the linear algebra or vector space interface.
virtual void set(const Vector &x)
Set where .
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 .
EKrylov StringToEKrylov(std::string s)
std::string NumberToString(T Number)
ESecant StringToESecant(std::string s)
std::string ECGFlagToString(ECGFlag cgf)