triSurfaceMesh.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-2021 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 "triSurfaceMesh.H"
27 #include "Random.H"
29 #include "EdgeMap.H"
30 #include "triSurfaceFields.H"
31 #include "Time.H"
32 #include "PatchTools.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(triSurfaceMesh, 0);
39  addToRunTimeSelectionTable(searchableSurface, triSurfaceMesh, dict);
40 }
41 
42 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
43 
44 Foam::fileName Foam::triSurfaceMesh::checkFile
45 (
46  const regIOobject& io,
47  const bool isGlobal
48 )
49 {
50  const fileName fName
51  (
52  isGlobal
53  ? io.globalFilePath(typeName)
54  : io.localFilePath(typeName)
55  );
56  if (fName.empty())
57  {
59  << "Cannot find triSurfaceMesh starting from "
60  << io.objectPath() << exit(FatalError);
61  }
62 
63  return fName;
64 }
65 
66 
67 Foam::fileName Foam::triSurfaceMesh::relativeFilePath
68 (
69  const regIOobject& io,
70  const fileName& f,
71  const bool isGlobal
72 )
73 {
74  fileName fName(f);
75  fName.expand();
76  if (!fName.isAbsolute())
77  {
78  // Is the specified file:
79  // - local to the cwd?
80  // - local to the case dir?
81  // - or just another name?
82  fName = fileHandler().filePath
83  (
84  isGlobal,
85  IOobject(io, fName),
87  );
88  }
89  return fName;
90 }
91 
92 
93 Foam::fileName Foam::triSurfaceMesh::checkFile
94 (
95  const regIOobject& io,
96  const dictionary& dict,
97  const bool isGlobal
98 )
99 {
100  fileName dictFName, fName;
101 
102  if (dict.readIfPresent("file", dictFName, false, false))
103  {
104  fName = relativeFilePath(io, dictFName, isGlobal);
105 
106  if (!exists(fName))
107  {
109  << "Cannot find triSurfaceMesh at " << io.path()/dictFName
110  << exit(FatalError);
111  }
112  }
113  else
114  {
115  fName =
116  (
117  isGlobal
118  ? io.globalFilePath(typeName)
119  : io.localFilePath(typeName)
120  );
121 
122  if (!exists(fName))
123  {
125  << "Cannot find triSurfaceMesh starting from "
126  << io.objectPath() << exit(FatalError);
127  }
128  }
129 
130  return fName;
131 }
132 
133 
134 bool Foam::triSurfaceMesh::addFaceToEdge
135 (
136  const edge& e,
137  EdgeMap<label>& facesPerEdge
138 )
139 {
140  EdgeMap<label>::iterator eFnd = facesPerEdge.find(e);
141  if (eFnd != facesPerEdge.end())
142  {
143  if (eFnd() == 2)
144  {
145  return false;
146  }
147  eFnd()++;
148  }
149  else
150  {
151  facesPerEdge.insert(e, 1);
152  }
153  return true;
154 }
155 
156 
157 bool Foam::triSurfaceMesh::isSurfaceClosed() const
158 {
159  const pointField& pts = triSurface::points();
160 
161  // Construct pointFaces. Let's hope surface has compact point
162  // numbering ...
163  labelListList pointFaces;
164  invertManyToMany(pts.size(), *this, pointFaces);
165 
166  // Loop over all faces surrounding point. Count edges emanating from point.
167  // Every edge should be used by two faces exactly.
168  // To prevent doing work twice per edge only look at edges to higher
169  // point
170  EdgeMap<label> facesPerEdge(100);
171  forAll(pointFaces, pointi)
172  {
173  const labelList& pFaces = pointFaces[pointi];
174 
175  facesPerEdge.clear();
176  forAll(pFaces, i)
177  {
178  const triSurface::FaceType& f = triSurface::operator[](pFaces[i]);
179  label fp = findIndex(f, pointi);
180 
181  // Something weird: if I expand the code of addFaceToEdge in both
182  // below instances it gives a segmentation violation on some
183  // surfaces. Compiler (4.3.2) problem?
184 
185 
186  // Forward edge
187  label nextPointi = f[f.fcIndex(fp)];
188 
189  if (nextPointi > pointi)
190  {
191  bool okFace = addFaceToEdge
192  (
193  edge(pointi, nextPointi),
194  facesPerEdge
195  );
196 
197  if (!okFace)
198  {
199  return false;
200  }
201  }
202  // Reverse edge
203  label prevPointi = f[f.rcIndex(fp)];
204 
205  if (prevPointi > pointi)
206  {
207  bool okFace = addFaceToEdge
208  (
209  edge(pointi, prevPointi),
210  facesPerEdge
211  );
212 
213  if (!okFace)
214  {
215  return false;
216  }
217  }
218  }
219 
220  // Check for any edges used only once.
221  forAllConstIter(EdgeMap<label>, facesPerEdge, iter)
222  {
223  if (iter() != 2)
224  {
225  return false;
226  }
227  }
228  }
229 
230  return true;
231 }
232 
233 
234 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
235 
237 :
238  searchableSurface(io),
240  (
241  IOobject
242  (
243  io.name(),
244  io.instance(),
245  io.local(),
246  io.db(),
247  io.readOpt(),
248  io.writeOpt(),
249  false // searchableSurface already registered under name
250  )
251  ),
252  triSurface(s),
254  minQuality_(-1),
255  surfaceClosed_(-1)
256 {
257  const pointField& pts = triSurface::points();
258 
259  bounds() = boundBox(pts);
260 }
261 
262 
264 :
265  // Find instance for triSurfaceMesh
266  searchableSurface(io),
267  // Reused found instance in objectRegistry
269  (
270  IOobject
271  (
272  io.name(),
274  io.local(),
275  io.db(),
276  io.readOpt(),
277  io.writeOpt(),
278  false // searchableSurface already registered under name
279  )
280  ),
281  triSurface(checkFile(static_cast<const searchableSurface&>(*this), true)),
282  triSurfaceRegionSearch(static_cast<const triSurface&>(*this)),
283  minQuality_(-1),
284  surfaceClosed_(-1)
285 {
286  const pointField& pts = triSurface::points();
287 
288  bounds() = boundBox(pts);
289 }
290 
291 
293 (
294  const IOobject& io,
295  const dictionary& dict
296 )
297 :
298  searchableSurface(io),
299  // Reused found instance in objectRegistry
301  (
302  IOobject
303  (
304  io.name(),
306  io.local(),
307  io.db(),
308  io.readOpt(),
309  io.writeOpt(),
310  false // searchableSurface already registered under name
311  )
312  ),
313  triSurface
314  (
315  checkFile(static_cast<const searchableSurface&>(*this), dict, true)
316  ),
317  triSurfaceRegionSearch(static_cast<const triSurface&>(*this), dict),
318  minQuality_(-1),
319  surfaceClosed_(-1)
320 {
321  // Reading from supplied file name instead of objectPath/filePath
322  if (dict.readIfPresent("file", fName_, false, false))
323  {
324  fName_ = relativeFilePath
325  (
326  static_cast<const searchableSurface&>(*this),
327  fName_,
328  true
329  );
330  }
331 
332  scalar scaleFactor = 0;
333 
334  // Allow rescaling of the surface points
335  // eg, CAD geometries are often done in millimeters
336  if (dict.readIfPresent("scale", scaleFactor) && scaleFactor > 0)
337  {
338  Info<< searchableSurface::name() << " : using scale " << scaleFactor
339  << endl;
340  triSurface::scalePoints(scaleFactor);
341  }
342 
343  const pointField& pts = triSurface::points();
344 
345  bounds() = boundBox(pts);
346 
347  // Have optional minimum quality for normal calculation
348  if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0)
349  {
351  << " : ignoring triangles with quality < "
352  << minQuality_ << " for normals calculation." << endl;
353  }
354 }
355 
356 
357 Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const bool isGlobal)
358 :
359  // Find instance for triSurfaceMesh
360  searchableSurface(io),
361  // Reused found instance in objectRegistry
363  (
364  IOobject
365  (
366  io.name(),
368  io.local(),
369  io.db(),
370  io.readOpt(),
371  io.writeOpt(),
372  false // searchableSurface already registered under name
373  )
374  ),
375  triSurface
376  (
377  checkFile(static_cast<const searchableSurface&>(*this), isGlobal)
378  ),
379  triSurfaceRegionSearch(static_cast<const triSurface&>(*this)),
380  minQuality_(-1),
381  surfaceClosed_(-1)
382 {
383  const pointField& pts = triSurface::points();
384 
385  bounds() = boundBox(pts);
386 }
387 
388 
390 (
391  const IOobject& io,
392  const dictionary& dict,
393  const bool isGlobal
394 )
395 :
396  searchableSurface(io),
397  // Reused found instance in objectRegistry
399  (
400  IOobject
401  (
402  io.name(),
404  io.local(),
405  io.db(),
406  io.readOpt(),
407  io.writeOpt(),
408  false // searchableSurface already registered under name
409  )
410  ),
411  triSurface
412  (
413  checkFile(static_cast<const searchableSurface&>(*this), dict, isGlobal)
414  ),
415  triSurfaceRegionSearch(static_cast<const triSurface&>(*this), dict),
416  minQuality_(-1),
417  surfaceClosed_(-1)
418 {
419  // Reading from supplied file name instead of objectPath/filePath
420  if (dict.readIfPresent("file", fName_, false, false))
421  {
422  fName_ = relativeFilePath
423  (
424  static_cast<const searchableSurface&>(*this),
425  fName_,
426  isGlobal
427  );
428  }
429 
430  scalar scaleFactor = 0;
431 
432  // Allow rescaling of the surface points
433  // eg, CAD geometries are often done in millimeters
434  if (dict.readIfPresent("scale", scaleFactor) && scaleFactor > 0)
435  {
436  Info<< searchableSurface::name() << " : using scale " << scaleFactor
437  << endl;
438  triSurface::scalePoints(scaleFactor);
439  }
440 
441  const pointField& pts = triSurface::points();
442 
443  bounds() = boundBox(pts);
444 
445  // Have optional minimum quality for normal calculation
446  if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0)
447  {
449  << " : ignoring triangles with quality < "
450  << minQuality_ << " for normals calculation." << endl;
451  }
452 }
453 
454 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
455 
457 {
458  clearOut();
459 }
460 
461 
463 {
465  edgeTree_.clear();
467 }
468 
469 
470 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
471 
473 {
474  tmp<pointField> tPts(new pointField(8));
475  pointField& pt = tPts.ref();
476 
477  // Use copy to calculate face centres so they don't get stored
479  (
482  ).faceCentres();
483 
484  return tPts;
485 }
486 
487 
489 (
490  pointField& centres,
491  scalarField& radiusSqr
492 ) const
493 {
494  centres = coordinates();
495  radiusSqr.setSize(size());
496  radiusSqr = 0.0;
497 
498  const pointField& pts = triSurface::points();
499 
500  forAll(*this, facei)
501  {
502  const labelledTri& f = triSurface::operator[](facei);
503  const point& fc = centres[facei];
504  forAll(f, fp)
505  {
506  const point& pt = pts[f[fp]];
507  radiusSqr[facei] = max(radiusSqr[facei], Foam::magSqr(fc-pt));
508  }
509  }
510 
511  // Add a bit to make sure all points are tested inside
512  radiusSqr += Foam::sqr(small);
513 }
514 
515 
517 {
518  return triSurface::points();
519 }
520 
521 
523 {
524  const indexedOctree<treeDataTriSurface>& octree = tree();
525 
526  labelList indices = octree.findBox(treeBoundBox(bb));
527 
528  return !indices.empty();
529 }
530 
531 
533 {
535  edgeTree_.clear();
536  triSurface::movePoints(newPoints);
537 }
538 
539 
542 {
543  if (edgeTree_.empty())
544  {
545  // Boundary edges
546  labelList bEdges
547  (
548  identity
549  (
550  nEdges()
551  -nInternalEdges()
552  )
553  + nInternalEdges()
554  );
555 
556  treeBoundBox bb(Zero, Zero);
557 
558  if (bEdges.size())
559  {
560  label nPoints;
562  (
563  *this,
564  bb,
565  nPoints
566  );
567 
568  // Slightly extended bb. Slightly off-centred just so on symmetric
569  // geometry there are less face/edge aligned items.
570  bb = bb.extend(1e-4);
571  }
572 
573  scalar oldTol = indexedOctree<treeDataEdge>::perturbTol();
575 
576  edgeTree_.reset
577  (
579  (
581  (
582  false, // cachebb
583  edges(), // edges
584  localPoints(), // points
585  bEdges // selected edges
586  ),
587  bb, // bb
588  maxTreeDepth(), // maxLevel
589  10, // leafsize
590  3.0 // duplicity
591  )
592  );
593 
595  }
596  return edgeTree_();
597 }
598 
599 
601 {
602  if (regions_.empty())
603  {
604  regions_.setSize(patches().size());
605  forAll(regions_, regionI)
606  {
607  regions_[regionI] = patches()[regionI].name();
608  }
609  }
610  return regions_;
611 }
612 
613 
615 {
616  if (surfaceClosed_ == -1)
617  {
618  if (isSurfaceClosed())
619  {
620  surfaceClosed_ = 1;
621  }
622  else
623  {
624  surfaceClosed_ = 0;
625  }
626  }
627 
628  return surfaceClosed_ == 1;
629 }
630 
631 
633 (
634  const pointField& samples,
635  const scalarField& nearestDistSqr,
636  List<pointIndexHit>& info
637 ) const
638 {
639  triSurfaceSearch::findNearest(samples, nearestDistSqr, info);
640 }
641 
642 
644 (
645  const pointField& samples,
646  const scalarField& nearestDistSqr,
647  const labelList& regionIndices,
648  List<pointIndexHit>& info
649 ) const
650 {
652  (
653  samples,
654  nearestDistSqr,
655  regionIndices,
656  info
657  );
658 }
659 
660 
662 (
663  const pointField& start,
664  const pointField& end,
665  List<pointIndexHit>& info
666 ) const
667 {
668  triSurfaceSearch::findLine(start, end, info);
669 }
670 
671 
673 (
674  const pointField& start,
675  const pointField& end,
676  List<pointIndexHit>& info
677 ) const
678 {
679  triSurfaceSearch::findLineAny(start, end, info);
680 }
681 
682 
684 (
685  const pointField& start,
686  const pointField& end,
688 ) const
689 {
690  triSurfaceSearch::findLineAll(start, end, info);
691 }
692 
693 
695 (
696  const List<pointIndexHit>& info,
697  labelList& region
698 ) const
699 {
700  region.setSize(info.size());
701  forAll(info, i)
702  {
703  if (info[i].hit())
704  {
705  region[i] = triSurface::operator[](info[i].index()).region();
706  }
707  else
708  {
709  region[i] = -1;
710  }
711  }
712 }
713 
714 
716 (
717  const List<pointIndexHit>& info,
718  vectorField& normal
719 ) const
720 {
721  const triSurface& s = *this;
722  const pointField& pts = s.points();
723 
724  normal.setSize(info.size());
725 
726  if (minQuality_ >= 0)
727  {
728  // Make sure we don't use triangles with low quality since
729  // normal is not reliable.
730 
731  const labelListList& faceFaces = s.faceFaces();
732 
733  forAll(info, i)
734  {
735  if (info[i].hit())
736  {
737  label facei = info[i].index();
738  normal[i] = s[facei].area(pts);
739 
740  scalar qual = s[facei].tri(pts).quality();
741 
742  if (qual < minQuality_)
743  {
744  // Search neighbouring triangles
745  const labelList& fFaces = faceFaces[facei];
746 
747  forAll(fFaces, j)
748  {
749  label nbrI = fFaces[j];
750  scalar nbrQual = s[nbrI].tri(pts).quality();
751  if (nbrQual > qual)
752  {
753  qual = nbrQual;
754  normal[i] = s[nbrI].area(pts);
755  }
756  }
757  }
758 
759  normal[i] /= mag(normal[i]) + vSmall;
760  }
761  else
762  {
763  // Set to what?
764  normal[i] = Zero;
765  }
766  }
767  }
768  else
769  {
770  forAll(info, i)
771  {
772  if (info[i].hit())
773  {
774  label facei = info[i].index();
775  // Cached:
776  // normal[i] = faceNormals()[facei];
777 
778  // Uncached
779  normal[i] = s[facei].normal(pts);
780  }
781  else
782  {
783  // Set to what?
784  normal[i] = Zero;
785  }
786  }
787  }
788 }
789 
790 
792 (
793  const pointField& points,
794  List<volumeType>& volType
795 ) const
796 {
797  volType.setSize(points.size());
798 
801 
802  forAll(points, pointi)
803  {
804  const point& pt = points[pointi];
805 
806  if (!tree().bb().contains(pt))
807  {
808  // Have to calculate directly as outside the octree
809  volType[pointi] = tree().shapes().getVolumeType(tree(), pt);
810  }
811  else
812  {
813  // - use cached volume type per each tree node
814  volType[pointi] = tree().getVolumeType(pt);
815  }
816  }
817 
819 }
820 
821 
823 {
825  (
827  (
828  IOobject
829  (
830  "values",
833  *this,
836  ),
837  *this,
838  dimless,
839  labelField(values)
840  )
841  );
842 
843  // Store field on triMesh
844  fldPtr.ptr()->store();
845 }
846 
847 
849 (
850  const List<pointIndexHit>& info,
851  labelList& values
852 ) const
853 {
854  if (foundObject<triSurfaceLabelField>("values"))
855  {
856  values.setSize(info.size());
857 
858  const triSurfaceLabelField& fld = lookupObject<triSurfaceLabelField>
859  (
860  "values"
861  );
862 
863  forAll(info, i)
864  {
865  if (info[i].hit())
866  {
867  values[i] = fld[info[i].index()];
868  }
869  }
870  }
871 }
872 
873 
875 (
879  const bool write
880 ) const
881 {
882  fileName fullPath;
883  if (fName_.size())
884  {
885  // Override file name
886 
887  fullPath = fName_;
888 
889  fullPath.expand();
890  if (!fullPath.isAbsolute())
891  {
892  // Add directory from regIOobject
893  fullPath = searchableSurface::objectPath().path()/fullPath;
894  }
895  }
896  else
897  {
898  fullPath = searchableSurface::objectPath();
899  }
900 
901  if (!mkDir(fullPath.path()))
902  {
903  return false;
904  }
905 
906  triSurface::write(fullPath);
907 
908  if (!isFile(fullPath))
909  {
910  return false;
911  }
912 
913  return true;
914 }
915 
916 
917 // ************************************************************************* //
virtual void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Find first intersection on segment from start to end.
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
label nPoints() const
Return number of points supporting patch faces.
label maxTreeDepth() const
Return max tree depth of octree.
Field< label > labelField
Specialisation of Field<T> for label.
Definition: labelField.H:49
bool exists(const fileName &, const bool checkVariants=true, const bool followLink=true)
Does the name exist (as directory or file) in the file system?
Definition: POSIX.C:520
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual tmp< pointField > coordinates() const
Get representative set of element coordinates.
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:325
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
const word & name() const
Return name.
Definition: IOobject.H:303
A class for handling file names.
Definition: fileName.H:79
triSurfaceMesh(const IOobject &, const triSurface &)
Construct from triSurface.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void clearOut()
Clear storage.
void findNearest(const pointField &samples, const scalarField &nearestDistSqr, const labelList &regionIndices, List< pointIndexHit > &info) const
Find the nearest point on the surface out of the regions.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void boundingSpheres(pointField &centres, scalarField &radiusSqr) const
Get bounding spheres (centre and radius squared). Any point.
T & operator[](const label)
Return element of UList.
Definition: UListI.H:167
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const Field< PointType > & faceCentres() const
Return face centres for patch.
triSurfaceRegionSearch(const triSurface &)
Construct from surface. Holds reference to surface!
Fields for triSurface.
virtual tmp< pointField > points() const
Get the points that define the surface.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
virtual ~triSurfaceMesh()
Destructor.
bool isFile(const fileName &, const bool checkVariants=true, const bool followLink=true)
Does the name exist as a file in the file system?
Definition: POSIX.C:555
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
virtual void findNearest(const pointField &sample, const scalarField &nearestDistSqr, List< pointIndexHit > &) const
virtual void setField(const labelList &values)
WIP. Store element-wise field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
static scalar & perturbTol()
Get the perturbation tolerance.
Helper class to search on triSurface. Creates an octree for each region of the surface and only searc...
virtual bool overlaps(const boundBox &bb) const
Does any part of the surface overlap the supplied bound box?
Holds data for octree to work on an edges subset.
Definition: treeDataEdge.H:53
label nInternalEdges() const
Number of internal edges.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:58
Base class of (analytical or triangulated) surface. Encapsulates all the search routines. WIP.
const dimensionSet dimless
virtual void scalePoints(const scalar)
Scale points. A non-positive factor is ignored.
Definition: triSurface.C:750
T * ptr()
Return object pointer for reuse.
Definition: autoPtrI.H:90
const Field< PointType > & localPoints() const
Return pointField of points in patch.
Macros for easy insertion into run-time selection tables.
void findNearest(const pointField &samples, const scalarField &nearestDistSqr, List< pointIndexHit > &info) const
virtual bool hasVolumeType() const
Whether supports volume type below. I.e. whether is closed.
virtual void movePoints(const pointField &)
Move points.
Definition: triSurface.C:737
virtual void getVolumeType(const pointField &, List< volumeType > &) const
Determine type (inside/outside/mixed) for point. unknown if.
objectRegistry(const Time &db, const label nIoObjects=128)
Construct the time objectRegistry given an initial estimate.
A list of faces which address into the list of points.
A List obtained as a section of another List.
Definition: SubList.H:53
const labelListList & faceFaces() const
Return face-face addressing.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool write=true) const
Write using given format, version and compression.
static void calcBounds(const PrimitivePatch< FaceList, PointField > &p, boundBox &bb, label &nPoints)
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &info) const
void findLineAll(const pointField &start, const pointField &end, List< List< pointIndexHit >> &info) const
Calculate all intersections from start to end.
bool isAbsolute() const
Return true if file name is absolute.
Definition: fileName.C:73
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
static const word & geometryDir()
Return the geometry directory name.
void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &info) const
const fileName & local() const
Definition: IOobject.H:400
virtual void findLineAll(const pointField &start, const pointField &end, List< List< pointIndexHit >> &) const
Get all intersections in order from start to end.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
streamFormat
Enumeration for the format of data in the stream.
Definition: IOstream.H:86
const Field< PointType > & points() const
Return reference to global points.
static const word null
An empty word.
Definition: word.H:77
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
word timeName
Definition: getTimeIndex.H:3
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const zero Zero
Definition: zero.H:97
virtual fileName filePath(const bool checkGlobal, const IOobject &, const word &typeName) const =0
Search for an object. checkGlobal : also check undecomposed case.
Triangle with additional region number.
Definition: labelledTri.H:57
const fileOperation & fileHandler()
Get current file handler.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
searchableSurface(const IOobject &io)
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:193
dimensioned< scalar > magSqr(const dimensioned< Type > &)
label nEdges() const
Return number of edges in patch.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const indexedOctree< treeDataEdge > & edgeTree() const
Demand driven construction of octree for boundary edges.
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:290
virtual void getRegion(const List< pointIndexHit > &, labelList &region) const
From a set of points and indices get the region.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void getNormal(const List< pointIndexHit > &, vectorField &normal) const
From a set of points and indices get the normal.
void setSize(const label)
Reset size of List.
Definition: List.C:281
const geometricSurfacePatchList & patches() const
Definition: triSurface.H:322
writeOption writeOpt() const
Definition: IOobject.H:363
const fileName & instance() const
Definition: IOobject.H:390
triSurface()
Construct null.
Definition: triSurface.C:534
std::remove_reference< ::Foam::List< labelledTri > >::type::value_type FaceType
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
string & expand(const bool allowEmpty=false)
Expand initial tildes and all occurrences of environment variables.
Definition: string.C:103
virtual label size() const
Range of local indices that can be returned.
virtual void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Return any intersection on segment from start to end.
virtual void movePoints(const pointField &)
Move points.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
scalar tolerance() const
Return tolerance to use in searches.
void invertManyToMany(const label len, const UList< InList > &, List< OutList > &)
Invert many-to-many.
Version number type.
Definition: IOstream.H:96
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
fileName path() const
Return directory path name (part before last /)
Definition: fileName.C:265
const boundBox & bounds() const
Return const reference to boundBox.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
A class for managing temporary objects.
Definition: PtrList.H:53
Registry of regIOobjects.
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:312
friend class iterator
Declare friendship with the iterator.
Definition: HashTable.H:194
Triangulated surface description with patch information.
Definition: triSurface.H:66
readOption readOpt() const
Definition: IOobject.H:353
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
label size() const
Return the number of elements in the UList.
Definition: ListI.H:171
fileName objectPath() const
Return complete path + object name.
Definition: IOobject.H:419
const indexedOctree< treeDataTriSurface > & tree() const
Demand driven construction of the octree.
virtual const wordList & regions() const
Names of regions.
treeBoundBox extend(const scalar s) const
Return asymetrically extended bounding box, with guaranteed.
Namespace for OpenFOAM.
virtual void getField(const List< pointIndexHit > &, labelList &) const
WIP. From a set of hits (points and.