Tempus Version of the Day
Time Integration
Loading...
Searching...
No Matches
Tempus_TimeEventRange_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_TimeEventRange_impl_hpp
10#define Tempus_TimeEventRange_impl_hpp
11
13
14namespace Tempus {
15
16template<class Scalar>
18 : start_ (0.0),
19 stop_ (0.0),
20 stride_(0.0),
21 numEvents_(1),
22 timeScale_(1.0),
23 relTol_(this->getDefaultTol()),
24 absTol_(this->getDefaultTol()),
25 landOnExactly_(true)
26{
27 this->setType("Range");
28 setRelTol(this->getDefaultTol()),
29 setTimeRange(0.0, 0.0, 0.0);
30 setLandOnExactly(true);
31 std::ostringstream oss;
32 oss << "TimeEventRange (" << start_<< "; " << stop_<< "; " << stride_<< ")";
33 this->setName(oss.str());
34}
35
36
37template<class Scalar>
39 Scalar start, Scalar stop, Scalar stride,
40 std::string name, bool landOnExactly, Scalar relTol)
41 : start_ (start),
42 stop_ (stop),
43 stride_(stride),
44 numEvents_(0),
45 timeScale_(std::max(std::abs(start_), std::abs(stop_))),
46 relTol_(relTol),
47 absTol_(relTol_*timeScale_),
48 landOnExactly_(landOnExactly)
49{
50 this->setType("Range");
51 if (name == "") {
52 std::ostringstream oss;
53 oss << "TimeEventRange (" << start << "; " << stop << "; " << stride << ")";
54 this->setName(oss.str());
55 } else {
56 this->setName(name);
57 }
58 setRelTol(relTol);
59 setTimeRange(start, stop, stride);
60 setLandOnExactly(landOnExactly);
61}
62
63
64template<class Scalar>
66 Scalar start, Scalar stop, int numEvents,
67 std::string name, bool landOnExactly, Scalar relTol)
68 : start_ (start),
69 stop_ (stop),
70 stride_(0.0),
71 numEvents_(numEvents),
72 timeScale_(std::max(std::abs(start_), std::abs(stop_))),
73 relTol_(relTol),
74 absTol_(relTol_*timeScale_),
75 landOnExactly_(landOnExactly)
76{
77 if (name == "") {
78 std::ostringstream oss;
79 oss << "TimeEventRange (" << start << "; " << stop << "; " << numEvents << ")";
80 this->setName(oss.str());
81 } else {
82 this->setName(name);
83 }
84 setRelTol(relTol);
85 setTimeRange(start, stop, numEvents);
86 setLandOnExactly(landOnExactly);
87}
88
89
90template<class Scalar>
92 Scalar start, Scalar stop, Scalar stride)
93{
94 start_ = start;
95 stop_ = stop;
96 if (stop_ < start_) {
97 Scalar tmp = start_;
98 start_ = stop_;
99 stop_ = tmp;
100 }
101 setTimeScale();
102 setTimeStride(stride);
103}
104
105
106template<class Scalar>
108 Scalar start, Scalar stop, int numEvents)
109{
110 start_ = start;
111 stop_ = stop;
112 if (stop_ < start_) {
113 Scalar tmp = start_;
114 start_ = stop_;
115 stop_ = tmp;
116 }
117 setTimeScale();
118 setNumEvents(numEvents);
119}
120
121
122template<class Scalar>
124{
125 start_ = start;
126 if (stop_ < start_) stop_ = start_;
127 setTimeScale();
128 setTimeStride(stride_); // Reset numEvents with the current stride.
129}
130
131
132template<class Scalar>
134{
135 stop_ = stop;
136 if (start_ > stop_) start_ = stop_;
137 setTimeScale();
138 setTimeStride(stride_); // Reset numEvents with the current stride.
139}
140
141
142template<class Scalar>
144{
145 timeScale_ = std::max(std::abs(start_), std::abs(stop_));
146 absTol_ = relTol_*timeScale_;
147
148 // Check if timeScale is near zero.
149 if ( approxZero(timeScale_, absTol_) ) {
150 timeScale_ = 1.0;
151 absTol_ = relTol_*timeScale_;
152 }
153}
154
155
156template<class Scalar>
158{
159 stride_ = Teuchos::ScalarTraits<Scalar>::magnitude(stride);
160 if ( approxEqualAbsTol(start_, stop_, absTol_) ) {
161 stride_ = 0.0;
162 numEvents_ = 1;
163 return;
164 }
165
166 if ((stride_ > stop_ - start_) || (stride_ < 2*absTol_)) {
167 stride_ = stop_ - start_;
168 }
169
170 numEvents_ = int((stop_+absTol_ - start_) / stride_) + 1;
171}
172
173
174template<class Scalar>
176{
177 numEvents_ = numEvents;
178 if (numEvents_ < 0) { // Do not use numEvents_ to set stride! Reset numEvents_ from stride_.
179 if (stride_ < 2 * absTol_) stride_ = 2*absTol_;
180 numEvents_ = int((stop_+absTol_ - start_) / stride_) + 1;
181 return;
182 } else if ( approxEqualAbsTol(start_, stop_, absTol_) ) {
183 numEvents_ = 1;
184 stride_ = 0.0;
185 stride_ = stop_ - start_;
186 } else {
187 if (numEvents_ < 2) numEvents_ = 2;
188 stride_ = (stop_ - start_)/Scalar(numEvents_-1);
189 }
190
191 // If stride_ is smaller than twice the absolute tolerance,
192 // the time steps cannot be differentiated!
193 if (stride_ <= 2 * absTol_) setTimeStride(2*absTol_);
194}
195
196
197template<class Scalar>
199{
200 relTol_ = std::abs(relTol);
201 setTimeScale();
202}
203
204
205template<class Scalar>
206bool TimeEventRange<Scalar>::isTime(Scalar time) const
207{
208 // Check if before first event.
209 if (time < start_-absTol_) return false;
210
211 // Check if after last event.
212 const Scalar timeOfLast = start_ + (numEvents_-1) * stride_;
213 if (time > timeOfLast+absTol_) return false;
214
215 int numStrides = 0;
216 if ( !approxZero(stride_, 2*absTol_) )
217 numStrides = (time - start_) / stride_;
218
219 numStrides = std::min(std::max(0, numStrides), int(numEvents_-1));
220 const Scalar leftBracket = start_ + numStrides * stride_;
221
222 // Check if close to left bracket.
223 if ( approxEqualAbsTol(time, leftBracket, absTol_) )
224 return true;
225
226 // Check if close to right bracket.
227 const Scalar rightBracket = leftBracket + stride_;
228 if ( approxEqualAbsTol(time, rightBracket, absTol_) )
229 return true;
230
231 return false;
232}
233
234
235template<class Scalar>
237{
238 return timeOfNextEvent(time) - time; // Neg. time indicates in the past.
239}
240
241
242template<class Scalar>
244{
245 // Check if before first event.
246 if (time < start_-absTol_) return start_;
247
248 const Scalar timeOfLast = start_ + (numEvents_-1) * stride_;
249 // Check if after or close to last event.
250 if (time > timeOfLast-absTol_) return std::numeric_limits<Scalar>::max();
251
252 int numStrides = 0;
253 if ( !approxZero(stride_, 2*absTol_) )
254 numStrides = int((time - start_) / stride_) + 1;
255 const Scalar timeEvent = start_ + numStrides*stride_;
256
257 // Check timeEvent is near time. If so, return the next event.
258 if ( approxEqualAbsTol(time, timeEvent, absTol_) )
259 return timeEvent + stride_;
260
261 return timeEvent;
262}
263
264
265template<class Scalar>
266bool TimeEventRange<Scalar>::eventInRange(Scalar time1, Scalar time2) const
267{
268 if (time1 > time2) {
269 Scalar tmp = time1;
270 time1 = time2;
271 time2 = tmp;
272 }
273
274 // Check if range is completely outside time events.
275 const Scalar timeOfLast = start_ + (numEvents_-1) * stride_;
276 if (time2 < start_-absTol_ || timeOfLast+absTol_ < time1) return false;
277
278 if ( approxZero(stride_) )
279 return (time1 < start_-absTol_ && start_-absTol_ <= time2);
280
281 const int strideJustBeforeTime1 = std::min(int(numEvents_-1),
282 std::max(int(0), int((time1 - start_ + absTol_) / stride_ - 0.5)));
283
284 const int strideJustAfterTime2 = std::max(int(0), std::min(int(numEvents_-1),
285 int((time2 - start_ + absTol_) / stride_ + 0.5)));
286
287 for ( int i = strideJustBeforeTime1; i <= strideJustAfterTime2; i++ ) {
288 const Scalar timeEvent = start_ + i * stride_;
289 if (time1 < timeEvent-absTol_ && timeEvent-absTol_ <= time2) return true;
290 }
291
292 return false;
293}
294
295
296template<class Scalar>
297void TimeEventRange<Scalar>::describe(Teuchos::FancyOStream &out,
298 const Teuchos::EVerbosityLevel verbLevel) const
299{
300 auto l_out = Teuchos::fancyOStream( out.getOStream() );
301 Teuchos::OSTab ostab(*l_out, 2, "TimeEventRange");
302 l_out->setOutputToRootOnly(0);
303
304 *l_out << "TimeEventRange:" << "\n"
305 << " name = " << this->getName() << "\n"
306 << " Type = " << this->getType() << "\n"
307 << " start_ = " << start_ << "\n"
308 << " stop_ = " << stop_ << "\n"
309 << " stride_ = " << stride_ << "\n"
310 << " numEvents_ = " << numEvents_ << "\n"
311 << " timeScale_ = " << timeScale_ << "\n"
312 << " relTol_ = " << relTol_ << "\n"
313 << " absTol_ = " << absTol_ << "\n"
314 << " landOnExactly_ = " << landOnExactly_ << std::endl;
315}
316
317
318template<class Scalar>
319Teuchos::RCP<const Teuchos::ParameterList>
321{
322 Teuchos::RCP<Teuchos::ParameterList> pl =
323 Teuchos::parameterList("Time Event Range");
324
325 pl->setName(this->getName());
326 pl->set("Name", this->getName());
327 pl->set("Type", this->getType());
328
329 pl->set("Start Time", getTimeStart(), "Start of time range");
330 pl->set("Stop Time", getTimeStop(), "Stop of time range");
331 pl->set("Stride Time", getTimeStride(), "Stride of time range");
332
333 if ( getTimeStride()*Scalar(getNumEvents()-1) - (getTimeStop()-getTimeStart()) < getAbsTol() )
334 pl->set("Number of Events", getNumEvents(),
335 "Number of events in time range. If specified, 'Stride Time' is reset.");
336
337 pl->set("Relative Tolerance", getRelTol(),
338 "Relative time tolerance for matching time events.");
339
340 pl->set("Land On Exactly", getLandOnExactly(),
341 "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.");
342
343 return pl;
344}
345
346
347// Nonmember constructors.
348// ------------------------------------------------------------------------
349
350template<class Scalar>
351Teuchos::RCP<TimeEventRange<Scalar> > createTimeEventRange(
352 Teuchos::RCP<Teuchos::ParameterList> pl)
353{
354 auto ter = Teuchos::rcp(new TimeEventRange<Scalar>());
355 if (pl == Teuchos::null) return ter; // Return default TimeEventRange.
356
357 TEUCHOS_TEST_FOR_EXCEPTION(
358 pl->get<std::string>("Type", "Range") != "Range",
359 std::logic_error,
360 "Error - Time Event Type != 'Range'. (='"
361 + pl->get<std::string>("Type")+"')\n");
362
363 auto validPL = *ter->getValidParameters();
364 bool numEventsFound = pl->isParameter("Number of Events");
365 if (!numEventsFound) validPL.remove("Number of Events");
366
367 pl->validateParametersAndSetDefaults(validPL);
368
369 ter->setName (pl->get("Name", "From createTimeEventRange"));
370 ter->setTimeStart (pl->get("Start Time", ter->getTimeStart()));
371 ter->setTimeStop (pl->get("Stop Time", ter->getTimeStop()));
372 ter->setTimeStride (pl->get("Stride Time", ter->getTimeStride()));
373 if (numEventsFound)
374 ter->setNumEvents (pl->get("Number of Events", ter->getNumEvents()));
375 ter->setRelTol (pl->get("Relative Tolerance", ter->getRelTol()));
376 ter->setLandOnExactly (pl->get("Land On Exactly", ter->getLandOnExactly()));
377
378 return ter;
379}
380
381
382} // namespace Tempus
383#endif // Tempus_TimeEventRange_impl_hpp
virtual void setType(std::string s)
virtual void setName(std::string name)
Set the name of the TimeEvent.
virtual Scalar getDefaultTol() const
Return the default tolerance used by TimeEvents.
TimeEventRange specifies a start, stop and stride time.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Return a valid ParameterList with current settings.
virtual void setRelTol(Scalar relTol)
Set the relative tolerance.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Describe member data.
virtual void setNumEvents(int numEvents)
Set the number of time events.
virtual void setTimeScale()
Set the time scale for the time events.
virtual void setLandOnExactly(bool LOE)
Set if the time event should be landed on exactly.
virtual void setTimeStop(Scalar stop)
Set the stop of the time range.
virtual void setTimeStride(Scalar stride)
Set the stride of the time range.
Scalar stride_
Stride of time range.
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.
Scalar start_
Start of time range.
TimeEventRange()
Default constructor.
virtual void setTimeStart(Scalar start)
Set the start of the time range.
virtual void setTimeRange(Scalar start, Scalar stop, Scalar stride)
Set the range of time events from start, stop and stride.
virtual bool isTime(Scalar time) const
Test if time is near an event (within tolerance).
Scalar stop_
Stop of time range.
virtual Scalar timeToNextEvent(Scalar time) const
How much time until the next event.
bool approxZero(Scalar value, Scalar tol=Teuchos::ScalarTraits< Scalar >::sfmin())
Test if value is approximately zero within tolerance.
Teuchos::RCP< TimeEventRange< Scalar > > createTimeEventRange(Teuchos::RCP< Teuchos::ParameterList > pList)
Nonmember Constructor via ParameterList.
bool approxEqualAbsTol(Scalar a, Scalar b, Scalar absTol)
Test if values are approximately equal within the absolute tolerance.