My Project
BlackoilAquiferModel.hpp
1 /*
2  File adapted from BlackoilWellModel.hpp
3 
4  Copyright 2017 TNO - Heat Transfer & Fluid Dynamics, Modelling & Optimization of the Subsurface
5  Copyright 2017 Statoil ASA.
6 
7  This file is part of the Open Porous Media project (OPM).
8 
9  OPM is free software: you can redistribute it and/or modify
10  it under the terms of the GNU General Public License as published by
11  the Free Software Foundation, either version 3 of the License, or
12  (at your option) any later version.
13 
14  OPM is distributed in the hope that it will be useful,
15  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  GNU General Public License for more details.
18 
19  You should have received a copy of the GNU General Public License
20  along with OPM. If not, see <http://www.gnu.org/licenses/>.
21 */
22 
23 
24 #ifndef OPM_BLACKOILAQUIFERMODEL_HEADER_INCLUDED
25 #define OPM_BLACKOILAQUIFERMODEL_HEADER_INCLUDED
26 
27 #include <ebos/eclbaseaquifermodel.hh>
28 
29 #include <opm/input/eclipse/EclipseState/Aquifer/Aquancon.hpp>
30 #include <opm/input/eclipse/EclipseState/Aquifer/AquiferCT.hpp>
31 #include <opm/input/eclipse/EclipseState/Aquifer/Aquifetp.hpp>
32 
33 #include <opm/output/data/Aquifer.hpp>
34 
35 #include <opm/simulators/aquifers/AquiferCarterTracy.hpp>
36 #include <opm/simulators/aquifers/AquiferFetkovich.hpp>
37 #include <opm/simulators/aquifers/AquiferNumerical.hpp>
38 
39 #include <opm/grid/CpGrid.hpp>
40 #include <opm/grid/polyhedralgrid.hh>
41 #if HAVE_DUNE_ALUGRID
42 #include <dune/alugrid/grid.hh>
43 #endif
44 
45 #include <opm/material/densead/Math.hpp>
46 
47 #include <vector>
48 #include <type_traits>
49 
50 namespace Opm
51 {
52 
53 template<class Grid>
55  : public std::bool_constant<false>
56 {};
57 
58 
59 template<>
60 class SupportsFaceTag<Dune::CpGrid>
61  : public std::bool_constant<true>
62 {};
63 
64 
65 template<>
66 class SupportsFaceTag<Dune::PolyhedralGrid<3, 3>>
67  : public std::bool_constant<true>
68 {};
69 
70 #if HAVE_DUNE_ALUGRID
71 template<>
72 class SupportsFaceTag<Dune::ALUGrid<3, 3, Dune::cube, Dune::nonconforming>>
73  : public std::bool_constant<true>
74 {};
75 #endif
76 
77 
79 template <typename TypeTag>
81 {
82  using Simulator = GetPropType<TypeTag, Properties::Simulator>;
83  using RateVector = GetPropType<TypeTag, Properties::RateVector>;
84 
85 
86 public:
87  explicit BlackoilAquiferModel(Simulator& simulator);
88 
89  void initialSolutionApplied();
90  void initFromRestart(const data::Aquifers& aquiferSoln);
91 
92  void beginEpisode();
93  void beginTimeStep();
94  void beginIteration();
95  // add the water rate due to aquifers to the source term.
96  template <class Context>
97  void addToSource(RateVector& rates, const Context& context, unsigned spaceIdx, unsigned timeIdx) const;
98  void endIteration();
99  void endTimeStep();
100  void endEpisode();
101 
102  data::Aquifers aquiferData() const;
103 
104  template <class Restarter>
105  void serialize(Restarter& res);
106 
107  template <class Restarter>
108  void deserialize(Restarter& res);
109 
110 protected:
111  // --------- Types ---------
112  using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
113  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
114 
117 
118  Simulator& simulator_;
119 
120  // TODO: probably we can use one variable to store both types of aquifers, because
121  // they share the base class
122  mutable std::vector<AquiferCarterTracy_object> aquifers_CarterTracy;
123  mutable std::vector<AquiferFetkovich_object> aquifers_Fetkovich;
124  std::vector<AquiferNumerical<TypeTag>> aquifers_numerical;
125 
126  // This initialization function is used to connect the parser objects with the ones needed by AquiferCarterTracy
127  void init();
128 
129  bool aquiferActive() const;
130  bool aquiferCarterTracyActive() const;
131  bool aquiferFetkovichActive() const;
132  bool aquiferNumericalActive() const;
133 };
134 
135 
136 } // namespace Opm
137 
138 #include "BlackoilAquiferModel_impl.hpp"
139 #endif
Definition: AquiferCarterTracy.hpp:38
Definition: AquiferFetkovich.hpp:37
Class for handling the blackoil well model.
Definition: BlackoilAquiferModel.hpp:81
Definition: BlackoilAquiferModel.hpp:56
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27