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#ifdef USE_POLYHEDRALGRID
41#include <opm/grid/polyhedralgrid.hh>
42#endif
43#if HAVE_DUNE_ALUGRID
44#include <dune/alugrid/grid.hh>
45#endif
46
47#include <opm/material/densead/Math.hpp>
48
49#include <vector>
50#include <type_traits>
51
52namespace Opm
53{
54
55template<class Grid>
57 : public std::bool_constant<false>
58{};
59
60
61template<>
62class SupportsFaceTag<Dune::CpGrid>
63 : public std::bool_constant<true>
64{};
65
66
67#ifdef USE_POLYHEDRALGRID
68template<>
69class SupportsFaceTag<Dune::PolyhedralGrid<3, 3>>
70 : public std::bool_constant<true>
71{};
72#endif
73
74#if HAVE_DUNE_ALUGRID
75template<>
76class SupportsFaceTag<Dune::ALUGrid<3, 3, Dune::cube, Dune::nonconforming>>
77 : public std::bool_constant<true>
78{};
79#endif
80
81
83template <typename TypeTag>
85{
86 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
87 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
88
89
90public:
91 explicit BlackoilAquiferModel(Simulator& simulator);
92
93 void initialSolutionApplied();
94 void initFromRestart(const data::Aquifers& aquiferSoln);
95
96 void beginEpisode();
97 void beginTimeStep();
98 void beginIteration();
99 // add the water rate due to aquifers to the source term.
100 template <class Context>
101 void addToSource(RateVector& rates, const Context& context, unsigned spaceIdx, unsigned timeIdx) const;
102 void addToSource(RateVector& rates, unsigned globalSpaceIdx, unsigned timeIdx) const;
103 void endIteration();
104 void endTimeStep();
105 void endEpisode();
106
107 data::Aquifers aquiferData() const;
108
109 template <class Restarter>
110 void serialize(Restarter& res);
111
112 template <class Restarter>
113 void deserialize(Restarter& res);
114
115protected:
116 // --------- Types ---------
117 using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
118 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
119
120 Simulator& simulator_;
121
122 std::vector<std::unique_ptr<AquiferInterface<TypeTag>>> aquifers;
123
124 // This initialization function is used to connect the parser objects with the ones needed by AquiferCarterTracy
125 void init();
126};
127
128
129} // namespace Opm
130
131#include "BlackoilAquiferModel_impl.hpp"
132
133#endif
Class for handling the blackoil well model.
Definition: BlackoilAquiferModel.hpp:85
Definition: BlackoilAquiferModel.hpp:58
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27