Xyce  6.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
N_NLS_DampedNewton.h
Go to the documentation of this file.
1 //-----------------------------------------------------------------------------
2 // Copyright Notice
3 //
4 // Copyright 2002 Sandia Corporation. Under the terms
5 // of Contract DE-AC04-94AL85000 with Sandia Corporation, the U.S.
6 // Government retains certain rights in this software.
7 //
8 // Xyce(TM) Parallel Electrical Simulator
9 // Copyright (C) 2002-2014 Sandia Corporation
10 //
11 // This program is free software: you can redistribute it and/or modify
12 // it under the terms of the GNU General Public License as published by
13 // the Free Software Foundation, either version 3 of the License, or
14 // (at your option) any later version.
15 //
16 // This program is distributed in the hope that it will be useful,
17 // but WITHOUT ANY WARRANTY; without even the implied warranty of
18 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU General Public License
22 // along with this program. If not, see <http://www.gnu.org/licenses/>.
23 //-----------------------------------------------------------------------------
24 
25 //-------------------------------------------------------------------------
26 // Filename : $RCSfile: N_NLS_DampedNewton.h,v $
27 //
28 // Purpose : Specification file for the implemenation of the Newton
29 // trust-region related methods.
30 //
31 // Special Notes :
32 //
33 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences.
34 //
35 // Creation Date : 04/28/00
36 //
37 // Revision Information:
38 // ---------------------
39 //
40 // Revision Number: $Revision: 1.77 $
41 //
42 // Revision Date : $Date: 2014/08/07 23:08:54 $
43 //
44 // Current Owner : $Author: dgbaur $
45 //-------------------------------------------------------------------------
46 
47 #ifndef Xyce_N_NLS_DampedNewton_h
48 #define Xyce_N_NLS_DampedNewton_h
49 
50 #include <N_IO_fwd.h>
51 #include <N_NLS_fwd.h>
52 #include <N_NLS_NonLinearSolver.h>
53 #include <N_NLS_NLParams.h>
54 #include <N_NLS_ParamMgr.h>
55 #include <N_NLS_ReturnCodes.h>
56 
57 namespace Xyce {
58 namespace Nonlinear {
59 
60 //-----------------------------------------------------------------------------
61 // Class : DampedNewton
62 // Purpose : This class implements a damped Newton nonlinear solver.
63 // Special Notes :
64 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
65 // Creation Date : 5/03/00
66 //-----------------------------------------------------------------------------
68 {
69 public:
70 
71  DampedNewton(N_IO_CmdParse & cp);
72  ~DampedNewton();
73 
74  bool setOptions(const N_UTL_OptionBlock& OB);
75  bool setTranOptions(const N_UTL_OptionBlock& OB);
76  bool setHBOptions(const N_UTL_OptionBlock& OB);
77 
78  bool initializeAll();
79 
80  int solve (NonLinearSolver * nlsTmpPtr = NULL);
81 
82  int takeFirstSolveStep (NonLinearSolver * nlsTmpPtr = NULL);
83  int takeOneSolveStep ();
84 
85  int getNumIterations() const;
86 #ifdef Xyce_DEBUG_NONLINEAR
87  int getDebugLevel() const;
88  bool getScreenOutputFlag () const;
89  double getDebugMinTime() const;
90  double getDebugMaxTime() const;
91  int getDebugMinTimeStep() const;
92  int getDebugMaxTimeStep() const;
93  bool getMMFormat () const;
94 #endif
95 
96  int getContinuationStep() const;
97  int getParameterNumber() const;
98 
99  bool isFirstContinuationParam() const;
100  bool isFirstSolveComplete() const;
101 
102  void setAnalysisMode(AnalysisMode mode);
103 
104  double getMaxNormF() const
105  { return maxNormRHS_; };
106 
107  int getMaxNormFindex() const
108  { return maxNormRHSindex_; };
109 
110 protected:
111 
112 private:
113 
114  void updateWeights_();
115 
116  void printHeader_(std::ostream &os);
117  void printFooter_(std::ostream &os);
118  void printStepInfo_(std::ostream &os, int step);
119 
120  bool rhs_();
121  bool newton_();
122 
123  void direction_();
124  void updateX_();
125  bool computeStepLength_();
126 
127  bool divide_();
128  bool backtrack_();
129  bool fullNewton_();
130  bool spiceNewton_();
131  bool bankRose_();
132  bool descent_();
133  bool simpleBacktrack_();
134  bool simpleBt_(double gsinit, double finit);
135 
136  double constrain_();
137  void setForcing_(const double);
138  void evalModNewton_();
139  int converged_();
140 
142 
143 public:
144 
145 protected:
146 
147 private:
148  //! Convergence Rate
149  double resConvRate_;
150 
151  //! Weighted convergence Rate
153 
154  //! Constraint object pointer
156 
157  // Current RHS norms:
158  double normRHS_;
159  double maxNormRHS_;
161 
162  double normRHS_init_; // used in calculating the relative norm.
163 
164  // Current dx norm:
165  double normDX_;
166 
167  // Current weighted dx norm:
168  double wtNormDX_;
169 
170  // Current relative RHS norm:
171  double normRHS_rel_;
172 
173  // Current solution norm:
174  double normSoln_;
175 
176  // step length:
177  double stepLength_;
178 
179  // backtracking upper and lower bounds and percentage change (0..1):
180  double BTUpper_;
181  double BTLower_;
182 
183  // constraint factor:
185 
186  // current nonlinear solver step:
187  unsigned nlStep_;
188 
189  // current Newton step:
190  unsigned newtonStep_;
191 
192  // current modified Newton step:
193  unsigned modNewtonStep_;
194 
195  // current steepest descent step:
196  unsigned descentStep_;
197 
198  // current line-search step:
199  unsigned searchStep_;
200 
201  //! Pointer to direction vector
202  /*! \todo Is this vector internal or external?? Is it allocated or just a pointer?? */
203  N_LAS_Vector* searchDirectionPtr_;
204 
205  // Needed for some reason
206  N_LAS_Vector* tmpVectorPtr_;
207 
209 
210  double delta_;
211 
212  // Flag for controlling the loading of the Jacobian matrix.
214 
215  // Flag for determining if this solver has been called
216  // before, and deltaXTol. The first time the solver is
217  // called we may want to tighten up the convergence tolerance
218  // a la Petzold, et al. These were static variables
219  // down in the "solve" function of this class. ERK.
220  bool firstTime;
222 
224 
225  double etaOld;
226  double nlResNormOld;
227  double tmpConvRate;
228 
229  int count;
230 };
231 
232 //---------------------------------------------------------------------------
233 // Function : DampedNewton::getNumIterations
234 //
235 // Return Type : Integer (current number of nonlinear iterations)
236 //---------------------------------------------------------------------------
238 {
239  return nlStep_;
240 }
241 
242 #ifdef Xyce_DEBUG_NONLINEAR
243 //---------------------------------------------------------------------------
244 // Function : DampedNewton::getDebugLevel
245 //
246 // Return Type : int
247 //---------------------------------------------------------------------------
248 inline int DampedNewton::getDebugLevel() const
249 {
250  return nlParams.getDebugLevel();
251 }
252 
253 //---------------------------------------------------------------------------
254 // Function : DampedNewton::getScreenOutputFlag
255 //
256 // Return Type : int
257 //---------------------------------------------------------------------------
258 inline bool DampedNewton::getScreenOutputFlag () const
259 {
260  return nlParams.getScreenOutputFlag ();
261 }
262 
263 //---------------------------------------------------------------------------
264 // Function : DampedNewton::getDebugMinTime
265 //
266 // Return Type : double
267 //---------------------------------------------------------------------------
268 inline double DampedNewton::getDebugMinTime() const
269 {
270  return nlParams.getDebugMinTime();
271 }
272 
273 //---------------------------------------------------------------------------
274 // Function : DampedNewton::getDebugMaxTime
275 //
276 // Return Type : double
277 //---------------------------------------------------------------------------
278 inline double DampedNewton::getDebugMaxTime() const
279 {
280  return nlParams.getDebugMaxTime();
281 }
282 
283 //---------------------------------------------------------------------------
284 // Function : DampedNewton::getDebugMinTimeStep
285 //
286 // Return Type : int
287 //---------------------------------------------------------------------------
288 inline int DampedNewton::getDebugMinTimeStep() const
289 {
290  return nlParams.getDebugMinTimeStep();
291 }
292 
293 //---------------------------------------------------------------------------
294 // Function : DampedNewton::getDebugMaxTimeStep
295 //
296 // Return Type : int
297 //---------------------------------------------------------------------------
298 inline int DampedNewton::getDebugMaxTimeStep() const
299 {
300  return nlParams.getDebugMaxTimeStep();
301 }
302 
303 //---------------------------------------------------------------------------
304 // Function : DampedNewton::getMMFormat
305 //
306 // Return Type : bool
307 //---------------------------------------------------------------------------
308 inline bool DampedNewton::getMMFormat () const
309 {
310  return nlParams.getMMFormat ();
311 }
312 
313 #endif // debug nonlin
314 
315 //---------------------------------------------------------------------------
316 // Function : DampedNewton::setAnalysisMode
317 //
318 // Purpose : Specify the analysis mode to be used by the nonlinear
319 // solver in the next call to solve(). This *may* affect
320 // the parameters used by the solver.
321 //
322 // See Also : setOptions, setTranOptions
323 //
324 // - Input Arguments -
325 //
326 // mode : Mode to be used in the next nonlinear solve.
327 //---------------------------------------------------------------------------
329 {
331 }
332 
333 //---------------------------------------------------------------------------
334 // Function : DampedNewton::resetCountersAndTimers_
335 // Purpose : Reset the counters and timers
336 //---------------------------------------------------------------------------
338 {
340  resConvRate_ = 0.0;
341  wtUpdateConvRate_ = 1.0;
342 }
343 
344 //---------------------------------------------------------------------------
345 // Function : DampedNewton::resetCountersAndTimers_
346 // Purpose : Reset the counters and timers
347 //---------------------------------------------------------------------------
349 {
350  return true;
351 }
352 
353 } // namespace Nonlinear
354 } // namespace Xyce
355 
357 
358 #endif
359