cellCellStencil.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) 2017-2021 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify i
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 "cellCellStencil.H"
30 #include "volFields.H"
31 #include "syncTools.H"
32 #include "globalIndex.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(cellCellStencil, 0);
39  defineRunTimeSelectionTable(cellCellStencil, mesh);
40 }
41 
42 const Foam::Enum
43 <
45 >
47 ({
48  { cellType::CALCULATED, "calculated" },
49  { cellType::INTERPOLATED, "interpolated" },
50  { cellType::HOLE, "hole" },
51 });
52 
53 
54 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
55 
56 Foam::cellCellStencil::cellCellStencil(const fvMesh& mesh)
57 :
58  mesh_(mesh),
59  nonInterpolatedFields_({"zoneID"})
60 {}
61 
62 
64 (
65  const fvMesh& mesh,
66  const dictionary& dict,
67  const bool update
68 )
69 {
70  DebugInFunction << "Constructing cellCellStencil" << endl;
71 
72  const word stencilType(dict.get<word>("method"));
73 
74  auto* ctorPtr = meshConstructorTable(stencilType);
75 
76  if (!ctorPtr)
77  {
79  (
80  dict,
81  "cellCellStencil",
82  stencilType,
83  *meshConstructorTablePtr_
84  ) << exit(FatalIOError);
85  }
86 
87  return autoPtr<cellCellStencil>(ctorPtr(mesh, dict, update));
88 }
89 
90 
91 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
92 
94 {}
95 
96 
97 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
98 
100 {
101  if (!mesh.foundObject<labelIOList>("zoneID"))
102  {
103  labelIOList* zoneIDPtr = new labelIOList
104  (
105  IOobject
106  (
107  "zoneID",
110  mesh,
113  ),
114  mesh.nCells()
115  );
116  labelIOList& zoneID = *zoneIDPtr;
117 
118  volScalarField volZoneID
119  (
120  IOobject
121  (
122  "zoneID",
123  mesh.time().findInstance(mesh.dbDir(), "zoneID"),
124  mesh,
127  false
128  ),
129  mesh
130  );
131  forAll(volZoneID, cellI)
132  {
133  zoneID[cellI] = label(volZoneID[cellI]);
134  }
135 
136  zoneIDPtr->store();
137  }
138  return mesh.lookupObject<labelIOList>("zoneID");
139 }
140 
141 
143 (
144  const label size,
145  const labelUList& lst
146 )
147 {
148  labelList count(size, Zero);
149  forAll(lst, i)
150  {
151  count[lst[i]]++;
152  }
154  return count;
155 }
156 
157 
159 {
160  return nonInterpolatedFields_;
161 }
162 
163 
165 {
166  return nonInterpolatedFields_;
167 }
168 
169 
171 {
172  forAll(slots, i)
173  {
174  if (slots[i] >= mesh_.nCells())
175  {
176  return false;
177  }
178  }
179  return true;
180 }
181 
182 
184 (
185  const globalIndex& gi,
186  const polyMesh& mesh,
187  const boolList& isValidCell,
188  const labelList& selectedCells,
189  labelListList& cellCells,
190  pointListList& cellCellCentres
191 )
192 {
193  // For selected cells determine the face neighbours (in global numbering)
194 
195  const pointField& cellCentres = mesh.cellCentres();
196  const labelList& faceOwner = mesh.faceOwner();
197  const labelList& faceNeighbour = mesh.faceNeighbour();
198  const cellList& cells = mesh.cells();
199 
200 
201  // 1. Determine global cell number on other side of coupled patches
202 
203  labelList globalCellIDs(identity(gi.localSize(), gi.localStart()));
204 
205  labelList nbrGlobalCellIDs;
207  (
208  mesh,
209  globalCellIDs,
210  nbrGlobalCellIDs
211  );
212  pointField nbrCellCentres;
214  (
215  mesh,
216  cellCentres,
217  nbrCellCentres
218  );
219 
220  boolList nbrIsValidCell;
222  (
223  mesh,
224  isValidCell,
225  nbrIsValidCell
226  );
227 
228 
229  // 2. Collect cell and all its neighbours
230 
231  cellCells.setSize(mesh.nCells());
232  cellCellCentres.setSize(cellCells.size());
233 
234  forAll(selectedCells, i)
235  {
236  label celli = selectedCells[i];
237 
238  const cell& cFaces = cells[celli];
239  labelList& stencil = cellCells[celli];
240  pointList& stencilPoints = cellCellCentres[celli];
241  stencil.setSize(cFaces.size()+1);
242  stencilPoints.setSize(stencil.size());
243  label compacti = 0;
244 
245  // First entry is cell itself
246  if (isValidCell[celli])
247  {
248  stencil[compacti] = globalCellIDs[celli];
249  stencilPoints[compacti++] = cellCentres[celli];
250  }
251 
252  // Other entries are cell neighbours
253  forAll(cFaces, i)
254  {
255  label facei = cFaces[i];
256  label bFacei = facei-mesh.nInternalFaces();
257  label own = faceOwner[facei];
258  label nbrCelli;
259  point nbrCc;
260  bool isValid = false;
261  if (bFacei >= 0)
262  {
263  nbrCelli = nbrGlobalCellIDs[bFacei];
264  nbrCc = nbrCellCentres[bFacei];
265  isValid = nbrIsValidCell[bFacei];
266  }
267  else
268  {
269  if (own != celli)
270  {
271  nbrCelli = gi.toGlobal(own);
272  nbrCc = cellCentres[own];
273  isValid = isValidCell[own];
274  }
275  else
276  {
277  label nei = faceNeighbour[facei];
278  nbrCelli = gi.toGlobal(nei);
279  nbrCc = cellCentres[nei];
280  isValid = isValidCell[nei];
281  }
282  }
283 
284  if (isValid)
285  {
286  SubList<label> current(stencil, compacti);
287  if (!current.found(nbrCelli))
288  {
289  stencil[compacti] = nbrCelli;
290  stencilPoints[compacti++] = nbrCc;
291  }
292  }
293  }
294  stencil.setSize(compacti);
295  stencilPoints.setSize(compacti);
296  }
297 }
298 
299 
300 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
volFields.H
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::Enum
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: IOstreamOption.H:57
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::cellCellStencil::localStencil
bool localStencil(const labelUList &) const
Helper: is stencil fully local.
Definition: cellCellStencil.C:170
Foam::globalIndex::localStart
label localStart() const
My local start.
Definition: globalIndexI.H:175
update
mesh update()
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::defineRunTimeSelectionTable
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
Foam::SubList
A List obtained as a section of another List.
Definition: SubList.H:54
Foam::polyMesh::dbDir
virtual const fileName & dbDir() const
Override the objectRegistry dbDir for a single-region case.
Definition: polyMesh.C:829
globalIndex.H
Foam::cellCellStencil::count
static labelList count(const label size, const labelUList &lst)
Count occurrences (in parallel)
Definition: cellCellStencil.C:143
Foam::cellCellStencil::nonInterpolatedFields
virtual const wordHashSet & nonInterpolatedFields() const
Return the names of any (stencil or mesh specific) fields that.
Definition: cellCellStencil.C:158
Foam::polyMesh::meshSubDir
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:321
Foam::primitiveMesh::cells
const cellList & cells() const
Definition: primitiveMeshCells.C:138
Foam::polyMesh::facesInstance
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:852
Foam::globalIndex::localSize
label localSize() const
My local size.
Definition: globalIndexI.H:187
Foam::FatalIOError
IOerror FatalIOError
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::dictionary::get
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:107
Foam::HashSet< word, Hash< word > >
Foam::labelIOList
IOList< label > labelIOList
Label container classes.
Definition: labelIOList.H:44
syncTools.H
Foam::objectRegistry::foundObject
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
Definition: objectRegistryTemplates.C:379
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::regIOobject::store
bool store()
Definition: regIOobjectI.H:37
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::primitiveMesh::nCells
label nCells() const noexcept
Number of mesh cells.
Definition: primitiveMeshI.H:96
FatalIOErrorInLookup
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
Definition: error.H:478
Foam::Field< vector >
Foam::cellCellStencil::globalCellCells
static void globalCellCells(const globalIndex &gi, const polyMesh &mesh, const boolList &isValidDonor, const labelList &selectedCells, labelListList &cellCells, pointListList &cellCellCentres)
Helper: create cell-cell addressing in global numbering.
Definition: cellCellStencil.C:184
DebugInFunction
#define DebugInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:388
Foam::List::setSize
void setSize(const label n)
Alias for resize()
Definition: List.H:222
Foam::cellCellStencil::zoneID
const labelIOList & zoneID() const
Helper: get reference to registered zoneID. Loads volScalarField.
Definition: cellCellStencil.H:210
Foam::polyMesh::faceOwner
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1107
Foam::objectRegistry::lookupObject
const Type & lookupObject(const word &name, const bool recursive=false) const
Definition: objectRegistryTemplates.C:434
Foam::cellCellStencil::~cellCellStencil
virtual ~cellCellStencil()
Destructor.
Definition: cellCellStencil.C:93
Foam::syncTools::swapBoundaryCellList
static void swapBoundaryCellList(const polyMesh &mesh, const UList< T > &cellData, List< T > &neighbourCellData)
Swap to obtain neighbour cell values for all boundary faces.
Definition: syncToolsTemplates.C:1392
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
zoneID
const labelIOList & zoneID
Definition: interpolatedFaces.H:22
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::globalIndex
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:68
Foam::Time::findInstance
word findInstance(const fileName &dir, const word &name=word::null, const IOobject::readOption rOpt=IOobject::MUST_READ, const word &stopInstance=word::null) const
Definition: Time.C:797
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::Pstream::listCombineGather
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
Definition: combineGatherScatter.C:290
Foam::autoPtr< Foam::cellCellStencil >
Foam::primitiveMesh::cellCentres
const vectorField & cellCentres() const
Definition: primitiveMeshCellCentresAndVols.C:84
Foam::BitOps::count
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
Definition: BitOps.H:77
Foam::Vector< scalar >
Foam::primitiveMesh::nInternalFaces
label nInternalFaces() const noexcept
Number of internal faces.
Definition: primitiveMeshI.H:78
Foam::List< label >
Foam::UList< label >
Foam::identity
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
Definition: labelList.C:38
Foam::IOList< label >
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:280
cells
const cellShapeList & cells
Definition: gmvOutputHeader.H:3
Foam::plusEqOp
Definition: ops.H:72
cellCellStencil.H
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::cellCellStencil::cellType
cellType
Definition: cellCellStencil.H:72
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::cellCellStencil::New
static autoPtr< cellCellStencil > New(const fvMesh &, const dictionary &dict, const bool update=true)
New function which constructs and returns pointer to a.
Definition: cellCellStencil.C:64
Foam::cellCellStencil::cellTypeNames_
static const Enum< cellType > cellTypeNames_
Mode type names.
Definition: cellCellStencil.H:85
Foam::cell
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:54
Foam::polyMesh::faceNeighbour
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1113
Foam::globalIndex::toGlobal
label toGlobal(const label i) const
From local to global index.
Definition: globalIndexI.H:240
Foam::IOobject::MUST_READ
Definition: IOobject.H:185