57void Foam::searchableSurfaceCollection::findNearest
60 scalarField& minDistSqr,
61 List<pointIndexHit>& nearestInfo,
62 labelList& nearestSurf
77 subGeom_[surfI].findNearest
81 transform_[surfI].localPosition(
samples),
90 if (hitInfo[pointi].hit())
94 point globalPt = transform_[surfI].globalPosition
98 hitInfo[pointi].rawPoint(),
105 if (distSqr < minDistSqr[pointi])
107 minDistSqr[pointi] = distSqr;
108 nearestInfo[pointi].setPoint(globalPt);
109 nearestInfo[pointi].setHit();
110 nearestInfo[pointi].setIndex
112 hitInfo[pointi].index()
113 + indexOffset_[surfI]
115 nearestSurf[pointi] = surfI;
125void Foam::searchableSurfaceCollection::sortHits
127 const List<pointIndexHit>& info,
128 List<List<pointIndexHit>>& surfInfo,
137 if (info[pointi].hit())
139 label index = info[pointi].index();
140 label surfI =
findLower(indexOffset_, index+1);
146 surfInfo.setSize(subGeom_.size());
148 infoMap.setSize(subGeom_.size());
152 surfInfo[surfI].setSize(nHits[surfI]);
153 infoMap[surfI].setSize(nHits[surfI]);
159 if (info[pointi].hit())
161 label index = info[pointi].index();
162 label surfI =
findLower(indexOffset_, index+1);
166 label localI = nHits[surfI]++;
170 info[pointi].rawPoint(),
171 index-indexOffset_[surfI]
173 infoMap[surfI][localI] = pointi;
188 instance_(
dict.size()),
190 transform_(
dict.size()),
191 subGeom_(
dict.size()),
192 mergeSubRegions_(
dict.get<
bool>(
"mergeSubRegions")),
193 indexOffset_(
dict.size()+1)
198 label startIndex = 0;
203 instance_[surfI] = dEntry.keyword();
214 coordinateSystem::typeName_(),
224 <<
"--> FOAM IOWarning :" <<
nl
225 <<
" Found [v1806] '"
226 << coordinateSystem::typeName_()
227 <<
"' entry within transform dictionary" <<
nl
239 const word subGeomName(sDict.
get<
word>(
"surface"));
248 if (
s.size() !=
s.globalSize())
251 <<
"Cannot use a distributed surface in a collection."
255 subGeom_.
set(surfI, &
s);
257 indexOffset_[surfI] = startIndex;
258 startIndex += subGeom_[surfI].
size();
260 Info<<
" instance : " << instance_[surfI] <<
endl;
261 Info<<
" surface : " <<
s.name() <<
endl;
262 Info<<
" scale : " << scale_[surfI] <<
endl;
263 Info<<
" transform: " << transform_[surfI] <<
endl;
268 indexOffset_[surfI] = startIndex;
272 transform_.setSize(surfI);
281 const boundBox& surfBb = subGeom_[surfI].bounds();
284 const point surfBbMin = transform_[surfI].globalPosition
292 const point surfBbMax = transform_[surfI].globalPosition
317 if (regions_.empty())
319 regionOffset_.
setSize(subGeom_.size());
324 regionOffset_[surfI] = allRegions.
size();
326 if (mergeSubRegions_)
333 const wordList& subRegions = subGeom_[surfI].regions();
341 regions_.transfer(allRegions);
349 return indexOffset_.last();
357 auto& ctrs = tctrs.ref();
368 ctrs[coordI++] = transform_[surfI].globalPosition
397 scalar maxScale =
cmptMax(scale_[surfI]);
401 subGeom_[surfI].boundingSpheres(subCentres, subRadiusSqr);
405 centres[coordI] = transform_[surfI].globalPosition
413 radiusSqr[coordI] = maxScale*subRadiusSqr[i];
428 nPoints += subGeom_[surfI].points()().size();
432 auto& pts = tpts.ref();
443 pts[
nPoints++] = transform_[surfI].globalPosition
458void Foam::searchableSurfaceCollection::findNearest
499 transform_[surfI].localPosition
509 transform_[surfI].localPosition
516 subGeom_[surfI].findLine(e0, e1, hitInfo);
520 if (hitInfo[pointi].hit())
523 nearest[pointi] = transform_[surfI].globalPosition
527 hitInfo[pointi].rawPoint(),
531 info[pointi] = hitInfo[pointi];
532 info[pointi].rawPoint() = nearest[pointi];
533 info[pointi].setIndex
535 hitInfo[pointi].index()
536 + indexOffset_[surfI]
548 if (info[pointi].hit())
550 vector n(end[pointi] - start[pointi]);
551 scalar magN =
mag(
n);
557 scalar
s = ((info[pointi].rawPoint()-start[pointi])&
n);
562 <<
"point:" << info[pointi]
564 <<
" outside vector "
565 <<
" start:" << start[pointi]
566 <<
" end:" << end[pointi]
584 findLine(start, end, info);
597 findLine(start, end, nearestInfo);
599 info.setSize(start.
size());
602 if (nearestInfo[pointi].hit())
604 info[pointi].setSize(1);
605 info[pointi][0] = nearestInfo[pointi];
609 info[pointi].
clear();
621 if (subGeom_.size() == 0)
623 else if (subGeom_.size() == 1)
625 if (mergeSubRegions_)
628 region = regionOffset_[0];
632 subGeom_[0].getRegion(info, region);
643 sortHits(info, surfInfo, infoMap);
650 if (mergeSubRegions_)
658 region[map[i]] = regionOffset_[surfI];
667 subGeom_[surfI].getRegion(surfInfo[surfI], surfRegion);
672 region[map[i]] = regionOffset_[surfI] + surfRegion[i];
686 if (subGeom_.size() == 0)
688 else if (subGeom_.size() == 1)
690 subGeom_[0].getNormal(info, normal);
700 sortHits(info, surfInfo, infoMap);
708 subGeom_[surfI].getNormal(surfInfo[surfI], surfNormal);
711 surfNormal = transform_[surfI].globalVector(surfNormal);
716 normal[map[i]] = surfNormal[i];
730 <<
"Volume type not supported for collection."
738 const bool keepNonLocal,
755 subGeom_[surfI].distribute
770 subGeom_[surfI].setField
777 subGeom_[surfI].size(),
792 if (subGeom_.size() == 0)
794 else if (subGeom_.size() == 1)
796 subGeom_[0].getField(info, values);
806 sortHits(info, surfInfo, infoMap);
812 subGeom_[surfI].getField(surfInfo[surfI], surfValues);
814 if (surfValues.
size())
817 values.setSize(info.
size());
822 values[map[i]] = surfValues[i];
Various functions to operate on Lists.
Macros for easy insertion into run-time selection tables.
#define addNamedToRunTimeSelectionTable(baseType, thisType, argNames, lookupName)
Add to construction table with 'lookupName' as the key.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
void append(const T &val)
Copy append an element to the end of this list.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const word & name() const noexcept
Return the object name.
const objectRegistry & db() const noexcept
Return the local objectRegistry.
void setSize(const label n)
Alias for resize()
void clear()
Clear the list, i.e. set size to zero.
A List obtained as a section of another List.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
void size(const label n)
Older name for setAddressableSize.
const T * set(const label i) const
void setSize(const label n)
Alias for resize()
label size() const noexcept
The number of elements in the list.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A bounding box defined in terms of min/max extrema points.
const point & min() const
Minimum describing the bounding box.
const point & max() const
Maximum describing the bounding box.
static const boundBox invertedBox
A large inverted boundBox: min/max == +/- ROOTVGREAT.
A Cartesian coordinate system.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
dictionary * findDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX)
Find and return a sub-dictionary pointer if present.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
A keyword and a list of tokens is an 'entry'.
static bool warnAboutAge(const int version) noexcept
Test if an age warning should be emitted.
Type & lookupObjectRef(const word &name, const bool recursive=false) const
Set of transformed searchableSurfaces. Does not do boolean operations so when meshing might find part...
virtual void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Find first intersection on segment from start to end.
virtual void distribute(const List< treeBoundBox > &, const bool keepNonLocal, autoPtr< mapDistribute > &faceMap, autoPtr< mapDistribute > &pointMap)
Set bounds of surface. Bounds currently set as list of.
virtual void findLineAll(const pointField &start, const pointField &end, List< List< pointIndexHit > > &) const
Get all intersections in order from start to end.
virtual void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Return any intersection on segment from start to end.
virtual label size() const
Range of local indices that can be returned.
virtual void boundingSpheres(pointField ¢res, scalarField &radiusSqr) const
Get bounding spheres (centre and radius squared), one per element.
virtual void getVolumeType(const pointField &points, List< volumeType > &volType) const
Determine type (inside/outside/mixed) for point.
virtual void setField(const labelList &values)
WIP. Store element-wise field.
virtual void getField(const List< pointIndexHit > &, labelList &) const
WIP. From a set of hits (points and.
virtual void getNormal(const List< pointIndexHit > &, vectorField &normal) const
From a set of points and indices get the normal.
virtual void getRegion(const List< pointIndexHit > &, labelList ®ion) const
From a set of points and indices get the region.
virtual const wordList & regions() const
Names of regions.
virtual ~searchableSurfaceCollection()
Destructor.
virtual tmp< pointField > coordinates() const
Get representative set of element coordinates.
virtual tmp< pointField > points() const
Get the points that define the surface.
Base class of (analytical or triangulated) surface. Encapsulates all the search routines....
virtual const boundBox & bounds() const
Return const reference to boundBox.
splitCell * master() const
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Foam::word regionName(Foam::polyMesh::defaultRegion)
PtrList< coordinateSystem > coordinates(solidRegions.size())
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
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))
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
List< label > labelList
A List of labels.
PointIndexHit< point > pointIndexHit
A PointIndexHit for 3D points.
dimensioned< Type > cmptDivide(const dimensioned< Type > &, const dimensioned< Type > &)
messageStream Info
Information stream (stdout output on master, null elsewhere)
label findLower(const ListType &input, const T &val, const label start, const ComparePredicate &comp)
vector point
Point is a vector.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
List< labelList > labelListList
A List of labelList.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManip< error > abort(error &err)
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.
scalarField samples(nIntervals, Zero)