MPLICcell.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::MPLICcell
26 
27 Description
28  Class performs geometric matching of volume fraction and calculates surface
29  interpolation of volume fraction field.
30 
31  Cut algorithms:
32  - Single cell cut
33  - Face-edge walk multiple cell cuts
34  - Tet decomposition cell cuts
35 
36 SourceFiles
37  MPLICcell.C
38 
39 \*---------------------------------------------------------------------------*/
40 
41 #ifndef MPLICcell_H
42 #define MPLICcell_H
43 
44 #include "MPLICface.H"
45 #include "MPLICcellStorage.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 /*---------------------------------------------------------------------------*\
53  Class MPLICcell Declaration
54 \*---------------------------------------------------------------------------*/
55 
56 class MPLICcell
57 {
58  // Private member data
59 
60  // Face related objects
61 
62  //- Select unweighted interpolation if true
63  // velocity flux corrected if false
64  const bool unweighted_;
65 
66  //- Select multi-cut if true or single-cut if false
67  const bool multiCut_;
68 
69  //- Surface interpolated alpha
70  DynamicList<scalar> alphaf_;
71 
72  //- Flux of alpha from point interpolated U
73  DynamicList<scalar> alphaPhiU_;
74 
75  //- Face flux from point interpolated U
76  DynamicList<scalar> phiU_;
77 
78 
79  // Geometry related objects
80 
81  //- Calculated volume fraction corresponding to the cut
82  scalar cutAlpha_;
83 
84  //- Cut surface area vector
85  vector cutSf_;
86 
87  //- Cut normal
88  vector cutNormal_;
89 
90  //- Face cutter
91  MPLICface faceCutter_;
92 
93  //- Sub cell volume
94  scalar subCellVolume_;
95 
96  //- Cut points for single cut
97  DynamicList<point> cutPoints_;
98 
99  //- Cut edge labels for single cut
100  DynamicList<label> cutEdges_;
101 
102  //- Submerged face areas
103  DynamicList<vector> subFaceAreas_;
104 
105  //- Submerged face areas
106  DynamicList<scalar> subFaceMagSf_;
107 
108  //- Submerged face centers
109  DynamicList<vector> subFaceCentres_;
110 
111  template<class Type>
112  class Vector4
113  :
114  public VectorSpace<Vector4<Type>, Type, 4>
115  {
116  public:
117 
118  // Constructors
119 
120  //- Construct null
121  inline Vector4()
122  {}
123 
124  inline Type a() const
125  {
126  return this->v_[0];
127  }
128 
129  inline Type b() const
130  {
131  return this->v_[1];
132  }
133 
134  inline Type c() const
135  {
136  return this->v_[2];
137  }
138 
139  inline Type d() const
140  {
141  return this->v_[3];
142  }
143 
144 
145  inline Type& a()
146  {
147  return this->v_[0];
148  }
149 
150  inline Type& b()
151  {
152  return this->v_[1];
153  }
154 
155  inline Type& c()
156  {
157  return this->v_[2];
158  }
159 
160  inline Type& d()
161  {
162  return this->v_[3];
163  }
164  };
165 
166  typedef Vector4<scalar> vector4;
167 
168  //- Four point alpha values for cubic polynomial fit
169  vector4 pCubicAlphas_;
170 
171  //- Four cell alpha values for cubic polynomial fit
172  vector4 cCubicAlphas_;
173 
174 
175  // Tetrahedron storage objects
176 
177  //- Tetrahedron alpha point values
178  scalarField tetPointsAlpha_;
179 
180  //- Tetrahedron velocity point values
181  vectorField tetPointsU_;
182 
183  //- Tetrahedron points
184  pointField tetPoints_;
185 
186  //- Tetrahedron surface area vectors
187  Vector4<vector> tetSf_;
188 
189  //- Tetrahedron face centres
190  Vector4<vector> tetCf_;
191 
192  //- Tetrahedron triangular faces addressing
193  const FixedList<face, 4> tetFaces_;
194 
195 
196  // Multicut addressing
197 
198  //- Is addressnig computed?
199  bool addressingCalculated_;
200 
201  //- Local edge faces
202  DynamicList<DynamicList<label>> localEdgeFaces_;
203 
204  //- Local face edges
205  DynamicList<DynamicList<label>> localFaceEdges_;
206 
207 
208  // Cell-point work arrays
209 
210  DynamicList<scalar> cellPointsAlpha_;
211 
212  DynamicList<scalar> pointsAlpha_;
213 
214 
215  // Private Member Functions
216 
217  //- Match cell volume ratio to phase fraction
218  // Returns:
219  // - +1: success
220  // - 0: fail
221  // - -1: base scheme
222  label calcMatchAlphaCutCell
223  (
224  const MPLICcellStorage& cellInfo,
225  const bool tetDecom = false
226  );
227 
228  //- Find in between which two point alpha values
229  // the target volume fraction lies
230  void findPointAlphaBounds
231  (
232  const MPLICcellStorage& cellInfo,
233  const bool tetDecom
234  );
235 
236  //- Calculate the two interior point alpha values for the cubic fit
237  void calcPointAlphaInterior
238  (
239  const MPLICcellStorage& cellInfo,
240  const bool tetDecom
241  );
242 
243  //- Solve the cubic fit
244  FixedList<scalar, 4> solveVanderMatrix() const;
245 
246  //- Identifying roots of cubic equation matching target volume fraction
247  void findRoots
248  (
249  const MPLICcellStorage& cellInfo,
250  const FixedList<scalar, 4>& coeffs,
251  const bool tetDecom
252  );
253 
254  //- Select simple cut or tet decomposition
255  scalar calcAlpha
256  (
257  const MPLICcellStorage& cellInfo,
258  const scalar target,
259  const bool tetDecom
260  );
261 
262  //- Calculate current sub-cell volume
263  void calcSubCellVolume();
264 
265  //- Calculate cell cut, volume and alpha
266  scalar calcCutCellVolumeAlpha
267  (
268  const MPLICcellStorage& cellInfo,
269  const scalar target
270  );
271 
272  //- Calculate cell cut, volume and alpha by tet decomposition
273  scalar calcTetCutCellVolumeAlpha
274  (
275  const MPLICcellStorage& cellInfo,
276  const scalar target
277  );
278 
279  //- Attempt single cut through the cell
280  // Returns:
281  // - 1: success
282  // - 0: fail
283  bool singleCutCell
284  (
285  const MPLICcellStorage& cellInfo,
286  const scalar target
287  );
288 
289  //- Attempt multiple cuts through the cell
290  bool multiCutCell
291  (
292  const MPLICcellStorage& cellInfo,
293  const scalar target
294  );
295 
296  //- Single cut through tet
297  bool cutTetCell
298  (
299  const scalar target,
300  const label faceOrig,
301  const bool ow
302  );
303 
304  //- Calculate local addressing for multi-cut
305  inline void calcAddressing
306  (
307  const MPLICcellStorage& cellInfo
308  );
309 
310  //- Append face area vectors and centers to cache
311  inline void appendSfCf
312  (
313  const vector& Sf,
314  const vector& Cf,
315  const scalar magSf,
316  const bool own = true
317  );
318 
319  //- Calculating surface vector of unordered edges
320  inline bool cutStatusCalcSf();
321 
322  //- Calculate cut surface area vector
323  inline vector calcCutSf() const;
324 
325  //- Calculating surface center of unordered edges
326  inline vector calcCutCf(const vector& cutSf) const;
327 
328  //- Calculate phiU
329  inline void phiU
330  (
331  const pointField& points,
332  const faceList& faces,
333  const labelList& cFaces,
334  const vectorField& pointsU
335  );
336 
337  //- Resize and set all the cell face fields to 0
338  inline void resetFaceFields(const label nFaces);
339 
340  //- Compute surface interpolation from area ratio
341  inline void calcAlphaf(const UIndirectList<scalar>& magSfs);
342 
343  //- Compute surface interpolation from flux ratio
344  inline void calcAlphaUf();
345 
346  //- Clear storage
347  inline void clear();
348 
349  //- Clear single cut storage
350  inline void clearOneCut();
351 
352 
353 public:
354 
355  // Constructors
356 
357  //- Construct for given interpolation and PLIC type
358  MPLICcell
359  (
360  const bool unweighted = true,
361  const bool multiCut = true
362  );
363 
364 
365  // Member functions
366 
367  //- Match cell volume ratios
368  bool matchAlpha
369  (
370  const MPLICcellStorage& cellInfo
371  );
372 
373  //- Return face volume fraction
374  inline const DynamicList<scalar>& alphaf() const;
375 
376  //- Return cut normal
377  inline const vector& cutNormal() const;
378 
379  //- Return volume fraction corresponding to the cut
380  inline scalar cutAlpha() const;
381 
382  //- Return sub-cell volume
383  inline scalar subCellVolume() const;
384 };
385 
386 
387 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
388 
389 } // End namespace Foam
390 
391 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
392 
393 #include "MPLICcellI.H"
394 
395 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
396 
397 #endif
398 
399 // ************************************************************************* //
tUEqn clear()
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:54
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
Templated vector space.
Definition: VectorSpace.H:53
Class that deals with cutting faces based on face point values and target value.
Definition: MPLICface.H:52
Holds information regarding type of cell. Used in inside/outside determination in cellClassification...
Definition: cellInfo.H:63
Provides local cell addressing for geometry and data for MPLIC class.
const dimensionedScalar c
Speed of light in a vacuum.
scalar subCellVolume() const
Return sub-cell volume.
Definition: MPLICcellI.H:281
const DynamicList< scalar > & alphaf() const
Return face volume fraction.
Definition: MPLICcellI.H:269
Class performs geometric matching of volume fraction and calculates surface interpolation of volume f...
Definition: MPLICcell.H:55
scalar cutAlpha() const
Return volume fraction corresponding to the cut.
Definition: MPLICcellI.H:287
const pointField & points
const vector & cutNormal() const
Return cut normal.
Definition: MPLICcellI.H:275
MPLICcell(const bool unweighted=true, const bool multiCut=true)
Construct for given interpolation and PLIC type.
Definition: MPLICcell.C:1002
bool matchAlpha(const MPLICcellStorage &cellInfo)
Match cell volume ratios.
Definition: MPLICcell.C:1027
Namespace for OpenFOAM.