Xyce  6.1
N_TIA_NoTimeIntegration.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_NoTimeIntegration.h,v $
27 //
28 // Purpose : This file defines the classes for the "no time integration"
29 // method.
30 //
31 // Special Notes :
32 //
33 // Creator : Eric Keiter
34 //
35 // Creation Date : 7/21/00
36 //
37 // Revision Information:
38 // ---------------------
39 //
40 // Revision Number: $Revision: 1.47.2.1 $
41 //
42 // Revision Date : $Date: 2015/04/02 18:20:18 $
43 //
44 // Current Owner : $Author: tvrusso $
45 //-----------------------------------------------------------------------------
46 
47 #ifndef Xyce_N_TIA_NO_TIME_INTEGRATION_H
48 #define Xyce_N_TIA_NO_TIME_INTEGRATION_H
49 
50 // ---------- Standard Includes ----------
51 #include <N_UTL_Math.h>
52 
53 // ---------- Xyce Includes ----------
54 #include <N_TIA_fwd.h>
55 #include <N_TIA_DataStore.h>
57 
58 #define MATRIX_FAILSAFE 1
59 #define NUM_LIMIT 1.0e-20
60 
61 namespace Xyce {
62 namespace TimeIntg {
63 
64 //-----------------------------------------------------------------------------
65 // Class : NoTimeIntegration
66 // Purpose : Class objects for use during DC analysis (applies only when
67 // all time derivatives are set to 0) (derived from
68 // TimeIntegrationMethod)
69 // Special Notes :
70 // Creator : Buddy Watts, SNL
71 // Creation Date : 6/01/00
72 //-----------------------------------------------------------------------------
74 {
75 public:
76  static const int type = 0;
77  static const char *name;
78 
79  static TimeIntegrationMethod *factory(const TIAParams &tia_params, StepErrorControl &step_error_control, DataStore &data_store);
80 
82  const TIAParams & tiaP,
83  StepErrorControl & secTmp,
84  DataStore & dsTmp);
85 
87 
88  const char *getName() const /* override */ {
89  return "None";
90  }
91 
92  // Predict solution at next time point (No Integration) /* override */.
93  void obtainPredictor() /* override */
94  {
96  }
97 
98  // Evaluate the predictor derivative formula (No Integration) /* override */.
99  void obtainPredictorDeriv() /* override */
100  {}
101 
102  // Evaluate the corrector derivative formula (No Integration) /* override */.
103  void obtainCorrectorDeriv() /* override */;
104 
105  // Compute an estimate of the error in the integration step (No Integration) /* override */.
106  void updateDerivsBlock(const std::list< IndexPair > & solGIDList,
107  const std::list< IndexPair > & staGIDList) /* override */
108  {}
109 
110  // Compute an estimate of the error in the integration step (No Integration) /* override */.
111  double computeErrorEstimate() const /* override */ { return 0.0; }
112 
113  // Interpolate solution approximation at prescribed time point (No
114  // Integration) /* override */.
115  bool interpolateSolution(double timepoint, Linear::Vector * tmpSolVectorPtr, std::vector<Linear::Vector*> & historyVec ) /* override */
116  {
117  return false;
118  }
119 
120  bool printOutputSolution(
121  Analysis::OutputMgrAdapter & outputManagerAdapter,
122  const TIAParams & tia_params,
123  const double time,
124  Linear::Vector * solnVecPtr,
125  const bool doNotInterpolate,
126  const std::vector<double> & outputInterpolationTimes,
127  bool skipPrintLineOutput) /* override */;
128 
129  bool saveOutputSolution(
130  Analysis::OutputMgrAdapter & outputManagerAdapter,
131  const TIAParams & tia_params,
132  Linear::Vector * solnVecPtr,
133  const double saveTime,
134  const bool doNotInterpolate) /* override */;
135 
136  // Computes the step adjustment (No Integration) /* override */.
137  double computeExpoStepAdjust(double stepadjust) /* override */ { return 0.0; }
138 
139  // Gets the time-integration order (No Integration) /* override */.
140  int getOrder() const /* override */
141  {
142  return 0;
143  }
144 
145  int getUsedOrder() const /* override */
146  {
147  return 0;
148  }
149 
150  int getNumberOfSteps() const /* override */
151  {
152  return 0;
153  }
154 
155  int getNscsco() const /* override */
156  {
157  return 0;
158  }
159 
160  int getMaxOrder() const /* override */
161  {
162  return -2;
163  }
164 
165  void getInitialQnorm(TwoLevelError & tle) const /* override */;
166  void getTwoLevelError(TwoLevelError & tle) const /* override */;
167 
168  double partialTimeDeriv() const /* override */
169  {
170 #ifdef MATRIX_FAILSAFE
171  return NUM_LIMIT;
172 #else
173  return 0.0;
174 #endif
175  }
176 
177  // Gets the leading coefficient for the specified time-integration method.
178  double getLeadingCoeff() const /* override */ { return leadingCoeff; }
179 
180  // sets the leading coefficient for the specified time-integration method.
181  void setLeadingCoeff(double & LC) /* override */ { leadingCoeff = LC; }
182 
183  // 03/08/04 erkeite: New functions necessary new-DAE:
184  // Evaluate residual for nonlinear solver
185  void obtainResidual() /* override */;
186 
187  // obtain residuals for transient sensitivity analysis
188  void obtainSensitivityResiduals() /* override */;
189 
190  // obtain final derivatives w.r.t. time for transient sensitivity analysis
191  void loadFinalSensitivityDerivatives() /* override */ {return;}
192 
193  // Evaluate Jacobian for nonlinear solver
194  void obtainJacobian() /* override */;
195 
196  // Apply Jacobian for nonlinear solver
197  void applyJacobian(const Linear::Vector& input, Linear::Vector& result) /* override */;
198 
199  void initialize(const TIAParams &tia_params) /* override */
200  {}
201 
202  void setTwoLevelTimeInfo(const TimeIntInfo & tiInfo) /* override */
203  {}
204 
205  void rejectStep(const TIAParams &tia_params) /* override */;
206  void completeStep(const TIAParams &tia_params) /* override */;
207 
208 private:
209  double alphas; ///< $\alpha_s$ fixed-leading coefficient of this BDF method
210  DataStore & ds; ///< Reference to the TIA data-store object.
211  StepErrorControl & sec; ///< Reference to step-error control object
212  double leadingCoeff; ///< Time-integration method leading coefficient value.
213 };
214 
215 } // namespace TimeIntg
216 } // namespace Xyce
217 
218 #endif // Xyce_N_TIA_NO_TIME_INTEGRATION_H
bool interpolateSolution(double timepoint, Linear::Vector *tmpSolVectorPtr, std::vector< Linear::Vector * > &historyVec)
void initialize(const TIAParams &tia_params)
void applyJacobian(const Linear::Vector &input, Linear::Vector &result)
StepErrorControl & sec
Reference to step-error control object.
Pure virtual class to augment a linear system.
double computeExpoStepAdjust(double stepadjust)
void completeStep(const TIAParams &tia_params)
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)
int getMaxOrder() const
Return max order of method (this should obey user option maxorder)
double alphas
$$ fixed-leading coefficient of this BDF method
NoTimeIntegration(const TIAParams &tiaP, StepErrorControl &secTmp, DataStore &dsTmp)
double leadingCoeff
Time-integration method leading coefficient value.
void updateDerivsBlock(const std::list< IndexPair > &solGIDList, const std::list< IndexPair > &staGIDList)
bool saveOutputSolution(Analysis::OutputMgrAdapter &outputManagerAdapter, const TIAParams &tia_params, Linear::Vector *solnVecPtr, const double saveTime, const bool doNotInterpolate)
void rejectStep(const TIAParams &tia_params)
DataStore & ds
Reference to the TIA data-store object.
#define NUM_LIMIT
void setTwoLevelTimeInfo(const TimeIntInfo &tiInfo)
void getTwoLevelError(TwoLevelError &tle) const
void getInitialQnorm(TwoLevelError &tle) const
static TimeIntegrationMethod * factory(const TIAParams &tia_params, StepErrorControl &step_error_control, DataStore &data_store)