My Project
Loading...
Searching...
No Matches
PiecewiseLinearTwoPhaseMaterialParams.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_PIECEWISE_LINEAR_TWO_PHASE_MATERIAL_PARAMS_HPP
28#define OPM_PIECEWISE_LINEAR_TWO_PHASE_MATERIAL_PARAMS_HPP
29
30#include <algorithm>
31#include <cassert>
32#include <vector>
33#include <stdexcept>
34#include <type_traits>
35
36#include <opm/common/ErrorMacros.hpp>
38#include <vector>
39
40#include <opm/common/utility/gpuDecorators.hpp>
41
42namespace Opm {
49template <class TraitsT, class VectorT = std::vector<typename TraitsT::Scalar>>
51{
52 using Scalar = typename TraitsT::Scalar;
53
54public:
55 using ValueVector = VectorT;
56
57 using Traits = TraitsT;
58
60 {
61 }
62
64 ValueVector pcwnSamples,
65 ValueVector SwKrwSamples,
66 ValueVector krwSamples,
67 ValueVector SwKrnSamples,
68 ValueVector krnSamples)
69 : SwPcwnSamples_(SwPcwnSamples)
70 , SwKrwSamples_(SwKrwSamples)
71 , SwKrnSamples_(SwKrnSamples)
72 , pcwnSamples_(pcwnSamples)
73 , krwSamples_(krwSamples)
74 , krnSamples_(krnSamples)
75 {
76 finalize();
77 }
78
83 void finalize()
84 {
85 EnsureFinalized ::finalize();
86
87 // revert the order of the sampling points if they were given
88 // in reverse direction.
89 if (SwPcwnSamples_.front() > SwPcwnSamples_.back())
90 swapOrderIfPossibleThrowOtherwise_(SwPcwnSamples_, pcwnSamples_);
91
92 if (SwKrwSamples_.front() > SwKrwSamples_.back())
93 swapOrderIfPossibleThrowOtherwise_(SwKrwSamples_, krwSamples_);
94
95 if (SwKrnSamples_.front() > SwKrnSamples_.back())
96 swapOrderIfPossibleThrowOtherwise_(SwKrnSamples_, krnSamples_);
97 }
98
102 OPM_HOST_DEVICE const ValueVector& SwKrwSamples() const
103 {
104 EnsureFinalized::check();
105 return SwKrwSamples_;
106 }
107
111 OPM_HOST_DEVICE const ValueVector& SwKrnSamples() const
112 {
113 EnsureFinalized::check();
114 return SwKrnSamples_;
115 }
116
120 OPM_HOST_DEVICE const ValueVector& SwPcwnSamples() const
121 {
122 EnsureFinalized::check();
123 return SwPcwnSamples_;
124 }
125
131 OPM_HOST_DEVICE const ValueVector& pcwnSamples() const
132 {
133 EnsureFinalized::check();
134 return pcwnSamples_;
135 }
136
142 template <class Container>
143 void setPcnwSamples(const Container& SwValues, const Container& values)
144 {
145 assert(SwValues.size() == values.size());
146
147 size_t n = SwValues.size();
148 SwPcwnSamples_.resize(n);
149 pcwnSamples_.resize(n);
150
151 std::copy(SwValues.begin(), SwValues.end(), SwPcwnSamples_.begin());
152 std::copy(values.begin(), values.end(), pcwnSamples_.begin());
153 }
154
161 OPM_HOST_DEVICE const ValueVector& krwSamples() const
162 {
163 EnsureFinalized::check();
164 return krwSamples_;
165 }
166
173 template <class Container>
174 void setKrwSamples(const Container& SwValues, const Container& values)
175 {
176 assert(SwValues.size() == values.size());
177
178 size_t n = SwValues.size();
179 SwKrwSamples_.resize(n);
180 krwSamples_.resize(n);
181
182 std::copy(SwValues.begin(), SwValues.end(), SwKrwSamples_.begin());
183 std::copy(values.begin(), values.end(), krwSamples_.begin());
184 }
185
192 OPM_HOST_DEVICE const ValueVector& krnSamples() const
193 {
194 EnsureFinalized::check();
195 return krnSamples_;
196 }
197
204 template <class Container>
205 void setKrnSamples(const Container& SwValues, const Container& values)
206 {
207 assert(SwValues.size() == values.size());
208
209 size_t n = SwValues.size();
210 SwKrnSamples_.resize(n);
211 krnSamples_.resize(n);
212
213 std::copy(SwValues.begin(), SwValues.end(), SwKrnSamples_.begin());
214 std::copy(values.begin(), values.end(), krnSamples_.begin());
215 }
216
217private:
218 void swapOrderIfPossibleThrowOtherwise_(ValueVector& swValues, ValueVector& values) const
219 {
220 // TODO: comparing saturation values to the actual values we sample from looks strange
221 // TODO: yet changing to swValues.back() breaks tests
222 if (swValues.front() > values.back()) {
223 // Reverting the order involves swapping which only works for non-consts.
224 // The const expr ensures we can create constant parameter views.
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;
228
229 std::swap(swValues[origSampleIdx], swValues[newSampleIdx]);
230 std::swap(values[origSampleIdx], values[newSampleIdx]);
231 }
232 }
233 else{
234 OPM_THROW(std::logic_error, "Saturation values in interpolation table provided in wrong order, but table is immutable");
235 }
236 }
237 }
238
239 ValueVector SwPcwnSamples_;
240 ValueVector SwKrwSamples_;
241 ValueVector SwKrnSamples_;
242 ValueVector pcwnSamples_;
243 ValueVector krwSamples_;
244 ValueVector krnSamples_;
245};
246} // namespace Opm
247
248namespace Opm::gpuistl{
249
256template <class TraitsT, class ContainerType, class ViewType>
257PiecewiseLinearTwoPhaseMaterialParams<TraitsT, ViewType> make_view(const PiecewiseLinearTwoPhaseMaterialParams<TraitsT, ContainerType>& params) {
258
259 using containedType = typename ContainerType::value_type;
260 using viewedTypeNoConst = typename std::remove_const_t<typename ViewType::value_type>;
261
262 static_assert(std::is_same_v<containedType, viewedTypeNoConst>);
263
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());
270
271 return PiecewiseLinearTwoPhaseMaterialParams<TraitsT, ViewType> (SwPcwnSamples,
272 pcwnSamples,
273 SwKrwSamples,
274 krwSamples,
275 SwKrnSamples,
276 krnSamples);
277}
278}
279
280#endif
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