Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Tempus_SolutionStateMetaData_impl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ****************************************************************************
3 // Tempus: Copyright (2017) Sandia Corporation
4 //
5 // Distributed under BSD 3-clause license (See accompanying file Copyright.txt)
6 // ****************************************************************************
7 // @HEADER
8 
9 #ifndef Tempus_SolutionStateMetaData_impl_hpp
10 #define Tempus_SolutionStateMetaData_impl_hpp
11 
12 namespace Tempus {
13 
14 template <class Scalar>
16  : time_(0.0),
17  iStep_(0),
18  dt_(0.0),
19  errorAbs_(0.0),
20  errorRel_(0.0),
21  errorRelNm1_(0.0),
22  errorRelNm2_(0.0),
23  order_(1),
24  nFailures_(0),
25  nRunningFailures_(0),
26  nConsecutiveFailures_(0),
27  tolRel_(1.0e-02),
28  tolAbs_(0.0),
29  xNormL2_(0.0),
30  dxNormL2Rel_(0.0),
31  dxNormL2Abs_(0.0),
32  computeNorms_(false),
33  solutionStatus_(WORKING),
34  output_(false),
35  outputScreen_(false),
36  isSynced_(true),
37  isInterpolated_(false),
38  accuracy_(0.0)
39 {
40 }
41 
42 template <class Scalar>
44  const Scalar time, const int iStep, const Scalar dt, const Scalar errorAbs,
45  const Scalar errorRel, const Scalar errorRelNm1, const Scalar errorRelNm2,
46  const int order, const int nFailures, const int nRunningFailures,
47  const int nConsecutiveFailures, const Scalar tolRel, const Scalar tolAbs,
48  const Scalar xNormL2, const Scalar dxNormL2Rel, const Scalar dxNormL2Abs,
49  const bool computeNorms, const Status solutionStatus, const bool output,
50  const bool outputScreen, const bool isSynced, const bool isInterpolated,
51  const Scalar accuracy)
52  : time_(time),
53  iStep_(iStep),
54  dt_(dt),
55  errorAbs_(errorAbs),
56  errorRel_(errorRel),
57  errorRelNm1_(errorRelNm1),
58  errorRelNm2_(errorRelNm2),
59  order_(order),
60  nFailures_(nFailures),
61  nRunningFailures_(nRunningFailures),
62  nConsecutiveFailures_(nConsecutiveFailures),
63  tolRel_(tolRel),
64  tolAbs_(tolAbs),
65  xNormL2_(xNormL2),
66  dxNormL2Rel_(dxNormL2Rel),
67  dxNormL2Abs_(dxNormL2Abs),
68  computeNorms_(computeNorms),
69  solutionStatus_(solutionStatus),
70  output_(output),
71  outputScreen_(outputScreen),
72  isSynced_(isSynced),
73  isInterpolated_(isInterpolated),
74  accuracy_(accuracy)
75 {
76 }
77 
78 template <class Scalar>
81  : time_(ssmd.time_),
82  iStep_(ssmd.iStep_),
83  dt_(ssmd.dt_),
84  errorAbs_(ssmd.errorAbs_),
85  errorRel_(ssmd.errorRel_),
86  errorRelNm1_(ssmd.errorRelNm1_),
87  errorRelNm2_(ssmd.errorRelNm2_),
88  order_(ssmd.order_),
89  nFailures_(ssmd.nFailures_),
90  nRunningFailures_(ssmd.nRunningFailures_),
91  nConsecutiveFailures_(ssmd.nConsecutiveFailures_),
92  tolRel_(ssmd.tolRel_),
93  tolAbs_(ssmd.tolAbs_),
94  dxNormL2Rel_(ssmd.dxNormL2Rel_),
95  dxNormL2Abs_(ssmd.dxNormL2Abs_),
96  computeNorms_(ssmd.computeNorms_),
97  solutionStatus_(ssmd.solutionStatus_),
98  output_(ssmd.output_),
99  outputScreen_(ssmd.outputScreen_),
100  isSynced_(ssmd.isSynced_),
101  isInterpolated_(ssmd.isInterpolated_),
102  accuracy_(ssmd.accuracy_)
103 {
104 }
105 
106 template <class Scalar>
109 {
112  time_, iStep_, dt_, errorAbs_, errorRel_, errorRelNm1_, errorRelNm2_,
113  order_, nFailures_, nRunningFailures_, nConsecutiveFailures_, tolRel_,
114  tolAbs_, xNormL2_, dxNormL2Rel_, dxNormL2Abs_, computeNorms_,
115  solutionStatus_, output_, outputScreen_, isSynced_, isInterpolated_,
116  accuracy_));
117 
118  return md;
119 }
120 
121 template <class Scalar>
123  const Teuchos::RCP<const SolutionStateMetaData<Scalar> >& ssmd)
124 {
125  time_ = ssmd->time_;
126  iStep_ = ssmd->iStep_;
127  dt_ = ssmd->dt_;
128  errorAbs_ = ssmd->errorAbs_;
129  errorRel_ = ssmd->errorRel_;
130  errorRelNm1_ = ssmd->errorRelNm1_;
131  errorRelNm2_ = ssmd->errorRelNm2_;
132  order_ = ssmd->order_;
133  nFailures_ = ssmd->nFailures_;
134  nRunningFailures_ = ssmd->nRunningFailures_;
135  nConsecutiveFailures_ = ssmd->nConsecutiveFailures_;
136  tolRel_ = ssmd->tolRel_, tolAbs_ = ssmd->tolAbs_, xNormL2_ = ssmd->xNormL2_,
137  dxNormL2Rel_ = ssmd->dxNormL2Rel_, dxNormL2Abs_ = ssmd->dxNormL2Abs_,
138  computeNorms_ = ssmd->computeNorms_, solutionStatus_ = ssmd->solutionStatus_;
139  output_ = ssmd->output_;
140  outputScreen_ = ssmd->outputScreen_;
141  isSynced_ = ssmd->isSynced_;
142  isInterpolated_ = ssmd->isInterpolated_;
143  accuracy_ = ssmd->accuracy_;
144 }
145 
146 template <class Scalar>
148 {
149  std::string name = "Tempus::SolutionStateMetaData";
150  return (name);
151 }
152 
153 template <class Scalar>
155  Teuchos::FancyOStream& out, const Teuchos::EVerbosityLevel verbLevel) const
156 {
157  auto l_out = Teuchos::fancyOStream(out.getOStream());
158  Teuchos::OSTab ostab(*l_out, 2, this->description());
159  l_out->setOutputToRootOnly(0);
160 
161  *l_out << "\n--- " << this->description() << " ---" << std::endl;
162 
163  if (verbLevel >= Teuchos::VERB_MEDIUM) {
164  *l_out << " time = " << time_ << std::endl
165  << " iStep = " << iStep_ << std::endl
166  << " dt = " << dt_ << std::endl
167  << " errorAbs = " << errorAbs_ << std::endl
168  << " errorRel = " << errorRel_ << std::endl
169  << " errorRelNm1 = " << errorRelNm1_ << std::endl
170  << " errorRelNm2 = " << errorRelNm2_ << std::endl
171  << " order = " << order_ << std::endl
172  << " nFailures = " << nFailures_ << std::endl
173  << " nRunningFailures = " << nRunningFailures_ << std::endl
174  << " nConsecutiveFailures = " << nConsecutiveFailures_ << std::endl
175  << " tolRel = " << tolRel_ << std::endl
176  << " tolAbs = " << tolAbs_ << std::endl
177  << " xNormL2 = " << xNormL2_ << std::endl
178  << " dxNormL2Rel = " << dxNormL2Rel_ << std::endl
179  << " dxNormL2Abs = " << dxNormL2Abs_ << std::endl
180  << " computeNorms = " << computeNorms_ << std::endl
181  << " solutionStatus = " << toString(solutionStatus_) << std::endl
182  << " output = " << output_ << std::endl
183  << " outputScreen = " << outputScreen_ << std::endl
184  << " isSynced = " << isSynced_ << std::endl
185  << " isInterpolated = " << isInterpolated_ << std::endl
186  << " accuracy = " << accuracy_ << std::endl;
187  }
188  *l_out << std::string(this->description().length() + 8, '-') << std::endl;
189 }
190 
191 } // namespace Tempus
192 #endif // Tempus_SolutionStateMetaData_impl_hpp
Teuchos::RCP< SolutionStateMetaData< Scalar > > clone() const
Clone constructor.
const std::string toString(const Status status)
Convert Status to string.
void copy(const Teuchos::RCP< const SolutionStateMetaData< Scalar > > &ssmd)
This is a deep copy.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Status
Status for the Integrator, the Stepper and the SolutionState.
RCP< std::basic_ostream< char_type, traits_type > > getOStream()