Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Tempus_TimeEventList_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_TimeEventList_impl_hpp
10 #define Tempus_TimeEventList_impl_hpp
11 
13 
14 
15 namespace Tempus {
16 
17 template<class Scalar>
19  : timeScale_(1.0), relTol_(this->getDefaultTol()),
20  absTol_(this->getDefaultTol()), landOnExactly_(true)
21 {
22  this->setType("List");
23  this->setName("TimeEventList");
24  setRelTol(this->getDefaultTol());
25  setTimeScale();
26  setLandOnExactly(true);
27 }
28 
29 
30 template<class Scalar>
32  std::vector<Scalar> timeList, std::string name,
33  bool landOnExactly, Scalar relTol)
34  : timeScale_(1.0), relTol_(this->getDefaultTol()),
35  absTol_(this->getDefaultTol()), landOnExactly_(true)
36 {
37  this->setType("List");
38  this->setName(name);
39  setRelTol(relTol);
40  setTimeScale();
41  setLandOnExactly(landOnExactly);
42  setTimeList(timeList);
43 }
44 
45 
46 template<class Scalar>
48 {
49  if (timeList_.size() == 0) {
50  timeScale_ = 1.0;
51  absTol_ = relTol_*timeScale_;
52  return;
53  }
54 
55  timeScale_ = std::max(std::abs(timeList_.front()),
56  std::abs(timeList_.back()));
57  absTol_ = relTol_*timeScale_;
58 
59  // Check if timeScale is near zero.
60  if ( approxZero(timeScale_, absTol_) ) {
61  timeScale_ = 1.0;
62  absTol_ = relTol_*timeScale_;
63  }
64 }
65 
66 
67 template<class Scalar>
68 void TimeEventList<Scalar>::setTimeList(std::vector<Scalar> timeList, bool sort)
69 {
70  timeList_ = timeList;
71  if (sort) {
72  std::sort(timeList_.begin(),timeList_.end());
73  timeList_.erase(std::unique(timeList_.begin(),
74  timeList_.end() ),
75  timeList_.end() );
76  }
77 }
78 
79 
80 template<class Scalar>
82 {
83  if (timeList_.size() == 0) {
84  timeList_.push_back(time);
85  setTimeScale();
86  return;
87  }
88 
89  auto it = std::upper_bound(timeList_.begin(), timeList_.end(), time);
90  if (timeList_.size() == 1) {
91  // Only add if a "different" time.
92  if (std::abs(timeList_.front()-time) >= absTol_)
93  timeList_.insert(it, time);
94  setTimeScale();
95  return;
96  }
97 
98  // Don't add if already in list. Check both ends.
99  if (it == timeList_.begin()) {
100  if (std::abs(timeList_.front()-time) >= absTol_)
101  timeList_.insert(it, time);
102  } else if (it == timeList_.end()) {
103  if (std::abs(timeList_.back()-time) >= absTol_)
104  timeList_.insert(it, time);
105  } else if (std::abs(*(it-1) - time) >= absTol_ &&
106  std::abs(*(it ) - time) >= absTol_) {
107  timeList_.insert(it, time);
108  }
109  setTimeScale();
110 }
111 
112 
113 template<class Scalar>
115 {
116  relTol_ = std::abs(relTol);
117  setTimeScale();
118 }
119 
120 
121 template<class Scalar>
122 bool TimeEventList<Scalar>::isTime(Scalar time) const
123 {
124  for (auto it = timeList_.begin(); it != timeList_.end(); ++it)
125  if ( approxEqualAbsTol(time, *it, absTol_) ) return true;
126 
127  return false;
128 }
129 
130 
131 template<class Scalar>
132 Scalar TimeEventList<Scalar>::timeToNextEvent(Scalar time) const
133 {
134  return timeOfNextEvent(time) - time; // Neg. indicating in the past.
135 }
136 
137 
138 template<class Scalar>
139 Scalar TimeEventList<Scalar>::timeOfNextEvent(Scalar time) const
140 {
141  if (timeList_.size() == 0) return this->getDefaultTime();
142 
143  // Check if before first event.
144  if (time < timeList_.front()-absTol_) return timeList_.front();
145 
146  // Check if after or close to last event.
147  if (time > timeList_.back()-absTol_)
148  return std::numeric_limits<Scalar>::max();
149 
150  typename std::vector<Scalar>::const_iterator it =
151  std::upper_bound(timeList_.begin(), timeList_.end(), time);
152  const Scalar timeEvent = *it;
153 
154  // Check timeEvent is near time. If so, return the next event.
155  if ( approxEqualAbsTol(time, timeEvent, absTol_) )
156  return *(it+1);
157 
158  return timeEvent;
159 }
160 
161 
162 template<class Scalar>
164  Scalar time1, Scalar time2) const
165 {
166  if (time1 > time2) {
167  Scalar tmp = time1;
168  time1 = time2;
169  time2 = tmp;
170  }
171 
172  if (timeList_.size() == 0) return false;
173 
174  for (auto it = timeList_.begin(); it != timeList_.end(); ++it)
175  if (time1+absTol_ < *it && *it < time2+absTol_) return true;
176 
177  return false;
178 }
179 
180 
181 template<class Scalar>
183  const Teuchos::EVerbosityLevel verbLevel) const
184 {
185  auto l_out = Teuchos::fancyOStream( out.getOStream() );
186  Teuchos::OSTab ostab(*l_out, 2, "TimeEventList");
187  l_out->setOutputToRootOnly(0);
188 
189  *l_out << "TimeEventList:" << "\n"
190  << " name = " << this->getName() << "\n"
191  << " Type = " << this->getType() << "\n"
192  << " timeScale_ = " << timeScale_ << "\n"
193  << " relTol_ = " << relTol_ << "\n"
194  << " absTol_ = " << absTol_ << "\n"
195  << " landOnExactly_ = " << landOnExactly_ << "\n"
196  << " timeList_ = ";
197  if (!timeList_.empty()) {
198  for (auto it = timeList_.begin(); it != timeList_.end()-1; ++it)
199  *l_out << *it << ", ";
200  *l_out << *(timeList_.end()-1) << std::endl;
201  } else {
202  *l_out << "<empty>" << std::endl;
203  }
204 }
205 
206 
207 template<class Scalar>
210 {
212  Teuchos::parameterList("Time Event List");
213 
214  pl->setName(this->getName());
215  pl->set("Name", this->getName());
216  pl->set("Type", this->getType());
217 
218  pl->set("Relative Tolerance", this->getRelTol(),
219  "Relative time tolerance for matching time events.");
220 
221  pl->set("Land On Exactly", this->getLandOnExactly(),
222  "Should these time events be landed on exactly, i.e, adjust the timestep to hit time event, versus stepping over and keeping the time step unchanged.");
223 
224  std::vector<Scalar> times = this->getTimeList();
225  std::ostringstream list;
226  if (!times.empty()) {
227  for(std::size_t i = 0; i < times.size()-1; ++i) list << times[i] << ", ";
228  list << times[times.size()-1];
229  }
230  pl->set<std::string>("Time List", list.str(),
231  "Comma deliminated list of times");
232 
233  return pl;
234 }
235 
236 
237 // Nonmember constructors.
238 // ------------------------------------------------------------------------
239 
240 template<class Scalar>
243 {
244  auto tel = Teuchos::rcp(new TimeEventList<Scalar>());
245  if (pl == Teuchos::null) return tel; // Return default TimeEventList.
246 
248  pl->get<std::string>("Type", "List") != "List",
249  std::logic_error,
250  "Error - Time Event Type != 'List'. (='"
251  + pl->get<std::string>("Type")+"')\n");
252 
253  pl->validateParametersAndSetDefaults(*tel->getValidParameters());
254 
255  tel->setName (pl->get("Name", "From createTimeEventList"));
256  tel->setRelTol (pl->get("Relative Tolerance", tel->getRelTol()));
257  tel->setLandOnExactly (pl->get("Land On Exactly", tel->getLandOnExactly()));
258 
259  std::vector<Scalar> timeList;
260  std::string str = pl->get<std::string>("Time List");
261  std::string delimiters(",");
262  // Skip delimiters at the beginning
263  std::string::size_type lastPos = str.find_first_not_of(delimiters, 0);
264  // Find the first delimiter
265  std::string::size_type pos = str.find_first_of(delimiters, lastPos);
266  while ((pos != std::string::npos) || (lastPos != std::string::npos)) {
267  // Found a token, add it to the vector
268  std::string token = str.substr(lastPos,pos-lastPos);
269  timeList.push_back(Scalar(std::stod(token)));
270  if(pos==std::string::npos) break;
271 
272  lastPos = str.find_first_not_of(delimiters, pos); // Skip delimiters
273  pos = str.find_first_of(delimiters, lastPos); // Find next delimiter
274  }
275  tel->setTimeList(timeList);
276 
277  return tel;
278 }
279 
280 
281 } // namespace Tempus
282 #endif // Tempus_TimeEventList_impl_hpp
virtual bool isTime(Scalar time) const
Test if time is near an event (within tolerance).
virtual void setType(std::string s)
virtual void setName(std::string name)
Set the name of the TimeEvent.
T & get(const std::string &name, T def_value)
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Return a valid ParameterList with current settings.
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual void addTime(Scalar time)
Add the time to event vector.
virtual Scalar timeOfNextEvent(Scalar time) const
Return the time of the next event following the input time.
virtual bool eventInRange(Scalar time1, Scalar time2) const
Test if an event occurs within the time range.
virtual Scalar getDefaultTol() const
Return the default tolerance used by TimeEvents.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
virtual void setTimeScale()
Set the time scale for the time events.
void validateParametersAndSetDefaults(ParameterList const &validParamList, int const depth=1000)
virtual void setRelTol(Scalar relTol)
Set the relative tolerance.
virtual void setTimeList(std::vector< Scalar > timeList, bool sort=true)
Set the list of time events.
TimeEventList()
Default constructor.
RCP< std::basic_ostream< char_type, traits_type > > getOStream()
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Describe member data.
bool approxEqualAbsTol(Scalar a, Scalar b, Scalar absTol)
Test if values are approximately equal within the absolute tolerance.
bool approxZero(Scalar value, Scalar tol=Teuchos::ScalarTraits< Scalar >::sfmin())
Test if value is approximately zero within tolerance.
virtual Scalar timeToNextEvent(Scalar time) const
How much time until the next event.
ParameterList & setName(const std::string &name)
virtual void setLandOnExactly(bool LOE)
Set if the time event should be landed on exactly.
TimeEventList specifies a list of time events.
Teuchos::RCP< TimeEventList< Scalar > > createTimeEventList(Teuchos::RCP< Teuchos::ParameterList > pList)
Nonmember Constructor via ParameterList.