60 Teuchos::RCP<Vector<Real> >
xnew_;
61 Teuchos::RCP<LineSearch<Real> >
btls_;
69 tol_ = parlist.sublist(
"Step").sublist(
"Line Search").sublist(
"Line-Search Method").get(
"Bracketing Tolerance",1.e-8);
77 btls_->initialize(x,s,g,obj,con);
80 void run( Real &alpha, Real &fval,
int &ls_neval,
int &ls_ngrad,
101 Real val_tr = obj.
value(*xnew_,tol);
105 if ( val_tr < val_tl ) {
106 if (
LineSearch<Real>::status(
LINESEARCH_BRENTS,ls_neval,ls_ngrad,tr,fval,gs,val_tr,x,s,obj,con) ) {
116 if ( val_tl < val_tr ) {
126 const Real gr = (1.0+sqrt(5.0))/2.0;
127 const Real inv_gr2 = 1.0/(gr*gr);
128 const Real goldinv = 1.0/(1.0+gr);
130 const Real max_extrap_factor = 100.0;
139 while ( val_tr > val_tl && itbt < 8 ) {
143 tr = goldinv * (tc + gr*tl);
146 val_tr = obj.
value(*xnew_,tol);
151 if ( val_tr > val_tl ) {
162 tc = tl + (gr-1.0)*(tr-tl);
165 val_tc = obj.
value(*xnew_,tol);
169 if ( val_tl <= val_tr && val_tl <= val_tc ) {
173 else if ( val_tc <= val_tr && val_tc <= val_tl ) {
182 if (
LineSearch<Real>::status(
LINESEARCH_BISECTION,ls_neval,ls_ngrad,t,fval,gs,val_t,x,s,obj,con) ) {
189 while ( val_tr >= val_tc && itb < 8 ) {
190 q = ( val_tr-val_tl ) * (tr - tc);
191 r = ( val_tr-val_tc ) * (tr - tl);
193 tmp = (tmp > tiny ? tmp : -tmp);
194 tm = tr - (q*(tr-tc) - r*(tr-tl))/(2.0*tmp);
196 tlim = tl + max_extrap_factor * (tc-tr);
198 if ( (tr-tm)*(tm-tc) > 0.0 ) {
201 val_tm = obj.
value(*xnew_,tol);
203 if ( val_tm < val_tc ) {
209 else if ( val_tm > val_tr) {
213 tm = tc + gr*(tc-tr);
216 val_tm = obj.
value(*xnew_,tol);
219 else if ( (tc - tm)*(tm -tlim) > 0.0 ) {
222 val_tm = obj.
value(*xnew_,tol);
224 if ( val_tm < val_tc ) {
231 tm = tc + gr*(tc-tr);
234 val_tm = obj.
value(*xnew_,tol);
238 else if ( (tm-tlim)*(tlim-tc) >= 0.0 ) {
242 val_tm = obj.
value(*xnew_,tol);
246 tm = tc + gr*(tc-tr);
249 val_tm = obj.
value(*xnew_,tol);
261 if ( val_tl <= val_tr && val_tl <= val_tc ) {
265 else if ( val_tc <= val_tr && val_tc <= val_tl ) {
274 if (
LineSearch<Real>::status(
LINESEARCH_BISECTION,ls_neval,ls_ngrad,t,fval,gs,val_t,x,s,obj,con) ) {
296 fw = (val_tl<val_tc ? val_tl : val_tc);
297 if ( fw == val_tl ) {
309 a = (tr < tc ? tr : tc);
310 b = (tr > tc ? tr : tc);
312 while ( !
LineSearch<Real>::status(
LINESEARCH_BRENTS,ls_neval,ls_ngrad,t,fval,gs,ft,x,s,obj,con)
313 && std::abs(t - tm) > tol_*(b-a) ) {
319 Real tol1 = tol_*std::abs(t) + tiny;
320 Real tol2 = 2.0*tol1;
322 if ( std::abs(e) > tol1 || it < 2 ) {
333 if ( std::abs(p) > std::abs(0.5*q*etemp) || p <= q*(a-t) || p >= q*(b-t) ) {
334 d = inv_gr2*(e=(t>=tm ? a-t : b-t));
339 if ( u-a < tol2 || b-u < tol2 ) {
340 d = ( tm-t > 0.0 ? std::abs(tol1) : -std::abs(tol1) );
345 d = inv_gr2*(e = (t>=tm ? a-t : b-t) );
347 u = (std::abs(d)>=tol1 ? t+d : t+(d>=0.0 ? std::abs(tol1) : -std::abs(tol1)));
350 fu = obj.
value(*xnew_,tol);
374 if ( fu <= fw || w == t ) {
380 else if ( fu <= fv || v == t || v == w ) {
391 btls_->run(alpha,fval,ls_neval,ls_ngrad,gs,s,x,obj,con);
Provides the interface to evaluate objective functions.
void updateIterate(Vector< Real > &xnew, const Vector< Real > &x, const Vector< Real > &s, Real alpha, BoundConstraint< Real > &con)
virtual Real getInitialAlpha(int &ls_neval, int &ls_ngrad, const Real fval, const Real gs, const Vector< Real > &x, const Vector< Real > &s, Objective< Real > &obj, BoundConstraint< Real > &con)
virtual Real value(const Vector< Real > &x, Real &tol)=0
Compute value.
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
Defines the linear algebra or vector space interface.
void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &con)
Provides interface for and implements line searches.
Implements a Brent's method line search.
Implements a simple back tracking line search.
void run(Real &alpha, Real &fval, int &ls_neval, int &ls_ngrad, const Real &gs, const Vector< Real > &s, const Vector< Real > &x, Objective< Real > &obj, BoundConstraint< Real > &con)
Brents(Teuchos::ParameterList &parlist)
Provides the interface to apply upper and lower bound constraints.
Teuchos::RCP< Vector< Real > > xnew_
Teuchos::RCP< LineSearch< Real > > btls_
virtual void update(const Vector< Real > &x, bool flag=true, int iter=-1)
Update objective function.
virtual void initialize(const Vector< Real > &x, const Vector< Real > &s, const Vector< Real > &g, Objective< Real > &obj, BoundConstraint< Real > &con)
static const double ROL_EPSILON
Platform-dependent machine epsilon.