3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
linteractionvolume.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 *****************************************************************************/
19
20#ifndef DUMUX_FVMPFALINTERACTIONVOLUME_HH
21#define DUMUX_FVMPFALINTERACTIONVOLUME_HH
22
29
30namespace Dumux
31{
32
34
39template<class TypeTag>
41{
42private:
46
47 enum
48 {
49 dim = GridView::dimension
50 };
51
52 using Element = typename GridView::template Codim<0>::Entity;
53 using ElementSeed = typename Grid::template Codim<0>::EntitySeed;
54
57 using PrimaryVariables = typename SolutionTypes::PrimaryVariables;
58
59 using DimVector = Dune::FieldVector<Scalar, dim>;
60 using FieldVectorVector = Dune::FieldVector<DimVector, dim>;
61 using IndexVector = Dune::FieldVector<int, dim>;
62 using BCTypeVector = std::vector<BoundaryTypes>;
63 using BCVector = std::vector<PrimaryVariables>;
64
65public:
67 {
68 inside = 1,
70 outside = -1
71 };
72
74 FVMPFALInteractionVolume(const Grid& grid)
75 : grid_(&grid),
76 stored_(false), normal_(FieldVectorVector(DimVector(0.0))),
77 facePos_(FieldVectorVector(DimVector(0.0))),
78 faceArea_(DimVector(0.0)),
79 faceType_(2*dim, inside),
80 indexOnElement_(IndexVector(0.0)),
81 elements_(2*dim),
82 centerVertexPos_(0),
83 elementNum_(0)
84 {
85 faceIndexOnSubVolume_[0][0] = 0;
86 faceIndexOnSubVolume_[0][1] = 3;
87 faceIndexOnSubVolume_[1][0] = 1;
88 faceIndexOnSubVolume_[1][1] = 0;
89 faceIndexOnSubVolume_[2][0] = 2;
90 faceIndexOnSubVolume_[2][1] = 1;
91 faceIndexOnSubVolume_[3][0] = 3;
92 faceIndexOnSubVolume_[3][1] = 2;
93 }
94
96 void reset()
97 {
98 stored_ = false;
99 elements_.clear();
100 elements_.resize(2*dim);
101 centerVertexPos_ = 0;
102 indexOnElement_ = IndexVector(0.0);
103 faceType_.clear();
104 faceType_.resize(2*dim, inside);
105 faceArea_ = DimVector(0.0);
106 facePos_ = FieldVectorVector(DimVector(0.0));
107 normal_ = FieldVectorVector(DimVector(0.0));
108 elementNum_ = 0;
109 boundaryTypes_.clear();
110 neumannValues_.clear();
111 dirichletValues_.clear();
112 }
113
116 {
117 stored_ = true;
118 }
119
121 bool isStored() const
122 {
123 return stored_;
124 }
125
127 void setCenterPosition(DimVector &centerVertexPos)
128 {
129 centerVertexPos_ = centerVertexPos;
130 }
131
133
137 void setSubVolumeElement(const Element& element, int subVolumeIdx)
138 {
139 elements_[subVolumeIdx].push_back(element.seed());
140 elementNum_++;
141 }
142
144
149 void setFacePosition(const DimVector& pos, int subVolumeIdx, int subVolumeFaceIdxInInside)
150 {
151 facePos_[subVolumeIdx][subVolumeFaceIdxInInside] = pos;
152 }
153
155
160 void setFaceArea(Scalar& faceArea, int subVolumeIdx, int subVolumeFaceIdxInInside)
161 {
162 faceArea_[subVolumeIdx][subVolumeFaceIdxInInside] = faceArea;
163 }
164
166
171 void setNormal(DimVector& normal, int subVolumeIdx, int subVolumeFaceIdxInInside)
172 {
173 normal_[subVolumeIdx][subVolumeFaceIdxInInside] = normal;
174 }
175
177
181 void setBoundary(BoundaryTypes& boundaryTypes, int subVolumeFaceIdx)
182 {
183 if (boundaryTypes_.size() == 0)
184 {
185 boundaryTypes_.resize(2 * dim);
186 }
187 boundaryTypes_[subVolumeFaceIdx] = boundaryTypes;
188 faceType_[subVolumeFaceIdx] = boundary;
189 }
190
192
195 void setOutsideFace(int subVolumeFaceIdx)
196 {
197 faceType_[subVolumeFaceIdx] = outside;
198 }
199
201
205 void setDirichletCondition(PrimaryVariables& condition, int subVolumeFaceIdx)
206 {
207 if (dirichletValues_.size() == 0)
208 {
209 dirichletValues_.resize(2 * dim);
210 }
211 dirichletValues_[subVolumeFaceIdx] = condition;
212 }
213
215
219 void setNeumannCondition(PrimaryVariables& condition, int subVolumeFaceIdx)
220 {
221 if (neumannValues_.size() == 0)
222 {
223 neumannValues_.resize(2 * dim);
224 }
225 neumannValues_[subVolumeFaceIdx] = condition;
226 }
227
229
234 void setIndexOnElement(int indexInInside, int subVolumeIdx, int subVolumeFaceIdxInInside)
235 {
236 indexOnElement_[subVolumeIdx][subVolumeFaceIdxInInside] = indexInInside;
237 }
238
240 DimVector& getCenterPosition()
241 {
242 return centerVertexPos_;
243 }
244
247 {
248 return elementNum_;
249 }
250
252
258 int getIndexOnElement(int subVolumeIdx, int subVolumeFaceIdx)
259 {
260 return indexOnElement_[subVolumeIdx][subVolumeFaceIdx];
261 }
262
264
270 int getFaceIndexFromSubVolume(int subVolumeIdx, int subVolumeFaceIdx)
271 {
272 return faceIndexOnSubVolume_[subVolumeIdx][subVolumeFaceIdx];
273 }
274
276
281 Element getSubVolumeElement(int subVolumeIdx)
282 {
283 return grid_->entity(elements_[subVolumeIdx][0]);
284 }
285
287
292 BoundaryTypes& getBoundaryType(int subVolumeFaceIdx)
293 {
294 return boundaryTypes_[subVolumeFaceIdx];
295 }
296
298
301 bool isOutsideFace(int subVolumeFaceIdx)
302 {
303 return (faceType_[subVolumeFaceIdx] == outside);
304 }
305
307
310 bool isInsideFace(int subVolumeFaceIdx)
311 {
312 return (faceType_[subVolumeFaceIdx] == inside);
313 }
314
317 {
318 for (unsigned int i = 0; i < faceType_.size(); i++)
319 {
320 if (isOutsideFace(i) || isBoundaryFace(i))
321 return false;
322 }
323 return true;
324 }
325
327
330 bool isBoundaryFace(int subVolumeFaceIdx)
331 {
332 return (faceType_[subVolumeFaceIdx] == boundary);
333 }
334
336
341 PrimaryVariables& getDirichletValues(int subVolumeFaceIdx)
342 {
343 return dirichletValues_[subVolumeFaceIdx];
344 }
345
347
352 PrimaryVariables& getNeumannValues(int subVolumeFaceIdx)
353 {
354 return neumannValues_[subVolumeFaceIdx];
355 }
356
358
364 DimVector& getNormal(int subVolumeIdx, int subVolumeFaceIdxInInside)
365 {
366 return normal_[subVolumeIdx][subVolumeFaceIdxInInside];
367 }
368
370
376 DimVector& getFacePosition(int subVolumeIdx, int subVolumeFaceIdxInInside)
377 {
378 return facePos_[subVolumeIdx][subVolumeFaceIdxInInside];
379 }
380
382
388 Scalar& getFaceArea(int subVolumeIdx, int subVolumeFaceIdxInInside)
389 {
390 return faceArea_[subVolumeIdx][subVolumeFaceIdxInInside];
391 }
392
395 {
396 std::cout<<"\nNumber of stored elements: "<<elementNum_<<"\n";
397 std::cout<<" center position: "<<centerVertexPos_<<"\n";
398 for (int i = 0; i < 2*dim; i++)
399 {
400 if (elements_[i].size() > 0)
401 {
402 std::cout<<"element "<<i<<":\n";
403 std::cout<<"element position: "<<elements_[i][0]->geometry().center()<<"\n";
404 std::cout<<"face indices on element: "<<indexOnElement_[i]<<"\n";
405 std::cout<<"face normals on element: "<<normal_[i]<<"\n";
406 std::cout<<"face areas on element: "<<faceArea_[i]<<"\n";
407 std::cout<<"face position on element: "<<facePos_[i]<<"\n";
408 std::cout<<"face type: "<<faceType_[i]<<"\n";
409 }
410 }
411 }
412
413private:
414 const Grid* grid_;
415 bool stored_;
416 Dune::FieldVector<FieldVectorVector, 2*dim> normal_;
417 Dune::FieldVector<FieldVectorVector, 2*dim> facePos_;
418 Dune::FieldVector<DimVector, 2*dim> faceArea_;
419 BCTypeVector boundaryTypes_;
420 std::vector<int> faceType_;
421 Dune::FieldVector<IndexVector, 2*dim> indexOnElement_;
422 Dune::FieldVector<IndexVector, 2*dim> faceIndexOnSubVolume_;
423 std::vector<std::vector<ElementSeed> > elements_;
424 BCVector neumannValues_;
425 BCVector dirichletValues_;
426 DimVector centerVertexPos_;
427 int elementNum_;
428};
429}
430#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 L-method that does not change with...
Definition: linteractionvolume.hh:41
bool isInnerVolume()
Returns true if the interaction volume is completely inside the model domain.
Definition: linteractionvolume.hh:316
void setSubVolumeElement(const Element &element, int subVolumeIdx)
Store an element of the interaction volume.
Definition: linteractionvolume.hh:137
int getElementNumber()
Get number of stored elements.
Definition: linteractionvolume.hh:246
void reset()
Delete stored information.
Definition: linteractionvolume.hh:96
void setFacePosition(const DimVector &pos, int subVolumeIdx, int subVolumeFaceIdxInInside)
Store position of the center of a flux face.
Definition: linteractionvolume.hh:149
void setBoundary(BoundaryTypes &boundaryTypes, int subVolumeFaceIdx)
Store boundary condtion types for a flux face.
Definition: linteractionvolume.hh:181
int getIndexOnElement(int subVolumeIdx, int subVolumeFaceIdx)
Map from local interaction volume numbering to numbering of the Dune reference element.
Definition: linteractionvolume.hh:258
void setIndexOnElement(int indexInInside, int subVolumeIdx, int subVolumeFaceIdxInInside)
Store map from local interaction volume numbering to numbering of the Dune reference element.
Definition: linteractionvolume.hh:234
bool isInsideFace(int subVolumeFaceIdx)
Returns true if an interaction volume flux face is inside the model domain and no boundary face.
Definition: linteractionvolume.hh:310
void setOutsideFace(int subVolumeFaceIdx)
Mark a flux face to be outside the model domain.
Definition: linteractionvolume.hh:195
DimVector & getNormal(int subVolumeIdx, int subVolumeFaceIdxInInside)
Get a flux face normal.
Definition: linteractionvolume.hh:364
Element getSubVolumeElement(int subVolumeIdx)
Get an element of the interaction volume.
Definition: linteractionvolume.hh:281
DimVector & getFacePosition(int subVolumeIdx, int subVolumeFaceIdxInInside)
Get position of the center of a flux face.
Definition: linteractionvolume.hh:376
void setStored()
Mark storage as completed.
Definition: linteractionvolume.hh:115
BoundaryTypes & getBoundaryType(int subVolumeFaceIdx)
Get boundary condtion types for a flux face.
Definition: linteractionvolume.hh:292
void printInteractionVolumeInfo()
Outputs stored information.
Definition: linteractionvolume.hh:394
Scalar & getFaceArea(int subVolumeIdx, int subVolumeFaceIdxInInside)
Get a flux face area.
Definition: linteractionvolume.hh:388
void setNormal(DimVector &normal, int subVolumeIdx, int subVolumeFaceIdxInInside)
Store a flux face normal.
Definition: linteractionvolume.hh:171
PrimaryVariables & getNeumannValues(int subVolumeFaceIdx)
Get the Neumann boundary condtions for a flux face.
Definition: linteractionvolume.hh:352
bool isBoundaryFace(int subVolumeFaceIdx)
Returns true if an interaction volume flux face is a boundary face.
Definition: linteractionvolume.hh:330
void setDirichletCondition(PrimaryVariables &condition, int subVolumeFaceIdx)
Store Dirichlet boundary condtions for a flux face.
Definition: linteractionvolume.hh:205
FVMPFALInteractionVolume(const Grid &grid)
Constructs a FVMPFALInteractionVolume object.
Definition: linteractionvolume.hh:74
FaceTypes
Definition: linteractionvolume.hh:67
@ inside
Definition: linteractionvolume.hh:68
@ outside
Definition: linteractionvolume.hh:70
@ boundary
Definition: linteractionvolume.hh:69
int getFaceIndexFromSubVolume(int subVolumeIdx, int subVolumeFaceIdx)
Map from local interaction volume numbering on element to numbering on interaction volume.
Definition: linteractionvolume.hh:270
void setCenterPosition(DimVector &centerVertexPos)
Store position of the center vertex of the interaction volume.
Definition: linteractionvolume.hh:127
void setNeumannCondition(PrimaryVariables &condition, int subVolumeFaceIdx)
Store Neumann boundary condtions for a flux face.
Definition: linteractionvolume.hh:219
bool isStored() const
Returns true if information has already been stored.
Definition: linteractionvolume.hh:121
DimVector & getCenterPosition()
Get position vector of central vertex.
Definition: linteractionvolume.hh:240
void setFaceArea(Scalar &faceArea, int subVolumeIdx, int subVolumeFaceIdxInInside)
Store a flux face area.
Definition: linteractionvolume.hh:160
PrimaryVariables & getDirichletValues(int subVolumeFaceIdx)
Get the Dirichlet boundary condtions for a flux face.
Definition: linteractionvolume.hh:341
bool isOutsideFace(int subVolumeFaceIdx)
Returns true if an interaction volume flux face is outside the model domain.
Definition: linteractionvolume.hh:301
Properties for a MPFA method.