33 void Foam::refinementFeatures::read
35 const objectRegistry& io,
36 const PtrList<dictionary>& featDicts
41 const dictionary& dict = featDicts[featI];
43 fileName featFileName(dict.lookup(
"file"));
52 "extendedFeatureEdgeMesh",
59 const fileName fName(extFeatObj.filePath());
68 Info<<
"Read extendedFeatureEdgeMesh " << extFeatObj.name()
70 eMeshPtr().writeStats(
Info);
73 set(featI,
new extendedFeatureEdgeMesh(extFeatObj, eMeshPtr()));
90 const fileName fName(featObj.filePath());
96 "refinementFeatures::read" 97 "(const objectRegistry&" 98 ", const PtrList<dictionary>&)",
100 ) <<
"Could not open " << featObj.objectPath()
107 const edgeMesh& eMesh = eMeshPtr();
109 Info<<
"Read edgeMesh " << featObj.name() <<
nl 111 eMesh.writeStats(
Info);
119 labelList oldToNew(eMesh.points().size(), -1);
120 DynamicField<point> newPoints(eMesh.points().size());
121 forAll(pointEdges, pointI)
123 if (pointEdges[pointI].
size() > 2)
125 oldToNew[pointI] = newPoints.size();
126 newPoints.append(eMesh.points()[pointI]);
131 label nFeatures = newPoints.size();
134 if (oldToNew[pointI] == -1)
136 oldToNew[pointI] = newPoints.size();
137 newPoints.append(eMesh.points()[pointI]);
142 const edgeList& edges = eMesh.edges();
146 const edge&
e = edges[edgeI];
147 newEdges[edgeI] = edge
158 extendedEdgeMesh eeMesh
170 List<extendedEdgeMesh::sideVolumeType>(0),
184 set(featI,
new extendedFeatureEdgeMesh(featObj, eeMesh));
187 const extendedEdgeMesh& eMesh =
operator[](featI);
189 if (dict.found(
"levels"))
191 List<Tuple2<scalar, label> > distLevels(dict[
"levels"]);
197 "refinementFeatures::read" 198 "(const objectRegistry&" 199 ", const PtrList<dictionary>&)" 200 ) <<
" : levels should be at least size 1" << endl
201 <<
"levels : " << dict[
"levels"]
205 distances_[featI].
setSize(distLevels.size());
206 levels_[featI].
setSize(distLevels.size());
210 distances_[featI][j] = distLevels[j].
first();
211 levels_[featI][j] = distLevels[j].second();
218 (distances_[featI][j] <= distances_[featI][j-1])
219 || (levels_[featI][j] > levels_[featI][j-1])
224 "refinementFeatures::read" 225 "(const objectRegistry&" 226 ", const PtrList<dictionary>&)" 227 ) <<
" : Refinement should be specified in order" 228 <<
" of increasing distance" 229 <<
" (and decreasing refinement level)." << endl
230 <<
"Distance:" << distances_[featI][j]
231 <<
" refinementLevel:" << levels_[featI][j]
244 Info<<
"Refinement level according to distance to " 245 << featFileName <<
" (" << eMesh.points().size() <<
" points, " 246 << eMesh.edges().size() <<
" edges)." <<
endl;
249 Info<<
" level " << levels_[featI][j]
250 <<
" for all cells within " << distances_[featI][j]
251 <<
" metre." <<
endl;
257 void Foam::refinementFeatures::buildTrees(
const label featI)
259 const extendedEdgeMesh& eMesh =
operator[](featI);
261 const edgeList& edges = eMesh.edges();
264 treeBoundBox bb(points);
272 bb.min() -=
point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
273 bb.max() +=
point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
278 new indexedOctree<treeDataEdge>
300 new indexedOctree<treeDataPoint>
302 treeDataPoint(points, featurePoints),
315 if (!regionEdgeTreesPtr_.valid())
317 regionEdgeTreesPtr_.reset
338 bb.
min() -=
point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
339 bb.
max() +=
point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
361 return regionEdgeTreesPtr_();
366 void Foam::refinementFeatures::findHigherLevel
384 label candidateI = 0;
390 if (levels[levelI] > maxLevel[pointI])
392 candidates[candidateI] = pt[pointI];
393 candidateMap[candidateI] = pointI;
394 candidateDistSqr[candidateI] =
sqr(distances[levelI]);
400 candidates.setSize(candidateI);
401 candidateMap.setSize(candidateI);
402 candidateDistSqr.setSize(candidateI);
408 forAll(candidates, candidateI)
410 nearInfo[candidateI] = tree.findNearest
412 candidates[candidateI],
413 candidateDistSqr[candidateI]
418 forAll(nearInfo, candidateI)
420 if (nearInfo[candidateI].hit())
426 mag(nearInfo[candidateI].hitPoint()-candidates[candidateI])
429 label pointI = candidateMap[candidateI];
432 maxLevel[pointI] = levels[minDistI+1];
447 distances_(featDicts.
size()),
448 levels_(featDicts.
size()),
449 edgeTrees_(featDicts.
size()),
450 pointTrees_(featDicts.
size())
555 if (tree.
shapes().size() > 0)
559 const point& sample = samples[sampleI];
562 if (nearInfo[sampleI].hit())
564 distSqr =
magSqr(nearInfo[sampleI].hitPoint()-sample);
568 distSqr = nearestDistSqr[sampleI];
575 nearFeature[sampleI] = featI;
584 const edge&
e = td.edges()[nearInfo[sampleI].index()];
585 nearNormal[sampleI] = e.
vec(td.points());
586 nearNormal[sampleI] /=
mag(nearNormal[sampleI])+VSMALL;
614 forAll(regionTrees, featI)
620 const point& sample = samples[sampleI];
623 if (nearInfo[sampleI].hit())
625 distSqr =
magSqr(nearInfo[sampleI].hitPoint()-sample);
629 distSqr = nearestDistSqr[sampleI];
633 pointIndexHit info = regionTree.findNearest(sample, distSqr);
639 nearFeature[sampleI] = featI;
647 const edge&
e = td.
edges()[nearInfo[sampleI].index()];
648 nearNormal[sampleI] = e.vec(td.
points());
649 nearNormal[sampleI] /=
mag(nearNormal[sampleI])+VSMALL;
722 forAll(pointTrees_, featI)
726 if (tree.
shapes().pointLabels().size() > 0)
730 const point& sample = samples[sampleI];
733 if (nearFeature[sampleI] != -1)
735 distSqr =
magSqr(nearInfo[sampleI].hitPoint()-sample);
739 distSqr = nearestDistSqr[sampleI];
746 nearFeature[sampleI] = featI;
760 void Foam::refinementFeatures::findHigherLevel
772 findHigherLevel(pt, featI, maxLevel);
779 scalar overallMax = -GREAT;
782 overallMax =
max(overallMax,
max(distances_[featI]));
const PtrList< indexedOctree< treeDataEdge > > & regionEdgeTrees() const
Simple random number generator.
cachedRandom rndGen(label(0),-1)
const point & min() const
Minimum describing the bounding box.
void findNearestPoint(const pointField &samples, const scalarField &nearestDistSqr, labelList &nearFeature, List< pointIndexHit > &nearInfo) const
Find nearest feature point. Sets.
label size() const
Return the number of elements in the PtrList.
bool set(const label) const
Is element set.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
vector point
Point is a vector.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Tuple2< scalar, label > nearInfo
Private class for finding nearest.
dimensioned< scalar > mag(const dimensioned< Type > &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const Type & shapes() const
Reference to shape.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const pointField & points() const
Return points.
Standard boundBox + extra functionality for use in octree.
void size(const label)
Override size to be inconsistent with allocated storage.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const edgeList & edges() const
PointIndexHit< point > pointIndexHit
static autoPtr< extendedEdgeMesh > New(const fileName &, const word &ext)
Select constructed from filename (explicit extension)
bool hit() const
Is there a hit.
const labelList & regionEdges() const
Return the feature edges which are on the boundary between.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
T & first()
Return the first element of the list.
label index() const
Return index.
vectorField pointField
pointField is a vectorField.
label readLabel(Istream &is)
const Point & hitPoint() const
Return hit point.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static bool canRead(const fileName &, const bool verbose=false)
Can we read this file format?
static autoPtr< edgeMesh > New(const fileName &, const word &ext)
Select constructed from filename (explicit extension)
Ostream & incrIndent(Ostream &os)
Increment the indent level.
const double e
Elementary charge.
void setSize(const label)
Reset size of List.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
treeBoundBox extend(Random &, const scalar s) const
Return slightly wider bounding box.
const pointField & points() const
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
refinementFeatures(const objectRegistry &io, const PtrList< dictionary > &featDicts)
Construct from description.
const List< scalarField > & distances() const
Per featureEdgeMesh the list of ranges.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
const point & max() const
Maximum describing the bounding box.
label findLower(const ListType &, typename ListType::const_reference, const label stary, const BinaryOp &bop)
Find last element < given value in sorted list and return index,.
vector vec(const pointField &) const
Return the vector (end - start)
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Non-pointer based hierarchical recursive searching.
Holds data for octree to work on an edges subset.
void findNearestRegionEdge(const pointField &samples, const scalarField &nearestDistSqr, labelList &nearFeature, List< pointIndexHit > &nearInfo, vectorField &nearNormal) const
Find nearest point on nearest region edge. Sets.
#define forAllReverse(list, i)
Registry of regIOobjects.
List< label > labelList
A List of labels.
scalar maxDistance() const
Highest distance of all features.
const labelListList & levels() const
Per featureEdgeMesh the list of level.
Description of feature edges and points.
void findNearestEdge(const pointField &samples, const scalarField &nearestDistSqr, labelList &nearFeature, List< pointIndexHit > &nearInfo, vectorField &nearNormal) const
Find nearest point on nearest feature edge. Sets.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
List< labelList > labelListList
A List of labelList.
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
const extendedFeatureEdgeMesh & operator[](const label) const
Return element const reference.
const edgeList & edges() const
Return edges.