Xyce  6.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
N_NLS_NOX_XyceTests.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-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_NOX_XyceTests.C,v $
27 //
28 // Purpose : Status test.
29 //
30 // Special Notes :
31 //
32 // Creator : Roger Pawlowski, SNL 9233
33 //
34 // Creation Date : 04/15/03
35 //
36 // Revision Information:
37 // ---------------------
38 //
39 // Revision Number: $Revision: 1.44 $
40 //
41 // Revision Date : $Date: 2014/02/24 23:49:25 $
42 //
43 // Current Owner : $Author: tvrusso $
44 //-------------------------------------------------------------------------
45 
46 #include <Xyce_config.h>
47 
48 
49 // ---------- Standard Includes ----------
50 
51 // ---------- Xyce Includes ----------
52 #include <N_UTL_Misc.h>
53 #include "N_NLS_NOX_XyceTests.h"
54 #include "N_NLS_NOX_Vector.h"
55 #include "N_LAS_Vector.h"
56 #include "N_LOA_Loader.h"
57 #include "NOX.H"
58 #include "NOX_Solver_LineSearchBased.H"
59 
60 // ---------- Namespaces ----------
61 
62 using namespace N_NLS_NOX;
63 
64 // ---------- Code ----------
65 
66 XyceTests::XyceTests(bool isTransient,
67  double normF,
68  double machPrec,
69  N_LAS_Vector** currSolVectorPtrPtr,
70  double epsilon_a,
71  double epsilon_r,
72  double tol,
73  int maxIters,
74  double convRate,
75  double relConvRate,
76  double maxConvRate,
77  double stagnationTol,
78  int maxBadSteps,
79  int checkDeviceConvergence,
80  double smallUpdateTol,
81  N_LOA_Loader* loader
82 #ifdef Xyce_NLS_MASKED_WRMS_NORMS
83  , bool nonTrivialDeviceMaskFlag,
84  N_LAS_Vector * maskVectorPtr
85 #endif
86  ) :
87 
88  status_(NOX::StatusTest::Unconverged),
89  returnTest_(0),
90  isTransient_(isTransient),
91  niters_(-1),
92  maxNormFindex_(-1),
93  maxNormF_(0.0),
94  requestedMaxNormF_(normF),
95  requestedMachPrecTol_(machPrec),
96  oldTimeStepVectorPtrPtr_(currSolVectorPtrPtr),
97  weightsVectorPtr_(0),
98  updateVectorPtr_(0),
99  tmpVectorPtr_(0),
100  epsilon_a_(epsilon_a),
101  epsilon_r_(epsilon_r),
102  tol_(tol),
103  weightedUpdate_(0.0),
104  maxIters_(maxIters),
105  requestedConvRate_(convRate),
106  currentConvRate_(1.0),
107  requestedRelativeConvRate_(relConvRate),
108  currentRelativeConvRate_(1.0),
109  normResidualInit_(1.0),
110  maxConvRate_(maxConvRate),
111  lastIteration_(-1),
112  badStepCount_(0),
113  minConvRate_(1.0),
114  stagnationTol_(stagnationTol),
115  maxBadSteps_(maxBadSteps),
116  xyceReturnCode_(0),
117  smallUpdateTol_(smallUpdateTol),
118  checkDeviceConvergence_(checkDeviceConvergence),
119  loaderPtr_(loader),
120  allDevicesConverged_(false),
121  innerDevicesConverged_(false)
122 #ifdef Xyce_NLS_MASKED_WRMS_NORMS
123  , deviceMaskFlag_( false ),
124  weightMaskVectorPtr_( maskVectorPtr)
125 #endif
126 {
127 
128 }
129 
131 {
132  if( weightsVectorPtr_ != 0 )
133  {
134  delete weightsVectorPtr_;
135  delete updateVectorPtr_;
136  delete tmpVectorPtr_;
137  }
138 }
139 
140 
141 NOX::StatusTest::StatusType XyceTests::
142 checkStatus(const NOX::Solver::Generic& problem,
143  NOX::StatusTest::CheckType checkType)
144 {
145  status_ = NOX::StatusTest::Unconverged;
146  xyceReturnCode_ = 0;
147  niters_ = problem.getNumIterations();
148 
149  returnTest_ = 0;
150 
151 
152 
153 
154  /*
155  problem.getSolutionGroup();
156  const NOX::Abstract::Vector& v = problem.getSolutionGroup().getX();
157  cout << "ROGER - norm = " << v.norm() << endl;
158  const N_NLS_NOX::Vector* testPtr = 0;
159  testPtr = dynamic_cast<const N_NLS_NOX::Vector*>(&v);
160  cout << "ROGER - ptr = " << testPtr << endl;
161 
162  if (testPtr == 0) {
163  const N_NLS_NOX::Vector* testPtr = 0;
164  testPtr = dynamic_cast<const N_NLS_NOX::Vector*>(&v);
165  }
166  //const N_NLS_NOX::Vector& v_nox = dynamic_cast<const N_NLS_NOX::Vector&>
167  //(problem.getSolutionGroup().getX());
168 
169  //const N_LAS_Vector& x__ = (dynamic_cast<const N_NLS_NOX::Vector&>
170  //(problem.getSolutionGroup().getX())).getNativeVectorRef();
171  exit(0);
172  */
173 
174 
175  // Get the current and previous solutions
176  const N_LAS_Vector& x = (dynamic_cast<const N_NLS_NOX::Vector&>
177  (problem.getSolutionGroup().getX())).getNativeVectorRef();
178  const N_LAS_Vector& oldX = (dynamic_cast<const N_NLS_NOX::Vector&>
179  (problem.getPreviousSolutionGroup().getX())).getNativeVectorRef();
180 
181  // Test0 - NaN/Inf checker
182  NOX::StatusTest::StatusType check = finiteTest_.checkStatus(problem,
183  checkType);
184  if (check == NOX::StatusTest::Failed)
185  {
186  status_ = check;
187  returnTest_ = 0;
188  xyceReturnCode_ = retCodes_.nanFail; // default: -6
189  return status_;
190  }
191 
192  // This test is for 2-level solves only.
193  // If the inner solve failed, then the whole thing needs to fail.
194  // If 2-level is not being used, this test doesn't do anything.
195  innerDevicesConverged_ = loaderPtr_->innerDevsConverged();
197  {
198  status_ = NOX::StatusTest::Failed;
199  returnTest_ = 9;
200  xyceReturnCode_ = retCodes_.innerSolveFailed; // default: -5
201  return status_;
202  }
203 
204  // Test 8 - Devices need to satisfy their own convergence criteria
206  {
207  allDevicesConverged_ = loaderPtr_->allDevsConverged();
209  {
210  status_ = NOX::StatusTest::Unconverged;
211  returnTest_ = 8;
212  xyceReturnCode_ = 0;
213  return status_;
214  }
215  }
216 
217  //Test 1 - Make sure the residual isn't too small, hardwired tolerances
218  maxNormF_ = problem.getSolutionGroup().getF()
219  .norm(NOX::Abstract::Vector::MaxNorm);
220 
221  const N_LAS_Vector& F = (dynamic_cast<const N_NLS_NOX::Vector&>
222  (problem.getSolutionGroup().getF())).getNativeVectorRef();
223 
224  std::vector<int> index(1, -1);
225  F.infNormIndex( &index[0] );
226  maxNormFindex_ = index[0];
227 
228  if ((maxNormF_ < requestedMaxNormF_) &&
230  {
231  status_ = NOX::StatusTest::Converged;
232  returnTest_ = 1;
233  xyceReturnCode_ = retCodes_.normTooSmall; // default: 1
234  return status_;
235  }
236 
237  // Test 2 - Normal convergence based on rhs residual (2a) and
238  // update norm (2b).
239 
240  // Copy into local reference
241  N_LAS_Vector& oldTimeStepX = **oldTimeStepVectorPtrPtr_;
242 
243  // Allocate space if necessary
244  if (weightsVectorPtr_ == 0)
245  {
246  weightsVectorPtr_ = new N_LAS_Vector(x);
247  // when we create weightsVectorPtr_ from the latest solution, there
248  // is a chance that one of the values will be zero. If this isn't
249  // the DC op step or the first iteration of a time step, then
250  // we'll end up dividing by zero in the wMaxNorm function below.
251  // So to be safe we'll just add epsilon_a_ on the when we create
252  // this vector.
253  for (int i=0; i< x.localLength() ; ++i )
254  {
255  (*weightsVectorPtr_)[i] += epsilon_a_;
256  }
257  updateVectorPtr_ = new N_LAS_Vector(x);
258  tmpVectorPtr_ = new N_LAS_Vector(x);
259  }
260 
261  // Local references
262  N_LAS_Vector& weights = *weightsVectorPtr_;
263  N_LAS_Vector& update = *updateVectorPtr_;
264  N_LAS_Vector& tmp = *tmpVectorPtr_;
265 
266  // Compute local portion of weights vector
267  // Weights are recomputed at each nonlinear iteration of a DC Op calc
268  // but only at the beginning of a transient nonlinear solve.
269  if ((!isTransient_) || (niters_ == 0))
270  {
271  int length = x.localLength();
272  for (int i = 0; i < length; ++i )
273  {
274  //update[i] = x[i] - oldX[i];
275  weights[i] =
276  epsilon_r_ * Xycemax(fabs(x[i]), fabs(oldTimeStepX[i])) + epsilon_a_;
277 #ifdef Xyce_NLS_MASKED_WRMS_NORMS
278  if( deviceMaskFlag_ && ((*weightMaskVectorPtr_)[i] == 0.0) )
279  {
280  weights[i] = N_UTL_MachineDependentParams::MachineBig();
281  }
282 #endif
283  }
284  }
285 
286 
287 
288  if (niters_ < 1)
289  {
290  weightedUpdate_ = 1.0;
291  }
292  else
293  {
294  // Next compute the update
295  update.update(1.0, x, -1.0, oldX, 0.0);
296 
297  // Compute final result
298 #ifdef Xyce_SPICE_NORMS
299  update.wMaxNorm(weights,tmp,&weightedUpdate_);
300 #else
301  update.wRMSNorm(weights,&weightedUpdate_);
302 #endif
303 
304  // RPP: If a line search is being used, we must account for any
305  // damping of the step length. Otherwise delta X could be small due
306  // the line search and not due to being close to a solution.
307  const NOX::Solver::LineSearchBased* test = 0;
308  test = dynamic_cast<const NOX::Solver::LineSearchBased*>(&problem);
309  if (test != 0)
310  {
311  weightedUpdate_ = weightedUpdate_/(test->getStepSize());
312  }
313 
314  // RPP: 11/11/2003 - Fix for Bug 354 - Xyce fails in DC Op calc
315  // but proceeds to transient.
316  // Check to see if WRMS is exactly zero. If so the linear solver
317  // has failed. Return a failure. This should be commented out for
318  // Broyden runs.
319  // RPP: 02/18/2004 - Causing premature failures if deltaX is ~1e-16.
320  // Adding a check to look at norm of dx
321  // RPP: 08/09/2004 - Still goes into transient from failed DC Op on
322  // freebsd platforms for Down_8-bit_03.cir (bug 354). Changed return
323  // code from +4 to -4. This forces TIA to assume failure.
324  if ((weightedUpdate_ == 0.0) && (niters_ > 0))
325  {
326  if (!(problem.getPreviousSolutionGroup().isNewton()))
327  {
328  NOX::Abstract::Group & tmpGrp =
329  (const_cast<NOX::Abstract::Group&>(problem.getPreviousSolutionGroup()));
330  Teuchos::ParameterList tmpParams;
331  tmpGrp.computeNewton (tmpParams);
332  }
333 
334  const N_LAS_Vector& dx = (dynamic_cast<const N_NLS_NOX::Vector&>
335  (problem.getPreviousSolutionGroup().getNewton())).getNativeVectorRef();
336 
337  double tmp = 0.0;
338  dx.lpNorm(2, &tmp );
339  if (tmp == 0.0)
340  {
341  status_ = NOX::StatusTest::Failed;
342  returnTest_ = 4;
343  xyceReturnCode_ = retCodes_.wrmsExactZero; // default: -4
344  return status_;
345  }
346  }
347 
349  {
350  status_ = NOX::StatusTest::Converged;
351  returnTest_ = 2;
353  return status_;
354  }
355  }
356 
357  // Test 3 - Near Convergence - Hit max iterations but residual
358  // and convergence rate indicate we may be near a converged solution.
359  // Therefore, let the time stepper decide whether or not the step is ok.
360  // Transient mode ONLY!
361  // NOTE: Convergence rates are based on the 2-Norm, not max norm!
362  if (niters_ > 0)
363  {
364  // ||F(x_current)|| / ||F(x_previous)||
365  currentConvRate_ = (problem.getSolutionGroup().getNormF()) /
366  (problem.getPreviousSolutionGroup().getNormF());
367 
368  // ||F(x)|| / ||F(x_init)||
369  currentRelativeConvRate_ = (problem.getSolutionGroup().getNormF()) /
371  }
372  else
373  {
374  currentConvRate_ = 1.0;
376  }
377 
378  if (isTransient_)
379  {
380  if (niters_ == 0)
381  {
382  normResidualInit_ = problem.getSolutionGroup().getNormF();
383  }
384 
385  // Test only if we hit the max number of iterations
386  if (niters_ >= maxIters_)
387  {
390  {
391  status_ = NOX::StatusTest::Converged;
392  xyceReturnCode_ = retCodes_.nearConvergence; // default: 3
393 
394  if (xyceReturnCode_ < 0)
395  status_ = NOX::StatusTest::Failed;
396  else
397  status_ = NOX::StatusTest::Converged;
398  }
399  else
400  {
401  status_ = NOX::StatusTest::Failed;
402  xyceReturnCode_ = retCodes_.tooManySteps; // default: -1
403  }
404 
405  returnTest_ = 3;
406  return status_;
407  }
408  } // end test 3
409 
410  // Test 4 - Update is too small
411  if ((niters_ > 0) && (weightedUpdate_ < smallUpdateTol_) && (niters_ < maxIters_))
412  {
413  if (isTransient_) // Let the time integrator determine convergence (+4)
414  {
415  xyceReturnCode_ = retCodes_.smallUpdate; // default: 4
416  status_ = NOX::StatusTest::Failed;
417  }
418  else // Steady state should always return a hard failure -4
419  {
420  xyceReturnCode_ = 0; // neither pass nor fail.
421  status_ = NOX::StatusTest::Unconverged;
422 
423  //xyceReturnCode_ = retCodes_.wrmsExactZero; // default: -4
424  //status_ = NOX::StatusTest::Failed;
425  }
426 
427  returnTest_ = 4;
428 
429  return status_;
430  }
431 
432  // Test 5 - Max nonlinear iterations (if transient, this will be checked
433  // in the NearConvergence test (#3)
434  if ((!isTransient_) && (niters_ >= maxIters_))
435  {
436  status_ = NOX::StatusTest::Failed;
437  returnTest_ = 5;
438  xyceReturnCode_ = retCodes_.tooManySteps; // default: -1
439  return status_;
440  }
441 
442  // Test 6 - update is too big
444  {
445  status_ = NOX::StatusTest::Failed;
446  returnTest_ = 6;
447  xyceReturnCode_ = retCodes_.updateTooBig; // default: -2
448  return status_;
449  }
450 
451  // Test 7 - Stall in the convergence rate. Transient mode ONLY!
452  if (isTransient_)
453  {
454  // First time through we don't do anything but reset the counters
455  if (niters_ == 0)
456  {
457  badStepCount_ = 0;
458  lastIteration_ = 0;
459  //minConvRate = 1.0; // Don't reset this. Xyce solver never does.
460  }
461 
462  // Make sure we have not already counted the last nonlinear iteration.
463  // This protects against multiple calls to checkStatus() in between
464  // nonlinear iterations.
465  bool isCounted = false;
466  if (niters_ == lastIteration_)
467  {
468  isCounted = true;
469  }
470  else
471  {
473  }
474 
475  // Set counter appropriately
476  if (!isCounted)
477  {
478  if (fabs(currentConvRate_ - 1.0) <= stagnationTol_)
479  {
480  if ((badStepCount_ == 0) || (currentConvRate_ < minConvRate_))
481  {
483  }
484  ++badStepCount_ ;
485  }
486  else
487  {
488  badStepCount_ = 0;
489  }
490  }
491 
492  if (badStepCount_ >= maxBadSteps_)
493  {
494  if ((currentRelativeConvRate_ <= 0.9) && (minConvRate_ <= 1.0))
495  {
496  status_ = NOX::StatusTest::Converged;
497  returnTest_ = 7;
498  xyceReturnCode_ = retCodes_.nearConvergence; // default: 3
499  // note - I'm not sure if this is
500  // really a near convergece test - but 3 is the code for it...
501 
502  if (xyceReturnCode_ < 0)
503  status_ = NOX::StatusTest::Failed;
504  else
505  status_ = NOX::StatusTest::Converged;
506  }
507  else
508  {
509  status_ = NOX::StatusTest::Failed;
510  returnTest_ = 7;
511  xyceReturnCode_ = retCodes_.stalled; // default: -3
512  }
513  }
514  }
515 
516  return status_;
517 }
518 
519 std::ostream& XyceTests::print(std::ostream& stream, int indent) const
520 {
521  // precision
522  int p = 5;
523 
524  for (int j = 0; j < indent; ++j )
525  stream << ' ';
526  stream << status_ << "by Test #" << returnTest_ << "\n";
527 
528  indent += 4;
529 
530  //for (int j = 0; j < indent; ++j )
531  // stream << ' ';
532  finiteTest_.print(stream, indent);
533 
535  for (int j = 0; j < indent; ++j )
536  stream << ' ';
537  stream << "8. Devices are Converged: ";
539  stream << "true" << "\n";
540  else
541  stream << "false" << "\n";
542  }
543 
544  for (int j = 0; j < indent; ++j )
545  stream << ' ';
546  stream << "1. Inf-Norm F too small" << "\n";
547 
548  for (int j = 0; j < indent; ++j )
549  stream << ' ';
550  stream << " Machine Precision: " << NOX::Utils::sciformat(maxNormF_, p)
551  << " < " << NOX::Utils::sciformat(requestedMachPrecTol_, p) << "\n";
552 
553  for (int j = 0; j < indent; ++j )
554  stream << ' ';
555  stream << " Requested Tolerance: " << NOX::Utils::sciformat(maxNormF_, p)
556  << " < " << NOX::Utils::sciformat(requestedMaxNormF_, p) << "\n";
557 
558  for (int j = 0; j < indent; ++j )
559  stream << ' ';
560  stream << "2. Normal Convergence" << "\n";
561 
562  for (int j = 0; j < indent; ++j )
563  stream << ' ';
564  stream << " Inf-Norm F: " << NOX::Utils::sciformat(maxNormF_, p)
565  << " < " << NOX::Utils::sciformat(requestedMaxNormF_, p) << "\n";
566 
567  for (int j = 0; j < indent; ++j )
568  stream << ' ';
569  stream << " Weighted Update: " << NOX::Utils::sciformat(weightedUpdate_, p)
570  << " < " << NOX::Utils::sciformat(tol_, p) << "\n";
571 
572  for (int j = 0; j < indent; ++j )
573  stream << ' ';
574  stream << "3. Near Convergence" << "\n";
575 
576  for (int j = 0; j < indent; ++j )
577  stream << ' ';
578  stream << " Max Iters: " << niters_
579  << " < " << maxIters_ << "\n";
580 
581  for (int j = 0; j < indent; ++j )
582  stream << ' ';
583  stream << " Convergence Rate: "
584  << NOX::Utils::sciformat(currentConvRate_, p)
585  << " < " << NOX::Utils::sciformat(requestedConvRate_, p) << "\n";
586 
587  for (int j = 0; j < indent; ++j )
588  stream << ' ';
589  stream << " Relative Convergence Rate: "
590  << NOX::Utils::sciformat(currentRelativeConvRate_, p)
591  << " < " << NOX::Utils::sciformat(requestedRelativeConvRate_, p)
592  << "\n";
593 
594  for (int j = 0; j < indent; ++j )
595  stream << ' ';
596  stream << "4. Small Weighted Update: "
597  << NOX::Utils::sciformat(weightedUpdate_, p)
598  << " < " << NOX::Utils::sciformat(smallUpdateTol_, p) << "\n";
599 
600  if (!isTransient_) {
601 
602  for (int j = 0; j < indent; ++j )
603  stream << ' ';
604  stream << "5. Maximum Iterations: "
605  << niters_
606  << " < " << maxIters_ << "\n";
607 
608  for (int j = 0; j < indent; ++j )
609  stream << ' ';
610  stream << "6. Large Conv Rate: "
611  << NOX::Utils::sciformat(currentConvRate_, p)
612  << " < " << NOX::Utils::sciformat(maxConvRate_, p) << "\n";
613  }
614 
615  for (int j = 0; j < indent; ++j )
616  stream << ' ';
617  stream << "7. Stagnation " << "\n";
618 
619  for (int j = 0; j < indent; ++j )
620  stream << ' ';
621  stream << " Bad Step Count: "
622  << badStepCount_ << " < " << maxBadSteps_ << "\n";
623 
624  for (int j = 0; j < indent; ++j )
625  stream << ' ';
626  stream << " Stagnation Tolerance: "
627  << NOX::Utils::sciformat(fabs(currentConvRate_ - 1.0), p)
628  << " < " << NOX::Utils::sciformat(stagnationTol_, p) << "\n";
629 
630  stream << std::endl;
631  return stream;
632 }
633 
635 {
636  return xyceReturnCode_;
637 }
638 
640 {
641  return maxNormF_;
642 }
643 
645 {
646  return maxNormFindex_;
647 }
648