setAndPointToFaceZone.C
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) 2024-2025 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 \*---------------------------------------------------------------------------*/
25 
26 #include "setAndPointToFaceZone.H"
27 #include "FaceCellWave.H"
28 #include "faceZoneSet.H"
29 #include "meshSearch.H"
30 #include "minData.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(setAndPointToFaceZone, 0);
39  (
40  topoSetSource,
41  setAndPointToFaceZone,
42  word
43  );
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const polyMesh& mesh,
52  const word& setName,
53  const vector& point
54 )
55 :
56  topoSetSource(mesh),
57  setName_(setName),
58  point_(point)
59 {}
60 
61 
63 (
64  const polyMesh& mesh,
65  const dictionary& dict
66 )
67 :
68  topoSetSource(mesh),
69  setName_(dict.lookup<word>("set")),
70  point_(dict.lookup<vector>("point", dimLength))
71 {}
72 
73 
74 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
75 
77 {}
78 
79 
80 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
81 
83 (
84  const topoSetSource::setAction action,
85  topoSet& set
86 ) const
87 {
88  if (!isA<faceZoneSet>(set))
89  {
91  << "Operation only allowed on a faceZoneSet."
92  << endl;
93  return;
94  }
95 
96  const meshSearch& searchEngine = meshSearch::New(mesh_);
97 
98  // Cast to get the zone
99  faceZoneSet& fzSet = refCast<faceZoneSet>(set);
100 
101  // Load the set
102  faceSet loadedSet(mesh_, setName_);
103 
104  // Allocate new topology
105  DynamicList<label> newAddressing(fzSet.addressing().size());
106  DynamicList<bool> newFlipMap(fzSet.flipMap().size());
107 
108  if ((action == topoSetSource::NEW) || (action == topoSetSource::ADD))
109  {
110  Info<< " Adding all faces from faceSet " << setName_
111  << " ..." << endl;
112 
113  // Construct a wave to transport a minimum index
114  List<minData> allFaceData(mesh_.nFaces());
115  List<minData> allCellData(mesh_.nCells());
116 
117  // Set the existing faces of the set to a higher index than that of
118  // uninitialised data
119  forAllConstIter(faceSet, loadedSet, iter)
120  {
121  allFaceData[iter.key()] = minData(-1);
122  }
123 
124  // Find the cell containing the given point and give that an even
125  // higher index
126  const label celli = searchEngine.findCell(point_);
127  if (returnReduce(celli, maxOp<label>()) == -1)
128  {
130  << "Point " << point_ << " was not found in the mesh"
131  << exit(FatalError);
132  }
133  if (celli != -1)
134  {
135  allCellData[celli] = minData(0);
136  }
137 
138  // Create seed faces as all those of the cell containing the given
139  // point. Give these the same index as the cell.
140  const label nSeedFaces = celli != -1 ? mesh_.cells()[celli].size() : 0;
141  labelList seedFaces(nSeedFaces);
142  List<minData> seedFaceData(nSeedFaces);
143  if (celli != -1)
144  {
145  forAll(mesh_.cells()[celli], cellFacei)
146  {
147  const label facei = mesh_.cells()[celli][cellFacei];
148  seedFaces[cellFacei] = facei;
149  seedFaceData = minData(0);
150  }
151  }
152 
153  // Wave from the point to the set faces
154  FaceCellWave<minData> wave
155  (
156  mesh_,
157  seedFaces,
158  seedFaceData,
159  allFaceData,
160  allCellData,
161  mesh().globalData().nTotalCells()+1
162  );
163 
164  // Initialise as a copy
165  newAddressing = fzSet.addressing();
166  newFlipMap = fzSet.flipMap();
167 
168  // Unpack. Each set face should be at the boundary between cells that
169  // were and were not visited by the wave. If it is the neighbour that
170  // was visited by the wave, then set the flip map.
171  forAllConstIter(faceSet, loadedSet, iter)
172  {
173  const label facei = iter.key();
174 
175  const label owni = mesh_.faceOwner()[facei];
176  const bool ownValid = allCellData[owni].valid(wave.data());
177 
178  if (mesh_.isInternalFace(facei))
179  {
180  const label nbri = mesh_.faceNeighbour()[facei];
181  const bool nbrValid = allCellData[nbri].valid(wave.data());
182 
183  if (ownValid && nbrValid)
184  {
186  << "Internal face #" << facei << " at "
187  << mesh_.faceCentres()[facei]
188  << " was visited from both sides by a wave from "
189  << point_ << exit(FatalError);
190  }
191 
192  if (!ownValid && !nbrValid)
193  {
195  << "Internal face #" << facei << " at "
196  << mesh_.faceCentres()[facei]
197  << " was not visited by a wave from "
198  << point_ << exit(FatalError);
199  }
200  }
201  else
202  {
203  if (!ownValid)
204  {
206  << "Boundary face #" << facei << " at "
207  << mesh_.faceCentres()[facei]
208  << " was not visited by a wave from "
209  << point_ << exit(FatalError);
210  }
211  }
212 
213  newAddressing.append(facei);
214  newFlipMap.append(!ownValid);
215  }
216  }
217  else if (action == topoSetSource::DELETE)
218  {
219  Info<< " Removing all faces from faceSet " << setName_
220  << " ..." << endl;
221 
222  // Remove everything in the set from the zone. Don't have to worry
223  // about computing the flips here, seeing as we are only removing
224  // faces.
225  DynamicList<label> newAddressing(fzSet.addressing().size());
226  DynamicList<bool> newFlipMap(fzSet.flipMap().size());
227  forAll(fzSet.addressing(), i)
228  {
229  if (!loadedSet.found(fzSet.addressing()[i]))
230  {
231  newAddressing.append(fzSet.addressing()[i]);
232  newFlipMap.append(fzSet.flipMap()[i]);
233  }
234  }
235  }
236 
237  // Transfer into the zone
238  fzSet.addressing().transfer(newAddressing);
239  fzSet.flipMap().transfer(newFlipMap);
240  fzSet.updateSet();
241 }
242 
243 
244 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:449
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:492
Macros for easy insertion into run-time selection tables.
static const meshSearch & New(const polyMesh &mesh, const pointInCellShapes=pointInCellShapes::tets)
Lookup or construct from mesh and cell decomposition option.
Definition: meshSearch.C:61
Motion of the mesh specified as a list of pointMeshMovers.
virtual ~setAndPointToFaceZone()
Destructor.
setAndPointToFaceZone(const polyMesh &mesh, const word &setName, const vector &normal)
Construct from components.
virtual void applyToSet(const topoSetSource::setAction action, topoSet &) const
setAction
Enumeration defining the valid actions.
Definition: topoSetSource.H:83
Foam::fvMesh mesh(Foam::IOobject(regionName, runTime.name(), runTime, Foam::IOobject::MUST_READ), false)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
#define WarningInFunction
Report a warning using Foam::Warning.
label wave(const fvMesh &mesh, const List< labelPair > &changedPatchAndFaces, const label nCorrections, GeometricField< scalar, GeoMesh > &distance, TrackingData &td, GeometricField< DataType, GeoMesh > &... data)
Wave distance (and maybe additional) data from faces. If nCorrections is.
const unitSet & lookup(const word &unitName)
Lookup and return the named unit from the table.
Definition: units.C:346
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:288
const dimensionSet & dimLength
Definition: dimensions.C:141
messageStream Info
vector point
Point is a vector.
Definition: point.H:41
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
error FatalError
defineTypeNameAndDebug(atmosphericBoundaryLayer, 0)
dictionary dict