Tempus Version of the Day
Time Integration
Loading...
Searching...
No Matches
Tempus_TimeEventListIndex_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_TimeEventListIndex_impl_hpp
10#define Tempus_TimeEventListIndex_impl_hpp
11
12
13namespace Tempus {
14
15template<class Scalar>
17{
18 this->setType("List Index");
19 this->setName("TimeEventListIndex");
20}
21
22
23template<class Scalar>
25 std::vector<int> indexList, std::string name)
26{
27 this->setType("List Index");
28 if (name == "" && !indexList.empty()) {
29 std::ostringstream oss;
30 oss << "TimeEventListIndex (" << indexList_.front() << ", ... ,"
31 << indexList_.back() << ")";
32 this->setName(oss.str());
33 } else {
34 this->setName(name);
35 }
36
37 this->setIndexList(indexList);
38}
39
40
41template<class Scalar>
42void TimeEventListIndex<Scalar>::setIndexList(std::vector<int> indexList, bool sort)
43{
44 indexList_ = indexList;
45 if (sort) {
46 std::sort(indexList_.begin(), indexList_.end());
47 indexList_.erase(std::unique(
48 indexList_.begin(), indexList_.end()), indexList_.end());
49 }
50}
51
52
53template<class Scalar>
55{
56 if (indexList_.size() == 0) {
57 indexList_.push_back(index);
58 return;
59 }
60
61 std::vector<int>::iterator it;
62 it = std::find(indexList_.begin(), indexList_.end(), index);
63 // Check if index is already in list.
64 if (it != indexList_.end()) return;
65
66 it = std::upper_bound(indexList_.begin(), indexList_.end(), index);
67 indexList_.insert(it, index);
68}
69
70
71template<class Scalar>
73{
74 return (std::find(indexList_.begin(), indexList_.end(), index) != indexList_.end() );
75}
76
77
78template<class Scalar>
80{
81 return indexOfNextEvent(index) - index; // Neg. indicating in the past.
82}
83
84
85template<class Scalar>
87{
88 if (indexList_.size() == 0) return this->getDefaultIndex();
89
90 // Check if before first event.
91 if (index < indexList_.front()) return indexList_.front();
92
93 // Check if after last event.
94 if (index >= indexList_.back()) return this->getDefaultIndex();
95
96 std::vector<int>::const_iterator it =
97 std::upper_bound(indexList_.begin(), indexList_.end(), index);
98
99 return int(*it);
100}
101
102
103template<class Scalar>
104bool TimeEventListIndex<Scalar>::eventInRangeIndex(int index1, int index2) const
105{
106 if (index1 > index2) {
107 int tmp = index1;
108 index1 = index2;
109 index2 = tmp;
110 }
111
112 if (indexList_.size() == 0) return false;
113
114 // Check if range is completely outside index events.
115 if (index2 < indexList_.front() || indexList_.back() < index1) return false;
116
117 Scalar indexEvent1 = indexOfNextEvent(index1);
118 Scalar indexEvent2 = indexOfNextEvent(index2);
119 // Check if the next index event is different for the two indices.
120 if (indexEvent1 != indexEvent2) return true;
121
122 // Check if indices bracket index event.
123 if (index1 < indexEvent1 && indexEvent1 <= index2) return true;
124
125 return false;
126}
127
128
129template<class Scalar>
130void TimeEventListIndex<Scalar>::describe(Teuchos::FancyOStream &out,
131 const Teuchos::EVerbosityLevel verbLevel) const
132{
133 auto l_out = Teuchos::fancyOStream( out.getOStream() );
134 Teuchos::OSTab ostab(*l_out, 2, "TimeEventListIndex");
135 l_out->setOutputToRootOnly(0);
136
137 *l_out << "TimeEventListIndex:" << "\n"
138 << " name = " << this->getName() << "\n"
139 << " Type = " << this->getType() << "\n"
140 << " IndexList_ = ";
141 if (!indexList_.empty()) {
142 for (auto it = indexList_.begin(); it != indexList_.end()-1; ++it)
143 *l_out << *it << ", ";
144 *l_out << *(indexList_.end()-1) << std::endl;
145 } else {
146 *l_out << "<empty>" << std::endl;
147 }
148}
149
150
151template<class Scalar>
152Teuchos::RCP<const Teuchos::ParameterList>
154{
155 Teuchos::RCP<Teuchos::ParameterList> pl =
156 Teuchos::parameterList("Time Event List Index");
157
158 pl->setName(this->getName());
159 pl->set("Name", this->getName());
160 pl->set("Type", this->getType());
161
162 std::ostringstream list;
163 if (!indexList_.empty()) {
164 for (std::size_t i = 0; i < indexList_.size()-1; ++i)
165 list << indexList_[i] << ", ";
166 list << indexList_[indexList_.size()-1];
167 }
168 pl->set<std::string>("Index List", list.str(),
169 "Comma deliminated list of indices");
170
171 return pl;
172}
173
174
175// Nonmember constructors.
176// ------------------------------------------------------------------------
177
178template<class Scalar>
179Teuchos::RCP<TimeEventListIndex<Scalar> > createTimeEventListIndex(
180 Teuchos::RCP<Teuchos::ParameterList> pl)
181{
182 auto teli = Teuchos::rcp(new TimeEventListIndex<Scalar>());
183 if (pl == Teuchos::null) return teli; // Return default TimeEventListIndex.
184
185 TEUCHOS_TEST_FOR_EXCEPTION(
186 pl->get<std::string>("Type", "List Index") != "List Index",
187 std::logic_error,
188 "Error - Time Event Type != 'List Index'. (='"
189 + pl->get<std::string>("Type")+"')\n");
190
191 pl->validateParametersAndSetDefaults(*teli->getValidParameters());
192
193 teli->setName (pl->get("Name", "From createTimeEventListIndex"));
194
195 std::vector<int> indexList;
196 indexList.clear();
197 std::string str = pl->get<std::string>("Index List");
198 std::string delimiters(",");
199 std::string::size_type lastPos = str.find_first_not_of(delimiters, 0);
200 std::string::size_type pos = str.find_first_of(delimiters, lastPos);
201 while ((pos != std::string::npos) || (lastPos != std::string::npos)) {
202 std::string token = str.substr(lastPos,pos-lastPos);
203 indexList.push_back(int(std::stoi(token)));
204 if(pos==std::string::npos) break;
205
206 lastPos = str.find_first_not_of(delimiters, pos);
207 pos = str.find_first_of(delimiters, lastPos);
208 }
209 teli->setIndexList(indexList);
210
211 return teli;
212}
213
214
215} // namespace Tempus
216#endif // Tempus_TimeEventListIndex_impl_hpp
virtual bool eventInRangeIndex(int index1, int index2) const
Test if an event occurs within the index range.
virtual void setIndexList(std::vector< int > indexList, bool sort=true)
Set the vector of event indices.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Describe member data.
virtual void addIndex(int index)
Add the index to event vector.
virtual bool isIndex(int index) const
Test if index is a time event.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Return a valid ParameterList with current settings.
virtual int indexToNextEvent(int index) const
How many indices until the next event.
virtual int indexOfNextEvent(int index) const
Return the index of the next event following the input index.
Teuchos::RCP< TimeEventListIndex< Scalar > > createTimeEventListIndex(Teuchos::RCP< Teuchos::ParameterList > pList)
Nonmember Constructor via ParameterList.