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>
43template <
class TraitsT,
class VectorT = std::vector<
typename TraitsT::Scalar>>
49 template <
class ViewType,
class TraitsT,
class ContainerType>
60template <
class TraitsT,
class VectorT = std::vector<
typename TraitsT::Scalar>>
63 using Scalar =
typename TraitsT::Scalar;
66 using ValueVector = VectorT;
68 using Traits = TraitsT;
87 if constexpr (std::is_same_v<ValueVector, std::vector<Scalar>>){
106 if (SwPcwnSamples_.front() > SwPcwnSamples_.back())
107 swapOrderIfPossibleThrowOtherwise_(SwPcwnSamples_, pcwnSamples_);
109 if (SwKrwSamples_.front() > SwKrwSamples_.back())
110 swapOrderIfPossibleThrowOtherwise_(SwKrwSamples_, krwSamples_);
112 if (SwKrnSamples_.front() > SwKrnSamples_.back())
113 swapOrderIfPossibleThrowOtherwise_(SwKrnSamples_, krnSamples_);
121 EnsureFinalized::check();
129 EnsureFinalized::check();
130 return SwKrwSamples_;
138 EnsureFinalized::check();
139 return SwKrnSamples_;
147 EnsureFinalized::check();
148 return SwPcwnSamples_;
158 EnsureFinalized::check();
167 template <
class Container>
170 assert(SwValues.size() == values.size());
172 size_t n = SwValues.size();
173 SwPcwnSamples_.resize(n);
174 pcwnSamples_.resize(n);
176 std::copy(SwValues.begin(), SwValues.end(), SwPcwnSamples_.begin());
177 std::copy(values.begin(), values.end(), pcwnSamples_.begin());
188 EnsureFinalized::check();
198 template <
class Container>
201 assert(SwValues.size() == values.size());
203 size_t n = SwValues.size();
204 SwKrwSamples_.resize(n);
205 krwSamples_.resize(n);
207 std::copy(SwValues.begin(), SwValues.end(), SwKrwSamples_.begin());
208 std::copy(values.begin(), values.end(), krwSamples_.begin());
219 EnsureFinalized::check();
229 template <
class Container>
232 assert(SwValues.size() == values.size());
234 size_t n = SwValues.size();
235 SwKrnSamples_.resize(n);
236 krnSamples_.resize(n);
238 std::copy(SwValues.begin(), SwValues.end(), SwKrnSamples_.begin());
239 std::copy(values.begin(), values.end(), krnSamples_.begin());
243 void swapOrderIfPossibleThrowOtherwise_(ValueVector& swValues, ValueVector& values)
const
247 if (swValues.front() > values.back()) {
250 if constexpr (!std::is_const_v<typename ValueVector::value_type> && !std::is_const_v<ValueVector>) {
251 for (
unsigned origSampleIdx = 0; origSampleIdx < swValues.size() / 2; ++origSampleIdx) {
252 size_t newSampleIdx = swValues.size() - origSampleIdx - 1;
254 std::swap(swValues[origSampleIdx], swValues[newSampleIdx]);
255 std::swap(values[origSampleIdx], values[newSampleIdx]);
259 OPM_THROW(std::logic_error,
"Saturation values in interpolation table provided in wrong order, but table is immutable");
264 template <
class ViewType,
class Traits,
class Container>
267 ValueVector SwPcwnSamples_;
268 ValueVector SwKrwSamples_;
269 ValueVector SwKrnSamples_;
270 ValueVector pcwnSamples_;
271 ValueVector krwSamples_;
272 ValueVector krnSamples_;
276namespace Opm::gpuistl{
284template <
class GPUContainerType,
class TraitsT>
288 params.checkFinalized();
290 auto SwPcwnSamples = GPUContainerType(params.SwPcwnSamples());
291 auto pcwnSamples = GPUContainerType(params.pcwnSamples());
292 auto SwKrwSamples = GPUContainerType(params.SwKrwSamples());
293 auto krwSamples = GPUContainerType(params.krwSamples());
294 auto SwKrnSamples = GPUContainerType(params.SwKrnSamples());
295 auto krnSamples = GPUContainerType(params.krnSamples());
311template <
class ViewType,
class TraitsT,
class ContainerType>
317 using ContainedType =
typename ViewType::value_type;
319 ViewType SwPcwnSamples = make_view<ContainedType>(params.SwPcwnSamples_);
320 ViewType pcwnSamples = make_view<ContainedType>(params.pcwnSamples_);
321 ViewType SwKrwSamples = make_view<ContainedType>(params.SwKrwSamples_);
322 ViewType krwSamples = make_view<ContainedType>(params.krwSamples_);
323 ViewType SwKrnSamples = make_view<ContainedType>(params.SwKrnSamples_);
324 ViewType krnSamples = make_view<ContainedType>(params.krnSamples_);
Default implementation for asserting finalization of parameter objects.
Default implementation for asserting finalization of parameter objects.
Definition EnsureFinalized.hpp:49
OPM_HOST_DEVICE void finalize()
Mark the object as finalized.
Definition EnsureFinalized.hpp:82
Specification of the material parameters for a two-phase material law which uses a table and piecewis...
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:62
OPM_HOST_DEVICE const ValueVector & SwKrwSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:127
OPM_HOST_DEVICE const ValueVector & krwSamples() const
Return the sampling points for the relative permeability curve of the wetting phase.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:186
void setPcnwSamples(const Container &SwValues, const Container &values)
Set the sampling points for the capillary pressure curve.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:168
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:230
OPM_HOST_DEVICE const ValueVector & SwKrnSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:136
void finalize()
Calculate all dependent quantities once the independent quantities of the parameter object have been ...
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:100
void setKrwSamples(const Container &SwValues, const Container &values)
Set the sampling points for the relative permeability curve of the wetting phase.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:199
OPM_HOST_DEVICE const ValueVector & pcwnSamples() const
Return the sampling points for the capillary pressure curve.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:156
void checkFinalized() const
Check if the parameter object has been finalized.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:119
OPM_HOST_DEVICE const ValueVector & krnSamples() const
Return the sampling points for the relative permeability curve of the non-wetting phase.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:217
OPM_HOST_DEVICE const ValueVector & SwPcwnSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:145
Definition PiecewiseLinearTwoPhaseMaterialParams.hpp:44
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30