remove some generic utility classes
moved to opm-common
This commit is contained in:
Arne Morten Kvarving
parent
5881cfe5da
commit
f4111ead14
@ -94,12 +94,9 @@ list (APPEND MAIN_SOURCE_FILES
|
|||||||
list (APPEND TEST_SOURCE_FILES
|
list (APPEND TEST_SOURCE_FILES
|
||||||
tests/test_compressedpropertyaccess.cpp
|
tests/test_compressedpropertyaccess.cpp
|
||||||
tests/test_dgbasis.cpp
|
tests/test_dgbasis.cpp
|
||||||
tests/test_cubic.cpp
|
|
||||||
tests/test_flowdiagnostics.cpp
|
tests/test_flowdiagnostics.cpp
|
||||||
tests/test_parallelistlinformation.cpp
|
tests/test_parallelistlinformation.cpp
|
||||||
tests/test_sparsevector.cpp
|
|
||||||
tests/test_velocityinterpolation.cpp
|
tests/test_velocityinterpolation.cpp
|
||||||
tests/test_uniformtablelinear.cpp
|
|
||||||
tests/test_wells.cpp
|
tests/test_wells.cpp
|
||||||
tests/test_wachspresscoord.cpp
|
tests/test_wachspresscoord.cpp
|
||||||
tests/test_linearsolver.cpp
|
tests/test_linearsolver.cpp
|
||||||
@ -191,7 +188,6 @@ list (APPEND PUBLIC_HEADER_FILES
|
|||||||
opm/core/linalg/LinearSolverPetsc.hpp
|
opm/core/linalg/LinearSolverPetsc.hpp
|
||||||
opm/core/linalg/LinearSolverUmfpack.hpp
|
opm/core/linalg/LinearSolverUmfpack.hpp
|
||||||
opm/core/linalg/ParallelIstlInformation.hpp
|
opm/core/linalg/ParallelIstlInformation.hpp
|
||||||
opm/core/linalg/blas_lapack.h
|
|
||||||
opm/core/linalg/call_umfpack.h
|
opm/core/linalg/call_umfpack.h
|
||||||
opm/core/linalg/sparse_sys.h
|
opm/core/linalg/sparse_sys.h
|
||||||
opm/core/pressure/CompressibleTpfa.hpp
|
opm/core/pressure/CompressibleTpfa.hpp
|
||||||
@ -254,10 +250,8 @@ list (APPEND PUBLIC_HEADER_FILES
|
|||||||
opm/core/utility/CompressedPropertyAccess.hpp
|
opm/core/utility/CompressedPropertyAccess.hpp
|
||||||
opm/core/utility/initHydroCarbonState.hpp
|
opm/core/utility/initHydroCarbonState.hpp
|
||||||
opm/core/utility/RegionMapping.hpp
|
opm/core/utility/RegionMapping.hpp
|
||||||
opm/core/utility/UniformTableLinear.hpp
|
|
||||||
opm/core/utility/VelocityInterpolation.hpp
|
opm/core/utility/VelocityInterpolation.hpp
|
||||||
opm/core/utility/WachspressCoord.hpp
|
opm/core/utility/WachspressCoord.hpp
|
||||||
opm/core/utility/buildUniformMonotoneTable.hpp
|
|
||||||
opm/core/utility/compressedToCartesian.hpp
|
opm/core/utility/compressedToCartesian.hpp
|
||||||
opm/core/utility/extractPvtTableIndex.hpp
|
opm/core/utility/extractPvtTableIndex.hpp
|
||||||
opm/core/utility/miscUtilities.hpp
|
opm/core/utility/miscUtilities.hpp
|
||||||
|
|||||||
@ -1,148 +0,0 @@
|
|||||||
/*
|
|
||||||
Copyright 2010 SINTEF ICT, Applied Mathematics.
|
|
||||||
|
|
||||||
This file is part of the Open Porous Media project (OPM).
|
|
||||||
|
|
||||||
OPM is free software: you can redistribute it and/or modify
|
|
||||||
it under the terms of the GNU General Public License as published by
|
|
||||||
the Free Software Foundation, either version 3 of the License, or
|
|
||||||
(at your option) any later version.
|
|
||||||
|
|
||||||
OPM is distributed in the hope that it will be useful,
|
|
||||||
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
||||||
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
||||||
GNU General Public License for more details.
|
|
||||||
|
|
||||||
You should have received a copy of the GNU General Public License
|
|
||||||
along with OPM. If not, see <http://www.gnu.org/licenses/>.
|
|
||||||
*/
|
|
||||||
|
|
||||||
#ifndef OPM_BLAS_LAPACK_HEADER_INCLUDED
|
|
||||||
#define OPM_BLAS_LAPACK_HEADER_INCLUDED
|
|
||||||
|
|
||||||
#ifdef __cplusplus
|
|
||||||
extern "C" {
|
|
||||||
#endif
|
|
||||||
|
|
||||||
#if defined(MATLAB_MEX_FILE) && MATLAB_MEX_FILE
|
|
||||||
#include <mex.h>
|
|
||||||
#undef MAT_SIZE_T
|
|
||||||
#define MAT_SIZE_T mwSignedIndex
|
|
||||||
#endif
|
|
||||||
|
|
||||||
#ifndef MAT_SIZE_T
|
|
||||||
#define MAT_SIZE_T int
|
|
||||||
#endif
|
|
||||||
|
|
||||||
|
|
||||||
/* C <- a1*op(A)*op(B) + a2*C where op(X) in {X, X.'} */
|
|
||||||
void dgemm_(const char *transA , const char *transB ,
|
|
||||||
const MAT_SIZE_T* m, const MAT_SIZE_T* n , const MAT_SIZE_T* k ,
|
|
||||||
const double* a1, const double* A , const MAT_SIZE_T* ldA,
|
|
||||||
const double* B, const MAT_SIZE_T* ldB,
|
|
||||||
const double* a2, double* C , const MAT_SIZE_T* ldC);
|
|
||||||
|
|
||||||
|
|
||||||
/* C <- a1*A*A' + a2*C *or* C <- a1*A'*A + a2*C */
|
|
||||||
void dsyrk_(const char *uplo, const char *trans,
|
|
||||||
const MAT_SIZE_T *n , const MAT_SIZE_T *k ,
|
|
||||||
const double *a1 , const double *A , const MAT_SIZE_T *ldA,
|
|
||||||
const double *a2 , double *C , const MAT_SIZE_T *ldC);
|
|
||||||
|
|
||||||
|
|
||||||
void dgeqrf_(const MAT_SIZE_T *m , const MAT_SIZE_T *n ,
|
|
||||||
double *A , const MAT_SIZE_T *ld ,
|
|
||||||
double *tau , double *work,
|
|
||||||
const MAT_SIZE_T *lwork, MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
|
|
||||||
void dorgqr_(const MAT_SIZE_T *m , const MAT_SIZE_T *n , const MAT_SIZE_T *k ,
|
|
||||||
double *A , const MAT_SIZE_T *ld , const double *tau,
|
|
||||||
double *work, const MAT_SIZE_T *lwork, MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* A <- LU(A) */
|
|
||||||
void dgetrf_(const MAT_SIZE_T *m , const MAT_SIZE_T *n ,
|
|
||||||
double *A , const MAT_SIZE_T *ld,
|
|
||||||
MAT_SIZE_T *ipiv, MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* B <- A \ B, when A is LU(A) from dgetrf() */
|
|
||||||
void dgetrs_(const char *trans, const MAT_SIZE_T *n,
|
|
||||||
const MAT_SIZE_T *nrhs ,
|
|
||||||
const double *A , const MAT_SIZE_T *lda,
|
|
||||||
const MAT_SIZE_T *ipiv , double *B,
|
|
||||||
const MAT_SIZE_T *ldb , MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* B <- A \ B, tridiagonal A with bands DL, D, DU */
|
|
||||||
void dgtsv_(const MAT_SIZE_T *n ,
|
|
||||||
const MAT_SIZE_T *nrhs ,
|
|
||||||
double *DL ,
|
|
||||||
double *D ,
|
|
||||||
double *DU ,
|
|
||||||
double *B ,
|
|
||||||
const MAT_SIZE_T *ldb ,
|
|
||||||
MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* B <- A \ B, band matrix A stored in AB with kl subdiagonals, ku superdiagonals */
|
|
||||||
void dgbsv_(const MAT_SIZE_T *n ,
|
|
||||||
const MAT_SIZE_T *kl ,
|
|
||||||
const MAT_SIZE_T *ku ,
|
|
||||||
const MAT_SIZE_T *nrhs ,
|
|
||||||
double *AB ,
|
|
||||||
const MAT_SIZE_T *ldab ,
|
|
||||||
MAT_SIZE_T *ipiv ,
|
|
||||||
double *B ,
|
|
||||||
const MAT_SIZE_T *ldb ,
|
|
||||||
MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* B <- A \ B, general solver */
|
|
||||||
void dgesv_(const MAT_SIZE_T *n,
|
|
||||||
const MAT_SIZE_T *nrhs ,
|
|
||||||
double *A ,
|
|
||||||
const MAT_SIZE_T *lda ,
|
|
||||||
MAT_SIZE_T *piv ,
|
|
||||||
double *B ,
|
|
||||||
const MAT_SIZE_T *ldb ,
|
|
||||||
MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* A <- chol(A) */
|
|
||||||
void dpotrf_(const char *uplo, const MAT_SIZE_T *n,
|
|
||||||
double *A , const MAT_SIZE_T *lda,
|
|
||||||
MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* B <- (A \ (A' \ B)), when A=DPOTRF(A_orig) */
|
|
||||||
void dpotrs_(const char *uplo, const MAT_SIZE_T *n , const MAT_SIZE_T *nrhs,
|
|
||||||
double *A , const MAT_SIZE_T *lda,
|
|
||||||
double *B , const MAT_SIZE_T *ldb, MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* A <- chol(A), packed format. */
|
|
||||||
void dpptrf_(const char *uplo, const MAT_SIZE_T *n,
|
|
||||||
double *Ap , MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* A <- (A \ (A' \ eye(n)) when A=DPPTRF(A_orig) (packed format). */
|
|
||||||
void dpptri_(const char *uplo, const MAT_SIZE_T *n,
|
|
||||||
double *Ap , MAT_SIZE_T *info);
|
|
||||||
|
|
||||||
/* y <- a1*op(A)*x + a2*y */
|
|
||||||
void dgemv_(const char *trans,
|
|
||||||
const MAT_SIZE_T *m , const MAT_SIZE_T *n,
|
|
||||||
const double *a1 , const double *A, const MAT_SIZE_T *ldA ,
|
|
||||||
const double *x, const MAT_SIZE_T *incX,
|
|
||||||
const double *a2 , double *y, const MAT_SIZE_T *incY);
|
|
||||||
|
|
||||||
|
|
||||||
/* y <- a*x + y */
|
|
||||||
void daxpy_(const MAT_SIZE_T *n, const double *a,
|
|
||||||
const double *x, const MAT_SIZE_T *incx,
|
|
||||||
double *y, const MAT_SIZE_T *incy);
|
|
||||||
|
|
||||||
/* s <- x' * y */
|
|
||||||
double ddot_(const MAT_SIZE_T *n, const double *x, const MAT_SIZE_T *incx,
|
|
||||||
const double *y, const MAT_SIZE_T *incy);
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
#ifdef __cplusplus
|
|
||||||
}
|
|
||||||
#endif
|
|
||||||
|
|
||||||
#endif /* OPM_BLAS_LAPACK_HEADER_INCLUDED */
|
|
||||||
@ -1,268 +0,0 @@
|
|||||||
/*
|
|
||||||
Copyright 2010 SINTEF ICT, Applied Mathematics.
|
|
||||||
|
|
||||||
This file is part of the Open Porous Media project (OPM).
|
|
||||||
|
|
||||||
OPM is free software: you can redistribute it and/or modify
|
|
||||||
it under the terms of the GNU General Public License as published by
|
|
||||||
the Free Software Foundation, either version 3 of the License, or
|
|
||||||
(at your option) any later version.
|
|
||||||
|
|
||||||
OPM is distributed in the hope that it will be useful,
|
|
||||||
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
||||||
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
||||||
GNU General Public License for more details.
|
|
||||||
|
|
||||||
You should have received a copy of the GNU General Public License
|
|
||||||
along with OPM. If not, see <http://www.gnu.org/licenses/>.
|
|
||||||
*/
|
|
||||||
|
|
||||||
#ifndef OPM_UNIFORMTABLELINEAR_HEADER_INCLUDED
|
|
||||||
#define OPM_UNIFORMTABLELINEAR_HEADER_INCLUDED
|
|
||||||
|
|
||||||
#include <cmath>
|
|
||||||
#include <exception>
|
|
||||||
#include <vector>
|
|
||||||
#include <utility>
|
|
||||||
#include <iostream>
|
|
||||||
|
|
||||||
#include <opm/common/ErrorMacros.hpp>
|
|
||||||
|
|
||||||
namespace Opm {
|
|
||||||
|
|
||||||
/// @brief This class uses linear interpolation to compute the value
|
|
||||||
/// (and its derivative) of a function f sampled at uniform points.
|
|
||||||
/// @tparam T the range type of the function (should be an algebraic ring type)
|
|
||||||
template<typename T>
|
|
||||||
class UniformTableLinear
|
|
||||||
{
|
|
||||||
public:
|
|
||||||
/// @brief Default constructor.
|
|
||||||
UniformTableLinear();
|
|
||||||
|
|
||||||
/// @brief Construct from vector of y-values.
|
|
||||||
/// @param xmin the x value corresponding to the first y value.
|
|
||||||
/// @param xmax the x value corresponding to the last y value.
|
|
||||||
/// @param y_values vector of range values.
|
|
||||||
UniformTableLinear(double xmin,
|
|
||||||
double xmax,
|
|
||||||
const std::vector<T>& y_values);
|
|
||||||
|
|
||||||
/// @brief Construct from array of y-values.
|
|
||||||
/// @param xmin the x value corresponding to the first y value.
|
|
||||||
/// @param xmax the x value corresponding to the last y value.
|
|
||||||
/// @param y_values array of range values.
|
|
||||||
/// @param num_y_values the number of values in y_values.
|
|
||||||
UniformTableLinear(double xmin,
|
|
||||||
double xmax,
|
|
||||||
const T* y_values,
|
|
||||||
int num_y_values);
|
|
||||||
|
|
||||||
/// @brief Get the domain.
|
|
||||||
/// @return the domain as a pair of doubles.
|
|
||||||
std::pair<double, double> domain();
|
|
||||||
|
|
||||||
/// @brief Rescale the domain.
|
|
||||||
/// @param new_domain the new domain as a pair of doubles.
|
|
||||||
void rescaleDomain(std::pair<double, double> new_domain);
|
|
||||||
|
|
||||||
/// @brief Evaluate the value at x.
|
|
||||||
/// @param x a domain value
|
|
||||||
/// @return f(x)
|
|
||||||
double operator()(const double x) const;
|
|
||||||
|
|
||||||
/// @brief Evaluate the derivative at x.
|
|
||||||
/// @param x a domain value
|
|
||||||
/// @return f'(x)
|
|
||||||
double derivative(const double x) const;
|
|
||||||
|
|
||||||
/// @brief Equality operator.
|
|
||||||
/// @param other another UniformTableLinear.
|
|
||||||
/// @return true if they are represented exactly alike.
|
|
||||||
bool operator==(const UniformTableLinear& other) const;
|
|
||||||
|
|
||||||
/// @brief Policies for how to behave when trying to evaluate outside the domain.
|
|
||||||
enum RangePolicy {Throw = 0, ClosestValue = 1, Extrapolate = 2};
|
|
||||||
|
|
||||||
/// @brief Sets the behavioural policy for evaluation to the left of the domain.
|
|
||||||
/// @param rp the policy
|
|
||||||
void setLeftPolicy(RangePolicy rp);
|
|
||||||
|
|
||||||
/// @brief Sets the behavioural policy for evaluation to the right of the domain.
|
|
||||||
/// @param rp the policy
|
|
||||||
void setRightPolicy(RangePolicy rp);
|
|
||||||
|
|
||||||
protected:
|
|
||||||
double xmin_;
|
|
||||||
double xmax_;
|
|
||||||
double xdelta_;
|
|
||||||
std::vector<T> y_values_;
|
|
||||||
RangePolicy left_;
|
|
||||||
RangePolicy right_;
|
|
||||||
template <typename U>
|
|
||||||
friend std::ostream& operator<<(std::ostream& os, const UniformTableLinear<U>& t);
|
|
||||||
};
|
|
||||||
|
|
||||||
|
|
||||||
// Member implementations.
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::UniformTableLinear()
|
|
||||||
: left_(ClosestValue), right_(ClosestValue)
|
|
||||||
{
|
|
||||||
}
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::UniformTableLinear(double xmin,
|
|
||||||
double xmax,
|
|
||||||
const std::vector<T>& y_values)
|
|
||||||
: xmin_(xmin), xmax_(xmax), y_values_(y_values),
|
|
||||||
left_(ClosestValue), right_(ClosestValue)
|
|
||||||
{
|
|
||||||
assert(xmax > xmin);
|
|
||||||
assert(y_values.size() > 1);
|
|
||||||
xdelta_ = (xmax - xmin)/(y_values.size() - 1);
|
|
||||||
}
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::UniformTableLinear(double xmin,
|
|
||||||
double xmax,
|
|
||||||
const T* y_values,
|
|
||||||
int num_y_values)
|
|
||||||
: xmin_(xmin), xmax_(xmax),
|
|
||||||
y_values_(y_values, y_values + num_y_values),
|
|
||||||
left_(ClosestValue), right_(ClosestValue)
|
|
||||||
{
|
|
||||||
assert(xmax > xmin);
|
|
||||||
assert(y_values_.size() > 1);
|
|
||||||
xdelta_ = (xmax - xmin)/(y_values_.size() - 1);
|
|
||||||
}
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline std::pair<double, double>
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::domain()
|
|
||||||
{
|
|
||||||
return std::make_pair(xmin_, xmax_);
|
|
||||||
}
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline void
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::rescaleDomain(std::pair<double, double> new_domain)
|
|
||||||
{
|
|
||||||
xmin_ = new_domain.first;
|
|
||||||
xmax_ = new_domain.second;
|
|
||||||
xdelta_ = (xmax_ - xmin_)/(y_values_.size() - 1);
|
|
||||||
}
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline double
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::operator()(const double xparam) const
|
|
||||||
{
|
|
||||||
// Implements ClosestValue policy.
|
|
||||||
double x = std::min(xparam, xmax_);
|
|
||||||
x = std::max(x, xmin_);
|
|
||||||
|
|
||||||
// Lookup is easy since we are uniform in x.
|
|
||||||
double pos = (x - xmin_)/xdelta_;
|
|
||||||
double posi = std::floor(pos);
|
|
||||||
int left = int(posi);
|
|
||||||
if (left == int(y_values_.size()) - 1) {
|
|
||||||
// We are at xmax_
|
|
||||||
return y_values_.back();
|
|
||||||
}
|
|
||||||
double w = pos - posi;
|
|
||||||
return (1.0 - w)*y_values_[left] + w*y_values_[left + 1];
|
|
||||||
}
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline double
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::derivative(const double xparam) const
|
|
||||||
{
|
|
||||||
// Implements derivative consistent
|
|
||||||
// with ClosestValue policy for function
|
|
||||||
double value;
|
|
||||||
if (xparam > xmax_ || xparam < xmin_) {
|
|
||||||
value = 0.0;
|
|
||||||
} else {
|
|
||||||
double x = std::min(xparam, xmax_);
|
|
||||||
x = std::max(x, xmin_);
|
|
||||||
// Lookup is easy since we are uniform in x.
|
|
||||||
double pos = (x - xmin_)/xdelta_;
|
|
||||||
double posi = std::floor(pos);
|
|
||||||
int left = int(posi);
|
|
||||||
if (left == int(y_values_.size()) - 1) {
|
|
||||||
// We are at xmax_
|
|
||||||
--left;
|
|
||||||
}
|
|
||||||
value = (y_values_[left + 1] - y_values_[left])/xdelta_;
|
|
||||||
}
|
|
||||||
return value;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline bool
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::operator==(const UniformTableLinear<T>& other) const
|
|
||||||
{
|
|
||||||
return xmin_ == other.xmin_
|
|
||||||
&& xdelta_ == other.xdelta_
|
|
||||||
&& y_values_ == other.y_values_
|
|
||||||
&& left_ == other.left_
|
|
||||||
&& right_ == other.right_;
|
|
||||||
}
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline void
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::setLeftPolicy(RangePolicy rp)
|
|
||||||
{
|
|
||||||
if (rp != ClosestValue) {
|
|
||||||
OPM_THROW(std::runtime_error, "Only ClosestValue RangePolicy implemented.");
|
|
||||||
}
|
|
||||||
left_ = rp;
|
|
||||||
}
|
|
||||||
|
|
||||||
template<typename T>
|
|
||||||
inline void
|
|
||||||
UniformTableLinear<T>
|
|
||||||
::setRightPolicy(RangePolicy rp)
|
|
||||||
{
|
|
||||||
if (rp != ClosestValue) {
|
|
||||||
OPM_THROW(std::runtime_error, "Only ClosestValue RangePolicy implemented.");
|
|
||||||
}
|
|
||||||
right_ = rp;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
template <typename T>
|
|
||||||
inline std::ostream& operator<<(std::ostream& os, const UniformTableLinear<T>& t)
|
|
||||||
{
|
|
||||||
int n = t.y_values_.size();
|
|
||||||
for (int i = 0; i < n; ++i) {
|
|
||||||
double f = double(i)/double(n - 1);
|
|
||||||
os << (1.0 - f)*t.xmin_ + f*t.xmax_
|
|
||||||
<< " " << t.y_values_[i] << '\n';
|
|
||||||
}
|
|
||||||
return os;
|
|
||||||
}
|
|
||||||
|
|
||||||
namespace utils
|
|
||||||
{
|
|
||||||
using Opm::UniformTableLinear;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
} // namespace Opm
|
|
||||||
|
|
||||||
#endif // OPM_UNIFORMTABLELINEAR_HEADER_INCLUDED
|
|
||||||
@ -1,50 +0,0 @@
|
|||||||
/*
|
|
||||||
Copyright 2010 SINTEF ICT, Applied Mathematics.
|
|
||||||
|
|
||||||
This file is part of the Open Porous Media project (OPM).
|
|
||||||
|
|
||||||
OPM is free software: you can redistribute it and/or modify
|
|
||||||
it under the terms of the GNU General Public License as published by
|
|
||||||
the Free Software Foundation, either version 3 of the License, or
|
|
||||||
(at your option) any later version.
|
|
||||||
|
|
||||||
OPM is distributed in the hope that it will be useful,
|
|
||||||
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
||||||
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
||||||
GNU General Public License for more details.
|
|
||||||
|
|
||||||
You should have received a copy of the GNU General Public License
|
|
||||||
along with OPM. If not, see <http://www.gnu.org/licenses/>.
|
|
||||||
*/
|
|
||||||
|
|
||||||
#ifndef OPM_BUILDUNIFORMMONOTONETABLE_HEADER_INCLUDED
|
|
||||||
#define OPM_BUILDUNIFORMMONOTONETABLE_HEADER_INCLUDED
|
|
||||||
|
|
||||||
#include <opm/core/utility/MonotCubicInterpolator.hpp>
|
|
||||||
#include <opm/core/utility/UniformTableLinear.hpp>
|
|
||||||
|
|
||||||
namespace Opm {
|
|
||||||
|
|
||||||
template <typename T>
|
|
||||||
void buildUniformMonotoneTable(const std::vector<double>& xv,
|
|
||||||
const std::vector<T>& yv,
|
|
||||||
const int samples,
|
|
||||||
UniformTableLinear<T>& table)
|
|
||||||
{
|
|
||||||
MonotCubicInterpolator interp(xv, yv);
|
|
||||||
std::vector<T> uniform_yv(samples);
|
|
||||||
double xmin = xv[0];
|
|
||||||
double xmax = xv.back();
|
|
||||||
for (int i = 0; i < samples; ++i) {
|
|
||||||
double w = double(i)/double(samples - 1);
|
|
||||||
double x = (1.0 - w)*xmin + w*xmax;
|
|
||||||
uniform_yv[i] = interp(x);
|
|
||||||
}
|
|
||||||
table = UniformTableLinear<T>(xmin, xmax, uniform_yv);
|
|
||||||
}
|
|
||||||
|
|
||||||
} // namespace Opm
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
#endif // OPM_BUILDUNIFORMMONOTONETABLE_HEADER_INCLUDED
|
|
||||||
@ -1,73 +0,0 @@
|
|||||||
//===========================================================================
|
|
||||||
//
|
|
||||||
// File: monotcubicinterpolator_test.cpp
|
|
||||||
//
|
|
||||||
// Created: Tue Dec 8 12:25:30 2009
|
|
||||||
//
|
|
||||||
// Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
|
|
||||||
// Bård Skaflestad <bard.skaflestad@sintef.no>
|
|
||||||
//
|
|
||||||
// $Date$
|
|
||||||
//
|
|
||||||
// $Revision$
|
|
||||||
//
|
|
||||||
//===========================================================================
|
|
||||||
|
|
||||||
/*
|
|
||||||
Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
|
|
||||||
Copyright 2009, 2010 Statoil ASA.
|
|
||||||
Portions Copyright 2013 Uni Research AS.
|
|
||||||
|
|
||||||
This file is part of The Open Reservoir Simulator Project (OpenRS).
|
|
||||||
|
|
||||||
OpenRS is free software: you can redistribute it and/or modify
|
|
||||||
it under the terms of the GNU General Public License as published by
|
|
||||||
the Free Software Foundation, either version 3 of the License, or
|
|
||||||
(at your option) any later version.
|
|
||||||
|
|
||||||
OpenRS is distributed in the hope that it will be useful,
|
|
||||||
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
||||||
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
||||||
GNU General Public License for more details.
|
|
||||||
|
|
||||||
You should have received a copy of the GNU General Public License
|
|
||||||
along with OpenRS. If not, see <http://www.gnu.org/licenses/>.
|
|
||||||
*/
|
|
||||||
|
|
||||||
#include "config.h"
|
|
||||||
|
|
||||||
/* --- Boost.Test boilerplate --- */
|
|
||||||
#if HAVE_DYNAMIC_BOOST_TEST
|
|
||||||
#define BOOST_TEST_DYN_LINK
|
|
||||||
#endif
|
|
||||||
|
|
||||||
#define NVERBOSE // Suppress own messages when throw()ing
|
|
||||||
|
|
||||||
#define BOOST_TEST_MODULE CubicTest
|
|
||||||
#include <boost/test/unit_test.hpp>
|
|
||||||
#include <boost/test/floating_point_comparison.hpp>
|
|
||||||
|
|
||||||
/* --- our own headers --- */
|
|
||||||
#include <opm/core/utility/MonotCubicInterpolator.hpp>
|
|
||||||
using namespace Opm;
|
|
||||||
|
|
||||||
BOOST_AUTO_TEST_SUITE ()
|
|
||||||
|
|
||||||
BOOST_AUTO_TEST_CASE (cubic)
|
|
||||||
{
|
|
||||||
const int num_v = 3;
|
|
||||||
double xv[num_v] = {0.0, 1.0, 2.0};
|
|
||||||
double fv[num_v] = {10.0, 21.0, 2.0};
|
|
||||||
std::vector<double> x(xv, xv + num_v);
|
|
||||||
std::vector<double> f(fv, fv + num_v);
|
|
||||||
MonotCubicInterpolator interp(x, f);
|
|
||||||
BOOST_REQUIRE_CLOSE (interp.evaluate(-1.0), 10., 0.00001);
|
|
||||||
BOOST_REQUIRE_CLOSE (interp.evaluate(0.0), 10., 0.00001);
|
|
||||||
BOOST_REQUIRE_CLOSE (interp.evaluate(0.0001), 10.0011, 0.00001);
|
|
||||||
BOOST_REQUIRE_CLOSE (interp.evaluate(0.5), 17.375, 0.00001);
|
|
||||||
BOOST_REQUIRE_CLOSE (interp.evaluate(1.0), 21., 0.00001);
|
|
||||||
BOOST_REQUIRE_CLOSE (interp.evaluate(2.0), 2., 0.00001);
|
|
||||||
BOOST_REQUIRE_CLOSE (interp.evaluate(4.0), 2., 0.00001);
|
|
||||||
}
|
|
||||||
|
|
||||||
BOOST_AUTO_TEST_SUITE_END()
|
|
||||||
@ -1,109 +0,0 @@
|
|||||||
//===========================================================================
|
|
||||||
//
|
|
||||||
// File: sparsevector_test.cpp
|
|
||||||
//
|
|
||||||
// Created: Mon Jun 29 21:00:53 2009
|
|
||||||
//
|
|
||||||
// Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
|
|
||||||
// Bård Skaflestad <bard.skaflestad@sintef.no>
|
|
||||||
//
|
|
||||||
// $Date$
|
|
||||||
//
|
|
||||||
// $Revision$
|
|
||||||
//
|
|
||||||
//===========================================================================
|
|
||||||
|
|
||||||
/*
|
|
||||||
Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
|
|
||||||
Copyright 2009, 2010 Statoil ASA.
|
|
||||||
|
|
||||||
This file is part of The Open Reservoir Simulator Project (OpenRS).
|
|
||||||
|
|
||||||
OpenRS is free software: you can redistribute it and/or modify
|
|
||||||
it under the terms of the GNU General Public License as published by
|
|
||||||
the Free Software Foundation, either version 3 of the License, or
|
|
||||||
(at your option) any later version.
|
|
||||||
|
|
||||||
OpenRS is distributed in the hope that it will be useful,
|
|
||||||
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
||||||
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
||||||
GNU General Public License for more details.
|
|
||||||
|
|
||||||
You should have received a copy of the GNU General Public License
|
|
||||||
along with OpenRS. If not, see <http://www.gnu.org/licenses/>.
|
|
||||||
*/
|
|
||||||
|
|
||||||
#include <config.h>
|
|
||||||
|
|
||||||
#if HAVE_DYNAMIC_BOOST_TEST
|
|
||||||
#define BOOST_TEST_DYN_LINK
|
|
||||||
#endif
|
|
||||||
#define NVERBOSE // to suppress our messages when throwing
|
|
||||||
|
|
||||||
#define BOOST_TEST_MODULE SparseVectorTest
|
|
||||||
#include <boost/test/unit_test.hpp>
|
|
||||||
|
|
||||||
#include <opm/core/utility/SparseVector.hpp>
|
|
||||||
|
|
||||||
using namespace Opm;
|
|
||||||
|
|
||||||
BOOST_AUTO_TEST_CASE(construction_and_queries)
|
|
||||||
{
|
|
||||||
const SparseVector<int> sv1;
|
|
||||||
BOOST_CHECK(sv1.empty());
|
|
||||||
BOOST_CHECK_EQUAL(sv1.size(), 0);
|
|
||||||
BOOST_CHECK_EQUAL(sv1.nonzeroSize(), 0);
|
|
||||||
|
|
||||||
const int size = 100;
|
|
||||||
const int num_elem = 9;
|
|
||||||
const int elem[num_elem] = { 9, 8, 7, 6, 5, 4, 3, 2, 1 };
|
|
||||||
const int indices[num_elem] = { 1, 2, 3, 5, 8, 13, 21, 34, 55 };
|
|
||||||
const SparseVector<int> sv2(size, elem, elem + num_elem, indices, indices + num_elem);
|
|
||||||
BOOST_CHECK(!sv2.empty());
|
|
||||||
BOOST_CHECK_EQUAL(sv2.size(), size);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.element(0), 0);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.element(1), 9);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.element(2), 8);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.element(3), 7);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.element(4), 0);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.element(5), 6);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.element(55), 1);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.element(99), 0);
|
|
||||||
BOOST_CHECK_EQUAL(sv2.nonzeroSize(), num_elem);
|
|
||||||
for (int i = 0; i < num_elem; ++i) {
|
|
||||||
BOOST_CHECK_EQUAL(sv2.nonzeroElement(i), elem[i]);
|
|
||||||
}
|
|
||||||
const SparseVector<int> sv2_again(size, elem, elem + num_elem, indices, indices + num_elem);
|
|
||||||
BOOST_CHECK(sv2 == sv2_again);
|
|
||||||
SparseVector<int> sv2_append(size, elem, elem + num_elem - 1, indices, indices + num_elem - 1);
|
|
||||||
BOOST_CHECK_EQUAL(sv2_append.nonzeroSize(), num_elem - 1);
|
|
||||||
sv2_append.addElement(elem[num_elem - 1], indices[num_elem - 1]);
|
|
||||||
BOOST_CHECK(sv2 == sv2_append);
|
|
||||||
SparseVector<int> sv2_append2(size);
|
|
||||||
for (int i = 0; i < num_elem; ++i) {
|
|
||||||
sv2_append2.addElement(elem[i], indices[i]);
|
|
||||||
}
|
|
||||||
BOOST_CHECK(sv2 == sv2_append2);
|
|
||||||
sv2_append2.clear();
|
|
||||||
SparseVector<int> sv_empty;
|
|
||||||
BOOST_CHECK(sv2_append2 == sv_empty);
|
|
||||||
|
|
||||||
// Tests that only run in debug mode.
|
|
||||||
#ifndef NDEBUG
|
|
||||||
// One element too few.
|
|
||||||
BOOST_CHECK_THROW(const SparseVector<int> sv3(size, elem, elem + num_elem - 1, indices, indices + num_elem), std::exception);
|
|
||||||
// One element too many.
|
|
||||||
BOOST_CHECK_THROW(const SparseVector<int> sv4(size, elem, elem + num_elem, indices, indices + num_elem - 1), std::exception);
|
|
||||||
// Indices out of range.
|
|
||||||
BOOST_CHECK_THROW(const SparseVector<int> sv5(4, elem, elem + num_elem, indices, indices + num_elem), std::exception);
|
|
||||||
// Indices not strictly increasing. Cheating by using the elements as indices.
|
|
||||||
BOOST_CHECK_THROW(const SparseVector<int> sv5(size, elem, elem + num_elem, elem, elem + num_elem), std::exception);
|
|
||||||
|
|
||||||
// Do not ask for out of range indices.
|
|
||||||
BOOST_CHECK_THROW(sv1.element(0), std::exception);
|
|
||||||
BOOST_CHECK_THROW(sv2.element(-1), std::exception);
|
|
||||||
BOOST_CHECK_THROW(sv2.element(sv2.size()), std::exception);
|
|
||||||
BOOST_CHECK_THROW(sv2.nonzeroElement(sv2.nonzeroSize()), std::exception);
|
|
||||||
#endif
|
|
||||||
}
|
|
||||||
|
|
||||||
@ -1,73 +0,0 @@
|
|||||||
/*
|
|
||||||
Copyright 2010 SINTEF ICT, Applied Mathematics.
|
|
||||||
|
|
||||||
This file is part of the Open Porous Media project (OPM).
|
|
||||||
|
|
||||||
OPM is free software: you can redistribute it and/or modify
|
|
||||||
it under the terms of the GNU General Public License as published by
|
|
||||||
the Free Software Foundation, either version 3 of the License, or
|
|
||||||
(at your option) any later version.
|
|
||||||
|
|
||||||
OPM is distributed in the hope that it will be useful,
|
|
||||||
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
||||||
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
||||||
GNU General Public License for more details.
|
|
||||||
|
|
||||||
You should have received a copy of the GNU General Public License
|
|
||||||
along with OPM. If not, see <http://www.gnu.org/licenses/>.
|
|
||||||
*/
|
|
||||||
#include <config.h>
|
|
||||||
|
|
||||||
#if defined(HAVE_DYNAMIC_BOOST_TEST)
|
|
||||||
#define BOOST_TEST_DYN_LINK
|
|
||||||
#endif
|
|
||||||
#define NVERBOSE // to suppress our messages when throwing
|
|
||||||
|
|
||||||
|
|
||||||
#define BOOST_TEST_MODULE UniformTableLinearTests
|
|
||||||
#include <boost/test/unit_test.hpp>
|
|
||||||
#include <opm/core/utility/UniformTableLinear.hpp>
|
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
BOOST_AUTO_TEST_CASE(table_operations)
|
|
||||||
{
|
|
||||||
// Make a simple table.
|
|
||||||
double yva[] = { 1.0, -1.0, 3.0, 4.0, 2.0 };
|
|
||||||
const int numvals = sizeof(yva)/sizeof(yva[0]);
|
|
||||||
std::vector<double> yv(yva, yva + numvals);
|
|
||||||
const double xmin = 1.0;
|
|
||||||
const double xdelta = 2.5;
|
|
||||||
const double xmax = xmin + (numvals - 1)*xdelta;
|
|
||||||
Opm::utils::UniformTableLinear<double> t1(xmin, xmax, yv);
|
|
||||||
Opm::utils::UniformTableLinear<double> t1_copy1(1.0, 11.0, yv);
|
|
||||||
Opm::utils::UniformTableLinear<double> t1_copy2(t1);
|
|
||||||
|
|
||||||
// Check equality.
|
|
||||||
BOOST_CHECK(t1 == t1_copy1);
|
|
||||||
BOOST_CHECK(t1 == t1_copy2);
|
|
||||||
|
|
||||||
// Check some evaluations.
|
|
||||||
for (int i = 0; i < numvals; ++i) {
|
|
||||||
BOOST_CHECK_EQUAL(t1(xmin + i*xdelta), yv[i]);
|
|
||||||
}
|
|
||||||
BOOST_CHECK_EQUAL(t1(2.25), 0.0);
|
|
||||||
BOOST_CHECK_EQUAL(t1(9.75), 3.0);
|
|
||||||
BOOST_CHECK_EQUAL(t1.derivative(9.75), -2.0/xdelta);
|
|
||||||
// Until we implement anything but the ClosestValue end policy, we only test that.
|
|
||||||
BOOST_CHECK_EQUAL(t1(xmin - 1.0), yv[0]);
|
|
||||||
BOOST_CHECK_EQUAL(t1(xmax + 1.0), yv.back());
|
|
||||||
|
|
||||||
// Domains.
|
|
||||||
BOOST_CHECK_EQUAL(t1.domain().first, xmin);
|
|
||||||
BOOST_CHECK_EQUAL(t1.domain().second, xmin + (numvals-1)*xdelta);
|
|
||||||
std::pair<double, double> new_domain(-100.0, 20.0);
|
|
||||||
t1.rescaleDomain(new_domain);
|
|
||||||
BOOST_CHECK_EQUAL(t1.domain().first, new_domain.first);
|
|
||||||
BOOST_CHECK_EQUAL(t1.domain().second, new_domain.second);
|
|
||||||
for (int i = 0; i < numvals; ++i) {
|
|
||||||
BOOST_CHECK_EQUAL(t1(-100.0 + i*120.0/(double(numvals - 1))), yv[i]);
|
|
||||||
}
|
|
||||||
BOOST_CHECK_EQUAL(t1(-85.0), 0.0);
|
|
||||||
BOOST_CHECK(std::fabs(t1.derivative(0.0) + 2.0/30.0) < 1e-14);
|
|
||||||
}
|
|
||||||
Loading…
Reference in New Issue
Block a user