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
42// forward declaration of the class so the function in the next namespace can be declared
43template <class TraitsT, class VectorT = std::vector<typename TraitsT::Scalar>>
45
46// declaration of make_view in correct namespace so friend function can be declared in the class
47namespace Opm::gpuistl
48{
49 template <class ViewType, class TraitsT, class ContainerType>
51}
52
53namespace Opm {
60template <class TraitsT, class VectorT = std::vector<typename TraitsT::Scalar>>
62{
63 using Scalar = typename TraitsT::Scalar;
64
65public:
66 using ValueVector = VectorT;
67
68 using Traits = TraitsT;
69
71 {
72 }
73
75 ValueVector pcwnSamples,
76 ValueVector SwKrwSamples,
77 ValueVector krwSamples,
78 ValueVector SwKrnSamples,
79 ValueVector krnSamples)
80 : SwPcwnSamples_(SwPcwnSamples)
81 , SwKrwSamples_(SwKrwSamples)
82 , SwKrnSamples_(SwKrnSamples)
83 , pcwnSamples_(pcwnSamples)
84 , krwSamples_(krwSamples)
85 , krnSamples_(krnSamples)
86 {
87 if constexpr (std::is_same_v<ValueVector, std::vector<Scalar>>){
88 finalize();
89 }
90 else{
91 // safe if we have a GPU type instantiated by copy_to_gpu or make_view
93 }
94 }
95
100 void finalize()
101 {
103
104 // revert the order of the sampling points if they were given
105 // in reverse direction.
106 if (SwPcwnSamples_.front() > SwPcwnSamples_.back())
107 swapOrderIfPossibleThrowOtherwise_(SwPcwnSamples_, pcwnSamples_);
108
109 if (SwKrwSamples_.front() > SwKrwSamples_.back())
110 swapOrderIfPossibleThrowOtherwise_(SwKrwSamples_, krwSamples_);
111
112 if (SwKrnSamples_.front() > SwKrnSamples_.back())
113 swapOrderIfPossibleThrowOtherwise_(SwKrnSamples_, krnSamples_);
114 }
115
119 void checkFinalized() const
120 {
121 EnsureFinalized::check();
122 }
123
127 OPM_HOST_DEVICE const ValueVector& SwKrwSamples() const
128 {
129 EnsureFinalized::check();
130 return SwKrwSamples_;
131 }
132
136 OPM_HOST_DEVICE const ValueVector& SwKrnSamples() const
137 {
138 EnsureFinalized::check();
139 return SwKrnSamples_;
140 }
141
145 OPM_HOST_DEVICE const ValueVector& SwPcwnSamples() const
146 {
147 EnsureFinalized::check();
148 return SwPcwnSamples_;
149 }
150
156 OPM_HOST_DEVICE const ValueVector& pcwnSamples() const
157 {
158 EnsureFinalized::check();
159 return pcwnSamples_;
160 }
161
167 template <class Container>
168 void setPcnwSamples(const Container& SwValues, const Container& values)
169 {
170 assert(SwValues.size() == values.size());
171
172 size_t n = SwValues.size();
173 SwPcwnSamples_.resize(n);
174 pcwnSamples_.resize(n);
175
176 std::copy(SwValues.begin(), SwValues.end(), SwPcwnSamples_.begin());
177 std::copy(values.begin(), values.end(), pcwnSamples_.begin());
178 }
179
186 OPM_HOST_DEVICE const ValueVector& krwSamples() const
187 {
188 EnsureFinalized::check();
189 return krwSamples_;
190 }
191
198 template <class Container>
199 void setKrwSamples(const Container& SwValues, const Container& values)
200 {
201 assert(SwValues.size() == values.size());
202
203 size_t n = SwValues.size();
204 SwKrwSamples_.resize(n);
205 krwSamples_.resize(n);
206
207 std::copy(SwValues.begin(), SwValues.end(), SwKrwSamples_.begin());
208 std::copy(values.begin(), values.end(), krwSamples_.begin());
209 }
210
217 OPM_HOST_DEVICE const ValueVector& krnSamples() const
218 {
219 EnsureFinalized::check();
220 return krnSamples_;
221 }
222
229 template <class Container>
230 void setKrnSamples(const Container& SwValues, const Container& values)
231 {
232 assert(SwValues.size() == values.size());
233
234 size_t n = SwValues.size();
235 SwKrnSamples_.resize(n);
236 krnSamples_.resize(n);
237
238 std::copy(SwValues.begin(), SwValues.end(), SwKrnSamples_.begin());
239 std::copy(values.begin(), values.end(), krnSamples_.begin());
240 }
241
242private:
243 void swapOrderIfPossibleThrowOtherwise_(ValueVector& swValues, ValueVector& values) const
244 {
245 // TODO: comparing saturation values to the actual values we sample from looks strange
246 // TODO: yet changing to swValues.back() breaks tests
247 if (swValues.front() > values.back()) {
248 // Reverting the order involves swapping which only works for non-consts.
249 // The const expr ensures we can create constant parameter views.
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;
253
254 std::swap(swValues[origSampleIdx], swValues[newSampleIdx]);
255 std::swap(values[origSampleIdx], values[newSampleIdx]);
256 }
257 }
258 else{
259 OPM_THROW(std::logic_error, "Saturation values in interpolation table provided in wrong order, but table is immutable");
260 }
261 }
262 }
263
264 template <class ViewType, class Traits, class Container>
266
267 ValueVector SwPcwnSamples_;
268 ValueVector SwKrwSamples_;
269 ValueVector SwKrnSamples_;
270 ValueVector pcwnSamples_;
271 ValueVector krwSamples_;
272 ValueVector krnSamples_;
273};
274} // namespace Opm
275
276namespace Opm::gpuistl{
277
284template <class GPUContainerType, class TraitsT>
286
287 // only create the GPU object if the CPU object is finalized
288 params.checkFinalized();
289
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());
296
298 pcwnSamples,
299 SwKrwSamples,
300 krwSamples,
301 SwKrnSamples,
302 krnSamples);
303}
304
311template <class ViewType, class TraitsT, class ContainerType>
313
314 // only create the GPU object if the CPU object is finalized
315 params.checkFinalized();
316
317 using ContainedType = typename ViewType::value_type;
318
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_);
325
327 pcwnSamples,
328 SwKrwSamples,
329 krwSamples,
330 SwKrnSamples,
331 krnSamples);
332}
333
334}
335
336#endif
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