Go to the documentation of this file.
44 Foam::sampledMeshedSurface::samplingSourceNames_
47 { samplingSource::insideCells,
"insideCells" },
48 { samplingSource::boundaryFaces,
"boundaryFaces" },
84 if (
y.first() <
x.first())
119 void Foam::sampledMeshedSurface::setZoneMap()
125 const auto& zones =
s.surfZones();
130 if (zoneIds_.empty() || zones.size() <= 1)
141 const label len =
min(zones[zonei].
size(), zoneIds_.size() - beg);
144 SubList<label>(zoneIds_, len, beg) = zonei;
151 const label len = (zoneIds_.size() - beg);
155 SubList<label>(zoneIds_, len, beg) =
max(0, zones.size()-1);
167 globalIndex globalCells(onBoundary() ?
mesh().nFaces() :
mesh().nCells());
171 const pointField& fc = surface_.faceCentres();
175 if (sampleSource_ ==
cells)
179 const indexedOctree<treeDataCell>& cellTree = meshSearcher.cellTree();
183 const point& pt = fc[facei];
189 nearest[facei].first() =
magSqr(info.hitPoint()-pt);
190 nearest[facei].second() = globalCells.toGlobal(info.index());
194 else if (sampleSource_ == insideCells)
198 const auto& cellTree = meshSearcher.cellTree();
202 const point& pt = fc[facei];
204 if (cellTree.bb().contains(pt))
206 const label index = cellTree.findInside(pt);
209 nearest[facei].first() = 0;
210 nearest[facei].second() = globalCells.toGlobal(index);
220 const auto& bndTree = meshSearcher.nonCoupledBoundaryTree();
224 const point& pt = fc[facei];
230 nearest[facei].first() =
magSqr(info.hitPoint()-pt);
231 nearest[facei].second() =
234 bndTree.shapes().faceLabels()[info.index()]
247 labelList cellOrFaceLabels(fc.size(), -1);
249 bitSet facesToSubset(fc.size());
253 const label index = nearest[facei].second();
259 else if (globalCells.isLocal(index))
261 cellOrFaceLabels[facei] = globalCells.toLocal(index);
262 facesToSubset.set(facei);
269 Pout<<
"Local out of faces:" << cellOrFaceLabels.size()
270 <<
" keeping:" << facesToSubset.count() <<
endl;
280 s = surface_.subsetMesh(facesToSubset, pointMap,
faceMap);
298 samplePoints_.resize(pointMap.size());
299 sampleElements_.resize(pointMap.size(), -1);
302 labelList pointToFace(std::move(pointMap));
306 const face&
f =
s[facei];
308 for (
const label labi :
f)
310 pointToFace[labi] = facei;
315 if (sampleSource_ ==
cells)
324 const label celli = cellOrFaceLabels[pointToFace[pointi]];
326 sampleElements_[pointi] = celli;
334 sampleElements_[pointi],
335 meshSearcher.decompMode()
339 samplePoints_[pointi] = pt;
345 scalar minDistSqr = VGREAT;
347 for (
const label facei :
mesh().
cells()[celli])
353 if (info.distance() < minDistSqr)
355 minDistSqr = info.distance();
356 samplePoints_[pointi] = info.rawPoint();
362 else if (sampleSource_ == insideCells)
371 const label celli = cellOrFaceLabels[pointToFace[pointi]];
373 sampleElements_[pointi] = celli;
374 samplePoints_[pointi] = pt;
387 const label facei = cellOrFaceLabels[pointToFace[pointi]];
389 sampleElements_[pointi] = facei;
390 samplePoints_[pointi] =
mesh().
faces()[facei].nearestPoint
409 sampleElements_.
transfer(cellOrFaceLabels);
410 samplePoints_.clear();
417 OFstream str(
mesh().time().
path()/
"surfaceToMesh.obj");
418 Info<<
"Dumping correspondence from local surface (points or faces)"
419 <<
" to mesh (cells or faces) to " << str.name() <<
endl;
431 forAll(samplePoints_, pointi)
439 label elemi = sampleElements_[pointi];
442 str <<
"l " << vertI-2 <<
' ' << vertI-1 <<
' ' << vertI <<
nl;
448 forAll(sampleElements_, triI)
453 label elemi = sampleElements_[triI];
456 str <<
"l " << vertI-1 <<
' ' << vertI <<
nl;
461 needsUpdate_ =
false;
472 const word& surfaceName,
478 surfaceName_(surfaceName),
484 sampleSource_(sampleSource),
504 meshedSurface::findFile
515 sampleSource_(samplingSourceNames_.get(
"source",
dict)),
524 includePatches.
uniq();
529 if (!includePatches.empty())
531 Info<<
"Subsetting surface " << surfaceName_
548 bitSet includeMap(surface_.size());
550 for (
const label zonei : zoneIndices)
556 if (includeMap.none())
559 <<
"Patch selection results in an empty surface"
560 <<
" - ignoring" <<
nl;
562 else if (!includeMap.all())
569 <<
"Bad surface subset (empty)"
570 <<
" - skip and hope for the best" <<
nl;
575 surface_.transfer(subSurf);
602 sampleElements_.clear();
603 samplePoints_.clear();
618 treeBoundBox bb(surface_.points(), surface_.meshPoints());
629 <<
"Surface " << surfaceName_
630 <<
" does not overlap bounding box of mesh " <<
mesh().
bounds()
634 const vector span(bb.span());
636 bb.
min() += (0.5-1
e-6)*span;
637 bb.
max() -= (0.5-1
e-6)*span;
642 const vector span(bb.span());
643 bb.
min() -= 0.5*span;
644 bb.
max() += 0.5*span;
652 return update(meshSearcher);
666 return update(meshSearcher);
675 return sampleOnFaces(sampler);
684 return sampleOnFaces(sampler);
693 return sampleOnFaces(sampler);
702 return sampleOnFaces(sampler);
711 return sampleOnFaces(sampler);
720 return sampleOnPoints(interpolator);
729 return sampleOnPoints(interpolator);
737 return sampleOnPoints(interpolator);
746 return sampleOnPoints(interpolator);
755 return sampleOnPoints(interpolator);
761 os <<
"meshedSurface: " <<
name() <<
" :"
762 <<
" surface:" << surfaceName_
763 <<
" faces:" << faces().size()
764 <<
" points:" <<
points().size()
765 <<
" zoneids:" << zoneIds().size();
int debug
Static debugging option.
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
virtual bool needsUpdate() const
Does the surface need an update?
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
A class for handling words, derived from Foam::string.
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search.
samplingSource
Types of sampling regions.
sampledMeshedSurface(const word &name, const polyMesh &mesh, const word &surfaceName, const samplingSource sampleSource)
Construct from components.
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;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
A class for managing temporary objects.
static IOobject selectReadIO(const word &name, const Time &runTime)
labelRange range() const
The start/size range of this zone.
Standard boundBox with extra functionality for use in octree.
const cellList & cells() const
Ostream & endl(Ostream &os)
Add newline and flush stream.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
virtual bool update()
Update the surface as required.
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.
Field< vector > vectorField
Specialisation of Field<T> for vector.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
Tuple2< scalar, label > nearInfo
Private class for finding nearest.
static const dictionary null
An empty dictionary, which is also the parent for all dictionaries.
bool intersect(const boundBox &bb)
Intersection (union) with the second box.
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
void transfer(List< T > &list)
addNamedToRunTimeSelectionTable(topoSetCellSource, badQualityToCell, word, badQuality)
Extract name (as a word) from an object, typically using its name() method.
An abstract class for surfaces with sampling.
labelList findMatching(const StringListType &input, const wordRes &whitelist, const wordRes &blacklist=wordRes(), AccessOp aop=noOp())
Return ids for items with matching names.
virtual void print(Ostream &os) const
Write.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
void resize(const label newSize)
Adjust allocated size of list.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
PointHit< point > pointHit
Macros for easy insertion into run-time selection tables.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
const boundBox & bounds() const
Return mesh bounding box.
const vectorField & cellCentres() const
virtual const faceList & faces() const
Return raw faces.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
PointIndexHit< point > pointIndexHit
FlatOutput< Container > flatOutput(const Container &obj, label len=0)
Global flatOutput function.
MeshedSurface< face > meshedSurface
A surface zone on a MeshedSurface.
const vectorField & faceCentres() const
const dimensionedScalar e
Elementary charge.
A List of wordRe with additional matching capabilities.
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
virtual void clearGeom() const
Additional cleanup when clearing the geometry.
void operator()(nearInfo &x, const nearInfo &y) const
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
const Time & time() const
Return the top-level database.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
virtual bool expire()
Mark the surface as needing an update.
virtual tmp< scalarField > sample(const interpolation< scalar > &sampler) const
Sample volume field onto surface faces.
A 2-tuple for storing two objects of dissimilar types. The container is similar in purpose to std::pa...
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
vector point
Point is a vector.
static wordRes uniq(const UList< wordRe > &input)
Return a wordRes with duplicate entries filtered out.
const word & constant() const
Return constant name.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
label size() const
The surface size is the number of faces.
bool interpolate() const
Interpolation to nodes requested for surface.
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
UIndirectList< label > labelUIndList
UIndirectList of labels.