Go to the documentation of this file.
51 void Foam::sampledIsoSurface::getIsoFields()
const
53 const fvMesh& fvm =
static_cast<const fvMesh&
>(
mesh());
63 <<
"Lookup volField " << isoField_ <<
endl;
65 storedVolFieldPtr_.clear();
72 <<
"Checking " << isoField_
73 <<
" for same time " << fvm.time().timeName() <<
endl;
77 storedVolFieldPtr_.empty()
78 || (fvm.time().timeName() != storedVolFieldPtr_().instance())
82 <<
"Reading volField " << isoField_
83 <<
" from time " << fvm.time().timeName() <<
endl;
88 fvm.time().timeName(),
97 storedVolFieldPtr_.reset
105 volFieldPtr_ = storedVolFieldPtr_.get();
110 <<
"Cannot find isosurface field " << isoField_
111 <<
" in database or directory " << vfHeader.path()
129 if (subMeshPtr_.empty())
131 const word pointFldName =
132 "volPointInterpolate_"
144 <<
"lookup pointField " << pointFldName <<
endl;
146 if (!pointFieldPtr_->upToDate(*volFieldPtr_))
149 <<
"updating pointField " << pointFldName <<
endl;
164 <<
"creating pointField " << pointFldName <<
endl;
185 storedVolFieldPtr_.reset
189 volFieldPtr_ = storedVolFieldPtr_.get();
194 <<
"volField " << volFieldPtr_->name()
195 <<
" min:" <<
min(*volFieldPtr_).value()
196 <<
" max:" <<
max(*volFieldPtr_).value() <<
nl
197 <<
"pointField " << pointFieldPtr_->name()
204 const fvMesh& subFvm = subMeshPtr_().subMesh();
213 <<
"Sub-mesh lookup volField " << isoField_ <<
endl;
215 storedVolSubFieldPtr_.clear();
220 <<
"Sub-setting volField " << isoField_ <<
endl;
222 storedVolSubFieldPtr_.reset
229 storedVolSubFieldPtr_->checkOut();
230 volSubFieldPtr_ = storedVolSubFieldPtr_.get();
236 const word pointFldName =
237 "volPointInterpolate_"
240 + volSubFieldPtr_->name()
249 <<
"Sub-mesh lookup pointField " << pointFldName <<
endl;
251 if (!pointFieldPtr_->upToDate(*volSubFieldPtr_))
254 <<
"Updating submesh pointField " << pointFldName <<
endl;
267 <<
"Interpolating submesh volField "
268 << volSubFieldPtr_->name()
269 <<
" to get submesh pointField " << pointFldName <<
endl;
285 storedVolSubFieldPtr_.reset
289 volSubFieldPtr_ = storedVolSubFieldPtr_.get();
295 << volSubFieldPtr_->name()
296 <<
" min:" <<
min(*volSubFieldPtr_).value()
297 <<
" max:" <<
max(*volSubFieldPtr_).value() <<
nl
299 << pointSubFieldPtr_->name()
306 bool Foam::sampledIsoSurface::updateGeometry()
const
308 const fvMesh& fvm =
static_cast<const fvMesh&
>(
mesh());
311 if (fvm.time().timeIndex() == prevTimeIndex_)
320 && subMeshPtr_.empty()
329 <<
"Allocating subset of size "
331 <<
" with exposed faces into patch "
339 mesh().cellZones().selection(zoneNames_),
346 prevTimeIndex_ = fvm.time().timeIndex();
355 if (subMeshPtr_.valid())
393 Pout<<
"sampledIsoSurface::updateGeometry() : constructed iso:"
395 <<
" filter : " << Switch(
bool(filter_)) <<
nl
396 <<
" average : " << Switch(average_) <<
nl
397 <<
" isoField : " << isoField_ <<
nl
398 <<
" isoValue : " << isoVal_ <<
nl;
399 if (subMeshPtr_.valid())
401 Pout<<
" zone size : " << subMeshPtr_().subMesh().nCells()
405 <<
" faces : " << surface().size() <<
nl
406 <<
" cut cells : " << surface().meshCells().size()
425 isoVal_(
dict.
get<scalar>(
"isoValue")),
441 storedVolFieldPtr_(
nullptr),
442 volFieldPtr_(
nullptr),
443 pointFieldPtr_(
nullptr)
448 <<
"Non-interpolated iso surface not supported since triangles"
455 zoneNames_.resize(1);
466 <<
"Cannot find patch " << exposedPatchName_
467 <<
" in which to put exposed faces." <<
endl
473 <<
"Restricting to cellZone(s) " <<
flatOutput(zoneNames_)
474 <<
" with exposed internal faces into patch "
475 << exposedPatchName_ <<
endl;
505 if (prevTimeIndex_ == -1)
518 return updateGeometry();
527 return sampleOnFaces(sampler);
536 return sampleOnFaces(sampler);
545 return sampleOnFaces(sampler);
554 return sampleOnFaces(sampler);
563 return sampleOnFaces(sampler);
572 return sampleOnPoints(interpolator);
581 return sampleOnPoints(interpolator);
589 return sampleOnPoints(interpolator);
598 return sampleOnPoints(interpolator);
607 return sampleOnPoints(interpolator);
613 os <<
"sampledIsoSurface: " <<
name() <<
" :"
614 <<
" field :" << isoField_
615 <<
" value :" << isoVal_;
int debug
Static debugging option.
const word & name() const
Return name.
Remove points from pyramid edges and face-diagonals.
A class for handling words, derived from Foam::string.
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.
static const boundBox invertedBox
A large inverted boundBox: min/max == +/- ROOTVGREAT.
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.
const cellZoneMesh & cellZones() const
Return cell zone mesh.
label findIndex(const keyType &key) const
Return zone index for the first match, return -1 if not found.
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
An Ostream wrapper for parallel output to std::cout.
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.
virtual void print(Ostream &) const
Write.
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.
wordList names() const
Return a list of patch names.
virtual tmp< scalarField > sample(const interpolation< scalar > &sampler) const
Sample volume field onto surface faces.
#define DebugInFunction
Report an information message using Foam::Info.
word name(const complex &c)
Return string representation of complex.
label findIndex(const ListType &input, typename ListType::const_reference val, const label start=0)
Deprecated(2017-10) search for first occurrence of the given element.
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,...
Macros for easy insertion into run-time selection tables.
Mesh data needed to do the Finite Volume discretisation.
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
FlatOutput< Container > flatOutput(const Container &obj, label len=0)
Global flatOutput function.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const dimensionedScalar e
Elementary charge.
static filterType getFilterType(const dictionary &dict, const isoSurfaceBase::filterType deflt)
Get 'regularise' as bool or enumeration.
const polyBoundaryMesh & patches
const Time & time() const
Return the top-level database.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
const polyMesh & mesh() const
Access to the underlying mesh.
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)
label timeIndex() const
Return current time index.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
bool interpolate() const
Interpolation to nodes requested for surface.
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
sampledIsoSurface(const word &name, const polyMesh &mesh, const dictionary &dict)
Construct from dictionary.