27#ifndef OPM_PIECEWISE_LINEAR_TWO_PHASE_MATERIAL_PARAMS_HPP
28#define OPM_PIECEWISE_LINEAR_TWO_PHASE_MATERIAL_PARAMS_HPP
36#include <opm/common/ErrorMacros.hpp>
40#include <opm/common/utility/gpuDecorators.hpp>
49template <
class TraitsT,
class VectorT = std::vector<
typename TraitsT::Scalar>>
52 using Scalar =
typename TraitsT::Scalar;
55 using ValueVector = VectorT;
57 using Traits = TraitsT;
85 EnsureFinalized ::finalize();
89 if (SwPcwnSamples_.front() > SwPcwnSamples_.back())
90 swapOrderIfPossibleThrowOtherwise_(SwPcwnSamples_, pcwnSamples_);
92 if (SwKrwSamples_.front() > SwKrwSamples_.back())
93 swapOrderIfPossibleThrowOtherwise_(SwKrwSamples_, krwSamples_);
95 if (SwKrnSamples_.front() > SwKrnSamples_.back())
96 swapOrderIfPossibleThrowOtherwise_(SwKrnSamples_, krnSamples_);
104 EnsureFinalized::check();
105 return SwKrwSamples_;
113 EnsureFinalized::check();
114 return SwKrnSamples_;
122 EnsureFinalized::check();
123 return SwPcwnSamples_;
133 EnsureFinalized::check();
142 template <
class Container>
145 assert(SwValues.size() == values.size());
147 size_t n = SwValues.size();
148 SwPcwnSamples_.resize(n);
149 pcwnSamples_.resize(n);
151 std::copy(SwValues.begin(), SwValues.end(), SwPcwnSamples_.begin());
152 std::copy(values.begin(), values.end(), pcwnSamples_.begin());
163 EnsureFinalized::check();
173 template <
class Container>
176 assert(SwValues.size() == values.size());
178 size_t n = SwValues.size();
179 SwKrwSamples_.resize(n);
180 krwSamples_.resize(n);
182 std::copy(SwValues.begin(), SwValues.end(), SwKrwSamples_.begin());
183 std::copy(values.begin(), values.end(), krwSamples_.begin());
194 EnsureFinalized::check();
204 template <
class Container>
207 assert(SwValues.size() == values.size());
209 size_t n = SwValues.size();
210 SwKrnSamples_.resize(n);
211 krnSamples_.resize(n);
213 std::copy(SwValues.begin(), SwValues.end(), SwKrnSamples_.begin());
214 std::copy(values.begin(), values.end(), krnSamples_.begin());
218 void swapOrderIfPossibleThrowOtherwise_(ValueVector& swValues, ValueVector& values)
const
222 if (swValues.front() > values.back()) {
225 if constexpr (!std::is_const_v<typename ValueVector::value_type> && !std::is_const_v<ValueVector>) {
226 for (
unsigned origSampleIdx = 0; origSampleIdx < swValues.size() / 2; ++origSampleIdx) {
227 size_t newSampleIdx = swValues.size() - origSampleIdx - 1;
229 std::swap(swValues[origSampleIdx], swValues[newSampleIdx]);
230 std::swap(values[origSampleIdx], values[newSampleIdx]);
234 OPM_THROW(std::logic_error,
"Saturation values in interpolation table provided in wrong order, but table is immutable");
239 ValueVector SwPcwnSamples_;
240 ValueVector SwKrwSamples_;
241 ValueVector SwKrnSamples_;
242 ValueVector pcwnSamples_;
243 ValueVector krwSamples_;
244 ValueVector krnSamples_;
248namespace Opm::gpuistl{
256template <
class TraitsT,
class ContainerType,
class ViewType>
257PiecewiseLinearTwoPhaseMaterialParams<TraitsT, ViewType> make_view(
const PiecewiseLinearTwoPhaseMaterialParams<TraitsT, ContainerType>& params) {
259 using containedType =
typename ContainerType::value_type;
260 using viewedTypeNoConst =
typename std::remove_const_t<typename ViewType::value_type>;
262 static_assert(std::is_same_v<containedType, viewedTypeNoConst>);
264 ViewType SwPcwnSamples = make_view<viewedTypeNoConst>(params.SwPcwnSamples());
265 ViewType pcwnSamples = make_view<viewedTypeNoConst>(params.pcwnSamples());
266 ViewType SwKrwSamples = make_view<viewedTypeNoConst>(params.SwKrwSamples());
267 ViewType krwSamples = make_view<viewedTypeNoConst>(params.krwSamples());
268 ViewType SwKrnSamples = make_view<viewedTypeNoConst>(params.SwKrnSamples());
269 ViewType krnSamples = make_view<viewedTypeNoConst>(params.krnSamples());
271 return PiecewiseLinearTwoPhaseMaterialParams<TraitsT, ViewType> (SwPcwnSamples,
Default implementation for asserting finalization of parameter objects.
Default implementation for asserting finalization of parameter objects.
Definition EnsureFinalized.hpp:49
Specification of the material parameters for a two-phase material law which uses a table and piecewis...
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:51
OPM_HOST_DEVICE const ValueVector & SwKrwSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:102
OPM_HOST_DEVICE const ValueVector & krwSamples() const
Return the sampling points for the relative permeability curve of the wetting phase.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:161
void setPcnwSamples(const Container &SwValues, const Container &values)
Set the sampling points for the capillary pressure curve.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:143
void setKrnSamples(const Container &SwValues, const Container &values)
Set the sampling points for the relative permeability curve of the non-wetting phase.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:205
OPM_HOST_DEVICE const ValueVector & SwKrnSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:111
void finalize()
Calculate all dependent quantities once the independent quantities of the parameter object have been ...
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:83
void setKrwSamples(const Container &SwValues, const Container &values)
Set the sampling points for the relative permeability curve of the wetting phase.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:174
OPM_HOST_DEVICE const ValueVector & pcwnSamples() const
Return the sampling points for the capillary pressure curve.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:131
OPM_HOST_DEVICE const ValueVector & krnSamples() const
Return the sampling points for the relative permeability curve of the non-wetting phase.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:192
OPM_HOST_DEVICE const ValueVector & SwPcwnSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:120
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30