Xyce  6.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
N_NLS_TwoLevelNewton.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_TwoLevelNewton.C,v $
27 //
28 // Purpose : Body for the two level Newton class.
29 //
30 // Special Notes :
31 //
32 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
33 //
34 // Creation Date : 10/20/02
35 //
36 // Revision Information:
37 // ---------------------
38 //
39 // Revision Number: $Revision: 1.135.2.1 $
40 //
41 // Revision Date : $Date $
42 //
43 // Current Owner : $Author: dgbaur $
44 //-------------------------------------------------------------------------
45 
46 #include <Xyce_config.h>
47 
48 
49 // ---------- Standard Includes ----------
50 #include <stdio.h>
51 #include <string>
52 #include <list>
53 
54 // ---------- Xyce Includes ----------
55 
56 #include <N_NLS_TwoLevelNewton.h>
57 #include <N_NLS_Manager.h>
58 #include <N_NLS_NOX_Interface.h>
59 #include <N_NLS_DampedNewton.h>
60 #include <N_NLS_ReturnCodes.h>
61 
62 #include <N_ANP_AnalysisManager.h>
63 #include <N_TIA_TimeIntInfo.h>
64 #include <N_LOA_Loader.h>
65 
66 #include <N_LAS_System.h>
67 #include <N_LAS_Vector.h>
68 #include <N_LAS_Matrix.h>
69 #include <N_LAS_Builder.h>
70 
71 #include <N_ERH_ErrorMgr.h>
72 #include <N_UTL_Xyce.h>
73 #include <N_UTL_Param.h>
74 #include <N_UTL_OptionBlock.h>
75 
76 #include <N_IO_OutputMgr.h>
77 
78 // ---------- Static Declarations ----------
79 
80 namespace Xyce {
81 namespace Nonlinear {
82 
83 //-----------------------------------------------------------------------------
84 // Function : TwoLevelNewton::TwoLevelNewton
85 // Purpose : constructor
86 // Special Notes :
87 // Scope : public
88 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
89 // Creation Date : 10/20/02
90 //-----------------------------------------------------------------------------
92  (bool noxFlag, bool noxFlagInner, N_IO_CmdParse & cp)
93  :
94  NonLinearSolver(cp),
95  twoLevelAlgorithm_(3),
96  twoLevelAlgorithmTran_(0),
97  maxOuterSteps_(20),
98  maxContSteps_(10),
99  contStep_(0),
100  setupOuterLoopParamsFlag_(false),
101  setupTranParamsFlag_(false),
102  externalAnalysisMode(DC_OP),
103  outerLoopActiveFlag_(true),
104  noxFlag_(noxFlag),
105  noxFlagInner_(noxFlagInner),//inner loop possibly different than outer loop.
106  numInterfaceNodesSetup_(false),
107  twoLevelCouplingMode_(FULL_PROBLEM),
108  savedRHSPtr_(0),
109  savedNextSolPtr_(0),
110  jdxpVectorPtr_(0),
111  numSubProblems_(0),
112  continuationType_(1),
113  innerLoopFailFatal_(true),
114  totalSolveFailFatal_(false),
115  doFullNewtonFinalEnforcement_(true),
116  nlsPassingPtr_(0),
117  firstDCOPFlag_(false),
118  increaseContScalar_(1.5),
119  decreaseContScalar_(0.2), // was 0.125
120  continuationCalledBefore_(false),
121  voltLimTol_(1.0e-6)
122 {
123  // allocate the "outer loop" and "inner loop" solvers.
124  if (noxFlag_)
125  {
126  nlsOuterPtr_ = new N_NLS_NOX::Interface(commandLine_);
127  }
128  else
129  {
130  nlsOuterPtr_ = new DampedNewton(commandLine_);
131  }
132 
133  if (noxFlagInner_)
134  {
135  nlsInnerPtr_ = new N_NLS_NOX::Interface(commandLine_);
136  }
137  else
138  {
139  nlsInnerPtr_ = new DampedNewton(commandLine_);
140  }
141 
142  nlsOuterPtr_->registerTwoLevelSolver(this);
143  nlsInnerPtr_->registerTwoLevelSolver(this);
144 
145 }
146 
147 //-----------------------------------------------------------------------------
148 // Function : TwoLevelNewton::~TwoLevelNewton
149 // Purpose : destructor
150 // Special Notes :
151 // Scope : public
152 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
153 // Creation Date : 10/20/02
154 //-----------------------------------------------------------------------------
155 
157 {
158  delete nlsOuterPtr_;
159  delete nlsInnerPtr_;
160 
161  if (savedRHSPtr_!=0) delete savedRHSPtr_;
162  if (savedNextSolPtr_!=0) delete savedNextSolPtr_;
163 }
164 
165 //-----------------------------------------------------------------------------
166 // Function : TwoLevelNewton::registerAnalysisManager
167 // Purpose :
168 // Special Notes :
169 // Scope : public
170 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
171 // Creation Date : 10/20/02
172 //-----------------------------------------------------------------------------
174  (N_ANP_AnalysisManager * tiaPtr_tmp)
175 {
176  bool bsuccess = true;
177  bool tmpBool = true;
178 
179  tiaPtr_ = tiaPtr_tmp;
180  tmpBool = nlsOuterPtr_->registerAnalysisManager(tiaPtr_);
181  bsuccess = bsuccess && tmpBool;
182 
183  tmpBool = nlsInnerPtr_->registerAnalysisManager(tiaPtr_);
184  bsuccess = bsuccess && tmpBool;
185 
186  return bsuccess;
187 }
188 
189 //-----------------------------------------------------------------------------
190 // Function : TwoLevelNewton::registerLinearSystem
191 // Purpose :
192 // Special Notes :
193 // Scope : public
194 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
195 // Creation Date : 10/20/02
196 //-----------------------------------------------------------------------------
198  (N_LAS_System * ptr)
199 {
200  bool bsuccess = true;
201  bool tmpBool = true;
202 
203  lasSysPtr_ = ptr;
204  tmpBool = nlsOuterPtr_->registerLinearSystem (ptr);
205  bsuccess = bsuccess && tmpBool;
206 
207  tmpBool = nlsInnerPtr_->registerLinearSystem (ptr);
208  bsuccess = bsuccess && tmpBool;
209 
210  return bsuccess;
211 }
212 
213 //-----------------------------------------------------------------------------
214 // Function : TwoLevelNewton::registerLoader
215 // Purpose :
216 // Special Notes :
217 // Scope : public
218 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
219 // Creation Date : 10/20/02
220 //-----------------------------------------------------------------------------
222  (N_LOA_Loader * loaderPtr_tmp)
223 {
224  bool bsuccess = true;
225  bool tmpBool = true;
226 
227  loaderPtr_ = loaderPtr_tmp;
228  tmpBool = nlsOuterPtr_->registerLoader (loaderPtr_tmp);
229  bsuccess = bsuccess && tmpBool;
230 
231  tmpBool = nlsInnerPtr_->registerLoader (loaderPtr_tmp);
232  bsuccess = bsuccess && tmpBool;
233 
234  return bsuccess;
235 }
236 
237 //-----------------------------------------------------------------------------
238 // Function : TwoLevelNewton::registerOutputMgr
239 // Purpose :
240 // Special Notes :
241 // Scope : public
242 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
243 // Creation Date : 9/23/03
244 //-----------------------------------------------------------------------------
245 bool TwoLevelNewton::registerOutputMgr (N_IO_OutputMgr * ptr)
246 {
247  bool bsuccess = true;
248  bool tmpBool = true;
249 
250  outMgrPtr_ = ptr;
251  tmpBool = nlsOuterPtr_->registerOutputMgr (ptr);
252  bsuccess = bsuccess && tmpBool;
253 
254  tmpBool = nlsInnerPtr_->registerOutputMgr (ptr);
255  bsuccess = bsuccess && tmpBool;
256 
257  return bsuccess;
258 }
259 
260 
261 //-----------------------------------------------------------------------------
262 // Function : TwoLevelNewton::getCouplingMode
263 // Purpose :
264 // Special Notes :
265 // Scope : public
266 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
267 // Creation Date : 12/05/02
268 //-----------------------------------------------------------------------------
270 {
271  return twoLevelCouplingMode_;
272 }
273 
274 //-----------------------------------------------------------------------------
275 // Function : TwoLevelNewton::getNumResidualLoads
276 // Purpose :
277 // Special Notes :
278 // Scope : public
279 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
280 // Creation Date : 10/20/02
281 //-----------------------------------------------------------------------------
283 {
284  int numResLoads = 0;
285  numResLoads += nlsOuterPtr_->getNumResidualLoads ();
286  numResLoads += numResidualLoads_;
287  return numResLoads;
288 }
289 
290 //-----------------------------------------------------------------------------
291 // Function : TwoLevelNewton::getNumJacobianLoads
292 // Purpose :
293 // Special Notes :
294 // Scope : public
295 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
296 // Creation Date : 10/20/02
297 //-----------------------------------------------------------------------------
299 {
300  int numJacLoads = 0;
301  numJacLoads += nlsOuterPtr_->getNumJacobianLoads ();
302  numJacLoads += numJacobianLoads_;
303  return numJacLoads;
304 }
305 
306 //-----------------------------------------------------------------------------
307 // Function : TwoLevelNewton::getNumLinearSolves
308 // Purpose :
309 // Special Notes :
310 // Scope : public
311 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
312 // Creation Date : 10/20/02
313 //-----------------------------------------------------------------------------
315 {
316  int numLinSolves = 0;
317  numLinSolves += nlsOuterPtr_->getNumLinearSolves ();
318  numLinSolves += numLinearSolves_;
319  return numLinSolves;
320 }
321 
322 //-----------------------------------------------------------------------------
323 // Function : TwoLevelNewton::getNumFailedLinearSolves
324 // Purpose :
325 // Special Notes :
326 // Scope : public
327 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
328 // Creation Date : 10/20/02
329 //-----------------------------------------------------------------------------
331 {
332  int numFLinSolves = 0;
333  numFLinSolves += nlsOuterPtr_->getNumFailedLinearSolves ();
334  numFLinSolves += numFailedLinearSolves_;
335  return numFLinSolves;
336 }
337 
338 //-----------------------------------------------------------------------------
339 // Function : TwoLevelNewton::getNumJacobianFactorizations
340 // Purpose :
341 // Special Notes :
342 // Scope : public
343 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
344 // Creation Date : 10/20/02
345 //-----------------------------------------------------------------------------
347 {
348  int numJacFact = 0;
349  numJacFact += nlsOuterPtr_->getNumJacobianFactorizations ();
350  numJacFact += numJacobianFactorizations_;
351  return numJacFact;
352 }
353 
354 //-----------------------------------------------------------------------------
355 // Function : TwoLevelNewton::getTotalNumLinearIters
356 // Purpose :
357 // Special Notes :
358 // Scope : public
359 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
360 // Creation Date : 10/20/02
361 //-----------------------------------------------------------------------------
363 {
364  int numLinIters = 0;
365  numLinIters += nlsOuterPtr_->getTotalNumLinearIters ();
366  numLinIters += totalNumLinearIters_;
367  return numLinIters ;
368 }
369 
370 //-----------------------------------------------------------------------------
371 // Function : TwoLevelNewton::getTotalLinearSolveTime
372 // Purpose :
373 // Special Notes :
374 // Scope : public
375 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
376 // Creation Date : 10/20/02
377 //-----------------------------------------------------------------------------
379 {
380  double totalLinSolveTime = 0.0;
381  totalLinSolveTime += nlsOuterPtr_->getTotalLinearSolveTime();
382  totalLinSolveTime += totalLinearSolveTime_;
383  return totalLinSolveTime;
384 }
385 
386 //-----------------------------------------------------------------------------
387 // Function : TwoLevelNewton::getTotalResidualLoadTime
388 // Purpose :
389 // Special Notes :
390 // Scope : public
391 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
392 // Creation Date : 10/20/02
393 //-----------------------------------------------------------------------------
395 {
396  double totalResLoadTime = 0.0;
397  totalResLoadTime += nlsOuterPtr_-> getTotalResidualLoadTime();
398  totalResLoadTime += totalResidualLoadTime_;
399  return totalResLoadTime;
400 }
401 
402 //-----------------------------------------------------------------------------
403 // Function : TwoLevelNewton::getTotalJacobianLoadTime
404 // Purpose :
405 // Special Notes :
406 // Scope : public
407 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
408 // Creation Date : 10/20/02
409 //-----------------------------------------------------------------------------
411 {
412  double totalJacLoadTime = 0.0;
413  totalJacLoadTime += nlsOuterPtr_->getTotalJacobianLoadTime();
414  totalJacLoadTime += totalJacobianLoadTime_;
415  return totalJacLoadTime;
416 }
417 
418 //-----------------------------------------------------------------------------
419 // Function : TwoLevelNewton::getNumIterations
420 // Purpose : Returns the current number of nonlinear iterations, for
421 // the solver currently being used (for two-level more than
422 // one solver is usually being invoked).
423 // Special Notes :
424 // Scope : private
425 // Creator : Eric R. Keiter, SNL, Compuational Sciences
426 // Creation Date : 10/24/02
427 //-----------------------------------------------------------------------------
429 {
430  int numIters = 0;
432  {
433  numIters = nlsOuterPtr_->getNumIterations ();
434  }
435  else
436  {
437  numIters = nlsInnerPtr_->getNumIterations ();
438  }
439 
440  return numIters;
441 }
442 
443 //-----------------------------------------------------------------------------
444 // Function : TwoLevelNewton::getMaxNormF
445 // Purpose :
446 // Special Notes :
447 // Scope : public
448 // Creator : Rich Schiek, Electrical and Memes Modeling
449 // Creation Date : 9/27/2009
450 //-----------------------------------------------------------------------------
452 {
453  double result = nlsInnerPtr_->getMaxNormF() + nlsOuterPtr_->getMaxNormF();
454  return result;
455 }
456 
457 //-----------------------------------------------------------------------------
458 // Function : TwoLevelNewton::getMaxNormFindex
459 // Purpose :
460 // Special Notes :
461 // Scope : public
462 // Creator : Rich Schiek, Electrical and Memes Modeling
463 // Creation Date : 9/27/2009
464 //-----------------------------------------------------------------------------
466 {
467  int indexInner = nlsInnerPtr_->getMaxNormFindex();
468  int indexOuter = nlsOuterPtr_->getMaxNormFindex();
469 
470  return indexInner; // usually the inner will be the one you want,
471  // but there should be more of a detailed comparison
472  // here just in case. FIX THIS
473 }
474 
475 
476 //-----------------------------------------------------------------------------
477 // Function : TwoLevelNewton::initializeAll
478 //
479 // Purpose : This serves the same purpose as the initializeAll
480 // function in the other solvers.
481 //
482 // Special Notes : Each of the solvers (inner, outer, and this one) will
483 // ultimately call the base class initialize all function,
484 // meaning that there are potentially 3 of each linear
485 // solver (3 iterative, 3 direct). This sort of makes
486 // sense for the inner and outer, as they may want
487 // different linear solvers. But, for the wrapper class,
488 // TwoLevelNewton, I'm only doing this so that it is
489 // possible to pass a solver pointer to Sensitivity,
490 // if needed.
491 //
492 // In the future, it may make sense to clean this up a
493 // litte.
494 //
495 // Scope : public
496 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
497 // Creation Date : 10/20/02
498 //-----------------------------------------------------------------------------
500 {
501  bool bsuccess = true;
502  bool tmpBool = true;
503 
504  tmpBool = nlsInnerPtr_->initializeAll ();
505  bsuccess = bsuccess && tmpBool;
506 
507  tmpBool= nlsOuterPtr_->initializeAll ();
508  bsuccess = bsuccess && tmpBool;
509 
510  tmpBool = NonLinearSolver::initializeAll();
511  bsuccess = bsuccess && tmpBool;
512 
513  savedRHSPtr_ = lasSysPtr_->builder().createVector ();
514  savedNextSolPtr_ = lasSysPtr_->builder().createVector ();
515  jdxpVectorPtr_ = lasSysPtr_->getJDXPVector ();
516 
517  // set up the return codes so that the "inner" solver is subject to
518  // greater restrictions than the outter solver.
519  ReturnCodes retCode;
520  retCode.nearConvergence = -3;
521  nlsInnerPtr_->setReturnCodes (retCode);
522 
523  return bsuccess;
524 }
525 
526 //-----------------------------------------------------------------------------
527 // Function : TwoLevelNewton::setOptions
528 //
529 // Purpose : This function processes the options set in the
530 // ".options NONLIN" line of the netlist.
531 //
532 // Special Notes : Mostly, these options are just passed on through to the
533 // "outer" nonlinear solver. However, since the outer
534 // solver is incremented one Newton step at a time, and
535 // the actual outer control loop sits here in the two level
536 // class, a few of the parameters are needed. In
537 // particular, the maximum number of Newton steps.
538 //
539 // Scope : public
540 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
541 // Creation Date : 10/20/02
542 //-----------------------------------------------------------------------------
543 bool TwoLevelNewton::setOptions(const N_UTL_OptionBlock & OB)
544 {
545  std::list<N_UTL_Param>::const_iterator it_tpL = OB.getParams().begin();
546  std::list<N_UTL_Param>::const_iterator end_tpL = OB.getParams().end();
547  for ( ; it_tpL != end_tpL; ++it_tpL)
548  {
549  if (it_tpL->uTag() == "MAXSTEP")
550  {
551  maxOuterSteps_ = static_cast<int>(it_tpL->getImmutableValue<int>());
552  }
553  }
554 
555  return nlsOuterPtr_->setOptions(OB);
556 }
557 
558 //-----------------------------------------------------------------------------
559 // Function : TwoLevelNewton::setTranOptions
560 // Purpose :
561 // Special Notes :
562 // Scope : public
563 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
564 // Creation Date : 10/20/02
565 //-----------------------------------------------------------------------------
566 bool TwoLevelNewton::setTranOptions(const N_UTL_OptionBlock & OB)
567 {
568  return nlsOuterPtr_->setTranOptions(OB);
569 }
570 
571 bool TwoLevelNewton::setHBOptions(const N_UTL_OptionBlock & OB)
572 {
573  return true;
574 }
575 
576 
577 //-----------------------------------------------------------------------------
578 // Function : TwoLevelNewton::setLocaOptions
579 // Purpose :
580 // Special Notes :
581 // Scope : public
582 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
583 // Creation Date : 09/10/03
584 //-----------------------------------------------------------------------------
585 bool TwoLevelNewton::setLocaOptions (const N_UTL_OptionBlock & OB)
586 {
587  outerLocaOptions_ = OB;
588  return nlsOuterPtr_->setLocaOptions(OB);
589 }
590 
591 //-----------------------------------------------------------------------------
592 // Function : TwoLevelNewton::setTwoLevelLocaOptions
593 // Purpose :
594 // Special Notes :
595 // Scope : public
596 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
597 // Creation Date : 09/10/03
598 //-----------------------------------------------------------------------------
599 bool TwoLevelNewton::setTwoLevelLocaOptions (const N_UTL_OptionBlock & OB)
600 {
601  innerLocaOptions_ = OB;
602  return nlsInnerPtr_->setLocaOptions(OB);
603 }
604 
605 
606 //-----------------------------------------------------------------------------
607 // Function : TwoLevelNewton::setTwoLevelOptions
608 // Purpose :
609 // Special Notes :
610 // Scope : public
611 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
612 // Creation Date : 10/20/02
613 //-----------------------------------------------------------------------------
614 bool TwoLevelNewton::setTwoLevelOptions(const N_UTL_OptionBlock & OB)
615 {
616  N_UTL_OptionBlock OBtmp;
617 
618  std::list<N_UTL_Param>::const_iterator it_tpL = OB.getParams().begin();
619  std::list<N_UTL_Param>::const_iterator end_tpL = OB.getParams().end();
620  for ( ; it_tpL != end_tpL; ++it_tpL)
621  {
622  if (it_tpL->uTag() == "ALGORITHM")
623  {
624  twoLevelAlgorithm_ = static_cast<int>(it_tpL->getImmutableValue<int>());
625  }
626  else if (it_tpL->uTag() == "NOX")
627  {
628  noxFlagInner_ = it_tpL->getImmutableValue<int>();
629  }
630  else if (it_tpL->uTag() == "MAXCONTSTEPS")
631  {
632  maxContSteps_ = static_cast<int>(it_tpL->getImmutableValue<int>());
633  }
634  else if (it_tpL->uTag() == "CONTINUATIONFLAG")
635  {
636  int tmp = static_cast<int>(it_tpL->getImmutableValue<int>());
637  continuationType_ = tmp;
638  }
639  else if (it_tpL->uTag() == "INNERFAIL")
640  {
641  int tmp = static_cast<int>(it_tpL->getImmutableValue<int>());
642  innerLoopFailFatal_ = (tmp!=0);
643  }
644  else if (it_tpL->uTag() == "EXITWITHFAILURE")
645  {
646  int tmp = static_cast<int>(it_tpL->getImmutableValue<int>());
647  totalSolveFailFatal_ = (tmp!=0);
648  }
649  else if (it_tpL->uTag() == "FULLNEWTONENFORCE")
650  {
651  int tmp = static_cast<int>(it_tpL->getImmutableValue<int>());
653  }
654  else if (it_tpL->uTag() == "CONPARAM")
655  {
656  paramNameList.push_back(it_tpL->stringValue());
657  }
658  else if (it_tpL->uTag() == "VOLTLIMTOL")
659  {
660  voltLimTol_ = static_cast<double>(it_tpL->getImmutableValue<double>());
661  }
662  else // anything that is not a "special" two-level param, push
663  // back on the the tmp param list.
664  {
665  OBtmp.getParams().push_back(*it_tpL);
666  }
667  }
668 
669  innerSolverOptions_ = OBtmp; // keep a copy of these.
670 
671  nlsInnerPtr_->setOptions(OBtmp);
672 
673  if (twoLevelAlgorithm_ < 0 || twoLevelAlgorithm_ > 5)
674  {
675  Xyce::lout() << "***** WARNING: Right now the only two-level algorithms are algorithm=0-5. Resetting to 0.\n\n" << std::endl;
676 
678  }
679 
680 #ifdef Xyce_VERBOSE_NONLINEAR
681  Xyce::dout() << "\n" << std::endl
682  << Xyce::section_divider << std::endl
683  << "\n***** 2-level Inner Loop Nonlinear solver options:\n" << std::endl
684  << "\talgorithm:\t\t\t" << twoLevelAlgorithm_ << std::endl
685  << "\toutersteps:\t\t\t" << maxOuterSteps_ << std::endl
686  << "\tmaxContSteps:\t\t\t" << maxContSteps_ << std::endl;
687 
688 #if 0
689  innerLoopParams.printParams(Xyce::dout());
690 #endif
691 
692  Xyce::dout() << "\n***** Done printing Inner Loop Params:\n" << std::endl
693  << Xyce::section_divider << std::endl
694  << "\n" << std::endl;
695 #endif
696 
697  // Now that the loop is done allocate the param val array:
698  paramFinalVal.resize(paramNameList.size(),0.0);
699  paramCurrentVal.resize(paramNameList.size(),0.0);
700 
701  return true;
702 }
703 
704 //-----------------------------------------------------------------------------
705 // Function : TwoLevelNewton::setTwoLevelTranOptions
706 // Purpose :
707 // Special Notes :
708 // Scope : public
709 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
710 // Creation Date : 10/20/02
711 //-----------------------------------------------------------------------------
712 bool TwoLevelNewton::setTwoLevelTranOptions(const N_UTL_OptionBlock & OB)
713 {
714  setupTranParamsFlag_ = true;
715 #ifdef Xyce_VERBOSE_NONLINEAR
716  Xyce::dout() << "In TwoLevelNewton::setTwoLevelTranOptions" << std::endl;
717 #endif
718 
719  N_UTL_OptionBlock OBtmp;
720 
721  std::list<N_UTL_Param>::const_iterator it_tpL = OB.getParams().begin();
722  std::list<N_UTL_Param>::const_iterator end_tpL = OB.getParams().end();
723  for ( ; it_tpL != end_tpL; ++it_tpL)
724  {
725  if (it_tpL->uTag() == "ALGORITHM")
726  {
727  twoLevelAlgorithmTran_ = static_cast<int>(it_tpL->getImmutableValue<int>());
728  }
729  else if ( it_tpL->uTag() == "MAXCONTSTEPS" )
730  {
731  maxContStepsTran_ = static_cast<int>(it_tpL->getImmutableValue<int>());
732  }
733  else // anything that is not a "special" two-level param, push
734  // back on the the tmp param list.
735  {
736  OBtmp.getParams().push_back(*it_tpL);
737  }
738  }
739 
741 
742  if (twoLevelAlgorithmTran_ < 0 || twoLevelAlgorithmTran_ > 3)
743  {
744  Xyce::lout() << "***** WARNING: Right now the only two-level algorithms are algorithm=0-3. Resetting to 0.\n\n" << std::endl;
745 
747  }
748 
749 #ifdef Xyce_VERBOSE_NONLINEAR
750 #if 0
751  Xyce::lout() << "\n" << std::endl
752  << Xyce::section_divider << std::endl
753  << "\n***** 2-level Transient Inner Loop Nonlinear solver options:\n" << std::endl;
754  innerLoopTranParams.printParams(Xuce::lout());
755 
756  Xyce::lout() << "\n***** Done printing Inner Loop Transient Params:\n" << std::endl
757  << Xyce::section_divider << std::endl
758  << "\n" << std::endl;
759 #endif
760 #endif
761 
762  return true;
763 }
764 
765 
766 //-----------------------------------------------------------------------------
767 // Function : TwoLevelNewton::setAnalysisMode
768 // Purpose : This function is slightly different than the
769 // function of the same name in NLParams.
770 //
771 // Special Notes :
772 // Scope : public
773 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
774 // Creation Date : 10/21/02
775 //-----------------------------------------------------------------------------
777 {
778  // this should be 0,1 or 2 (DC_OP, DC_SWEEP, or TRANSIENT)
779 #ifdef Xyce_VERBOSE_NONLINEAR
780  Xyce::dout() << std::endl;
781  Xyce::dout() << "Setting the externalAnalysisMode = " << mode << std::endl;
782 #endif
783  externalAnalysisMode = mode;
784 
787 }
788 
789 //-----------------------------------------------------------------------------
790 // Function : TwoLevelNewton::registerPetraOptions
791 // Purpose : see header file
792 // Special Notes : see header file
793 // Scope : public
794 // Creator : Robert Hoekstra, SNL, Parallel Computational Sciences
795 // Creation Date : 11/9/00
796 //-----------------------------------------------------------------------------
797 
798 bool TwoLevelNewton::setPetraOptions(const N_UTL_OptionBlock & OB)
799 {
800  bool bsuccess = true;
801  bool tmpBool = true;
802  tmpBool = nlsOuterPtr_->setPetraOptions(OB);
803  bsuccess = bsuccess && tmpBool;
804 
805  tmpBool = nlsInnerPtr_->setPetraOptions(OB);
806  bsuccess = bsuccess && tmpBool;
807 
808  tmpBool = NonLinearSolver::setPetraOptions(OB);
809  bsuccess = bsuccess && tmpBool;
810 
811  return bsuccess;
812 }
813 
814 
815 //-----------------------------------------------------------------------------
816 // Function : TwoLevelNewton::printStepInfo_
817 // Purpose : Print out 2-level Newton step information.
818 // Special Notes :
819 // Scope : public
820 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
821 // Creation Date : 10/21/02
822 //-----------------------------------------------------------------------------
824  (int step, int success, TwoLevelNewtonMode solveType)
825 {
826  char tmpChar[128];
827  if (solveType==FULL_PROBLEM)
828  {
829  Xyce::lout() << "---------- 2LNiter: " << step << "\t" << success << "\tFULL PROBLEM --------------------------------" << std::endl;
830  }
831  else if (solveType==INNER_PROBLEM)
832  {
833  Xyce::lout() << "---------- 2LNiter: " << step << "\t" << success << "\tINNER PROBLEM ----------------------------" << std::endl;
834  }
835  else
836  {
837  Xyce::lout() << "---------- 2LNiter: " << step << "\t" << success << "\tOUTER PROBLEM ----------------------------" << std::endl;
838  }
839 }
840 
841 //-----------------------------------------------------------------------------
842 // Function : TwoLevelNewton::zeroInnerLoopStatistics_
843 // Purpose :
844 // Special Notes :
845 // Scope : public
846 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
847 // Creation Date : 10/21/02
848 //-----------------------------------------------------------------------------
850 {
851  numResidualLoads_ = 0;
852  numJacobianLoads_ = 0;
853  numLinearSolves_ = 0;
857  totalLinearSolveTime_ = 0.0;
860 }
861 
862 //-----------------------------------------------------------------------------
863 // Function : TwoLevelNewton::calcInnerLoopStatistics_
864 // Purpose :
865 // Special Notes :
866 // Scope : public
867 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
868 // Creation Date : 10/21/02
869 //-----------------------------------------------------------------------------
871 {
881 }
882 
883 //-----------------------------------------------------------------------------
884 // Function : TwoLevelNewton::calcOuterLoopStatistics_
885 // Purpose :
886 // Special Notes :
887 // Scope : public
888 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
889 // Creation Date : 10/21/02
890 //-----------------------------------------------------------------------------
892 {
902 }
903 
904 
905 //-----------------------------------------------------------------------------
906 // Function : TwoLevelNewton::algorithm0_
907 // Purpose : This algorithm is the full newton algorithm.
908 //
909 // Special Notes : The main thing that is different about running this
910 // function, rather than just not using the 2-level class
911 // at all, is that here the conductance and capacitance
912 // calculations are performed at the end of the solve.
913 //
914 // Scope : private
915 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
916 // Creation Date : 05/07/03
917 //-----------------------------------------------------------------------------
919 {
920  int status = -1;
921 
922 #ifdef Xyce_VERBOSE_NONLINEAR
923  Xyce::dout() << std::endl << "Running algorithm 0:" << std::endl;
924 #endif
925 
926  // set up the local algorithm variable:
927  int algorithm = twoLevelAlgorithm_;
928  if ( externalAnalysisMode ==2) algorithm = twoLevelAlgorithmTran_;
929 
930  // This step is neccessary in case we've just switched algorithms.
931  if ( algorithm == 0) twoLevelCouplingMode_=FULL_PROBLEM;
932 
933  status = nlsOuterPtr_->solve ();
934 
935  // Now do conductance/capacitance extractions, if it makes sense.
936  if (!firstDCOPFlag_)
937  {
939  }
940 
941  // Do this so that I can test the ConductanceExtractor class.
943  loaderPtr_->loadJacobian ();
944 
945  return status;
946 }
947 
948 //-----------------------------------------------------------------------------
949 // Function : TwoLevelNewton::algorithm1_
950 // Purpose : This is the first two-level algorithm implemented.
951 // The idea is to have nested Newton solves.
952 //
953 // Inner loop = PDE-device only
954 // Outer loop = Full Newton.
955 //
956 // Special Notes : Doesn't work with NOX yet...
957 //
958 // This will exit immediately if the inner loop fails.
959 // The reason for this is that the point of 2-level
960 // Newton is to solve problems that are difficult
961 // as a result of two very different problem types
962 // being coupled together.
963 //
964 // Ideally, for the two-level Newton, you have a
965 // circuit problem that can be easily solved on its
966 // own, and a PDE-device problem that can easily be
967 // solved on its own, and the only problem is that
968 // they can't easily be solved together. Two level
969 // Newton separates them, and applies different
970 // (hopefully ideal) solver settings to each.
971 //
972 // If the PDE problem
973 // (the inner loop) can't be solved with its
974 // current settings, even in stand-alone mode, 2-level
975 // Newton won't help.
976 //
977 // Scope : private
978 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
979 // Creation Date : 10/22/02
980 //-----------------------------------------------------------------------------
982 {
983  int status = -1;
984 
985 #ifdef Xyce_VERBOSE_NONLINEAR
986  Xyce::dout() << std::endl << "Running algorithm 1:" << std::endl;
987 #endif
988 
989  bool firstOuterStepTaken = false;
990  nlsPassingPtr_ = 0;
991 
992  if (status < 0)
993  {
994  int twoLevelStep;
995  for (twoLevelStep=0;twoLevelStep<maxOuterSteps_; ++twoLevelStep)
996  {
997  // Device Only:
999  outerLoopActiveFlag_= false;
1000  int statInner = nlsInnerPtr_->solve (nlsPassingPtr_);
1001  nlsPassingPtr_ = 0;
1003 
1004 #ifdef Xyce_VERBOSE_NONLINEAR
1005  // Output the nonlinear solver step information:
1006  printStepInfo_(twoLevelStep+1,statInner, twoLevelCouplingMode_);
1007 #endif
1008 
1009  // If we can't even get the inner loop to converge,
1010  // just give up.
1011  if (statInner < 0)
1012  {
1013  break;
1014  }
1015 
1016  // Full Problem:
1018  outerLoopActiveFlag_= true;
1019  if (firstOuterStepTaken)
1020  {
1021  status = nlsOuterPtr_->takeOneSolveStep ();
1022  }
1023  else
1024  {
1025  firstOuterStepTaken = true;
1027  }
1028 
1030 
1031 #ifdef Xyce_VERBOSE_NONLINEAR
1032  // Output the nonlinear solver step information:
1033  printStepInfo_(twoLevelStep+1,status, twoLevelCouplingMode_);
1034 #endif
1035  // exit?
1036  // Use the total "full Newton" error to evaluate this step.
1037  if (status >= 0)
1038  {
1039  break;
1040  }
1041 
1042  } // end of for loop
1043  } // end of status if statement
1044 
1045 #ifdef Xyce_VERBOSE_NONLINEAR
1046  if (status >=0)
1047  {
1048  Xyce::dout() << "TWO LEVEL Newton succeeded!" << std::endl;
1049  }
1050 #endif
1051 
1052  return status;
1053 }
1054 
1055 //-----------------------------------------------------------------------------
1056 // Function : TwoLevelNewton::algorithm2_
1057 // Purpose : Similar to algorithm1, but this one allows for the
1058 // inner PDE Newton solve to be gradually stepped up to
1059 // the circuit imposed boundary conditions. Essentially,
1060 // the inner loop Newton solve is a continuation method.
1061 //
1062 // This could almost be called a "3-level" Newton, as there
1063 // are now 3 nested loops.
1064 //
1065 // Special Notes : not finished yet...
1066 //
1067 // Scope : private
1068 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1069 // Creation Date : 10/22/02
1070 //-----------------------------------------------------------------------------
1072 {
1073  int status = -1;
1074  int statInner = 0;
1075  bool statusFull = false;
1076  bool firstOuterStepTaken = false;
1077 
1078 #ifdef Xyce_VERBOSE_NONLINEAR
1079  Xyce::dout() << std::endl << "Running algorithm 2:" << std::endl;
1080 #endif
1081 
1082  // Start out with a single step of the full problem:
1083  // Full Problem:
1085  firstOuterStepTaken = true;
1086  outerLoopActiveFlag_= true;
1088 
1089 #ifdef Xyce_VERBOSE_NONLINEAR
1090  // Output the nonlinear solver step information:
1092 #endif
1093 
1094  nlsPassingPtr_ = 0;
1095  if (status <= 0)
1096  {
1097  int twoLevelStep;
1098  for (twoLevelStep=0;twoLevelStep<maxOuterSteps_; ++twoLevelStep)
1099  {
1100  // Device Only:
1102  outerLoopActiveFlag_= false;
1103 
1104  if (continuationType_ == 1)
1105  {
1106  statInner = continuationLoop_ ();
1107  }
1108  else if (continuationType_ == 2)
1109  {
1110  statInner = locaLoop_ ();
1111  }
1112  else
1113  {
1114  statInner = nlsInnerPtr_->solve (nlsPassingPtr_);
1115  nlsPassingPtr_ = 0;
1117  }
1118 
1119 #ifdef Xyce_VERBOSE_NONLINEAR
1120  // Output the nonlinear solver step information:
1121  printStepInfo_(twoLevelStep+1,statInner,twoLevelCouplingMode_);
1122 #endif
1123 
1124  if (innerLoopFailFatal_)
1125  {
1126  // If we can't even get the inner loop to converge,
1127  // just give up.
1128  if (statInner <= 0)
1129  {
1130  break;
1131  }
1132  }
1133 
1134  // Full Problem:
1136  outerLoopActiveFlag_= true;
1137 
1138  if (firstOuterStepTaken)
1139  {
1140  status = nlsOuterPtr_->takeOneSolveStep ();
1141  }
1142  else
1143  {
1144  firstOuterStepTaken = true;
1146  }
1148 
1149 #ifdef Xyce_VERBOSE_NONLINEAR
1150  // Output the nonlinear solver step information:
1151  printStepInfo_(twoLevelStep+1,status,twoLevelCouplingMode_);
1152 #endif
1153 
1154  // exit?
1155  // Use the total "full Newton" error to evaluate this step.
1156  // If the ckt considers itself converged at this point, do a
1157  // final reality check by taking a full Newton step.
1158  statusFull = false;
1159 
1160  if (status > 0 && statInner>0) statusFull = true;
1161 
1162  if (statusFull) break;
1163 
1164  } // end of for loop
1165  } // end of status if statement
1166 
1167 #ifdef Xyce_VERBOSE_NONLINEAR
1168  if (status >0 && statInner > 0)
1169  {
1170  Xyce::dout() << "TWO LEVEL Newton succeeded!" << std::endl;
1171  }
1172 #endif
1173 
1174  return status;
1175 }
1176 
1177 //-----------------------------------------------------------------------------
1178 // Function : TwoLevelNewton::algorithm3_
1179 //
1180 // Purpose : This is another 2-level algorithm. For this
1181 // algorithm, the de-coupling is greater than for
1182 // algorithms 1 and 2.
1183 //
1184 // The circuit Newton steps do *not* include any
1185 // PDE device elements, and never do. The PDE
1186 // devices are replaced by approximated conductances.
1187 //
1188 // Special Notes : This is the best algorithm.
1189 //
1190 // Scope : private
1191 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1192 // Creation Date : 10/25/02
1193 //-----------------------------------------------------------------------------
1195 {
1196  int status = -1;
1197  int statInner = 0;
1198  bool statusFull = false;
1199  bool firstOuterStepTaken = false;
1200 
1201 #ifdef Xyce_VERBOSE_NONLINEAR
1202  Xyce::dout() << std::endl << "Running algorithm 3:" << std::endl;
1203 #endif
1204 
1205  nlsPassingPtr_ = 0;
1206 
1207  int twoLevelStep;
1208  if (status <= 0)
1209  {
1210  for (twoLevelStep=0;twoLevelStep<maxOuterSteps_; ++twoLevelStep)
1211  {
1212  // Device Only:
1213  // Solve device problem with continuation:
1215  outerLoopActiveFlag_= false;
1216  statInner = 0;
1217 
1218  if (continuationType_ == 1)
1219  {
1220  statInner = continuationLoop_ ();
1221  }
1222  else if (continuationType_ == 2)
1223  {
1224  statInner = locaLoop_ ();
1225  }
1226  else
1227  {
1228  statInner = nlsInnerPtr_->solve (nlsPassingPtr_);
1229  nlsPassingPtr_ = 0;
1231  }
1232 
1233 #ifdef Xyce_VERBOSE_NONLINEAR
1234  // Output the nonlinear solver step information:
1235  printStepInfo_(twoLevelStep+1,statInner,twoLevelCouplingMode_);
1236 #endif
1237 
1238  if (innerLoopFailFatal_)
1239  {
1240  // If we can't even get the inner loop to converge,
1241  // just give up.
1242  if (statInner <= 0)
1243  {
1244  break;
1245  }
1246  }
1247 
1248  // Now do all the conductance extractions:
1249  calcCouplingTerms_ ();
1250 
1251  // ckt-only problem:
1253  outerLoopActiveFlag_= true;
1254  if (firstOuterStepTaken)
1255  {
1256  status = nlsOuterPtr_->takeOneSolveStep ();
1257  }
1258  else
1259  {
1260  firstOuterStepTaken = true;
1262  }
1263 
1264  if (noxFlag_)
1265  {
1267  }
1268 
1269 #ifdef Xyce_VERBOSE_NONLINEAR
1270  // Output the nonlinear solver step information:
1271  printStepInfo_(twoLevelStep+1,status,twoLevelCouplingMode_);
1272 #endif
1273 
1274  // check if voltage limiting is still active by getting the norm of the
1275  // jdxp vector:
1276  double twoNormJDXP_ = 0.0;
1277  jdxpVectorPtr_->lpNorm(2, &twoNormJDXP_);
1278 
1279  bool voltLimStat = (twoNormJDXP_ <= voltLimTol_);
1280 
1281 #ifdef Xyce_VERBOSE_NONLINEAR
1282  Xyce::dout() << std::endl;
1283  Xyce::dout() << " 2-norm of voltage limiting vector: " << twoNormJDXP_ << std::endl;
1284 #endif
1285 
1286  // exit?
1287  statusFull = false;
1288  if (status > 0 && statInner>0 && voltLimStat) statusFull = true;
1289 
1290  if (statusFull) break;
1291 
1292  } // end of outer steps "for" loop
1293 
1294  } // end of status if statement
1295 
1296  // Do a final few "full Newton" steps, to get a final consistency between
1297  // the two solvers.
1298  int statFinal = 2;
1299  if (doFullNewtonFinalEnforcement_ && statusFull)
1300  {
1302  statFinal = nlsOuterPtr_->solve ();
1303 #ifdef Xyce_VERBOSE_NONLINEAR
1304  // Output the nonlinear solver step information:
1305  ++twoLevelStep;
1306  printStepInfo_(twoLevelStep+1,status,twoLevelCouplingMode_);
1307 #endif
1308  }
1309 
1310  // Do one final conductance/capacitance extraction:
1311  calcCouplingTerms_ ();
1312 
1313 #ifdef Xyce_VERBOSE_NONLINEAR
1314  if (status >0 && statInner > 0 && statFinal > 0)
1315  {
1316  Xyce::dout() << "TWO LEVEL Newton succeeded!" << std::endl;
1317  }
1318 #endif
1319 
1320  return status;
1321 }
1322 
1323 //-----------------------------------------------------------------------------
1324 // Function : TwoLevelNewton::algorithm4_
1325 //
1326 // Purpose : This is mostly a test algorithm.
1327 //
1328 // This is not a two-level class algorithm, but rather a
1329 // continuation algorithm to test out
1330 // some of the new continuation capabilities.
1331 //
1332 // Mostly, I'm putting this here to test out the new
1333 // "setParam" function in the device package, which is
1334 // needed by LOCA.
1335 //
1336 // Special Notes :
1337 // Scope : private
1338 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1339 // Creation Date : 07/26/03
1340 //-----------------------------------------------------------------------------
1342 {
1343  int status = -1;
1344  int statNL;
1345  nlsPassingPtr_ = 0;
1346  bool continuationLoop = false;
1347  bool successBool;
1348  int contMaxTmp = 100; // this is just a guess for number of steps.
1349  int stepsLeft = contMaxTmp;
1350 
1351 #ifdef Xyce_VERBOSE_NONLINEAR
1352  Xyce::dout() << std::endl << "Running algorithm 4:" << std::endl;
1353  Xyce::dout() << std::endl << "Initial continuation steps: ";
1354  Xyce::dout() << contMaxTmp << std::endl;
1355 #endif
1356 
1357  // The solves should be FULL newton solves.
1359 
1360  // These values assume that any parameter will start at zero and ramp up
1361  // to a final value. The final value is the value that is set for this
1362  // parameter in the netlist. So, for example, if a voltage source is set
1363  // to 5.0 volts in the netlist:
1364  //
1365  // vid 4 0 5.00
1366  //
1367  // and one of the continuation parameters is vid:v0, then the final value
1368  // for vid:v0 is 5.00.
1369 
1370  double initVal = 0.0;
1371  double currVal = initVal;
1372  double prevVal = initVal;
1373  double stepSizeEst;
1374 
1375  std::vector<std::string>::iterator iter;
1376  std::vector<std::string>::iterator begin = paramNameList.begin ();
1377  std::vector<std::string>::iterator end = paramNameList.end ();
1378 
1379  std::vector<double>::iterator iterFinalVal;
1380  std::vector<double>::iterator beginFinalVal = paramFinalVal.begin ();
1381  std::vector<double>::iterator endFinalVal = paramFinalVal.end ();
1382 
1383  std::vector<double>::iterator iterCurrentVal;
1384  std::vector<double>::iterator beginCurrentVal = paramCurrentVal.begin ();
1385  std::vector<double>::iterator endCurrentVal = paramCurrentVal.end ();
1386 
1387  // Get the final values for each param., then set each to zero.
1388  for (iter=begin, iterFinalVal=beginFinalVal, iterCurrentVal=beginCurrentVal;
1389  iter!=end;
1390  ++iter, ++iterFinalVal, ++iterCurrentVal)
1391  {
1392  //*iterFinalVal = loaderPtr_->getParam(*iter);
1393  *iterFinalVal = 1.0;
1394  loaderPtr_->setParam(*iter, 0.0);
1395  }
1396 
1397  // now, loop over each parameter, and do a continuation loop for each.
1398  // For each param, the loop will adjust the param from 0.0 to the final
1399  // value. Once each parameter has achieved a final value, the solve is
1400  // finished.
1401  for (iter=begin, iterFinalVal=beginFinalVal;
1402  iter!=end;
1403  ++iter, ++iterFinalVal)
1404  {
1405  // get the initial stepsize:
1406  stepSizeEst = (*iterFinalVal-0.0)/(static_cast<double>(contMaxTmp));
1407  currVal = 0.0;
1408  prevVal = 0.0;
1409  contStep_=1;
1410 
1411 #ifdef Xyce_VERBOSE_NONLINEAR
1412  Xyce::dout() << "Parameter = " << *iter;
1413  Xyce::dout() << " finalVal = " << *iterFinalVal << std::endl;
1414 #endif
1415 
1416  // begin continuation Loop for the current parameter, *iter:
1417  bool continuationLoopFinished = false;
1418  int numTotalFailures = 0;
1419  while (!continuationLoopFinished)
1420  {
1421  bool stepFinished = false;
1422  int numFailures = 0;
1423 
1424  while(!stepFinished)
1425  {
1426  if (stepSizeEst != 0.0)
1427  {
1428  stepsLeft = static_cast<int>((*iterFinalVal-currVal)/stepSizeEst) + 1;
1429  }
1430  else
1431  {
1432  stepsLeft = 1;
1433  }
1434 
1435 #ifdef Xyce_VERBOSE_NONLINEAR
1436  Xyce::dout() << std::endl << "Continuation Step: " << contStep_;
1437  Xyce::dout() << " Estimated Remaining Steps: " << stepsLeft;
1438  Xyce::dout() << " " << *iter;
1439  Xyce::dout() << std::endl;
1440  Xyce::dout() << "currVal= " << currVal;
1441  Xyce::dout() << " prevVal= " << prevVal;
1442  Xyce::dout() << " step= " << stepSizeEst;
1443  Xyce::dout() << std::endl;
1444 #endif
1445 
1446  if (stepsLeft < 0)
1447  {
1448  std::string tmp = "Continuation step estimate broken. Exiting\n";
1449  N_ERH_ErrorMgr::report(N_ERH_ErrorMgr::USR_FATAL_0, tmp);
1450  }
1451 
1452  // save a copy of solution:
1453  (*savedNextSolPtr_) = (**nextSolVectorPtrPtr_);
1454 
1455  // set the continuation-parameter:
1456  loaderPtr_->setParam (*iter, currVal);
1457  *iterCurrentVal = currVal;
1458 
1459  // perform the nonlinear solve:
1460  statNL = nlsOuterPtr_->solve (nlsPassingPtr_);
1461 
1462  int stepsTaken = nlsOuterPtr_->getNumIterations ();
1463  int stepsNotTaken = abs(maxOuterSteps_ - stepsTaken);
1464 
1465  nlsPassingPtr_ = 0;
1466 
1467  // Add to the stats:
1469 
1470  // Earlier (in initializeAll), return codes should have been set,
1471  // hopefully correctly...
1472  successBool = (statNL > 0);
1473 
1474  if (successBool) // success!
1475  {
1476  stepFinished = true;
1477 
1478  //if (numFailures <= 0 && stepsNotTaken > 7)
1479  if (numFailures <= 0)
1480  {
1481  stepSizeEst *= increaseContScalar_;
1482  }
1483 
1484  --numFailures;
1485  if (numFailures < 0) numFailures = 0;
1486 
1487  prevVal = currVal;
1488  currVal += stepSizeEst;
1489 
1490  if ( (*iterFinalVal >= 0 && currVal > *iterFinalVal) ||
1491  (*iterFinalVal < 0 && currVal < *iterFinalVal) )
1492  {
1493  currVal = *iterFinalVal;
1494  stepSizeEst = currVal - prevVal;
1495  }
1496 
1497 #ifdef Xyce_DEBUG_NONLINEAR
1498  Xyce::dout() << "\nRight before outputHOMOTOPY:" << std::endl;
1499  for (int ieric=0;ieric<paramNameList.size();++ieric)
1500  {
1501  Xyce::dout() << paramNameList[ieric] << "\t";
1502  Xyce::dout() << paramCurrentVal[ieric] << std::endl;
1503  }
1504 #endif
1505 
1506  outMgrPtr_->outputHomotopy(pdsMgrPtr_->getPDSComm()->comm(), paramNameList, paramCurrentVal, **nextSolVectorPtrPtr_);
1507  }
1508  else // failure!
1509  {
1510  stepSizeEst *= decreaseContScalar_;
1511 
1512  // restore the solution:
1513  **nextSolVectorPtrPtr_ = (*savedNextSolPtr_);
1514 
1515  ++numFailures;
1516  ++numTotalFailures;
1517 
1518  currVal = prevVal + stepSizeEst;
1519  }
1520  } // end of stepFinished while loop
1521 
1522  if (!successBool) // failure...
1523  {
1524  break;
1525  }
1526 
1527  ++contStep_;
1528 
1529  continuationLoopFinished =
1530  ( (*iterFinalVal >= 0 && prevVal >= *iterFinalVal) ||
1531  (*iterFinalVal < 0 && prevVal <= *iterFinalVal) );
1532 
1533  } // end of continuation loop.
1534 
1535 #ifdef Xyce_VERBOSE_NONLINEAR
1536  Xyce::dout() << "currVal= " << currVal;
1537  Xyce::dout() << " prevVal= " << prevVal;
1538  Xyce::dout() << std::endl;
1539  Xyce::dout() << std::endl;
1540  Xyce::dout() << "Total number of failures = " << numTotalFailures << std::endl;
1541  Xyce::dout() << "Number of actual steps = " << contStep_-1;
1542  Xyce::dout() << std::endl;
1543 #endif
1544 
1545  } // end of parameter loop.
1546 
1547  status = statNL;
1548 
1549  return status;
1550 }
1551 
1552 //-----------------------------------------------------------------------------
1553 // Function : TwoLevelNewton::algorithm5_
1554 //
1555 // Purpose : Same as algorithm0, but calls the inner solver rather
1556 // than the outer.
1557 //
1558 // Special Notes : This is best used if you want to use LOCA most of the
1559 // time, but not on the nonlinear poisson solve
1560 // (firstDCOPstep).
1561 //
1562 // Scope : private
1563 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1564 // Creation Date : 05/07/03
1565 //-----------------------------------------------------------------------------
1567 {
1568  int status = -1;
1569 
1570 #ifdef Xyce_VERBOSE_NONLINEAR
1571  Xyce::dout() << std::endl << "Running algorithm 5:" << std::endl;
1572 #endif
1573 
1574  // set up the local algorithm variable:
1575  int algorithm = twoLevelAlgorithm_;
1576  if ( externalAnalysisMode ==2) algorithm = twoLevelAlgorithmTran_;
1577 
1579 
1580  status = nlsInnerPtr_->solve ();
1581 
1582  // Now do conductance/capacitance extractions, if it makes sense.
1583  //if (!firstDCOPFlag_)
1584  //{
1585  //calcCouplingTerms_ ();
1586  //}
1587 
1588  return status;
1589 }
1590 
1591 
1592 //-----------------------------------------------------------------------------
1593 // Function : TwoLevelNewton::solve
1594 // Purpose :
1595 // Special Notes : Doesn't work with NOX yet...
1596 // Scope : public
1597 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1598 // Creation Date : 10/20/02
1599 //-----------------------------------------------------------------------------
1601 {
1602  int status = -1;
1603 
1604  firstDCOPFlag_ = false;
1605 
1606  // Is this the first step of a double DCOP step?
1607  // If so, then just do a simple full Newton solve, no matter
1608  // what algorithm is specified.
1609  if (tiaPtr_ == 0)
1610  {
1611  Xyce::dout() << std::endl;
1612  Xyce::dout() << "tiaPtr is null. exiting" << std::endl; exit(0);
1613  }
1614 
1615  N_TIA_TimeIntInfo tiInfo;
1616  tiaPtr_->getTimeIntInfo(tiInfo);
1617  bool doubleDCOPEnable = tiInfo.doubleDCOPEnabled;
1618  int tiaMode = tiInfo.timeIntMode;
1619  int ddcopStep = tiInfo.doubleDCOPStep;
1620 
1621  if (doubleDCOPEnable && (tiaMode == 0) && (ddcopStep==0) )
1622  firstDCOPFlag_ = true;
1623 
1625 
1626  // set up the local algorithm variable:
1627  int algorithm = twoLevelAlgorithm_;
1628  if ( externalAnalysisMode ==2) algorithm = twoLevelAlgorithmTran_;
1629 
1630  // Some intitial setup:
1632  {
1633  numInterfaceNodesSetup_ = true;
1634  loaderPtr_->getNumInterfaceNodes (numInterfaceNodes_);
1636 
1637 #ifdef Xyce_VERBOSE_NONLINEAR
1638  Xyce::dout() << std::endl;
1639  Xyce::dout() << "numSubProblems_ = " << numSubProblems_ << std::endl;
1640 #endif
1641  }
1642 
1643  // Algorithm 0:
1644  // Full Newton, same as if two level was never called.
1645  if (algorithm == 0 || firstDCOPFlag_)
1646  {
1647  status = algorithm0_();
1648  }
1649  // else if algorithm 1:
1650  // outter loop is full Newton, inner loop PDE device only.
1651  else if (algorithm == 1)
1652  {
1653  status = algorithm1_ ();
1654  }
1655  // same as algorithm 1, but with continuation applied to the inner loop.
1656  else if (algorithm == 2)
1657  {
1658  status = algorithm2_ ();
1659  }
1660  // outter loop is circuit only, inner loop PDE device only, with
1661  // continuation.
1662  else if (algorithm == 3)
1663  {
1664  status = algorithm3_ ();
1665  }
1666  else if (algorithm == 4)
1667  {
1668  status = algorithm4_ ();
1669  }
1670  else if (algorithm == 5)
1671  {
1672  status = algorithm5_ ();
1673  }
1674  else
1675  {
1676  std::string tmp =
1677  "Two-Level Newton Algorithm set to invalid number.\n";
1678  N_ERH_ErrorMgr::report(N_ERH_ErrorMgr::USR_FATAL_0, tmp);
1679  }
1680 
1681  // Sometimes, when trying to debug, you want the code to exit the first
1682  // time the solver bombs. The totalSolveFailFlag gives that option.
1683  // Right before exiting, however, it instructs the device package to
1684  // output all its stuff. (tecplot files, text files, etc.)
1685  if (totalSolveFailFatal_ && status<=0)
1686  {
1687 #ifndef Xyce_PARALLEL_MPI
1688  loaderPtr_->output();
1689 #endif
1690 
1691  std::string tmp =
1692  "Two-Level Newton Algorithm failed to converge. Exiting.\n";
1693  N_ERH_ErrorMgr::report(N_ERH_ErrorMgr::USR_FATAL_0, tmp);
1694  }
1695 
1696  return status;
1697 }
1698 
1699 //-----------------------------------------------------------------------------
1700 // Function : TwoLevelNewton::calcCouplingTerms_
1701 //
1702 // Purpose : The purpose of this function is to manage the extraction
1703 // of conductances (coupling terms) from each PDE device.
1704 // These conductances may be neccessary if the outter Newton
1705 // loop is to be a ckt-only loop. In that case, the
1706 // PDE devices are replaced by much smaller conductance
1707 // based models.
1708 //
1709 // Special Notes : Performing a conductance calculation on a PDE device
1710 // requires a chain rule calculation, which includes several
1711 // linear solves of the Jacobian. It is because of the
1712 // required linear solve that some of the work is done
1713 // up here in the nonlinear solver, instead of having all
1714 // of it be done by the individual PDE devices.
1715 //
1716 // There may be a better solution later, but for the time
1717 // being I didn't want linear solves to be performed down
1718 // in the device package.
1719 //
1720 // Scope : public
1721 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1722 // Creation Date : 12/03/02
1723 //-----------------------------------------------------------------------------
1725 {
1726  bool bsuccess = true;
1727  bool tmpBool = true;
1728 
1729  char filename1[256];
1730 
1731  for (int ich = 0; ich < 256; ++ich)
1732  { filename1[ich] = 0; }
1733 
1734  // If the current mode is not "INNER_PROBLEM", then we need to re-load
1735  // the matrix, as the ckt part of the matrix needs to be kept out of
1736  // this calculation. Also, re-load Jacobian to make sure it is
1737  // completely up-to-date with the current solution.
1740  {
1742  }
1743  loaderPtr_->loadJacobian ();
1744 
1745  sprintf(filename1,"%s","tmpJac.txt");
1746  N_LAS_Matrix *A = lasSysPtr_->getJacobianMatrix();
1747  A->writeToFile(filename1);
1748 
1749  // save a copy of the RHS, as it is going to get mangled.
1750  N_LAS_Vector *rhsVecPtr = lasSysPtr_->getRHSVector();
1751  N_LAS_Vector *newtVecPtr = nlsOuterPtr_->NewtonVectorPtr_;
1752  *savedRHSPtr_ = *rhsVecPtr;
1753 
1754  // Loop over each sub-problem (PDE device).
1755  // Within each sub-problem, loop over each coupling term. (electrode)
1756  for (int iSubProblem=0; iSubProblem<numSubProblems_; ++iSubProblem)
1757  {
1758  int iCouple;
1759  int numCoupleTerms = numInterfaceNodes_[iSubProblem];
1760 
1761 #ifdef Xyce_VERBOSE_NONLINEAR
1762  Xyce::dout() << "\n numCoupleTerms = " << numCoupleTerms << std::endl;
1763 #endif
1764 
1765  for (iCouple=0;iCouple<numCoupleTerms;++iCouple)
1766  {
1767  // first zero out the RHS vector
1768  rhsVecPtr->putScalar(0.0);
1769 
1770  // load RHS vector with dFdV.
1771  tmpBool = loaderPtr_->loadCouplingRHS (iSubProblem, iCouple, rhsVecPtr );
1772  bsuccess = bsuccess && tmpBool;
1773 
1774  sprintf(filename1,"dfdv%02d.txt", iCouple);
1775  rhsVecPtr->writeToFile(filename1);
1776 
1777  // solve linear system to get dVdX.
1778  tmpBool = nlsOuterPtr_->newton_();
1779  bsuccess = bsuccess && tmpBool;
1780  numLinearSolves_ += 1;
1781 
1782  sprintf(filename1,"dvdx%02d.txt", iCouple);
1783  newtVecPtr->writeToFile(filename1);
1784 
1785  // copy the newton vector (result of linear solve) into the
1786  // RHSvector. Doing this b/c device package knows about rhs, but
1787  // doesn't know about newton.
1788  *rhsVecPtr = *newtVecPtr;
1789 
1790  // instruct each PDE device to finish the conductance calc.
1791  tmpBool = loaderPtr_->calcCouplingTerms (iSubProblem, iCouple, rhsVecPtr);
1792  bsuccess = bsuccess && tmpBool;
1793  }
1794  }
1795 
1796  // now restore everything, just in case:
1797  *rhsVecPtr = *savedRHSPtr_;
1798  twoLevelCouplingMode_ = savedMode;
1799 
1800  return bsuccess;
1801 }
1802 
1803 //-----------------------------------------------------------------------------
1804 // Function : TwoLevelNewton::continuationLoop_
1805 //
1806 // Purpose : Boundary condition continuation loop:
1807 //
1808 // Special Notes : This continuation loop allows for variable step size.
1809 //
1810 // Scope : public
1811 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1812 // Creation Date : 12/06/02
1813 //-----------------------------------------------------------------------------
1815 {
1816  int statInner;
1817  bool successBool;
1818  int contMaxTmp;
1819 
1820  // Instruct the PDE devices to set up the incremental
1821  // boundary conditions
1822  int suggestedSteps = 10; // this is a guess to the number of steps.
1823  suggestedSteps = loaderPtr_->enablePDEContinuation ();
1824 
1825  if (suggestedSteps < 1) suggestedSteps = 1;
1826  contMaxTmp = suggestedSteps;
1827 
1828  double stepSizeEst = 1.0/(static_cast<double>(contMaxTmp));
1829  double currentAlpha = 0.0;
1830  double previousAlpha = 0.0;
1831  int stepsLeft = contMaxTmp;
1832 
1833  // If the continuation loop has never been called before, then leave the
1834  // initial alpha at zero, because we have no solution to start from,
1835  // probably. If it has been called before, go ahead and take the first
1836  // step, as it should just be one step beyond an already obtained
1837  // solution.
1839  {
1840  currentAlpha = stepSizeEst;
1841  }
1843 
1844  contStep_=1;
1845  bool continuationLoopFinished = false;
1846 
1847  int numTotalFailures = 0;
1848 
1849  while (!continuationLoopFinished)
1850  {
1851  bool stepFinished = false;
1852  int numFailures = 0;
1853 
1854  while(!stepFinished)
1855  {
1856 
1857  stepsLeft = static_cast<int>((1.0-currentAlpha)/stepSizeEst) + 1;
1858 
1859 #ifdef Xyce_VERBOSE_NONLINEAR
1860  Xyce::dout() << std::endl << "Continuation Step: " << contStep_;
1861  Xyce::dout() << " Estimated Remaining Steps: " << stepsLeft;
1862  Xyce::dout() << std::endl;
1863  Xyce::dout() << "current alpha = " << currentAlpha;
1864  Xyce::dout() << " prev. alpha = " << previousAlpha;
1865  Xyce::dout() << " step = " << stepSizeEst;
1866  Xyce::dout() << std::endl;
1867 #endif
1868  if (stepsLeft < 0)
1869  {
1870  std::string tmp = "Continuation step estimate broken. Exiting\n";
1871  N_ERH_ErrorMgr::report(N_ERH_ErrorMgr::USR_FATAL_0, tmp);
1872  }
1873 
1874  // save a copy of solution:
1875  (*savedNextSolPtr_) = (**nextSolVectorPtrPtr_);
1876 
1877  std::string paramName = "pdealpha";
1878  loaderPtr_->setParam (paramName, currentAlpha);
1879 
1880  // perform the nonlinear solve:
1881  statInner = nlsInnerPtr_->solve (nlsPassingPtr_);
1882  nlsPassingPtr_ = 0;
1884 
1885  // Earlier (in initializeAll), return codes for the inner solver
1886  // should have been set so that "nearConvergence" is not considered
1887  // an adequate success. In other words, it returns a negative
1888  // number, not positive.
1889 
1890 #ifdef Xyce_DEBUG_NONLINEAR
1891  Xyce::dout() << "Status of inner loop solve: " << statInner << std::endl;
1892 #endif
1893  successBool = (statInner > 0);
1894 
1895  if (successBool) // success!
1896  {
1897  stepFinished = true;
1898 
1899  if (numFailures <= 0)
1900  {
1901  stepSizeEst *= increaseContScalar_;
1902  }
1903 
1904  --numFailures;
1905  if (numFailures < 0) numFailures = 0;
1906 
1907  previousAlpha = currentAlpha;
1908  currentAlpha += stepSizeEst;
1909 
1910  if (currentAlpha > 1.0)
1911  {
1912  currentAlpha = 1.0;
1913  stepSizeEst = currentAlpha - previousAlpha;
1914  }
1915  }
1916  else // failure!
1917  {
1918  stepSizeEst *= decreaseContScalar_;
1919 
1920  // restore the solution:
1921  (**nextSolVectorPtrPtr_) = (*savedNextSolPtr_);
1922 
1923  ++numFailures;
1924  ++numTotalFailures;
1925 
1926  currentAlpha = previousAlpha + stepSizeEst;
1927  }
1928  } // end of stepFinished while loop
1929 
1930  if (!successBool) // failure...
1931  {
1932  break;
1933  }
1934 
1935  ++contStep_;
1936  continuationLoopFinished = (previousAlpha >= 1.0);
1937 
1938  } // end of continuation loop.
1939 
1940 #ifdef Xyce_VERBOSE_NONLINEAR
1941  Xyce::dout() << "current alpha = " << currentAlpha;
1942  Xyce::dout() << " previous alpha = " << previousAlpha;
1943  Xyce::dout() << std::endl;
1944  Xyce::dout() << std::endl;
1945  Xyce::dout() << "Total number of failures = " << numTotalFailures << std::endl;
1946  Xyce::dout() << "Number of actual steps = " << contStep_-1;
1947  Xyce::dout() << std::endl;
1948 #endif
1949 
1950  loaderPtr_->disablePDEContinuation ();
1951  contStep_=0;
1952 
1953  return statInner;
1954 }
1955 
1956 //-----------------------------------------------------------------------------
1957 // Function : TwoLevelNewton::locaLoop_
1958 //
1959 // Purpose : Boundary condition continuation (using LOCA) loop:
1960 //
1961 // Special Notes : This continuation loop allows for variable step size.
1962 //
1963 // Scope : public
1964 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1965 // Creation Date : 2/26/04
1966 //-----------------------------------------------------------------------------
1968 {
1969  int statInner;
1970 
1971  int suggestedSteps = loaderPtr_->enablePDEContinuation ();
1972 
1973  // this ifdef'd block of code was an attempt to adjust loca parameters
1974  // in-situ. If the voltage change is really, really small, then there
1975  // isn't much point in doing a continuation, so I wanted to have to
1976  // option to on-the-fly switch to straight Newton.
1977  //
1978  // Unfortunately, it doesn't work at the moment - the second time it gets
1979  // called, everything goes haywire. ERK 2/25/04.
1980 #if 0
1981  Xyce::dout() << "suggested steps are: " << suggestedSteps << std::endl;
1982 
1983  std::list<N_UTL_Param>::iterator it_tpL;
1984  for (it_tpL = innerSolverOptions_.params.begin();
1985  it_tpL != innerSolverOptions_.params.end();
1986  ++it_tpL)
1987  {
1988  ExtendedString tmpTag = it_tpL->tag ();
1989  tmpTag.toUpper ();
1990 
1991  Xyce::dout() << "tmpTag = " << tmpTag << std::endl;
1992 
1993  if (tmpTag == "CONTINUATION")
1994  {
1995  if (suggestedSteps <= 1) // then just do one newton step.
1996  {
1997  Xyce::dout() << "Setting the solver type to 0" << std::endl;
1998  suggestedSteps = 1;
1999  it_tpL->setVal(0);
2000  }
2001  else
2002  {
2003  it_tpL->setVal(1);
2004  }
2005  }
2006  }
2007 
2010 #endif
2011  statInner = nlsInnerPtr_->solve (nlsPassingPtr_);
2012 
2013  nlsPassingPtr_ = 0;
2015 
2016  loaderPtr_->disablePDEContinuation ();
2017 
2018  return statInner;
2019 }
2020 
2021 //-----------------------------------------------------------------------------
2022 // Function : TwoLevelNewton::enableSensitivity
2023 //
2024 // Purpose : This re-sets the code for a sensitivity calculation.
2025 // Mainly, it loads the jacobian and rhs vectors in FULL
2026 // newton mode, if neccessary.
2027 //
2028 // Special Notes : Since the implementation of the "fullNewtonEnforce"
2029 // parameter, this has become less crucial. When this option
2030 // is enabled (the default), then the last loads prior to this
2031 // function being called were full newton loads, so there is
2032 // no need to re-do them.
2033 //
2034 // Scope : public
2035 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
2036 // Creation Date : 04/21/03
2037 //-----------------------------------------------------------------------------
2039 {
2040 #ifdef Xyce_VERBOSE_NONLINEAR
2041  Xyce::dout() << std::endl;
2042  Xyce::dout() << Xyce::section_divider << std::endl;
2043  Xyce::dout() << "TwoLevelNewton::enableSensitivity " << std::endl;
2044 #endif
2045  bool bsuccess = true;
2046  bool tmpBool = true;
2047 
2048  // if the last rhs/jacobian load was not done with FULL_PROBLEM, then these
2049  // things need to be re-loaded.
2050  //if (twoLevelCouplingMode_ != FULL_PROBLEM)
2051  //{
2053  tmpBool = NonLinearSolver::rhs_ (); bsuccess = bsuccess && tmpBool;
2054  tmpBool = NonLinearSolver::jacobian_ (); bsuccess = bsuccess && tmpBool;
2055 
2056 #ifdef Xyce_VERBOSE_NONLINEAR
2057  // print out the norm info for this "full newton" residual:
2058  double maxNormRHS_=0, twoNormRHS_ = 0.0;
2059  rhsVectorPtr_->infNorm(&maxNormRHS_);
2060  rhsVectorPtr_->lpNorm(2, &twoNormRHS_);
2061  Xyce::dout().width(21); Xyce::dout().precision(13); Xyce::dout().setf(std::ios::scientific);
2062  Xyce::dout() << std::endl;
2063  Xyce::dout() << "Max. norm of full Newton RHS: " << maxNormRHS_ << std::endl;
2064  Xyce::dout() << " 2-norm of full Newton RHS: " << twoNormRHS_ << std::endl;
2065  Xyce::dout() << Xyce::section_divider << std::endl;
2066 #endif
2067  //}
2068 
2069 #ifdef Xyce_DEBUG_NONLINEAR
2070  static int callsSens = 0;
2071  char filename1[256]; for (int ich = 0; ich < 256; ++ich) filename1[ich] = 0;
2072  char filename2[256]; for (int ich = 0; ich < 256; ++ich) filename2[ich] = 0;
2073 
2074  sprintf(filename1, "matrixTmp%d.txt",callsSens);
2075  N_LAS_Matrix *A = lasSysPtr_->getJacobianMatrix();
2076  A->writeToFile(filename1);
2077 
2078 
2079  N_LAS_Vector *b = lasSysPtr_->getRHSVector();
2080  sprintf(filename2, "rhsTmp%d.txt", callsSens);
2081  int size, i;
2082  FILE *fp1;
2083 
2084 #ifndef Xyce_PARALLEL_MPI
2085  fp1 = fopen(filename2,"w");
2086  size = lasSysPtr_->getSolutionSize();
2087  for (i=0;i<size;++i)
2088  {
2089  double output = b->getElementByGlobalIndex(i);
2090  fprintf(fp1,"%25.18e\n",output);
2091  }
2092  fclose(fp1);
2093 #endif
2094 
2095  N_LAS_Vector *x = (*nextSolVectorPtrPtr_);
2096  sprintf(filename2, "solTmp%d.txt", callsSens);
2097 
2098 #ifndef Xyce_PARALLEL_MPI
2099  fp1 = fopen(filename2,"w");
2100  size = lasSysPtr_->getSolutionSize();
2101  for (i=0;i<size;++i)
2102  {
2103  double output = x->getElementByGlobalIndex(i);
2104  fprintf(fp1,"%25.18e\n",output);
2105  }
2106  fclose(fp1);
2107 #endif
2108 
2109  ++callsSens;
2110 #endif // if 1
2111 
2112  return bsuccess;
2113 }
2114 
2115 //-----------------------------------------------------------------------------
2116 // Dummy function since homotopy doesn't work with 2-level. (for the outer loop)
2117 //-----------------------------------------------------------------------------
2119 {
2120  return true;
2121 }
2122 
2123 //-----------------------------------------------------------------------------
2124 // Dummy function ... for now. This probably needs to get revamped.
2125 //-----------------------------------------------------------------------------
2127 {
2128  return 0;
2129 }
2130 
2131 //-----------------------------------------------------------------------------
2132 // Dummy function ... for now. This probably needs to get revamped.
2133 //-----------------------------------------------------------------------------
2135 {
2136  return 0;
2137 }
2138 
2139 } // namespace Nonlinear
2140 } // namespace Xyce