sampledTriSurfaceMesh.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 "sampledTriSurfaceMesh.H"
27 #include "meshSearch.H"
28 #include "treeDataCell.H"
29 #include "treeDataFace.H"
30 #include "meshTools.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace sampledSurfaces
38 {
39  defineTypeNameAndDebug(triSurfaceMesh, 0);
40  addToRunTimeSelectionTable(sampledSurface, triSurfaceMesh, word);
41 
42  //- Private class for finding nearest
43  // Comprising:
44  // - global index
45  // - sqr(distance)
47 
49  {
50 
51  public:
52 
53  void operator()(nearInfo& x, const nearInfo& y) const
54  {
55  if (y.first() < x.first())
56  {
57  x = y;
58  }
59  }
60  };
61 }
62 }
63 
64 template<>
65 const char* Foam::NamedEnum
66 <
68  3
69 >::names[] =
70 {
71  "cells",
72  "insideCells",
73  "boundaryFaces"
74 };
75 
76 const Foam::NamedEnum
77 <
79  3
80 > Foam::sampledSurfaces::triSurfaceMesh::samplingSourceNames_;
81 
82 
83 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
84 
86 Foam::sampledSurfaces::triSurfaceMesh::nonCoupledboundaryTree() const
87 {
88  // Variant of meshSearch::boundaryTree() that only does non-coupled
89  // boundary faces.
90 
91  if (!boundaryTreePtr_.valid())
92  {
93  // all non-coupled boundary faces (not just walls)
94  const polyBoundaryMesh& patches = mesh().boundaryMesh();
95 
96  labelList bndFaces(mesh().nFaces()-mesh().nInternalFaces());
97  label bndI = 0;
98  forAll(patches, patchi)
99  {
100  const polyPatch& pp = patches[patchi];
101  if (!pp.coupled())
102  {
103  forAll(pp, i)
104  {
105  bndFaces[bndI++] = pp.start()+i;
106  }
107  }
108  }
109  bndFaces.setSize(bndI);
110 
111 
112  treeBoundBox overallBb(mesh().points());
113  overallBb = overallBb.extend(1e-4);
114 
115  boundaryTreePtr_.reset
116  (
118  (
119  treeDataFace // all information needed to search faces
120  (
121  false, // do not cache bb
122  mesh(),
123  bndFaces // boundary faces only
124  ),
125  overallBb, // overall search domain
126  8, // maxLevel
127  10, // leafsize
128  3.0 // duplicity
129  )
130  );
131  }
132 
133  return boundaryTreePtr_();
134 }
135 
136 
138 (
139  const meshSearch& meshSearcher
140 )
141 {
142  // Find the cells the triangles of the surface are in.
143  // Does approximation by looking at the face centres only
144  const pointField& fc = surface_.faceCentres();
145 
146  List<nearInfo> nearest(fc.size());
147 
148  // Global numbering for cells/faces - only used to uniquely identify local
149  // elements
150  globalIndex globalCells
151  (
152  (sampleSource_ == cells || sampleSource_ == insideCells)
153  ? mesh().nCells()
154  : mesh().nFaces()
155  );
156 
157  forAll(nearest, i)
158  {
159  nearest[i].first() = great;
160  nearest[i].second() = labelMax;
161  }
162 
163  if (sampleSource_ == cells)
164  {
165  // Search for nearest cell
166 
167  const indexedOctree<treeDataCell>& cellTree = meshSearcher.cellTree();
168 
169  forAll(fc, triI)
170  {
171  pointIndexHit nearInfo = cellTree.findNearest
172  (
173  fc[triI],
174  sqr(great)
175  );
176  if (nearInfo.hit())
177  {
178  nearest[triI].first() = magSqr(nearInfo.hitPoint()-fc[triI]);
179  nearest[triI].second() = globalCells.toGlobal(nearInfo.index());
180  }
181  }
182  }
183  else if (sampleSource_ == insideCells)
184  {
185  // Search for cell containing point
186 
187  const indexedOctree<treeDataCell>& cellTree = meshSearcher.cellTree();
188 
189  forAll(fc, triI)
190  {
191  if (cellTree.bb().contains(fc[triI]))
192  {
193  label index = cellTree.findInside(fc[triI]);
194  if (index != -1)
195  {
196  nearest[triI].first() = 0.0;
197  nearest[triI].second() = globalCells.toGlobal(index);
198  }
199  }
200  }
201  }
202  else
203  {
204  // Search on all non-coupled boundary faces
205 
206  const indexedOctree<treeDataFace>& bTree = nonCoupledboundaryTree();
207 
208  forAll(fc, triI)
209  {
210  pointIndexHit nearInfo = bTree.findNearest
211  (
212  fc[triI],
213  sqr(great)
214  );
215  if (nearInfo.hit())
216  {
217  nearest[triI].first() = magSqr(nearInfo.hitPoint()-fc[triI]);
218  nearest[triI].second() = globalCells.toGlobal
219  (
220  bTree.shapes().faceLabels()[nearInfo.index()]
221  );
222  }
223  }
224  }
225 
226 
227  // See which processor has the nearest. Mark and subset
228  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
229 
230  Pstream::listCombineGather(nearest, nearestEqOp());
232 
233  labelList cellOrFaceLabels(fc.size(), -1);
234 
235  label nFound = 0;
236  forAll(nearest, triI)
237  {
238  if (nearest[triI].second() == labelMax)
239  {
240  // Not found on any processor. How to map?
241  }
242  else if (globalCells.isLocal(nearest[triI].second()))
243  {
244  cellOrFaceLabels[triI] = globalCells.toLocal
245  (
246  nearest[triI].second()
247  );
248  nFound++;
249  }
250  }
251 
252 
253  if (debug)
254  {
255  Pout<< "Local out of faces:" << cellOrFaceLabels.size()
256  << " keeping:" << nFound << endl;
257  }
258 
259  // Now subset the surface. Do not use triSurface::subsetMesh since requires
260  // original surface to be in compact numbering.
261 
262  const triSurface& s = surface_;
263 
264  // Compact to original triangle
265  labelList faceMap(s.size());
266  // Compact to original points
267  labelList pointMap(s.points().size());
268  // From original point to compact points
269  labelList reversePointMap(s.points().size(), -1);
270 
271  {
272  label newPointi = 0;
273  label newFacei = 0;
274 
275  forAll(s, facei)
276  {
277  if (cellOrFaceLabels[facei] != -1)
278  {
279  faceMap[newFacei++] = facei;
280 
281  const triSurface::FaceType& f = s[facei];
282  forAll(f, fp)
283  {
284  if (reversePointMap[f[fp]] == -1)
285  {
286  pointMap[newPointi] = f[fp];
287  reversePointMap[f[fp]] = newPointi++;
288  }
289  }
290  }
291  }
292  faceMap.setSize(newFacei);
293  pointMap.setSize(newPointi);
294  }
295 
296 
297  // Subset cellOrFaceLabels
298  cellOrFaceLabels = UIndirectList<label>(cellOrFaceLabels, faceMap)();
299 
300  // Store any face per point (without using pointFaces())
301  labelList pointToFace(pointMap.size());
302 
303  // Create faces and points for subsetted surface
304  faceList& faces = this->storedFaces();
305  faces.setSize(faceMap.size());
306  forAll(faceMap, i)
307  {
308  const triFace& f = s[faceMap[i]];
309  triFace newF
310  (
311  reversePointMap[f[0]],
312  reversePointMap[f[1]],
313  reversePointMap[f[2]]
314  );
315  faces[i] = newF.triFaceFace();
316 
317  forAll(newF, fp)
318  {
319  pointToFace[newF[fp]] = i;
320  }
321  }
322 
323  this->storedPoints() = pointField(s.points(), pointMap);
324 
325  if (debug)
326  {
327  print(Pout);
328  Pout<< endl;
329  }
330 
331 
332 
333  // Collect the samplePoints and sampleElements
334  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
335 
337  {
338  samplePoints_.setSize(pointMap.size());
339  sampleElements_.setSize(pointMap.size(), -1);
340 
341  if (sampleSource_ == cells)
342  {
343  // samplePoints_ : per surface point a location inside the cell
344  // sampleElements_ : per surface point the cell
345 
346  forAll(points(), pointi)
347  {
348  const point& pt = points()[pointi];
349  label celli = cellOrFaceLabels[pointToFace[pointi]];
350  sampleElements_[pointi] = celli;
351 
352  // Check if point inside cell
353  if
354  (
355  mesh().pointInCell
356  (
357  pt,
358  sampleElements_[pointi],
359  meshSearcher.decompMode()
360  )
361  )
362  {
363  samplePoints_[pointi] = pt;
364  }
365  else
366  {
367  // Find nearest point on faces of cell
368  const cell& cFaces = mesh().cells()[celli];
369 
370  scalar minDistSqr = vGreat;
371 
372  forAll(cFaces, i)
373  {
374  const face& f = mesh().faces()[cFaces[i]];
375  pointHit info = f.nearestPoint(pt, mesh().points());
376  if (info.distance() < minDistSqr)
377  {
378  minDistSqr = info.distance();
379  samplePoints_[pointi] = info.rawPoint();
380  }
381  }
382  }
383  }
384  }
385  else if (sampleSource_ == insideCells)
386  {
387  // samplePoints_ : per surface point a location inside the cell
388  // sampleElements_ : per surface point the cell
389 
390  forAll(points(), pointi)
391  {
392  const point& pt = points()[pointi];
393  label celli = cellOrFaceLabels[pointToFace[pointi]];
394  sampleElements_[pointi] = celli;
395  samplePoints_[pointi] = pt;
396  }
397  }
398  else
399  {
400  // samplePoints_ : per surface point a location on the boundary
401  // sampleElements_ : per surface point the boundary face containing
402  // the location
403 
404  forAll(points(), pointi)
405  {
406  const point& pt = points()[pointi];
407  label facei = cellOrFaceLabels[pointToFace[pointi]];
408  sampleElements_[pointi] = facei;
409  samplePoints_[pointi] = mesh().faces()[facei].nearestPoint
410  (
411  pt,
412  mesh().points()
413  ).rawPoint();
414  }
415  }
416  }
417  else
418  {
419  // if sampleSource_ == cells:
420  // samplePoints_ : n/a
421  // sampleElements_ : per surface triangle the cell
422  // if sampleSource_ == insideCells:
423  // samplePoints_ : n/a
424  // sampleElements_ : -1 or per surface triangle the cell
425  // else:
426  // samplePoints_ : n/a
427  // sampleElements_ : per surface triangle the boundary face
428  samplePoints_.clear();
429  sampleElements_.transfer(cellOrFaceLabels);
430  }
431 
432 
433  if (debug)
434  {
435  this->clearOut();
436  OFstream str(mesh().time().path()/"surfaceToMesh.obj");
437  Info<< "Dumping correspondence from local surface (points or faces)"
438  << " to mesh (cells or faces) to " << str.name() << endl;
439  label vertI = 0;
440 
442  {
443  if (sampleSource_ == cells || sampleSource_ == insideCells)
444  {
445  forAll(samplePoints_, pointi)
446  {
447  meshTools::writeOBJ(str, points()[pointi]);
448  vertI++;
449 
450  meshTools::writeOBJ(str, samplePoints_[pointi]);
451  vertI++;
452 
453  label celli = sampleElements_[pointi];
454  meshTools::writeOBJ(str, mesh().cellCentres()[celli]);
455  vertI++;
456  str << "l " << vertI-2 << ' ' << vertI-1 << ' ' << vertI
457  << nl;
458  }
459  }
460  else
461  {
462  forAll(samplePoints_, pointi)
463  {
464  meshTools::writeOBJ(str, points()[pointi]);
465  vertI++;
466 
467  meshTools::writeOBJ(str, samplePoints_[pointi]);
468  vertI++;
469 
470  label facei = sampleElements_[pointi];
471  meshTools::writeOBJ(str, mesh().faceCentres()[facei]);
472  vertI++;
473  str << "l " << vertI-2 << ' ' << vertI-1 << ' ' << vertI
474  << nl;
475  }
476  }
477  }
478  else
479  {
480  if (sampleSource_ == cells || sampleSource_ == insideCells)
481  {
482  forAll(sampleElements_, triI)
483  {
484  meshTools::writeOBJ(str, faceCentres()[triI]);
485  vertI++;
486 
487  label celli = sampleElements_[triI];
488  meshTools::writeOBJ(str, mesh().cellCentres()[celli]);
489  vertI++;
490  str << "l " << vertI-1 << ' ' << vertI << nl;
491  }
492  }
493  else
494  {
495  forAll(sampleElements_, triI)
496  {
497  meshTools::writeOBJ(str, faceCentres()[triI]);
498  vertI++;
499 
500  label facei = sampleElements_[triI];
501  meshTools::writeOBJ(str, mesh().faceCentres()[facei]);
502  vertI++;
503  str << "l " << vertI-1 << ' ' << vertI << nl;
504  }
505  }
506  }
507  }
508 
509 
510  needsUpdate_ = false;
511  return true;
512 }
513 
514 
515 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
516 
518 (
519  const word& name,
520  const polyMesh& mesh,
521  const word& surfaceName,
522  const samplingSource sampleSource
523 )
524 :
525  sampledSurface(name, mesh),
526  surface_
527  (
528  IOobject
529  (
530  surfaceName,
531  mesh.time().constant(),
533  mesh,
536  false
537  )
538  ),
539  sampleSource_(sampleSource),
540  needsUpdate_(true),
541  sampleElements_(0),
542  samplePoints_(0)
543 {}
544 
545 
547 (
548  const word& name,
549  const polyMesh& mesh,
550  const dictionary& dict
551 )
552 :
553  sampledSurface(name, mesh, dict),
554  surface_
555  (
556  IOobject
557  (
558  dict.lookup("surface"),
559  mesh.time().constant(),
561  mesh,
564  false
565  )
566  ),
567  sampleSource_(samplingSourceNames_[dict.lookup("source")]),
568  needsUpdate_(true),
569  sampleElements_(0),
570  samplePoints_(0)
571 {}
572 
573 
575 (
576  const word& name,
577  const polyMesh& mesh,
578  const triSurface& surface,
579  const word& sampleSourceName
580 )
581 :
582  sampledSurface(name, mesh),
583  surface_
584  (
585  IOobject
586  (
587  name,
588  mesh.time().constant(),
590  mesh,
593  false
594  ),
595  surface
596  ),
597  sampleSource_(samplingSourceNames_[sampleSourceName]),
598  needsUpdate_(true),
599  sampleElements_(0),
600  samplePoints_(0)
601 {}
602 
603 
604 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
605 
607 {}
608 
609 
610 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
611 
613 {
614  return needsUpdate_;
615 }
616 
617 
619 {
620  // already marked as expired
621  if (needsUpdate_)
622  {
623  return false;
624  }
625 
628 
629  boundaryTreePtr_.clear();
630  sampleElements_.clear();
631  samplePoints_.clear();
632 
633  needsUpdate_ = true;
634  return true;
635 }
636 
637 
639 {
640  if (!needsUpdate_)
641  {
642  return false;
643  }
644 
645  // Mesh search engine, no triangulation of faces.
646  meshSearch meshSearcher(mesh(), polyMesh::FACE_PLANES);
647 
648  return update(meshSearcher);
649 }
650 
651 
653 {
654  if (!needsUpdate_)
655  {
656  return false;
657  }
658 
659  // Mesh search engine on subset, no triangulation of faces.
660  meshSearch meshSearcher(mesh(), bb, polyMesh::FACE_PLANES);
661 
662  return update(meshSearcher);
663 }
664 
665 
668 (
669  const volScalarField& vField
670 ) const
671 {
672  return sampleField(vField);
673 }
674 
675 
678 (
679  const volVectorField& vField
680 ) const
681 {
682  return sampleField(vField);
683 }
684 
687 (
688  const volSphericalTensorField& vField
689 ) const
690 {
691  return sampleField(vField);
692 }
693 
694 
697 (
698  const volSymmTensorField& vField
699 ) const
700 {
701  return sampleField(vField);
702 }
703 
704 
707 (
708  const volTensorField& vField
709 ) const
710 {
711  return sampleField(vField);
712 }
713 
714 
717 (
718  const interpolation<scalar>& interpolator
719 ) const
720 {
721  return interpolateField(interpolator);
722 }
723 
724 
727 (
728  const interpolation<vector>& interpolator
729 ) const
730 {
731  return interpolateField(interpolator);
732 }
733 
736 (
737  const interpolation<sphericalTensor>& interpolator
738 ) const
739 {
740  return interpolateField(interpolator);
741 }
742 
743 
746 (
747  const interpolation<symmTensor>& interpolator
748 ) const
749 {
750  return interpolateField(interpolator);
751 }
752 
753 
756 (
757  const interpolation<tensor>& interpolator
758 ) const
759 {
760  return interpolateField(interpolator);
761 }
762 
763 
765 {
766  os << "triSurfaceMesh: " << name() << " :"
767  << " surface:" << surface_.objectRegistry::name()
768  << " faces:" << faces().size()
769  << " points:" << points().size();
770 }
771 
772 
773 // ************************************************************************* //
const fvPatchList & patches
polyMesh::cellDecomposition decompMode() const
Definition: meshSearch.H:187
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search...
Definition: meshSearch.H:57
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
tUEqn clear()
const indexedOctree< treeDataCell > & cellTree() const
Get (demand driven) reference to octree holding all cells.
Definition: meshSearch.C:656
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
An abstract class for surfaces with sampling.
A 2-tuple for storing two objects of different types.
Definition: HashTable.H:65
Encapsulation of data needed to search for faces.
Definition: treeDataFace.H:58
defineTypeNameAndDebug(distanceSurface, 0)
triSurfaceMesh(const word &name, const polyMesh &mesh, const word &surfaceName, const samplingSource sampleSource)
Construct from components.
const Type1 & first() const
Return first.
Definition: Tuple2.H:116
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Output to file stream.
Definition: OFstream.H:82
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
bool interpolate() const
Interpolation requested for surface.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
virtual bool update()
Update the surface as required.
Tuple2< scalar, label > nearInfo
Private class for finding nearest.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
fvMesh & mesh
const fileName & name() const
Return the name of the stream.
Definition: OFstream.H:120
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
Macros for easy insertion into run-time selection tables.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
scalar y
const Point & hitPoint() const
Return hit point.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
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
virtual void print(Ostream &) const
Write.
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
Definition: polyPatch.H:316
virtual void clearGeom() const
const cellShapeList & cells
const pointField & points
static const word & geometryDir()
Return the geometry directory name.
A class for handling words, derived from string.
Definition: word.H:59
A triangular face using a FixedList of labels corresponding to mesh vertices.
Definition: triFace.H:68
const word & constant() const
Return constant name.
Definition: TimePaths.H:123
bool hit() const
Is there a hit.
const Field< PointType > & points() const
Return reference to global points.
const Type & shapes() const
Reference to shape.
static const label labelMax
Definition: label.H:62
const treeBoundBox & bb() const
Top bounding box.
Foam::polyBoundaryMesh.
virtual bool expire()
Mark the surface as needing an update.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
const Point & rawPoint() const
Return point with no checking.
Definition: PointHit.H:158
static const char nl
Definition: Ostream.H:260
label findInside(const point &) const
Find shape containing point. Only implemented for certain.
const Time & time() const
Return time.
labelList f(nPoints)
addToRunTimeSelectionTable(sampledSurface, distanceSurface, word)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
labelList second(const UList< labelPair > &p)
Definition: patchToPatch.C:48
void setSize(const label)
Reset size of List.
Definition: List.C:281
label patchi
std::remove_reference< ::Foam::List< labelledTri > >::type::value_type FaceType
pointHit nearestPoint(const point &p, const pointField &) const
Return nearest point to face.
label newPointi
Definition: readKivaGrid.H:501
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:57
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
A List with indirect addressing.
Definition: fvMatrix.H:106
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:306
virtual bool needsUpdate() const
Does the surface need an update?
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
void operator()(nearInfo &x, const nearInfo &y) const
messageStream Info
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
A class for managing temporary objects.
Definition: PtrList.H:53
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Triangulated surface description with patch information.
Definition: triSurface.H:66
A topoSetSource to select faces based on use of points.
Definition: pointToFace.H:49
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
label index() const
Return index.
virtual tmp< scalarField > sample(const volScalarField &) const
Sample field on surface.
treeBoundBox extend(const scalar s) const
Return asymetrically extended bounding box, with guaranteed.
face triFaceFace() const
Return triangle as a face.
Definition: triFaceI.H:140
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864