3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
ointeractionvolume.hh
Go to the documentation of this file.
1/*****************************************************************************
2 * See the file COPYING for full copying permissions. *
3 * *
4 * This program is free software: you can redistribute it and/or modify *
5 * it under the terms of the GNU General Public License as published by *
6 * the Free Software Foundation, either version 3 of the License, or *
7 * (at your option) any later version. *
8 * *
9 * This program is distributed in the hope that it will be useful, *
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
12 * GNU General Public License for more details. *
13 * *
14 * You should have received a copy of the GNU General Public License *
15 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
16 *****************************************************************************/
17#ifndef DUMUX_FVMPFAOINTERACTIONVOLUME_HH
18#define DUMUX_FVMPFAOINTERACTIONVOLUME_HH
19
26
27namespace Dumux
28{
29
31
36template<class TypeTag>
38{
39private:
43
44 enum
45 {
46 dim = GridView::dimension
47 };
48
49 using Element = typename GridView::template Codim<0>::Entity;
50 using ElementSeed = typename Grid::template Codim<0>::EntitySeed;
51
54 using PrimaryVariables = typename SolutionTypes::PrimaryVariables;
55
56 using DimVector = Dune::FieldVector<Scalar, dim>;
57 using DimMatrix = Dune::FieldMatrix<Scalar, dim, dim>;
58 using FieldVectorVector = Dune::FieldVector<DimVector, dim>;
59 using IndexVector = Dune::FieldVector<int, dim>;
60 using BCTypeVector = std::vector<BoundaryTypes>;
61 using BCVector = std::vector<PrimaryVariables>;
62
63public:
65 {
66 inside = 1,
68 outside = -1
69 };
70
72 FVMPFAOInteractionVolume(const Grid& grid)
73 : grid_(&grid),
74 stored_(false), permTimesNu_(FieldVectorVector(DimVector(0.0))),
75 nu_(FieldVectorVector(DimVector(0.0))), normal_(FieldVectorVector(DimVector(0.0))),
76 faceArea_(DimVector(0.0)), dF_(0.0),
77 faceType_(2*dim, inside),
78 indexOnElement_(IndexVector(0.0)),
79 elements_(2*dim)
80 {
81 faceIndexOnSubVolume_[0][0] = 0;
82 faceIndexOnSubVolume_[0][1] = 3;
83 faceIndexOnSubVolume_[1][0] = 1;
84 faceIndexOnSubVolume_[1][1] = 0;
85 faceIndexOnSubVolume_[2][0] = 2;
86 faceIndexOnSubVolume_[2][1] = 1;
87 faceIndexOnSubVolume_[3][0] = 3;
88 faceIndexOnSubVolume_[3][1] = 2;
89 }
90
92 void setStored()
93 {
94 stored_ = true;
95 }
96
98 bool isStored() const
99 {
100 return stored_;
101 }
102
104
108 void setSubVolumeElement(const Element& element, int subVolumeIdx)
109 {
110 elements_[subVolumeIdx].push_back(element.seed());
111 }
112
114
118 void setDF(Scalar dF, int subVolumeIdx)
119 {
120 dF_[subVolumeIdx] = dF;
121 }
122
124
129 void setFaceArea(Scalar& faceArea, int subVolumeIdx, int subVolumeFaceIdxInInside)
130 {
131 faceArea_[subVolumeIdx][subVolumeFaceIdxInInside] = faceArea;
132 }
133
135
141 void setPermTimesNu(DimVector& nu, DimMatrix& perm, int subVolumeIdx, int subVolumeFaceIdxInInside)
142 {
143 permTimesNu_[subVolumeIdx][subVolumeFaceIdxInInside] = 0;
144
145 perm.mv(nu, permTimesNu_[subVolumeIdx][subVolumeFaceIdxInInside]);
146
147 nu_[subVolumeIdx][subVolumeFaceIdxInInside] = nu;
148 }
149
151
156 void setNormal(DimVector& normal, int subVolumeIdx, int subVolumeFaceIdxInInside)
157 {
158 normal_[subVolumeIdx][subVolumeFaceIdxInInside] = normal;
159 }
160
162
166 void setBoundary(BoundaryTypes& boundaryTypes, int subVolumeFaceIdx)
167 {
168// std::cout<<"old BCType = "<<boundaryType_[subVolumeFaceIdx]<<"\n";
169// std::cout<<"new BCType = "<<boundaryType<<"\n";
170 if (boundaryTypes_.size() == 0)
171 {
172 boundaryTypes_.resize(2 * dim);
173 }
174 boundaryTypes_[subVolumeFaceIdx] = boundaryTypes;
175 faceType_[subVolumeFaceIdx] = boundary;
176 }
177
179
182 void setOutsideFace(int subVolumeFaceIdx)
183 {
184 faceType_[subVolumeFaceIdx] = outside;
185 }
186
188
192 void setDirichletCondition(PrimaryVariables& condition, int subVolumeFaceIdx)
193 {
194 if (dirichletValues_.size() == 0)
195 {
196 dirichletValues_.resize(2 * dim);
197 }
198 dirichletValues_[subVolumeFaceIdx] = condition;
199 }
200
202
206 void setNeumannCondition(PrimaryVariables& condition, int subVolumeFaceIdx)
207 {
208 if (neumannValues_.size() == 0)
209 {
210 neumannValues_.resize(2 * dim);
211 }
212 neumannValues_[subVolumeFaceIdx] = condition;
213 }
214
216
221 void setIndexOnElement(int indexInInside, int subVolumeIdx, int subVolumeFaceIdxInInside)
222 {
223 indexOnElement_[subVolumeIdx][subVolumeFaceIdxInInside] = indexInInside;
224 }
225
227
233 int getIndexOnElement(int subVolumeIdx, int subVolumeFaceIdx)
234 {
235 return indexOnElement_[subVolumeIdx][subVolumeFaceIdx];
236 }
237
239
245 int getFaceIndexFromSubVolume(int subVolumeIdx, int subVolumeFaceIdx)
246 {
247 return faceIndexOnSubVolume_[subVolumeIdx][subVolumeFaceIdx];
248 }
249
251
256 Element getSubVolumeElement(int subVolumeIdx)
257 {
258 return grid_->entity(elements_[subVolumeIdx][0]);
259 }
260
262
267 BoundaryTypes& getBoundaryType(int subVolumeFaceIdx)
268 {
269 return boundaryTypes_[subVolumeFaceIdx];
270 }
271
273
276 bool isOutsideFace(int subVolumeFaceIdx)
277 {
278 return (faceType_[subVolumeFaceIdx] == outside);
279 }
280
282
285 bool isInsideFace(int subVolumeFaceIdx)
286 {
287 return (faceType_[subVolumeFaceIdx] == inside);
288 }
289
292 {
293 for (unsigned int i = 0; i < faceType_.size(); i++)
294 {
295 if (isOutsideFace(i) || isBoundaryFace(i))
296 return false;
297 }
298 return true;
299 }
300
302
305 bool isBoundaryFace(int subVolumeFaceIdx)
306 {
307 return (faceType_[subVolumeFaceIdx] == boundary);
308 }
309
311
316 PrimaryVariables& getDirichletValues(int subVolumeFaceIdx)
317 {
318 return dirichletValues_[subVolumeFaceIdx];
319 }
320
322
327 PrimaryVariables& getNeumannValues(int subVolumeFaceIdx)
328 {
329 return neumannValues_[subVolumeFaceIdx];
330 }
331
333
339 DimVector& getNormal(int subVolumeIdx, int subVolumeFaceIdxInInside)
340 {
341 return normal_[subVolumeIdx][subVolumeFaceIdxInInside];
342 }
343
345
351 Scalar& getFaceArea(int subVolumeIdx, int subVolumeFaceIdxInInside)
352 {
353 return faceArea_[subVolumeIdx][subVolumeFaceIdxInInside];
354 }
355
357
363 DimVector& getNu(int subVolumeIdx, int subVolumeFaceIdxInInside)
364 {
365 return nu_[subVolumeIdx][subVolumeFaceIdxInInside];
366 }
367
369
378 Scalar getNtkNu(int subVolumeIdx, int subVolumeFaceIdxInInsideN, int subVolumeFaceIdxInInsideNu) const
379 {
380 return normal_[subVolumeIdx][subVolumeFaceIdxInInsideN] * permTimesNu_[subVolumeIdx][subVolumeFaceIdxInInsideNu];
381 }
382
384
394 Scalar getNtkrkNu(Scalar& relPerm, int subVolumeIdx, int subVolumeFaceIdxInInsideN,
395 int subVolumeFaceIdxInInsideNu) const
396 {
397 DimVector krKNu(permTimesNu_[subVolumeIdx][subVolumeFaceIdxInInsideNu]);
398
399 krKNu *= relPerm;
400
401 return normal_[subVolumeIdx][subVolumeFaceIdxInInsideN] * krKNu;
402 }
403
406
415 Scalar getNtkNu_df(int subVolumeIdx, int subVolumeFaceIdxInInsideN, int subVolumeFaceIdxInInsideNu) const
416 {
417 return faceArea_[subVolumeIdx][subVolumeFaceIdxInInsideN]
418 * getNtkNu(subVolumeIdx, subVolumeFaceIdxInInsideN, subVolumeFaceIdxInInsideNu)
419 / dF_[subVolumeIdx];
420 }
421
424
434 Scalar getNtkrkNu_df(Scalar& relPerm, int subVolumeIdx, int subVolumeFaceIdxInInsideN, int subVolumeFaceIdxInInsideNu) const
435 {
436 return faceArea_[subVolumeIdx][subVolumeFaceIdxInInsideN]
437 * getNtkrkNu(relPerm, subVolumeIdx, subVolumeFaceIdxInInsideN, subVolumeFaceIdxInInsideNu)
438 / dF_[subVolumeIdx];
439 }
440
441private:
442 const Grid* grid_;
443 bool stored_;
444 Dune::FieldVector<FieldVectorVector, 2*dim> permTimesNu_;
445 Dune::FieldVector<FieldVectorVector, 2*dim> nu_;
446 Dune::FieldVector<FieldVectorVector, 2*dim> normal_;
447 Dune::FieldVector<DimVector, 2*dim> faceArea_;
448 Dune::FieldVector<Scalar, 2*dim> dF_;
449 BCTypeVector boundaryTypes_;
450 std::vector<int> faceType_;
451 Dune::FieldVector<IndexVector, 2*dim> indexOnElement_;
452 Dune::FieldVector<IndexVector, 2*dim> faceIndexOnSubVolume_;
453 std::vector<std::vector<ElementSeed> > elements_;
454 BCVector neumannValues_;
455 BCVector dirichletValues_;
456};
457}
458#endif
Vector normal(const Vector &v)
Create a vector normal to the given one (v is expected to be non-zero)
Definition: geometry/normal.hh:36
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
Class including the information of an interaction volume of a MPFA O-method that does not change with...
Definition: ointeractionvolume.hh:38
DimVector & getNu(int subVolumeIdx, int subVolumeFaceIdxInInside)
Get vector used for the transmissiblity calculation.
Definition: ointeractionvolume.hh:363
bool isInnerVolume()
Returns true if the interaction volume is completely inside the model domain.
Definition: ointeractionvolume.hh:291
FVMPFAOInteractionVolume(const Grid &grid)
Constructs a FVMPFAOInteractionVolume object.
Definition: ointeractionvolume.hh:72
bool isInsideFace(int subVolumeFaceIdx)
Returns true if an interaction volume flux face is inside the model domain and no boundary face.
Definition: ointeractionvolume.hh:285
void setIndexOnElement(int indexInInside, int subVolumeIdx, int subVolumeFaceIdxInInside)
Store map from local interaction volume numbering to numbering of the Dune reference element.
Definition: ointeractionvolume.hh:221
bool isOutsideFace(int subVolumeFaceIdx)
Returns true if an interaction volume flux face is outside the model domain.
Definition: ointeractionvolume.hh:276
Scalar getNtkrkNu(Scalar &relPerm, int subVolumeIdx, int subVolumeFaceIdxInInsideN, int subVolumeFaceIdxInInsideNu) const
Get for the transmissiblity calculation.
Definition: ointeractionvolume.hh:394
BoundaryTypes & getBoundaryType(int subVolumeFaceIdx)
Get boundary condtion types for a flux face.
Definition: ointeractionvolume.hh:267
void setDirichletCondition(PrimaryVariables &condition, int subVolumeFaceIdx)
Store Dirichlet boundary condtions for a flux face.
Definition: ointeractionvolume.hh:192
void setSubVolumeElement(const Element &element, int subVolumeIdx)
Store an element of the interaction volume.
Definition: ointeractionvolume.hh:108
bool isBoundaryFace(int subVolumeFaceIdx)
Returns true if an interaction volume flux face is a boundary face.
Definition: ointeractionvolume.hh:305
PrimaryVariables & getDirichletValues(int subVolumeFaceIdx)
Get the Dirichlet boundary condtions for a flux face.
Definition: ointeractionvolume.hh:316
void setOutsideFace(int subVolumeFaceIdx)
Mark a flux face to be outside the model domain.
Definition: ointeractionvolume.hh:182
void setBoundary(BoundaryTypes &boundaryTypes, int subVolumeFaceIdx)
Store boundary condtion types for a flux face.
Definition: ointeractionvolume.hh:166
void setNeumannCondition(PrimaryVariables &condition, int subVolumeFaceIdx)
Store Neumann boundary condtions for a flux face.
Definition: ointeractionvolume.hh:206
Scalar & getFaceArea(int subVolumeIdx, int subVolumeFaceIdxInInside)
Get a flux face area.
Definition: ointeractionvolume.hh:351
Scalar getNtkrkNu_df(Scalar &relPerm, int subVolumeIdx, int subVolumeFaceIdxInInsideN, int subVolumeFaceIdxInInsideNu) const
Definition: ointeractionvolume.hh:434
void setStored()
Mark storage as completed.
Definition: ointeractionvolume.hh:92
Scalar getNtkNu_df(int subVolumeIdx, int subVolumeFaceIdxInInsideN, int subVolumeFaceIdxInInsideNu) const
Definition: ointeractionvolume.hh:415
PrimaryVariables & getNeumannValues(int subVolumeFaceIdx)
Get the Neumann boundary condtions for a flux face.
Definition: ointeractionvolume.hh:327
Scalar getNtkNu(int subVolumeIdx, int subVolumeFaceIdxInInsideN, int subVolumeFaceIdxInInsideNu) const
Get for the transmissiblity calculation.
Definition: ointeractionvolume.hh:378
DimVector & getNormal(int subVolumeIdx, int subVolumeFaceIdxInInside)
Get a flux face normal.
Definition: ointeractionvolume.hh:339
void setNormal(DimVector &normal, int subVolumeIdx, int subVolumeFaceIdxInInside)
Store a flux face normal.
Definition: ointeractionvolume.hh:156
void setFaceArea(Scalar &faceArea, int subVolumeIdx, int subVolumeFaceIdxInInside)
Store a flux face area.
Definition: ointeractionvolume.hh:129
int getFaceIndexFromSubVolume(int subVolumeIdx, int subVolumeFaceIdx)
Map from local interaction volume numbering on element to numbering on interaction volume.
Definition: ointeractionvolume.hh:245
bool isStored() const
Returns true if information has already been stored.
Definition: ointeractionvolume.hh:98
void setDF(Scalar dF, int subVolumeIdx)
Store the for the transmissiblity calculation.
Definition: ointeractionvolume.hh:118
FaceTypes
Definition: ointeractionvolume.hh:65
@ outside
Definition: ointeractionvolume.hh:68
@ boundary
Definition: ointeractionvolume.hh:67
@ inside
Definition: ointeractionvolume.hh:66
void setPermTimesNu(DimVector &nu, DimMatrix &perm, int subVolumeIdx, int subVolumeFaceIdxInInside)
Store for the transmissiblity calculation.
Definition: ointeractionvolume.hh:141
Element getSubVolumeElement(int subVolumeIdx)
Get an element of the interaction volume.
Definition: ointeractionvolume.hh:256
int getIndexOnElement(int subVolumeIdx, int subVolumeFaceIdx)
Map from local interaction volume numbering to numbering of the Dune reference element.
Definition: ointeractionvolume.hh:233
Properties for a MPFA method.