sampledCuttingPlane.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2016-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "sampledCuttingPlane.H"
30 #include "dictionary.H"
31 #include "fvMesh.H"
32 #include "volFields.H"
33 #include "volPointInterpolation.H"
35 #include "PtrList.H"
36 
37 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41  defineTypeNameAndDebug(sampledCuttingPlane, 0);
43  (
44  sampledSurface,
45  sampledCuttingPlane,
46  word,
47  cuttingPlane
48  );
49 }
50 
51 
52 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
53 
54 void Foam::sampledCuttingPlane::checkBoundsIntersection
55 (
56  const plane& pln,
57  const boundBox& meshBb
58 ) const
59 {
60  // Verify specified bounding box
61  const boundBox& clipBb = isoParams_.getClipBounds();
62 
63  if (clipBb.valid())
64  {
65  // Bounding box does not overlap with (global) mesh!
66  if (!clipBb.overlaps(meshBb))
67  {
69  << nl
70  << name() << " : "
71  << "Bounds " << clipBb
72  << " do not overlap the mesh bounding box " << meshBb
73  << nl << endl;
74  }
75 
76  // Plane does not intersect the bounding box
77  if (!clipBb.intersects(pln))
78  {
80  << nl
81  << name() << " : "
82  << "Plane "<< pln << " does not intersect the bounds "
83  << clipBb
84  << nl << endl;
85  }
86  }
87 
88  // Plane does not intersect the (global) mesh!
89  if (!meshBb.intersects(pln))
90  {
92  << nl
93  << name() << " : "
94  << "Plane "<< pln << " does not intersect the mesh bounds "
95  << meshBb
96  << nl << endl;
97  }
98 }
99 
100 
101 void Foam::sampledCuttingPlane::setDistanceFields(const plane& pln)
102 {
103  volScalarField& cellDistance = cellDistancePtr_();
104 
105  // Get mesh from volField,
106  // so automatically submesh or baseMesh
107 
108  const fvMesh& mesh = cellDistance.mesh();
109 
110  // Distance to cell centres
111  // ~~~~~~~~~~~~~~~~~~~~~~~~
112 
113  // Internal field
114  {
115  const auto& cc = mesh.cellCentres();
116  scalarField& fld = cellDistance.primitiveFieldRef();
117 
118  forAll(cc, i)
119  {
120  fld[i] = pln.signedDistance(cc[i]);
121  }
122  }
123 
124  // Patch fields
125  {
126  volScalarField::Boundary& cellDistanceBf =
127  cellDistance.boundaryFieldRef();
128 
129  forAll(cellDistanceBf, patchi)
130  {
131  if
132  (
133  isA<emptyFvPatchScalarField>
134  (
135  cellDistanceBf[patchi]
136  )
137  )
138  {
139  cellDistanceBf.set
140  (
141  patchi,
142  new calculatedFvPatchScalarField
143  (
144  mesh.boundary()[patchi],
145  cellDistance
146  )
147  );
148 
149  const polyPatch& pp = mesh.boundary()[patchi].patch();
150  pointField::subField cc = pp.patchSlice(mesh.faceCentres());
151 
152  fvPatchScalarField& fld = cellDistanceBf[patchi];
153  fld.setSize(pp.size());
154  forAll(fld, i)
155  {
156  fld[i] = pln.signedDistance(cc[i]);
157  }
158  }
159  else
160  {
161  // Other side cell centres?
162  const pointField& cc = mesh.C().boundaryField()[patchi];
163  fvPatchScalarField& fld = cellDistanceBf[patchi];
164 
165  forAll(fld, i)
166  {
167  fld[i] = pln.signedDistance(cc[i]);
168  }
169  }
170  }
171  }
172 
173 
174  // On processor patches the mesh.C() will already be the cell centre
175  // on the opposite side so no need to swap cellDistance.
176 
177  // Distance to points
178  pointDistance_.resize(mesh.nPoints());
179  {
180  const pointField& pts = mesh.points();
181 
182  forAll(pointDistance_, i)
183  {
184  pointDistance_[i] = pln.signedDistance(pts[i]);
185  }
186  }
187 }
188 
189 
190 void Foam::sampledCuttingPlane::combineSurfaces
191 (
192  PtrList<isoSurfaceBase>& isoSurfPtrs
193 )
194 {
195  isoSurfacePtr_.reset(nullptr);
196 
197  // Already checked previously for ALGO_POINT, but do it again
198  // - ALGO_POINT still needs fields (for interpolate)
199  // The others can do straight transfer
200  if
201  (
202  isoParams_.algorithm() == isoSurfaceParams::ALGO_POINT
203  && isoSurfPtrs.size() == 1
204  )
205  {
206  // Shift ownership from list to autoPtr
207  isoSurfacePtr_.reset(isoSurfPtrs.release(0));
208  }
209  else if (isoSurfPtrs.size() == 1)
210  {
211  autoPtr<isoSurfaceBase> ptr(isoSurfPtrs.release(0));
212  auto& surf = *ptr;
213 
214  surface_.transfer(static_cast<meshedSurface&>(surf));
215  meshCells_.transfer(surf.meshCells());
216  }
217  else
218  {
219  // Combine faces with point offsets
220  //
221  // Note: use points().size() from surface, not nPoints()
222  // since there may be uncompacted dangling nodes
223 
224  label nFaces = 0, nPoints = 0;
225 
226  for (const auto& surf : isoSurfPtrs)
227  {
228  nFaces += surf.size();
229  nPoints += surf.points().size();
230  }
231 
232  faceList newFaces(nFaces);
233  pointField newPoints(nPoints);
234  meshCells_.resize(nFaces);
235 
236  surfZoneList newZones(isoSurfPtrs.size());
237 
238  nFaces = 0;
239  nPoints = 0;
240  forAll(isoSurfPtrs, surfi)
241  {
242  autoPtr<isoSurfaceBase> ptr(isoSurfPtrs.release(surfi));
243  auto& surf = *ptr;
244 
245  SubList<face> subFaces(newFaces, surf.size(), nFaces);
246  SubList<point> subPoints(newPoints, surf.points().size(), nPoints);
247  SubList<label> subCells(meshCells_, surf.size(), nFaces);
248 
249  newZones[surfi] = surfZone
250  (
252  subFaces.size(), // size
253  nFaces, // start
254  surfi // index
255  );
256 
257  subFaces = surf.surfFaces();
258  subPoints = surf.points();
259  subCells = surf.meshCells();
260 
261  if (nPoints)
262  {
263  for (face& f : subFaces)
264  {
265  for (label& pointi : f)
266  {
267  pointi += nPoints;
268  }
269  }
270  }
271 
272  nFaces += subFaces.size();
273  nPoints += subPoints.size();
274  }
275 
276  meshedSurface combined
277  (
278  std::move(newPoints),
279  std::move(newFaces),
280  newZones
281  );
282 
283  surface_.transfer(combined);
284  }
285 
286  // Addressing into the full mesh
287  if (subMeshPtr_ && meshCells_.size())
288  {
289  meshCells_ =
290  UIndirectList<label>(subMeshPtr_->cellMap(), meshCells_);
291  }
292 }
293 
294 
295 void Foam::sampledCuttingPlane::createGeometry()
296 {
297  if (debug)
298  {
299  Pout<< "sampledCuttingPlane::createGeometry :updating geometry."
300  << endl;
301  }
302 
303  // Clear any previously stored topologies
304  surface_.clear();
305  meshCells_.clear();
306  isoSurfacePtr_.reset(nullptr);
307 
308  // Clear derived data
310 
311  // Clear any stored fields
312  pointDistance_.clear();
313  cellDistancePtr_.clear();
314 
315  const bool hasCellZones =
316  (-1 != mesh().cellZones().findIndex(zoneNames_));
317 
318  const fvMesh& fvm = static_cast<const fvMesh&>(this->mesh());
319 
320  // Geometry
321  if
322  (
323  simpleSubMesh_
324  && isoParams_.algorithm() != isoSurfaceParams::ALGO_POINT
325  )
326  {
327  subMeshPtr_.reset(nullptr);
328 
329  // Handle cell zones as inverse (blocked) selection
330  if (!ignoreCellsPtr_)
331  {
332  ignoreCellsPtr_.reset(new bitSet);
333 
334  if (hasCellZones)
335  {
336  bitSet select(mesh().cellZones().selection(zoneNames_));
337 
338  if (select.any() && !select.all())
339  {
340  // From selection to blocking
341  select.flip();
342 
343  *ignoreCellsPtr_ = std::move(select);
344  }
345  }
346  }
347  }
348  else
349  {
350  // A standard subMesh treatment
351 
352  if (ignoreCellsPtr_)
353  {
354  ignoreCellsPtr_->clearStorage();
355  }
356  else
357  {
358  ignoreCellsPtr_.reset(new bitSet);
359  }
360 
361  // Get sub-mesh if any
362  if (!subMeshPtr_ && hasCellZones)
363  {
364  const label exposedPatchi =
365  mesh().boundaryMesh().findPatchID(exposedPatchName_);
366 
367  bitSet cellsToSelect(mesh().cellZones().selection(zoneNames_));
368 
369  DebugInfo
370  << "Allocating subset of size "
371  << cellsToSelect.count() << " with exposed faces into patch "
372  << exposedPatchi << endl;
373 
374 
375  // If we will use a fvMeshSubset so can apply bounds as well to make
376  // the initial selection smaller.
377 
378  const boundBox& clipBb = isoParams_.getClipBounds();
379  if (clipBb.valid() && cellsToSelect.any())
380  {
381  const auto& cellCentres = fvm.C();
382 
383  for (const label celli : cellsToSelect)
384  {
385  const point& cc = cellCentres[celli];
386 
387  if (!clipBb.contains(cc))
388  {
389  cellsToSelect.unset(celli);
390  }
391  }
392 
393  DebugInfo
394  << "Bounded subset of size "
395  << cellsToSelect.count() << endl;
396  }
397 
398  subMeshPtr_.reset
399  (
400  new fvMeshSubset(fvm, cellsToSelect, exposedPatchi)
401  );
402  }
403  }
404 
405 
406  // Select either the submesh or the underlying mesh
407  const fvMesh& mesh =
408  (
409  subMeshPtr_
410  ? subMeshPtr_->subMesh()
411  : fvm
412  );
413 
414  checkBoundsIntersection(plane_, mesh.bounds());
415 
416 
417  // Distance to cell centres
418  // ~~~~~~~~~~~~~~~~~~~~~~~~
419 
420  cellDistancePtr_.reset
421  (
422  new volScalarField
423  (
424  IOobject
425  (
426  "cellDistance",
427  mesh.time().timeName(),
428  mesh.time(),
431  false
432  ),
433  mesh,
435  )
436  );
437  const volScalarField& cellDistance = cellDistancePtr_();
438 
439  setDistanceFields(plane_);
440 
441  if (debug)
442  {
443  Pout<< "Writing cell distance:" << cellDistance.objectPath() << endl;
444  cellDistance.write();
445  pointScalarField pointDist
446  (
447  IOobject
448  (
449  "pointDistance",
450  mesh.time().timeName(),
451  mesh.time(),
454  false
455  ),
458  );
459  pointDist.primitiveFieldRef() = pointDistance_;
460 
461  Pout<< "Writing point distance:" << pointDist.objectPath() << endl;
462  pointDist.write();
463  }
464 
465 
466  // Create surfaces for each offset
467 
468  PtrList<isoSurfaceBase> isoSurfPtrs(offsets_.size());
469 
470  forAll(offsets_, surfi)
471  {
472  isoSurfPtrs.set
473  (
474  surfi,
476  (
477  isoParams_,
478  cellDistance,
479  pointDistance_,
480  offsets_[surfi],
481  *ignoreCellsPtr_
482  )
483  );
484  }
485 
486  // And flatten
487  combineSurfaces(isoSurfPtrs);
488 
489 
490  // Discard fields if not required by an iso-surface
491  if (!isoSurfacePtr_)
492  {
493  cellDistancePtr_.reset(nullptr);
494  pointDistance_.clear();
495  }
496 
497  if (debug)
498  {
499  print(Pout, debug);
500  Pout<< endl;
501  }
502 }
503 
504 
505 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
506 
508 (
509  const word& name,
510  const polyMesh& mesh,
511  const dictionary& dict
512 )
513 :
515  plane_(dict),
516  offsets_(),
517  isoParams_
518  (
519  dict,
522  ),
523  average_(dict.getOrDefault("average", false)),
524  simpleSubMesh_(dict.getOrDefault("simpleSubMesh", false)),
525  zoneNames_(),
526  exposedPatchName_(),
527  needsUpdate_(true),
528 
529  surface_(),
530  meshCells_(),
531  isoSurfacePtr_(nullptr),
532 
533  subMeshPtr_(nullptr),
534  ignoreCellsPtr_(nullptr),
535  cellDistancePtr_(nullptr),
536  pointDistance_()
537 {
538  dict.readIfPresent("offsets", offsets_);
539 
540  if (offsets_.empty())
541  {
542  offsets_.resize(1);
543  offsets_.first() = Zero;
544  }
545 
546  if (offsets_.size() > 1)
547  {
548  const label nOrig = offsets_.size();
549 
550  inplaceUniqueSort(offsets_);
551 
552  if (nOrig != offsets_.size())
553  {
555  << "Removed non-unique offsets" << nl;
556  }
557  }
558 
559  if (isoParams_.algorithm() == isoSurfaceParams::ALGO_POINT)
560  {
561  // Not possible for ALGO_POINT
562  simpleSubMesh_ = false;
563 
564  // Not possible for ALGO_POINT
565  if (offsets_.size() > 1)
566  {
568  << "Multiple offsets with iso-surface (point) not supported"
569  << " since needs original interpolators." << nl
570  << exit(FatalIOError);
571  }
572  }
573 
574 
575  // Zones
576 
577  if (!dict.readIfPresent("zones", zoneNames_) && dict.found("zone"))
578  {
579  zoneNames_.resize(1);
580  dict.readEntry("zone", zoneNames_.first());
581  }
582 
583  if (-1 != mesh.cellZones().findIndex(zoneNames_))
584  {
585  dict.readIfPresent("exposedPatchName", exposedPatchName_);
586 
587  DebugInfo
588  << "Restricting to cellZone(s) " << flatOutput(zoneNames_)
589  << " with exposed internal faces into patch "
590  << mesh.boundaryMesh().findPatchID(exposedPatchName_) << endl;
591  }
592 }
593 
594 
595 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
596 
598 {
599  return needsUpdate_;
600 }
601 
602 
604 {
605  if (debug)
606  {
607  Pout<< "sampledCuttingPlane::expire :"
608  << " needsUpdate:" << needsUpdate_ << endl;
609  }
610 
611  surface_.clear();
612  meshCells_.clear();
613  isoSurfacePtr_.reset(nullptr);
614 
615  // Clear derived data
617 
618  // Already marked as expired
619  if (needsUpdate_)
620  {
621  return false;
622  }
623 
624  needsUpdate_ = true;
625  return true;
626 }
627 
628 
630 {
631  if (debug)
632  {
633  Pout<< "sampledCuttingPlane::update :"
634  << " needsUpdate:" << needsUpdate_ << endl;
635  }
636 
637  if (!needsUpdate_)
638  {
639  return false;
640  }
641 
642  createGeometry();
643 
644  needsUpdate_ = false;
645  return true;
646 }
647 
648 
651 (
652  const interpolation<scalar>& sampler
653 ) const
654 {
655  return sampleOnFaces(sampler);
656 }
657 
658 
661 (
662  const interpolation<vector>& sampler
663 ) const
664 {
665  return sampleOnFaces(sampler);
666 }
667 
668 
671 (
672  const interpolation<sphericalTensor>& sampler
673 ) const
674 {
675  return sampleOnFaces(sampler);
676 }
677 
678 
681 (
682  const interpolation<symmTensor>& sampler
683 ) const
684 {
685  return sampleOnFaces(sampler);
686 }
687 
688 
691 (
692  const interpolation<tensor>& sampler
693 ) const
694 {
695  return sampleOnFaces(sampler);
696 }
697 
698 
701 (
702  const interpolation<scalar>& interpolator
703 ) const
704 {
705  return sampleOnPoints(interpolator);
706 }
707 
708 
711 (
712  const interpolation<vector>& interpolator
713 ) const
714 {
715  return sampleOnPoints(interpolator);
716 }
717 
718 
721 (
722  const interpolation<sphericalTensor>& interpolator
723 ) const
724 {
725  return sampleOnPoints(interpolator);
726 }
727 
728 
731 (
732  const interpolation<symmTensor>& interpolator
733 ) const
734 {
735  return sampleOnPoints(interpolator);
736 }
737 
738 
741 (
742  const interpolation<tensor>& interpolator
743 ) const
744 {
745  return sampleOnPoints(interpolator);
746 }
747 
748 
750 {
751  os << "sampledCuttingPlane: " << name() << " :"
752  << " plane:" << plane_
753  << " offsets:" << flatOutput(offsets_);
754 
755  if (level)
756  {
757  os << " faces:" << faces().size()
758  << " points:" << points().size();
759  }
760 }
761 
762 
763 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::fvPatchScalarField
fvPatchField< scalar > fvPatchScalarField
Definition: fvPatchFieldsFwd.H:40
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
volFields.H
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Foam::polyMesh::points
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1069
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
Foam::isoSurfaceParams::filterType::DIAGCELL
Remove pyramid edge points, face-diagonals.
Foam::isoSurfaceBase::New
static autoPtr< isoSurfaceBase > New(const isoSurfaceParams &params, const volScalarField &cellValues, const scalarField &pointValues, const scalar iso, const bitSet &ignoreCells=bitSet())
Create for specified algorithm type.
Definition: isoSurfaceBaseNew.C:37
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::ZoneMesh::findIndex
label findIndex(const wordRe &key) const
Zone index for the first match, return -1 if not found.
Definition: ZoneMesh.C:491
Foam::dimLength
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::dictionary::found
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Definition: dictionaryI.H:87
Foam::MeshObject< polyMesh, UpdateableMeshObject, pointMesh >::New
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
Definition: MeshObject.C:48
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::inplaceUniqueSort
void inplaceUniqueSort(ListType &input)
Inplace sorting and removal of duplicates.
Definition: ListOpsTemplates.C:468
Foam::FatalIOError
IOerror FatalIOError
Foam::polyMesh::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:444
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
meshBb
List< treeBoundBox > meshBb(1, treeBoundBox(boundBox(coarseMesh.points(), false)).extend(rndGen, 1e-3))
Foam::Pout
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::primitiveMesh::nPoints
label nPoints() const noexcept
Number of mesh points.
Definition: primitiveMeshI.H:37
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::sampledCuttingPlane::print
virtual void print(Ostream &os, int level=0) const
Print information.
Definition: sampledCuttingPlane.C:749
nPoints
label nPoints
Definition: gmvOutputHeader.H:2
sampledCuttingPlane.H
Foam::sampledCuttingPlane::needsUpdate
virtual bool needsUpdate() const
Does the surface need an update?
Definition: sampledCuttingPlane.C:597
Foam::sampledSurface::interpolate
bool interpolate() const noexcept
Same as isPointData()
Definition: sampledSurface.H:598
Foam::isoSurfaceParams::getClipBounds
const boundBox & getClipBounds() const noexcept
Get optional clipping bounding box.
Definition: isoSurfaceParams.H:261
Foam::polyMesh::cellZones
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
Definition: polyMesh.H:492
Foam::addNamedToRunTimeSelectionTable
addNamedToRunTimeSelectionTable(topoSetCellSource, badQualityToCell, word, badQuality)
Foam::dictionary::readEntry
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
Definition: dictionaryTemplates.C:302
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::sampledCuttingPlane::update
virtual bool update()
Update the surface as required.
Definition: sampledCuttingPlane.C:629
Foam::sampledSurface
An abstract class for surfaces with sampling.
Definition: sampledSurface.H:121
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Foam::polyBoundaryMesh::findPatchID
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
Definition: polyBoundaryMesh.C:765
fld
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputLagrangian.H:23
Foam::fvMesh::C
const volVectorField & C() const
Return cell centres as volVectorField.
Definition: fvMeshGeometry.C:341
Foam::interpolation< scalar >
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::sampledSurface::name
const word & name() const noexcept
Name of surface.
Definition: sampledSurface.H:322
Foam::isoSurfaceParams::ALGO_POINT
Definition: isoSurfaceParams.H:119
Foam::Ostream::write
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
fvMesh.H
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::flatOutput
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
Definition: FlatOutput.H:216
Foam::sampledCuttingPlane::sampledCuttingPlane
sampledCuttingPlane(const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.
Definition: sampledCuttingPlane.C:508
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::GeometricField::primitiveFieldRef
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
Definition: GeometricField.C:766
Foam::isoSurfaceParams::ALGO_TOPO
Definition: isoSurfaceParams.H:117
volPointInterpolation.H
Foam::polyMesh::bounds
const boundBox & bounds() const
Return mesh bounding box.
Definition: polyMesh.H:450
Foam::fvMesh::boundary
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:685
Foam::surfZoneIdentifier::defaultName
static word defaultName(const label n=-1)
Default zone name: "zone" or "zoneN".
Definition: surfZoneIdentifier.H:83
Foam::primitiveMesh::cellCentres
const vectorField & cellCentres() const
Definition: primitiveMeshCellCentresAndVols.C:84
DebugInfo
#define DebugInfo
Report an information message using Foam::Info.
Definition: messageStream.H:382
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::meshedSurface
MeshedSurface< face > meshedSurface
Definition: MeshedSurfacesFwd.H:41
f
labelList f(nPoints)
Foam::faceList
List< face > faceList
A List of faces.
Definition: faceListFwd.H:47
Foam::primitiveMesh::faceCentres
const vectorField & faceCentres() const
Definition: primitiveMeshFaceCentresAndAreas.C:77
Foam::BitOps::print
Ostream & print(Ostream &os, UIntType value, char off='0', char on='1')
Print 0/1 bits in the (unsigned) integral type.
Definition: BitOps.H:199
points
const pointField & points
Definition: gmvOutputHeader.H:1
dictionary.H
Foam::primitiveMesh::reset
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
Definition: primitiveMesh.C:207
Foam::sampledSurface::clearGeom
virtual void clearGeom() const
Additional cleanup when clearing the geometry.
Definition: sampledSurface.C:53
Foam::surfZoneList
List< surfZone > surfZoneList
Definition: surfZoneList.H:47
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
PtrList.H
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:280
Foam::Field::subField
SubField< Type > subField
Declare type of subField.
Definition: Field.H:89
FatalIOErrorInFunction
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
Foam::sampledCuttingPlane::sample
virtual tmp< scalarField > sample(const interpolation< scalar > &sampler) const
Sample volume field onto surface faces.
Definition: sampledCuttingPlane.C:651
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
IOWarningInFunction
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
Definition: messageStream.H:340
Foam::point
vector point
Point is a vector.
Definition: point.H:43
Foam::dictionary::getOrDefault
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:148
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::sampledCuttingPlane::expire
virtual bool expire()
Mark the surface as needing an update.
Definition: sampledCuttingPlane.C:603
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:328
Foam::dictionary::readIfPresent
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:405
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::pointScalarField
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
Definition: pointFieldsFwd.H:56