Go to the documentation of this file.
45 Foam::sampledTriSurfaceMesh::samplingSourceNames_
48 { samplingSource::insideCells,
"insideCells" },
49 { samplingSource::boundaryFaces,
"boundaryFaces" },
59 sampledTriSurfaceMesh,
75 if (
y.first() <
x.first())
101 const surfZone& zn = zoneLst[zonei];
112 Foam::sampledTriSurfaceMesh::nonCoupledboundaryTree()
const
117 if (!boundaryTreePtr_.valid())
130 bndFaces[bndI++] = pp.start()+i;
134 bndFaces.setSize(bndI);
140 overallBb = overallBb.extend(
rndGen, 1
e-4);
141 overallBb.min() -=
point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
142 overallBb.max() +=
point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
144 boundaryTreePtr_.reset
146 new indexedOctree<treeDataFace>
162 return *boundaryTreePtr_;
170 const pointField& fc = surface_.faceCentres();
172 List<nearInfo> nearest(fc.size());
176 globalIndex globalCells(onBoundary() ?
mesh().nFaces() :
mesh().nCells());
180 near.first() = GREAT;
184 if (sampleSource_ ==
cells)
188 const indexedOctree<treeDataCell>& cellTree = meshSearcher.cellTree();
200 nearest[triI].second() = globalCells.toGlobal(
nearInfo.index());
204 else if (sampleSource_ == insideCells)
208 const indexedOctree<treeDataCell>& cellTree = meshSearcher.cellTree();
212 if (cellTree.bb().contains(fc[triI]))
214 const label index = cellTree.findInside(fc[triI]);
217 nearest[triI].first() = 0.0;
218 nearest[triI].second() = globalCells.toGlobal(index);
230 const indexedOctree<treeDataFace>& bTree = nonCoupledboundaryTree();
242 nearest[triI].second() = globalCells.toGlobal
244 bTree.shapes().faceLabels()[
nearInfo.index()]
257 labelList cellOrFaceLabels(fc.size(), -1);
262 if (nearest[triI].second() ==
labelMax)
266 else if (globalCells.isLocal(nearest[triI].second()))
268 cellOrFaceLabels[triI] = globalCells.toLocal
270 nearest[triI].second()
279 Pout<<
"Local out of faces:" << cellOrFaceLabels.size()
280 <<
" keeping:" << nFound <<
endl;
286 const triSurface&
s = surface_;
293 labelList reversePointMap(
s.points().size(), -1);
302 Map<label> zoneSizes;
321 zoneSizes.set(patchi, 0);
332 if (cellOrFaceLabels[facei] != -1)
335 const label regionid =
f.region();
337 auto fnd = zoneSizes.find(regionid);
345 zoneSizes.insert(regionid, 1);
351 zoneIds_[newFacei] = regionid;
355 for (
const label labi :
f)
357 if (reversePointMap[labi] == -1)
368 zoneIds_.setSize(newFacei);
381 const label regionid = iter.key();
384 auto fnd = zoneNames.cfind(regionid);
394 zoneLst[zoneI] = surfZone
418 const label zonei = zoneIds_[facei];
419 label sortedFacei = zoneLst[zonei].start() + zoneLst[zonei].size()++;
420 sortedFaceMap[sortedFacei] =
faceMap[facei];
424 setZoneMap(zoneLst, zoneIds_);
427 faceMap.transfer(sortedFaceMap);
441 faceList& surfFaces = this->storedFaces();
442 surfFaces.setSize(
faceMap.size());
444 this->storedZones().transfer(zoneLst);
448 const labelledTri& origF =
s[
faceMap[facei]];
449 face&
f = surfFaces[facei];
453 reversePointMap[origF[0]],
454 reversePointMap[origF[1]],
455 reversePointMap[origF[2]]
458 for (
const label labi :
f)
460 pointToFace[labi] = facei;
464 this->storedPoints() =
pointField(
s.points(), pointMap);
479 samplePoints_.setSize(pointMap.size());
480 sampleElements_.setSize(pointMap.size(), -1);
482 if (sampleSource_ ==
cells)
490 label celli = cellOrFaceLabels[pointToFace[pointi]];
491 sampleElements_[pointi] = celli;
499 sampleElements_[pointi],
500 meshSearcher.decompMode()
504 samplePoints_[pointi] = pt;
509 const cell& cFaces =
mesh().
cells()[celli];
511 scalar minDistSqr = VGREAT;
517 if (info.distance() < minDistSqr)
519 minDistSqr = info.distance();
520 samplePoints_[pointi] = info.rawPoint();
526 else if (sampleSource_ == insideCells)
534 label celli = cellOrFaceLabels[pointToFace[pointi]];
535 sampleElements_[pointi] = celli;
536 samplePoints_[pointi] = pt;
548 label facei = cellOrFaceLabels[pointToFace[pointi]];
549 sampleElements_[pointi] = facei;
550 samplePoints_[pointi] =
mesh().
faces()[facei].nearestPoint
569 sampleElements_.
transfer(cellOrFaceLabels);
570 samplePoints_.clear();
577 OFstream str(
mesh().time().
path()/
"surfaceToMesh.obj");
578 Info<<
"Dumping correspondence from local surface (points or faces)"
579 <<
" to mesh (cells or faces) to " << str.name() <<
endl;
591 forAll(samplePoints_, pointi)
599 label elemi = sampleElements_[pointi];
602 str <<
"l " << vertI-2 <<
' ' << vertI-1 <<
' ' << vertI <<
nl;
608 forAll(sampleElements_, triI)
613 label elemi = sampleElements_[triI];
616 str <<
"l " << vertI-1 <<
' ' << vertI <<
nl;
621 needsUpdate_ =
false;
632 const word& surfaceName,
651 sampleSource_(sampleSource),
684 sampleSource_(samplingSourceNames_.get(
"source",
dict)),
699 const word& sampleSourceName
718 sampleSource_(samplingSourceNames_[sampleSourceName]),
754 originalIds_.clear();
755 boundaryTreePtr_.clear();
756 sampleElements_.clear();
757 samplePoints_.clear();
774 surface_.triSurface::points(),
775 surface_.triSurface::meshPoints()
787 <<
"Surface " << surface_.searchableSurface::name()
788 <<
" does not overlap bounding box of mesh " <<
mesh().
bounds()
792 const vector span(bb.span());
794 bb.
min() += (0.5-1
e-6)*span;
795 bb.
max() -= (0.5-1
e-6)*span;
800 const vector span(bb.span());
801 bb.
min() -= 0.5*span;
802 bb.
max() += 0.5*span;
810 return update(meshSearcher);
824 return update(meshSearcher);
833 return sampleOnFaces(sampler);
842 return sampleOnFaces(sampler);
851 return sampleOnFaces(sampler);
860 return sampleOnFaces(sampler);
869 return sampleOnFaces(sampler);
878 return sampleOnPoints(interpolator);
887 return sampleOnPoints(interpolator);
895 return sampleOnPoints(interpolator);
904 return sampleOnPoints(interpolator);
913 return sampleOnPoints(interpolator);
919 os <<
"sampledTriSurfaceMesh: " <<
name() <<
" :"
920 <<
" surface:" << surface_.objectRegistry::name()
921 <<
" faces:" << faces().size()
922 <<
" points:" <<
points().size()
923 <<
" zoneids:" << zoneIds().size();
virtual const pointField & points() const
Points of surface.
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
label start() const
Return start label of this zone in the face list.
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,...
Tuple2< scalar, label > nearInfo
Private class for finding nearest.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
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.
virtual bool needsUpdate() const
Does the surface need an update?
sampledTriSurfaceMesh(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 polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
virtual bool update()
Update the surface as required.
A class for managing temporary objects.
virtual tmp< scalarField > sample(const interpolation< scalar > &sampler) const
Sample volume field onto surface faces.
samplingSource
Types of communications.
A List obtained as a section of another List.
Standard boundBox with extra functionality for use in octree.
const cellList & cells() const
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
T lookupOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
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)
bool intersect(const boundBox &bb)
Intersection (union) with the second box.
List< geometricSurfacePatch > geometricSurfacePatchList
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Triangulated surface description with patch information.
messageStream Info
Information stream (uses stdout - output is on the master only)
A patch is a list of labels that address the faces in the global face list.
word name(const complex &c)
Return string representation of complex.
void transfer(List< T > &list)
An abstract class for surfaces with sampling.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
label size() const
Return size of this zone in the face list.
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
static void setZoneMap(const surfZoneList &zoneLst, labelList &zoneIds)
Set new zoneIds list based on the surfZoneList information.
List< face > faceList
A List of faces.
A surface zone on a MeshedSurface.
label ListType::const_reference const label start
const vectorField & faceCentres() const
Ostream & print(Ostream &os, UIntType value, char off='0', char on='1')
Print 0/1 bits in the (unsigned) integral type.
const dimensionedScalar e
Elementary charge.
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
const polyBoundaryMesh & patches
virtual void clearGeom() const
Additional cleanup when clearing the geometry.
void operator()(nearInfo &x, const nearInfo &y) const
List< surfZone > surfZoneList
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.
virtual ~sampledTriSurfaceMesh()
Destructor.
static word printf(const char *fmt, const PrimitiveType &val)
Use a printf-style formatter for a primitive.
const Time & time() const
Return the top-level database.
const polyMesh & mesh() const
Access to the underlying mesh.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
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.
void setSize(const label newSize)
Alias for resize(const label)
const word & constant() const
Return constant name.
virtual bool expire()
Mark the surface as needing an update.
bool interpolate() const
Interpolation to nodes requested for surface.
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
virtual void print(Ostream &os) const
Write.
UIndirectList< label > labelUIndList
UIndirectList of labels.