lineFace.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) 2011-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 "lineFace.H"
27 #include "meshSearch.H"
28 #include "meshBoundarySearch.H"
29 #include "sampledSetCloud.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace sampledSets
37 {
40 }
41 }
42 
43 
44 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
45 
46 void Foam::sampledSets::lineFace::calcSamples
47 (
48  const polyMesh& mesh,
49  const vector& start,
50  const vector& end,
51  const label storeFaces,
52  const bool storeCells,
53  DynamicList<point>& samplingPositions,
54  DynamicList<scalar>& samplingDistances,
55  DynamicList<label>& samplingSegments,
56  DynamicList<label>& samplingCells,
57  DynamicList<label>& samplingFaces
58 )
59 {
60  const meshSearch& searchEngine =
62  const meshBoundarySearch& boundarySearchEngine =
64 
65  // Get all candidates for starting the tracks
66  List<DynamicList<label>> procCandidateCells(Pstream::nProcs());
67  List<DynamicList<scalar>> procCandidateTs(Pstream::nProcs());
68  {
69  const label startCelli = searchEngine.findCell(start);
70  if (startCelli != -1)
71  {
72  procCandidateCells[Pstream::myProcNo()].append(startCelli);
73  procCandidateTs[Pstream::myProcNo()].append(0);
74  }
75 
76  const label endCelli = searchEngine.findCell(end);
77  if (endCelli != -1)
78  {
79  procCandidateCells[Pstream::myProcNo()].append(endCelli);
80  procCandidateTs[Pstream::myProcNo()].append(1);
81  }
82 
83  const List<pointIndexHit> bHits =
84  boundarySearchEngine.intersections(start, end);
85  forAll(bHits, bHiti)
86  {
87  for (label bHitj = bHiti + 1; bHitj < bHits.size(); ++ bHitj)
88  {
89  const point midP =
90  (bHits[bHiti].hitPoint() + bHits[bHitj].hitPoint())/2;
91 
92  const label midCelli = searchEngine.findCell(midP);
93  const scalar midT = mag(midP - start)/mag(end - start);
94 
95  if (midCelli != -1)
96  {
97  procCandidateCells[Pstream::myProcNo()].append(midCelli);
98  procCandidateTs[Pstream::myProcNo()].append(midT);
99  }
100  }
101  }
102  }
103  Pstream::gatherList(procCandidateCells);
104  Pstream::scatterList(procCandidateCells);
105  Pstream::gatherList(procCandidateTs);
106  Pstream::scatterList(procCandidateTs);
107 
108  // Tracking data
109  const List<point> startEnd({start, end});
110  const List<point> endStart({end, start});
111  DynamicList<point> halfSegmentPositions;
112  DynamicList<scalar> halfSegmentDistances;
113  DynamicList<label> halfSegmentCells;
114  DynamicList<label> halfSegmentFaces;
115 
116  // Create a cloud with which to track segments
117  sampledSetCloud particles
118  (
119  mesh,
122  );
123 
124  // Create each segment in turn
125  label segmenti = 0;
126  label nLocateBoundaryHits = 0;
127  forAll(procCandidateCells, proci)
128  {
129  forAll(procCandidateCells[proci], candidatei)
130  {
131  const label celli = procCandidateCells[proci][candidatei];
132 
133  if (celli == -1) continue;
134 
135  const scalar t = procCandidateTs[proci][candidatei];
136  const point p = (1 - t)*start + t*end;
137 
138  // Track in both directions to form parts of the segment either
139  // side of the candidate
140  Pair<scalar> segmentT;
141  forAll(segmentT, i)
142  {
144  (
145  particles,
146  i == 0 ? endStart : startEnd,
147  false,
148  storeFaces,
149  storeCells,
150  halfSegmentPositions,
151  halfSegmentDistances,
152  halfSegmentCells,
153  halfSegmentFaces
154  );
155 
156  particles.clear();
157  if (proci == Pstream::myProcNo())
158  {
159  particles.addParticle
160  (
162  (
163  searchEngine,
164  p,
165  celli,
166  nLocateBoundaryHits,
167  0,
168  i == 0 ? t : 1 - t,
169  0
170  )
171  );
172  }
173 
174  particles.move(particles, tdBwd);
175 
176  segmentT[i] =
178  (
179  particles.size()
180  ? mag(particles.first()->position(mesh) - start)
181  /mag(end - start)
182  : i == 0 ? vGreat : -vGreat,
183  [i](const scalar a, const scalar b)
184  {
185  return (i == 0) == (a < b) ? a : b;
186  }
187  );
188 
189  if (i == 0)
190  {
191  reverse(halfSegmentPositions);
192  reverse(halfSegmentDistances);
193  reverse(halfSegmentCells);
194  reverse(halfSegmentFaces);
195  }
196 
197  const label n = halfSegmentPositions.size();
198 
199  samplingPositions.append(halfSegmentPositions);
200  samplingSegments.append(labelList(n, segmenti));
201  samplingCells.append(halfSegmentCells);
202  samplingFaces.append(halfSegmentFaces);
203 
204  halfSegmentPositions.clear();
205  halfSegmentDistances.clear();
206  halfSegmentCells.clear();
207  halfSegmentFaces.clear();
208 
209  // If storing cells we need to store the starting cells between
210  // the tracks
211  if (proci == Pstream::myProcNo() && i == 0 && storeCells)
212  {
213  particle trackBwd
214  (
215  searchEngine,
216  p,
217  celli,
218  nLocateBoundaryHits
219  );
220  particle trackFwd(trackBwd);
221  trackBwd.trackToFace(mesh, start - p, 0);
222  trackFwd.trackToFace(mesh, end - p, 0);
223  if (trackBwd.onFace() && trackFwd.onFace())
224  {
225  samplingPositions.append
226  (
227  (
228  trackBwd.position(mesh)
229  + trackFwd.position(mesh)
230  )/2
231  );
232  samplingSegments.append(segmenti);
233  samplingCells.append(celli);
234  samplingFaces.append(-1);
235  }
236  }
237  }
238 
239  // Disable all candidates that fall within the bounds of the
240  // computed segment
241  forAll(procCandidateCells, procj)
242  {
243  forAll(procCandidateCells[procj], candidatej)
244  {
245  const label cellj = procCandidateCells[procj][candidatej];
246 
247  if (cellj == -1) continue;
248 
249  const scalar t = procCandidateTs[procj][candidatej];
250 
251  if
252  (
253  t > segmentT.first() - rootSmall
254  && t < segmentT.second() + rootSmall
255  )
256  {
257  procCandidateCells[procj][candidatej] = -1;
258  procCandidateTs[procj][candidatej] = NaN;
259  }
260  }
261  }
262 
263  // Move onto the next segment
264  segmenti ++;
265  }
266  }
267 
268  // Set the distances
269  samplingDistances.resize(samplingPositions.size());
270  forAll(samplingPositions, i)
271  {
272  samplingDistances[i] = mag(samplingPositions[i] - start);
273  }
274 }
275 
276 
277 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
278 
279 bool Foam::sampledSets::lineFace::calcSamples
280 (
281  DynamicList<point>& samplingPositions,
282  DynamicList<scalar>& samplingDistances,
283  DynamicList<label>& samplingSegments,
284  DynamicList<label>& samplingCells,
285  DynamicList<label>& samplingFaces
286 ) const
287 {
288  calcSamples
289  (
290  mesh(),
291  start_,
292  end_,
293  1,
294  false,
295  samplingPositions,
296  samplingDistances,
297  samplingSegments,
298  samplingCells,
299  samplingFaces
300  );
301 
302  // This set is ordered. Distances have been created.
303  return true;
304 }
305 
306 
307 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
308 
310 (
311  const word& name,
312  const polyMesh& mesh,
313  const dictionary& dict
314 )
315 :
317  start_(dict.lookup("start")),
318  end_(dict.lookup("end"))
319 {}
320 
321 
322 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
323 
325 {}
326 
327 
328 // ************************************************************************* //
label n
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:449
Macros for easy insertion into run-time selection tables.
static meshBoundarySearch & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:78
void resize(const label)
Alter the addressed list size.
Definition: DynamicListI.H:216
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
void clear()
Clear the addressed list, i.e. set the size to zero.
Definition: DynamicListI.H:236
Template class for intrusive linked lists.
Definition: ILList.H:67
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
const Type & second() const
Return second.
Definition: PairI.H:121
const Type & first() const
Return first.
Definition: PairI.H:107
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
T * first()
Return the first entry.
Definition: UILList.H:109
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
A list of keywords followed by any number of values (e.g. words and numbers) or sub-dictionaries.
Definition: dictionary.H:162
label size() const
Return the number of particles in the cloud.
Definition: Cloud.H:191
void addParticle(ParticleType *pPtr)
Transfer particle to cloud.
Definition: Cloud.C:227
void move(TrackCloudType &cloud, typename ParticleType::trackingData &td)
Move the particles.
Definition: Cloud.C:292
Mesh object that implements searches within the local boundary faces.
List< pointIndexHit > intersections(const point &pStart, const point &pEnd) const
Find all intersections of the boundary with the line.
Mesh object that implements searches within the local cells and faces.
Definition: meshSearch.H:59
label findCell(const point &p, const pointInCellShapes=pointInCellShapes::tets) const
Find the cell containing the given point.
Definition: meshSearch.C:173
static const meshSearch & New(const polyMesh &mesh, const pointInCellShapes=pointInCellShapes::tets)
Lookup or construct from mesh and cell decomposition option.
Definition: meshSearch.C:61
Base particle class.
Definition: particle.H:85
scalar trackToFace(const polyMesh &mesh, const vector &displacement, const scalar fraction)
As particle::track, but stops when a face is hit.
Definition: particle.C:203
bool onFace() const
Is the particle on a face?
Definition: particleI.H:135
vector position(const polyMesh &mesh) const
Return current particle position.
Definition: particleI.H:159
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:78
A Cloud of sampledSet particles.
Particle for generating line-type sampled sets.
Holds list of sampling points which is filled at construction time. Various implementations of this b...
Definition: sampledSet.H:64
const polyMesh & mesh() const
Access the mesh.
Definition: sampledSetI.H:36
Face-intersections along a line.
Definition: lineFace.H:92
lineFace(const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.
Definition: lineFace.C:310
virtual ~lineFace()
Destructor.
Definition: lineFace.C:324
Set of sets to sample. Call sampledSets.write() to sample&write files.
A class for handling words, derived from string.
Definition: word.H:63
Foam::fvMesh mesh(Foam::IOobject(regionName, runTime.name(), runTime, Foam::IOobject::MUST_READ), false)
volScalarField & b
Definition: createFields.H:27
defineTypeNameAndDebug(arcUniform, 0)
addToRunTimeSelectionTable(sampledSet, arcUniform, word)
const unitSet & lookup(const word &unitName)
Lookup and return the named unit from the table.
Definition: units.C:346
Namespace for OpenFOAM.
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
String typeName(const std::type_info &info)
Return the un-mangled name given the standard type info.
void reverse(UList< T > &, const label n)
Definition: UListI.H:334
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
word name(const LagrangianState state)
Return a string representation of a Lagrangian state enumeration.
tmp< DimensionedField< scalar, GeoMesh, Field > > mag(const DimensionedField< Type, GeoMesh, PrimitiveField > &df)
dictionary dict
volScalarField & p