3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
gridadaptionindicator.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_GRIDADAPTIONINDICATOR2P_HH
25#define DUMUX_GRIDADAPTIONINDICATOR2P_HH
26
30
31namespace Dumux {
32
39template<class TypeTag>
41{
42private:
46 using Element = typename GridView::Traits::template Codim<0>::Entity;
47
49 using ScalarSolutionType = typename SolutionTypes::ScalarSolution;
50 using ElementMapper = typename SolutionTypes::ElementMapper;
51
53
54 enum
55 {
56 sw = Indices::saturationW,
57 sn = Indices::saturationNw
58 };
59 enum
60 {
61 wPhaseIdx = Indices::wPhaseIdx,
62 nPhaseIdx = Indices::nPhaseIdx
63 };
64
65public:
72 {
73 // prepare an indicator for refinement
74 if(indicatorVector_.size() != problem_.variables().cellDataGlobal().size())
75 {
76 indicatorVector_.resize(problem_.variables().cellDataGlobal().size());
77 }
78 indicatorVector_ = -1e100;
79
80 Scalar globalMax = -1e100;
81 Scalar globalMin = 1e100;
82
83 // 1) calculate Indicator -> min, maxvalues
84 // loop over all leaf-elements
85 for (const auto& element : elements(problem_.gridView()))
86 {
87 // calculate minimum and maximum saturation
88 // index of the current leaf-elements
89 int globalIdxI = problem_.variables().index(element);
90
91 Scalar satI = 0.0;
92 switch (saturationType_)
93 {
94 case sw:
95 satI = problem_.variables().cellData(globalIdxI).saturation(wPhaseIdx);
96 break;
97 case sn:
98 satI = problem_.variables().cellData(globalIdxI).saturation(nPhaseIdx);
99 break;
100 }
101
102 using std::min;
103 globalMin = min(satI, globalMin);
104 using std::max;
105 globalMax = max(satI, globalMax);
106
107 // calculate refinement indicator in all cells
108 for (const auto& intersection : intersections(problem_.gridView(), element))
109 {
110 // Only consider internal intersections
111 if (intersection.neighbor())
112 {
113 // Access neighbor
114 auto outside = intersection.outside();
115 int globalIdxJ = problem_.variables().index(outside);
116
117 // Visit intersection only once
118 if (element.level() > outside.level() || (element.level() == outside.level() && globalIdxI < globalIdxJ))
119 {
120 Scalar satJ = 0.;
121 switch (saturationType_)
122 {
123 case sw:
124 satJ = problem_.variables().cellData(globalIdxJ).saturation(wPhaseIdx);
125 break;
126 case sn:
127 satJ = problem_.variables().cellData(globalIdxJ).saturation(nPhaseIdx);
128 break;
129 }
130
131 using std::abs;
132 Scalar localdelta = abs(satI - satJ);
133 indicatorVector_[globalIdxI][0] = max(indicatorVector_[globalIdxI][0], localdelta);
134 indicatorVector_[globalIdxJ][0] = max(indicatorVector_[globalIdxJ][0], localdelta);
135 }
136 }
137 }
138 }
139
140 Scalar globaldelta = globalMax - globalMin;
141
142 refineBound_ = refinetol_*globaldelta;
143 coarsenBound_ = coarsentol_*globaldelta;
144
145#if HAVE_MPI
146 // communicate updated values
147 using DataHandle = VectorCommDataHandleEqual<ElementMapper, ScalarSolutionType, 0/*elementCodim*/>;
148 DataHandle dataHandle(problem_.elementMapper(), indicatorVector_);
149 problem_.gridView().template communicate<DataHandle>(dataHandle,
150 Dune::InteriorBorder_All_Interface,
151 Dune::ForwardCommunication);
152
153 using std::max;
154 refineBound_ = problem_.gridView().comm().max(refineBound_);
155 coarsenBound_ = problem_.gridView().comm().max(coarsenBound_);
156
157#endif
158 }
159
167 bool refine(const Element& element)
168 {
169 return (indicatorVector_[problem_.elementMapper().index(element)] > refineBound_);
170 }
171
179 bool coarsen(const Element& element)
180 {
181 return (indicatorVector_[problem_.elementMapper().index(element)] < coarsenBound_);
182 }
183
186 void init()
187 {
188 refineBound_ = 0.;
189 coarsenBound_ = 0.;
190 };
191
202 GridAdaptionIndicator2P (Problem& problem):
203 problem_(problem)
204 {
205 refinetol_ = getParam<Scalar>("GridAdapt.RefineTolerance");
206 coarsentol_ = getParam<Scalar>("GridAdapt.CoarsenTolerance");
207 }
208
209protected:
210 Problem& problem_;
215 ScalarSolutionType indicatorVector_;
216 static const int saturationType_ = getPropValue<TypeTag, Properties::SaturationFormulation>();
217};
218} // end namespace Dumux
219
220#endif
Contains a class to exchange entries of a vector.
Base file for properties related to sequential IMPET algorithms.
Definition: adapt.hh:29
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type GetProp
get the type of a property (equivalent to old macro GET_PROP(...))
Definition: propertysystem.hh:140
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
A data handle class to exchange entries of a vector.
Definition: vectorcommdatahandle.hh:78
Class defining a standard, saturation dependent indicator for grid adaption.
Definition: gridadaptionindicator.hh:41
Problem & problem_
Definition: gridadaptionindicator.hh:210
GridAdaptionIndicator2P(Problem &problem)
Constructs a GridAdaptionIndicator instance.
Definition: gridadaptionindicator.hh:202
Scalar refinetol_
Definition: gridadaptionindicator.hh:211
bool coarsen(const Element &element)
Indicator function for marking of grid cells for coarsening.
Definition: gridadaptionindicator.hh:179
Scalar coarsentol_
Definition: gridadaptionindicator.hh:212
static const int saturationType_
Definition: gridadaptionindicator.hh:216
Scalar refineBound_
Definition: gridadaptionindicator.hh:213
void calculateIndicator()
Calculates the indicator used for refinement/coarsening for each grid cell.
Definition: gridadaptionindicator.hh:71
ScalarSolutionType indicatorVector_
Definition: gridadaptionindicator.hh:215
void init()
Initializes the adaption indicator class.
Definition: gridadaptionindicator.hh:186
Scalar coarsenBound_
Definition: gridadaptionindicator.hh:214
bool refine(const Element &element)
Indicator function for marking of grid cells for refinement.
Definition: gridadaptionindicator.hh:167
Defines the properties required for (immiscible) two-phase sequential models.