fvbasediscretizationfemadapt.hh
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*/
28#ifndef EWOMS_FV_BASE_DISCRETIZATION_FEMADAPT_HH
29#define EWOMS_FV_BASE_DISCRETIZATION_FEMADAPT_HH
30
32
33#include <dune/fem/space/common/adaptationmanager.hh>
34#include <dune/fem/space/common/restrictprolongtuple.hh>
35#include <dune/fem/function/blockvectorfunction.hh>
36#include <dune/fem/misc/capabilities.hh>
37
38namespace Opm {
39
40template<class TypeTag>
41class FvBaseDiscretizationFemAdapt;
42
43namespace Properties {
44
45template<class TypeTag>
46struct BaseDiscretizationType<TypeTag,TTag::FvBaseDiscretization> {
48};
49
50template<class TypeTag>
51struct DiscreteFunction<TypeTag, TTag::FvBaseDiscretization> {
54 using type = Dune::Fem::ISTLBlockVectorDiscreteFunction<DiscreteFunctionSpace, PrimaryVariables>;
55};
56
57} // namespace Properties
58
65template <class TypeTag>
67{
73
74 static constexpr unsigned historySize = getPropValue<TypeTag, Properties::TimeDiscHistorySize>();
75
77 // discrete function storing solution data
78 using DiscreteFunction = Dune::Fem::ISTLBlockVectorDiscreteFunction<DiscreteFunctionSpace, PrimaryVariables>;
79
80 // problem restriction and prolongation operator for adaptation
81 using ProblemRestrictProlongOperator = typename Problem::RestrictProlongOperator;
82
83 // discrete function restriction and prolongation operator for adaptation
84 using DiscreteFunctionRestrictProlong = Dune::Fem::RestrictProlongDefault<DiscreteFunction>;
85 using RestrictProlong
86 = Dune::Fem::RestrictProlongTuple<DiscreteFunctionRestrictProlong, ProblemRestrictProlongOperator>;
87
88 // adaptation classes
89 using AdaptationManager = Dune::Fem::AdaptationManager<Grid, RestrictProlong>;
90
91public:
92 template<class Serializer>
94 template<class SolutionType>
95 static void serializeOp(Serializer& serializer,
96 SolutionType& solution)
97 {
98 for (auto& sol : solution) {
99 serializer(sol->blockVector());
100 }
101 }
102 };
103
104 FvBaseDiscretizationFemAdapt(Simulator& simulator)
105 : ParentType(simulator)
106 , space_(simulator.vanguard().gridPart())
107 {
108 if (this->enableGridAdaptation_ && !Dune::Fem::Capabilities::isLocallyAdaptive<Grid>::v) {
109 throw std::invalid_argument("Grid adaptation enabled, but chosen Grid is not capable"
110 " of adaptivity");
111 }
112
113 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
114 this->solution_[timeIdx] = std::make_unique<DiscreteFunction>("solution", space_);
115 }
116 }
117
119 {
120 // adapt the grid if enabled and if all dependencies are available
121 // adaptation is only done if markForGridAdaptation returns true
122 if (this->enableGridAdaptation_) {
123 // check if problem allows for adaptation and cells were marked
124 if (this->simulator_.problem().markForGridAdaptation()) {
125 // adapt the grid and load balance if necessary
126 adaptationManager().adapt();
127
128 // if the grid has potentially changed, we need to re-create the
129 // supporting data structures.
130 this->elementMapper_.update(this->gridView_);
131 this->vertexMapper_.update(this->gridView_);
132 this->resetLinearizer();
133 // this is a bit hacky because it supposes that Problem::finishInit()
134 // works fine multiple times in a row.
135 //
136 // TODO: move this to Problem::gridChanged()
137 this->finishInit();
138
139 // notify the problem that the grid has changed
140 //
141 // TODO: come up with a mechanism to access the unadapted data structures
142 // outside of the problem (i.e., grid, mappers, solutions)
143 this->simulator_.problem().gridChanged();
144
145 // notify the modules for visualization output
146 auto outIt = this->outputModules_.begin();
147 auto outEndIt = this->outputModules_.end();
148 for (; outIt != outEndIt; ++outIt)
149 (*outIt)->allocBuffers();
150 }
151 }
152 }
153
154 AdaptationManager& adaptationManager()
155 {
156 if (!adaptationManager_) {
157 // create adaptation objects here, because when doing so in constructor
158 // problem is not yet intialized, aka seg fault
159 restrictProlong_ = std::make_unique<RestrictProlong>(DiscreteFunctionRestrictProlong(*(this->solution_[/*timeIdx=*/0])),
160 this->simulator_.problem().restrictProlongOperator());
161 adaptationManager_ = std::make_unique<AdaptationManager>(this->simulator_.vanguard().grid(), *restrictProlong_);
162 }
163 return *adaptationManager_;
164 }
165
166private:
167 DiscreteFunctionSpace space_;
168 std::unique_ptr<RestrictProlong> restrictProlong_;
169 std::unique_ptr<AdaptationManager> adaptationManager_;
170};
171
172} // namespace Opm
173
174#endif
The base class for the finite volume discretization schemes.
Definition: fvbasediscretizationfemadapt.hh:67
AdaptationManager & adaptationManager()
Definition: fvbasediscretizationfemadapt.hh:154
FvBaseDiscretizationFemAdapt(Simulator &simulator)
Definition: fvbasediscretizationfemadapt.hh:104
void adaptGrid()
Definition: fvbasediscretizationfemadapt.hh:118
The base class for the finite volume discretization schemes without adaptation.
Definition: fvbasediscretization.hh:1947
The base class for the finite volume discretization schemes.
Definition: fvbasediscretization.hh:294
void finishInit()
Apply the initial conditions to the model.
Definition: fvbasediscretization.hh:473
std::list< BaseOutputModule< TypeTag > * > outputModules_
Definition: fvbasediscretization.hh:1926
std::array< std::unique_ptr< DiscreteFunction >, historySize > solution_
Definition: fvbasediscretization.hh:1924
bool enableGridAdaptation_
Definition: fvbasediscretization.hh:1934
ElementMapper elementMapper_
Definition: fvbasediscretization.hh:1899
GridView gridView_
Definition: fvbasediscretization.hh:1896
const SolutionVector & solution(unsigned timeIdx) const
Reference to the solution at a given history index as a block vector.
Definition: fvbasediscretization.hh:1137
Simulator & simulator_
Definition: fvbasediscretization.hh:1893
VertexMapper vertexMapper_
Definition: fvbasediscretization.hh:1900
void resetLinearizer()
Resets the Jacobian matrix linearizer, so that the boundary types can be altered.
Definition: fvbasediscretization.hh:1527
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:235
Definition: fvbasediscretizationfemadapt.hh:93
static void serializeOp(Serializer &serializer, SolutionType &solution)
Definition: fvbasediscretizationfemadapt.hh:95
GetPropType< TypeTag, Properties::DiscreteFunctionSpace > DiscreteFunctionSpace
Definition: fvbasediscretizationfemadapt.hh:52
GetPropType< TypeTag, Properties::PrimaryVariables > PrimaryVariables
Definition: fvbasediscretizationfemadapt.hh:53
typename BaseDiscretization::BlockVectorWrapper type
Definition: fvbasediscretization.hh:279
Definition: fvbaseproperties.hh:77