Xyce  6.1
N_TIA_TimeIntegrationMethods.h
Go to the documentation of this file.
1 //-----------------------------------------------------------------------------
2 // Copyright Notice
3 //
4 // Copyright 2002 Sandia Corporation. Under the terms
5 // of Contract DE-AC04-94AL85000 with Sandia Corporation, the U.S.
6 // Government retains certain rights in this software.
7 //
8 // Xyce(TM) Parallel Electrical Simulator
9 // Copyright (C) 2002-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_TIA_TimeIntegrationMethods.h,v $
27 //
28 // Purpose : This file defines the classes for the time integration
29 // methods --- the "interface base class" along with the
30 // accompanying integration methods classes which can be
31 // used in the time integration algorithm.
32 //
33 // Special Notes :
34 //
35 // Creator : Buddy Watts
36 //
37 // Creation Date : 6/1/00
38 //
39 // Revision Information:
40 // ---------------------
41 //
42 // Revision Number: $Revision: 1.74.2.1 $
43 //
44 // Revision Date : $Date: 2015/04/02 18:20:18 $
45 //
46 // Current Owner : $Author: tvrusso $
47 //-----------------------------------------------------------------------------
48 
49 #ifndef Xyce_N_TIA_TimeIntegrationMethods_h
50 #define Xyce_N_TIA_TimeIntegrationMethods_h
51 
52 // ---------- Standard Includes ----------
53 #include <list>
54 
55 // ---------- Xyce Includes ----------
56 #include <N_ANP_fwd.h>
57 #include <N_IO_fwd.h>
58 #include <N_TIA_fwd.h>
59 #include <N_LAS_fwd.h>
60 
61 #include <N_UTL_IndexPair.h>
62 
63 namespace Xyce {
64 namespace TimeIntg {
65 
66 //-----------------------------------------------------------------------------
67 // Class : TimeIntegrationMethod
68 //
69 // Purpose : This is the integration methods interface class, from which
70 // specific integration methods (such as BDF15, trap, etc) are
71 // derrived.
72 //
73 // Special Notes :
74 // Creator : Buddy Watts, SNL
75 // Creation Date : 6/01/00
76 //-----------------------------------------------------------------------------
78 {
79 protected:
81  {}
82 
83 public:
85  {}
86 
87  virtual const char *getName() const = 0;
88 
89  // Predict solution at next time point.
90  virtual void obtainPredictor() = 0;
91 
92  // Evaluate the predictor derivative formula.
93  virtual void obtainPredictorDeriv() = 0;
94 
95  // Evaluate the corrector derivative formula.
96  virtual void obtainCorrectorDeriv() = 0;
97 
98  // Compute an estimate of the error in the integration step.
99  virtual void updateDerivsBlock(
100  const std::list<IndexPair> & solGIDList,
101  const std::list<IndexPair> & staGIDList) = 0;
102 
103  // Compute an estimate of the error in the integration step.
104  virtual double computeErrorEstimate() const = 0;
105 
106  // Interpolate solution, state or store approximation at prescribed time point.
107  virtual bool interpolateSolution(double timepoint, Linear::Vector * tmpSolVectorPtr, std::vector<Linear::Vector*> & historyVec) = 0;
108 
109  // Print output using interpolation when order is high
110  virtual bool printOutputSolution(
111  Analysis::OutputMgrAdapter & outputManagerAdapter,
112  const TIAParams & tia_params,
113  const double time,
114  Linear::Vector * solnVecPtr,
115  const bool doNotInterpolate,
116  const std::vector<double> & outputInterpolationTimes,
117  bool skipPrintLineOutput) = 0;
118 
119  // Print MPDE output using local interpolation methods
121  Analysis::OutputMgrAdapter & outputManagerAdapter,
122  const double time,
123  Linear::Vector * solnVecPtr,
124  const std::vector<double> & fastTimes)
125  {
126  return false;
127  }
128 
129  // Print WaMPDE output using local interpolation methods
131  Analysis::OutputMgrAdapter & outputManagerAdapter,
132  const double time,
133  Linear::Vector * solnVecPtr,
134  const std::vector<double> & fastTimes,
135  const int phiGID)
136  {
137  return false;
138  }
139 
140  // .SAVE output using interpolation when order is high
141  virtual bool saveOutputSolution(
142  Analysis::OutputMgrAdapter & outputManagerAdapter,
143  const TIAParams & tia_params,
144  Linear::Vector * solnVecPtr,
145  const double saveTime,
146  const bool doNotInterpolate) = 0;
147 
148  // Computes the step adjustment.
149  virtual double computeExpoStepAdjust(double stepadjust) = 0;
150 
151  // Gets the time-integration order.
152  virtual int getOrder() const = 0;
153  virtual int getNumberOfSteps() const = 0;
154  virtual int getUsedOrder() const = 0;
155  virtual int getNscsco() const = 0;
156  virtual int getMaxOrder() const = 0; ///< Return max order of method (this should obey user option maxorder)
157 
158  virtual void getInitialQnorm(TwoLevelError & tle) const = 0;
159  virtual void getTwoLevelError(TwoLevelError & tle) const = 0;
160 
161  // This is for new-DAE.
162  virtual void updateStateDeriv ()
163  {}
164 
165  // calculates dQ/dt component of store vector and adds it to store vector
166  virtual void updateLeadCurrent ()
167  {}
168 
169  // calculates dQ/dt component of lead current vector and adds it to the lead current vector
170  virtual void updateLeadCurrentVec ()
171  {}
172 
173 
174  // Returns the current partial time derivative for either the solution or
175  // state vector.
176  virtual double partialTimeDeriv() const = 0;
177 
178  // Gets the leading coefficient for the specified time-integration method.
179  virtual double getLeadingCoeff() const = 0;
180 
181  // sets the leading coefficient for the specified time-integration method.
182  virtual void setLeadingCoeff(double & LC) = 0;
183 
184  // Evaluate residual for nonlinear solver
185  virtual void obtainResidual() = 0;
186 
187  // obtain residuals for transient sensitivity analysis
188  virtual void obtainSensitivityResiduals() = 0;
189 
190  // obtain final derivatives w.r.t. time for transient sensitivity analysis
191  virtual void loadFinalSensitivityDerivatives() = 0;
192 
193  // Evaluate Jacobian for nonlinear solver
194  virtual void obtainJacobian()
195  {}
196 
197  // Apply Jacobian for nonlinear solver
198  virtual void applyJacobian(const Linear::Vector& input, Linear::Vector& result)
199  {}
200 
201  // Update history array after a successful step
202  virtual void updateHistory()
203  {}
204 
205  // Restore history array after a failed step
206  virtual void restoreHistory()
207  {}
208 
209  // Update method coefficients
210  virtual void updateCoeffs()
211  {}
212 
213  // Initialize method with initial solution & step-size
214  virtual void initialize(const TIAParams &tia_params) = 0;
215 
216  // setup 2-level data.
217  virtual void setTwoLevelTimeInfo(const TimeIntInfo & tiInfo) = 0;
218 
219  // Reject a step (this turns back history and chooses new order & step-size)
220  virtual void rejectStep(const TIAParams &tia_params) = 0;
221 
222  // Reject a step, but only restore the history, as the new step will be
223  // imposed by an outside program.
224  virtual void rejectStepForHabanero()
225  {}
226 
227  // Complete a step (this updates history and chooses new order & step-size)
228  virtual void completeStep(const TIAParams &tia_params) = 0;
229 };
230 
231 } // namespace TimeIntg
232 } // namespace Xyce
233 
234 #endif // Xyce_N_TIA_TimeIntegrationMethods_h
virtual void updateDerivsBlock(const std::list< IndexPair > &solGIDList, const std::list< IndexPair > &staGIDList)=0
virtual const char * getName() const =0
virtual bool printMPDEOutputSolution(Analysis::OutputMgrAdapter &outputManagerAdapter, const double time, Linear::Vector *solnVecPtr, const std::vector< double > &fastTimes)
Pure virtual class to augment a linear system.
virtual void getTwoLevelError(TwoLevelError &tle) const =0
virtual int getNumberOfSteps() const =0
virtual double getLeadingCoeff() const =0
virtual bool interpolateSolution(double timepoint, Linear::Vector *tmpSolVectorPtr, std::vector< Linear::Vector * > &historyVec)=0
virtual void initialize(const TIAParams &tia_params)=0
virtual bool printWaMPDEOutputSolution(Analysis::OutputMgrAdapter &outputManagerAdapter, const double time, Linear::Vector *solnVecPtr, const std::vector< double > &fastTimes, const int phiGID)
virtual void completeStep(const TIAParams &tia_params)=0
virtual int getUsedOrder() const =0
virtual double computeErrorEstimate() const =0
virtual void applyJacobian(const Linear::Vector &input, Linear::Vector &result)
virtual void setLeadingCoeff(double &LC)=0
virtual double partialTimeDeriv() const =0
virtual void getInitialQnorm(TwoLevelError &tle) const =0
virtual bool printOutputSolution(Analysis::OutputMgrAdapter &outputManagerAdapter, const TIAParams &tia_params, const double time, Linear::Vector *solnVecPtr, const bool doNotInterpolate, const std::vector< double > &outputInterpolationTimes, bool skipPrintLineOutput)=0
virtual double computeExpoStepAdjust(double stepadjust)=0
virtual void setTwoLevelTimeInfo(const TimeIntInfo &tiInfo)=0
virtual int getNscsco() const =0
virtual int getMaxOrder() const =0
Return max order of method (this should obey user option maxorder)
virtual void loadFinalSensitivityDerivatives()=0
virtual bool saveOutputSolution(Analysis::OutputMgrAdapter &outputManagerAdapter, const TIAParams &tia_params, Linear::Vector *solnVecPtr, const double saveTime, const bool doNotInterpolate)=0
virtual void rejectStep(const TIAParams &tia_params)=0