libsurf
Programmer's Documentation

pentagrow.h (r6227/r5385)
1 
2 /* Copyright (C) 2015 David Eller <david@larosterna.com>
3  *
4  * Commercial License Usage
5  * Licensees holding valid commercial licenses may use this file in accordance
6  * with the terms contained in their respective non-exclusive license agreement.
7  * For further information contact david@larosterna.com .
8  *
9  * GNU General Public License Usage
10  * Alternatively, this file may be used under the terms of the GNU General
11  * Public License version 3.0 as published by the Free Software Foundation and
12  * appearing in the file gpl.txt included in the packaging of this file.
13  */
14 
15 #ifndef SURF_PENTAGROW_H
16 #define SURF_PENTAGROW_H
17 
18 #include <genua/trimesh.h>
19 #include <genua/mxmesh.h>
20 #include <genua/ndpointtree.h>
21 #include <genua/logger.h>
22 #include <boost/atomic.hpp>
23 
44 class PentaGrow : public MxMesh, public Logger
45 {
46 public:
47 
49  PentaGrow() {}
50 
52  PentaGrow(const TriMesh & m);
53 
55  void configure(const ConfigParser &cfg);
56 
58  static int maximumTagValue() {
59  return std::numeric_limits<int>::max();
60  }
61 
63  static uint maximumTriangleCount() {
64  return std::numeric_limits<int>::max();
65  }
66 
68  void generateShell(int hiter, int niter, int ncrititer,
69  int laplaceiter);
70 
72  size_t nWallNodes() const {return mwall.nvertices();}
73 
75  uint extrude(bool curvedGrowth);
76 
78  void adaptWall(const DVector<uint> &faceTags);
79 
81  void readTets(const std::string &basename);
82 
84  void writeTetgen(const std::string &fname, const TriMesh & farf,
85  const PointList<3> & holes, const TriMesh &refr = TriMesh(),
86  Real nearBoxEdge = 0.0);
87 
89  void shrink();
90 
92  size_t countNegativeVolumes(std::ostream &msg);
93 
95  void envelopeBounds(Vct3 &plo, Vct3 &phi) const;
96 
98  bool ellipsoidEncloses(const Vct3 &ctr, const Vct3 &hax) const;
99 
101  void envelopeEdgeStats(Real &lmean, Real &lmax) const;
102 
104  Vector prismQualitySumCos(const std::string &fname,
105  uint isection, uint nbin=NotFound) const;
106 
107 #ifdef HAVE_NLOPT
108 
110  void initializeBounds(double *x, double *lbound, double *ubound);
111 
113  double inversionConstraint(const double *x, double *grad) const;
114 
116  double intersectionConstraint(const double *x, double *grad) const;
117 
119  double qualityObjective(const double *x, double *grad) const;
120 
122  void optimizeEnvelope();
123 
124 #endif // HAVE_NLOPT
125 
126  // write outermost layer to file (debugging)
127  void writeShell(const std::string & fname);
128 
129  // test connectivity
130  void debugConnect();
131 
132 private:
133 
134  enum VertexCategory { Undefined = 0,
135  Concave = 1,
136  Convex = 2,
137  Conical = 4,
138  Corner = 8,
139  Ridge = 16,
140  StrongCurvature = 32,
141  Sharp = 64,
142  Saddle = (Concave | Convex), // == 3
143  ConcaveCorner = (Corner | Concave), // == 9
144  ConvexCorner = (Corner | Convex), // == 10
145  SaddleCorner = (Corner | Saddle), // == 11
146  ConeDipp = (Concave | Corner | Conical), // == 13
147  ConeTip = (Convex | Corner | Conical), // == 14
148  BluntCorner = (Corner|Saddle|Conical), // == 15
149  Trench = (Ridge | Concave), // == 17
150  ConvexEdge = (Ridge | Convex), // == 18
151  Wedge = (Ridge | Convex | Concave), // == 19
152  RidgeConeTip = (Convex | Conical | Ridge), // == 22
153  LeadingEdgeIntersection = (Trench | StrongCurvature), // == 49
154  TrailingEdgeIntersection = (SaddleCorner | Sharp), // == 83
155  CriticalCorner = 512,
156  Flat = 1024,
157  Anything = -1};
158 
160  void extractWall(const MxMesh & gm);
161 
163  void classify();
164 
166  void adjustRidgeNormals();
167 
169  bool isClass(size_t i, int cat) const {
170  return ( i < vtype.size() ) and ( vtype[i] == cat );
171  }
172 
174  bool hasClass(size_t i, int cat) const {
175  return ( i < vtype.size() ) and ( (vtype[i] & cat) == cat );
176  }
177 
179  Real convexity(const Vct3 &p1, const Vct3 &n1,
180  const Vct3 &p2, const Vct3 &n2)
181  {
182  return arg( p2-p1, n1 ) + arg( p1-p2, n2 ) - PI;
183  }
184 
186  Real convexity(uint i1, uint i2)
187  {
188  return convexity( mwall.vertex(i1), wfn[i1],
189  mwall.vertex(i2), wfn[i2] );
190  }
191 
193  void smoothThickness(Vector &lyt, int niter) const;
194 
196  void edgeLengthStats(uint k, Real &lmean, Real &lmax, Real &lmin) const;
197 
199  void prismPattern(Real rhfirst, Real rhlast, Vector &xpp) const;
200 
202  Vct3 projectToWall(const TriMesh & mout,
203  const Vct3 & pout, uint inear) const;
204 
206  Vct3 findWallVertex(const TriMesh &oldShell,
207  const TriMesh &newShell, uint niShell) const;
208 
210  void rebuildTree();
211 
213  void findNeighbors(const Vct3 &p, Real r, Indices &neighbors) const {
214  assert(nodeTree.npoints() == vout.size());
215  neighbors.clear();
216  nodeTree.find( Vct3f(p), r, neighbors );
217  }
218 
220  bool collisions(Indices &colliding, uint iwall, Real safety,
221  Real nrmdev=cos(rad(60.)), Real fnrmdev=cos(rad(120.))) const;
222 
224  int collisions(uint iwall, Real safety = 1.2, Real nrmdev = cos(rad(60.)),
225  Real fnrmdev = cos(rad(120.))) const;
226 
228  uint32_t extractSectionTag(uint32_t tag) const {
229  return id2section[tag];
230  }
231 
233  uint32_t extractElementTag(uint32_t tag) const {
234  return id2index[tag];
235  }
236 
238  void edgeMap(ConnectMap &map) const;
239 
241  template <class Container>
242  void smooth(const ConnectMap &map, Container &c) const
243  {
244  typedef typename Container::value_type ItemType;
245  Container b(c);
246  const size_t n = map.size();
247  ConnectMap::const_iterator itr, last;
248 
249 //#pragma omp parallel for
250  for (size_t i=0; i<n; ++i) {
251  ItemType sum;
252  sum = 0.0;
253  int nsum = 0;
254  last = map.end(i);
255  for (itr = map.begin(i); itr != last; ++itr) {
256  sum += c[*itr];
257  ++nsum;
258  }
259  if (nsum > 0)
260  b[i] = 0.5*c[i] + (0.5/nsum) * sum;
261  }
262 
263  c.swap(b);
264  }
265 
267  template <class Container, class WritePred, class ReadPred>
268  void smooth(const ConnectMap &map, Container &c,
269  const WritePred &writeNode,
270  const ReadPred &readNode) const
271  {
272  typedef typename Container::value_type ItemType;
273  Container b(c);
274  const size_t n = map.size();
275  ConnectMap::const_iterator itr, last;
276 
277 //#pragma omp parallel for
278  for (size_t i=0; i<n; ++i) {
279  if (not writeNode(i))
280  continue;
281  ItemType sum;
282  sum = 0.0;
283  int nsum = 0;
284  last = map.end(i);
285  for (itr = map.begin(i); itr != last; ++itr) {
286  if (readNode(*itr)) {
287  sum += c[*itr];
288  ++nsum;
289  }
290  }
291  if (nsum > 0)
292  b[i] = 0.5*c[i] + (0.5/nsum) * sum;
293  }
294 
295  c.swap(b);
296  }
297 
299  void smoothShellNodes(const Vector &lyt, int niter, Real omega=0.5);
300 
302  Vct3 nbBarycenter(const PointList<3> &pts, size_t k) const;
303 
305  int untangle(Vector &lyt, int niter, Real permitted_etwist);
306 
308  void unwarp(int niter, Real permitted_angle);
309 
311  void uncollide(int niter, Real safety,
312  Real retraction, Real limitphi, Real limitphif);
313 
315  uint uncollideVertex(uint i, Vector &lyt, Real safety,
316  Real retraction, Real cphi, Real cphif);
317 
319  void retractNeighbors(const Indices &afv, Vector &lyt, int ring=3);
320 
322  void extrudeVertex(int ivx, int nl, Real hi,
323  bool curvedGrowth, PointGrid<3> &grid) const;
324 
326  void smoothWallTransition(int niter);
327 
329  void mergeNeighbors(Indices &idx) const;
330 
332  size_t untangleGrid(PointGrid<3> &grid);
333 
335  uint appendPrismLayer(const PointGrid<3> &grid);
336 
338  void updateShellNormals();
339 
341  void centerGridNodes(uint niter, PointGrid<3> &grid) const;
342 
344  void centerGridNodesPass(const PointGrid<3> &cgrid, PointGrid<3> &grid) const;
345 
347  void findEnvelopeNeighbors(Indices &interfaceNodes,
348  Indices &nearTetNodes) const;
349 
350 private:
351 
354  public:
355  UncollideTask(PentaGrow &pg, Vector &lyt,
356  Real safety, Real retract, Real cphi, Real cphif)
357  : m_pg(pg), m_lyt(lyt), m_safety(safety), m_retract(retract),
358  m_cphi(cphi), m_cphif(cphif) { m_ncol.store(0); }
359  void operator() (uint a, uint b) {
360  int sum = 0;
361  for (uint i=a; i<b; ++i) {
362  uint n = m_pg.uncollideVertex(i, m_lyt, m_safety, m_retract, m_cphi,
363  m_cphif);
364  sum += n;
365  if ( n > 0 )
366  m_afv.push_back(i);
367  }
368  m_ncol += sum; // atomic
369  }
370  const Indices &affected() const{return m_afv;}
371  uint ncollisions() const {return m_ncol.load();}
372  private:
373  PentaGrow &m_pg;
374  Vector &m_lyt;
375  Real m_safety, m_retract, m_cphi, m_cphif;
376  boost::atomic<int> m_ncol;
377  Indices m_afv;
378  };
379 
382  public:
383  ExtrusionTask(const PentaGrow &pg, PointGrid<3> &grid,
384  uint nl, Real hi, bool curved)
385  : m_pg(pg), m_grid(grid), m_hi(hi), m_nl(nl), m_curved(curved) {}
386  void operator() (uint a, uint b) {
387  for (uint i=a; i<b; ++i)
388  m_pg.extrudeVertex(i, m_nl, m_hi, m_curved, m_grid);
389  }
390  private:
391  const PentaGrow &m_pg;
392  PointGrid<3> &m_grid;
393  Real m_hi;
394  uint m_nl;
395  bool m_curved;
396  };
397 
398 private:
399 
402 
405 
408 
411 
414 
416  Vector targetHeight;
417 
420 
423 
426 
429 
432 
434  Indices wallTags, farTags;
435 
438 
440  std::vector<bool> gridBaseTangled;
441 
444 
447 
450 
453 
455  Real firstCellHeight, maxRelHeight, maxAbsHeight, maxExpansionFactor;
456 
459 
462 
465 
468 
471 };
472 
473 #endif // PENTAGROW_H
void smoothWallTransition(int niter)
optionally initialize, then distribute wall normal transition parameters
Definition: pentagrow.cpp:1441
PointList< 3 > envNormals
vertex normals for outer layer
Definition: pentagrow.h:410
Hybrid prismatic mesh generation.
Definition: pentagrow.h:44
uint size() const
void findNeighbors(const Vct3 &p, Real r, Indices &neighbors) const
find all outer-layer nodes closer than r to p
Definition: pentagrow.h:213
Real cosSharpAngle
cosine of angle limit for classification as sharp (default 150deg)
Definition: pentagrow.h:449
bool ellipsoidEncloses(const Vct3 &ctr, const Vct3 &hax) const
check whether axis-aligned ellipsoid encloses all envelope vertices
Definition: pentagrow.cpp:2253
size_t countNegativeVolumes(std::ostream &msg)
check volume elements in final mesh for positive volume
Definition: pentagrow.cpp:2140
void prismPattern(Real rhfirst, Real rhlast, Vector &xpp) const
determine suitable normalized pattern for prism heights
Definition: pentagrow.cpp:1575
Vct3 findWallVertex(const TriMesh &oldShell, const TriMesh &newShell, uint niShell) const
determine corresponding wall mesh vertex
Definition: pentagrow.cpp:1752
void readTets(const std::string &basename)
read tetgen result and collect face tags
Definition: pentagrow.cpp:2011
std::vector< bool > gridBaseTangled
tags surface nodes which resulted in tangled grid nodes
Definition: pentagrow.h:440
PentaGrow()
empty object
Definition: pentagrow.h:49
Vector invGrowthExponent
exponent factor for curved growth direction (0.0 -&gt; straight)
Definition: pentagrow.h:419
void updateShellNormals()
determine vertex normals for envelope
Definition: pentagrow.cpp:2330
bool chattyOptimization
whether to log function values during optimization
Definition: pentagrow.h:467
void edgeLengthStats(uint k, Real &lmean, Real &lmax, Real &lmin) const
compute wall mesh edge length statistics around node k
Definition: pentagrow.cpp:1559
Vector prismQualitySumCos(const std::string &fname, uint isection, uint nbin=NotFound) const
compute prism quality histogram and write to file
Definition: pentagrow.cpp:2282
Real convexity(const Vct3 &p1, const Vct3 &n1, const Vct3 &p2, const Vct3 &n2)
test for convexity, returns positive values for convex features
Definition: pentagrow.h:179
const Vct3 & vertex(uint i) const
void smooth(const ConnectMap &map, Container &c) const
generalized nodal smoothing.
Definition: pentagrow.h:242
uint extrude(bool curvedGrowth)
extrude between wall and envelope, return prism layer mesh section index
Definition: pentagrow.cpp:1317
bool isClass(size_t i, int cat) const
test whether a vertex is exactly of a certain category
Definition: pentagrow.h:169
void smoothShellNodes(const Vector &lyt, int niter, Real omega=0.5)
Laplace smoothing of outer shell node coordinates.
Definition: pentagrow.cpp:1272
static uint maximumTriangleCount()
maximum permitted number of boundary triangles
Definition: pentagrow.h:63
bool hasClass(size_t i, int cat) const
test whether a vertex has at least a certain category
Definition: pentagrow.h:174
void mergeNeighbors(Indices &idx) const
augment a set of vertices with its direct neighbors
Definition: pentagrow.cpp:2018
NDPointTree< 3, float > nodeTree
search tree for nodes in the outer layer
Definition: pentagrow.h:437
uint npoints() const
Real maxOptimizationTime
maximum time to be used by numerical optimization (default 30 seconds)
Definition: pentagrow.h:458
Vct3 projectToWall(const TriMesh &mout, const Vct3 &pout, uint inear) const
determine wall point for outer mesh point pout
Definition: pentagrow.cpp:1780
static int maximumTagValue()
maximum permitted section tag value
Definition: pentagrow.h:58
bool find(const NDPoint &pt, FloatType r, Indices &fnd) const
void configure(const ConfigParser &cfg)
set configuration options
Definition: pentagrow.cpp:40
const_iterator begin(uint ir) const
size_t size() const
void extrudeVertex(int ivx, int nl, Real hi, bool curvedGrowth, PointGrid< 3 > &grid) const
extrude a single vertex (nucleaus for parallelization)
Definition: pentagrow.cpp:1375
Real defaultInvGrowthExp
growth exponent factor : make this larger for improved wall-normality
Definition: pentagrow.h:461
void smooth(const ConnectMap &map, Container &c, const WritePred &writeNode, const ReadPred &readNode) const
generalized nodal smoothing.
Definition: pentagrow.h:268
Real cosFeatureAngle
Feature angle for geometrical identification.
Definition: pentagrow.h:443
Vector targetHeight
target height values for NLOPT
Definition: pentagrow.h:416
DVector< int > etype
integer flag indicating mesh edge category
Definition: pentagrow.h:425
Real cosconcave
Angle for concave/convex identification.
Definition: pentagrow.h:446
size_t nWallNodes() const
number of wall nodes
Definition: pentagrow.h:72
bool collisions(Indices &colliding, uint iwall, Real safety, Real nrmdev=cos(rad(60.)), Real fnrmdev=cos(rad(120.))) const
find collision candidates using normal criterion
Definition: pentagrow.cpp:1079
void shrink()
reduce memory footprint by erasing all working data (only raw mesh left)
Definition: pentagrow.cpp:2236
uint size() const
void envelopeBounds(Vct3 &plo, Vct3 &phi) const
compute bounding box of wall mesh
Definition: pentagrow.cpp:2245
PointList< 3 > vout
outermost layer
Definition: pentagrow.h:407
uint nvertices() const
Task for parallel resolution of indirect collisions.
Definition: pentagrow.h:353
const_iterator end(uint ir) const
int untangle(Vector &lyt, int niter, Real permitted_etwist)
reduce edge twist
Definition: pentagrow.cpp:914
void writeTetgen(const std::string &fname, const TriMesh &farf, const PointList< 3 > &holes, const TriMesh &refr=TriMesh(), Real nearBoxEdge=0.0)
export boundaries to tetgen smesh file
Definition: pentagrow.cpp:1899
DVector< uint > id2index
maps triangle id passed to tetgen to original triangle index
Definition: pentagrow.h:428
MxMeshSection farfieldSection
farfield section generated by adaptWall
Definition: pentagrow.h:452
void smoothThickness(Vector &lyt, int niter) const
smooth thickness of prismatic layers
Definition: pentagrow.cpp:1198
PointList< 3 > fudir
local wall coordinate system for optimization
Definition: pentagrow.h:413
void classify()
classify and rank vertices
Definition: pentagrow.cpp:89
void edgeMap(ConnectMap &map) const
build vertex-to-vertex connectivity for the surface
Definition: pentagrow.cpp:70
void adjustRidgeNormals()
determine initial wall normals
Definition: pentagrow.cpp:316
uint appendPrismLayer(const PointGrid< 3 > &grid)
add pentahedral elements to mesh, return mesh section index
Definition: pentagrow.cpp:1480
Vct3 nbBarycenter(const PointList< 3 > &pts, size_t k) const
compute barycenter of local neighborhood of node k
Definition: pentagrow.cpp:1299
bool attemptGridUntangling
whether to attamept grid untangling or not (default - yes)
Definition: pentagrow.h:470
void findEnvelopeNeighbors(Indices &interfaceNodes, Indices &nearTetNodes) const
find nodes which are part of prismatic and tetrahedral region
Definition: pentagrow.cpp:2405
size_t untangleGrid(PointGrid< 3 > &grid)
attempt to untangle remaining tangled pentahedra
Definition: pentagrow.cpp:2038
void retractNeighbors(const Indices &afv, Vector &lyt, int ring=3)
ring-2 smoothing of affected vertices in untangle/unwarp/uncollide
Definition: pentagrow.cpp:776
Real firstCellHeight
configuration parameters
Definition: pentagrow.h:455
int numPrismLayers
configuration parameters
Definition: pentagrow.h:464
uint32_t extractSectionTag(uint32_t tag) const
extract section tag from tetgen tag
Definition: pentagrow.h:228
void envelopeEdgeStats(Real &lmean, Real &lmax) const
used to suggest near-field refinement factor: envelope edge lengths
Definition: pentagrow.cpp:2266
void unwarp(int niter, Real permitted_angle)
reduce penta warp
Definition: pentagrow.cpp:810
PointList< 3 > wfn
smoothed wall normals
Definition: pentagrow.h:404
void adaptWall(const DVector< uint > &faceTags)
adapt wall from refined outer shell (from tetgen)
Definition: pentagrow.cpp:1587
DVector< uint > id2section
maps triangle id passed to tetgen to section index
Definition: pentagrow.h:431
Indices wallTags
mesh tags which contain wall boundary
Definition: pentagrow.h:434
void rebuildTree()
rebuild search tree using current set of outer-layer vertices
Definition: pentagrow.cpp:2222
Indices::const_iterator const_iterator
TriMesh mwall
wall mesh (must be watertight)
Definition: pentagrow.h:401
void generateShell(int hiter, int niter, int ncrititer, int laplaceiter)
generate the outermost layer
Definition: pentagrow.cpp:399
Real convexity(uint i1, uint i2)
test for convexity, returns positive values for convex features
Definition: pentagrow.h:186
void centerGridNodes(uint niter, PointGrid< 3 > &grid) const
move grid vertices to the barycenter of their neighborhood
Definition: pentagrow.cpp:2338
uint uncollideVertex(uint i, Vector &lyt, Real safety, Real retraction, Real cphi, Real cphif)
uncollide a single vertex (nucleus function for parallelization)
Definition: pentagrow.cpp:1058
DVector< int > vtype
integer flag indicating vertex category
Definition: pentagrow.h:422
void uncollide(int niter, Real safety, Real retraction, Real limitphi, Real limitphif)
resolve indirect collisions
Definition: pentagrow.cpp:995
Task for parallel extrusion of prismatic mesh grid.
Definition: pentagrow.h:381
void extractWall(const MxMesh &gm)
extract wall from global mesh
Definition: pentagrow.cpp:1873
uint32_t extractElementTag(uint32_t tag) const
extract element tag from tetgen tag
Definition: pentagrow.h:233
void centerGridNodesPass(const PointGrid< 3 > &cgrid, PointGrid< 3 > &grid) const
move grid vertices to the barycenter of their neighborhood, single pass
Definition: pentagrow.cpp:2376
Generated on Wed Jan 19 2022 03:03:15 for libsurf by   doxygen 1.8.5