27 #ifndef OPM_PIECEWISE_LINEAR_TWO_PHASE_MATERIAL_PARAMS_HPP
28 #define OPM_PIECEWISE_LINEAR_TWO_PHASE_MATERIAL_PARAMS_HPP
43 template<
class TraitsT>
46 using Scalar =
typename TraitsT::Scalar;
49 using ValueVector = std::vector<Scalar>;
51 using Traits = TraitsT;
67 if (SwPcwnSamples_.front() > SwPcwnSamples_.back())
68 swapOrder_(SwPcwnSamples_, pcwnSamples_);
70 if (SwKrwSamples_.front() > SwKrwSamples_.back())
71 swapOrder_(SwKrwSamples_, krwSamples_);
74 if (SwKrnSamples_.front() > SwKrnSamples_.back())
75 swapOrder_(SwKrnSamples_, krnSamples_);
83 { EnsureFinalized::check();
return SwKrwSamples_; }
89 { EnsureFinalized::check();
return SwKrnSamples_; }
95 { EnsureFinalized::check();
return SwPcwnSamples_; }
103 { EnsureFinalized::check();
return pcwnSamples_; }
110 template <
class Container>
113 assert(SwValues.size() == values.size());
115 size_t n = SwValues.size();
116 SwPcwnSamples_.resize(n);
117 pcwnSamples_.resize(n);
119 std::copy(SwValues.begin(), SwValues.end(), SwPcwnSamples_.begin());
120 std::copy(values.begin(), values.end(), pcwnSamples_.begin());
130 { EnsureFinalized::check();
return krwSamples_; }
138 template <
class Container>
141 assert(SwValues.size() == values.size());
143 size_t n = SwValues.size();
144 SwKrwSamples_.resize(n);
145 krwSamples_.resize(n);
147 std::copy(SwValues.begin(), SwValues.end(), SwKrwSamples_.begin());
148 std::copy(values.begin(), values.end(), krwSamples_.begin());
158 { EnsureFinalized::check();
return krnSamples_; }
166 template <
class Container>
169 assert(SwValues.size() == values.size());
171 size_t n = SwValues.size();
172 SwKrnSamples_.resize(n);
173 krnSamples_.resize(n);
175 std::copy(SwValues.begin(), SwValues.end(), SwKrnSamples_.begin());
176 std::copy(values.begin(), values.end(), krnSamples_.begin());
180 void swapOrder_(ValueVector& swValues, ValueVector& values)
const
182 if (swValues.front() > values.back()) {
183 for (
unsigned origSampleIdx = 0;
184 origSampleIdx < swValues.size() / 2;
187 size_t newSampleIdx = swValues.size() - origSampleIdx - 1;
189 std::swap(swValues[origSampleIdx], swValues[newSampleIdx]);
190 std::swap(values[origSampleIdx], values[newSampleIdx]);
195 ValueVector SwPcwnSamples_;
196 ValueVector SwKrwSamples_;
197 ValueVector SwKrnSamples_;
198 ValueVector pcwnSamples_;
199 ValueVector krwSamples_;
200 ValueVector krnSamples_;
Default implementation for asserting finalization of parameter objects.
Default implementation for asserting finalization of parameter objects.
Definition: EnsureFinalized.hpp:47
void finalize()
Mark the object as finalized.
Definition: EnsureFinalized.hpp:75
Specification of the material parameters for a two-phase material law which uses a table and piecewis...
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:45
const ValueVector & SwPcwnSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:94
void finalize()
Calculate all dependent quantities once the independent quantities of the parameter object have been ...
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:61
void setKrwSamples(const Container &SwValues, const Container &values)
Set the sampling points for the relative permeability curve of the wetting phase.
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:139
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:167
const ValueVector & SwKrnSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:88
const ValueVector & krnSamples() const
Return the sampling points for the relative permeability curve of the non-wetting phase.
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:157
void setPcnwSamples(const Container &SwValues, const Container &values)
Set the sampling points for the capillary pressure curve.
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:111
const ValueVector & krwSamples() const
Return the sampling points for the relative permeability curve of the wetting phase.
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:129
const ValueVector & SwKrwSamples() const
Return the wetting-phase saturation values of all sampling points.
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:82
const ValueVector & pcnwSamples() const
Return the sampling points for the capillary pressure curve.
Definition: PiecewiseLinearTwoPhaseMaterialParams.hpp:102