3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
c/sequential/celldataadaptive.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_ELEMENTDATA2P2C_ADAPTIVE_HH
25#define DUMUX_ELEMENTDATA2P2C_ADAPTIVE_HH
26
27#include <dune/common/float_cmp.hh>
28#include <dune/grid/utility/persistentcontainer.hh>
29
32
33namespace Dumux {
43template<class TypeTag>
45{
46private:
49 using Grid = typename GridView::Grid;
52
53 enum
54 {
55 dim = GridView::dimension
56 };
57
59
60 enum
61 {
62 wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx,
63 wCompIdx = Indices::wCompIdx, nCompIdx = Indices::nCompIdx
64 };
65 enum
66 {
67 numPhases = getPropValue<TypeTag, Properties::NumPhases>()
68 };
69 using Element = typename GridView::Traits::template Codim<0>::Entity;
70
72 static constexpr int pressureType = getPropValue<TypeTag, Properties::PressureFormulation>();
73 int upwindError_[numPhases];
74
75public:
88 {
89 Dune::FieldVector<Scalar,2> totalConcentration_;
90 Dune::FieldVector<Scalar,2> pressure_;
92 Dune::FieldVector<Scalar,3> volumeDerivatives_;
93 Scalar cellVolume;
96 int count;
99 {
101 pressure_ = 0.;
103 cellVolume = 0.;
104 subdomain=-1;
105 count = 0;
106 isRefined = false;
107 }
108 };
109
110
113 {
114 for (int i = 0; i < numPhases;i++)
115 upwindError_[i] = 0;
116 }
117
126 void storeAdaptionValues(AdaptedValues& adaptedValues, const Element& element)
127 {
128 adaptedValues.totalConcentration_[wCompIdx]= this->massConcentration(wCompIdx);
129 adaptedValues.totalConcentration_[nCompIdx]= this->massConcentration(nCompIdx);
130 adaptedValues.pressure_[wPhaseIdx] = this->pressure(wPhaseIdx);
131 adaptedValues.pressure_[nPhaseIdx] = this->pressure(nPhaseIdx);
132 adaptedValues.volumeDerivatives_[wCompIdx] = this->dv(wPhaseIdx);
133 adaptedValues.volumeDerivatives_[nCompIdx] = this->dv(nPhaseIdx);
134 adaptedValues.volumeDerivatives_[2] = this->dv_dp();
135 adaptedValues.cellVolume = element.geometry().volume();
136 adaptedValues.subdomain = this->subdomain();
137 adaptedValues.fluxData_=this->fluxData();
138 }
139
150 static void storeAdaptionValues(AdaptedValues& adaptedValues,
151 AdaptedValues& adaptedValuesFather,
152 const Element& fatherElement)
153 {
154 adaptedValuesFather.totalConcentration_[wCompIdx]
155 += adaptedValues.cellVolume* adaptedValues.totalConcentration_[wCompIdx];
156 adaptedValuesFather.totalConcentration_[nCompIdx]
157 += adaptedValues.cellVolume* adaptedValues.totalConcentration_[nCompIdx];
158 // if all cells are summed up, re-convert mass into total concentrations
159 Scalar fatherVolume = fatherElement.geometry().volume();
160 if(adaptedValuesFather.count == 1 << dim)
161 {
162 adaptedValuesFather.totalConcentration_[wCompIdx] /= fatherVolume;
163 adaptedValuesFather.totalConcentration_[nCompIdx] /= fatherVolume;
164 }
165 adaptedValuesFather.cellVolume = fatherVolume;
166
167 adaptedValuesFather.pressure_[wPhaseIdx] += adaptedValues.pressure_[wPhaseIdx] / adaptedValues.count;
168 adaptedValuesFather.pressure_[nPhaseIdx] += adaptedValues.pressure_[nPhaseIdx] / adaptedValues.count;
169 // apply maximum complexity for new cell
170 using std::max;
171 adaptedValuesFather.subdomain = max(adaptedValuesFather.subdomain, adaptedValues.subdomain);
172 }
173
185 void setAdaptionValues(AdaptedValues& adaptedValues, const Element& element)
186 {
187 // in new cells, there is a cellData object, but not yet a fluidstate.
188 // To write the adapted values in this fluidstate, we have to enshure that
189 // it was created. We therefore specify the type of FluidState that should be stored.
190 this->setSubdomainAndFluidStateType(adaptedValues.subdomain);
191 this->setMassConcentration(wCompIdx,
192 adaptedValues.totalConcentration_[wCompIdx]);
193 this->setMassConcentration(nCompIdx,
194 adaptedValues.totalConcentration_[nCompIdx]);
195 this->setPressure(wPhaseIdx,
196 adaptedValues.pressure_[wPhaseIdx] / adaptedValues.count);
197 this->setPressure(nPhaseIdx,
198 adaptedValues.pressure_[nPhaseIdx] / adaptedValues.count);
199
200 //copy flux directions
201 this->fluxData()=adaptedValues.fluxData_;
202
203 //indicate if cell was just refined in this TS
204 this->wasRefined()= adaptedValues.isRefined;
205 if(adaptedValues.isRefined)
206 {
207 this->dv(wPhaseIdx) = adaptedValues.volumeDerivatives_[wCompIdx];
208 this->dv(nPhaseIdx) = adaptedValues.volumeDerivatives_[nCompIdx];
209 this->dv_dp() = adaptedValues.volumeDerivatives_[2];
210 this->volumeDerivativesAvailable(true);
211 }
212 else
213 this->volumeDerivativesAvailable(false); // recalculate volume derivatives
214 }
215
227 static void reconstructAdaptionValues(Dune::PersistentContainer<Grid, AdaptedValues>& adaptationMap,
228 const Element& father, const Element& son, const Problem& problem)
229 {
230 // acess father and son
231 AdaptedValues& adaptedValues = adaptationMap[son];
232 AdaptedValues& adaptedValuesFather = adaptationMap[father];
233
234 adaptedValues.subdomain = adaptedValuesFather.subdomain;
235 adaptedValues.totalConcentration_[wCompIdx]
236 = adaptedValuesFather.totalConcentration_[wCompIdx] / adaptedValuesFather.count;
237 adaptedValues.totalConcentration_[nCompIdx]
238 = adaptedValuesFather.totalConcentration_[nCompIdx] / adaptedValuesFather.count;
239
240 // Introduce a hydrostatic pressure distribution inside the father cell
241 Scalar gTimesHeight = problem.gravity()
242 * (son.geometry().center() - father.geometry().center());
243 gTimesHeight *= (adaptedValuesFather.totalConcentration_[wCompIdx]+adaptedValuesFather.totalConcentration_[nCompIdx])/
244 problem.spatialParams().porosity(son);
245// int globalIdxSon = problem.variables().index(son);
246
247
248 // recalculate new Primary variable and store pc (i.e. other pressure)
249 if(pressureType == wPhaseIdx)
250 {
251 // recompute pressure and pc
252 Scalar pressure = adaptedValuesFather.pressure_[wPhaseIdx] / adaptedValuesFather.count;
253 Scalar pc = adaptedValuesFather.pressure_[nPhaseIdx] / adaptedValuesFather.count
254 - pressure;
255
256 adaptedValues.pressure_[wPhaseIdx]
257 = pressure + gTimesHeight ;
258 adaptedValues.pressure_[nPhaseIdx]
259 = pc + adaptedValues.pressure_[wPhaseIdx];
260 }
261 else
262 {
263 // recompute pressure and pc
264 Scalar pressure = adaptedValuesFather.pressure_[nPhaseIdx] / adaptedValuesFather.count;
265 Scalar pc = pressure
266 - adaptedValuesFather.pressure_[wPhaseIdx] / adaptedValuesFather.count;
267
268 adaptedValues.pressure_[nPhaseIdx]
269 = pressure + gTimesHeight ;
270 adaptedValues.pressure_[wPhaseIdx]
271 = adaptedValues.pressure_[nPhaseIdx] - pc;
272 }
273 // copy volume derivatives and compressibility, and indicate that cell was refined
274 adaptedValues.volumeDerivatives_[wCompIdx] = adaptedValuesFather.volumeDerivatives_[wCompIdx];
275 adaptedValues.volumeDerivatives_[nCompIdx] = adaptedValuesFather.volumeDerivatives_[nCompIdx];
276 adaptedValues.volumeDerivatives_[2] = adaptedValuesFather.volumeDerivatives_[2];
277 adaptedValues.isRefined = true;
278 }
279};
280
281} // end namespace Dumux
282#endif
Storage container for discretized data for multi-physics models.
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
void volumeDerivativesAvailable(bool value)
Specifies if volume derivatives are computed and available.
Definition: 2p2c/sequential/celldata.hh:365
bool & wasRefined()
Indicates if current cell was refined at this time step.
Definition: 2p2c/sequential/celldata.hh:374
const Scalar massConcentration(int compIdx) const
Returns the total mass concentration of a component .
Definition: 2p2c/sequential/celldata.hh:145
void setMassConcentration(int compIdx, Scalar value)
Sets the total mass concentration of a component .
Definition: 2p2c/sequential/celldata.hh:162
Scalar & dv_dp()
Return the derivative of specific volume w.r.t. pressure For details, see description of FVPressureCo...
Definition: 2p2c/sequential/celldata.hh:242
FluxData & fluxData()
Acess to flux data, representing information living on the intersections.
Definition: 2p2c/sequential/celldata.hh:101
Scalar & dv(int compIdx)
Return the derivative of spec. volume w.r.t. change of mass For details, see description of FVPressur...
Definition: 2p2c/sequential/celldata.hh:258
Class including the data of a grid cell needed if an adaptive grid is used.
Definition: c/sequential/celldataadaptive.hh:45
CellData2P2CAdaptive()
Constructs an adaptive CellData object.
Definition: c/sequential/celldataadaptive.hh:112
void setAdaptionValues(AdaptedValues &adaptedValues, const Element &element)
Set adapted values in CellData This methods stores reconstructed values into the cellData object,...
Definition: c/sequential/celldataadaptive.hh:185
static void storeAdaptionValues(AdaptedValues &adaptedValues, AdaptedValues &adaptedValuesFather, const Element &fatherElement)
Adds cell information to father element for possible averaging / coarsening Sum up the adaptedValues ...
Definition: c/sequential/celldataadaptive.hh:150
void storeAdaptionValues(AdaptedValues &adaptedValues, const Element &element)
Stores leaf cell primary variables to transfer to new indexing Stores values to be adapted from the c...
Definition: c/sequential/celldataadaptive.hh:126
static void reconstructAdaptionValues(Dune::PersistentContainer< Grid, AdaptedValues > &adaptationMap, const Element &father, const Element &son, const Problem &problem)
Reconstructs sons entries from data of father cell Reconstructs an new solution from a father cell in...
Definition: c/sequential/celldataadaptive.hh:227
A container for all necessary variables to map an old solution to a new grid If the primary variables...
Definition: c/sequential/celldataadaptive.hh:88
FluxData2P2C< TypeTag > fluxData_
Information on old flux direction.
Definition: c/sequential/celldataadaptive.hh:94
Dune::FieldVector< Scalar, 2 > totalConcentration_
Transport primary variables.
Definition: c/sequential/celldataadaptive.hh:89
Dune::FieldVector< Scalar, 2 > pressure_
Definition: c/sequential/celldataadaptive.hh:90
Dune::FieldVector< Scalar, 3 > volumeDerivatives_
Storage for volume derivatives, as transport estimate on old pressure field is not correct after refi...
Definition: c/sequential/celldataadaptive.hh:92
int subdomain
subdomain
Definition: c/sequential/celldataadaptive.hh:95
int count
counts the number of cells averaged
Definition: c/sequential/celldataadaptive.hh:96
Scalar cellVolume
Cell volume for transformation of volume-specific primary variables.
Definition: c/sequential/celldataadaptive.hh:93
AdaptedValues()
Definition: c/sequential/celldataadaptive.hh:98
int isRefined
Indicates if cell is refined.
Definition: c/sequential/celldataadaptive.hh:97
Storage container for discretized data for multiphysics models.
Definition: celldatamultiphysics.hh:45
Scalar pressure(int phaseIdx)
DOC ME!
Definition: celldatamultiphysics.hh:79
void setPressure(int phaseIdx, Scalar value)
DOC ME!
Definition: celldatamultiphysics.hh:99
int & subdomain()
Return subdomain information.
Definition: celldatamultiphysics.hh:118
void setSubdomainAndFluidStateType(int index)
Specify subdomain information and fluidStateType.
Definition: celldatamultiphysics.hh:136
Class including the variables and data of discretized data of the constitutive relations.
Definition: 2p2c/sequential/fluxdata.hh:43
Storage container for discretized data of the constitutive relations for one element.