3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
linteractionvolume3d.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_FVMPFAL3DINTERACTIONVOLUME_HH
21#define DUMUX_FVMPFAL3DINTERACTIONVOLUME_HH
22
23#include "properties.hh"
24
31namespace Dumux
32{
33
35// Mapper for local interaction volume indices (see doc/docextra/3dmpfa).
36class IndexTranslator
37{
38public:
39 enum
40 {
41 subVolumeTotalNum = 8,
42 fluxFacesTotalNum = 12,
43 fluxFacesNumOnSubVolume = 3,
44 fluxEdgesTotalNum = 6,
45 edgesNumOnFluxFace = 2
46 };
47
48 static int getFaceIndexFromSubVolume(int subVolumeIdx, int subVolumeFaceIdx)
49 {
50 return faceIndexOnSubVolume_[subVolumeIdx][subVolumeFaceIdx];
51 }
52
53 static int getEdgeIndexFromSubVolumeFace(int subVolumeIdx, int subVolumeFaceIdx, int subVolumeEdgeIdx)
54 {
55 return edgeIndexOnSubVolumeFace_[subVolumeIdx][subVolumeFaceIdx][subVolumeEdgeIdx];
56 }
57
58 static int getFaceIndexFromElements(int elem1Idx, int elem2Idx)
59 {
60 return faceIndexFromElements_[elem1Idx][elem2Idx];
61 }
62
63
64private:
65 static const int faceIndexOnSubVolume_[subVolumeTotalNum][fluxFacesNumOnSubVolume];
66 static const int edgeIndexOnSubVolumeFace_[subVolumeTotalNum][fluxFacesNumOnSubVolume][edgesNumOnFluxFace];
67 static const int faceIndexFromElements_[subVolumeTotalNum][subVolumeTotalNum];
68};
69
70const int IndexTranslator::faceIndexOnSubVolume_[subVolumeTotalNum][fluxFacesNumOnSubVolume] =
71{
72 {0, 3, 8},
73 {1, 0, 9},
74 {3, 2, 11},
75 {2, 1, 10},
76 {8, 4, 7},
77 {9, 5, 4},
78 {11, 7, 6},
79 {10, 6, 5}
80};
81
82const int IndexTranslator::edgeIndexOnSubVolumeFace_[subVolumeTotalNum][fluxFacesNumOnSubVolume][edgesNumOnFluxFace] =
83{
84 {{2, 0},{0, 5},{5, 2}},
85 {{3, 0},{0, 2},{2, 3}},
86 {{5, 0},{0, 4},{4, 5}},
87 {{4, 0},{0, 3},{3, 4}},
88 {{2, 5},{1, 2},{5, 1}},
89 {{3, 2},{1, 3},{2, 1}},
90 {{5, 4},{1, 5},{4, 1}},
91 {{4, 3},{1, 4},{3, 1}}
92};
93
94const int IndexTranslator::faceIndexFromElements_[subVolumeTotalNum][subVolumeTotalNum] =
95{
96 {-1, 0, 3, -1, 8, -1, -1, -1},
97 {0, -1, -1, 1, -1, 9, -1, -1},
98 {3, -1, -1, 2, -1, -1, 11, -1},
99 {-1, 1, 2, -1, -1, -1, -1, 10},
100 {8, -1, -1, -1, -1, 4, 7, -1},
101 {-1, 9, -1, -1, 4, -1, -1, 5},
102 {-1, -1, 11, -1, 7, -1, -1, 6},
103 {-1, -1, -1, 10, -1, 5, 6, -1}
104};
106
108
113template<class TypeTag>
115{
116private:
121
122 enum
123 {
124 dim = GridView::dimension,
125 dimWorld = GridView::dimensionworld,
126 };
127
128 using Element = typename GridView::template Codim<0>::Entity;
129 using ElementSeed = typename Grid::template Codim<0>::EntitySeed;
130
133 using PrimaryVariables = typename SolutionTypes::PrimaryVariables;
134
135 using DimVector = Dune::FieldVector<Scalar, dim>;
136 using FieldVectorVector = Dune::FieldVector<DimVector, dim>;
137 using FieldVectorVector2 = Dune::FieldVector<DimVector, 2>;
138 using FieldVectorVectorVector = Dune::FieldVector<FieldVectorVector2, dim>;
139 using IndexVector = Dune::FieldVector<int, dim>;
140 using BCTypeVector = std::vector<BoundaryTypes>;
141 using BCVector = std::vector<PrimaryVariables>;
142
143public:
144 //An interaction volume around a vertex includes in general 12 flux faces (see doc/docextra/3dmpfa).
145 //If the vertex is on the boundary these flux faces can be inside the domain on the boundary or outside the domain.
147 {
151 };
152
153 enum
154 {
155 subVolumeTotalNum = IndexTranslator::subVolumeTotalNum,
156 fluxFacesTotalNum = IndexTranslator::fluxFacesTotalNum,
157 fluxEdgesTotalNum = IndexTranslator::fluxEdgesTotalNum
158 };
159
161
164 : grid_(&grid)
165 , normal_(FieldVectorVector(DimVector(0.0)))
166 , facePos_(DimVector(0.0))
167 , edgePos_((DimVector(0.0)))
168 , faceArea_(0.0)
169 , faceType_(fluxFacesTotalNum, inside)
170 , indexOnElement_(IndexVector(0.0))
171 , elements_(subVolumeTotalNum)
172 , centerVertexPos_(0)
173 , elementNum_(0)
174 {}
175
177 void reset()
178 {
179 elements_.clear();
180 elements_.resize(subVolumeTotalNum);
181 centerVertexPos_ = 0;
182 indexOnElement_ = IndexVector(0.0);
183 faceType_.clear();
184 faceType_.resize(fluxFacesTotalNum, inside);
185 faceArea_ = 0.0;
186 facePos_ = DimVector(0.0);
187 edgePos_ = DimVector(0.0);
188 normal_ = FieldVectorVector(DimVector(0.0));
189 elementNum_ = 0;
190 boundaryTypes_.clear();
191 neumannValues_.clear();
192 dirichletValues_.clear();
193 }
194
196 /*
197 * \param centerVertexPos Position of the central vertex
198 */
199 void setCenterPosition(const DimVector &centerVertexPos)
200 {
201 centerVertexPos_ = centerVertexPos;
202 }
204
208 void setSubVolumeElement(const Element& element, int subVolumeIdx)
209 {
210 if (!hasSubVolumeElement(subVolumeIdx))
211 {
212 elements_[subVolumeIdx].push_back(element.seed());
213 elementNum_++;
214 }
215 else
216 {
217 elements_[subVolumeIdx].insert(elements_[subVolumeIdx].begin(), element.seed());
218 }
219 }
220
222
226 void setFacePosition(const DimVector& pos, int fIdx)
227 {
228 facePos_[fIdx] = pos;
229 }
230
232
236 void setEdgePosition(const DimVector& pos, int edgeIdx)
237 {
238 edgePos_[edgeIdx] = pos;
239 }
240
242
246 void setFaceArea(Scalar faceArea, int fIdx)
247 {
248 faceArea_[fIdx] = faceArea;
249 }
250
252
257 void setNormal(DimVector& normal, int subVolumeIdx, int subVolumeFaceIdx)
258 {
259 normal_[subVolumeIdx][subVolumeFaceIdx] = normal;
260 }
261
263
267 void setBoundary(BoundaryTypes& boundaryTypes, int subVolumeFaceIdx)
268 {
269 if (boundaryTypes_.size() == 0)
270 {
271 boundaryTypes_.resize(fluxFacesTotalNum);
272 }
273 boundaryTypes_[subVolumeFaceIdx] = boundaryTypes;
274 faceType_[subVolumeFaceIdx] = boundary;
275 }
276
278
281 void setOutsideFace(int subVolumeFaceIdx)
282 {
283 faceType_[subVolumeFaceIdx] = outside;
284 }
285
287
291 void setDirichletCondition(PrimaryVariables& condition, int subVolumeFaceIdx)
292 {
293 if (dirichletValues_.size() == 0)
294 {
295 dirichletValues_.resize(fluxFacesTotalNum, PrimaryVariables(0.0));
296 }
297 dirichletValues_[subVolumeFaceIdx] = condition;
298 }
299
301
305 void setNeumannCondition(PrimaryVariables& condition, int subVolumeFaceIdx)
306 {
307 if (neumannValues_.size() == 0)
308 {
309 neumannValues_.resize(fluxFacesTotalNum, PrimaryVariables(0.0));
310 }
311 neumannValues_[subVolumeFaceIdx] = condition;
312 }
313
315
320 void setIndexOnElement(int indexInInside, int subVolumeIdx, int subVolumeFaceIdx)
321 {
322 indexOnElement_[subVolumeIdx][subVolumeFaceIdx] = indexInInside;
323 }
324
326 DimVector& getCenterPosition()
327 {
328 return centerVertexPos_;
329 }
330
332 /*
333 * \param subVolumeIdx The local element index in the interaction volume
334 * \param subVolumeFaceIdx The local face index in the interaction volume element
335 *
336 * \return The interaction volume face index
337 */
338 int getFaceIndexFromSubVolume(int subVolumeIdx, int subVolumeFaceIdx)
339 {
340 return IndexTranslator::getFaceIndexFromSubVolume(subVolumeIdx, subVolumeFaceIdx);
341 }
342
344 /*
345 * \param subVolumeIdx The local element index in the interaction volume
346 * \param subVolumeFaceIdx The local face index in the interaction volume element
347 * \param subVolumeEdgeIdx The local edge index in the interaction volume element
348 *
349 * \return The interaction volume edge index
350 */
351 int getEdgeIndexFromSubVolumeFace(int subVolumeIdx, int subVolumeFaceIdx, int subVolumeEdgeIdx)
352 {
353 return IndexTranslator::getEdgeIndexFromSubVolumeFace(subVolumeIdx, subVolumeFaceIdx, subVolumeEdgeIdx);
354 }
355
358 {
359 return elementNum_;
360 }
361
363
369 int getIndexOnElement(int subVolumeIdx, int subVolumeFaceIdx)
370 {
371 return indexOnElement_[subVolumeIdx][subVolumeFaceIdx];
372 }
373
375
380 Element getSubVolumeElement(int subVolumeIdx)
381 {
382 if (hasSubVolumeElement(subVolumeIdx))
383 return grid_->entity(elements_[subVolumeIdx][0]);
384 else
385 {
386 std::cout<<"Problems when calling getSubVolumeElement("<<subVolumeIdx<<")\n";
388 DUNE_THROW(Dune::RangeError, "element not in interaction volume!");
389 }
390 }
391
393
398 bool hasSubVolumeElement(int subVolumeIdx)
399 {
400 return (elements_[subVolumeIdx].size() > 0);
401 }
402
404
409 BoundaryTypes& getBoundaryType(int subVolumeFaceIdx)
410 {
411 return boundaryTypes_[subVolumeFaceIdx];
412 }
413
415
420 bool isOutsideFace(int subVolumeFaceIdx)
421 {
422 return (faceType_[subVolumeFaceIdx] == outside);
423 }
424
426
431 bool isInsideFace(int subVolumeFaceIdx)
432 {
433 return (faceType_[subVolumeFaceIdx] == inside);
434 }
435
437
442 bool isBoundaryFace(int subVolumeFaceIdx)
443 {
444 return (faceType_[subVolumeFaceIdx] == boundary);
445 }
446
448
452 {
453 for (unsigned int i = 0; i < faceType_.size(); i++)
454 {
455 if (isOutsideFace(i) || isBoundaryFace(i))
456 return false;
457 }
458 return true;
459 }
460
462
466 {
467 return (boundaryTypes_.size() > 0);
468 }
469
471
476 PrimaryVariables& getDirichletValues(int subVolumeFaceIdx)
477 {
478 return dirichletValues_[subVolumeFaceIdx];
479 }
480
482
487 PrimaryVariables& getNeumannValues(int subVolumeFaceIdx)
488 {
489 return neumannValues_[subVolumeFaceIdx];
490 }
491
493
499 DimVector& getNormal(int subVolumeIdx, int subVolumeFaceIdx)
500 {
501 return normal_[subVolumeIdx][subVolumeFaceIdx];
502 }
503
505
511 DimVector& getFacePosition(int subVolumeIdx, int subVolumeFaceIdx)
512 {
513 return facePos_[IndexTranslator::getFaceIndexFromSubVolume(subVolumeIdx, subVolumeFaceIdx)];
514 }
515
517
524 DimVector& getEdgePosition(int subVolumeIdx, int subVolumeFaceIdx, int subVolumeEdgeIdx)
525 {
526 return edgePos_[IndexTranslator::getEdgeIndexFromSubVolumeFace(subVolumeIdx, subVolumeFaceIdx, subVolumeEdgeIdx)];
527 }
528
530
536 Scalar& getFaceArea(int subVolumeIdx, int subVolumeFaceIdx)
537 {
538 return faceArea_[IndexTranslator::getFaceIndexFromSubVolume(subVolumeIdx, subVolumeFaceIdx)];
539 }
540
542
547 DimVector& getFacePosition(int fIdx)
548 {
549 return facePos_[fIdx];
550 }
551
553
558 DimVector& getEdgePosition(int edgeIdx)
559 {
560 return edgePos_[edgeIdx];
561 }
562
564
569 Scalar& getFaceArea(int fIdx)
570 {
571 return faceArea_[fIdx];
572 }
573
576 {
577 std::cout<<"\nNumber of stored elements: "<<elementNum_<<"\n";
578 std::cout<<" center position: "<<centerVertexPos_<<"\n";
579 for (int i = 0; i < subVolumeTotalNum; i++)
580 {
581 if (elements_[i].size() > 0)
582 {
583 std::cout<<"element "<<i<<":\n";
584 std::cout<<"element level: "<<grid_->entity(elements_[i][0]).level()<<"\n";
585 std::cout<<"element position: "<<grid_->entity(elements_[i][0]).geometry().center()<<"\n";
586 std::cout<<"element volume: "<<grid_->entity(elements_[i][0]).geometry().volume()<<"\n";
587 std::cout<<"face indices on element: "<<indexOnElement_[i]<<"\n";
588 std::cout<<"face normals on element: "<<normal_[i]<<"\n";
589 std::cout<<"face areas on element: ";
590 for (int j = 0; j < 3; j++)
591 {
592 std::cout<<getFaceArea(i, j)<<" ";
593 }
594 std::cout<<"\n";
595 std::cout<<"face position on element: ";
596 for (int j = 0; j < 3; j++)
597 {
598 std::cout<<getFacePosition(i, j)<<" ";
599 }
600 std::cout<<"\n";
601 std::cout<<"edgePos:";
602 for (int j = 0; j < 3; j++)
603 {
604 std::cout<<" ";
605 for (int k = 0; k < 2; k++)
606 {
607 std::cout<<getEdgePosition(i, j, k) <<" ";
608 }
609 }
610 std::cout<<"\n";
611 }
612 }
613 }
614
615// void printInteractionVolumeInfoToFile(std::ofstream& dataFile)
616// {
617// dataFile<<"\nNumber of stored elements: "<<elementNum_<<"\n";
618// dataFile<<" center position: "<<centerVertexPos_<<"\n";
619// for (int i = 0; i < subVolumeTotalNum; i++)
620// {
621// if (elements_[i].size() > 0)
622// {
623// dataFile<<"element "<<i<<":\n";
624// dataFile<<"element level: "<<elements_[i][0]->level()<<"\n";
625// dataFile<<"element position: "<<elements_[i][0]->geometry().center()<<"\n";
626// dataFile<<"element volume: "<<elements_[i][0]->geometry().volume()<<"\n";
627// dataFile<<"face indices on element: "<<indexOnElement_[i]<<"\n";
628// dataFile<<"face normals on element: "<<normal_[i]<<"\n";
629// dataFile<<"face areas on element: ";
630// for (int j = 0; j < 3; j++)
631// {
632// dataFile<<getFaceArea(i, j)<<" ";
633// }
634// dataFile<<"\n";
635// dataFile<<"face position on element: ";
636// for (int j = 0; j < 3; j++)
637// {
638// dataFile<<getFacePosition(i, j)<<" ";
639// }
640// dataFile<<"\n";
641// dataFile<<"edgePos:";
642// for (int j = 0; j < 3; j++)
643// {
644// dataFile<<" ";
645// for (int k = 0; k < 2; k++)
646// {
647// dataFile<<getEdgePosition(i, j, k) <<" ";
648// }
649// }
650// dataFile<<"\n";
651//
652// }
653// dataFile<<"face types: ";
654// for (int i = 0; i < fluxFacesTotalNum; i++)
655// {
656// dataFile<<faceType_[i]<<" ";
657// }
658// dataFile<<"\n";
659// if (isBoundaryInteractionVolume())
660// {
661// dataFile<<"Boundaries:\n";
662// dataFile<<"dirichlet: ";
663// for (int i = 0; i < fluxFacesTotalNum; i++)
664// {
665// dataFile<<boundaryTypes_[i].isDirichlet(0)<<" ";
666// dataFile<<boundaryTypes_[i].isDirichlet(1)<<" ";
667// }
668// dataFile<<"\n";
669// dataFile<<"neumann: ";
670// for (int i = 0; i < fluxFacesTotalNum; i++)
671// {
672// dataFile<<boundaryTypes_[i].isNeumann(0)<<" ";
673// dataFile<<boundaryTypes_[i].isNeumann(1)<<" ";
674// }
675// dataFile<<"\n";
676// dataFile<<"values: ";
677// for (int i = 0; i < fluxFacesTotalNum; i++)
678// {
679// if (dirichletValues_.size() > 0)
680// dataFile<<dirichletValues_[i]<<" ";
681// if (neumannValues_.size() > 0)
682// dataFile<<neumannValues_[i]<<" ";
683// }
684// dataFile<<"\n";
685// }
686// }
687
688private:
689 const Grid* grid_;
690 Dune::FieldVector<FieldVectorVector, subVolumeTotalNum> normal_;
691 Dune::FieldVector<DimVector, fluxFacesTotalNum> facePos_;
692 Dune::FieldVector<DimVector, fluxEdgesTotalNum> edgePos_;
693 Dune::FieldVector<Scalar, fluxFacesTotalNum> faceArea_;
694 BCTypeVector boundaryTypes_;
695 std::vector<int> faceType_;
696 Dune::FieldVector<IndexVector, subVolumeTotalNum> indexOnElement_;
697 std::vector<std::vector<ElementSeed> > elements_;
698 BCVector neumannValues_;
699 BCVector dirichletValues_;
700 DimVector centerVertexPos_;
701 int elementNum_;
702};
703}
704#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 a 3d interaction volume of a MPFA L-method that does not change wi...
Definition: linteractionvolume3d.hh:115
int getIndexOnElement(int subVolumeIdx, int subVolumeFaceIdx)
The local dune face index dependent on the local interaction volume indices.
Definition: linteractionvolume3d.hh:369
DimVector & getEdgePosition(int edgeIdx)
The position of the edge center.
Definition: linteractionvolume3d.hh:558
bool isInsideFace(int subVolumeFaceIdx)
Check if an interaction volume face is inside the model domain (for boundary vertices)
Definition: linteractionvolume3d.hh:431
void setEdgePosition(const DimVector &pos, int edgeIdx)
Store the center of the edges.
Definition: linteractionvolume3d.hh:236
void setSubVolumeElement(const Element &element, int subVolumeIdx)
Store a dune element as a sub volume element.
Definition: linteractionvolume3d.hh:208
bool isBoundaryFace(int subVolumeFaceIdx)
Check if an interaction volume face is a model domain boundary (for boundary vertices)
Definition: linteractionvolume3d.hh:442
BoundaryTypes & getBoundaryType(int subVolumeFaceIdx)
The boundary types.
Definition: linteractionvolume3d.hh:409
void setNormal(DimVector &normal, int subVolumeIdx, int subVolumeFaceIdx)
Store the normals.
Definition: linteractionvolume3d.hh:257
PrimaryVariables & getDirichletValues(int subVolumeFaceIdx)
The Dirichlet boundary values.
Definition: linteractionvolume3d.hh:476
Element getSubVolumeElement(int subVolumeIdx)
The dune element of the interaction volume sub-volume.
Definition: linteractionvolume3d.hh:380
void setBoundary(BoundaryTypes &boundaryTypes, int subVolumeFaceIdx)
Store the types of boundary conditions.
Definition: linteractionvolume3d.hh:267
PrimaryVariables & getNeumannValues(int subVolumeFaceIdx)
The Neumann boundary values.
Definition: linteractionvolume3d.hh:487
void setIndexOnElement(int indexInInside, int subVolumeIdx, int subVolumeFaceIdx)
Store the local dune face index dependent on the local interaction volume indices.
Definition: linteractionvolume3d.hh:320
DimVector & getEdgePosition(int subVolumeIdx, int subVolumeFaceIdx, int subVolumeEdgeIdx)
The position of the edge center.
Definition: linteractionvolume3d.hh:524
void setFaceArea(Scalar faceArea, int fIdx)
Store the flux face areas.
Definition: linteractionvolume3d.hh:246
void setNeumannCondition(PrimaryVariables &condition, int subVolumeFaceIdx)
Store Neumann boundary conditions.
Definition: linteractionvolume3d.hh:305
bool isInnerVolume()
Check if an interaction volume is not located around a boundary vertex.
Definition: linteractionvolume3d.hh:451
DimVector & getNormal(int subVolumeIdx, int subVolumeFaceIdx)
The face normal.
Definition: linteractionvolume3d.hh:499
int getFaceIndexFromSubVolume(int subVolumeIdx, int subVolumeFaceIdx)
The local interaction volume face index dependent on the local sub volume indices.
Definition: linteractionvolume3d.hh:338
void printInteractionVolumeInfo()
Print the stored interaction volume data.
Definition: linteractionvolume3d.hh:575
DimVector & getCenterPosition()
The position of the central vertex.
Definition: linteractionvolume3d.hh:326
Scalar & getFaceArea(int fIdx)
The interaction volume flux face area.
Definition: linteractionvolume3d.hh:569
bool hasSubVolumeElement(int subVolumeIdx)
Check if a dune element is stored for an interaction volume sub-volume.
Definition: linteractionvolume3d.hh:398
bool isOutsideFace(int subVolumeFaceIdx)
Check if an interaction volume face is outside the model domain (for boundary vertices)
Definition: linteractionvolume3d.hh:420
int getEdgeIndexFromSubVolumeFace(int subVolumeIdx, int subVolumeFaceIdx, int subVolumeEdgeIdx)
The local interaction volume edge index dependent on the local sub volume indices.
Definition: linteractionvolume3d.hh:351
DimVector & getFacePosition(int fIdx)
The position of the face center.
Definition: linteractionvolume3d.hh:547
void setDirichletCondition(PrimaryVariables &condition, int subVolumeFaceIdx)
Store Dirichlet boundary conditions.
Definition: linteractionvolume3d.hh:291
FaceTypes
Definition: linteractionvolume3d.hh:147
@ outside
Flux face is outside the model domain.
Definition: linteractionvolume3d.hh:150
@ inside
Flux face is inside the model domain.
Definition: linteractionvolume3d.hh:148
@ boundary
Flux face is a boundary face.
Definition: linteractionvolume3d.hh:149
Scalar & getFaceArea(int subVolumeIdx, int subVolumeFaceIdx)
The interaction volume flux face area.
Definition: linteractionvolume3d.hh:536
int getElementNumber()
Number of stored dune elements.
Definition: linteractionvolume3d.hh:357
FvMpfaL3dInteractionVolume(const Grid &grid)
Constructs a FvMpfaL3dInteractionVolume object.
Definition: linteractionvolume3d.hh:163
bool isBoundaryInteractionVolume()
Check if an interaction volume is located around a boundary vertex.
Definition: linteractionvolume3d.hh:465
DimVector & getFacePosition(int subVolumeIdx, int subVolumeFaceIdx)
The position of the face center.
Definition: linteractionvolume3d.hh:511
void setCenterPosition(const DimVector &centerVertexPos)
Store position of the central vertex.
Definition: linteractionvolume3d.hh:199
void reset()
Reset the interaction volume (deletes stored data)
Definition: linteractionvolume3d.hh:177
void setOutsideFace(int subVolumeFaceIdx)
Define a flux face to be outside the model domain (for boundary vertices)
Definition: linteractionvolume3d.hh:281
@ fluxFacesTotalNum
Number of flux faces in the interaction volume.
Definition: linteractionvolume3d.hh:156
@ subVolumeTotalNum
Number of sub-volumes in the interaction volume.
Definition: linteractionvolume3d.hh:155
@ fluxEdgesTotalNum
Number of edges in the interaction volume.
Definition: linteractionvolume3d.hh:157
void setFacePosition(const DimVector &pos, int fIdx)
Store the position of a flux face.
Definition: linteractionvolume3d.hh:226
Base file for properties related to sequential models.