Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 
12 
13 namespace Tempus {
14 
15 template<class Scalar>
17  : start_(this->getDefaultTime()),
18  stop_ (this->getDefaultTime()),
19  stride_(0.0),
20  numEvents_(1),
21  timeScale_(1.0),
22  relTol_(1.0e-14),
23  absTol_(1.0e-14),
24  landOnExactly_(true)
25 {
26  this->setName("TimeEventRange");
28 }
29 
30 
31 template<class Scalar>
33  std::string name, Scalar start, Scalar stop, Scalar stride,
34  Scalar relTol, bool landOnExactly)
35  : start_(start),
36  stop_ (stop),
37  stride_(stride),
38  numEvents_(std::abs(stop-start)/stride+1),
39  timeScale_(start),
40  relTol_(relTol),
41  absTol_(relTol*start),
42  landOnExactly_(landOnExactly)
43 {
44  this->setName(name);
45  setTimeRange(start, stop, stride);
46 }
47 
48 
49 template<class Scalar>
51  std::string name, Scalar start, Scalar stop, int numEvents,
52  Scalar relTol, bool landOnExactly)
53  : start_(start),
54  stop_ (stop),
55  stride_(std::abs(stop-start)/(numEvents-1)),
56  numEvents_(numEvents),
57  timeScale_(start),
58  relTol_(relTol),
59  absTol_(relTol*start),
60  landOnExactly_(landOnExactly)
61 {
62  this->setName(name);
63  setTimeRange(start, stop, numEvents);
64  setRelTol(relTol);
65  setLandOnExactly(landOnExactly);
66 }
67 
68 
69 template<class Scalar>
71 {
72  start_ = start;
73  if (stop_ < start_) stop_ = start_;
74  setTimeStride(stride_); // Reset numEvents with the current stride.
75  setTimeScale();
76 }
77 
78 
79 template<class Scalar>
81 {
82  stop_ = stop;
83  if (start_ > stop_) start_ = stop_;
84  setTimeStride(stride_); // Reset numEvents with the current stride.
85  setTimeScale();
86 }
87 
88 
89 template<class Scalar>
91 {
92  timeScale_ = std::max(std::abs(start_), std::abs(stop_));
93  absTol_ = relTol_*timeScale_;
94 
95  // Check if timeScale is near zero.
96  if ((-absTol_ <= timeScale_ ) && (timeScale_ <= absTol_)) {
97  timeScale_ = 1.0;
98  absTol_ = relTol_*timeScale_;
99  }
100 }
101 
102 
103 template<class Scalar>
105 {
106  stride_ = stride;
107  if ((start_ >= stop_-absTol_) && (start_ <= stop_+absTol_)) {
108  stride_ = 0.0;
109  numEvents_ = 1;
110  return;
111  }
112 
113  if ((stride_ > stop_ - start_) || (stride_ < 2*absTol_)) {
114  stride_ = stop_ - start_;
115  }
116 
117  numEvents_ = int((stop_+absTol_ - start_) / stride_) + 1;
118 }
119 
120 
121 template<class Scalar>
123 {
124  numEvents_ = numEvents;
125  if (numEvents_ < 1) numEvents_ = 1;
126  stride_ = (stop_ - start_)/Scalar(numEvents_-1);
127 
128  if (stride_ < 2 * absTol_) {
129  setTimeStride(2*absTol_);
130  }
131 }
132 
133 
134 template<class Scalar>
136 {
137  relTol_ = std::abs(relTol);
138  absTol_ = relTol_*timeScale_;
139 }
140 
141 
142 template<class Scalar>
143 bool TimeEventRange<Scalar>::isTime(Scalar time) const
144 {
145  return (std::abs(timeToNextEvent(time)) <= absTol_);
146 }
147 
148 
149 template<class Scalar>
151 {
152  return timeOfNextEvent(time) - time; // Neg. time indicates in the past.
153 }
154 
155 
156 template<class Scalar>
158 {
159  // Check if before or close to first event.
160  if (start_ >= time-absTol_) return start_;
161 
162  const Scalar timeOfLast = start_ + (numEvents_-1) * stride_;
163  // Check if after or close to last event.
164  if (timeOfLast <= time+absTol_) return timeOfLast;
165 
166  const int numStrides = (time - start_) / stride_;
167  const Scalar timeOfNext = start_ + numStrides * stride_;
168 
169  // Check if close to left-side time event
170  if (timeOfNext > time-absTol_ &&
171  timeOfNext < time+absTol_) return timeOfNext;
172 
173  // Otherwise it is the next event.
174  return timeOfNext + stride_;
175 }
176 
177 
178 template<class Scalar>
179 bool TimeEventRange<Scalar>::eventInRange(Scalar time1, Scalar time2) const
180 {
181  if (time1 > time2) {
182  Scalar tmp = time1;
183  time1 = time2;
184  time2 = tmp;
185  }
186 
187  const Scalar timeOfLast = start_ + (numEvents_-1) * stride_;
188  // Check if range is completely outside time events.
189  if (time2+absTol_ < start_ || timeOfLast < time1-absTol_) return false;
190 
191  Scalar timeEvent1 = timeOfNextEvent(time1);
192  Scalar timeEvent2 = timeOfNextEvent(time2);
193  // Check if the next time event is different for the two times.
194  if (timeEvent1 != timeEvent2) return true;
195 
196  // Check if times bracket time event.
197  if (time1-absTol_ <= timeEvent1 && timeEvent1 <= time2+absTol_) return true;
198 
199  return false;
200 }
201 
202 
203 template<class Scalar>
205 {
208  *out << "TimeEventRange:" << "\n"
209  << "name = " << this->getName() << "\n"
210  << "start_ = " << start_ << "\n"
211  << "stop_ = " << stop_ << "\n"
212  << "stride_ = " << stride_ << "\n"
213  << "numEvents_ = " << numEvents_ << "\n"
214  << "timeScale_ = " << timeScale_ << "\n"
215  << "relTol_ = " << relTol_ << "\n"
216  << "absTol_ = " << absTol_ << "\n"
217  << "landOnExactly_ = " << landOnExactly_ << std::endl;
218 }
219 
220 
221 } // namespace Tempus
222 #endif // Tempus_TimeEventRange_impl_hpp
virtual void describe() const
Describe member data.
virtual void setRelTol(Scalar relTol)
Scalar stride_
Stride of time range.
virtual void setNumEvents(int numEvents)
virtual void setName(std::string name)
virtual void setLandOnExactly(bool LOE)
Scalar start_
Start of time range.
virtual bool eventInRange(Scalar time1, Scalar time2) const
Test if an event occurs within the time range.
virtual Scalar timeToNextEvent(Scalar time) const
How much time until the next event. Negative indicating the last event is in the past.
static RCP< FancyOStream > getDefaultOStream()
virtual void setTimeStop(Scalar stop)
virtual void setTimeRange(Scalar start, Scalar stop, Scalar stride)
virtual void setTimeStart(Scalar start)
Scalar stop_
Stop of time range.
virtual void setTimeStride(Scalar stride)
virtual bool isTime(Scalar time) const
Test if time is near a TimeEvent (within tolerance).
TimeEventRange()
Default constructor.
virtual Scalar timeOfNextEvent(Scalar time) const
Time of the next event. Negative indicating the last event is in the past.