61 std::shared_ptr<const GridGeometry> gridGeometry,
62 std::shared_ptr<const GridVariables> gridVariables)
64 , gridGeometry_(gridGeometry)
65 , gridVariables_(gridVariables)
66 , minLevel_(
getParamFromGroup<std::size_t>(problem->paramGroup(),
"Adaptive.MinLevel"))
67 , maxLevel_(
getParamFromGroup<std::size_t>(problem->paramGroup(),
"Adaptive.MaxLevel"))
68 , refineAtDirichletBC_(
getParamFromGroup<bool>(problem->paramGroup(),
"Adaptive.RefineAtDirichletBC", true))
69 , refineAtFluxBC_(
getParamFromGroup<bool>(problem->paramGroup(),
"Adaptive.RefineAtFluxBC", true))
70 , refineAtSource_(
getParamFromGroup<bool>(problem->paramGroup(),
"Adaptive.RefineAtSource", true))
71 , eps_(
getParamFromGroup<Scalar>(problem->paramGroup(),
"Adaptive.BCRefinementThreshold", 1e-10))
132 indicatorVector_.assign(gridGeometry_->gridView().size(0),
false);
135 auto fvGeometry =
localView(*gridGeometry_);
136 auto elemVolVars =
localView(gridVariables_->curGridVolVars());
138 auto elemFluxVarsCache =
localView(gridVariables_->gridFluxVarsCache());
140 for (
const auto& element : elements(gridGeometry_->gridView()))
142 const auto eIdx = gridGeometry_->elementMapper().index(element);
145 if (element.level() < minLevel_)
147 indicatorVector_[eIdx] =
true;
152 if (!refineAtSource_ && !refineAtFluxBC_ && !refineAtDirichletBC_)
156 if (element.level() == maxLevel_)
160 fvGeometry.bind(element);
161 elemVolVars.bind(element, fvGeometry, sol);
162 elemFluxVarsCache.bind(element, fvGeometry, elemVolVars);
168 for (
const auto& scv : scvs(fvGeometry))
170 auto source = problem_->source(element, fvGeometry, elemVolVars, scv);
171 auto pointSource = problem_->scvPointSources(element, fvGeometry, elemVolVars, scv);
172 if (source.infinity_norm() + pointSource.infinity_norm() > eps_)
174 indicatorVector_[eIdx] =
true;
181 if (!indicatorVector_[eIdx]
182 && element.hasBoundaryIntersections()
183 && (refineAtDirichletBC_ || refineAtFluxBC_))
188 for (
const auto& scvf : scvfs(fvGeometry))
191 if (!scvf.boundary())
194 const auto bcTypes = problem_->boundaryTypes(element, scvf);
196 if(bcTypes.hasOnlyDirichlet() && refineAtDirichletBC_)
198 indicatorVector_[eIdx] =
true;
203 else if(refineAtFluxBC_)
205 const auto fluxes = problem_->neumann(element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
206 if (fluxes.infinity_norm() > eps_)
208 indicatorVector_[eIdx] =
true;
219 std::vector<BoundaryTypes> bcTypes(fvGeometry.numScv());
222 for (
const auto& scv : scvs(fvGeometry))
224 bcTypes[scv.localDofIndex()] = problem_->boundaryTypes(element, scv);
225 if (refineAtDirichletBC_ && bcTypes[scv.localDofIndex()].hasDirichlet())
227 indicatorVector_[eIdx] =
true;
233 if (!indicatorVector_[eIdx] && refineAtFluxBC_)
235 for (
const auto& scvf : scvfs(fvGeometry))
238 if (scvf.boundary() && bcTypes[scvf.insideScvIdx()].hasNeumann())
240 const auto fluxes = problem_->neumann(element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
241 if (fluxes.infinity_norm() > eps_)
243 indicatorVector_[eIdx] =
true;