3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
initializationindicator.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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program 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 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program 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 this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_GRIDADAPTINITIALIZATIONINDICATOR_HH
25#define DUMUX_GRIDADAPTINITIALIZATIONINDICATOR_HH
26
27#include <memory>
28
29#include <dune/geometry/type.hh>
34
35namespace Dumux {
36
44template<class TypeTag>
46{
50 using Element = typename GridView::Traits::template Codim<0>::Entity;
51
54
56
57public:
58
72 GridAdaptInitializationIndicator(std::shared_ptr<const Problem> problem,
73 std::shared_ptr<const GridGeometry> gridGeometry,
74 std::shared_ptr<const GridVariables> gridVariables)
75 : problem_(problem)
76 , gridGeometry_(gridGeometry)
77 , gridVariables_(gridVariables)
78 , minLevel_(getParamFromGroup<int>(problem->paramGroup(), "Adaptive.MinLevel"))
79 , maxLevel_(getParamFromGroup<int>(problem->paramGroup(), "Adaptive.MaxLevel"))
80 , refineAtDirichletBC_(getParamFromGroup<bool>(problem->paramGroup(), "Adaptive.RefineAtDirichletBC", true))
81 , refineAtFluxBC_(getParamFromGroup<bool>(problem->paramGroup(), "Adaptive.RefineAtFluxBC", true))
82 , refineAtSource_(getParamFromGroup<bool>(problem->paramGroup(), "Adaptive.RefineAtSource", true))
83 , eps_(getParamFromGroup<Scalar>(problem->paramGroup(), "Adaptive.BCRefinementThreshold", 1e-10))
84 {}
85
89 void setMinLevel(std::size_t minLevel)
90 {
91 minLevel_ = minLevel;
92 }
93
97 void setMaxLevel(std::size_t maxLevel)
98 {
99 maxLevel_ = maxLevel;
100 }
101
105 void setLevels(std::size_t minLevel, std::size_t maxLevel)
106 {
107 minLevel_ = minLevel;
108 maxLevel_ = maxLevel;
109 }
110
115 {
116 refineAtDirichletBC_ = refine;
117 }
118
123 {
124 refineAtSource_ = refine;
125 }
126
131 {
132 refineAtFluxBC_ = refine;
133 }
134
140 template<class SolutionVector>
141 void calculate(const SolutionVector& sol)
142 {
144 indicatorVector_.assign(gridGeometry_->gridView().size(0), false);
145
146 for (const auto& element : elements(gridGeometry_->gridView()))
147 {
148 const auto eIdx = gridGeometry_->elementMapper().index(element);
149
151 if (element.level() < minLevel_)
152 {
153 indicatorVector_[eIdx] = true;
154 continue; // proceed to the next element
155 }
156
157 // If refinement at sources/BCs etc is deactivated, skip the rest
158 if (!refineAtSource_ && !refineAtFluxBC_ && !refineAtDirichletBC_)
159 continue;
160
161 // if the element is already on the maximum permissive level, skip rest
162 if (element.level() == maxLevel_)
163 continue;
164
165 // get the fvGeometry and elementVolVars needed for the bc and source interfaces
166 auto fvGeometry = localView(*gridGeometry_);
167 fvGeometry.bind(element);
168
169 auto elemVolVars = localView(gridVariables_->curGridVolVars());
170 elemVolVars.bind(element, fvGeometry, sol);
171
172 // elemFluxVarsCache for neumann interface
173 auto elemFluxVarsCache = localView(gridVariables_->gridFluxVarsCache());
174 elemFluxVarsCache.bind(element, fvGeometry, elemVolVars);
175
177 if (refineAtSource_)
178 {
179 for (const auto& scv : scvs(fvGeometry))
180 {
181 auto source = problem_->source(element, fvGeometry, elemVolVars, scv);
182 auto pointSource = problem_->scvPointSources(element, fvGeometry, elemVolVars, scv);
183 if (source.infinity_norm() + pointSource.infinity_norm() > eps_)
184 {
185 indicatorVector_[eIdx] = true;
186 break; // element is marked, escape scv loop
187 }
188 }
189 }
190
192 if (!indicatorVector_[eIdx] // proceed if element is not already marked
193 && element.hasBoundaryIntersections() // proceed if element is on boundary
194 && (refineAtDirichletBC_ || refineAtFluxBC_)) // proceed if boundary refinement is active
195 {
196 // cell-centered schemes
197 if (!isBox)
198 {
199 for (const auto& scvf : scvfs(fvGeometry))
200 {
201 // skip non-boundary scvfs
202 if (!scvf.boundary())
203 continue;
204
205 const auto bcTypes = problem_->boundaryTypes(element, scvf);
206 // We assume pure BCs, mixed boundary types are not allowed anyway!
207 if(bcTypes.hasOnlyDirichlet() && refineAtDirichletBC_)
208 {
209 indicatorVector_[eIdx] = true;
210 break; // element is marked, escape scvf loop
211 }
212
213 // we are on a pure Neumann boundary
214 else if(refineAtFluxBC_)
215 {
216 const auto fluxes = problem_->neumann(element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
217 if (fluxes.infinity_norm() > eps_)
218 {
219 indicatorVector_[eIdx] = true;
220 break; // element is marked, escape scvf loop
221 }
222 }
223 }
224 }
225 // box-scheme
226 else
227 {
228 // container to store bcTypes
229 using BoundaryTypes = typename ProblemTraits<Problem>::BoundaryTypes;
230 std::vector<BoundaryTypes> bcTypes(fvGeometry.numScv());
231
232 // Get bcTypes and maybe mark for refinement on Dirichlet boundaries
233 for (const auto& scv : scvs(fvGeometry))
234 {
235 bcTypes[scv.localDofIndex()] = problem_->boundaryTypes(element, scv);
236 if (refineAtDirichletBC_ && bcTypes[scv.localDofIndex()].hasDirichlet())
237 {
238 indicatorVector_[eIdx] = true;
239 break; // element is marked, escape scv loop
240 }
241 }
242
243 // If element hasn't been marked because of Dirichlet BCS, check Neumann BCs
244 if (!indicatorVector_[eIdx] && refineAtFluxBC_)
245 {
246 for (const auto& scvf : scvfs(fvGeometry))
247 {
249 if (scvf.boundary() && bcTypes[scvf.insideScvIdx()].hasNeumann())
250 {
251 const auto fluxes = problem_->neumann(element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
252 if (fluxes.infinity_norm() > eps_)
253 {
254 indicatorVector_[eIdx] = true;
255 break; // element is marked, escape scvf loop
256 }
257 }
258 }
259 }
260 }
261 }
262 }
263 }
264
276 int operator() (const Element& element) const
277 {
278 if (indicatorVector_[gridGeometry_->elementMapper().index(element)])
279 return 1;
280 return 0;
281 }
282
283private:
284 std::shared_ptr<const Problem> problem_;
285 std::shared_ptr<const GridGeometry> gridGeometry_;
286 std::shared_ptr<const GridVariables> gridVariables_;
287 std::vector<bool> indicatorVector_;
288
289 int minLevel_;
290 int maxLevel_;
291 bool refineAtDirichletBC_;
292 bool refineAtFluxBC_;
293 bool refineAtSource_;
294 Scalar eps_;
295};
296
297} // end namespace Dumux
298
299#endif
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
The available discretization methods in Dumux.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition: parameters.hh:374
Definition: adapt.hh:29
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:149
Class defining an initialization indicator for grid adaption. Refines at sources and boundaries....
Definition: initializationindicator.hh:46
void setRefinementAtFluxBoundaries(bool refine)
Function to set if refinement at sources is to be used.
Definition: initializationindicator.hh:130
void setLevels(std::size_t minLevel, std::size_t maxLevel)
Function to set the minumum/maximum allowed levels.
Definition: initializationindicator.hh:105
void setMaxLevel(std::size_t maxLevel)
Function to set the maximum allowed level.
Definition: initializationindicator.hh:97
void calculate(const SolutionVector &sol)
Calculates the indicator used for refinement/coarsening for each grid cell.
Definition: initializationindicator.hh:141
void setRefinementAtDirichletBC(bool refine)
Function to set if refinement at Dirichlet boundaries is to be used.
Definition: initializationindicator.hh:114
int operator()(const Element &element) const
function call operator to return mark
Definition: initializationindicator.hh:276
GridAdaptInitializationIndicator(std::shared_ptr< const Problem > problem, std::shared_ptr< const GridGeometry > gridGeometry, std::shared_ptr< const GridVariables > gridVariables)
Constructor.
Definition: initializationindicator.hh:72
bool refine(const Element &element)
Indicator function for marking of grid cells for refinement.
Definition: gridadaptinitializationindicator.hh:316
int maxLevel()
Definition: gridadaptinitializationindicator.hh:346
void setRefinementAtSources(bool refine)
Function to set if refinement at sources is to be used.
Definition: initializationindicator.hh:122
void setMinLevel(std::size_t minLevel)
Function to set the minimum allowed level.
Definition: initializationindicator.hh:89
typename Detail::template ProblemTraits< Problem, GridGeometry::discMethod >::BoundaryTypes BoundaryTypes
Definition: common/typetraits/problem.hh:46
Declares all properties used in Dumux.
Type traits for problem classes.