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
48{
49
50class Geometry; // Forward declaration for typedef.
51typedef std::shared_ptr<Geometry> GeometrySharedPtr;
52typedef std::vector<GeometrySharedPtr> GeometryVector;
53typedef std::unordered_set<GeometrySharedPtr> GeometrySet;
54typedef std::shared_ptr<GeometryVector> GeometryVectorSharedPtr;
55
56class PointGeom;
57typedef std::shared_ptr<PointGeom> PointGeomSharedPtr;
58
59class Geometry1D;
60class Geometry2D;
61typedef std::shared_ptr<Geometry1D> Geometry1DSharedPtr;
62typedef std::shared_ptr<Geometry2D> Geometry2DSharedPtr;
63
64struct Curve;
65typedef std::shared_ptr<Curve> CurveSharedPtr;
66typedef std::unordered_map<int, CurveSharedPtr> CurveMap;
68
69/// \brief Less than operator to sort Geometry objects by global id when sorting
70/// STL containers.
72 const std::shared_ptr<Geometry> &lhs, const std::shared_ptr<Geometry> &rhs);
73
75 const std::shared_ptr<Geometry> &lhs, const std::shared_ptr<Geometry> &rhs);
76
77/// Base class for shape geometry information
79{
80public:
83 SPATIAL_DOMAINS_EXPORT virtual ~Geometry() = default;
84
85 //---------------------------------------
86 // Helper functions
87 //---------------------------------------
88
89 SPATIAL_DOMAINS_EXPORT inline int GetCoordim() const;
90 SPATIAL_DOMAINS_EXPORT inline void SetCoordim(int coordim);
91
97
98 //---------------------------------------
99 // Set and get ID
100 //---------------------------------------
101 SPATIAL_DOMAINS_EXPORT inline int GetGlobalID(void) const;
102 SPATIAL_DOMAINS_EXPORT inline void SetGlobalID(int globalid);
103
104 //---------------------------------------
105 // Vertex, edge and face access
106 //---------------------------------------
107 SPATIAL_DOMAINS_EXPORT int GetVid(int i) const;
108 SPATIAL_DOMAINS_EXPORT int GetEid(int i) const;
109 SPATIAL_DOMAINS_EXPORT int GetFid(int i) const;
110 SPATIAL_DOMAINS_EXPORT inline int GetTid(int i) const;
115 const int i) const;
117 const int i) const;
118 SPATIAL_DOMAINS_EXPORT inline int GetNumVerts() const;
119 SPATIAL_DOMAINS_EXPORT inline int GetNumEdges() const;
120 SPATIAL_DOMAINS_EXPORT inline int GetNumFaces() const;
121 SPATIAL_DOMAINS_EXPORT inline int GetShapeDim() const;
122
123 //---------------------------------------
124 // \chi mapping access
125 //---------------------------------------
127 const;
129 const int i) const;
130 SPATIAL_DOMAINS_EXPORT inline void FillGeom();
131
132 //---------------------------------------
133 // Point lookups
134 //---------------------------------------
135 SPATIAL_DOMAINS_EXPORT std::array<NekDouble, 6> GetBoundingBox();
137
139 const Array<OneD, const NekDouble> &gloCoord, NekDouble tol = 0.0);
141 const Array<OneD, const NekDouble> &gloCoord,
142 Array<OneD, NekDouble> &locCoord, NekDouble tol);
144 const Array<OneD, const NekDouble> &gloCoord,
145 Array<OneD, NekDouble> &locCoord, NekDouble tol, NekDouble &dist);
147 const Array<OneD, const NekDouble> &coords,
148 Array<OneD, NekDouble> &Lcoords);
150 const int i, const Array<OneD, const NekDouble> &Lcoord);
152 const Array<OneD, const NekDouble> &gloCoord);
154 const Array<OneD, const NekDouble> &gloCoord);
156 Array<OneD, NekDouble> &locCoord,
157 NekDouble tol = std::numeric_limits<NekDouble>::epsilon());
160
161 //---------------------------------------
162 // Misc. helper functions
163 //---------------------------------------
164 SPATIAL_DOMAINS_EXPORT inline int GetVertexEdgeMap(int i, int j) const;
165 SPATIAL_DOMAINS_EXPORT inline int GetVertexFaceMap(int i, int j) const;
166 SPATIAL_DOMAINS_EXPORT inline int GetEdgeFaceMap(int i, int j) const;
168 int j) const;
169 SPATIAL_DOMAINS_EXPORT inline int GetDir(const int i,
170 const int j = 0) const;
171
172 SPATIAL_DOMAINS_EXPORT inline void Reset(CurveMap &curvedEdges,
173 CurveMap &curvedFaces);
174 SPATIAL_DOMAINS_EXPORT inline void ResetNonRecursive(CurveMap &curvedEdges,
175 CurveMap &curvedFaces);
176
177 SPATIAL_DOMAINS_EXPORT inline void Setup();
178
179 /// Handles generation of geometry factors.
180 void GenGeomFactors();
181
182protected:
184 GeomFactorsSharedPtr geomFactor);
186
187 /// Coordinate dimension of this geometry object.
189 /// Geometric factors.
191 /// State of the geometric factors
193 /// \f$\chi\f$ mapping containing isoparametric transformation.
195 /// Enumeration to dictate whether coefficients are filled.
197 /// Wether or not the setup routines have been run
199 /// Type of geometry.
201 /// Type of shape.
203 /// Global ID
205 /// Array containing expansion coefficients of @p m_xmap
207 /// Array containing bounding box
212
213 //---------------------------------------
214 // Helper functions
215 //---------------------------------------
216 virtual PointGeomSharedPtr v_GetVertex(int i) const = 0;
217 virtual Geometry1DSharedPtr v_GetEdge(int i) const;
218 virtual Geometry2DSharedPtr v_GetFace(int i) const;
219 virtual StdRegions::Orientation v_GetEorient(const int i) const;
220 virtual StdRegions::Orientation v_GetForient(const int i) const;
221 virtual int v_GetNumVerts() const;
222 virtual int v_GetNumEdges() const;
223 virtual int v_GetNumFaces() const;
224 virtual int v_GetShapeDim() const;
225
227 virtual void v_FillGeom();
228
229 virtual bool v_ContainsPoint(const Array<OneD, const NekDouble> &gloCoord,
230 Array<OneD, NekDouble> &locCoord,
231 NekDouble tol, NekDouble &dist);
232 virtual int v_AllLeftCheck(const Array<OneD, const NekDouble> &gloCoord);
233
234 virtual NekDouble v_GetCoord(const int i,
235 const Array<OneD, const NekDouble> &Lcoord);
237 Array<OneD, NekDouble> &Lcoords);
240
241 virtual int v_GetVertexEdgeMap(int i, int j) const;
242 virtual int v_GetVertexFaceMap(int i, int j) const;
243 virtual int v_GetEdgeFaceMap(int i, int j) const;
244 virtual int v_GetEdgeNormalToFaceVert(const int i, const int j) const;
245 virtual int v_GetDir(const int faceidx, const int facedir) const;
246
247 virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces);
248 virtual void v_Setup();
249 virtual void v_GenGeomFactors() = 0;
250
251 inline void SetUpCoeffs(const int nCoeffs);
252 virtual void v_CalculateInverseIsoParam();
253}; // class Geometry
254
255/**
256 * @brief Unary function that constructs a hash of a Geometry object, based on
257 * the vertex IDs.
258 */
260{
261 std::size_t operator()(GeometrySharedPtr const &p) const
262 {
263 int i;
264 size_t seed = 0;
265 int nVert = p->GetNumVerts();
266 std::vector<unsigned int> ids(nVert);
267
268 for (i = 0; i < nVert; ++i)
269 {
270 ids[i] = p->GetVid(i);
271 }
272 std::sort(ids.begin(), ids.end());
273 hash_range(seed, ids.begin(), ids.end());
274
275 return seed;
276 }
277};
278
279/**
280 * @brief Return the coordinate dimension of this object (i.e. the dimension of
281 * the space in which this object is embedded).
282 */
283inline int Geometry::GetCoordim() const
284{
285 return m_coordim;
286}
287
288/**
289 * @brief Sets the coordinate dimension of this object (i.e. the dimension of
290 * the space in which this object is embedded).
291 */
292inline void Geometry::SetCoordim(int dim)
293{
294 m_coordim = dim;
295}
296
297/**
298 * @brief Get the geometric factors for this object, generating them if
299 * required.
300 */
302{
305}
306
307/**
308 * @brief Get the geometric factors for this object.
309 */
311{
312 return m_geomFactors;
313}
314
315/**
316 * @brief Get the geometric shape type of this object.
317 */
319{
320 return m_shapeType;
321}
322
323/**
324 * @brief Get the ID of this object.
325 */
326inline int Geometry::GetGlobalID(void) const
327{
328 return m_globalID;
329}
330
331/**
332 * @brief Set the ID of this object.
333 */
334inline void Geometry::SetGlobalID(int globalid)
335{
336 m_globalID = globalid;
337}
338
339/**
340 * @brief Get the ID of trace @p i of this object.
341 *
342 * The trace element is the facet one dimension lower than the object; for
343 * example, a quadrilateral has four trace segments forming its boundary.
344 */
345inline int Geometry::GetTid(int i) const
346{
347 const int nDim = GetShapeDim();
348 return nDim == 1 ? GetVid(i)
349 : nDim == 2 ? GetEid(i)
350 : nDim == 3 ? GetFid(i)
351 : 0;
352}
353
354/**
355 * @brief Returns vertex @p i of this object.
356 */
358{
359 return v_GetVertex(i);
360}
361
362/**
363 * @brief Returns edge @p i of this object.
364 */
366{
367 return v_GetEdge(i);
368}
369
370/**
371 * @brief Returns face @p i of this object.
372 */
374{
375 return v_GetFace(i);
376}
377
378/**
379 * @brief Returns the orientation of edge @p i with respect to the ordering of
380 * edges in the standard element.
381 */
383{
384 return v_GetEorient(i);
385}
386
387/**
388 * @brief Returns the orientation of face @p i with respect to the ordering of
389 * faces in the standard element.
390 */
392{
393 return v_GetForient(i);
394}
395
396/**
397 * @brief Get the number of vertices of this object.
398 */
399inline int Geometry::GetNumVerts() const
400{
401 return v_GetNumVerts();
402}
403
404/**
405 * @brief Get the number of edges of this object.
406 */
407inline int Geometry::GetNumEdges() const
408{
409 return v_GetNumEdges();
410}
411
412/**
413 * @brief Get the number of faces of this object.
414 */
415inline int Geometry::GetNumFaces() const
416{
417 return v_GetNumFaces();
418}
419
420/**
421 * @brief Get the object's shape dimension.
422 *
423 * For example, a segment is one dimensional and quadrilateral is two
424 * dimensional.
425 */
426inline int Geometry::GetShapeDim() const
427{
428 return v_GetShapeDim();
429}
430
431/**
432 * @brief Return the mapping object Geometry::m_xmap that represents the
433 * coordinate transformation from standard element to physical element.
434 */
436{
437 return v_GetXmap();
438}
439
440/**
441 * @brief Return the coefficients of the transformation Geometry::m_xmap in
442 * coordinate direction @p i.
443 */
445 const int i) const
446{
447 return m_coeffs[i];
448}
449
450/**
451 * @brief Populate the coordinate mapping Geometry::m_coeffs information from
452 * any children geometry elements.
453 *
454 * @see v_FillGeom()
455 */
457{
458 v_FillGeom();
459}
460
461/**
462 * @brief Determine whether an element contains a particular Cartesian
463 * coordinate \f$(x,y,z)\f$.
464 *
465 * @see Geometry::ContainsPoint
466 */
468 const Array<OneD, const NekDouble> &gloCoord, NekDouble tol)
469{
470 Array<OneD, NekDouble> locCoord(GetCoordim(), 0.0);
471 NekDouble dist;
472 return v_ContainsPoint(gloCoord, locCoord, tol, dist);
473}
474
475/**
476 * @brief Determine whether an element contains a particular Cartesian
477 * coordinate \f$(x,y,z)\f$.
478 *
479 * @see Geometry::ContainsPoint
480 */
482 const Array<OneD, const NekDouble> &gloCoord,
483 Array<OneD, NekDouble> &locCoord, NekDouble tol)
484{
485 NekDouble dist;
486 return v_ContainsPoint(gloCoord, locCoord, tol, dist);
487}
488
489/**
490 * @brief Determine whether an element contains a particular Cartesian
491 * coordinate \f$\vec{x} = (x,y,z)\f$.
492 *
493 * For curvilinear and non-affine elements (i.e. where the Jacobian varies as a
494 * function of the standard element coordinates), this is a non-linear
495 * optimisation problem that requires the use of a Newton iteration. Note
496 * therefore that this can be an expensive operation.
497 *
498 * The parameter @p tol which is by default 0, can be used to expand the
499 * coordinate range of the standard element from \f$[-1,1]^d\f$ to
500 * \f$[-1-\epsilon,1+\epsilon\f$ to handle challenging edge cases. The function
501 * also returns the local coordinates corresponding to @p gloCoord that can be
502 * used to speed up subsequent searches.
503 *
504 * @param gloCoord The coordinate \f$ (x,y,z) \f$.
505 * @param locCoord On exit, this is the local coordinate \f$\vec{\xi}\f$ such
506 * that \f$\chi(\vec{\xi}) = \vec{x}\f$.
507 * @param tol The tolerance used to dictate the bounding box of the
508 * standard coordinates \f$\vec{\xi}\f$.
509 * @param dist On exit, returns the minimum distance between @p gloCoord
510 * and the quadrature points inside the element.
511 *
512 * @return `true` if the coordinate @p gloCoord is contained in the element;
513 * false otherwise.
514 *
515 * @see Geometry::GetLocCoords.
516 */
518 const Array<OneD, const NekDouble> &gloCoord,
519 Array<OneD, NekDouble> &locCoord, NekDouble tol, NekDouble &dist)
520{
521 return v_ContainsPoint(gloCoord, locCoord, tol, dist);
522}
523
524/**
525 * @brief Determine the local collapsed coordinates that correspond to a
526 * given Cartesian coordinate for this geometry object.
527 *
528 * For curvilinear and non-affine elements (i.e. where the Jacobian varies as a
529 * function of the standard element coordinates), this is a non-linear
530 * optimisation problem that requires the use of a Newton iteration. Note
531 * therefore that this can be an expensive operation.
532 *
533 * Note that, clearly, the provided Cartesian coordinate lie outside the
534 * element. The function therefore returns the minimum distance from some
535 * position in the element to . @p Lcoords will also be constrained to fit
536 * within the range \f$[-1,1]^d\f$ where \f$ d \f$ is the dimension of the
537 * element.
538 *
539 * @param coords Input Cartesian global coordinates
540 * @param Lcoords Corresponding local coordinates
541 *
542 * @return Distance between obtained coordinates and provided ones.
543 */
546{
547 return v_GetLocCoords(coords, Lcoords);
548}
549
550/**
551 * @brief Given local collapsed coordinate @p Lcoord, return the value of
552 * physical coordinate in direction @p i.
553 */
554inline NekDouble Geometry::GetCoord(const int i,
555 const Array<OneD, const NekDouble> &Lcoord)
556{
557 return v_GetCoord(i, Lcoord);
558}
559
562{
563 return v_FindDistance(xs, xi);
564}
565
566/**
567 * @brief Returns the standard element edge IDs that are connected to a given
568 * vertex.
569 *
570 * For example, on a prism, vertex 0 is connnected to edges 0, 3, and 4;
571 * `GetVertexEdgeMap(0,j)` would therefore return the values 0, 1 and 4
572 * respectively. We assume that @p j runs between 0 and 2 inclusive, which is
573 * true for every 3D element asides from the pyramid.
574 *
575 * This function is used in the construction of the low-energy preconditioner.
576 *
577 * @param i The vertex to query connectivity for.
578 * @param j The local edge index between 0 and 2 connected to this element.
579 *
580 * @todo Expand to work with pyramid elements.
581 * @see MultiRegions::PreconditionerLowEnergy
582 */
583inline int Geometry::GetVertexEdgeMap(int i, int j) const
584{
585 return v_GetVertexEdgeMap(i, j);
586}
587
588/**
589 * @brief Returns the standard element face IDs that are connected to a given
590 * vertex.
591 *
592 * For example, on a hexahedron, vertex 0 is connnected to faces 0, 1, and 4;
593 * `GetVertexFaceMap(0,j)` would therefore return the values 0, 1 and 4
594 * respectively. We assume that @p j runs between 0 and 2 inclusive, which is
595 * true for every 3D element asides from the pyramid.
596 *
597 * This is used in the construction of the low-energy preconditioner.
598 *
599 * @param i The vertex to query connectivity for.
600 * @param j The local face index between 0 and 2 connected to this element.
601 *
602 * @todo Expand to work with pyramid elements.
603 * @see MultiRegions::PreconditionerLowEnergy
604 */
605inline int Geometry::GetVertexFaceMap(int i, int j) const
606{
607 return v_GetVertexFaceMap(i, j);
608}
609
610/**
611 * @brief Returns the standard element edge IDs that are connected to a given
612 * face.
613 *
614 * For example, on a prism, edge 0 is connnected to faces 0 and 1;
615 * `GetEdgeFaceMap(0,j)` would therefore return the values 0 and 1
616 * respectively. We assume that @p j runs between 0 and 1 inclusive, since every
617 * face is connected to precisely two faces for all 3D elements.
618 *
619 * This function is used in the construction of the low-energy preconditioner.
620 *
621 * @param i The edge to query connectivity for.
622 * @param j The local face index between 0 and 1 connected to this element.
623 *
624 * @see MultiRegions::PreconditionerLowEnergy
625 */
626inline int Geometry::GetEdgeFaceMap(int i, int j) const
627{
628 return v_GetEdgeFaceMap(i, j);
629}
630
631/**
632 * @brief Returns the standard lement edge IDs that are normal to a given face
633 * vertex.
634 *
635 * For example, on a hexahedron, on face 0 at vertices 0,1,2,3 the
636 * edges normal to that face are 4,5,6,7, ; so
637 * `GetEdgeNormalToFaceVert(0,j)` would therefore return the values 4,
638 * 5, 6 and 7 respectively. We assume that @p j runs between 0 and 3
639 * inclusive on a quadrilateral face and between 0 and 2 inclusive on
640 * a triangular face.
641 *
642 * This is used to help set up a length scale normal to an face
643 *
644 * @param i The face to query for the normal edge
645 * @param j The local vertex index between 0 and nverts on this face
646 *
647 */
648inline int Geometry::GetEdgeNormalToFaceVert(int i, int j) const
649{
650 return v_GetEdgeNormalToFaceVert(i, j);
651}
652
653/**
654 * @brief Returns the element coordinate direction corresponding to a given face
655 * coordinate direction
656 */
657inline int Geometry::GetDir(const int faceidx, const int facedir) const
658{
659 return v_GetDir(faceidx, facedir);
660}
661
662/**
663 * @brief Reset this geometry object: unset the current state, zero
664 * Geometry::m_coeffs and remove allocated GeomFactors.
665 */
666inline void Geometry::Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
667{
668 v_Reset(curvedEdges, curvedFaces);
669}
670
671/**
672 * @brief Reset this geometry object non-recursively: unset the current state,
673 * zero Geometry::m_coeffs and remove allocated GeomFactors.
674 */
675inline void Geometry::ResetNonRecursive(CurveMap &curvedEdges,
676 CurveMap &curvedFaces)
677{
678 Geometry::v_Reset(curvedEdges, curvedFaces);
679}
680
681inline void Geometry::Setup()
682{
683 v_Setup();
684}
685
686/**
687 * @brief Generate the geometric factors (i.e. derivatives of \f$\chi\f$) and
688 * related metrics.
689 *
690 * @see SpatialDomains::GeomFactors
691 */
693{
694 return v_GenGeomFactors();
695}
696
697/**
698 * @brief Initialise the Geometry::m_coeffs array.
699 */
700inline void Geometry::SetUpCoeffs(const int nCoeffs)
701{
703
704 for (int i = 0; i < m_coordim; ++i)
705 {
706 m_coeffs[i] = Array<OneD, NekDouble>(nCoeffs, 0.0);
707 }
708}
709
710} // namespace Nektar::SpatialDomains
711
712#endif // NEKTAR_SPATIALDOMAINS_GEOMETRY_H
#define SPATIAL_DOMAINS_EXPORT
1D geometry information
Definition: Geometry1D.h:54
2D geometry information
Definition: Geometry2D.h:65
Base class for shape geometry information.
Definition: Geometry.h:79
virtual int v_GetNumEdges() const
Get the number of edges of this object.
Definition: Geometry.cpp:180
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:169
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:554
LibUtilities::ShapeType GetShapeType(void)
Get the geometric shape type of this object.
Definition: Geometry.h:318
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:331
Array< OneD, NekDouble > m_boundingBox
Array containing bounding box.
Definition: Geometry.h:208
GeomFactorsSharedPtr GetRefGeomFactors(const Array< OneD, const LibUtilities::BasisSharedPtr > &tbasis)
bool m_setupState
Wether or not the setup routines have been run.
Definition: Geometry.h:198
PointGeomSharedPtr GetVertex(int i) const
Returns vertex i of this object.
Definition: Geometry.h:357
GeomState m_state
Enumeration to dictate whether coefficients are filled.
Definition: Geometry.h:196
virtual void v_CalculateInverseIsoParam()
Definition: Geometry.cpp:209
void SetUpCoeffs(const int nCoeffs)
Initialise the Geometry::m_coeffs array.
Definition: Geometry.h:700
virtual int v_AllLeftCheck(const Array< OneD, const NekDouble > &gloCoord)
Definition: Geometry.cpp:203
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:544
GeomType m_geomType
Type of geometry.
Definition: Geometry.h:200
int PreliminaryCheck(const Array< OneD, const NekDouble > &gloCoord)
A fast and robust check if a given global coord is outside of a deformed element. For regular element...
Definition: Geometry.cpp:481
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:309
Geometry1DSharedPtr GetEdge(int i) const
Returns edge i of this object.
Definition: Geometry.h:365
void SetGlobalID(int globalid)
Set the ID of this object.
Definition: Geometry.h:334
int GetNumFaces() const
Get the number of faces of this object.
Definition: Geometry.h:415
int GetShapeDim() const
Get the object's shape dimension.
Definition: Geometry.h:426
virtual void v_FillGeom()
Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements.
Definition: Geometry.cpp:355
Geometry()
Default constructor.
Definition: Geometry.cpp:50
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:343
virtual Geometry2DSharedPtr v_GetFace(int i) const
Returns face i of this object.
Definition: Geometry.cpp:138
int GetVid(int i) const
Get the ID of vertex i of this object.
Definition: Geometry.cpp:104
virtual int v_GetNumFaces() const
Get the number of faces of this object.
Definition: Geometry.cpp:188
virtual StdRegions::StdExpansionSharedPtr v_GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition: Geometry.cpp:218
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
Definition: Geometry.h:209
Array< OneD, Array< OneD, NekDouble > > m_invIsoParam
Definition: Geometry.h:210
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:364
bool ClampLocCoords(Array< OneD, NekDouble > &locCoord, NekDouble tol=std::numeric_limits< NekDouble >::epsilon())
Clamp local coords to be within standard regions [-1, 1]^dim.
Definition: Geometry.cpp:530
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:444
int GetGlobalID(void) const
Get the ID of this object.
Definition: Geometry.h:326
int GetCoordim() const
Return the coordinate dimension of this object (i.e. the dimension of the space in which this object ...
Definition: Geometry.h:283
int GetFid(int i) const
Get the ID of face i of this object.
Definition: Geometry.cpp:120
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:298
LibUtilities::ShapeType m_shapeType
Type of shape.
Definition: Geometry.h:202
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:657
void FillGeom()
Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements.
Definition: Geometry.h:456
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
Definition: Geometry.h:206
int GetVertexEdgeMap(int i, int j) const
Returns the standard element edge IDs that are connected to a given vertex.
Definition: Geometry.h:583
GeomState m_geomFactorsState
State of the geometric factors.
Definition: Geometry.h:192
void ResetNonRecursive(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object non-recursively: unset the current state, zero Geometry::m_coeffs and remo...
Definition: Geometry.h:675
virtual PointGeomSharedPtr v_GetVertex(int i) const =0
GeomFactorsSharedPtr GetMetricInfo()
Get the geometric factors for this object.
Definition: Geometry.h:310
StdRegions::StdExpansionSharedPtr m_xmap
mapping containing isoparametric transformation.
Definition: Geometry.h:194
virtual Geometry1DSharedPtr v_GetEdge(int i) const
Returns edge i of this object.
Definition: Geometry.cpp:128
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:158
StdRegions::StdExpansionSharedPtr GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition: Geometry.h:435
NekDouble FindDistance(const Array< OneD, const NekDouble > &xs, Array< OneD, NekDouble > &xi)
Definition: Geometry.h:560
void GenGeomFactors()
Handles generation of geometry factors.
Definition: Geometry.h:692
static GeomFactorsVector m_regGeomFactorsManager
Definition: Geometry.h:185
int GetNumEdges() const
Get the number of edges of this object.
Definition: Geometry.h:407
std::array< NekDouble, 6 > GetBoundingBox()
Generates the bounding box for the element.
Definition: Geometry.cpp:391
int GetNumVerts() const
Get the number of vertices of this object.
Definition: Geometry.h:399
Geometry2DSharedPtr GetFace(int i) const
Returns face i of this object.
Definition: Geometry.h:373
int GetVertexFaceMap(int i, int j) const
Returns the standard element face IDs that are connected to a given vertex.
Definition: Geometry.h:605
bool ContainsPoint(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
Determine whether an element contains a particular Cartesian coordinate .
Definition: Geometry.h:467
int GetEdgeFaceMap(int i, int j) const
Returns the standard element edge IDs that are connected to a given face.
Definition: Geometry.h:626
bool MinMaxCheck(const Array< OneD, const NekDouble > &gloCoord)
Check if given global coord is within the BoundingBox of the element.
Definition: Geometry.cpp:507
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:276
int GetEdgeNormalToFaceVert(int i, int j) const
Returns the standard lement edge IDs that are normal to a given face vertex.
Definition: Geometry.h:648
virtual NekDouble v_FindDistance(const Array< OneD, const NekDouble > &xs, Array< OneD, NekDouble > &xi)
Definition: Geometry.cpp:264
GeomFactorsSharedPtr GetGeomFactors()
Get the geometric factors for this object, generating them if required.
Definition: Geometry.h:301
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:81
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:382
virtual int v_GetNumVerts() const
Get the number of vertices of this object.
Definition: Geometry.cpp:148
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
Definition: Geometry.h:190
int m_coordim
Coordinate dimension of this geometry object.
Definition: Geometry.h:188
int GetTid(int i) const
Get the ID of trace i of this object.
Definition: Geometry.h:345
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:320
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:292
virtual int v_GetShapeDim() const
Get the object's shape dimension.
Definition: Geometry.cpp:196
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:666
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:287
virtual bool v_ContainsPoint(const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &dist)
Determine whether an element contains a particular Cartesian coordinate .
Definition: Geometry.cpp:229
int GetEid(int i) const
Get the ID of edge i of this object.
Definition: Geometry.cpp:112
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:391
std::vector< GeomFactorsSharedPtr > GeomFactorsVector
A vector of GeomFactor pointers.
Definition: GeomFactors.h:62
bool GlobalIdEquality(const std::shared_ptr< Geometry > &lhs, const std::shared_ptr< Geometry > &rhs)
Definition: Geometry.cpp:95
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
Definition: GeomFactors.h:60
std::shared_ptr< GeometryVector > GeometryVectorSharedPtr
Definition: Geometry.h:54
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:89
std::unordered_set< GeometrySharedPtr > GeometrySet
Definition: Geometry.h:53
std::shared_ptr< Curve > CurveSharedPtr
Definition: Curve.hpp:58
std::unordered_map< int, CurveSharedPtr > CurveMap
Definition: Curve.hpp:59
GeomType
Indicates the type of element geometry.
static CurveMap NullCurveMap
Definition: Geometry.h:67
std::shared_ptr< PointGeom > PointGeomSharedPtr
Definition: Geometry.h:57
std::shared_ptr< Geometry2D > Geometry2DSharedPtr
Definition: Geometry.h:62
GeomState
Indicates if the geometric information for an element has been populated.
std::shared_ptr< Geometry > GeometrySharedPtr
Definition: Geometry.h:51
std::shared_ptr< Geometry1D > Geometry1DSharedPtr
Definition: Geometry.h:61
std::vector< GeometrySharedPtr > GeometryVector
Definition: Geometry.h:52
std::shared_ptr< StdExpansion > StdExpansionSharedPtr
std::size_t hash_range(Iter first, Iter last)
Definition: HashUtils.hpp:64
double NekDouble
Unary function that constructs a hash of a Geometry object, based on the vertex IDs.
Definition: Geometry.h:260
std::size_t operator()(GeometrySharedPtr const &p) const
Definition: Geometry.h:261