ROL
ROL_TrustRegionUtilities.hpp
Go to the documentation of this file.
1// @HEADER
2// ************************************************************************
3//
4// Rapid Optimization Library (ROL) Package
5// Copyright (2014) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact lead developers:
38// Drew Kouri (dpkouri@sandia.gov) and
39// Denis Ridzal (dridzal@sandia.gov)
40//
41// ************************************************************************
42// @HEADER
43
44#ifndef ROL_TRUSTREGIONUTILITIES_U_H
45#define ROL_TRUSTREGIONUTILITIES_U_H
46
48
49namespace ROL {
50namespace TRUtils {
51
71
72inline std::string ETRFlagToString(ETRFlag trf) {
73 std::string retString;
74 switch(trf) {
75 case SUCCESS:
76 retString = "Both actual and predicted reductions are positive (success)";
77 break;
78 case POSPREDNEG:
79 retString = "Actual reduction is positive and predicted reduction is negative (impossible)";
80 break;
81 case NPOSPREDPOS:
82 retString = "Actual reduction is nonpositive and predicted reduction is positive";
83 break;
84 case NPOSPREDNEG:
85 retString = "Actual reduction is nonpositive and predicted reduction is negative (impossible)";
86 break;
87 case TRNAN:
88 retString = "Actual and/or predicted reduction is a NaN";
89 break;
90 case QMINSUFDEC:
91 retString = "Subproblem solution did not produce sufficient decrease";
92 break;
93 default:
94 retString = "INVALID ETRFlag";
95 }
96 return retString;
97}
98
99template<typename Real>
100inline Real initialRadius(int &nfval,
101 const Vector<Real> &x,
102 const Vector<Real> &g,
103 Vector<Real> &Bg,
104 const Real fx,
105 const Real gnorm,
106 Objective<Real> &obj,
108 const Real delMax,
109 std::ostream &outStream,
110 const bool print = false) {
111 const Real zero(0), half(0.5), one(1), two(2), three(3), six(6);
112 const Real eps(ROL_EPSILON<Real>());
113 Real del(ROL_INF<Real>());
114 Ptr<Vector<Real>> xcp = x.clone();
115 model.setData(obj,x,g);
116 Real htol = std::sqrt(eps);
117 model.hessVec(Bg,g.dual(),x,htol);
118 Real gBg = Bg.dot(g);
119 Real alpha = one;
120 if ( gBg > eps ) {
121 alpha = gnorm*gnorm/gBg;
122 }
123 // Evaluate the objective function at the Cauchy point
124 xcp->set(g.dual());
125 xcp->scale(-alpha);
126 //Real gs = xcp->dot(g.dual());
127 Real gs = xcp->apply(g);
128 xcp->plus(x);
129 obj.update(*xcp,UpdateType::Temp);
130 Real ftol = static_cast<Real>(0.1)*ROL_OVERFLOW<Real>();
131 Real fnew = obj.value(*xcp,ftol); // MUST DO SOMETHING HERE WITH FTOL
132 nfval++;
133 // Perform cubic interpolation to determine initial trust region radius
134 Real a = fnew - fx - gs - half*alpha*alpha*gBg;
135 if ( std::abs(a) < eps ) {
136 // a = 0 implies the objective is quadratic in the negative gradient direction
137 del = std::min(alpha*gnorm,delMax);
138 }
139 else {
140 Real b = half*alpha*alpha*gBg;
141 Real c = gs;
142 if ( b*b-three*a*c > eps ) {
143 // There is at least one critical point
144 Real t1 = (-b-std::sqrt(b*b-three*a*c))/(three*a);
145 Real t2 = (-b+std::sqrt(b*b-three*a*c))/(three*a);
146 if ( six*a*t1 + two*b > zero ) {
147 // t1 is the minimizer
148 del = std::min(t1*alpha*gnorm,delMax);
149 }
150 else {
151 // t2 is the minimizer
152 del = std::min(t2*alpha*gnorm,delMax);
153 }
154 }
155 else {
156 del = std::min(alpha*gnorm,delMax);
157 }
158 }
159 if (del <= eps*gnorm) {
160 del = one;
161 }
163 if ( print ) {
164 outStream << " In TrustRegionUtilities::initialRadius" << std::endl;
165 outStream << " Initial radius: " << del << std::endl;
166 }
167 return del;
168}
169
170template<typename Real>
171inline void analyzeRatio(Real &rho,
172 ETRFlag &flag,
173 const Real fold,
174 const Real ftrial,
175 const Real pRed,
176 const Real epsi,
177 std::ostream &outStream = std::cout,
178 const bool print = false) {
179 const Real zero(0), one(1);
180 Real eps = epsi*std::max(one,fold);
181 Real aRed = fold - ftrial;
182 Real aRed_safe = aRed + eps, pRed_safe = pRed + eps;
183 if (((std::abs(aRed_safe) < epsi) && (std::abs(pRed_safe) < epsi)) || aRed == pRed) {
184 rho = one;
185 flag = SUCCESS;
186 }
187 else if ( std::isnan(aRed_safe) || std::isnan(pRed_safe) ) {
188 rho = -one;
189 flag = TRNAN;
190 }
191 else {
192 rho = aRed_safe/pRed_safe;
193 if (pRed_safe < zero && aRed_safe > zero) {
194 flag = POSPREDNEG;
195 }
196 else if (aRed_safe <= zero && pRed_safe > zero) {
197 flag = NPOSPREDPOS;
198 }
199 else if (aRed_safe <= zero && pRed_safe < zero) {
200 flag = NPOSPREDNEG;
201 }
202 else {
203 flag = SUCCESS;
204 }
205 }
206 if ( print ) {
207 outStream << " In TrustRegionUtilities::analyzeRatio" << std::endl;
208 outStream << " Current objective function value: " << fold << std::endl;
209 outStream << " New objective function value: " << ftrial << std::endl;
210 outStream << " Actual reduction: " << aRed << std::endl;
211 outStream << " Predicted reduction: " << pRed << std::endl;
212 outStream << " Safeguard: " << epsi << std::endl;
213 outStream << " Actual reduction with safeguard: " << aRed_safe << std::endl;
214 outStream << " Predicted reduction with safeguard: " << pRed_safe << std::endl;
215 outStream << " Ratio of actual and predicted reduction: " << rho << std::endl;
216 outStream << " Trust-region flag: " << flag << std::endl;
217 outStream << std::endl;
218 }
219}
220
221template<typename Real>
222inline Real interpolateRadius(const Vector<Real> &g,
223 const Vector<Real> &s,
224 const Real snorm,
225 const Real pRed,
226 const Real fold,
227 const Real ftrial,
228 const Real del,
229 const Real gamma0,
230 const Real gamma1,
231 const Real eta2,
232 std::ostream &outStream = std::cout,
233 const bool print = false) {
234 const Real one(1);
235 //Real gs = g.dot(s.dual());
236 Real gs = g.apply(s);
237 Real modelVal = fold - pRed;
238 Real theta = (one-eta2)*gs/((one-eta2)*(fold+gs)+eta2*modelVal-ftrial);
239 if ( print ) {
240 outStream << " In TrustRegionUtilities::interpolateRadius" << std::endl;
241 outStream << " Interpolation model value: " << modelVal << std::endl;
242 outStream << " Interpolation step length: " << theta << std::endl;
243 outStream << std::endl;
244 }
245 return std::min(gamma1*std::min(snorm,del),std::max(gamma0,theta)*del);
246}
247
248} // namespace TrustRegion
249} // namespace ROL
250
251
252#endif
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0 zero)()
Provides the interface to evaluate objective functions.
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.
Provides the interface to evaluate trust-region model functions.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &s, Real &tol) override
Apply Hessian approximation to vector.
virtual void setData(Objective< Real > &obj, const Vector< Real > &x, const Vector< Real > &g)
Defines the linear algebra or vector space interface.
virtual Real apply(const Vector< Real > &x) const
Apply to a dual vector. This is equivalent to the call .
virtual const Vector & dual() const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis,...
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
virtual Real dot(const Vector &x) const =0
Compute where .
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)
std::string ETRFlagToString(ETRFlag trf)
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)
Real initialRadius(int &nfval, const Vector< Real > &x, const Vector< Real > &g, Vector< Real > &Bg, const Real fx, const Real gnorm, Objective< Real > &obj, TrustRegionModel_U< Real > &model, const Real delMax, std::ostream &outStream, const bool print=false)