Xyce  6.1
N_NLS_NLParams.C
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-2015 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_NLParams.C,v $
27 //
28 // Purpose :
29 //
30 // Special Notes :
31 //
32 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
33 //
34 // Creation Date : 07/13/01
35 //
36 // Revision Information:
37 // ---------------------
38 //
39 // Revision Number: $Revision: 1.83.2.2 $
40 //
41 // Revision Date : $Date: 2015/04/02 18:20:17 $
42 //
43 // Current Owner : $Author $
44 //-------------------------------------------------------------------------
45 
46 #include <Xyce_config.h>
47 
48 #include <N_ERH_ErrorMgr.h>
49 #include <N_IO_CmdParse.h>
50 #include <N_NLS_DampedNewton.h>
51 #include <N_NLS_NLParams.h>
52 #include <N_NLS_NonLinearSolver.h>
53 #include <N_NLS_TwoLevelNewton.h>
54 #include <N_UTL_FeatureTest.h>
55 #include <N_UTL_OptionBlock.h>
56 #include <N_UTL_Param.h>
57 
58 namespace Xyce {
59 namespace Nonlinear {
60 
61 //-----------------------------------------------------------------------------
62 // Function : NLParams::NLParams
63 // Purpose : constructor
64 // Special Notes :
65 // Scope : public
66 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
67 // Creation Date : 5/09/00
68 //-----------------------------------------------------------------------------
69 NLParams::NLParams(AnalysisMode mode, const IO::CmdParse & cp)
70  : commandLine_(&cp),
71  modeToggled_(true),
72  printParamsFlag_(true),
73  debugLevel_(1),
74  debugMinTimeStep_(0),
75  debugMaxTimeStep_(Util::MachineDependentParams::IntMax()),
76  debugMinTime_(0.0),
77  debugMaxTime_(Util::MachineDependentParams::DoubleMax()),
78  screenOutputFlag_(false),
79  maskingFlag_(false),
80  analysisMode_(mode)
81 {
82  // Set default update norm tolerance
84 
85  // Set default small update norm tolerance
87 
88  // Set default residual norm tolerance
89  resetRHSTol();
90 
91  // Set default absolute tolerance value for use in weighted norm
92  resetAbsTol();
93 
94  // Set default relative tolerance value for use in weighted norm
95  resetRelTol();
96 
98 
100  resetDirection();
101  resetNLStrategy();
106  resetNormLevel();
107  resetLinearOpt();
112 
113  // Set the default parameters for transient, if the specified mode
114  // is TRANSIENT.
115  // The defaults set in the NLParams constructor are for DC_OP,
116  // so the transient ones should be reset here.
117  if (mode==TRANSIENT)
118  {
121  setMaxSearchStep(2);
122  setMaxNewtonStep(20);
123  setDeltaXTol(0.33);
124  setForcingFlag(false);
125  setAbsTol(1.0e-06);
126  setRelTol(1.0e-02);
127  setRHSTol(1.0e-02);
128  }
129 
130  if (mode == HB_MODE)
131  {
134  setAbsTol(1.0e-9);
135  setRHSTol(1.0e-4);
136  }
137 
138 
140 }
141 
142 //-----------------------------------------------------------------------------
143 // Function : NLParams::NLParams
144 // Purpose : copy constructor
145 // Special Notes :
146 // Scope : public
147 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
148 // Creation Date : 5/09/00
149 //-----------------------------------------------------------------------------
151  : commandLine_(right.commandLine_),
152  absTol_(right.absTol_),
153  relTol_(right.relTol_),
154  INForcingFlag_(right.INForcingFlag_),
155  eta_(right.eta_),
156  normLevel_(right.normLevel_),
157  linearOptimization_(right.linearOptimization_),
158  modeToggled_(right.modeToggled_),
159  printParamsFlag_(right.printParamsFlag_),
160  analysisMode_(right.analysisMode_),
161  maxNewtonStep_(right.maxNewtonStep_),
162  maxSearchStep_(right.maxSearchStep_),
163  nlStrategy_(right.nlStrategy_),
164  searchMethod_(right.searchMethod_),
165  direction_(right.direction_),
166  deltaXTol_(right.deltaXTol_),
167  RHSTol_(right.RHSTol_),
168  constraintBT_(right.constraintBT_),
169  globalBTMax_(right.globalBTMax_),
170  globalBTMin_(right.globalBTMin_),
171  globalBTChange_(right.globalBTChange_),
172  debugLevel_(right.debugLevel_),
173  debugMinTimeStep_(right.debugMinTimeStep_),
174  debugMaxTimeStep_(right.debugMaxTimeStep_),
175  debugMinTime_(right.debugMinTime_),
176  debugMaxTime_(right.debugMaxTime_)
177 {
178 }
179 
180 //-----------------------------------------------------------------------------
181 // Function : NLParams::~NLParams
182 // Purpose : destructor
183 // Special Notes :
184 // Scope : public
185 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
186 // Creation Date : 5/09/00
187 //-----------------------------------------------------------------------------
189 {
190 // delete lasSolverPtr;
191 }
192 
193 //-----------------------------------------------------------------------------
194 // Function : NLParams::setOptions
195 // Purpose : This function takes an .options statement option block,
196 // and uses it to set the various parameters of
197 // the NLParams class.
198 //
199 // Special Notes : This was originally in DampedNewton, but it makes
200 // more sense to have it here.
201 //
202 // Scope : public
203 // Creator : Eric R. Keiter, SNL, Computational Sciences
204 // Creation Date : 10/20/02
205 //-----------------------------------------------------------------------------
206 
207 bool NLParams::setOptions (const Util::OptionBlock & OB)
208 {
209  for (Util::ParamList::const_iterator it_tpL = OB.begin();
210  it_tpL != OB.end(); ++ it_tpL)
211  {
212  if (it_tpL->uTag() == "ABSTOL")
213  {
214  setAbsTol(it_tpL->getImmutableValue<double>());
215  }
216  else if (it_tpL->uTag() == "RELTOL")
217  {
218  setRelTol(it_tpL->getImmutableValue<double>());
219  }
220  else if (it_tpL->uTag() == "DELTAXTOL")
221  {
222  setDeltaXTol(it_tpL->getImmutableValue<double>());
223  }
224  else if (it_tpL->uTag() == "SMALLUPDATETOL")
225  {
226  setSmallUpdateTol(it_tpL->getImmutableValue<double>());
227  }
228  else if (it_tpL->uTag() == "ENFORCEDEVICECONV")
229  {
230  setEnforceDeviceConvFlag(it_tpL->getImmutableValue<double>());
231  }
232  else if (it_tpL->uTag() == "RHSTOL")
233  {
234  setRHSTol(it_tpL->getImmutableValue<double>());
235  }
236  else if (it_tpL->uTag() == "MAXSTEP")
237  {
238  setMaxNewtonStep(it_tpL->getImmutableValue<int>());
239  }
240  else if (it_tpL->uTag() == "LINOPT")
241  {
242  setLinearOpt(it_tpL->getImmutableValue<int>());
243  }
244  else if (it_tpL->uTag() == "CONSTRAINTBT")
245  {
246  setConstraintBT(it_tpL->getImmutableValue<int>());
247  }
248  else if (it_tpL->uTag() == "CONSTRAINTMAX")
249  {
250  setGlobalBTMax(it_tpL->getImmutableValue<double>());
251  }
252  else if (it_tpL->uTag() == "CONSTRAINTMIN")
253  {
254  setGlobalBTMin(it_tpL->getImmutableValue<double>());
255  }
256  else if (it_tpL->uTag() == "CONSTRAINTCHANGE")
257  {
258  setGlobalBTChange(it_tpL->getImmutableValue<double>());
259  }
260  else if (it_tpL->uTag() == "NLSTRATEGY")
261  {
262  setNLStrategy(it_tpL->getImmutableValue<int>());
263  }
264  else if (it_tpL->uTag() == "SEARCHMETHOD")
265  {
266  setSearchMethod(it_tpL->getImmutableValue<int>());
267  }
268  else if (it_tpL->uTag() == "MAXSEARCHSTEP")
269  {
270  setMaxSearchStep(it_tpL->getImmutableValue<int>());
271  }
272  else if (it_tpL->uTag() == "IN_FORCING")
273  {
274  setForcingFlag(it_tpL->getImmutableValue<int>());
275  }
276  else if (it_tpL->uTag() == "NORMLVL")
277  {
278  setNormLevel(it_tpL->getImmutableValue<int>());
279  }
280  else if (it_tpL->uTag() == "NOX")
281  {
282  // do nothing.
283  }
284  else if (it_tpL->uTag() == "MATRIXMARKET")
285  {
286  if (DEBUG_NONLINEAR)
287  {
288  setMMFormat (static_cast<bool>(it_tpL->getImmutableValue<double>()));
289  }
290  }
291  else if (it_tpL->uTag() == "DEBUGLEVEL")
292  {
293  if (DEBUG_NONLINEAR)
294  {
295  setDebugLevel(it_tpL->getImmutableValue<int>());
296  }
297  }
298  else if (it_tpL->uTag() == "DEBUGMINTIMESTEP")
299  {
300  if (DEBUG_NONLINEAR)
301  {
302  setDebugMinTimeStep(it_tpL->getImmutableValue<int>());
303  }
304  }
305  else if (it_tpL->uTag() == "DEBUGMAXTIMESTEP")
306  {
307  if (DEBUG_NONLINEAR)
308  {
309  setDebugMaxTimeStep(it_tpL->getImmutableValue<int>());
310  }
311  }
312  else if (it_tpL->uTag() == "DEBUGMINTIME")
313  {
314  if (DEBUG_NONLINEAR)
315  {
316  setDebugMinTime(it_tpL->getImmutableValue<double>());
317  }
318  }
319  else if (it_tpL->uTag() == "DEBUGMAXTIME")
320  {
321  if (DEBUG_NONLINEAR)
322  {
323  setDebugMaxTime(it_tpL->getImmutableValue<double>());
324  }
325  }
326  else if (it_tpL->uTag() == "SCREENOUTPUT")
327  {
328  if (DEBUG_NONLINEAR)
329  {
330  setScreenOutputFlag (static_cast<bool>(it_tpL->getImmutableValue<double>()));
331  }
332  }
333  else if (NLS_MASKED_WRMS_NORMS && it_tpL->uTag() == "USEMASKING")
334  {
335  setMaskingFlag(static_cast<bool>(it_tpL->getImmutableValue<double>()));
336  }
337  else
338  {
339  std::string tmp = it_tpL->uTag() +
340  " is not a recognized nonlinear solver option.\n";
341  N_ERH_ErrorMgr::report (N_ERH_ErrorMgr::USR_FATAL_0, tmp);
342  }
343  }
344 
346 
347  return true;
348 }
349 
350 //-----------------------------------------------------------------------------
351 // Function : NLParams::printParams
352 // Purpose : Print out the nonlinear solver parameter values.
353 // Special Notes :
354 // Scope : public
355 // Creator : Scott A. Hutchinson, SNL, Computational Sciences
356 // Creation Date : 01/16/01
357 //-----------------------------------------------------------------------------
358 void NLParams::printParams(std::ostream &os)
359 {
360  os << "\n" << std::endl
361  << Xyce::section_divider << std::endl;
362  os << "\n***** Nonlinear solver options:\n" << std::endl
363  << "\tabsTol:\t\t\t" << getAbsTol()
364  << "\trelTol:\t\t\t" << getRelTol()
365  << "\tdeltaXTol (weighted):\t" << getDeltaXTol()
366  << "\tRHSTol:\t\t\t" << getRHSTol()
367  << "\tSmall Update Tol:\t" << getSmallUpdateTol()
368  << "\tmax NL Steps:\t\t" << getMaxNewtonStep();
369 
370  if (analysisMode_ == DC_OP)
371  os << "\tAnalysis Mode:\t\t" << analysisMode_ << "\t(DC Op)" << std::endl;
372  else if (analysisMode_ == DC_SWEEP)
373  os << "\tAnalysis Mode:\t\t" << analysisMode_ << "\t(DC Sweep)" << std::endl;
374  else if (analysisMode_ == TRANSIENT)
375  os << "\tAnalysis Mode:\t\t" << analysisMode_ << "\t(Transient)" << std::endl;
376 
377  NLStrategy strategy = getNLStrategy();
378  if (strategy == NEWTON)
379  os << "\tNL Strategy:\t\t" << strategy << "\t(None => Newton)" << std::endl;
380  else if (strategy == GRADIENT)
381  os << "\tNL Strategy:\t\t" << strategy << "\t(Gradient)" << std::endl;
382  else if (strategy == NEWTON_GRADIENT)
383  os << "\tNL Strategy:\t\t" << strategy << "\t(Newton/Gradient)" << std::endl;
384  else if (strategy == MOD_NEWTON)
385  os << "\tNL Strategy:\t\t" << strategy << "\t(Modified-Newton)" << std::endl;
386  else if (strategy == MOD_NEWTON_GRADIENT)
387  os << "\tNL Strategy:\t\t" << strategy << "\t(Modified-Newton/Gradient)" << std::endl;
388 
389  LineSearchMethod searchMethod = getSearchMethod();
390  if (searchMethod == FULL)
391  os << "\tsearch method:\t\t" << searchMethod << "\t(None => Full Newton Steps)" << std::endl;
392 
393  else if (searchMethod == DIVIDE)
394  os << "\tsearch method:\t\t" << searchMethod << "\t(Divide)" << std::endl;
395 
396  else if (searchMethod == BACKTRACK)
397  os << "\tsearch method:\t\t" << searchMethod << "\t(Backtrack)" << std::endl;
398 
399  else if (searchMethod == SIMPLE_BACKTRACK)
400  os << "\tsearch method:\t\t" << searchMethod << "\t(Simple Backtrack)" << std::endl;
401 
402  else if (searchMethod == BANK_ROSE)
403  os << "\tsearch method:\t\t" << searchMethod << "\t(Bank and Rose Algorithm)" << std::endl;
404 
405  else if (searchMethod == DESCENT)
406  os << "\tsearch method:\t\t" << searchMethod << "\t(Line Search)" << std::endl;
407 
408  os << "\tmax search steps:\t" << getMaxSearchStep()
409  << "\tinexact-Newton forcing:\t" << getForcingFlag()
410  << "\tnorm level:\t\t" << getNormLevel()
411  << "\tlinear optimization:\t" << getLinearOpt()
412  << "\tconstraint backtrack:\t" << getConstraintBT();
413 
414  if (DEBUG_NONLINEAR)
415  {
416  os << "\tdebugLevel:\t\t" << getDebugLevel ()
417  << "\tdebugMinTimeStep:\t" << getDebugMinTimeStep ()
418  << "\tdebugMaxTimeStep:\t" << getDebugMaxTimeStep ();
419  }
420 
421  os << Xyce::section_divider << "\n" << std::endl;
422 }
423 
424 //-----------------------------------------------------------------------------
425 // Function : NLParams::operator=
426 // Purpose : "=" operator.
427 // Special Notes :
428 // Scope : public
429 // Creator : Eric R. Keiter, SNL, Computational Sciences
430 // Creation Date : 09/05/01
431 //-----------------------------------------------------------------------------
433 {
434  if (this != &right)
435  {
436  commandLine_ = right.commandLine_;
437  nlStrategy_ = right.nlStrategy_;
439  direction_ = right.direction_;
440  deltaXTol_ = right.deltaXTol_;
441  RHSTol_ = right.RHSTol_;
442  absTol_ = right.absTol_;
443  relTol_ = right.relTol_;
447  eta_ = right.eta_;
448  normLevel_ = right.normLevel_;
449 
451 
453 
455  globalBTMax_ = right.globalBTMax_;
456  globalBTMin_ = right.globalBTMin_;
458 
459  // Debug output options:
460  debugLevel_ = right.debugLevel_;
465  }
466 
467  return *this;
468 }
469 
470 //-----------------------------------------------------------------------------
471 // Function : NLParams::setCmdLineOptions
472 // Purpose :
473 // Special Notes :
474 // Scope : public
475 // Creator : Eric R. Keiter, SNL, Computational Sciences
476 // Creation Date : 05/24/05
477 //-----------------------------------------------------------------------------
478 
480 {
481  // set (or override) debug levels based on command line options
482  if (DEBUG_NONLINEAR && commandLine_->argExists( "-ndl" ))
483  setDebugLevel(commandLine_->getArgumentIntValue("-ndl", getDebugLevel()));
484 
485  return true;
486 }
487 
488 } // namespace Nonlinear
489 } // namespace Xyce
490 
void setForcingFlag(bool flag)
void setSearchMethod(LineSearchMethod method)
void setGlobalBTMin(double value)
Pure virtual class to augment a linear system.
double getSmallUpdateTol() const
void setDebugMinTime(double value)
unsigned getMaxNewtonStep() const
void setMaxSearchStep(unsigned maxSearchStep)
void setSmallUpdateTol(double Tolerance)
void setScreenOutputFlag(bool bval)
void setEnforceDeviceConvFlag(bool flag)
void setNLStrategy(NLStrategy strategy)
void setConstraintBT(bool flag)
NLParams(AnalysisMode mode, const IO::CmdParse &cp)
void setDebugMinTimeStep(int value)
void setGlobalBTChange(double value)
LineSearchMethod getSearchMethod() const
void setDeltaXTol(double Tolerance)
NLParams & operator=(const NLParams &right)
LineSearchMethod searchMethod_
void setRelTol(double Tolerance)
void setDebugLevel(int value)
void setMaskingFlag(bool bval)
void setDebugMaxTimeStep(int value)
void setNormLevel(int level)
void setAbsTol(double Tolerance)
void setGlobalBTMax(double value)
const IO::CmdParse * commandLine_
void setMMFormat(bool value)
double getDeltaXTol() const
void printParams(std::ostream &os)
void setLinearOpt(bool flag)
void setMaxNewtonStep(unsigned maxNewtonStep)
NLStrategy getNLStrategy() const
bool setOptions(const Util::OptionBlock &OB)
void setRHSTol(double Tolerance)
void setDebugMaxTime(double value)
unsigned getMaxSearchStep() const