NonuniformTableLinear.hpp
1 /*
2  Copyright 2009, 2010, 2011, 2012 SINTEF ICT, Applied Mathematics.
3  Copyright 2009, 2010, 2011, 2012 Statoil ASA.
4 
5  This file is part of the Open Porous Media project (OPM).
6 
7  OPM is free software: you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  OPM is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with OPM. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 #ifndef OPM_NONUNIFORMTABLELINEAR_HEADER_INCLUDED
22 #define OPM_NONUNIFORMTABLELINEAR_HEADER_INCLUDED
23 
24 #include <cmath>
25 #include <exception>
26 #include <vector>
27 #include <utility>
28 
29 #include <opm/common/ErrorMacros.hpp>
30 #include <opm/parser/eclipse/EclipseState/Tables/TableColumn.hpp>
31 #include <opm/core/utility/linearInterpolation.hpp>
32 
33 
34 namespace Opm
35 {
36 
37 
43  template<typename T>
45  {
46  public:
49 
53  NonuniformTableLinear(const std::vector<double>& x_values,
54  const std::vector<T>& y_values);
55 
56  NonuniformTableLinear(const TableColumn& x_column,
57  const std::vector<T>& y_values);
58 
59  NonuniformTableLinear(const TableColumn& x_column,
60  const TableColumn& y_column);
61 
64  std::pair<double, double> domain();
65 
68  void rescaleDomain(std::pair<double, double> new_domain);
69 
73  double operator()(const double x) const;
74 
78  double derivative(const double x) const;
79 
83  double inverse(const double y) const;
84 
88  bool operator==(const NonuniformTableLinear& other) const;
89 
90  protected:
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_;
95  };
96 
97 
98  // A utility function
105  template <typename FI>
106  bool isNondecreasing(const FI beg, const FI end)
107  {
108  if (beg == end) return true;
109  FI it = beg;
110  ++it;
111  FI prev = beg;
112  for (; it != end; ++it, ++prev) {
113  if (*it < *prev) {
114  return false;
115  }
116  }
117  return true;
118  }
119 
120 
121 
122  // Member implementations.
123 
124  template<typename T>
125  inline
128  {
129  }
130 
131 
132  template<typename T>
133  inline
135  ::NonuniformTableLinear(const TableColumn& x_column,
136  const TableColumn& y_column)
137  : x_values_( x_column.begin() , x_column.end()),
138  y_values_( y_column.begin() , y_column.end())
139  {
140  assert(isNondecreasing(x_values_.begin(), x_values_.end()));
141  }
142 
143 
144  template<typename T>
145  inline
147  ::NonuniformTableLinear(const TableColumn& x_column,
148  const std::vector<T>& y_values)
149  : x_values_( x_column.begin() , x_column.end()),
150  y_values_(y_values)
151  {
152  assert(isNondecreasing(x_values_.begin(), x_values_.end()));
153  }
154 
155 
156 
157  template<typename T>
158  inline
160  ::NonuniformTableLinear(const std::vector<double>& x_values,
161  const std::vector<T>& y_values)
162  : x_values_(x_values), y_values_(y_values)
163  {
164  assert(isNondecreasing(x_values.begin(), x_values.end()));
165  }
166 
167  template<typename T>
168  inline std::pair<double, double>
171  {
172  return std::make_pair(x_values_[0], x_values_.back());
173  }
174 
175  template<typename T>
176  inline void
178  ::rescaleDomain(std::pair<double, double> new_domain)
179  {
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;
184  // x in [a, b] -> x in [c, d]
185  for (int i = 0; i < int(x_values_.size()); ++i) {
186  x_values_[i] = (x_values_[i] - a)*(d - c)/(b - a) + c;
187  }
188  }
189 
190  template<typename T>
191  inline double
193  ::operator()(const double x) const
194  {
195  return Opm::linearInterpolation(x_values_, y_values_, x);
196  }
197 
198  template<typename T>
199  inline double
201  ::derivative(const double x) const
202  {
203  return Opm::linearInterpolationDerivative(x_values_, y_values_, x);
204  }
205 
206  template<typename T>
207  inline double
209  ::inverse(const double y) const
210  {
211  if (y_values_.front() < y_values_.back()) {
212  return Opm::linearInterpolation(y_values_, x_values_, y);
213  } else {
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());
220  }
221  return Opm::linearInterpolation(y_values_reversed_, x_values_reversed_, y);
222  }
223  }
224 
225  template<typename T>
226  inline bool
229  {
230  return x_values_ == other.x_values_
231  && y_values_ == other.y_values_;
232  }
233 
234 } // namespace Opm
235 
236 #endif // OPM_NONUNIFORMTABLELINEAR_HEADER_INCLUDED
This class uses linear interpolation to compute the value (and its derivative) of a function f sample...
Definition: NonuniformTableLinear.hpp:44
bool operator==(const NonuniformTableLinear &other) const
Equality operator.
Definition: NonuniformTableLinear.hpp:228
double operator()(const double x) const
Evaluate the value at x.
Definition: NonuniformTableLinear.hpp:193
double inverse(const double y) const
Evaluate the inverse at y.
Definition: NonuniformTableLinear.hpp:209
double derivative(const double x) const
Evaluate the derivative at x.
Definition: NonuniformTableLinear.hpp:201
Definition: AnisotropicEikonal.cpp:446
void rescaleDomain(std::pair< double, double > new_domain)
Rescale the domain.
Definition: NonuniformTableLinear.hpp:178
NonuniformTableLinear()
Default constructor.
Definition: NonuniformTableLinear.hpp:127
std::pair< double, double > domain()
Get the domain.
Definition: NonuniformTableLinear.hpp:170