reconstructionSchemes.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2019-2020 DLR
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "reconstructionSchemes.H"
29 #include "cutCellPLIC.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
37 }
38 
39 
40 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
41 
42 bool Foam::reconstructionSchemes::alreadyReconstructed(bool forceUpdate) const
43 {
44  const Time& runTime = alpha1_.mesh().time();
45 
46  label& curTimeIndex = timeIndexAndIter_.first();
47  label& curIter = timeIndexAndIter_.second();
48 
49  // Reset timeIndex and curIter
50  if (curTimeIndex < runTime.timeIndex())
51  {
52  curTimeIndex = runTime.timeIndex();
53  curIter = 0;
54  return false;
55  }
56 
57  if (forceUpdate)
58  {
59  curIter = 0;
60  return false;
61  }
62 
63  // Always reconstruct when subcycling
64  if (runTime.subCycling() != 0)
65  {
66  return false;
67  }
68 
69  ++curIter;
70  if (curIter > 1)
71  {
72  return true;
73  }
74 
75  return false;
76 }
77 
78 
79 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
80 
82 (
83  const word& type,
85  const surfaceScalarField& phi,
86  const volVectorField& U,
87  const dictionary& dict
88 )
89 :
90  IOdictionary
91  (
92  IOobject
93  (
94  "reconstructionScheme",
95  alpha1.time().constant(),
96  alpha1.db(),
97  IOobject::NO_READ,
98  IOobject::NO_WRITE
99  )
100  ),
101  reconstructionSchemesCoeffs_(dict),
102  alpha1_(alpha1),
103  phi_(phi),
104  U_(U),
105  normal_
106  (
107  IOobject
108  (
109  IOobject::groupName("interfaceNormal", alpha1.group()),
110  alpha1_.mesh().time().timeName(),
111  alpha1_.mesh(),
112  IOobject::NO_READ,
113  dict.getOrDefault("writeFields",false)
114  ? IOobject::AUTO_WRITE
115  : IOobject::NO_WRITE
116  ),
117  alpha1_.mesh(),
119  ),
120  centre_
121  (
122  IOobject
123  (
124  IOobject::groupName("interfaceCentre", alpha1.group()),
125  alpha1_.mesh().time().timeName(),
126  alpha1_.mesh(),
127  IOobject::NO_READ,
128  dict.getOrDefault("writeFields",false)
129  ? IOobject::AUTO_WRITE
130  : IOobject::NO_WRITE
131  ),
132  alpha1_.mesh(),
134  ),
135  interfaceCell_(alpha1_.mesh().nCells(), false),
136  interfaceLabels_(0.2*alpha1_.mesh().nCells()),
137  timeIndexAndIter_(0, 0)
138 {}
139 
140 
141 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
142 
144 {
145  reconstruct(false);
146  const fvMesh& mesh = centre_.mesh();
147 
148  cutCellPLIC cellCut(mesh);
149 
150  DynamicList<point> dynPts;
151  DynamicList<face> dynFaces(mesh.nCells()*0.1);
152  bitSet interfaceCellAddressing(mesh.nCells());
153 
154  forAll(interfaceCell_, celli)
155  {
156  if (interfaceCell_[celli])
157  {
158  if (mag(normal_[celli]) != 0)
159  {
160  interfaceCellAddressing.set(celli);
161  vector n = -normal_[celli]/mag(normal_[celli]);
162 
163  scalar cutVal = (centre_[celli] - mesh.C()[celli]) & n;
164  cellCut.calcSubCell(celli, cutVal, n);
165 
166  // cellCut.facePoints() are ordered and not connected
167  // to the other face
168  // append face with the starting label: dynPts.size()
169  dynFaces.append
170  (
171  face(identity(cellCut.facePoints().size(), dynPts.size()))
172  );
173  dynPts.append(cellCut.facePoints());
174  }
175  }
176  }
177 
178 
179  labelList meshCells(interfaceCellAddressing.sortedToc());
180 
181  // Transfer to mesh storage
182  pointField pts(std::move(dynPts));
183  faceList faces(std::move(dynFaces));
184 
185  return interface(std::move(pts), std::move(faces), std::move(meshCells));
186 }
187 
188 
189 // ************************************************************************* //
Original code supplied by Henning Scheufler, DLR (2019)
const T & first() const noexcept
Access the first element.
Definition: Pair.H:137
Pair< label > timeIndexAndIter_
Store timeindex/iteration to avoid multiple reconstruction.
dictionary dict
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:116
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:68
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
engineTime & runTime
interfaceProperties interface(alpha1, U, thermo->transportPropertiesDict())
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:82
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:69
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:421
constexpr const char *const group
Group name for atomic constants.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:81
reconstructionSchemes(const reconstructionSchemes &)=delete
No copy construct.
List< face > faceList
List of faces.
Definition: faceListFwd.H:39
word timeName
Definition: getTimeIndex.H:3
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: POSIX.C:799
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:38
dynamicFvMesh & mesh
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
Definition: labelLists.C:44
interface surface()
Generated interface surface points/faces.
void append(const T &val)
Copy append an element to the end of this list.
Definition: DynamicList.H:584
Class for cutting a cell, cellI, of an fvMesh, mesh_, at its intersection with an surface defined by ...
Definition: cutCellPLIC.H:67
volScalarField & alpha1_
Reference to the VoF Field.
const Mesh & mesh() const noexcept
Return mesh.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
defineTypeNameAndDebug(combustionModel, 0)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> reconstruct(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
U
Definition: pEqn.H:72
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Definition: bitSet.H:59
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
label nCells() const noexcept
Number of mesh cells.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
bool alreadyReconstructed(bool forceUpdate=true) const
Is the interface already up-to-date?
label n
const volVectorField & C() const
Return cell centres as volVectorField.
List< label > labelList
A List of labels.
Definition: List.H:62
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const T & second() const noexcept
Access the second element.
Definition: Pair.H:147
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.
const pointField & pts
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:127
const volScalarField & alpha1