78 void simpleMarkFeatures
82 const scalar featureAngle,
83 const bool concaveMultiCells,
84 const bool doNotPreserveFaceZones,
92 const scalar minCos =
Foam::cos(featureAngle);
113 label meshEdgeI = meshEdges[edgeI];
114 featureEdgeSet.
insert(meshEdgeI);
115 singleCellFeaturePointSet.
insert(mesh.
edges()[meshEdgeI][0]);
116 singleCellFeaturePointSet.
insert(mesh.
edges()[meshEdgeI][1]);
141 allBoundary.checkPointManifold(
false, &singleCellFeaturePointSet);
145 const labelList& meshPoints = allBoundary.meshPoints();
149 const labelList& eFaces = edgeFaces[edgeI];
151 if (eFaces.
size() > 2)
153 const edge&
e = allBoundary.edges()[edgeI];
160 singleCellFeaturePointSet.
insert(meshPoints[
e[0]]);
161 singleCellFeaturePointSet.
insert(meshPoints[
e[1]]);
168 const labelList& eFaces = edgeFaces[edgeI];
170 if (eFaces.
size() == 2)
172 label f0 = eFaces[0];
173 label f1 = eFaces[1];
176 const vector& n0 = allBoundary.faceNormals()[f0];
177 const vector& n1 = allBoundary.faceNormals()[f1];
179 if ((n0 & n1) < minCos)
181 const edge&
e = allBoundary.edges()[edgeI];
182 label v0 = meshPoints[
e[0]];
183 label v1 = meshPoints[
e[1]];
186 featureEdgeSet.
insert(meshEdgeI);
192 allBoundary[f1].centre(allBoundary.points())
193 - allBoundary[f0].centre(allBoundary.points())
196 if (concaveMultiCells && (c1c0 & n0) > small)
199 Info<<
"Detected concave feature edge:" << edgeI
200 <<
" cos:" << (c1c0 & n0)
202 << allBoundary.points()[v0]
203 << allBoundary.points()[v1]
206 singleCellFeaturePointSet.
erase(v0);
207 multiCellFeaturePointSet.
insert(v0);
208 singleCellFeaturePointSet.
erase(v1);
209 multiCellFeaturePointSet.
insert(v1);
214 if (!multiCellFeaturePointSet.
found(v0))
216 singleCellFeaturePointSet.
insert(v0);
218 if (!multiCellFeaturePointSet.
found(v1))
220 singleCellFeaturePointSet.
insert(v1);
236 featureFaceSet.insert(facei);
242 if (doNotPreserveFaceZones)
244 if (faceZones.
size() > 0)
247 <<
"Detected " << faceZones.
size()
248 <<
" faceZones. These will not be preserved."
254 if (faceZones.
size() > 0)
256 Info<<
"Detected " << faceZones.
size()
257 <<
" faceZones. Preserving these by marking their"
258 <<
" points, edges and faces as features." <<
endl;
263 const faceZone& fz = faceZones[zoneI];
265 Info<<
"Inserting all faces in faceZone " << fz.
name()
266 <<
" as features." <<
endl;
274 featureFaceSet.insert(facei);
279 singleCellFeaturePointSet.
erase(
f[fp]);
280 multiCellFeaturePointSet.
insert(
f[fp]);
283 featureEdgeSet.
insert(fEdges[fp]);
290 featureFaces = featureFaceSet.toc();
291 featureEdges = featureEdgeSet.
toc();
292 singleCellFeaturePoints = singleCellFeaturePointSet.
toc();
293 multiCellFeaturePoints = multiCellFeaturePointSet.
toc();
303 const labelList& singleCellFeaturePoints,
309 Info<<
"Dumping centres of featureFaces to obj file " << str.name()
318 Info<<
"Dumping featureEdges to obj file " << str.name() <<
endl;
323 const edge&
e = mesh.
edges()[featureEdges[i]];
328 str<<
"l " << vertI-1 <<
' ' << vertI <<
nl;
332 OFstream str(
"singleCellFeaturePoints.obj");
333 Info<<
"Dumping featurePoints that become a single cell to obj file "
334 << str.name() <<
endl;
335 forAll(singleCellFeaturePoints, i)
341 OFstream str(
"multiCellFeaturePoints.obj");
342 Info<<
"Dumping featurePoints that become multiple cells to obj file "
343 << str.name() <<
endl;
344 forAll(multiCellFeaturePoints, i)
352 int main(
int argc,
char *argv[])
361 "have multiple faces in between cells"
366 "split cells on concave boundary edges into multiple cells"
370 "doNotPreserveFaceZones",
371 "disable the default behaviour of preserving faceZones by having"
372 " multiple faces in between cells"
391 isBoundaryEdge.
set(fEdges[i], 1);
396 const scalar minCos =
Foam::cos(featureAngle);
399 <<
"minCos :" << minCos <<
endl
406 Info<<
"Splitting all internal faces to create multiple faces"
407 <<
" between two cells." <<
nl
414 "doNotPreserveFaceZones"
417 if (concaveMultiCells)
419 Info<<
"Generating multiple cells for points on concave feature edges."
427 if (doNotPreserveFaceZones)
448 doNotPreserveFaceZones,
452 singleCellFeaturePoints,
453 multiCellFeaturePoints
474 singleCellFeaturePoints,
475 multiCellFeaturePoints
487 dualMaker.setRefinement
492 singleCellFeaturePoints,
493 multiCellFeaturePoints,
512 Info<<
"Writing dual mesh to " << runTime.name() <<
endl;
#define forAll(list, i)
Loop across all elements in list.
bool insert(const Key &key)
Insert a new entry.
bool erase(const iterator &)
Erase a hashedEntry specified by given iterator.
List< Key > toc() const
Return the table of contents.
bool found(const Key &) const
Return true if hashedEntry is found in table.
void size(const label)
Override size to be inconsistent with allocated storage.
void set(const PackedList< 1 > &)
Set specified bits.
label nEdges() const
Return number of edges in patch.
label nInternalEdges() const
Number of internal edges.
A List obtained as a section of another List.
label size() const
Return the number of elements in the UPtrList.
void clear()
Clear the zones.
const word & name() const
Return name.
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
bool optionFound(const word &opt) const
Return true if the named option is found.
static void noParallel()
Remove the parallel options.
static SLList< string > validArgs
A list of valid (mandatory) arguments.
T argRead(const label index) const
Read a value from the argument at index.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
A subset of mesh faces organised as a primitive patch.
A face is a list of labels corresponding to mesh vertices.
Creates dual of polyMesh. Every point becomes a cell (or multiple cells for feature points),...
Mesh consisting of general polyhedral cells.
const pointZoneList & pointZones() const
Return point zones.
const cellZoneList & cellZones() const
Return cell zones.
virtual const faceList & faces() const
Return raw faces.
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
const fileName & pointsInstance() const
Return the current instance directory for points.
const faceZoneList & faceZones() const
Return face zones.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const pointField & points() const
Return raw points.
void setInstance(const fileName &)
Set the instance for mesh files.
A patch is a list of labels that address the faces in the global face list.
const labelList & meshEdges() const
Return global edge index for local edges.
Direct mesh changes based on v1.3 polyTopoChange syntax.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
const vectorField & faceCentres() const
label nInternalFaces() const
const labelListList & faceEdges() const
virtual bool write(const bool write=true) const
Write using setting from DB.
A class for handling words, derived from string.
int main(int argc, char *argv[])
const fvPatchList & patches
#define WarningInFunction
Report a warning using Foam::Warning.
scalar radToDeg(const scalar rad)
Convert radians to degrees.
scalar degToRad(const scalar deg)
Convert degrees to radians.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
labelList identityMap(const label len)
Create identity map (map[i] == i) of given length.
dimensionedScalar cos(const dimensionedScalar &ds)
Foam::argList args(argc, argv)