My Project
MonotCubicInterpolator.hpp
1 /* -*-C++-*- */
2 
3 #ifndef _MONOTCUBICINTERPOLATOR_H
4 #define _MONOTCUBICINTERPOLATOR_H
5 
6 #include <vector>
7 #include <map>
8 #include <string>
9 
10 /*
11  MonotCubicInterpolator
12  Copyright (C) 2006 Statoil ASA
13 
14  This program is free software; you can redistribute it and/or
15  modify it under the terms of the GNU General Public License
16  as published by the Free Software Foundation; either version 2
17  of the License, or (at your option) any later version.
18 
19  This program is distributed in the hope that it will be useful,
20  but WITHOUT ANY WARRANTY; without even the implied warranty of
21  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22  GNU General Public License for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with this program; if not, write to the Free Software
26  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
27 */
28 
29 namespace Opm
30 {
31 
63  public:
64 
74  explicit MonotCubicInterpolator(const std::string & datafilename)
75  {
76  if (!read(datafilename)) {
77  throw("Unable to constuct MonotCubicInterpolator from file.") ;
78  } ;
79  }
80 
81 
95  explicit MonotCubicInterpolator(const char* datafilename)
96  {
97  if (!read(std::string(datafilename))) {
98  throw("Unable to constuct MonotCubicInterpolator from file.") ;
99  } ;
100  }
101 
102 
111  MonotCubicInterpolator(const char* datafilename, int xColumn, int fColumn)
112  {
113  if (!read(std::string(datafilename),xColumn,fColumn)) {
114  throw("Unable to constuct MonotCubicInterpolator from file.") ;
115  } ;
116  }
117 
126  MonotCubicInterpolator(const std::string & datafilename, int xColumn, int fColumn)
127  {
128  if (!read(datafilename,xColumn,fColumn)) {
129  throw("Unable to constuct MonotCubicInterpolator from file.") ;
130  } ;
131  }
132 
141  MonotCubicInterpolator(const std::vector<double> & x ,
142  const std::vector<double> & f);
143 
151 
152 
153 
167  bool read(const std::string & datafilename) {
168  return read(datafilename,1,2) ;
169  }
170 
179  bool read(const std::string & datafilename, int xColumn, int fColumn) ;
180 
181 
182 
194  double operator () (double x) const { return evaluate(x) ; }
195 
207  double evaluate(double x) const;
208 
228  double evaluate(double x, double & errorestimate_output ) const ;
229 
236  std::pair<double,double> getMinimumX() const {
237  // Easy since the data is sorted on x:
238  return *data.begin();
239  }
240 
247  std::pair<double,double> getMaximumX() const {
248  // Easy since the data is sorted on x:
249  return *data.rbegin();
250  }
251 
258  std::pair<double,double> getMaximumF() const ;
259 
266  std::pair<double,double> getMinimumF() const ;
267 
268 
276  std::vector<double> get_xVector() const ;
277 
286  std::vector<double> get_fVector() const ;
287 
293  void scaleData(double factor);
294 
303 
304  /* Use cached value if it can be trusted */
305  if (strictlyMonotoneCached) {
306  return strictlyMonotone;
307  }
308  else {
309  computeInternalFunctionData();
310  return strictlyMonotone;
311  }
312  }
313 
319  bool isMonotone() const {
320  if (monotoneCached) {
321  return monotone;
322  }
323  else {
324  computeInternalFunctionData();
325  return monotone;
326  }
327  }
336 
337  /* Use cached value if it can be trusted */
338  if (strictlyMonotoneCached) {
339  return (strictlyMonotone && strictlyIncreasing);
340  }
341  else {
342  computeInternalFunctionData();
343  return (strictlyMonotone && strictlyIncreasing);
344  }
345  }
346 
352  bool isMonotoneIncreasing() const {
353  if (monotoneCached) {
354  return (monotone && increasing);
355  }
356  else {
357  computeInternalFunctionData();
358  return (monotone && increasing);
359  }
360  }
369 
370  /* Use cached value if it can be trusted */
371  if (strictlyMonotoneCached) {
372  return (strictlyMonotone && strictlyDecreasing);
373  }
374  else {
375  computeInternalFunctionData();
376  return (strictlyMonotone && strictlyDecreasing);
377  }
378  }
379 
385  bool isMonotoneDecreasing() const {
386  if (monotoneCached) {
387  return (monotone && decreasing);
388  }
389  else {
390  computeInternalFunctionData();
391  return (monotone && decreasing);
392  }
393  }
394 
395 
396 
410  void addPair(double newx, double newf);
411 
424  std::pair<double,double> getMissingX() const;
425 
431  std::string toString() const;
432 
436  int getSize() const {
437  return data.size();
438  }
439 
459  void chopFlatEndpoints(const double);
460 
466  chopFlatEndpoints(1e-14);
467  }
468 
483  void shrinkFlatAreas(const double);
484 
490  shrinkFlatAreas(1e-14);
491  }
492 
493 
494 
495 private:
496 
497  // Data structure to store x- and f-values
498  std::map<double, double> data;
499 
500  // Data structure to store x- and d-values
501  mutable std::map<double, double> ddata;
502 
503 
504  // Storage containers for precomputed interpolation data
505  // std::vector<double> dvalues; // derivatives in Hermite interpolation.
506 
507  // Flag to determine whether the boolean strictlyMonotone can be
508  // trusted.
509  mutable bool strictlyMonotoneCached;
510  mutable bool monotoneCached; /* only monotone, not stricly montone */
511 
512  mutable bool strictlyMonotone;
513  mutable bool monotone;
514 
515  // if strictlyMonotone is true (and can be trusted), the two next are meaningful
516  mutable bool strictlyDecreasing;
517  mutable bool strictlyIncreasing;
518  mutable bool decreasing;
519  mutable bool increasing;
520 
521 
522  /* Hermite basis functions, t \in [0,1] ,
523  notation from:
524  http://en.wikipedia.org/w/index.php?title=Cubic_Hermite_spline&oldid=84495502
525  */
526 
527  double H00(double t) const {
528  return 2*t*t*t - 3*t*t + 1;
529  }
530  double H10(double t) const {
531  return t*t*t - 2*t*t + t;
532  }
533  double H01(double t) const {
534  return -2*t*t*t + 3*t*t;
535  }
536  double H11(double t) const {
537  return t*t*t - t*t;
538  }
539 
540 
541  void computeInternalFunctionData() const ;
542 
550  void computeSimpleDerivatives() const ;
551 
552 
559  void adjustDerivativesForMonotoneness() const ;
560 
570  bool isMonotoneCoeff(double alpha, double beta) const {
571  if ((alpha*alpha + beta*beta) <= 9) {
572  return true;
573  } else {
574  return false;
575  }
576  }
577 
578 };
579 
580 
581 } // namespace Opm
582 
583 #endif
Class to represent a one-dimensional function f with single-valued argument x.
Definition: MonotCubicInterpolator.hpp:62
double evaluate(double x, double &errorestimate_output) const
bool isStrictlyDecreasing()
Determines if the current function-value-data is strictly decreasing.
Definition: MonotCubicInterpolator.hpp:368
bool read(const std::string &datafilename)
Definition: MonotCubicInterpolator.hpp:167
std::vector< double > get_xVector() const
Provide a copy of the x-data as a vector.
std::string toString() const
Constructs a string containing the data in a table.
std::pair< double, double > getMaximumF() const
Maximum f-value, returns both x and f in a pair.
std::pair< double, double > getMissingX() const
Returns an x-value that is believed to yield the best improvement in global accuracy for the interpol...
bool isStrictlyMonotone()
Determines if the current function-value-data is strictly monotone.
Definition: MonotCubicInterpolator.hpp:302
std::pair< double, double > getMinimumF() const
Minimum f-value, returns both x and f in a pair.
int getSize() const
Definition: MonotCubicInterpolator.hpp:436
void addPair(double newx, double newf)
void chopFlatEndpoints(const double)
Checks if the function curve is flat at the endpoints, chop off endpoint data points if that is the c...
MonotCubicInterpolator(const std::string &datafilename)
Definition: MonotCubicInterpolator.hpp:74
bool isMonotoneDecreasing() const
Determines if the current function-value-data is monotone and decreasing.
Definition: MonotCubicInterpolator.hpp:385
MonotCubicInterpolator(const char *datafilename, int xColumn, int fColumn)
Definition: MonotCubicInterpolator.hpp:111
std::pair< double, double > getMinimumX() const
Minimum x-value, returns both x and f in a pair.
Definition: MonotCubicInterpolator.hpp:236
MonotCubicInterpolator()
No input, an empty function object is created.
Definition: MonotCubicInterpolator.hpp:150
MonotCubicInterpolator(const char *datafilename)
Definition: MonotCubicInterpolator.hpp:95
void shrinkFlatAreas()
Wrapper function for shrinkFlatAreas(const double) providing a default epsilon parameter.
Definition: MonotCubicInterpolator.hpp:489
std::vector< double > get_fVector() const
Provide a copy of tghe function data as a vector.
bool isMonotoneIncreasing() const
Determines if the current function-value-data is monotone and increasing.
Definition: MonotCubicInterpolator.hpp:352
void scaleData(double factor)
MonotCubicInterpolator(const std::vector< double > &x, const std::vector< double > &f)
bool read(const std::string &datafilename, int xColumn, int fColumn)
bool isStrictlyIncreasing()
Determines if the current function-value-data is strictly increasing.
Definition: MonotCubicInterpolator.hpp:335
std::pair< double, double > getMaximumX() const
Maximum x-value, returns both x and f in a pair.
Definition: MonotCubicInterpolator.hpp:247
double operator()(double x) const
Definition: MonotCubicInterpolator.hpp:194
double evaluate(double x) const
MonotCubicInterpolator(const std::string &datafilename, int xColumn, int fColumn)
Definition: MonotCubicInterpolator.hpp:126
void shrinkFlatAreas(const double)
If function is monotone, but not strictly monotone, this function will remove datapoints from interva...
bool isMonotone() const
Determines if the current function-value-data is monotone.
Definition: MonotCubicInterpolator.hpp:319
void chopFlatEndpoints()
Wrapper function for chopFlatEndpoints(const double) providing a default epsilon parameter.
Definition: MonotCubicInterpolator.hpp:465
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:29