diff --git a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C
index 60ad087b70a3d68239124c8a3960f79b3a4b5eec..40ad36bbc3a4e5b185388fbd006c220853d1694b 100644
--- a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C
+++ b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C
@@ -912,41 +912,6 @@ Foam::distributedTriSurfaceMesh::independentlyDistributedBbs
}
-void Foam::distributedTriSurfaceMesh::calcBounds
-(
- boundBox& bb,
- label& nPoints
-) const
-{
- // Unfortunately nPoints constructs meshPoints() so do compact version
- // ourselves
-
- PackedBoolList pointIsUsed(points().size());
-
- nPoints = 0;
- bb.min() = point(VGREAT, VGREAT, VGREAT);
- bb.max() = point(-VGREAT, -VGREAT, -VGREAT);
-
- const triSurface& s = static_cast<const triSurface&>(*this);
-
- forAll(s, triI)
- {
- const labelledTri& f = s[triI];
-
- forAll(f, fp)
- {
- label pointI = f[fp];
- if (pointIsUsed.set(pointI, 1u))
- {
- bb.min() = ::Foam::min(bb.min(), points()[pointI]);
- bb.max() = ::Foam::max(bb.max(), points()[pointI]);
- nPoints++;
- }
- }
- }
-}
-
-
// Does any part of triangle overlap bb.
bool Foam::distributedTriSurfaceMesh::overlaps
(
@@ -1935,20 +1900,7 @@ void Foam::distributedTriSurfaceMesh::getNormal
{
if (!Pstream::parRun())
{
- normal.setSize(info.size());
-
- forAll(info, i)
- {
- if (info[i].hit())
- {
- normal[i] = faceNormals()[info[i].index()];
- }
- else
- {
- // Set to what?
- normal[i] = vector::zero;
- }
- }
+ triSurfaceMesh::getNormal(info, normal);
return;
}
@@ -2006,70 +1958,64 @@ void Foam::distributedTriSurfaceMesh::getNormal
void Foam::distributedTriSurfaceMesh::getField
(
- const word& fieldName,
const List<pointIndexHit>& info,
labelList& values
) const
{
- const triSurfaceLabelField& fld = lookupObject<triSurfaceLabelField>
- (
- fieldName
- );
-
-
if (!Pstream::parRun())
{
- values.setSize(info.size());
- forAll(info, i)
- {
- if (info[i].hit())
- {
- values[i] = fld[info[i].index()];
- }
- }
+ triSurfaceMesh::getField(info, values);
return;
}
+ if (foundObject<triSurfaceLabelField>("values"))
+ {
+ const triSurfaceLabelField& fld = lookupObject<triSurfaceLabelField>
+ (
+ "values"
+ );
- // Get query data (= local index of triangle)
- // ~~~~~~~~~~~~~~
- labelList triangleIndex(info.size());
- autoPtr<mapDistribute> mapPtr
- (
- calcLocalQueries
+ // Get query data (= local index of triangle)
+ // ~~~~~~~~~~~~~~
+
+ labelList triangleIndex(info.size());
+ autoPtr<mapDistribute> mapPtr
(
- info,
- triangleIndex
- )
- );
- const mapDistribute& map = mapPtr();
+ calcLocalQueries
+ (
+ info,
+ triangleIndex
+ )
+ );
+ const mapDistribute& map = mapPtr();
- // Do my tests
- // ~~~~~~~~~~~
+ // Do my tests
+ // ~~~~~~~~~~~
- values.setSize(triangleIndex.size());
+ values.setSize(triangleIndex.size());
- forAll(triangleIndex, i)
- {
- label triI = triangleIndex[i];
- values[i] = fld[triI];
- }
+ forAll(triangleIndex, i)
+ {
+ label triI = triangleIndex[i];
+ values[i] = fld[triI];
+ }
- // Send back results
- // ~~~~~~~~~~~~~~~~~
+ // Send back results
+ // ~~~~~~~~~~~~~~~~~
- map.distribute
- (
- Pstream::nonBlocking,
- List<labelPair>(0),
- info.size(),
- map.constructMap(), // what to send
- map.subMap(), // what to receive
- values
- );
+ map.distribute
+ (
+ Pstream::nonBlocking,
+ List<labelPair>(0),
+ info.size(),
+ map.constructMap(), // what to send
+ map.subMap(), // what to receive
+ values
+ );
+ }
}
diff --git a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H
index dba72d2abf4b1491be861f6e5dc825533d9b1427..577d12e4c3544fad2b73e8b4343d5c2f80682262 100644
--- a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H
+++ b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H
@@ -217,9 +217,6 @@ private:
const triSurface&
);
- //- Calculate surface bounding box
- void calcBounds(boundBox& bb, label& nPoints) const;
-
//- Does any part of triangle overlap bb.
static bool overlaps
(
@@ -418,7 +415,7 @@ public:
// Should really be split into a routine to determine decomposition
// and one that does actual distribution but determining
// decomposition with duplicate triangle merging requires
- // same amoun as work as actual distribution.
+ // same amount as work as actual distribution.
virtual void distribute
(
const List<treeBoundBox>&,
@@ -430,14 +427,9 @@ public:
// Other
- //- Specific to triSurfaceMesh: from a set of hits (points and
+ //- WIP. From a set of hits (points and
// indices) get the specified field. Misses do not get set.
- virtual void getField
- (
- const word& fieldName,
- const List<pointIndexHit>&,
- labelList& values
- ) const;
+ virtual void getField(const List<pointIndexHit>&, labelList&) const;
//- Subset the part of surface that is overlapping bounds.
static triSurface overlappingSurface
diff --git a/src/meshTools/searchableSurface/searchableBox.C b/src/meshTools/searchableSurface/searchableBox.C
index a4cc593562a6460ac80f4998d7dc570c4651da71..cb854c3772552b9602a50531321c80a07268beac 100644
--- a/src/meshTools/searchableSurface/searchableBox.C
+++ b/src/meshTools/searchableSurface/searchableBox.C
@@ -229,6 +229,21 @@ const Foam::wordList& Foam::searchableBox::regions() const
}
+Foam::pointField Foam::searchableBox::coordinates() const
+{
+ pointField ctrs(6);
+
+ const pointField pts = treeBoundBox::points();
+ const faceList& fcs = treeBoundBox::faces;
+
+ forAll(fcs, i)
+ {
+ ctrs[i] = fcs[i].centre(pts);
+ }
+ return ctrs;
+}
+
+
Foam::pointIndexHit Foam::searchableBox::findNearest
(
const point& sample,
diff --git a/src/meshTools/searchableSurface/searchableBox.H b/src/meshTools/searchableSurface/searchableBox.H
index 9d1f2e8217ef3db841ebf9d1c8b4f5796c34b866..77e0ecd98d8c33b06a2d710eee10f47b1f8de70c 100644
--- a/src/meshTools/searchableSurface/searchableBox.H
+++ b/src/meshTools/searchableSurface/searchableBox.H
@@ -127,6 +127,9 @@ public:
return 6;
}
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const;
// Single point queries.
diff --git a/src/meshTools/searchableSurface/searchableCylinder.C b/src/meshTools/searchableSurface/searchableCylinder.C
index fcaa42194aab42ce1f06f9e199f215548b92c627..0df222ccb0fc929bb882c8c9b1fafa1535097ec2 100644
--- a/src/meshTools/searchableSurface/searchableCylinder.C
+++ b/src/meshTools/searchableSurface/searchableCylinder.C
@@ -40,6 +40,14 @@ addToRunTimeSelectionTable(searchableSurface, searchableCylinder, dict);
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+Foam::pointField Foam::searchableCylinder::coordinates() const
+{
+ pointField ctrs(1, 0.5*(point1_ + point2_));
+
+ return ctrs;
+}
+
+
Foam::pointIndexHit Foam::searchableCylinder::findNearest
(
const point& sample,
diff --git a/src/meshTools/searchableSurface/searchableCylinder.H b/src/meshTools/searchableSurface/searchableCylinder.H
index cae0f058db2863db95ff479381d207c8ad9e57cf..98b4d91e4117eceb782b45a976928820a211c759 100644
--- a/src/meshTools/searchableSurface/searchableCylinder.H
+++ b/src/meshTools/searchableSurface/searchableCylinder.H
@@ -150,6 +150,10 @@ public:
return 1;
}
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const;
+
// Multiple point queries.
diff --git a/src/meshTools/searchableSurface/searchablePlane.H b/src/meshTools/searchableSurface/searchablePlane.H
index 87777e3209bc6763b1e8cd155eb55ce4f3f7c4e1..c0b58cc8ac2d0c086deca4a4dda306756a295623 100644
--- a/src/meshTools/searchableSurface/searchablePlane.H
+++ b/src/meshTools/searchableSurface/searchablePlane.H
@@ -26,7 +26,7 @@ Class
Foam::searchablePlane
Description
- Searching on plane. See plane.H
+ Searching on (infinite) plane. See plane.H
SourceFiles
searchablePlane.C
@@ -122,6 +122,14 @@ public:
return 1;
}
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const
+ {
+ //notImplemented("searchablePlane::coordinates()")
+ return pointField(1, refPoint());
+ }
+
// Multiple point queries.
diff --git a/src/meshTools/searchableSurface/searchablePlate.H b/src/meshTools/searchableSurface/searchablePlate.H
index 7bc900f554bd3dc7a281873cf4b850077d5cf695..b05414bb6cb0e2dd764b732ff0e34f97a59539db 100644
--- a/src/meshTools/searchableSurface/searchablePlate.H
+++ b/src/meshTools/searchableSurface/searchablePlate.H
@@ -145,6 +145,13 @@ public:
return 1;
}
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const
+ {
+ return pointField(1, origin_);
+ }
+
// Multiple point queries.
diff --git a/src/meshTools/searchableSurface/searchableSphere.H b/src/meshTools/searchableSurface/searchableSphere.H
index 297890b9db1b6e4568db5a65ca6d2055fd25777e..fb4aeb22f81c9d3c37afc69da4c549880af66362 100644
--- a/src/meshTools/searchableSurface/searchableSphere.H
+++ b/src/meshTools/searchableSurface/searchableSphere.H
@@ -133,6 +133,13 @@ public:
return 1;
}
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const
+ {
+ return pointField(1, centre_);
+ }
+
// Multiple point queries.
diff --git a/src/meshTools/searchableSurface/searchableSurface.H b/src/meshTools/searchableSurface/searchableSurface.H
index 35b9dc9fb5a1ff9f3ba054492f7130ab2beee4a2..daf8b76daa7131344ef19d86263cda4883c1d6c9 100644
--- a/src/meshTools/searchableSurface/searchableSurface.H
+++ b/src/meshTools/searchableSurface/searchableSurface.H
@@ -190,6 +190,10 @@ public:
return size();
}
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const = 0;
+
// Single point queries.
@@ -319,6 +323,18 @@ public:
)
{}
+ //- WIP. Store element-wise field.
+ virtual void setField(const labelList& values)
+ {}
+
+ //- WIP. From a set of hits (points and
+ // indices) get the specified field. Misses do not get set. Return
+ // empty field if not supported.
+ virtual void getField(const List<pointIndexHit>&, labelList& values)
+ const
+ {
+ values.clear();
+ }
};
diff --git a/src/meshTools/searchableSurface/searchableSurfaceCollection.C b/src/meshTools/searchableSurface/searchableSurfaceCollection.C
index ba1dd32988e47a2aa53eb001e1014ad1686f83d2..9279087c1908b2298e61c15b71215646104ced6e 100644
--- a/src/meshTools/searchableSurface/searchableSurfaceCollection.C
+++ b/src/meshTools/searchableSurface/searchableSurfaceCollection.C
@@ -101,7 +101,11 @@ void Foam::searchableSurfaceCollection::findNearest
minDistSqr[pointI] = distSqr;
nearestInfo[pointI].setPoint(globalPt);
nearestInfo[pointI].setHit();
- nearestInfo[pointI].setIndex(hitInfo[pointI].index());
+ nearestInfo[pointI].setIndex
+ (
+ hitInfo[pointI].index()
+ + indexOffset_[surfI]
+ );
nearestSurf[pointI] = surfI;
}
}
@@ -110,6 +114,62 @@ void Foam::searchableSurfaceCollection::findNearest
}
+// Sort hits into per-surface bins. Misses are rejected. Maintains map back
+// to position
+void Foam::searchableSurfaceCollection::sortHits
+(
+ const List<pointIndexHit>& info,
+ List<List<pointIndexHit> >& surfInfo,
+ labelListList& infoMap
+) const
+{
+ // Count hits per surface.
+ labelList nHits(subGeom_.size(), 0);
+
+ forAll(info, pointI)
+ {
+ if (info[pointI].hit())
+ {
+ label index = info[pointI].index();
+ label surfI = findLower(indexOffset_, index+1);
+ nHits[surfI]++;
+ }
+ }
+
+ // Per surface the hit
+ surfInfo.setSize(subGeom_.size());
+ // Per surface the original position
+ infoMap.setSize(subGeom_.size());
+
+ forAll(surfInfo, surfI)
+ {
+ surfInfo[surfI].setSize(nHits[surfI]);
+ infoMap[surfI].setSize(nHits[surfI]);
+ }
+ nHits = 0;
+
+ forAll(info, pointI)
+ {
+ if (info[pointI].hit())
+ {
+ label index = info[pointI].index();
+ label surfI = findLower(indexOffset_, index+1);
+
+ // Store for correct surface and adapt indices back to local
+ // ones
+ label localI = nHits[surfI]++;
+ surfInfo[surfI][localI] = pointIndexHit
+ (
+ info[pointI].hit(),
+ info[pointI].rawPoint(),
+ index-indexOffset_[surfI]
+ );
+ infoMap[surfI][localI] = pointI;
+ }
+ }
+}
+
+
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::searchableSurfaceCollection::searchableSurfaceCollection
@@ -123,11 +183,13 @@ Foam::searchableSurfaceCollection::searchableSurfaceCollection
scale_(dict.size()),
transform_(dict.size()),
subGeom_(dict.size()),
- mergeSubRegions_(dict.lookup("mergeSubRegions"))
+ mergeSubRegions_(dict.lookup("mergeSubRegions")),
+ indexOffset_(dict.size()+1)
{
Info<< "SearchableCollection : " << name() << endl;
label surfI = 0;
+ label startIndex = 0;
forAllConstIter(dictionary, dict, iter)
{
if (dict.isDict(iter().keyword()))
@@ -153,8 +215,24 @@ Foam::searchableSurfaceCollection::searchableSurfaceCollection
const searchableSurface& s =
io.db().lookupObject<searchableSurface>(subGeomName);
+ // I don't know yet how to handle the globalSize combined with
+ // regionOffset. Would cause non-consecutive indices locally
+ // if all indices offset by globalSize() of the local region...
+ if (s.size() != s.globalSize())
+ {
+ FatalErrorIn
+ (
+ "searchableSurfaceCollection::searchableSurfaceCollection"
+ "(const IOobject&, const dictionary&)"
+ ) << "Cannot use a distributed surface in a collection."
+ << exit(FatalError);
+ }
+
subGeom_.set(surfI, &const_cast<searchableSurface&>(s));
+ indexOffset_[surfI] = startIndex;
+ startIndex += subGeom_[surfI].size();
+
Info<< " instance : " << instance_[surfI] << endl;
Info<< " surface : " << s.name() << endl;
Info<< " scale : " << scale_[surfI] << endl;
@@ -163,10 +241,13 @@ Foam::searchableSurfaceCollection::searchableSurfaceCollection
surfI++;
}
}
+ indexOffset_[surfI] = startIndex;
+
instance_.setSize(surfI);
scale_.setSize(surfI);
transform_.setSize(surfI);
subGeom_.setSize(surfI);
+ indexOffset_.setSize(surfI+1);
}
@@ -212,12 +293,36 @@ const Foam::wordList& Foam::searchableSurfaceCollection::regions() const
Foam::label Foam::searchableSurfaceCollection::size() const
{
- label n = 0;
+ return indexOffset_[indexOffset_.size()-1];
+}
+
+
+Foam::pointField Foam::searchableSurfaceCollection::coordinates() const
+{
+ // Get overall size
+ pointField coords(size());
+
+ // Append individual coordinates
+ label coordI = 0;
+
forAll(subGeom_, surfI)
{
- n += subGeom_[surfI].size();
+ const pointField subCoords = subGeom_[surfI].coordinates();
+
+ forAll(subCoords, i)
+ {
+ coords[coordI++] = transform_[surfI].globalPosition
+ (
+ cmptMultiply
+ (
+ subCoords[i],
+ scale_[surfI]
+ )
+ );
+ }
}
- return n;
+
+ return coords;
}
@@ -296,6 +401,11 @@ void Foam::searchableSurfaceCollection::findLine
);
info[pointI] = hitInfo[pointI];
info[pointI].rawPoint() = nearest[pointI];
+ info[pointI].setIndex
+ (
+ hitInfo[pointI].index()
+ + indexOffset_[surfI]
+ );
}
}
}
@@ -397,82 +507,42 @@ void Foam::searchableSurfaceCollection::getRegion
}
else
{
- region.setSize(info.size());
- region = -1;
+ // Multiple surfaces. Sort by surface.
- // Which region did point come from. Retest for now to see which
- // surface it originates from - crap solution! Should use global indices
- // in index inside pointIndexHit to do this better.
+ // Per surface the hit
+ List<List<pointIndexHit> > surfInfo;
+ // Per surface the original position
+ List<List<label> > infoMap;
+ sortHits(info, surfInfo, infoMap);
- pointField samples(info.size());
- forAll(info, pointI)
- {
- if (info[pointI].hit())
- {
- samples[pointI] = info[pointI].hitPoint();
- }
- else
- {
- samples[pointI] = vector::zero;
- }
- }
- //scalarField minDistSqr(info.size(), SMALL);
- scalarField minDistSqr(info.size(), GREAT);
+ region.setSize(info.size());
+ region = -1;
- labelList nearestSurf;
- List<pointIndexHit> nearestInfo;
- findNearest
- (
- samples,
- minDistSqr,
- nearestInfo,
- nearestSurf
- );
+ // Do region tests
- // Check
+ if (mergeSubRegions_)
{
- forAll(info, pointI)
+ // Actually no need for surfInfo. Just take region for surface.
+ forAll(infoMap, surfI)
{
- if (info[pointI].hit() && nearestSurf[pointI] == -1)
+ const labelList& map = infoMap[surfI];
+ forAll(map, i)
{
- FatalErrorIn
- (
- "searchableSurfaceCollection::getRegion(..)"
- ) << "pointI:" << pointI
- << " sample:" << samples[pointI]
- << " nearest:" << nearestInfo[pointI]
- << " nearestsurf:" << nearestSurf[pointI]
- << abort(FatalError);
+ region[map[i]] = regionOffset_[surfI];
}
}
}
-
- forAll(subGeom_, surfI)
+ else
{
- // Collect points from my surface
- labelList indices(findIndices(nearestSurf, surfI));
-
- if (mergeSubRegions_)
- {
- forAll(indices, i)
- {
- region[indices[i]] = regionOffset_[surfI];
- }
- }
- else
+ forAll(infoMap, surfI)
{
labelList surfRegion;
- subGeom_[surfI].getRegion
- (
- List<pointIndexHit>
- (
- UIndirectList<pointIndexHit>(info, indices)
- ),
- surfRegion
- );
- forAll(indices, i)
+ subGeom_[surfI].getRegion(surfInfo[surfI], surfRegion);
+
+ const labelList& map = infoMap[surfI];
+ forAll(map, i)
{
- region[indices[i]] = regionOffset_[surfI] + surfRegion[i];
+ region[map[i]] = regionOffset_[surfI] + surfRegion[i];
}
}
}
@@ -494,52 +564,26 @@ void Foam::searchableSurfaceCollection::getNormal
}
else
{
- normal.setSize(info.size());
+ // Multiple surfaces. Sort by surface.
- // See above - crap retest to find surface point originates from.
- pointField samples(info.size());
- forAll(info, pointI)
- {
- if (info[pointI].hit())
- {
- samples[pointI] = info[pointI].hitPoint();
- }
- else
- {
- samples[pointI] = vector::zero;
- }
- }
- //scalarField minDistSqr(info.size(), SMALL);
- scalarField minDistSqr(info.size(), GREAT);
-
- labelList nearestSurf;
- List<pointIndexHit> nearestInfo;
- findNearest
- (
- samples,
- minDistSqr,
- nearestInfo,
- nearestSurf
- );
+ // Per surface the hit
+ List<List<pointIndexHit> > surfInfo;
+ // Per surface the original position
+ List<List<label> > infoMap;
+ sortHits(info, surfInfo, infoMap);
+ normal.setSize(info.size());
- forAll(subGeom_, surfI)
+ // Do region tests
+ forAll(surfInfo, surfI)
{
- // Collect points from my surface
- labelList indices(findIndices(nearestSurf, surfI));
-
vectorField surfNormal;
- subGeom_[surfI].getNormal
- (
- List<pointIndexHit>
- (
- UIndirectList<pointIndexHit>(info, indices)
- ),
- surfNormal
- );
- forAll(indices, i)
+ subGeom_[surfI].getNormal(surfInfo[surfI], surfNormal);
+
+ const labelList& map = infoMap[surfI];
+ forAll(map, i)
{
- normal[indices[i]] = surfNormal[i];
+ normal[map[i]] = surfNormal[i];
}
}
}
@@ -561,4 +605,99 @@ void Foam::searchableSurfaceCollection::getVolumeType
}
+void Foam::searchableSurfaceCollection::distribute
+(
+ const List<treeBoundBox>& bbs,
+ const bool keepNonLocal,
+ autoPtr<mapDistribute>& faceMap,
+ autoPtr<mapDistribute>& pointMap
+)
+{
+ forAll(subGeom_, surfI)
+ {
+ // Note:Tranform the bounding boxes? Something like
+ // pointField bbPoints =
+ // cmptDivide
+ // (
+ // transform_[surfI].localPosition
+ // (
+ // bbs[i].points()
+ // ),
+ // scale_[surfI]
+ // );
+ // treeBoundBox newBb(bbPoints);
+
+ // Note: what to do with faceMap, pointMap from multiple surfaces?
+ subGeom_[surfI].distribute
+ (
+ bbs,
+ keepNonLocal,
+ faceMap,
+ pointMap
+ );
+ }
+}
+
+
+void Foam::searchableSurfaceCollection::setField(const labelList& values)
+{
+ forAll(subGeom_, surfI)
+ {
+ subGeom_[surfI].setField
+ (
+ static_cast<const labelList&>
+ (
+ SubList<label>
+ (
+ values,
+ subGeom_[surfI].size(),
+ indexOffset_[surfI]
+ )
+ )
+ );
+ }
+}
+
+
+void Foam::searchableSurfaceCollection::getField
+(
+ const List<pointIndexHit>& info,
+ labelList& values
+) const
+{
+ if (subGeom_.size() == 0)
+ {}
+ else if (subGeom_.size() == 1)
+ {
+ subGeom_[0].getField(info, values);
+ }
+ else
+ {
+ // Multiple surfaces. Sort by surface.
+
+ // Per surface the hit
+ List<List<pointIndexHit> > surfInfo;
+ // Per surface the original position
+ List<List<label> > infoMap;
+ sortHits(info, surfInfo, infoMap);
+
+ values.setSize(info.size());
+ //?Misses do not get set? values = 0;
+
+ // Do surface tests
+ forAll(surfInfo, surfI)
+ {
+ labelList surfValues;
+ subGeom_[surfI].getField(surfInfo[surfI], surfValues);
+
+ const labelList& map = infoMap[surfI];
+ forAll(map, i)
+ {
+ values[map[i]] = surfValues[i];
+ }
+ }
+ }
+}
+
+
// ************************************************************************* //
diff --git a/src/meshTools/searchableSurface/searchableSurfaceCollection.H b/src/meshTools/searchableSurface/searchableSurfaceCollection.H
index 3f7ee84baa4012d2482007fa75434a8d1dbc4ca2..e0de60827831f51f04c4cdb232fe5c18c9659d02 100644
--- a/src/meshTools/searchableSurface/searchableSurfaceCollection.H
+++ b/src/meshTools/searchableSurface/searchableSurfaceCollection.H
@@ -77,6 +77,10 @@ private:
Switch mergeSubRegions_;
+ //- offsets for indices coming from different surfaces
+ // (sized with size() of each surface)
+ labelList indexOffset_;
+
//- Region names
mutable wordList regions_;
//- From individual regions to collection regions
@@ -95,6 +99,15 @@ private:
labelList& nearestSurf
) const;
+ //- Sort hits into per-surface bins. Misses are rejected.
+ // Maintains map back to position
+ void sortHits
+ (
+ const List<pointIndexHit>& info,
+ List<List<pointIndexHit> >& surfInfo,
+ labelListList& infoMap
+ ) const;
+
//- Disallow default bitwise copy construct
searchableSurfaceCollection(const searchableSurfaceCollection&);
@@ -161,6 +174,10 @@ public:
//- Range of local indices that can be returned.
virtual label size() const;
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const;
+
// Multiple point queries.
@@ -215,6 +232,27 @@ public:
List<volumeType>&
) const;
+ // Other
+
+ //- Set bounds of surface. Bounds currently set as list of
+ // bounding boxes. The bounds are hints to the surface as for
+ // the range of queries it can expect. faceMap/pointMap can be
+ // set if the surface has done any redistribution.
+ virtual void distribute
+ (
+ const List<treeBoundBox>&,
+ const bool keepNonLocal,
+ autoPtr<mapDistribute>& faceMap,
+ autoPtr<mapDistribute>& pointMap
+ );
+
+ //- WIP. Store element-wise field.
+ virtual void setField(const labelList& values);
+
+ //- WIP. From a set of hits (points and
+ // indices) get the specified field. Misses do not get set. Return
+ // empty field if not supported.
+ virtual void getField(const List<pointIndexHit>&, labelList&) const;
// regIOobject implementation
diff --git a/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H b/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H
index 312c74b1d77b66ce6cb4bff37b6ccbac15068b2e..64db7ee91ceabb0343e6be0d68c6f40814fc96d2 100644
--- a/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H
+++ b/src/meshTools/searchableSurface/searchableSurfaceWithGaps.H
@@ -151,6 +151,13 @@ public:
return surface().size();
}
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const
+ {
+ return surface().coordinates();
+ }
+
// Multiple point queries.
@@ -225,6 +232,41 @@ public:
}
+ // Other
+
+ //- Set bounds of surface. Bounds currently set as list of
+ // bounding boxes. The bounds are hints to the surface as for
+ // the range of queries it can expect. faceMap/pointMap can be
+ // set if the surface has done any redistribution.
+ virtual void distribute
+ (
+ const List<treeBoundBox>& bbs,
+ const bool keepNonLocal,
+ autoPtr<mapDistribute>& faceMap,
+ autoPtr<mapDistribute>& pointMap
+ )
+ {
+ subGeom_[0].distribute(bbs, keepNonLocal, faceMap, pointMap);
+ }
+
+ //- WIP. Store element-wise field.
+ virtual void setField(const labelList& values)
+ {
+ subGeom_[0].setField(values);
+ }
+
+ //- WIP. From a set of hits (points and
+ // indices) get the specified field. Misses do not get set. Return
+ // empty field if not supported.
+ virtual void getField
+ (
+ const List<pointIndexHit>& info,
+ labelList& values
+ ) const
+ {
+ surface().getField(info, values);
+ }
+
// regIOobject implementation
bool writeData(Ostream& os) const
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.C b/src/meshTools/searchableSurface/triSurfaceMesh.C
index 11d687aa8b6f75d86f4b191d89d9f681f8382ea9..08d367e9746252451c5ddd9cb95e3c663235e3b9 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.C
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.C
@@ -30,6 +30,7 @@ License
#include "EdgeMap.H"
#include "triSurfaceFields.H"
#include "Time.H"
+#include "PackedBoolList.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -281,6 +282,36 @@ void Foam::triSurfaceMesh::getNextIntersections
}
+void Foam::triSurfaceMesh::calcBounds(boundBox& bb, label& nPoints) const
+{
+ // Unfortunately nPoints constructs meshPoints() so do compact version
+ // ourselves
+
+ const triSurface& s = static_cast<const triSurface&>(*this);
+
+ PackedBoolList pointIsUsed(points().size());
+
+ nPoints = 0;
+ bb = boundBox::invertedBox;
+
+ forAll(s, triI)
+ {
+ const labelledTri& f = s[triI];
+
+ forAll(f, fp)
+ {
+ label pointI = f[fp];
+ if (pointIsUsed.set(pointI, 1u))
+ {
+ bb.min() = ::Foam::min(bb.min(), points()[pointI]);
+ bb.max() = ::Foam::max(bb.max(), points()[pointI]);
+ nPoints++;
+ }
+ }
+ }
+}
+
+
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
@@ -301,6 +332,7 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
),
triSurface(s),
tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
+ maxTreeDepth_(10),
surfaceClosed_(-1)
{}
@@ -344,6 +376,7 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io)
)
),
tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
+ maxTreeDepth_(10),
surfaceClosed_(-1)
{}
@@ -390,6 +423,7 @@ Foam::triSurfaceMesh::triSurfaceMesh
)
),
tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
+ maxTreeDepth_(10),
surfaceClosed_(-1)
{
scalar scaleFactor = 0;
@@ -410,6 +444,14 @@ Foam::triSurfaceMesh::triSurfaceMesh
Info<< searchableSurface::name() << " : using intersection tolerance "
<< tolerance_ << endl;
}
+
+
+ // Have optional non-standard tree-depth to limit storage.
+ if (dict.readIfPresent("maxTreeDepth", maxTreeDepth_) && maxTreeDepth_ > 0)
+ {
+ Info<< searchableSurface::name() << " : using maximum tree depth "
+ << maxTreeDepth_ << endl;
+ }
}
@@ -431,6 +473,17 @@ void Foam::triSurfaceMesh::clearOut()
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+Foam::pointField Foam::triSurfaceMesh::coordinates() const
+{
+ // Use copy to calculate face centres so they don't get stored
+ return PrimitivePatch<labelledTri, SubList, const pointField&>
+ (
+ SubList<labelledTri>(*this, triSurface::size()),
+ triSurface::points()
+ ).faceCentres();
+}
+
+
void Foam::triSurfaceMesh::movePoints(const pointField& newPoints)
{
tree_.clear();
@@ -444,15 +497,28 @@ const Foam::indexedOctree<Foam::treeDataTriSurface>&
{
if (tree_.empty())
{
+ // Calculate bb without constructing local point numbering.
+ treeBoundBox bb;
+ label nPoints;
+ calcBounds(bb, nPoints);
+
+ if (nPoints != points().size())
+ {
+ WarningIn("triSurfaceMesh::tree() const")
+ << "Surface " << searchableSurface::name()
+ << " does not have compact point numbering."
+ << " Of " << points().size() << " only " << nPoints
+ << " are used. This might give problems in some routines."
+ << endl;
+ }
+
+
// Random number generator. Bit dodgy since not exactly random ;-)
Random rndGen(65431);
// Slightly extended bb. Slightly off-centred just so on symmetric
// geometry there are less face/edge aligned items.
- treeBoundBox bb
- (
- treeBoundBox(points(), meshPoints()).extend(rndGen, 1E-4)
- );
+ bb = bb.extend(rndGen, 1E-4);
bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
@@ -465,9 +531,9 @@ const Foam::indexedOctree<Foam::treeDataTriSurface>&
(
treeDataTriSurface(*this),
bb,
- 10, // maxLevel
- 10, // leafsize
- 3.0 // duplicity
+ maxTreeDepth_, // maxLevel
+ 10, // leafsize
+ 3.0 // duplicity
)
);
@@ -494,15 +560,17 @@ const Foam::indexedOctree<Foam::treeDataEdge>&
+ nInternalEdges()
);
+ treeBoundBox bb;
+ label nPoints;
+ calcBounds(bb, nPoints);
+
// Random number generator. Bit dodgy since not exactly random ;-)
Random rndGen(65431);
// Slightly extended bb. Slightly off-centred just so on symmetric
// geometry there are less face/edge aligned items.
- treeBoundBox bb
- (
- treeBoundBox(points(), meshPoints()).extend(rndGen, 1E-4)
- );
+
+ bb = bb.extend(rndGen, 1E-4);
bb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
bb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
@@ -517,10 +585,10 @@ const Foam::indexedOctree<Foam::treeDataEdge>&
localPoints(), // points
bEdges // selected edges
),
- bb, // bb
- 8, // maxLevel
- 10, // leafsize
- 3.0 // duplicity
+ bb, // bb
+ maxTreeDepth_, // maxLevel
+ 10, // leafsize
+ 3.0 // duplicity
)
);
}
@@ -754,24 +822,53 @@ void Foam::triSurfaceMesh::getNormal
}
+void Foam::triSurfaceMesh::setField(const labelList& values)
+{
+ autoPtr<triSurfaceLabelField> fldPtr
+ (
+ new triSurfaceLabelField
+ (
+ IOobject
+ (
+ "values",
+ objectRegistry::time().timeName(), // instance
+ "triSurface", // local
+ *this,
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ *this,
+ dimless,
+ labelField(values)
+ )
+ );
+
+ // Store field on triMesh
+ fldPtr.ptr()->store();
+}
+
+
void Foam::triSurfaceMesh::getField
(
- const word& fieldName,
const List<pointIndexHit>& info,
labelList& values
) const
{
- const triSurfaceLabelField& fld = lookupObject<triSurfaceLabelField>
- (
- fieldName
- );
-
- values.setSize(info.size());
- forAll(info, i)
+ if (foundObject<triSurfaceLabelField>("values"))
{
- if (info[i].hit())
+ values.setSize(info.size());
+
+ const triSurfaceLabelField& fld = lookupObject<triSurfaceLabelField>
+ (
+ "values"
+ );
+
+ forAll(info, i)
{
- values[i] = fld[info[i].index()];
+ if (info[i].hit())
+ {
+ values[i] = fld[info[i].index()];
+ }
}
}
}
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.H b/src/meshTools/searchableSurface/triSurfaceMesh.H
index e93b83ecdbb9a04046639c79b666edef979e0c72..dea7c642b896ab21f8670c2e6218b3a491e79c63 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.H
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.H
@@ -71,6 +71,9 @@ private:
//- Optional tolerance to use in searches
scalar tolerance_;
+ //- Optional max tree depth of octree
+ label maxTreeDepth_;
+
//- Search tree (triangles)
mutable autoPtr<indexedOctree<treeDataTriSurface> > tree_;
@@ -129,6 +132,11 @@ private:
void operator=(const triSurfaceMesh&);
+protected:
+
+ //- Calculate (number of)used points and their bounding box
+ void calcBounds(boundBox& bb, label& nPoints) const;
+
public:
//- Runtime type information
@@ -185,6 +193,10 @@ public:
return triSurface::size();
}
+ //- Get representative set of element coordinates
+ // Usually the element centres (should be of length size()).
+ virtual pointField coordinates() const;
+
virtual void findNearest
(
const pointField& sample,
@@ -236,6 +248,8 @@ public:
List<volumeType>&
) const;
+ // Other
+
//- Set bounds of surface. Bounds currently set as list of
// bounding boxes. The bounds are hints to the surface as for
// the range of queries it can expect. faceMap/pointMap can be
@@ -249,16 +263,12 @@ public:
)
{}
- // Other
+ //- WIP. Store element-wise field.
+ virtual void setField(const labelList& values);
- //- Specific to triSurfaceMesh: from a set of hits (points and
+ //- WIP. From a set of hits (points and
// indices) get the specified field. Misses do not get set.
- virtual void getField
- (
- const word& fieldName,
- const List<pointIndexHit>&,
- labelList& values
- ) const;
+ virtual void getField(const List<pointIndexHit>&, labelList&) const;
// regIOobject implementation