21 #ifndef OPM_NONUNIFORMTABLELINEAR_HEADER_INCLUDED
22 #define OPM_NONUNIFORMTABLELINEAR_HEADER_INCLUDED
29 #include <opm/common/ErrorMacros.hpp>
30 #include <opm/parser/eclipse/EclipseState/Tables/TableColumn.hpp>
31 #include <opm/core/utility/linearInterpolation.hpp>
54 const std::vector<T>& y_values);
57 const std::vector<T>& y_values);
60 const TableColumn& y_column);
64 std::pair<double, double>
domain();
83 double inverse(
const double y)
const;
91 std::vector<double> x_values_;
92 std::vector<T> y_values_;
93 mutable std::vector<T> x_values_reversed_;
94 mutable std::vector<T> y_values_reversed_;
105 template <
typename FI>
106 bool isNondecreasing(
const FI beg,
const FI end)
108 if (beg == end)
return true;
112 for (; it != end; ++it, ++prev) {
136 const TableColumn& y_column)
137 : x_values_( x_column.begin() , x_column.end()),
138 y_values_( y_column.begin() , y_column.end())
140 assert(isNondecreasing(x_values_.begin(), x_values_.end()));
148 const std::vector<T>& y_values)
149 : x_values_( x_column.begin() , x_column.end()),
152 assert(isNondecreasing(x_values_.begin(), x_values_.end()));
161 const std::vector<T>& y_values)
162 : x_values_(x_values), y_values_(y_values)
164 assert(isNondecreasing(x_values.begin(), x_values.end()));
168 inline std::pair<double, double>
172 return std::make_pair(x_values_[0], x_values_.back());
180 const double a = x_values_[0];
181 const double b = x_values_.back();
182 const double c = new_domain.first;
183 const double d = new_domain.second;
185 for (
int i = 0; i < int(x_values_.size()); ++i) {
186 x_values_[i] = (x_values_[i] - a)*(d - c)/(b - a) + c;
195 return Opm::linearInterpolation(x_values_, y_values_, x);
203 return Opm::linearInterpolationDerivative(x_values_, y_values_, x);
211 if (y_values_.front() < y_values_.back()) {
212 return Opm::linearInterpolation(y_values_, x_values_, y);
214 if (y_values_reversed_.empty()) {
215 y_values_reversed_ = y_values_;
216 std::reverse(y_values_reversed_.begin(), y_values_reversed_.end());
217 assert(isNondecreasing(y_values_reversed_.begin(), y_values_reversed_.end()));
218 x_values_reversed_ = x_values_;
219 std::reverse(x_values_reversed_.begin(), x_values_reversed_.end());
221 return Opm::linearInterpolation(y_values_reversed_, x_values_reversed_, y);
230 return x_values_ == other.x_values_
231 && y_values_ == other.y_values_;
236 #endif // OPM_NONUNIFORMTABLELINEAR_HEADER_INCLUDED