ROL
ROL_NonlinearCG.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_NONLINEARCG_H
45#define ROL_NONLINEARCG_H
46
72#include "ROL_Types.hpp"
73
74namespace ROL {
75
76template<class Real>
77struct NonlinearCGState {
78 std::vector<ROL::Ptr<Vector<Real> > > grad; // Gradient Storage
79 std::vector<ROL::Ptr<Vector<Real> > > pstep; // Step Storage
80 int iter; // Nonlinear-CG Iteration Counter
81 int restart; // Reinitialize every 'restart' iterations
82 ENonlinearCG nlcg_type; // Nonlinear-CG Type
83};
84
85template<class Real>
86class NonlinearCG {
87private:
88
89 ROL::Ptr<NonlinearCGState<Real> > state_; // Nonlinear-CG State
90
91 ROL::Ptr<Vector<Real> > y_;
92 ROL::Ptr<Vector<Real> > yd_;
93
94public:
95
96 virtual ~NonlinearCG() {}
97
98 // Constructor
99 NonlinearCG(ENonlinearCG type, int restart = 100) {
100 state_ = ROL::makePtr<NonlinearCGState<Real>>();
101 state_->iter = 0;
102 state_->grad.resize(1);
103 state_->pstep.resize(1);
104 ROL_TEST_FOR_EXCEPTION(!(isValidNonlinearCG(type)),
105 std::invalid_argument,
106 ">>> ERROR (ROL_NonlinearCG.hpp): Invalid nonlinear CG type in constructor!");
107 state_->nlcg_type = type;
108 ROL_TEST_FOR_EXCEPTION((restart < 1),
109 std::invalid_argument,
110 ">>> ERROR (ROL_NonlinearCG.hpp): Non-positive restart integer in constructor!");
111 state_->restart = restart;
112 }
113
114 ROL::Ptr<NonlinearCGState<Real> >& get_state() { return this->state_; }
115
116 // Run one step of nonlinear CG.
117 virtual void run( Vector<Real> &s , const Vector<Real> &g, const Vector<Real> &x, Objective<Real> &obj ) {
118 Real one(1);
119 // Initialize vector storage
120 if ( state_->iter == 0 ) {
121 if ( state_->nlcg_type != NONLINEARCG_FLETCHER_REEVES &&
122 state_->nlcg_type != NONLINEARCG_FLETCHER_CONJDESC ) {
123 y_ = g.clone();
124 }
125 if ( state_->nlcg_type == NONLINEARCG_HAGER_ZHANG ||
126 state_->nlcg_type == NONLINEARCG_OREN_LUENBERGER ) {
127 yd_ = g.clone();
128 }
129 }
130
131 s.set(g.dual());
132
133 if ((state_->iter % state_->restart) != 0) {
134 Real beta(0), zero(0);
135 switch(state_->nlcg_type) {
136
138 y_->set(g);
139 y_->axpy(-one, *(state_->grad[0]));
140 beta = - g.dot(*y_) / (state_->pstep[0]->dot(y_->dual()));
141 beta = std::max(beta, zero);
142 break;
143 }
144
146 beta = g.dot(g) / (state_->grad[0])->dot(*(state_->grad[0]));
147 break;
148 }
149
150 case NONLINEARCG_DANIEL: {
151 Real htol(0);
152 obj.hessVec( *y_, *(state_->pstep[0]), x, htol );
153 beta = - g.dot(*y_) / (state_->pstep[0])->dot(y_->dual());
154 beta = std::max(beta, zero);
155 break;
156 }
157
159 y_->set(g);
160 y_->axpy(-one, *(state_->grad[0]));
161 beta = g.dot(*y_) / (state_->grad[0])->dot(*(state_->grad[0]));
162 beta = std::max(beta, zero);
163 break;
164 }
165
167 beta = g.dot(g) / (state_->pstep[0])->dot((state_->grad[0])->dual());
168 break;
169 }
170
172 y_->set(g);
173 y_->axpy(-one, *(state_->grad[0]));
174 beta = g.dot(*y_) / (state_->pstep[0])->dot((state_->grad[0])->dual());
175 //beta = std::max(beta, 0.0); // Is this needed? May need research.
176 break;
177 }
178
180 y_->set(g);
181 y_->axpy(-one, *(state_->grad[0]));
182 beta = - g.dot(g) / (state_->pstep[0])->dot(y_->dual());
183 break;
184 }
185
187 Real eta_0(1e-2), two(2);
188 y_->set(g);
189 y_->axpy(-one, *(state_->grad[0]));
190 yd_->set(*y_);
191 Real mult = two * ( y_->dot(*y_) / (state_->pstep[0])->dot(y_->dual()) );
192 yd_->axpy(-mult, (state_->pstep[0])->dual());
193 beta = - yd_->dot(g) / (state_->pstep[0])->dot(y_->dual());
194 Real eta = -one / ((state_->pstep[0])->norm()*std::min(eta_0,(state_->grad[0])->norm()));
195 beta = std::max(beta, eta);
196 break;
197 }
198
200 Real eta_0(1e-2);
201 y_->set(g);
202 y_->axpy(-one, *(state_->grad[0]));
203 yd_->set(*y_);
204 Real mult = ( y_->dot(*y_) / (state_->pstep[0])->dot(y_->dual()) );
205 yd_->axpy(-mult, (state_->pstep[0])->dual());
206 beta = - yd_->dot(g) / (state_->pstep[0])->dot(y_->dual());
207 Real eta = -one / ((state_->pstep[0])->norm()*std::min(eta_0,(state_->grad[0])->norm()));
208 beta = std::max(beta, eta);
209 break;
210 }
211
212 default:
213 ROL_TEST_FOR_EXCEPTION(!(isValidNonlinearCG(state_->nlcg_type)),
214 std::invalid_argument,
215 ">>> ERROR (ROL_NonlinearCG.hpp): Invalid nonlinear CG type in the 'run' method!");
216 }
217
218 s.axpy(beta, *(state_->pstep[0]));
219 }
220
221 // Update storage.
222 if (state_->iter == 0) {
223 (state_->grad[0]) = g.clone();
224 (state_->pstep[0]) = s.clone();
225 }
226 (state_->grad[0])->set(g);
227 (state_->pstep[0])->set(s);
228 state_->iter++;
229 }
230
231};
232
233}
234
235#endif
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Contains definitions of custom data types in ROL.
ENonlinearCG
Definition: ROL_Types.hpp:566
@ NONLINEARCG_POLAK_RIBIERE
Definition: ROL_Types.hpp:570
@ NONLINEARCG_HESTENES_STIEFEL
Definition: ROL_Types.hpp:567
@ NONLINEARCG_DAI_YUAN
Definition: ROL_Types.hpp:573
@ NONLINEARCG_OREN_LUENBERGER
Definition: ROL_Types.hpp:575
@ NONLINEARCG_DANIEL
Definition: ROL_Types.hpp:569
@ NONLINEARCG_FLETCHER_REEVES
Definition: ROL_Types.hpp:568
@ NONLINEARCG_HAGER_ZHANG
Definition: ROL_Types.hpp:574
@ NONLINEARCG_LIU_STOREY
Definition: ROL_Types.hpp:572
@ NONLINEARCG_FLETCHER_CONJDESC
Definition: ROL_Types.hpp:571
int isValidNonlinearCG(ENonlinearCG s)
Verifies validity of a NonlinearCG enum.
Definition: ROL_Types.hpp:604