triSurfaceMesh.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 | Copyright (C) 2011-2017 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 Foam::fileName Foam::triSurfaceMesh::checkFile
93 (
94  const regIOobject& io,
95  const dictionary& dict,
96  const bool isGlobal
97 )
98 {
99  fileName fName;
100  if (dict.readIfPresent("file", fName, false, false))
101  {
102  fName = relativeFilePath(io, fName, isGlobal);
103 
104  if (!exists(fName))
105  {
107  << "Cannot find triSurfaceMesh at " << fName
108  << exit(FatalError);
109  }
110  }
111  else
112  {
113  fName =
114  (
115  isGlobal
116  ? io.globalFilePath(typeName)
117  : io.localFilePath(typeName)
118  );
119 
120  if (!exists(fName))
121  {
123  << "Cannot find triSurfaceMesh starting from "
124  << io.objectPath() << exit(FatalError);
125  }
126  }
127 
128  return fName;
129 }
130 
131 
132 bool Foam::triSurfaceMesh::addFaceToEdge
133 (
134  const edge& e,
135  EdgeMap<label>& facesPerEdge
136 )
137 {
138  EdgeMap<label>::iterator eFnd = facesPerEdge.find(e);
139  if (eFnd != facesPerEdge.end())
140  {
141  if (eFnd() == 2)
142  {
143  return false;
144  }
145  eFnd()++;
146  }
147  else
148  {
149  facesPerEdge.insert(e, 1);
150  }
151  return true;
152 }
153 
154 
155 bool Foam::triSurfaceMesh::isSurfaceClosed() const
156 {
157  const pointField& pts = triSurface::points();
158 
159  // Construct pointFaces. Let's hope surface has compact point
160  // numbering ...
161  labelListList pointFaces;
162  invertManyToMany(pts.size(), *this, pointFaces);
163 
164  // Loop over all faces surrounding point. Count edges emanating from point.
165  // Every edge should be used by two faces exactly.
166  // To prevent doing work twice per edge only look at edges to higher
167  // point
168  EdgeMap<label> facesPerEdge(100);
169  forAll(pointFaces, pointi)
170  {
171  const labelList& pFaces = pointFaces[pointi];
172 
173  facesPerEdge.clear();
174  forAll(pFaces, i)
175  {
176  const triSurface::FaceType& f = triSurface::operator[](pFaces[i]);
177  label fp = findIndex(f, pointi);
178 
179  // Something weird: if I expand the code of addFaceToEdge in both
180  // below instances it gives a segmentation violation on some
181  // surfaces. Compiler (4.3.2) problem?
182 
183 
184  // Forward edge
185  label nextPointi = f[f.fcIndex(fp)];
186 
187  if (nextPointi > pointi)
188  {
189  bool okFace = addFaceToEdge
190  (
191  edge(pointi, nextPointi),
192  facesPerEdge
193  );
194 
195  if (!okFace)
196  {
197  return false;
198  }
199  }
200  // Reverse edge
201  label prevPointi = f[f.rcIndex(fp)];
202 
203  if (prevPointi > pointi)
204  {
205  bool okFace = addFaceToEdge
206  (
207  edge(pointi, prevPointi),
208  facesPerEdge
209  );
210 
211  if (!okFace)
212  {
213  return false;
214  }
215  }
216  }
217 
218  // Check for any edges used only once.
219  forAllConstIter(EdgeMap<label>, facesPerEdge, iter)
220  {
221  if (iter() != 2)
222  {
223  return false;
224  }
225  }
226  }
227 
228  return true;
229 }
230 
231 
232 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
233 
234 Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
235 :
236  searchableSurface(io),
238  (
239  IOobject
240  (
241  io.name(),
242  io.instance(),
243  io.local(),
244  io.db(),
245  io.readOpt(),
246  io.writeOpt(),
247  false // searchableSurface already registered under name
248  )
249  ),
250  triSurface(s),
252  minQuality_(-1),
253  surfaceClosed_(-1)
254 {
255  const pointField& pts = triSurface::points();
256 
257  bounds() = boundBox(pts);
258 }
259 
260 
261 Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io)
262 :
263  // Find instance for triSurfaceMesh
264  searchableSurface(io),
265  // Reused found instance in objectRegistry
267  (
268  IOobject
269  (
270  io.name(),
272  io.local(),
273  io.db(),
274  io.readOpt(),
275  io.writeOpt(),
276  false // searchableSurface already registered under name
277  )
278  ),
279  triSurface(checkFile(static_cast<const searchableSurface&>(*this), true)),
280  triSurfaceRegionSearch(static_cast<const triSurface&>(*this)),
281  minQuality_(-1),
282  surfaceClosed_(-1)
283 {
284  const pointField& pts = triSurface::points();
285 
286  bounds() = boundBox(pts);
287 }
288 
289 
290 Foam::triSurfaceMesh::triSurfaceMesh
291 (
292  const IOobject& io,
293  const dictionary& dict
294 )
295 :
296  searchableSurface(io),
297  // Reused found instance in objectRegistry
299  (
300  IOobject
301  (
302  io.name(),
304  io.local(),
305  io.db(),
306  io.readOpt(),
307  io.writeOpt(),
308  false // searchableSurface already registered under name
309  )
310  ),
311  triSurface
312  (
313  checkFile(static_cast<const searchableSurface&>(*this), dict, true)
314  ),
315  triSurfaceRegionSearch(static_cast<const triSurface&>(*this), dict),
316  minQuality_(-1),
317  surfaceClosed_(-1)
318 {
319  // Reading from supplied file name instead of objectPath/filePath
320  if (dict.readIfPresent("file", fName_, false, false))
321  {
322  fName_ = relativeFilePath
323  (
324  static_cast<const searchableSurface&>(*this),
325  fName_,
326  true
327  );
328  }
329 
330  scalar scaleFactor = 0;
331 
332  // Allow rescaling of the surface points
333  // eg, CAD geometries are often done in millimeters
334  if (dict.readIfPresent("scale", scaleFactor) && scaleFactor > 0)
335  {
336  Info<< searchableSurface::name() << " : using scale " << scaleFactor
337  << endl;
338  triSurface::scalePoints(scaleFactor);
339  }
340 
341  const pointField& pts = triSurface::points();
342 
343  bounds() = boundBox(pts);
344 
345  // Have optional minimum quality for normal calculation
346  if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0)
347  {
349  << " : ignoring triangles with quality < "
350  << minQuality_ << " for normals calculation." << endl;
351  }
352 }
353 
354 
355 Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const bool isGlobal)
356 :
357  // Find instance for triSurfaceMesh
358  searchableSurface(io),
359  // Reused found instance in objectRegistry
361  (
362  IOobject
363  (
364  io.name(),
366  io.local(),
367  io.db(),
368  io.readOpt(),
369  io.writeOpt(),
370  false // searchableSurface already registered under name
371  )
372  ),
373  triSurface
374  (
375  checkFile(static_cast<const searchableSurface&>(*this), isGlobal)
376  ),
377  triSurfaceRegionSearch(static_cast<const triSurface&>(*this)),
378  minQuality_(-1),
379  surfaceClosed_(-1)
380 {
381  const pointField& pts = triSurface::points();
382 
383  bounds() = boundBox(pts);
384 }
385 
386 
387 Foam::triSurfaceMesh::triSurfaceMesh
388 (
389  const IOobject& io,
390  const dictionary& dict,
391  const bool isGlobal
392 )
393 :
394  searchableSurface(io),
395  // Reused found instance in objectRegistry
397  (
398  IOobject
399  (
400  io.name(),
402  io.local(),
403  io.db(),
404  io.readOpt(),
405  io.writeOpt(),
406  false // searchableSurface already registered under name
407  )
408  ),
409  triSurface
410  (
411  checkFile(static_cast<const searchableSurface&>(*this), dict, isGlobal)
412  ),
413  triSurfaceRegionSearch(static_cast<const triSurface&>(*this), dict),
414  minQuality_(-1),
415  surfaceClosed_(-1)
416 {
417  // Reading from supplied file name instead of objectPath/filePath
418  if (dict.readIfPresent("file", fName_, false, false))
419  {
420  fName_ = relativeFilePath
421  (
422  static_cast<const searchableSurface&>(*this),
423  fName_,
424  isGlobal
425  );
426  }
427 
428  scalar scaleFactor = 0;
429 
430  // Allow rescaling of the surface points
431  // eg, CAD geometries are often done in millimeters
432  if (dict.readIfPresent("scale", scaleFactor) && scaleFactor > 0)
433  {
434  Info<< searchableSurface::name() << " : using scale " << scaleFactor
435  << endl;
436  triSurface::scalePoints(scaleFactor);
437  }
438 
439  const pointField& pts = triSurface::points();
440 
441  bounds() = boundBox(pts);
442 
443  // Have optional minimum quality for normal calculation
444  if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0)
445  {
447  << " : ignoring triangles with quality < "
448  << minQuality_ << " for normals calculation." << endl;
449  }
450 }
451 
452 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
453 
455 {
456  clearOut();
457 }
458 
459 
461 {
463  edgeTree_.clear();
465 }
466 
467 
468 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
469 
471 {
472  tmp<pointField> tPts(new pointField(8));
473  pointField& pt = tPts.ref();
474 
475  // Use copy to calculate face centres so they don't get stored
477  (
480  ).faceCentres();
481 
482  return tPts;
483 }
484 
485 
487 (
488  pointField& centres,
489  scalarField& radiusSqr
490 ) const
491 {
492  centres = coordinates();
493  radiusSqr.setSize(size());
494  radiusSqr = 0.0;
495 
496  const pointField& pts = triSurface::points();
497 
498  forAll(*this, facei)
499  {
500  const labelledTri& f = triSurface::operator[](facei);
501  const point& fc = centres[facei];
502  forAll(f, fp)
503  {
504  const point& pt = pts[f[fp]];
505  radiusSqr[facei] = max(radiusSqr[facei], Foam::magSqr(fc-pt));
506  }
507  }
508 
509  // Add a bit to make sure all points are tested inside
510  radiusSqr += Foam::sqr(SMALL);
511 }
512 
513 
515 {
516  return triSurface::points();
517 }
518 
519 
521 {
522  const indexedOctree<treeDataTriSurface>& octree = tree();
523 
524  labelList indices = octree.findBox(treeBoundBox(bb));
525 
526  return !indices.empty();
527 }
528 
529 
531 {
533  edgeTree_.clear();
534  triSurface::movePoints(newPoints);
535 }
536 
537 
540 {
541  if (edgeTree_.empty())
542  {
543  // Boundary edges
544  labelList bEdges
545  (
546  identity
547  (
548  nEdges()
549  -nInternalEdges()
550  )
551  + nInternalEdges()
552  );
553 
554  treeBoundBox bb(Zero, Zero);
555 
556  if (bEdges.size())
557  {
558  label nPoints;
560  (
561  *this,
562  bb,
563  nPoints
564  );
565 
566  // Random number generator. Bit dodgy since not exactly random ;-)
567  Random rndGen(65431);
568 
569  // Slightly extended bb. Slightly off-centred just so on symmetric
570  // geometry there are less face/edge aligned items.
571 
572  bb = bb.extend(rndGen, 1e-4);
573  bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
574  bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
575  }
576 
577  scalar oldTol = indexedOctree<treeDataEdge>::perturbTol();
579 
580  edgeTree_.reset
581  (
583  (
585  (
586  false, // cachebb
587  edges(), // edges
588  localPoints(), // points
589  bEdges // selected edges
590  ),
591  bb, // bb
592  maxTreeDepth(), // maxLevel
593  10, // leafsize
594  3.0 // duplicity
595  )
596  );
597 
599  }
600  return edgeTree_();
601 }
602 
603 
605 {
606  if (regions_.empty())
607  {
608  regions_.setSize(patches().size());
609  forAll(regions_, regionI)
610  {
611  regions_[regionI] = patches()[regionI].name();
612  }
613  }
614  return regions_;
615 }
616 
617 
619 {
620  if (surfaceClosed_ == -1)
621  {
622  if (isSurfaceClosed())
623  {
624  surfaceClosed_ = 1;
625  }
626  else
627  {
628  surfaceClosed_ = 0;
629  }
630  }
631 
632  return surfaceClosed_ == 1;
633 }
634 
635 
637 (
638  const pointField& samples,
639  const scalarField& nearestDistSqr,
640  List<pointIndexHit>& info
641 ) const
642 {
643  triSurfaceSearch::findNearest(samples, nearestDistSqr, info);
644 }
645 
646 
648 (
649  const pointField& samples,
650  const scalarField& nearestDistSqr,
651  const labelList& regionIndices,
652  List<pointIndexHit>& info
653 ) const
654 {
656  (
657  samples,
658  nearestDistSqr,
659  regionIndices,
660  info
661  );
662 }
663 
664 
666 (
667  const pointField& start,
668  const pointField& end,
669  List<pointIndexHit>& info
670 ) const
671 {
672  triSurfaceSearch::findLine(start, end, info);
673 }
674 
675 
677 (
678  const pointField& start,
679  const pointField& end,
680  List<pointIndexHit>& info
681 ) const
682 {
683  triSurfaceSearch::findLineAny(start, end, info);
684 }
685 
686 
688 (
689  const pointField& start,
690  const pointField& end,
692 ) const
693 {
694  triSurfaceSearch::findLineAll(start, end, info);
695 }
696 
697 
699 (
700  const List<pointIndexHit>& info,
701  labelList& region
702 ) const
703 {
704  region.setSize(info.size());
705  forAll(info, i)
706  {
707  if (info[i].hit())
708  {
709  region[i] = triSurface::operator[](info[i].index()).region();
710  }
711  else
712  {
713  region[i] = -1;
714  }
715  }
716 }
717 
718 
720 (
721  const List<pointIndexHit>& info,
722  vectorField& normal
723 ) const
724 {
725  const triSurface& s = *this;
726  const pointField& pts = s.points();
727 
728  normal.setSize(info.size());
729 
730  if (minQuality_ >= 0)
731  {
732  // Make sure we don't use triangles with low quality since
733  // normal is not reliable.
734 
735  const labelListList& faceFaces = s.faceFaces();
736 
737  forAll(info, i)
738  {
739  if (info[i].hit())
740  {
741  label facei = info[i].index();
742  normal[i] = s[facei].normal(pts);
743 
744  scalar qual = s[facei].tri(pts).quality();
745 
746  if (qual < minQuality_)
747  {
748  // Search neighbouring triangles
749  const labelList& fFaces = faceFaces[facei];
750 
751  forAll(fFaces, j)
752  {
753  label nbrI = fFaces[j];
754  scalar nbrQual = s[nbrI].tri(pts).quality();
755  if (nbrQual > qual)
756  {
757  qual = nbrQual;
758  normal[i] = s[nbrI].normal(pts);
759  }
760  }
761  }
762 
763  normal[i] /= mag(normal[i]) + VSMALL;
764  }
765  else
766  {
767  // Set to what?
768  normal[i] = Zero;
769  }
770  }
771  }
772  else
773  {
774  forAll(info, i)
775  {
776  if (info[i].hit())
777  {
778  label facei = info[i].index();
779  // Cached:
780  //normal[i] = faceNormals()[facei];
781 
782  // Uncached
783  normal[i] = s[facei].normal(pts);
784  normal[i] /= mag(normal[i]) + VSMALL;
785  }
786  else
787  {
788  // Set to what?
789  normal[i] = Zero;
790  }
791  }
792  }
793 }
794 
795 
797 {
799  (
801  (
802  IOobject
803  (
804  "values",
805  objectRegistry::time().timeName(), // instance
806  "triSurface", // local
807  *this,
810  ),
811  *this,
812  dimless,
813  labelField(values)
814  )
815  );
816 
817  // Store field on triMesh
818  fldPtr.ptr()->store();
819 }
820 
821 
823 (
824  const List<pointIndexHit>& info,
825  labelList& values
826 ) const
827 {
828  if (foundObject<triSurfaceLabelField>("values"))
829  {
830  values.setSize(info.size());
831 
832  const triSurfaceLabelField& fld = lookupObject<triSurfaceLabelField>
833  (
834  "values"
835  );
836 
837  forAll(info, i)
838  {
839  if (info[i].hit())
840  {
841  values[i] = fld[info[i].index()];
842  }
843  }
844  }
845 }
846 
847 
849 (
850  const pointField& points,
851  List<volumeType>& volType
852 ) const
853 {
854  volType.setSize(points.size());
855 
858 
859  forAll(points, pointi)
860  {
861  const point& pt = points[pointi];
862 
863  if (!tree().bb().contains(pt))
864  {
865  // Have to calculate directly as outside the octree
866  volType[pointi] = tree().shapes().getVolumeType(tree(), pt);
867  }
868  else
869  {
870  // - use cached volume type per each tree node
871  volType[pointi] = tree().getVolumeType(pt);
872  }
873  }
874 
876 }
877 
878 
880 (
884  const bool valid
885 ) const
886 {
887  fileName fullPath;
888  if (fName_.size())
889  {
890  // Override file name
891 
892  fullPath = fName_;
893 
894  fullPath.expand();
895  if (!fullPath.isAbsolute())
896  {
897  // Add directory from regIOobject
898  fullPath = searchableSurface::objectPath().path()/fullPath;
899  }
900  }
901  else
902  {
903  fullPath = searchableSurface::objectPath();
904  }
905 
906  if (!mkDir(fullPath.path()))
907  {
908  return false;
909  }
910 
911  triSurface::write(fullPath);
912 
913  if (!isFile(fullPath))
914  {
915  return false;
916  }
917 
918  return true;
919 }
920 
921 
922 // ************************************************************************* //
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
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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:313
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:291
A class for handling file names.
Definition: fileName.H:69
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:137
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:319
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:174
virtual ~triSurfaceMesh()
Destructor.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
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:253
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
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
const Field< point > & faceCentres() const
Return face centres for patch.
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.
virtual void scalePoints(const scalar)
Scale points. A non-positive factor is ignored.
Definition: triSurface.C:803
T * ptr()
Return object pointer for reuse.
Definition: autoPtrI.H:90
Macros for easy insertion into run-time selection tables.
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool valid) const
Write using given format, version and compression.
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:790
virtual void getVolumeType(const pointField &, List< volumeType > &) const
Determine type (inside/outside/mixed) for point. unknown if.
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))
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:58
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:124
void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &info) const
const fileName & local() const
Definition: IOobject.H:396
virtual void findLineAll(const pointField &start, const pointField &end, List< List< pointIndexHit >> &) const
Get all intersections in order from start to end.
cachedRandom rndGen(label(0), -1)
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
static void calcBounds(const PrimitivePatch< Face, FaceList, PointField, PointType > &p, boundBox &bb, label &nPoints)
const Field< point > & 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:91
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
Simple random number generator.
Definition: Random.H:49
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:193
bool isFile(const fileName &, const bool checkGzip=true, const bool followLink=true)
Does the name exist as a FILE in the file system?
Definition: POSIX.C:551
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const Field< point > & localPoints() const
Return pointField of points in patch.
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)
treeBoundBox extend(Random &, const scalar s) const
Return slightly wider bounding box.
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence 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:297
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 point & max() const
Maximum describing the bounding box.
Definition: boundBoxI.H:60
const geometricSurfacePatchList & patches() const
Definition: triSurface.H:314
writeOption writeOpt() const
Definition: IOobject.H:363
const fileName & instance() const
Definition: IOobject.H:386
triSurface()
Construct null.
Definition: triSurface.C:596
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
vector point
Point is a vector.
Definition: point.H:41
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.H:47
string & expand(const bool allowEmpty=false)
Expand initial tildes and all occurences of environment variables.
Definition: string.C:95
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:250
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
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
A class for managing temporary objects.
Definition: PtrList.H:53
Registry of regIOobjects.
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:331
friend class iterator
Declare friendship with the iterator.
Definition: HashTable.H:194
Triangulated surface description with patch information.
Definition: triSurface.H:65
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:170
fileName objectPath() const
Return complete path + object name.
Definition: IOobject.H:412
const indexedOctree< treeDataTriSurface > & tree() const
Demand driven construction of the octree.
virtual const wordList & regions() const
Names of regions.
bool exists(const fileName &, const bool checkGzip=true, const bool followLink=true)
Does the name exist (as DIRECTORY or FILE) in the file system?
Definition: POSIX.C:517
Namespace for OpenFOAM.
virtual void getField(const List< pointIndexHit > &, labelList &) const
WIP. From a set of hits (points and.