Go to the documentation of this file.
53 void Foam::sampledIsoSurface::getIsoFields()
const
55 const fvMesh& fvm =
static_cast<const fvMesh&
>(
mesh());
65 <<
"Lookup volField " << isoField_ <<
endl;
67 storedVolFieldPtr_.clear();
74 <<
"Checking " << isoField_
75 <<
" for same time " << fvm.time().timeName() <<
endl;
80 || (fvm.time().timeName() != storedVolFieldPtr_().instance())
84 <<
"Reading volField " << isoField_
85 <<
" from time " << fvm.time().timeName() <<
endl;
90 fvm.time().timeName(),
99 storedVolFieldPtr_.reset
107 volFieldPtr_ = storedVolFieldPtr_.get();
112 <<
"Cannot find isosurface field " << isoField_
113 <<
" in database or directory " << vfHeader.path()
133 const word pointFldName =
134 "volPointInterpolate_"
146 <<
"lookup pointField " << pointFldName <<
endl;
148 if (!pointFieldPtr_->upToDate(*volFieldPtr_))
151 <<
"updating pointField " << pointFldName <<
endl;
166 <<
"creating pointField " << pointFldName <<
endl;
187 storedVolFieldPtr_.reset
191 volFieldPtr_ = storedVolFieldPtr_.get();
196 <<
"volField " << volFieldPtr_->name()
197 <<
" min:" <<
min(*volFieldPtr_).value()
198 <<
" max:" <<
max(*volFieldPtr_).value() <<
nl
199 <<
"pointField " << pointFieldPtr_->name()
206 const fvMesh& subFvm = subMeshPtr_().subMesh();
215 <<
"Sub-mesh lookup volField " << isoField_ <<
endl;
217 storedVolSubFieldPtr_.clear();
222 <<
"Sub-setting volField " << isoField_ <<
endl;
224 storedVolSubFieldPtr_.reset
231 storedVolSubFieldPtr_->checkOut();
232 volSubFieldPtr_ = storedVolSubFieldPtr_.get();
238 const word pointFldName =
239 "volPointInterpolate_"
242 + volSubFieldPtr_->name()
251 <<
"Sub-mesh lookup pointField " << pointFldName <<
endl;
253 if (!pointFieldPtr_->upToDate(*volSubFieldPtr_))
256 <<
"Updating submesh pointField " << pointFldName <<
endl;
269 <<
"Interpolating submesh volField "
270 << volSubFieldPtr_->name()
271 <<
" to get submesh pointField " << pointFldName <<
endl;
287 storedVolSubFieldPtr_.reset
291 volSubFieldPtr_ = storedVolSubFieldPtr_.get();
297 << volSubFieldPtr_->name()
298 <<
" min:" <<
min(*volSubFieldPtr_).value()
299 <<
" max:" <<
max(*volSubFieldPtr_).value() <<
nl
301 << pointSubFieldPtr_->name()
308 void Foam::sampledIsoSurface::combineSurfaces
310 PtrList<isoSurfaceBase>& isoSurfPtrs
313 isoSurfacePtr_.reset(
nullptr);
321 && isoSurfPtrs.size() == 1
325 isoSurfacePtr_.reset(isoSurfPtrs.release(0));
327 else if (isoSurfPtrs.size() == 1)
329 autoPtr<isoSurfaceBase> ptr(isoSurfPtrs.release(0));
333 meshCells_.transfer(surf.meshCells());
344 for (
const auto& surf : isoSurfPtrs)
346 nFaces += surf.size();
347 nPoints += surf.points().size();
352 meshCells_.resize(nFaces);
358 forAll(isoSurfPtrs, surfi)
360 autoPtr<isoSurfaceBase> ptr(isoSurfPtrs.release(surfi));
363 SubList<face> subFaces(newFaces, surf.size(), nFaces);
364 SubList<point> subPoints(newPoints, surf.points().size(),
nPoints);
365 SubList<label> subCells(meshCells_, surf.size(), nFaces);
367 newZones[surfi] = surfZone
375 subFaces = surf.surfFaces();
376 subPoints = surf.points();
377 subCells = surf.meshCells();
381 for (face&
f : subFaces)
383 for (label& pointi :
f)
390 nFaces += subFaces.size();
396 std::move(newPoints),
401 surface_.transfer(combined);
405 if (subMeshPtr_ && meshCells_.size())
408 UIndirectList<label>(subMeshPtr_->cellMap(), meshCells_);
413 bool Foam::sampledIsoSurface::updateGeometry()
const
415 const fvMesh& fvm =
static_cast<const fvMesh&
>(
mesh());
418 if (fvm.time().timeIndex() == prevTimeIndex_)
423 prevTimeIndex_ = fvm.time().timeIndex();
428 isoSurfacePtr_.reset(
nullptr);
433 const bool hasCellZones =
443 subMeshPtr_.reset(
nullptr);
446 if (!ignoreCellsPtr_)
448 ignoreCellsPtr_.reset(
new bitSet);
452 bitSet select(
mesh().cellZones().selection(zoneNames_));
454 if (select.any() && !select.all())
459 *ignoreCellsPtr_ = std::move(select);
470 ignoreCellsPtr_->clearStorage();
474 ignoreCellsPtr_.reset(
new bitSet);
478 if (!subMeshPtr_ && hasCellZones)
480 const label exposedPatchi =
484 <<
"Allocating subset of size "
486 <<
" with exposed faces into patch "
487 << exposedPatchi <<
endl;
494 mesh().cellZones().selection(zoneNames_),
505 refPtr<volScalarField> tvolFld(*volFieldPtr_);
506 refPtr<pointScalarField> tpointFld(*pointFieldPtr_);
510 tvolFld.cref(*volSubFieldPtr_);
511 tpointFld.cref(*pointSubFieldPtr_);
517 PtrList<isoSurfaceBase> isoSurfPtrs(isoValues_.size());
528 tpointFld().primitiveField(),
536 const_cast<sampledIsoSurface&
>(*this)
537 .combineSurfaces(isoSurfPtrs);
553 meshCells_ = UIndirectList<label>(meshCells_,
faceMap)();
558 Pout<<
"isoSurface::updateGeometry() : constructed iso:" <<
nl
559 <<
" field : " << isoField_ <<
nl
561 <<
" average : " << Switch(average_) <<
nl
563 << Switch(
bool(isoParams_.filter())) <<
nl
564 <<
" bounds : " << isoParams_.getClipBounds() <<
nl;
567 Pout<<
" zone size : "
568 << subMeshPtr_->subMesh().nCells() <<
nl;
571 <<
" faces : " << surface().size() <<
nl
572 <<
" cut cells : " << meshCells().size()
593 isoParams_(
dict, params),
603 isoSurfacePtr_(
nullptr),
605 subMeshPtr_(
nullptr),
606 ignoreCellsPtr_(
nullptr),
608 storedVolFieldPtr_(
nullptr),
609 volFieldPtr_(
nullptr),
610 pointFieldPtr_(
nullptr),
612 storedVolSubFieldPtr_(
nullptr),
613 volSubFieldPtr_(
nullptr),
614 pointSubFieldPtr_(
nullptr)
619 isoParams_.algorithm(params.
algorithm());
626 isoValues_.resize(1);
630 if (isoValues_.empty())
633 <<
"No isoValue or isoValues specified." <<
nl
637 if (isoValues_.size() > 1)
639 const label nOrig = isoValues_.size();
643 if (nOrig != isoValues_.size())
646 <<
"Removed non-unique isoValues" <<
nl;
653 simpleSubMesh_ =
false;
659 if (isoValues_.size() > 1)
662 <<
"Multiple values on iso-surface (point) not supported"
663 <<
" since needs original interpolators." <<
nl
677 <<
"Cannot triangulate without a regularise filter" <<
nl
687 zoneNames_.resize(1);
696 <<
"Restricting to cellZone(s) " <<
flatOutput(zoneNames_)
697 <<
" with exposed internal faces into patch "
734 isoSurfacePtr_.reset(
nullptr);
735 subMeshPtr_.reset(
nullptr);
741 if (prevTimeIndex_ == -1)
754 return updateGeometry();
764 return sampleOnFaces(sampler);
774 return sampleOnFaces(sampler);
784 return sampleOnFaces(sampler);
794 return sampleOnFaces(sampler);
804 return sampleOnFaces(sampler);
814 return sampleOnPoints(interpolator);
824 return sampleOnPoints(interpolator);
833 return sampleOnPoints(interpolator);
843 return sampleOnPoints(interpolator);
853 return sampleOnPoints(interpolator);
859 os <<
"isoSurface: " <<
name() <<
" :";
860 isoParams_.print(
os);
861 os <<
" field:" << isoField_
int debug
Static debugging option.
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
static autoPtr< isoSurfaceBase > New(const isoSurfaceParams ¶ms, const volScalarField &cellValues, const scalarField &pointValues, const scalar iso, const bitSet &ignoreCells=bitSet())
Create for specified algorithm type.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
A class for handling words, derived from Foam::string.
label findIndex(const wordRe &key) const
Zone index for the first match, return -1 if not found.
virtual bool needsUpdate() const
Does the surface need an update?
A class for managing temporary objects.
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
algorithmType algorithm() const noexcept
Get current algorithm.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
static const volPointInterpolation & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
void inplaceUniqueSort(ListType &input)
Inplace sorting and removal of duplicates.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
virtual bool update()
Update the surface as required.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
tmp< GeometricField< Type, pointPatchField, pointMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &) const
Interpolate volField using inverse distance weighting.
static tmp< GeometricField< Type, fvPatchField, volMesh > > pointAverage(const GeometricField< Type, pointPatchField, pointMesh > &pfld)
Create cell values by averaging the point values.
unsigned int count(const bool on=true) const
Count number of bits set.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
label timeIndex() const noexcept
Return current time index.
bool interpolate() const noexcept
Same as isPointData()
virtual tmp< scalarField > sample(const interpolation< scalar > &sampler) const
Sample volume field onto surface faces.
#define DebugInFunction
Report an information message using Foam::Info.
Use current 'standard' algorithm.
sampledIsoSurface(const isoSurfaceParams ¶ms, const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
A sampledSurface defined by a surface of iso value. It only recalculates the iso-surface if time chan...
addNamedToRunTimeSelectionTable(topoSetCellSource, badQualityToCell, word, badQuality)
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
An abstract class for surfaces with sampling.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual bool expire()
Mark the surface as needing an update.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Preferences for controlling iso-surface algorithms.
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
Mesh data needed to do the Finite Volume discretisation.
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
errorManipArg< error, int > exit(error &err, const int errNo=1)
static word defaultName(const label n=-1)
Default zone name: "zone" or "zoneN".
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
MeshedSurface< face > meshedSurface
List< face > faceList
A List of faces.
virtual void print(Ostream &os, int level=0) const
Print information.
const polyMesh & mesh() const noexcept
Access to the underlying mesh.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
virtual void clearGeom() const
Additional cleanup when clearing the geometry.
List< surfZone > surfZoneList
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const Time & time() const
Return the top-level database.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
bitSet selection(const labelUList &zoneIds) const
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Type gMin(const FieldField< Field, Type > &f)
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
virtual ~sampledIsoSurface()
Destructor.
defineTypeNameAndDebug(combustionModel, 0)
Type gMax(const FieldField< Field, Type > &f)
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField