Xyce  6.1
N_TIA_WorkingIntegrationMethod.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_WorkingIntegrationMethod.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.12 $
43 //
44 // Revision Date : $Date: 2015/06/05 20:49:43 $
45 //
46 // Current Owner : $Author: dgbaur $
47 //-----------------------------------------------------------------------------
48 
49 #ifndef Xyce_N_TIA_WorkingIntegrationMethods_h
50 #define Xyce_N_TIA_WorkingIntegrationMethods_h
51 
52 // ---------- Standard Includes ----------
53 #include <iostream>
54 #include <N_UTL_Math.h>
55 #include <list>
56 
57 // ---------- Xyce Includes ----------
58 #include <N_ANP_fwd.h>
59 #include <N_IO_fwd.h>
60 #include <N_LAS_fwd.h>
61 #include <N_TIA_fwd.h>
62 #include <N_UTL_fwd.h>
63 
64 #include <N_UTL_IndexPair.h>
65 
66 #include <N_UTL_Stats.h>
67 
68 namespace Xyce {
69 namespace TimeIntg {
70 
71 typedef TimeIntegrationMethod *(*Factory)(const TIAParams &tia_params, StepErrorControl &step_error_control, DataStore &data_store);
72 
74 
75 void registerFactory(int type, const char *name, Factory factory);
76 
77 template <class T>
78 inline
79 void
81 {
82  registerFactory(T::type, T::name, T::factory);
83 }
84 
85 
86 const char *getTimeIntegrationName(int type);
87 
88 TimeIntegrationMethod *createTimeIntegrationMethod(int type, const TIAParams & tia_params, StepErrorControl & step_error_control, DataStore & data_store);
89 
90 //-----------------------------------------------------------------------------
91 // Class : WorkingIntegrationMethod
92 // Purpose : This class provides a way for obtaining a specific
93 // working integration method and its associated data items.
94 // Special Notes :
95 // Creator : Buddy Watts, SNL
96 // Creation Date : 6/01/00
97 //-----------------------------------------------------------------------------
99 {
100 public:
101 // WorkingIntegrationMethod( );
102  WorkingIntegrationMethod(Stats::Stat parent_stat);
103 
104  virtual ~WorkingIntegrationMethod();
105 
106  // Method which creates a time-integration method.
108  int type,
109  const TIAParams & tia_params,
110  StepErrorControl & step_error_control,
111  DataStore & data_store);
112 
113  // accessors to the integration method object functions:
115 
116  double partialTimeDeriv() const;
117  void obtainPredictor();
118  void obtainPredictorDeriv();
119  void obtainCorrectorDeriv();
120  void updateDerivsBlock( const std::list<IndexPair> & solGIDList, const std::list<IndexPair> & staGIDList);
121 
122  int getOrder() const;
123  int getUsedOrder() const;
124  int getNumberOfSteps() const;
125  int getNscsco() const;
126  void getInitialQnorm (TwoLevelError & tle) const;
127  void getTwoLevelError(TwoLevelError & tle) const;
128  void setTwoLevelTimeInfo();
129  void updateCoeffs();
130  void rejectStepForHabanero ();
131  void initialize(const TIAParams &tia_params);
132  void completeStep(const TIAParams &tia_params);
133  void rejectStep(const TIAParams &tia_params);
134  double computeErrorEstimate() const;
135  void updateStateDeriv ();
136  void updateLeadCurrent ();
137  void updateLeadCurrentVec ();
138  void obtainResidual();
141  void obtainJacobian();
142  void applyJacobian(const Linear::Vector& input, Linear::Vector& result);
143 
145  Analysis::OutputMgrAdapter & outputManagerAdapter,
146  const double time,
147  Linear::Vector * solnVecPtr,
148  const std::vector<double> & fastTimes );
149 
151  Analysis::OutputMgrAdapter & outputManagerAdapter,
152  const double time,
153  Linear::Vector * solnVecPtr,
154  const std::vector<double> & fastTimes,
155  const int phiGID );
156 
157  bool printOutputSolution(
158  Analysis::OutputMgrAdapter & outputManagerAdapter,
159  const TIAParams & tia_params,
160  const double time,
161  Linear::Vector * solnVecPtr,
162  const bool doNotInterpolate,
163  const std::vector<double> & outputInterpolationTimes,
164  bool skipPrintLineOutput );
165 
166  bool saveOutputSolution(
167  Parallel::Machine comm,
168  IO::InitialConditionsManager & initial_conditions_manager,
169  const NodeNameMap & node_name_map,
170  const TIAParams & tia_params,
171  Linear::Vector * solnVecPtr,
172  const double saveTime,
173  const bool doNotInterpolate);
174 
175 private:
176  TimeIntegrationMethod * timeIntegrationMethod_; ///< Pointer to the integration method.
178  double jacLimit;
179 
180  Stats::Stat timeIntegratorStat_;
181  Stats::Stat predictorStat_;
182  Stats::Stat completeStepStat_;
183  Stats::Stat rejectStepStat_;
184  Stats::Stat updateCoefStat_;
185  Stats::Stat residualStat_;
186  Stats::Stat jacobianStat_;
187  Stats::Stat initializeStat_;
188 
189  Stats::Stat updateLeadStat_;
190 };
191 
192 } // namespace TimeIntg
193 } // namespace Xyce
194 
195 #endif // Xyce_N_TIA_WorkingIntegrationMethods_h
Pure virtual class to augment a linear system.
void updateDerivsBlock(const std::list< IndexPair > &solGIDList, const std::list< IndexPair > &staGIDList)
TimeIntegrationMethod * timeIntegrationMethod_
Pointer to the integration method.
const char * getTimeIntegrationName(int type)
TimeIntegrationMethod * createTimeIntegrationMethod(int type, const TIAParams &tia_params, StepErrorControl &step_error_control, DataStore &data_store)
void registerFactory(int type, const char *name, Factory factory)
bool printMPDEOutputSolution(Analysis::OutputMgrAdapter &outputManagerAdapter, const double time, Linear::Vector *solnVecPtr, const std::vector< double > &fastTimes)
TimeIntegrationMethod *(* Factory)(const TIAParams &tia_params, StepErrorControl &step_error_control, DataStore &data_store)
void createTimeIntegMethod(int type, const TIAParams &tia_params, StepErrorControl &step_error_control, DataStore &data_store)
bool saveOutputSolution(Parallel::Machine comm, IO::InitialConditionsManager &initial_conditions_manager, const NodeNameMap &node_name_map, const TIAParams &tia_params, Linear::Vector *solnVecPtr, const double saveTime, const bool doNotInterpolate)
void applyJacobian(const Linear::Vector &input, Linear::Vector &result)
bool printWaMPDEOutputSolution(Analysis::OutputMgrAdapter &outputManagerAdapter, const double time, Linear::Vector *solnVecPtr, const std::vector< double > &fastTimes, const int phiGID)
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)