Nektar++
Geometry.h
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // File: Geometry.h
4 //
5 // For more information, please see: http://www.nektar.info/
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // Permission is hereby granted, free of charge, to any person obtaining a
14 // copy of this software and associated documentation files (the "Software"),
15 // to deal in the Software without restriction, including without limitation
16 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
17 // and/or sell copies of the Software, and to permit persons to whom the
18 // Software is furnished to do so, subject to the following conditions:
19 //
20 // The above copyright notice and this permission notice shall be included
21 // in all copies or substantial portions of the Software.
22 //
23 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29 // DEALINGS IN THE SOFTWARE.
30 //
31 // Description: This file contains the base class specification for the
32 // Geometry class.
33 //
34 //
35 ////////////////////////////////////////////////////////////////////////////////
36 
37 #ifndef NEKTAR_SPATIALDOMAINS_GEOMETRY_H
38 #define NEKTAR_SPATIALDOMAINS_GEOMETRY_H
39 
43 
44 #include <array>
45 #include <unordered_map>
46 
47 namespace Nektar
48 {
49 
50 namespace SpatialDomains
51 {
52 
53 class Geometry; // Forward declaration for typedef.
54 typedef std::shared_ptr<Geometry> GeometrySharedPtr;
55 typedef std::vector<GeometrySharedPtr> GeometryVector;
56 typedef std::unordered_set<GeometrySharedPtr> GeometrySet;
57 typedef std::shared_ptr<GeometryVector> GeometryVectorSharedPtr;
58 
59 class PointGeom;
60 typedef std::shared_ptr<PointGeom> PointGeomSharedPtr;
61 
62 class Geometry1D;
63 class Geometry2D;
64 typedef std::shared_ptr<Geometry1D> Geometry1DSharedPtr;
65 typedef std::shared_ptr<Geometry2D> Geometry2DSharedPtr;
66 
67 struct Curve;
68 typedef std::shared_ptr<Curve> CurveSharedPtr;
69 typedef std::unordered_map<int, CurveSharedPtr> CurveMap;
71 
72 /// \brief Less than operator to sort Geometry objects by global id when sorting
73 /// STL containers.
75  const std::shared_ptr<Geometry> &lhs, const std::shared_ptr<Geometry> &rhs);
76 
78  const std::shared_ptr<Geometry> &lhs, const std::shared_ptr<Geometry> &rhs);
79 
80 /// Base class for shape geometry information
81 class Geometry
82 {
83 public:
85  SPATIAL_DOMAINS_EXPORT Geometry(int coordim);
86 
88 
89  //---------------------------------------
90  // Helper functions
91  //---------------------------------------
92 
93  SPATIAL_DOMAINS_EXPORT inline int GetCoordim() const;
94  SPATIAL_DOMAINS_EXPORT inline void SetCoordim(int coordim);
95 
101 
102  //---------------------------------------
103  // Set and get ID
104  //---------------------------------------
105  SPATIAL_DOMAINS_EXPORT inline int GetGlobalID(void) const;
106  SPATIAL_DOMAINS_EXPORT inline void SetGlobalID(int globalid);
107 
108  //---------------------------------------
109  // Vertex, edge and face access
110  //---------------------------------------
111  SPATIAL_DOMAINS_EXPORT int GetVid(int i) const;
112  SPATIAL_DOMAINS_EXPORT int GetEid(int i) const;
113  SPATIAL_DOMAINS_EXPORT int GetFid(int i) const;
114  SPATIAL_DOMAINS_EXPORT inline int GetTid(int i) const;
119  const int i) const;
121  const int i) const;
122  SPATIAL_DOMAINS_EXPORT inline int GetNumVerts() const;
123  SPATIAL_DOMAINS_EXPORT inline int GetNumEdges() const;
124  SPATIAL_DOMAINS_EXPORT inline int GetNumFaces() const;
125  SPATIAL_DOMAINS_EXPORT inline int GetShapeDim() const;
126 
127  //---------------------------------------
128  // \chi mapping access
129  //---------------------------------------
131  const;
133  const int i) const;
134  SPATIAL_DOMAINS_EXPORT inline void FillGeom();
135 
136  //---------------------------------------
137  // Point lookups
138  //---------------------------------------
139  SPATIAL_DOMAINS_EXPORT std::array<NekDouble, 6> GetBoundingBox();
140 
142  const Array<OneD, const NekDouble> &gloCoord, NekDouble tol = 0.0);
144  const Array<OneD, const NekDouble> &gloCoord,
145  Array<OneD, NekDouble> &locCoord, NekDouble tol);
147  const Array<OneD, const NekDouble> &gloCoord,
148  Array<OneD, NekDouble> &locCoord, NekDouble tol, NekDouble &dist);
150  const Array<OneD, const NekDouble> &coords,
151  Array<OneD, NekDouble> &Lcoords);
153  const int i, const Array<OneD, const NekDouble> &Lcoord);
155  const Array<OneD, const NekDouble> &gloCoord);
157  NekDouble tol);
158 
159  //---------------------------------------
160  // Misc. helper functions
161  //---------------------------------------
162  SPATIAL_DOMAINS_EXPORT inline int GetVertexEdgeMap(int i, int j) const;
163  SPATIAL_DOMAINS_EXPORT inline int GetVertexFaceMap(int i, int j) const;
164  SPATIAL_DOMAINS_EXPORT inline int GetEdgeFaceMap(int i, int j) const;
166  int j) const;
167  SPATIAL_DOMAINS_EXPORT inline int GetDir(const int i,
168  const int j = 0) const;
169 
170  SPATIAL_DOMAINS_EXPORT inline void Reset(CurveMap &curvedEdges,
171  CurveMap &curvedFaces);
172  SPATIAL_DOMAINS_EXPORT inline void Setup();
173 
174 protected:
176  GeomFactorsSharedPtr geomFactor);
178 
179  /// Coordinate dimension of this geometry object.
181  /// Geometric factors.
183  /// State of the geometric factors
185  /// \f$\chi\f$ mapping containing isoparametric transformation.
187  /// Enumeration to dictate whether coefficients are filled.
189  /// Wether or not the setup routines have been run
191  /// Type of geometry.
193  /// Type of shape.
195  /// Global ID
197  /// Array containing expansion coefficients of @p m_xmap
199  /// Array containing bounding box
201 
202  /// Handles generation of geometry factors.
203  void GenGeomFactors();
204 
205  //---------------------------------------
206  // Helper functions
207  //---------------------------------------
208  virtual PointGeomSharedPtr v_GetVertex(int i) const = 0;
209  virtual Geometry1DSharedPtr v_GetEdge(int i) const;
210  virtual Geometry2DSharedPtr v_GetFace(int i) const;
211  virtual StdRegions::Orientation v_GetEorient(const int i) const;
212  virtual StdRegions::Orientation v_GetForient(const int i) const;
213  virtual int v_GetNumVerts() const;
214  virtual int v_GetNumEdges() const;
215  virtual int v_GetNumFaces() const;
216  virtual int v_GetShapeDim() const;
217 
219  virtual void v_FillGeom();
220 
221  virtual bool v_ContainsPoint(const Array<OneD, const NekDouble> &gloCoord,
222  Array<OneD, NekDouble> &locCoord,
223  NekDouble tol, NekDouble &dist);
224 
225  virtual NekDouble v_GetCoord(const int i,
226  const Array<OneD, const NekDouble> &Lcoord);
228  Array<OneD, NekDouble> &Lcoords);
229 
230  virtual int v_GetVertexEdgeMap(int i, int j) const;
231  virtual int v_GetVertexFaceMap(int i, int j) const;
232  virtual int v_GetEdgeFaceMap(int i, int j) const;
233  virtual int v_GetEdgeNormalToFaceVert(const int i, const int j) const;
234  virtual int v_GetDir(const int faceidx, const int facedir) const;
235 
236  virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces);
237  virtual void v_Setup();
238  virtual void v_GenGeomFactors() = 0;
239 
240  inline void SetUpCoeffs(const int nCoeffs);
241 }; // class Geometry
242 
243 /**
244  * @brief Unary function that constructs a hash of a Geometry object, based on
245  * the vertex IDs.
246  */
247 struct GeometryHash : std::unary_function<GeometrySharedPtr, std::size_t>
248 {
249  std::size_t operator()(GeometrySharedPtr const &p) const
250  {
251  int i;
252  size_t seed = 0;
253  int nVert = p->GetNumVerts();
254  std::vector<unsigned int> ids(nVert);
255 
256  for (i = 0; i < nVert; ++i)
257  {
258  ids[i] = p->GetVid(i);
259  }
260  std::sort(ids.begin(), ids.end());
261  hash_range(seed, ids.begin(), ids.end());
262 
263  return seed;
264  }
265 };
266 
267 /**
268  * @brief Return the coordinate dimension of this object (i.e. the dimension of
269  * the space in which this object is embedded).
270  */
271 inline int Geometry::GetCoordim() const
272 {
273  return m_coordim;
274 }
275 
276 /**
277  * @brief Sets the coordinate dimension of this object (i.e. the dimension of
278  * the space in which this object is embedded).
279  */
280 inline void Geometry::SetCoordim(int dim)
281 {
282  m_coordim = dim;
283 }
284 
285 /**
286  * @brief Get the geometric factors for this object, generating them if
287  * required.
288  */
290 {
291  GenGeomFactors();
293 }
294 
295 /**
296  * @brief Get the geometric factors for this object.
297  */
299 {
300  return m_geomFactors;
301 }
302 
303 /**
304  * @brief Get the geometric shape type of this object.
305  */
307 {
308  return m_shapeType;
309 }
310 
311 /**
312  * @brief Get the ID of this object.
313  */
314 inline int Geometry::GetGlobalID(void) const
315 {
316  return m_globalID;
317 }
318 
319 /**
320  * @brief Set the ID of this object.
321  */
322 inline void Geometry::SetGlobalID(int globalid)
323 {
324  m_globalID = globalid;
325 }
326 
327 /**
328  * @brief Get the ID of trace @p i of this object.
329  *
330  * The trace element is the facet one dimension lower than the object; for
331  * example, a quadrilateral has four trace segments forming its boundary.
332  */
333 inline int Geometry::GetTid(int i) const
334 {
335  const int nDim = GetShapeDim();
336  return nDim == 1 ? GetVid(i)
337  : nDim == 2 ? GetEid(i)
338  : nDim == 3 ? GetFid(i)
339  : 0;
340 }
341 
342 /**
343  * @brief Returns vertex @p i of this object.
344  */
346 {
347  return v_GetVertex(i);
348 }
349 
350 /**
351  * @brief Returns edge @p i of this object.
352  */
354 {
355  return v_GetEdge(i);
356 }
357 
358 /**
359  * @brief Returns face @p i of this object.
360  */
362 {
363  return v_GetFace(i);
364 }
365 
366 /**
367  * @brief Returns the orientation of edge @p i with respect to the ordering of
368  * edges in the standard element.
369  */
371 {
372  return v_GetEorient(i);
373 }
374 
375 /**
376  * @brief Returns the orientation of face @p i with respect to the ordering of
377  * faces in the standard element.
378  */
380 {
381  return v_GetForient(i);
382 }
383 
384 /**
385  * @brief Get the number of vertices of this object.
386  */
387 inline int Geometry::GetNumVerts() const
388 {
389  return v_GetNumVerts();
390 }
391 
392 /**
393  * @brief Get the number of edges of this object.
394  */
395 inline int Geometry::GetNumEdges() const
396 {
397  return v_GetNumEdges();
398 }
399 
400 /**
401  * @brief Get the number of faces of this object.
402  */
403 inline int Geometry::GetNumFaces() const
404 {
405  return v_GetNumFaces();
406 }
407 
408 /**
409  * @brief Get the object's shape dimension.
410  *
411  * For example, a segment is one dimensional and quadrilateral is two
412  * dimensional.
413  */
414 inline int Geometry::GetShapeDim() const
415 {
416  return v_GetShapeDim();
417 }
418 
419 /**
420  * @brief Return the mapping object Geometry::m_xmap that represents the
421  * coordinate transformation from standard element to physical element.
422  */
424 {
425  return v_GetXmap();
426 }
427 
428 /**
429  * @brief Return the coefficients of the transformation Geometry::m_xmap in
430  * coordinate direction @p i.
431  */
433  const int i) const
434 {
435  return m_coeffs[i];
436 }
437 
438 /**
439  * @brief Populate the coordinate mapping Geometry::m_coeffs information from
440  * any children geometry elements.
441  *
442  * @see v_FillGeom()
443  */
444 inline void Geometry::FillGeom()
445 {
446  v_FillGeom();
447 }
448 
449 /**
450  * @brief Determine whether an element contains a particular Cartesian
451  * coordinate \f$(x,y,z)\f$.
452  *
453  * @see Geometry::ContainsPoint
454  */
456  const Array<OneD, const NekDouble> &gloCoord, NekDouble tol)
457 {
458  Array<OneD, NekDouble> locCoord(GetCoordim(), 0.0);
459  NekDouble dist;
460  return v_ContainsPoint(gloCoord, locCoord, tol, dist);
461 }
462 
463 /**
464  * @brief Determine whether an element contains a particular Cartesian
465  * coordinate \f$(x,y,z)\f$.
466  *
467  * @see Geometry::ContainsPoint
468  */
470  const Array<OneD, const NekDouble> &gloCoord,
471  Array<OneD, NekDouble> &locCoord, NekDouble tol)
472 {
473  NekDouble dist;
474  return v_ContainsPoint(gloCoord, locCoord, tol, dist);
475 }
476 
477 /**
478  * @brief Determine whether an element contains a particular Cartesian
479  * coordinate \f$\vec{x} = (x,y,z)\f$.
480  *
481  * For curvilinear and non-affine elements (i.e. where the Jacobian varies as a
482  * function of the standard element coordinates), this is a non-linear
483  * optimisation problem that requires the use of a Newton iteration. Note
484  * therefore that this can be an expensive operation.
485  *
486  * The parameter @p tol which is by default 0, can be used to expand the
487  * coordinate range of the standard element from \f$[-1,1]^d\f$ to
488  * \f$[-1-\epsilon,1+\epsilon\f$ to handle challenging edge cases. The function
489  * also returns the local coordinates corresponding to @p gloCoord that can be
490  * used to speed up subsequent searches.
491  *
492  * @param gloCoord The coordinate \f$ (x,y,z) \f$.
493  * @param locCoord On exit, this is the local coordinate \f$\vec{\xi}\f$ such
494  * that \f$\chi(\vec{\xi}) = \vec{x}\f$.
495  * @param tol The tolerance used to dictate the bounding box of the
496  * standard coordinates \f$\vec{\xi}\f$.
497  * @param dist On exit, returns the minimum distance between @p gloCoord
498  * and the quadrature points inside the element.
499  *
500  * @return `true` if the coordinate @p gloCoord is contained in the element;
501  * false otherwise.
502  *
503  * @see Geometry::GetLocCoords.
504  */
506  const Array<OneD, const NekDouble> &gloCoord,
507  Array<OneD, NekDouble> &locCoord, NekDouble tol, NekDouble &dist)
508 {
509  return v_ContainsPoint(gloCoord, locCoord, tol, dist);
510 }
511 
512 /**
513  * @brief Determine the local collapsed coordinates that correspond to a
514  * given Cartesian coordinate for this geometry object.
515  *
516  * For curvilinear and non-affine elements (i.e. where the Jacobian varies as a
517  * function of the standard element coordinates), this is a non-linear
518  * optimisation problem that requires the use of a Newton iteration. Note
519  * therefore that this can be an expensive operation.
520  *
521  * Note that, clearly, the provided Cartesian coordinate lie outside the
522  * element. The function therefore returns the minimum distance from some
523  * position in the element to . @p Lcoords will also be constrained to fit
524  * within the range \f$[-1,1]^d\f$ where \f$ d \f$ is the dimension of the
525  * element.
526  *
527  * @param coords Input Cartesian global coordinates
528  * @param Lcoords Corresponding local coordinates
529  *
530  * @return Distance between obtained coordinates and provided ones.
531  */
534 {
535  return v_GetLocCoords(coords, Lcoords);
536 }
537 
538 /**
539  * @brief Given local collapsed coordinate @p Lcoord, return the value of
540  * physical coordinate in direction @p i.
541  */
542 inline NekDouble Geometry::GetCoord(const int i,
543  const Array<OneD, const NekDouble> &Lcoord)
544 {
545  return v_GetCoord(i, Lcoord);
546 }
547 
548 /**
549  * @brief Returns the standard element edge IDs that are connected to a given
550  * vertex.
551  *
552  * For example, on a prism, vertex 0 is connnected to edges 0, 3, and 4;
553  * `GetVertexEdgeMap(0,j)` would therefore return the values 0, 1 and 4
554  * respectively. We assume that @p j runs between 0 and 2 inclusive, which is
555  * true for every 3D element asides from the pyramid.
556  *
557  * This function is used in the construction of the low-energy preconditioner.
558  *
559  * @param i The vertex to query connectivity for.
560  * @param j The local edge index between 0 and 2 connected to this element.
561  *
562  * @todo Expand to work with pyramid elements.
563  * @see MultiRegions::PreconditionerLowEnergy
564  */
565 inline int Geometry::GetVertexEdgeMap(int i, int j) const
566 {
567  return v_GetVertexEdgeMap(i, j);
568 }
569 
570 /**
571  * @brief Returns the standard element face IDs that are connected to a given
572  * vertex.
573  *
574  * For example, on a hexahedron, vertex 0 is connnected to faces 0, 1, and 4;
575  * `GetVertexFaceMap(0,j)` would therefore return the values 0, 1 and 4
576  * respectively. We assume that @p j runs between 0 and 2 inclusive, which is
577  * true for every 3D element asides from the pyramid.
578  *
579  * This is used in the construction of the low-energy preconditioner.
580  *
581  * @param i The vertex to query connectivity for.
582  * @param j The local face index between 0 and 2 connected to this element.
583  *
584  * @todo Expand to work with pyramid elements.
585  * @see MultiRegions::PreconditionerLowEnergy
586  */
587 inline int Geometry::GetVertexFaceMap(int i, int j) const
588 {
589  return v_GetVertexFaceMap(i, j);
590 }
591 
592 /**
593  * @brief Returns the standard element edge IDs that are connected to a given
594  * face.
595  *
596  * For example, on a prism, edge 0 is connnected to faces 0 and 1;
597  * `GetEdgeFaceMap(0,j)` would therefore return the values 0 and 1
598  * respectively. We assume that @p j runs between 0 and 1 inclusive, since every
599  * face is connected to precisely two faces for all 3D elements.
600  *
601  * This function is used in the construction of the low-energy preconditioner.
602  *
603  * @param i The edge to query connectivity for.
604  * @param j The local face index between 0 and 1 connected to this element.
605  *
606  * @see MultiRegions::PreconditionerLowEnergy
607  */
608 inline int Geometry::GetEdgeFaceMap(int i, int j) const
609 {
610  return v_GetEdgeFaceMap(i, j);
611 }
612 
613 /**
614  * @brief Returns the standard lement edge IDs that are normal to a given face
615  * vertex.
616  *
617  * For example, on a hexahedron, on face 0 at vertices 0,1,2,3 the
618  * edges normal to that face are 4,5,6,7, ; so
619  * `GetEdgeNormalToFaceVert(0,j)` would therefore return the values 4,
620  * 5, 6 and 7 respectively. We assume that @p j runs between 0 and 3
621  * inclusive on a quadrilateral face and between 0 and 2 inclusive on
622  * a triangular face.
623  *
624  * This is used to help set up a length scale normal to an face
625  *
626  * @param i The face to query for the normal edge
627  * @param j The local vertex index between 0 and nverts on this face
628  *
629  */
630 inline int Geometry::GetEdgeNormalToFaceVert(int i, int j) const
631 {
632  return v_GetEdgeNormalToFaceVert(i, j);
633 }
634 
635 /**
636  * @brief Returns the element coordinate direction corresponding to a given face
637  * coordinate direction
638  */
639 inline int Geometry::GetDir(const int faceidx, const int facedir) const
640 {
641  return v_GetDir(faceidx, facedir);
642 }
643 
644 /**
645  * @brief Reset this geometry object: unset the current state, zero
646  * Geometry::m_coeffs and remove allocated GeomFactors.
647  */
648 inline void Geometry::Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
649 {
650  v_Reset(curvedEdges, curvedFaces);
651 }
652 inline void Geometry::Setup()
653 {
654  v_Setup();
655 }
656 
657 /**
658  * @brief Generate the geometric factors (i.e. derivatives of \f$\chi\f$) and
659  * related metrics.
660  *
661  * @see SpatialDomains::GeomFactors
662  */
664 {
665  return v_GenGeomFactors();
666 }
667 
668 /**
669  * @brief Initialise the Geometry::m_coeffs array.
670  */
671 inline void Geometry::SetUpCoeffs(const int nCoeffs)
672 {
674 
675  for (int i = 0; i < m_coordim; ++i)
676  {
677  m_coeffs[i] = Array<OneD, NekDouble>(nCoeffs, 0.0);
678  }
679 }
680 
681 } // namespace SpatialDomains
682 } // namespace Nektar
683 
684 #endif // NEKTAR_SPATIALDOMAINS_GEOMETRY_H
#define SPATIAL_DOMAINS_EXPORT
1D geometry information
Definition: Geometry1D.h:55
2D geometry information
Definition: Geometry2D.h:69
Base class for shape geometry information.
Definition: Geometry.h:82
virtual int v_GetNumEdges() const
Get the number of edges of this object.
Definition: Geometry.cpp:217
virtual StdRegions::Orientation v_GetForient(const int i) const
Returns the orientation of face i with respect to the ordering of faces in the standard element.
Definition: Geometry.cpp:206
NekDouble GetCoord(const int i, const Array< OneD, const NekDouble > &Lcoord)
Given local collapsed coordinate Lcoord, return the value of physical coordinate in direction i.
Definition: Geometry.h:542
LibUtilities::ShapeType GetShapeType(void)
Get the geometric shape type of this object.
Definition: Geometry.h:306
virtual NekDouble v_GetCoord(const int i, const Array< OneD, const NekDouble > &Lcoord)
Given local collapsed coordinate Lcoord, return the value of physical coordinate in direction i.
Definition: Geometry.cpp:334
Array< OneD, NekDouble > m_boundingBox
Array containing bounding box.
Definition: Geometry.h:200
GeomFactorsSharedPtr GetRefGeomFactors(const Array< OneD, const LibUtilities::BasisSharedPtr > &tbasis)
bool m_setupState
Wether or not the setup routines have been run.
Definition: Geometry.h:190
PointGeomSharedPtr GetVertex(int i) const
Returns vertex i of this object.
Definition: Geometry.h:345
GeomState m_state
Enumeration to dictate whether coefficients are filled.
Definition: Geometry.h:188
bool ClampLocCoords(Array< OneD, NekDouble > &locCoord, NekDouble tol)
Clamp local coords to be within standard regions [-1, 1]^dim.
Definition: Geometry.cpp:496
void SetUpCoeffs(const int nCoeffs)
Initialise the Geometry::m_coeffs array.
Definition: Geometry.h:671
NekDouble GetLocCoords(const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this ge...
Definition: Geometry.h:532
GeomType m_geomType
Type of geometry.
Definition: Geometry.h:192
virtual int v_GetEdgeNormalToFaceVert(const int i, const int j) const
Returns the standard lement edge IDs that are normal to a given face vertex.
Definition: Geometry.cpp:312
Geometry1DSharedPtr GetEdge(int i) const
Returns edge i of this object.
Definition: Geometry.h:353
void SetGlobalID(int globalid)
Set the ID of this object.
Definition: Geometry.h:322
int GetNumFaces() const
Get the number of faces of this object.
Definition: Geometry.h:403
int GetShapeDim() const
Get the object's shape dimension.
Definition: Geometry.h:414
virtual void v_FillGeom()
Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements.
Definition: Geometry.cpp:358
Geometry()
Default constructor.
Definition: Geometry.cpp:54
virtual NekDouble v_GetLocCoords(const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this ge...
Definition: Geometry.cpp:346
virtual Geometry2DSharedPtr v_GetFace(int i) const
Returns face i of this object.
Definition: Geometry.cpp:174
int GetVid(int i) const
Get the ID of vertex i of this object.
Definition: Geometry.cpp:139
virtual int v_GetNumFaces() const
Get the number of faces of this object.
Definition: Geometry.cpp:225
virtual StdRegions::StdExpansionSharedPtr v_GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition: Geometry.cpp:243
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated Geo...
Definition: Geometry.cpp:367
const Array< OneD, const NekDouble > & GetCoeffs(const int i) const
Return the coefficients of the transformation Geometry::m_xmap in coordinate direction i.
Definition: Geometry.h:432
int GetGlobalID(void) const
Get the ID of this object.
Definition: Geometry.h:314
int GetCoordim() const
Return the coordinate dimension of this object (i.e. the dimension of the space in which this object ...
Definition: Geometry.h:271
int GetFid(int i) const
Get the ID of face i of this object.
Definition: Geometry.cpp:155
virtual int v_GetEdgeFaceMap(int i, int j) const
Returns the standard element edge IDs that are connected to a given face.
Definition: Geometry.cpp:301
LibUtilities::ShapeType m_shapeType
Type of shape.
Definition: Geometry.h:194
int GetDir(const int i, const int j=0) const
Returns the element coordinate direction corresponding to a given face coordinate direction.
Definition: Geometry.h:639
void FillGeom()
Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements.
Definition: Geometry.h:444
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
Definition: Geometry.h:198
virtual ~Geometry()
Default destructor.
Definition: Geometry.cpp:74
int GetVertexEdgeMap(int i, int j) const
Returns the standard element edge IDs that are connected to a given vertex.
Definition: Geometry.h:565
GeomState m_geomFactorsState
State of the geometric factors.
Definition: Geometry.h:184
virtual PointGeomSharedPtr v_GetVertex(int i) const =0
GeomFactorsSharedPtr GetMetricInfo()
Get the geometric factors for this object.
Definition: Geometry.h:298
StdRegions::StdExpansionSharedPtr m_xmap
mapping containing isoparametric transformation.
Definition: Geometry.h:186
virtual Geometry1DSharedPtr v_GetEdge(int i) const
Returns edge i of this object.
Definition: Geometry.cpp:163
virtual StdRegions::Orientation v_GetEorient(const int i) const
Returns the orientation of edge i with respect to the ordering of edges in the standard element.
Definition: Geometry.cpp:195
StdRegions::StdExpansionSharedPtr GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition: Geometry.h:423
void GenGeomFactors()
Handles generation of geometry factors.
Definition: Geometry.h:663
static GeomFactorsVector m_regGeomFactorsManager
Definition: Geometry.h:177
int GetNumEdges() const
Get the number of edges of this object.
Definition: Geometry.h:395
std::array< NekDouble, 6 > GetBoundingBox()
Generates the bounding box for the element.
Definition: Geometry.cpp:394
int GetNumVerts() const
Get the number of vertices of this object.
Definition: Geometry.h:387
Geometry2DSharedPtr GetFace(int i) const
Returns face i of this object.
Definition: Geometry.h:361
int GetVertexFaceMap(int i, int j) const
Returns the standard element face IDs that are connected to a given vertex.
Definition: Geometry.h:587
bool ContainsPoint(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
Determine whether an element contains a particular Cartesian coordinate .
Definition: Geometry.h:455
int GetEdgeFaceMap(int i, int j) const
Returns the standard element edge IDs that are connected to a given face.
Definition: Geometry.h:608
bool MinMaxCheck(const Array< OneD, const NekDouble > &gloCoord)
Check if given global coord is within the BoundingBox of the element.
Definition: Geometry.cpp:473
virtual int v_GetVertexEdgeMap(int i, int j) const
Returns the standard element edge IDs that are connected to a given vertex.
Definition: Geometry.cpp:279
int GetEdgeNormalToFaceVert(int i, int j) const
Returns the standard lement edge IDs that are normal to a given face vertex.
Definition: Geometry.h:630
GeomFactorsSharedPtr GetGeomFactors()
Get the geometric factors for this object, generating them if required.
Definition: Geometry.h:289
static GeomFactorsSharedPtr ValidateRegGeomFactor(GeomFactorsSharedPtr geomFactor)
Check to see if a geometric factor has already been created that contains the same regular informatio...
Definition: Geometry.cpp:92
StdRegions::Orientation GetEorient(const int i) const
Returns the orientation of edge i with respect to the ordering of edges in the standard element.
Definition: Geometry.h:370
virtual int v_GetNumVerts() const
Get the number of vertices of this object.
Definition: Geometry.cpp:185
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
Definition: Geometry.h:182
int m_coordim
Coordinate dimension of this geometry object.
Definition: Geometry.h:180
int GetTid(int i) const
Get the ID of trace i of this object.
Definition: Geometry.h:333
virtual int v_GetDir(const int faceidx, const int facedir) const
Returns the element coordinate direction corresponding to a given face coordinate direction.
Definition: Geometry.cpp:323
void SetCoordim(int coordim)
Sets the coordinate dimension of this object (i.e. the dimension of the space in which this object is...
Definition: Geometry.h:280
virtual int v_GetShapeDim() const
Get the object's shape dimension.
Definition: Geometry.cpp:233
void Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated Geo...
Definition: Geometry.h:648
virtual int v_GetVertexFaceMap(int i, int j) const
Returns the standard element face IDs that are connected to a given vertex.
Definition: Geometry.cpp:290
virtual bool v_ContainsPoint(const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &dist)
Definition: Geometry.cpp:253
int GetEid(int i) const
Get the ID of edge i of this object.
Definition: Geometry.cpp:147
StdRegions::Orientation GetForient(const int i) const
Returns the orientation of face i with respect to the ordering of faces in the standard element.
Definition: Geometry.h:379
std::vector< GeomFactorsSharedPtr > GeomFactorsVector
A vector of GeomFactor pointers.
Definition: GeomFactors.h:64
bool GlobalIdEquality(const std::shared_ptr< Geometry > &lhs, const std::shared_ptr< Geometry > &rhs)
Definition: Geometry.cpp:130
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
Definition: GeomFactors.h:62
std::shared_ptr< GeometryVector > GeometryVectorSharedPtr
Definition: Geometry.h:57
bool SortByGlobalId(const std::shared_ptr< Geometry > &lhs, const std::shared_ptr< Geometry > &rhs)
Less than operator to sort Geometry objects by global id when sorting STL containers.
Definition: Geometry.cpp:124
std::unordered_set< GeometrySharedPtr > GeometrySet
Definition: Geometry.h:56
std::shared_ptr< Curve > CurveSharedPtr
Definition: Curve.hpp:60
std::unordered_map< int, CurveSharedPtr > CurveMap
Definition: Curve.hpp:61
GeomType
Indicates the type of element geometry.
static CurveMap NullCurveMap
Definition: Geometry.h:70
std::shared_ptr< PointGeom > PointGeomSharedPtr
Definition: Geometry.h:59
std::shared_ptr< Geometry2D > Geometry2DSharedPtr
Definition: Geometry.h:65
GeomState
Indicates if the geometric information for an element has been populated.
std::shared_ptr< Geometry > GeometrySharedPtr
Definition: Geometry.h:53
std::shared_ptr< Geometry1D > Geometry1DSharedPtr
Definition: Geometry.h:63
std::vector< GeometrySharedPtr > GeometryVector
Definition: Geometry.h:55
std::shared_ptr< StdExpansion > StdExpansionSharedPtr
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:1
std::size_t hash_range(Iter first, Iter last)
Definition: HashUtils.hpp:68
double NekDouble
Unary function that constructs a hash of a Geometry object, based on the vertex IDs.
Definition: Geometry.h:248
std::size_t operator()(GeometrySharedPtr const &p) const
Definition: Geometry.h:249