My Project
Loading...
Searching...
No Matches
UniformTableLinear.hpp
1/*
2 Copyright 2010 SINTEF ICT, Applied Mathematics.
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19
20#ifndef OPM_UNIFORMTABLELINEAR_HEADER_INCLUDED
21#define OPM_UNIFORMTABLELINEAR_HEADER_INCLUDED
22
23#include <cmath>
24#include <exception>
25#include <vector>
26#include <utility>
27#include <ostream>
28
29#include <opm/common/ErrorMacros.hpp>
30
31namespace Opm {
32
36 template<typename T>
38 {
39 public:
42
47 UniformTableLinear(double xmin,
48 double xmax,
49 const std::vector<T>& y_values);
50
56 UniformTableLinear(double xmin,
57 double xmax,
58 const T* y_values,
59 int num_y_values);
60
63 std::pair<double, double> domain();
64
67 void rescaleDomain(std::pair<double, double> new_domain);
68
72 double operator()(const double x) const;
73
77 double derivative(const double x) const;
78
82 bool operator==(const UniformTableLinear& other) const;
83
85 enum RangePolicy {Throw = 0, ClosestValue = 1, Extrapolate = 2};
86
90
94
95 protected:
96 double xmin_;
97 double xmax_;
98 double xdelta_;
99 std::vector<T> y_values_;
100 RangePolicy left_;
101 RangePolicy right_;
102 template <typename U>
103 friend std::ostream& operator<<(std::ostream& os, const UniformTableLinear<U>& t);
104 };
105
106
107 // Member implementations.
108
109 template<typename T>
110 inline
111 UniformTableLinear<T>
112 ::UniformTableLinear()
113 : left_(ClosestValue), right_(ClosestValue)
114 {
115 }
116
117 template<typename T>
118 inline
121 double xmax,
122 const std::vector<T>& y_values)
123 : xmin_(xmin), xmax_(xmax), y_values_(y_values),
124 left_(ClosestValue), right_(ClosestValue)
125 {
126 assert(xmax > xmin);
127 assert(y_values.size() > 1);
128 xdelta_ = (xmax - xmin)/(y_values.size() - 1);
129 }
130
131 template<typename T>
132 inline
135 double xmax,
136 const T* y_values,
137 int num_y_values)
138 : xmin_(xmin), xmax_(xmax),
139 y_values_(y_values, y_values + num_y_values),
140 left_(ClosestValue), right_(ClosestValue)
141 {
142 assert(xmax > xmin);
143 assert(y_values_.size() > 1);
144 xdelta_ = (xmax - xmin)/(y_values_.size() - 1);
145 }
146
147 template<typename T>
148 inline std::pair<double, double>
151 {
152 return std::make_pair(xmin_, xmax_);
153 }
154
155 template<typename T>
156 inline void
158 ::rescaleDomain(std::pair<double, double> new_domain)
159 {
160 xmin_ = new_domain.first;
161 xmax_ = new_domain.second;
162 xdelta_ = (xmax_ - xmin_)/(y_values_.size() - 1);
163 }
164
165 template<typename T>
166 inline double
168 ::operator()(const double xparam) const
169 {
170 // Implements ClosestValue policy.
171 double x = std::min(xparam, xmax_);
172 x = std::max(x, xmin_);
173
174 // Lookup is easy since we are uniform in x.
175 double pos = (x - xmin_)/xdelta_;
176 double posi = std::floor(pos);
177 int left = int(posi);
178 if (left == int(y_values_.size()) - 1) {
179 // We are at xmax_
180 return y_values_.back();
181 }
182 double w = pos - posi;
183 return (1.0 - w)*y_values_[left] + w*y_values_[left + 1];
184 }
185
186 template<typename T>
187 inline double
189 ::derivative(const double xparam) const
190 {
191 // Implements derivative consistent
192 // with ClosestValue policy for function
193 double value;
194 if (xparam > xmax_ || xparam < xmin_) {
195 value = 0.0;
196 } else {
197 double x = std::min(xparam, xmax_);
198 x = std::max(x, xmin_);
199 // Lookup is easy since we are uniform in x.
200 double pos = (x - xmin_)/xdelta_;
201 double posi = std::floor(pos);
202 int left = int(posi);
203 if (left == int(y_values_.size()) - 1) {
204 // We are at xmax_
205 --left;
206 }
207 value = (y_values_[left + 1] - y_values_[left])/xdelta_;
208 }
209 return value;
210 }
211
212
213 template<typename T>
214 inline bool
217 {
218 return xmin_ == other.xmin_
219 && xdelta_ == other.xdelta_
220 && y_values_ == other.y_values_
221 && left_ == other.left_
222 && right_ == other.right_;
223 }
224
225 template<typename T>
226 inline void
229 {
230 if (rp != ClosestValue) {
231 OPM_THROW(std::runtime_error, "Only ClosestValue RangePolicy implemented.");
232 }
233 left_ = rp;
234 }
235
236 template<typename T>
237 inline void
240 {
241 if (rp != ClosestValue) {
242 OPM_THROW(std::runtime_error, "Only ClosestValue RangePolicy implemented.");
243 }
244 right_ = rp;
245 }
246
247
248 template <typename T>
249 inline std::ostream& operator<<(std::ostream& os, const UniformTableLinear<T>& t)
250 {
251 int n = t.y_values_.size();
252 for (int i = 0; i < n; ++i) {
253 double f = double(i)/double(n - 1);
254 os << (1.0 - f)*t.xmin_ + f*t.xmax_
255 << " " << t.y_values_[i] << '\n';
256 }
257 return os;
258 }
259
260 namespace utils
261 {
263 }
264
265
266} // namespace Opm
267
268#endif // OPM_UNIFORMTABLELINEAR_HEADER_INCLUDED
This class uses linear interpolation to compute the value (and its derivative) of a function f sample...
Definition UniformTableLinear.hpp:38
void rescaleDomain(std::pair< double, double > new_domain)
Rescale the domain.
Definition UniformTableLinear.hpp:158
double operator()(const double x) const
Evaluate the value at x.
Definition UniformTableLinear.hpp:168
UniformTableLinear()
Default constructor.
Definition UniformTableLinear.hpp:112
double derivative(const double x) const
Evaluate the derivative at x.
Definition UniformTableLinear.hpp:189
void setRightPolicy(RangePolicy rp)
Sets the behavioural policy for evaluation to the right of the domain.
Definition UniformTableLinear.hpp:239
RangePolicy
Policies for how to behave when trying to evaluate outside the domain.
Definition UniformTableLinear.hpp:85
std::pair< double, double > domain()
Get the domain.
Definition UniformTableLinear.hpp:150
void setLeftPolicy(RangePolicy rp)
Sets the behavioural policy for evaluation to the left of the domain.
Definition UniformTableLinear.hpp:228
bool operator==(const UniformTableLinear &other) const
Equality operator.
Definition UniformTableLinear.hpp:216
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30