Nektar++
Loading...
Searching...
No Matches
Public Member Functions | Protected Member Functions | Static Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
Nektar::SpatialDomains::Geometry Class Referenceabstract

Base class for shape geometry information. More...

#include <Geometry.h>

Inheritance diagram for Nektar::SpatialDomains::Geometry:
[legend]

Public Member Functions

 Geometry ()
 Default constructor.
 
 Geometry (int coordim)
 Constructor when supplied a coordinate dimension.
 
virtual ~Geometry ()=default
 
int GetCoordim () const
 Return the coordinate dimension of this object (i.e. the dimension of the space in which this object is embedded).
 
void SetCoordim (int coordim)
 Sets the coordinate dimension of this object (i.e. the dimension of the space in which this object is embedded).
 
GeomFactorsSharedPtr GetGeomFactors ()
 Get the geometric factors for this object, generating them if required.
 
GeomFactorsSharedPtr GetRefGeomFactors (const Array< OneD, const LibUtilities::BasisSharedPtr > &tbasis)
 
GeomFactorsSharedPtr GetMetricInfo ()
 Get the geometric factors for this object.
 
LibUtilities::ShapeType GetShapeType (void)
 Get the geometric shape type of this object.
 
int GetGlobalID (void) const
 Get the ID of this object.
 
void SetGlobalID (int globalid)
 Set the ID of this object.
 
int GetVid (int i) const
 Returns global id of vertex i of this object.
 
int GetEid (int i) const
 Get the ID of edge i of this object.
 
int GetFid (int i) const
 Get the ID of face i of this object.
 
int GetTid (int i) const
 Get the ID of trace i of this object.
 
PointGeomGetVertex (int i) const
 Returns vertex i of this object.
 
Geometry1DGetEdge (int i) const
 Returns edge i of this object.
 
Geometry2DGetFace (int i) const
 Returns face i of this object.
 
StdRegions::Orientation GetEorient (const int i) const
 Returns the orientation of edge i with respect to the ordering of edges in the standard element.
 
StdRegions::Orientation GetForient (const int i) const
 Returns the orientation of face i with respect to the ordering of faces in the standard element.
 
int GetNumVerts () const
 Get the number of vertices of this object.
 
int GetNumEdges () const
 Get the number of edges of this object.
 
int GetNumFaces () const
 Get the number of faces of this object.
 
int GetShapeDim () const
 Get the object's shape dimension.
 
StdRegions::StdExpansionSharedPtr GetXmap () const
 Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standard element to physical element.
 
const Array< OneD, const NekDouble > & GetCoeffs (const int i) const
 Return the coefficients of the transformation Geometry::m_xmap in coordinate direction i.
 
void FillGeom ()
 Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements.
 
std::array< NekDouble, 6 > GetBoundingBox ()
 Generates the bounding box for the element.
 
void ClearBoundingBox ()
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
 Determine whether an element contains a particular Cartesian coordinate \((x,y,z)\).
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol)
 Determine whether an element contains a particular Cartesian coordinate \((x,y,z)\).
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &dist)
 Determine whether an element contains a particular Cartesian coordinate \(\vec{x} = (x,y,z)\).
 
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 geometry object.
 
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.
 
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 elements, this check is unnecessary.
 
bool MinMaxCheck (const Array< OneD, const NekDouble > &gloCoord)
 Check if given global coord is within the BoundingBox of the element.
 
bool ClampLocCoords (Array< OneD, NekDouble > &locCoord, NekDouble tol=std::numeric_limits< NekDouble >::epsilon())
 Clamp local coords to be within standard regions [-1, 1]^dim.
 
NekDouble FindDistance (const Array< OneD, const NekDouble > &xs, Array< OneD, NekDouble > &xi)
 
int GetVertexEdgeMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given vertex.
 
int GetVertexFaceMap (int i, int j) const
 Returns the standard element face IDs that are connected to a given vertex.
 
int GetEdgeFaceMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given face.
 
int GetEdgeNormalToFaceVert (int i, int j) const
 Returns the standard lement edge IDs that are normal to a given face vertex.
 
int GetDir (const int i, const int j=0) const
 Returns the element coordinate direction corresponding to a given face coordinate direction.
 
void Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors.
 
void ResetNonRecursive (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object non-recursively: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors.
 
void Setup ()
 
void GenGeomFactors ()
 Handles generation of geometry factors.
 

Protected Member Functions

virtual int v_GetVid (int i) const
 Get the ID of vertex i of this object.
 
virtual PointGeomv_GetVertex (int i) const
 Returns vertex i of this object.
 
virtual Geometry1Dv_GetEdge (int i) const
 Returns edge i of this object.
 
virtual Geometry2Dv_GetFace (int i) const
 Returns face i of this object.
 
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.
 
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.
 
virtual int v_GetNumVerts () const
 Get the number of vertices of this object.
 
virtual int v_GetNumEdges () const
 Get the number of edges of this object.
 
virtual int v_GetNumFaces () const
 Get the number of faces of this object.
 
virtual int v_GetShapeDim () const
 Get the object's shape dimension.
 
virtual StdRegions::StdExpansionSharedPtr v_GetXmap () const
 Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standard element to physical element.
 
virtual void v_FillGeom ()
 Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements.
 
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 \(\vec{x} = (x,y,z)\).
 
virtual int v_AllLeftCheck (const Array< OneD, const NekDouble > &gloCoord)
 
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.
 
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 geometry object.
 
virtual NekDouble v_FindDistance (const Array< OneD, const NekDouble > &xs, Array< OneD, NekDouble > &xi)
 
virtual int v_GetVertexEdgeMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given vertex.
 
virtual int v_GetVertexFaceMap (int i, int j) const
 Returns the standard element face IDs that are connected to a given vertex.
 
virtual int v_GetEdgeFaceMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given face.
 
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.
 
virtual int v_GetDir (const int faceidx, const int facedir) const
 Returns the element coordinate direction corresponding to a given face coordinate direction.
 
virtual void v_Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors.
 
virtual void v_Setup ()
 
virtual void v_GenGeomFactors ()=0
 
void SetUpCoeffs (const int nCoeffs)
 Initialise the Geometry::m_coeffs array.
 
virtual void v_CalculateInverseIsoParam ()
 

Static Protected Member Functions

static GeomFactorsSharedPtr ValidateRegGeomFactor (GeomFactorsSharedPtr geomFactor)
 Check to see if a geometric factor has already been created that contains the same regular information.
 

Protected Attributes

int m_coordim
 Coordinate dimension of this geometry object.
 
GeomFactorsSharedPtr m_geomFactors
 Geometric factors.
 
GeomState m_geomFactorsState
 State of the geometric factors.
 
StdRegions::StdExpansionSharedPtr m_xmap
 \(\chi\) mapping containing isoparametric transformation.
 
GeomState m_state
 Enumeration to dictate whether coefficients are filled.
 
bool m_setupState
 Wether or not the setup routines have been run.
 
GeomType m_geomType
 Type of geometry.
 
LibUtilities::ShapeType m_shapeType
 Type of shape.
 
int m_globalID
 Global ID.
 
Array< OneD, Array< OneD, NekDouble > > m_coeffs
 Array containing expansion coefficients of m_xmap.
 
Array< OneD, NekDoublem_boundingBox
 Array containing bounding box.
 
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
 
Array< OneD, Array< OneD, NekDouble > > m_invIsoParam
 
int m_straightEdge
 

Static Protected Attributes

static GeomFactorsVector m_regGeomFactorsManager
 

Detailed Description

Base class for shape geometry information.

Definition at line 73 of file Geometry.h.

Constructor & Destructor Documentation

◆ Geometry() [1/2]

Nektar::SpatialDomains::Geometry::Geometry ( )

Default constructor.

Definition at line 50 of file Geometry.cpp.

54{
55}
bool m_setupState
Wether or not the setup routines have been run.
Definition Geometry.h:193
GeomState m_state
Enumeration to dictate whether coefficients are filled.
Definition Geometry.h:191
LibUtilities::ShapeType m_shapeType
Type of shape.
Definition Geometry.h:197
GeomState m_geomFactorsState
State of the geometric factors.
Definition Geometry.h:187
int m_coordim
Coordinate dimension of this geometry object.
Definition Geometry.h:183
@ eNotFilled
Geometric information has not been generated.

◆ Geometry() [2/2]

Nektar::SpatialDomains::Geometry::Geometry ( int  coordim)

Constructor when supplied a coordinate dimension.

Definition at line 60 of file Geometry.cpp.

◆ ~Geometry()

virtual Nektar::SpatialDomains::Geometry::~Geometry ( )
virtualdefault

Member Function Documentation

◆ ClampLocCoords()

bool Nektar::SpatialDomains::Geometry::ClampLocCoords ( Array< OneD, NekDouble > &  locCoord,
NekDouble  tol = std::numeric_limits<NekDouble>::epsilon() 
)

Clamp local coords to be within standard regions [-1, 1]^dim.

Parameters
LcoordsCorresponding local coordinates

Definition at line 536 of file Geometry.cpp.

537{
538 // Validation checks
539 ASSERTL1(locCoord.size() >= GetShapeDim(),
540 "Expects local coordinates to be same or "
541 "larger than shape dimension.");
542
543 // If out of range clamp locCoord to be within [-1,1]^dim
544 // since any larger value will be very oscillatory if
545 // called by 'returnNearestElmt' option in
546 // ExpList::GetExpIndex
547 bool clamp = false;
548 for (int i = 0; i < GetShapeDim(); ++i)
549 {
550 if (!std::isfinite(locCoord[i]))
551 {
552 locCoord[i] = 0.;
553 clamp = true;
554 }
555 else if (locCoord[i] < -(1. + tol))
556 {
557 locCoord[i] = -(1. + tol);
558 clamp = true;
559 }
560 else if (locCoord[i] > (1. + tol))
561 {
562 locCoord[i] = 1. + tol;
563 clamp = true;
564 }
565 }
566 return clamp;
567}
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
int GetShapeDim() const
Get the object's shape dimension.
Definition Geometry.h:430

References ASSERTL1, and GetShapeDim().

Referenced by Nektar::SpatialDomains::Geometry2D::NewtonIterationForLocCoord(), Nektar::SpatialDomains::Geometry3D::NewtonIterationForLocCoord(), Nektar::FieldUtils::ProcessWallNormalData::NewtonIterForLocCoordOnBndElmt(), v_ContainsPoint(), Nektar::SpatialDomains::Geometry2D::v_FindDistance(), Nektar::SpatialDomains::SegGeom::v_FindDistance(), Nektar::SpatialDomains::Geometry2D::v_GetLocCoords(), and Nektar::SpatialDomains::Geometry3D::v_GetLocCoords().

◆ ClearBoundingBox()

void Nektar::SpatialDomains::Geometry::ClearBoundingBox ( )

Definition at line 473 of file Geometry.cpp.

474{
475 m_boundingBox = {};
476}
Array< OneD, NekDouble > m_boundingBox
Array containing bounding box.
Definition Geometry.h:203

References m_boundingBox.

◆ ContainsPoint() [1/3]

bool Nektar::SpatialDomains::Geometry::ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
Array< OneD, NekDouble > &  locCoord,
NekDouble  tol 
)
inline

Determine whether an element contains a particular Cartesian coordinate \((x,y,z)\).

See also
Geometry::ContainsPoint

Definition at line 485 of file Geometry.h.

488{
489 NekDouble dist;
490 return v_ContainsPoint(gloCoord, locCoord, tol, dist);
491}
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:237

References v_ContainsPoint().

◆ ContainsPoint() [2/3]

bool Nektar::SpatialDomains::Geometry::ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
Array< OneD, NekDouble > &  locCoord,
NekDouble  tol,
NekDouble dist 
)
inline

Determine whether an element contains a particular Cartesian coordinate \(\vec{x} = (x,y,z)\).

For curvilinear and non-affine elements (i.e. where the Jacobian varies as a function of the standard element coordinates), this is a non-linear optimisation problem that requires the use of a Newton iteration. Note therefore that this can be an expensive operation.

The parameter tol which is by default 0, can be used to expand the coordinate range of the standard element from \([-1,1]^d\) to \([-1-\epsilon,1+\epsilon\) to handle challenging edge cases. The function also returns the local coordinates corresponding to gloCoord that can be used to speed up subsequent searches.

Parameters
gloCoordThe coordinate \( (x,y,z) \).
locCoordOn exit, this is the local coordinate \(\vec{\xi}\) such that \(\chi(\vec{\xi}) = \vec{x}\).
tolThe tolerance used to dictate the bounding box of the standard coordinates \(\vec{\xi}\).
distOn exit, returns the minimum distance between gloCoord and the quadrature points inside the element.
Returns
true if the coordinate gloCoord is contained in the element; false otherwise.
See also
Geometry::GetLocCoords.

Definition at line 521 of file Geometry.h.

524{
525 return v_ContainsPoint(gloCoord, locCoord, tol, dist);
526}

References v_ContainsPoint().

◆ ContainsPoint() [3/3]

bool Nektar::SpatialDomains::Geometry::ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
NekDouble  tol = 0.0 
)
inline

Determine whether an element contains a particular Cartesian coordinate \((x,y,z)\).

See also
Geometry::ContainsPoint

Definition at line 471 of file Geometry.h.

473{
474 Array<OneD, NekDouble> locCoord(GetCoordim(), 0.0);
475 NekDouble dist;
476 return v_ContainsPoint(gloCoord, locCoord, tol, dist);
477}
int GetCoordim() const
Return the coordinate dimension of this object (i.e. the dimension of the space in which this object ...
Definition Geometry.h:279

References GetCoordim(), and v_ContainsPoint().

Referenced by Geometry_ContainsPoint().

◆ FillGeom()

void Nektar::SpatialDomains::Geometry::FillGeom ( )
inline

Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements.

See also
v_FillGeom()

Definition at line 460 of file Geometry.h.

461{
462 v_FillGeom();
463}
virtual void v_FillGeom()
Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements.
Definition Geometry.cpp:363

References v_FillGeom().

Referenced by export_Geometry(), Nektar::LocalRegions::NodalTriExp::v_GetCoord(), Nektar::LocalRegions::HexExp::v_GetCoord(), Nektar::LocalRegions::PrismExp::v_GetCoord(), Nektar::LocalRegions::QuadExp::v_GetCoord(), Nektar::LocalRegions::SegExp::v_GetCoord(), Nektar::LocalRegions::TriExp::v_GetCoord(), and Nektar::LocalRegions::Expansion::v_GetCoords().

◆ FindDistance()

NekDouble Nektar::SpatialDomains::Geometry::FindDistance ( const Array< OneD, const NekDouble > &  xs,
Array< OneD, NekDouble > &  xi 
)
inline

Definition at line 564 of file Geometry.h.

566{
567 return v_FindDistance(xs, xi);
568}
virtual NekDouble v_FindDistance(const Array< OneD, const NekDouble > &xs, Array< OneD, NekDouble > &xi)
Definition Geometry.cpp:272

References v_FindDistance().

◆ GenGeomFactors()

void Nektar::SpatialDomains::Geometry::GenGeomFactors ( )
inline

Handles generation of geometry factors.

Generate the geometric factors (i.e. derivatives of \(\chi\)) and related metrics.

See also
SpatialDomains::GeomFactors

Definition at line 696 of file Geometry.h.

697{
698 return v_GenGeomFactors();
699}

References v_GenGeomFactors().

Referenced by GetGeomFactors().

◆ GetBoundingBox()

std::array< NekDouble, 6 > Nektar::SpatialDomains::Geometry::GetBoundingBox ( )

Generates the bounding box for the element.

For regular elements, the vertices are sufficient to define the extent of the bounding box. For non-regular elements, the extremes of the quadrature point coordinates are used. A 10% margin is added around this computed region to account for convex hull elements where the true extent of the element may extend slightly beyond the quadrature points.

Definition at line 399 of file Geometry.cpp.

400{
401 if (m_boundingBox.size() == 6)
402 {
403 return {{m_boundingBox[0], m_boundingBox[1], m_boundingBox[2],
405 }
406 // NekDouble minx, miny, minz, maxx, maxy, maxz;
407 Array<OneD, NekDouble> min(3), max(3);
408
409 // Always get vertexes min/max
410 PointGeom *p = GetVertex(0);
411 Array<OneD, NekDouble> x(3, 0.0);
412 p->GetCoords(x[0], x[1], x[2]);
413 for (int j = 0; j < 3; ++j)
414 {
415 min[j] = x[j];
416 max[j] = x[j];
417 }
418 for (int i = 1; i < GetNumVerts(); ++i)
419 {
420 p = GetVertex(i);
421 p->GetCoords(x[0], x[1], x[2]);
422 for (int j = 0; j < 3; ++j)
423 {
424 min[j] = (x[j] < min[j] ? x[j] : min[j]);
425 max[j] = (x[j] > max[j] ? x[j] : max[j]);
426 }
427 }
428 // If element is deformed loop over quadrature points
430 if (GetGeomFactors()->GetGtype() != eRegular)
431 {
432 marginFactor = 0.1;
433 const int nq = GetXmap()->GetTotPoints();
434 Array<OneD, Array<OneD, NekDouble>> xvec(3);
435 for (int j = 0; j < 3; ++j)
436 {
437 xvec[j] = Array<OneD, NekDouble>(nq, 0.0);
438 }
439 for (int j = 0; j < GetCoordim(); ++j)
440 {
441 GetXmap()->BwdTrans(m_coeffs[j], xvec[j]);
442 }
443 for (int j = 0; j < 3; ++j)
444 {
445 for (int i = 0; i < nq; ++i)
446 {
447 min[j] = (xvec[j][i] < min[j] ? xvec[j][i] : min[j]);
448 max[j] = (xvec[j][i] > max[j] ? xvec[j][i] : max[j]);
449 }
450 }
451 }
452 // Add margin to bounding box, in order to
453 // return the nearest element
454 for (int j = 0; j < 3; ++j)
455 {
456 NekDouble margin =
457 marginFactor * (max[j] - min[j]) + NekConstants::kFindDistanceMin;
458 min[j] -= margin;
459 max[j] += margin;
460 }
461
462 // save bounding box
463 m_boundingBox = Array<OneD, NekDouble>(6);
464 for (int j = 0; j < 3; ++j)
465 {
466 m_boundingBox[j] = min[j];
467 m_boundingBox[j + 3] = max[j];
468 }
469 // Return bounding box
470 return {{min[0], min[1], min[2], max[0], max[1], max[2]}};
471}
PointGeom * GetVertex(int i) const
Returns vertex i of this object.
Definition Geometry.h:361
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
Definition Geometry.h:201
StdRegions::StdExpansionSharedPtr GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition Geometry.h:439
int GetNumVerts() const
Get the number of vertices of this object.
Definition Geometry.h:403
GeomFactorsSharedPtr GetGeomFactors()
Get the geometric factors for this object, generating them if required.
Definition Geometry.h:297
static const NekDouble kGeomFactorsTol
static const NekDouble kFindDistanceMin
@ eRegular
Geometry is straight-sided with constant geometric factors.
std::vector< double > p(NPUPPER)
scalarT< T > max(scalarT< T > lhs, scalarT< T > rhs)
Definition scalar.hpp:305
scalarT< T > min(scalarT< T > lhs, scalarT< T > rhs)
Definition scalar.hpp:300

References Nektar::SpatialDomains::eRegular, GetCoordim(), GetGeomFactors(), GetNumVerts(), GetVertex(), GetXmap(), Nektar::NekConstants::kFindDistanceMin, Nektar::NekConstants::kGeomFactorsTol, m_boundingBox, m_coeffs, tinysimd::max(), and tinysimd::min().

Referenced by Nektar::SpatialDomains::MeshGraph::GeomRTree::InsertGeom(), and MinMaxCheck().

◆ GetCoeffs()

const Array< OneD, const NekDouble > & Nektar::SpatialDomains::Geometry::GetCoeffs ( const int  i) const
inline

◆ GetCoord()

NekDouble Nektar::SpatialDomains::Geometry::GetCoord ( const int  i,
const Array< OneD, const NekDouble > &  Lcoord 
)
inline

Given local collapsed coordinate Lcoord, return the value of physical coordinate in direction i.

Definition at line 558 of file Geometry.h.

560{
561 return v_GetCoord(i, Lcoord);
562}
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:339

References v_GetCoord().

Referenced by Nektar::SpatialDomains::Geometry2D::v_FindDistance(), Nektar::SpatialDomains::SegGeom::v_FindDistance(), Nektar::LocalRegions::NodalTriExp::v_GetCoord(), Nektar::LocalRegions::HexExp::v_GetCoord(), Nektar::LocalRegions::PrismExp::v_GetCoord(), Nektar::LocalRegions::PyrExp::v_GetCoord(), Nektar::LocalRegions::QuadExp::v_GetCoord(), Nektar::LocalRegions::SegExp::v_GetCoord(), Nektar::LocalRegions::TetExp::v_GetCoord(), and Nektar::LocalRegions::TriExp::v_GetCoord().

◆ GetCoordim()

int Nektar::SpatialDomains::Geometry::GetCoordim ( ) const
inline

Return the coordinate dimension of this object (i.e. the dimension of the space in which this object is embedded).

Definition at line 279 of file Geometry.h.

280{
281 return m_coordim;
282}

References m_coordim.

Referenced by Nektar::SpatialDomains::MeshGraph::CheckRange(), ContainsPoint(), Nektar::LocalRegions::SegExp::CreateMatrix(), Nektar::QuadCollectionTests::CreateSegGeom(), Nektar::TriCollectionTests::CreateSegGeom(), export_Geometry(), GetBoundingBox(), Nektar::LocalRegions::NodalTriExp::v_AlignVectorToCollapsedDir(), Nektar::LocalRegions::QuadExp::v_AlignVectorToCollapsedDir(), Nektar::LocalRegions::TriExp::v_AlignVectorToCollapsedDir(), Nektar::SpatialDomains::SegGeom::v_FindDistance(), Nektar::LocalRegions::NodalTriExp::v_GetCoord(), Nektar::LocalRegions::HexExp::v_GetCoord(), Nektar::LocalRegions::PrismExp::v_GetCoord(), Nektar::LocalRegions::PyrExp::v_GetCoord(), Nektar::LocalRegions::QuadExp::v_GetCoord(), Nektar::LocalRegions::SegExp::v_GetCoord(), Nektar::LocalRegions::TetExp::v_GetCoord(), Nektar::LocalRegions::TriExp::v_GetCoord(), Nektar::LocalRegions::Expansion::v_GetCoordim(), Nektar::LocalRegions::Expansion::v_GetCoords(), Nektar::LocalRegions::SegExp::v_HelmholtzMatrixOp(), Nektar::LocalRegions::SegExp::v_IProductWRTDerivBase(), Nektar::LocalRegions::QuadExp::v_IProductWRTDerivBase_SumFac(), Nektar::LocalRegions::SegExp::v_LaplacianMatrixOp(), Nektar::LocalRegions::SegExp::v_PhysDeriv_n(), Nektar::LocalRegions::SegExp::v_PhysDeriv_s(), Nektar::LocalRegions::QuadExp::v_PhysDirectionalDeriv(), and Nektar::LocalRegions::TriExp::v_PhysDirectionalDeriv().

◆ GetDir()

int Nektar::SpatialDomains::Geometry::GetDir ( const int  i,
const int  j = 0 
) const
inline

Returns the element coordinate direction corresponding to a given face coordinate direction.

Definition at line 661 of file Geometry.h.

662{
663 return v_GetDir(faceidx, facedir);
664}
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:328

References v_GetDir().

Referenced by Nektar::LocalRegions::Expansion3D::v_GetTracePhysVals(), Nektar::LocalRegions::Expansion2D::v_TraceNormLen(), and Nektar::LocalRegions::Expansion3D::v_TraceNormLen().

◆ GetEdge()

Geometry1D * Nektar::SpatialDomains::Geometry::GetEdge ( int  i) const
inline

◆ GetEdgeFaceMap()

int Nektar::SpatialDomains::Geometry::GetEdgeFaceMap ( int  i,
int  j 
) const
inline

Returns the standard element edge IDs that are connected to a given face.

For example, on a prism, edge 0 is connnected to faces 0 and 1; GetEdgeFaceMap(0,j) would therefore return the values 0 and 1 respectively. We assume that j runs between 0 and 1 inclusive, since every face is connected to precisely two faces for all 3D elements.

This function is used in the construction of the low-energy preconditioner.

Parameters
iThe edge to query connectivity for.
jThe local face index between 0 and 1 connected to this element.
See also
MultiRegions::PreconditionerLowEnergy

Definition at line 630 of file Geometry.h.

631{
632 return v_GetEdgeFaceMap(i, j);
633}
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:306

References v_GetEdgeFaceMap().

Referenced by Nektar::LocalRegions::Expansion3D::v_BuildTransformationMatrix().

◆ GetEdgeNormalToFaceVert()

int Nektar::SpatialDomains::Geometry::GetEdgeNormalToFaceVert ( int  i,
int  j 
) const
inline

Returns the standard lement edge IDs that are normal to a given face vertex.

For example, on a hexahedron, on face 0 at vertices 0,1,2,3 the edges normal to that face are 4,5,6,7, ; so GetEdgeNormalToFaceVert(0,j) would therefore return the values 4, 5, 6 and 7 respectively. We assume that j runs between 0 and 3 inclusive on a quadrilateral face and between 0 and 2 inclusive on a triangular face.

This is used to help set up a length scale normal to an face

Parameters
iThe face to query for the normal edge
jThe local vertex index between 0 and nverts on this face

Definition at line 652 of file Geometry.h.

653{
654 return v_GetEdgeNormalToFaceVert(i, j);
655}
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:317

References v_GetEdgeNormalToFaceVert().

Referenced by Nektar::LocalRegions::Expansion3D::v_TraceNormLen().

◆ GetEid()

int Nektar::SpatialDomains::Geometry::GetEid ( int  i) const

◆ GetEorient()

StdRegions::Orientation Nektar::SpatialDomains::Geometry::GetEorient ( const int  i) const
inline

Returns the orientation of edge i with respect to the ordering of edges in the standard element.

Definition at line 386 of file Geometry.h.

387{
388 return v_GetEorient(i);
389}
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:166

References v_GetEorient().

Referenced by export_Geometry(), Nektar::MultiRegions::DisContField::FindPeriodicTraces(), Nektar::LocalRegions::Expansion3D::GetEdgeInverseBoundaryMap(), Nektar::LocalRegions::Expansion2D::GetTraceInverseBoundaryMap(), Nektar::GlobalMapping::UpdateGeometry(), Nektar::LocalRegions::QuadExp::v_GetTraceOrient(), and Nektar::LocalRegions::TriExp::v_GetTraceOrient().

◆ GetFace()

Geometry2D * Nektar::SpatialDomains::Geometry::GetFace ( int  i) const
inline

Returns face i of this object.

Definition at line 377 of file Geometry.h.

378{
379 return v_GetFace(i);
380}
virtual Geometry2D * v_GetFace(int i) const
Returns face i of this object.
Definition Geometry.cpp:146

References v_GetFace().

Referenced by export_Geometry(), GetFid(), Nektar::SpatialDomains::MeshGraph::PopulateFaceToElMap(), Nektar::GlobalMapping::UpdateGeometry(), Nektar::LocalRegions::Expansion3D::v_GenTraceExp(), and Nektar::LocalRegions::Expansion3D::v_TraceNormLen().

◆ GetFid()

int Nektar::SpatialDomains::Geometry::GetFid ( int  i) const

◆ GetForient()

StdRegions::Orientation Nektar::SpatialDomains::Geometry::GetForient ( const int  i) const
inline

Returns the orientation of face i with respect to the ordering of faces in the standard element.

Definition at line 395 of file Geometry.h.

396{
397 return v_GetForient(i);
398}
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:177

References v_GetForient().

Referenced by export_Geometry(), and Nektar::LocalRegions::Expansion3D::v_GetTraceOrient().

◆ GetGeomFactors()

GeomFactorsSharedPtr Nektar::SpatialDomains::Geometry::GetGeomFactors ( )
inline

Get the geometric factors for this object, generating them if required.

Definition at line 297 of file Geometry.h.

298{
301}
void GenGeomFactors()
Handles generation of geometry factors.
Definition Geometry.h:696
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
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
Definition Geometry.h:185

References GenGeomFactors(), m_geomFactors, and ValidateRegGeomFactor().

Referenced by Geometry_GenGeomFactors(), Geometry_IsValid(), GetBoundingBox(), Nektar::FieldUtils::ProcessQualityMetric::GetQ(), and Nektar::LocalRegions::Expansion::Reset().

◆ GetGlobalID()

int Nektar::SpatialDomains::Geometry::GetGlobalID ( void  ) const
inline

Get the ID of this object.

Definition at line 322 of file Geometry.h.

323{
324 return m_globalID;
325}

References m_globalID.

Referenced by Nektar::SolverUtils::CouplingCwipi::AddElementsToMesh(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), CheckTetRotation(), Nektar::LocalRegions::Expansion::Expansion(), export_Geometry(), GenerateMapEidsv1v2(), Nektar::SpatialDomains::SegGeom::GenerateOneSpaceDimGeom(), Nektar::SpatialDomains::SegGeom::GetEdgeOrientation(), GetEid(), Nektar::SpatialDomains::MeshGraph::GetElementsFromEdge(), Nektar::SpatialDomains::MeshGraph::GetElementsFromFace(), Nektar::SpatialDomains::MeshGraph::GetExpansionInfo(), GetFid(), Nektar::FieldUtils::ProcessQualityMetric::GetQ(), Nektar::SpatialDomains::GlobalIdEquality(), Nektar::SpatialDomains::MeshGraph::GeomRTree::InsertGeom(), main(), Nektar::SpatialDomains::Geometry3D::NewtonIterationForLocCoord(), Nektar::SpatialDomains::Geometry2D::NewtonIterationForLocCoord(), Nektar::SpatialDomains::Geometry3D::NewtonIterationForLocCoord(), Nektar::FieldUtils::ProcessWallNormalData::NewtonIterForLocCoordOnBndElmt(), Nektar::SpatialDomains::MeshGraph::PopulateFaceToElMap(), Nektar::SpatialDomains::MeshGraph::PRefinementElmts(), Nektar::PulseWaveSystem::SetUpDomainInterfaceBCs(), Nektar::SpatialDomains::HexGeom::SetUpEdgeOrientation(), Nektar::SpatialDomains::PrismGeom::SetUpEdgeOrientation(), Nektar::SpatialDomains::PyrGeom::SetUpEdgeOrientation(), Nektar::SpatialDomains::TetGeom::SetUpEdgeOrientation(), Nektar::SpatialDomains::HexGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::PrismGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::PyrGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::TetGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::HexGeom::SetUpLocalEdges(), Nektar::SpatialDomains::PrismGeom::SetUpLocalEdges(), Nektar::SpatialDomains::PyrGeom::SetUpLocalEdges(), Nektar::SpatialDomains::TetGeom::SetUpLocalEdges(), Nektar::SpatialDomains::HexGeom::SetUpLocalVertices(), Nektar::SpatialDomains::PrismGeom::SetUpLocalVertices(), Nektar::SpatialDomains::PyrGeom::SetUpLocalVertices(), Nektar::SpatialDomains::TetGeom::SetUpLocalVertices(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalDGMap(), Nektar::SpatialDomains::SortByGlobalId(), Nektar::GlobalMapping::UpdateGeometry(), Nektar::LocalRegions::Expansion2D::v_AddEdgeNormBoundaryInt(), Nektar::LocalRegions::Expansion3D::v_AddFaceNormBoundaryInt(), Nektar::SpatialDomains::QuadGeom::v_AllLeftCheck(), Nektar::SpatialDomains::TriGeom::v_AllLeftCheck(), v_GetVid(), Nektar::SpatialDomains::WriteVert(), and Nektar::SpatialDomains::ZoneBase::ZoneBase().

◆ GetLocCoords()

NekDouble Nektar::SpatialDomains::Geometry::GetLocCoords ( const Array< OneD, const NekDouble > &  coords,
Array< OneD, NekDouble > &  Lcoords 
)
inline

Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this geometry object.

For curvilinear and non-affine elements (i.e. where the Jacobian varies as a function of the standard element coordinates), this is a non-linear optimisation problem that requires the use of a Newton iteration. Note therefore that this can be an expensive operation.

Note that, clearly, the provided Cartesian coordinate lie outside the element. The function therefore returns the minimum distance from some position in the element to . Lcoords will also be constrained to fit within the range \([-1,1]^d\) where \( d \) is the dimension of the element.

Parameters
coordsInput Cartesian global coordinates
LcoordsCorresponding local coordinates
Returns
Distance between obtained coordinates and provided ones.

Definition at line 548 of file Geometry.h.

550{
551 return v_GetLocCoords(coords, Lcoords);
552}
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:351

References v_GetLocCoords().

Referenced by v_ContainsPoint(), Nektar::SpatialDomains::Geometry2D::v_FindDistance(), Nektar::SpatialDomains::SegGeom::v_FindDistance(), Nektar::LocalRegions::HexExp::v_PhysEvalFirstDeriv(), Nektar::LocalRegions::PrismExp::v_PhysEvalFirstDeriv(), Nektar::LocalRegions::PyrExp::v_PhysEvalFirstDeriv(), Nektar::LocalRegions::QuadExp::v_PhysEvalFirstDeriv(), Nektar::LocalRegions::SegExp::v_PhysEvalFirstDeriv(), Nektar::LocalRegions::TetExp::v_PhysEvalFirstDeriv(), Nektar::LocalRegions::TriExp::v_PhysEvalFirstDeriv(), Nektar::LocalRegions::SegExp::v_PhysEvalFirstSecondDeriv(), Nektar::LocalRegions::NodalTriExp::v_PhysEvaluate(), Nektar::LocalRegions::PrismExp::v_PhysEvaluate(), Nektar::LocalRegions::PyrExp::v_PhysEvaluate(), Nektar::LocalRegions::QuadExp::v_PhysEvaluate(), Nektar::LocalRegions::SegExp::v_PhysEvaluate(), Nektar::LocalRegions::TriExp::v_PhysEvaluate(), Nektar::LocalRegions::HexExp::v_PhysEvaluate(), and Nektar::LocalRegions::TetExp::v_PhysEvaluate().

◆ GetMetricInfo()

GeomFactorsSharedPtr Nektar::SpatialDomains::Geometry::GetMetricInfo ( )
inline

◆ GetNumEdges()

int Nektar::SpatialDomains::Geometry::GetNumEdges ( ) const
inline

Get the number of edges of this object.

Definition at line 411 of file Geometry.h.

412{
413 return v_GetNumEdges();
414}
virtual int v_GetNumEdges() const
Get the number of edges of this object.
Definition Geometry.cpp:188

References v_GetNumEdges().

Referenced by Nektar::SpatialDomains::MeshGraph::CheckRange(), export_Geometry(), Nektar::SpatialDomains::MeshGraph::GetElementsFromEdge(), and Nektar::SpatialDomains::MeshGraphIOXml::WriteXMLGeometry().

◆ GetNumFaces()

int Nektar::SpatialDomains::Geometry::GetNumFaces ( ) const
inline

Get the number of faces of this object.

Definition at line 419 of file Geometry.h.

420{
421 return v_GetNumFaces();
422}
virtual int v_GetNumFaces() const
Get the number of faces of this object.
Definition Geometry.cpp:196

References v_GetNumFaces().

Referenced by Nektar::SpatialDomains::MeshGraph::CheckRange(), export_Geometry(), and Nektar::SpatialDomains::MeshGraphIOXml::WriteXMLGeometry().

◆ GetNumVerts()

int Nektar::SpatialDomains::Geometry::GetNumVerts ( ) const
inline

◆ GetRefGeomFactors()

GeomFactorsSharedPtr Nektar::SpatialDomains::Geometry::GetRefGeomFactors ( const Array< OneD, const LibUtilities::BasisSharedPtr > &  tbasis)

◆ GetShapeDim()

int Nektar::SpatialDomains::Geometry::GetShapeDim ( ) const
inline

Get the object's shape dimension.

For example, a segment is one dimensional and quadrilateral is two dimensional.

Definition at line 430 of file Geometry.h.

431{
432 return v_GetShapeDim();
433}
virtual int v_GetShapeDim() const
Get the object's shape dimension.
Definition Geometry.cpp:204

References v_GetShapeDim().

Referenced by ClampLocCoords(), export_Geometry(), Nektar::SpatialDomains::MeshGraph::GetCompositeString(), Nektar::SpatialDomains::MeshGraphIO::GetCompositeString(), GetTid(), and v_ContainsPoint().

◆ GetShapeType()

LibUtilities::ShapeType Nektar::SpatialDomains::Geometry::GetShapeType ( void  )
inline

◆ GetTid()

int Nektar::SpatialDomains::Geometry::GetTid ( int  i) const
inline

Get the ID of trace i of this object.

The trace element is the facet one dimension lower than the object; for example, a quadrilateral has four trace segments forming its boundary.

Definition at line 341 of file Geometry.h.

342{
343 const int nDim = GetShapeDim();
344 return nDim == 1 ? GetVid(i)
345 : nDim == 2 ? GetEid(i)
346 : nDim == 3 ? GetFid(i)
347 : 0;
348}
int GetVid(int i) const
Returns global id of vertex i of this object.
Definition Geometry.h:353
int GetFid(int i) const
Get the ID of face i of this object.
Definition Geometry.cpp:118
int GetEid(int i) const
Get the ID of edge i of this object.
Definition Geometry.cpp:110

References GetEid(), GetFid(), GetShapeDim(), and GetVid().

Referenced by export_Geometry().

◆ GetVertex()

PointGeom * Nektar::SpatialDomains::Geometry::GetVertex ( int  i) const
inline

Returns vertex i of this object.

Definition at line 361 of file Geometry.h.

362{
363 return v_GetVertex(i);
364}
virtual PointGeom * v_GetVertex(int i) const
Returns vertex i of this object.
Definition Geometry.cpp:126

References v_GetVertex().

Referenced by Nektar::SolverUtils::CouplingCwipi::AddElementsToMesh(), Nektar::SpatialDomains::MeshGraph::CheckRange(), Nektar::SpatialDomains::MeshGraph::CheckRange(), export_Geometry(), Nektar::MultiRegions::DisContField::FindPeriodicTraces(), GenerateMapEidsv1v2(), GetBoundingBox(), Nektar::SpatialDomains::SegGeom::GetEdgeOrientation(), GetNewVertexLocation(), Nektar::MappingIdealToRef(), Nektar::FieldUtils::MappingIdealToRef(), Nektar::SolverUtils::DriverArnoldi::MaskInit(), Nektar::SpatialDomains::MeshGraph::PRefinementElmts(), Nektar::VariableConverter::SetElmtMinHP(), Nektar::PulseWaveSystem::SetUpDomainInterfaceBCs(), Nektar::SpatialDomains::HexGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::PrismGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::PyrGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::TetGeom::SetUpFaceOrientation(), Nektar::GlobalMapping::UpdateGeometry(), v_GetVid(), Nektar::DiffusionLDGNS::v_InitObject(), Nektar::FieldUtils::ProcessPowerSpectrum::v_Process(), Nektar::LocalRegions::Expansion1D::v_TraceNormLen(), Nektar::LocalRegions::Expansion2D::v_TraceNormLen(), Nektar::LocalRegions::Expansion3D::v_TraceNormLen(), Nektar::FieldUtils::ProcessForceDecompose::VolumeIntegrateForce(), and Nektar::SpatialDomains::ZoneBase::ZoneBase().

◆ GetVertexEdgeMap()

int Nektar::SpatialDomains::Geometry::GetVertexEdgeMap ( int  i,
int  j 
) const
inline

Returns the standard element edge IDs that are connected to a given vertex.

For example, on a prism, vertex 0 is connnected to edges 0, 3, and 4; GetVertexEdgeMap(0,j) would therefore return the values 0, 1 and 4 respectively. We assume that j runs between 0 and 2 inclusive, which is true for every 3D element asides from the pyramid.

This function is used in the construction of the low-energy preconditioner.

Parameters
iThe vertex to query connectivity for.
jThe local edge index between 0 and 2 connected to this element.
Todo:
Expand to work with pyramid elements.
See also
MultiRegions::PreconditionerLowEnergy

Definition at line 587 of file Geometry.h.

588{
589 return v_GetVertexEdgeMap(i, j);
590}
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:284

References v_GetVertexEdgeMap().

Referenced by Nektar::LocalRegions::Expansion3D::v_BuildTransformationMatrix().

◆ GetVertexFaceMap()

int Nektar::SpatialDomains::Geometry::GetVertexFaceMap ( int  i,
int  j 
) const
inline

Returns the standard element face IDs that are connected to a given vertex.

For example, on a hexahedron, vertex 0 is connnected to faces 0, 1, and 4; GetVertexFaceMap(0,j) would therefore return the values 0, 1 and 4 respectively. We assume that j runs between 0 and 2 inclusive, which is true for every 3D element asides from the pyramid.

This is used in the construction of the low-energy preconditioner.

Parameters
iThe vertex to query connectivity for.
jThe local face index between 0 and 2 connected to this element.
Todo:
Expand to work with pyramid elements.
See also
MultiRegions::PreconditionerLowEnergy

Definition at line 609 of file Geometry.h.

610{
611 return v_GetVertexFaceMap(i, j);
612}
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:295

References v_GetVertexFaceMap().

Referenced by Nektar::LocalRegions::Expansion3D::v_BuildTransformationMatrix().

◆ GetVid()

int Nektar::SpatialDomains::Geometry::GetVid ( int  i) const
inline

◆ GetXmap()

StdRegions::StdExpansionSharedPtr Nektar::SpatialDomains::Geometry::GetXmap ( ) const
inline

Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standard element to physical element.

Definition at line 439 of file Geometry.h.

440{
441 return v_GetXmap();
442}
virtual StdRegions::StdExpansionSharedPtr v_GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition Geometry.cpp:226

References v_GetXmap().

Referenced by Nektar::FieldUtils::ProcessWallNormalData::BisectionForLocCoordOnBndElmt(), Nektar::FieldUtils::ProcessWallNormalData::BndElmtContainsPoint(), export_Geometry(), Nektar::FieldUtils::ProcessBodyFittedVelocity::GenPntwiseBodyFittedCoordSys(), GetBoundingBox(), Nektar::FieldUtils::ProcessWallNormalData::GetNormals(), Nektar::FieldUtils::ProcessQualityMetric::GetQ(), Nektar::FieldUtils::ProcessBodyFittedVelocity::LocCoordForNearestPntOnBndElmt(), Nektar::FieldUtils::ProcessBodyFittedVelocity::LocCoordForNearestPntOnBndElmt_2D(), Nektar::FieldUtils::ProcessWallNormalData::NewtonIterForLocCoordOnBndElmt(), Nektar::SpatialDomains::HexGeom::SetUpXmap(), Nektar::SpatialDomains::PrismGeom::SetUpXmap(), Nektar::SpatialDomains::PyrGeom::SetUpXmap(), Nektar::SpatialDomains::QuadGeom::SetUpXmap(), Nektar::SpatialDomains::TetGeom::SetUpXmap(), Nektar::SpatialDomains::TriGeom::SetUpXmap(), Nektar::SpatialDomains::QuadGeom::v_AllLeftCheck(), Nektar::SpatialDomains::TriGeom::v_AllLeftCheck(), Nektar::SpatialDomains::HexGeom::v_FillGeom(), Nektar::SpatialDomains::PrismGeom::v_FillGeom(), Nektar::SpatialDomains::PyrGeom::v_FillGeom(), Nektar::SpatialDomains::QuadGeom::v_FillGeom(), Nektar::SpatialDomains::TetGeom::v_FillGeom(), Nektar::SpatialDomains::TriGeom::v_FillGeom(), Nektar::LocalRegions::Expansion::v_GetCoords(), and Nektar::FieldUtils::ProcessWallNormalData::v_Process().

◆ MinMaxCheck()

bool Nektar::SpatialDomains::Geometry::MinMaxCheck ( const Array< OneD, const NekDouble > &  gloCoord)

Check if given global coord is within the BoundingBox of the element.

Parameters
coordsInput Cartesian global coordinates
Returns
True if within distance or False otherwise.

Definition at line 513 of file Geometry.cpp.

514{
515 // Validation checks
516 ASSERTL1(gloCoord.size() >= m_coordim,
517 "Expects number of global coordinates supplied to be greater than "
518 "or equal to the mesh dimension.");
519
520 std::array<NekDouble, 6> minMax = GetBoundingBox();
521 for (int i = 0; i < m_coordim; ++i)
522 {
523 if ((gloCoord[i] < minMax[i]) || (gloCoord[i] > minMax[i + 3]))
524 {
525 return false;
526 }
527 }
528 return true;
529}
std::array< NekDouble, 6 > GetBoundingBox()
Generates the bounding box for the element.
Definition Geometry.cpp:399

References ASSERTL1, GetBoundingBox(), and m_coordim.

Referenced by PreliminaryCheck().

◆ PreliminaryCheck()

int Nektar::SpatialDomains::Geometry::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 elements, this check is unnecessary.

Parameters
coordsInput Cartesian global coordinates
Returns
1 is inside of the element. 0 maybe inside -1 outside of the element

Definition at line 488 of file Geometry.cpp.

489{
490 // bounding box check
491 if (!MinMaxCheck(gloCoord))
492 {
493 return -1;
494 }
495
496 // regular element check
497 if (GetMetricInfo()->GetGtype() == eRegular)
498 {
499 return 0;
500 }
501
502 // All left check for straight edges/plane surfaces
503 return v_AllLeftCheck(gloCoord);
504}
virtual int v_AllLeftCheck(const Array< OneD, const NekDouble > &gloCoord)
Definition Geometry.cpp:211
GeomFactorsSharedPtr GetMetricInfo()
Get the geometric factors for this object.
Definition Geometry.h:306
bool MinMaxCheck(const Array< OneD, const NekDouble > &gloCoord)
Check if given global coord is within the BoundingBox of the element.
Definition Geometry.cpp:513

References Nektar::SpatialDomains::eRegular, GetMetricInfo(), MinMaxCheck(), and v_AllLeftCheck().

Referenced by v_ContainsPoint().

◆ Reset()

void Nektar::SpatialDomains::Geometry::Reset ( CurveMap curvedEdges,
CurveMap curvedFaces 
)
inline

Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors.

Definition at line 670 of file Geometry.h.

671{
672 v_Reset(curvedEdges, curvedFaces);
673}
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:372

References v_Reset().

◆ ResetNonRecursive()

void Nektar::SpatialDomains::Geometry::ResetNonRecursive ( CurveMap curvedEdges,
CurveMap curvedFaces 
)
inline

Reset this geometry object non-recursively: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors.

Definition at line 679 of file Geometry.h.

681{
682 Geometry::v_Reset(curvedEdges, curvedFaces);
683}

References v_Reset().

◆ SetCoordim()

void Nektar::SpatialDomains::Geometry::SetCoordim ( int  coordim)
inline

Sets the coordinate dimension of this object (i.e. the dimension of the space in which this object is embedded).

Definition at line 288 of file Geometry.h.

289{
290 m_coordim = dim;
291}

References m_coordim.

◆ SetGlobalID()

void Nektar::SpatialDomains::Geometry::SetGlobalID ( int  globalid)
inline

Set the ID of this object.

Definition at line 330 of file Geometry.h.

331{
332 m_globalID = globalid;
333}

References m_globalID.

Referenced by export_Geometry().

◆ Setup()

void Nektar::SpatialDomains::Geometry::Setup ( )
inline

Definition at line 685 of file Geometry.h.

686{
687 v_Setup();
688}

References v_Setup().

Referenced by export_Geometry().

◆ SetUpCoeffs()

void Nektar::SpatialDomains::Geometry::SetUpCoeffs ( const int  nCoeffs)
inlineprotected

◆ v_AllLeftCheck()

int Nektar::SpatialDomains::Geometry::v_AllLeftCheck ( const Array< OneD, const NekDouble > &  gloCoord)
protectedvirtual

Reimplemented in Nektar::SpatialDomains::Geometry3D, Nektar::SpatialDomains::QuadGeom, and Nektar::SpatialDomains::TriGeom.

Definition at line 211 of file Geometry.cpp.

213{
214 return 0;
215}

Referenced by PreliminaryCheck().

◆ v_CalculateInverseIsoParam()

void Nektar::SpatialDomains::Geometry::v_CalculateInverseIsoParam ( )
protectedvirtual

Reimplemented in Nektar::SpatialDomains::Geometry2D, and Nektar::SpatialDomains::Geometry3D.

Definition at line 217 of file Geometry.cpp.

218{
220 "This function is only valid for shape type geometries");
221}
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...

References Nektar::ErrorUtil::efatal, and NEKERROR.

◆ v_ContainsPoint()

bool Nektar::SpatialDomains::Geometry::v_ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
Array< OneD, NekDouble > &  locCoord,
NekDouble  tol,
NekDouble dist 
)
protectedvirtual

Determine whether an element contains a particular Cartesian coordinate \(\vec{x} = (x,y,z)\).

For curvilinear and non-affine elements (i.e. where the Jacobian varies as a function of the standard element coordinates), this is a non-linear optimisation problem that requires the use of a Newton iteration. Note therefore that this can be an expensive operation.

The parameter tol which is by default 0, can be used to expand the coordinate range of the standard element from \([-1,1]^d\) to \([-1-\epsilon,1+\epsilon\) to handle challenging edge cases. The function also returns the local coordinates corresponding to gloCoord that can be used to speed up subsequent searches.

Parameters
gloCoordThe coordinate \( (x,y,z) \).
locCoordOn exit, this is the local coordinate \(\vec{\xi}\) such that \(\chi(\vec{\xi}) = \vec{x}\).
tolThe tolerance used to dictate the bounding box of the standard coordinates \(\vec{\xi}\).
distOn exit, returns the minimum distance between gloCoord and the quadrature points inside the element.
Returns
true if the coordinate gloCoord is contained in the element; false otherwise.
See also
Geometry::GetLocCoords. dist is assigned value for curved elements

Reimplemented in Nektar::SpatialDomains::Geometry0D.

Definition at line 237 of file Geometry.cpp.

240{
241 int inside = PreliminaryCheck(gloCoord);
242 if (inside == -1)
243 {
244 dist = std::numeric_limits<double>::max();
245 return false;
246 }
247 dist = GetLocCoords(gloCoord, locCoord);
248 if (inside == 1)
249 {
250 dist = 0.;
251 return true;
252 }
253 else
254 {
255 Array<OneD, NekDouble> eta(GetShapeDim(), 0.);
256 m_xmap->LocCoordToLocCollapsed(locCoord, eta);
257 if (ClampLocCoords(eta, tol))
258 {
259 if (GetMetricInfo()->GetGtype() == eRegular)
260 {
261 dist = std::numeric_limits<double>::max();
262 }
263 return false;
264 }
265 return 3 != m_coordim ||
268 dist <= tol;
269 }
270}
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:548
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:488
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:536
StdRegions::StdExpansionSharedPtr m_xmap
mapping containing isoparametric transformation.
Definition Geometry.h:189

References ClampLocCoords(), Nektar::LibUtilities::eQuadrilateral, Nektar::SpatialDomains::eRegular, Nektar::LibUtilities::eTriangle, GetLocCoords(), GetMetricInfo(), GetShapeDim(), m_coordim, m_shapeType, m_xmap, and PreliminaryCheck().

Referenced by ContainsPoint(), ContainsPoint(), and ContainsPoint().

◆ v_FillGeom()

void Nektar::SpatialDomains::Geometry::v_FillGeom ( )
protectedvirtual

◆ v_FindDistance()

NekDouble Nektar::SpatialDomains::Geometry::v_FindDistance ( const Array< OneD, const NekDouble > &  xs,
Array< OneD, NekDouble > &  xi 
)
protectedvirtual

Reimplemented in Nektar::SpatialDomains::Geometry2D, and Nektar::SpatialDomains::SegGeom.

Definition at line 272 of file Geometry.cpp.

275{
277 "This function has not been defined for this geometry");
278 return false;
279}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by FindDistance().

◆ v_GenGeomFactors()

virtual void Nektar::SpatialDomains::Geometry::v_GenGeomFactors ( )
protectedpure virtual

◆ v_GetCoord()

NekDouble Nektar::SpatialDomains::Geometry::v_GetCoord ( const int  i,
const Array< OneD, const NekDouble > &  Lcoord 
)
protectedvirtual

Given local collapsed coordinate Lcoord, return the value of physical coordinate in direction i.

Reimplemented in Nektar::SpatialDomains::Geometry3D, Nektar::SpatialDomains::QuadGeom, Nektar::SpatialDomains::SegGeom, and Nektar::SpatialDomains::TriGeom.

Definition at line 339 of file Geometry.cpp.

342{
344 "This function is only valid for expansion type geometries");
345 return 0.0;
346}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetCoord().

◆ v_GetDir()

int Nektar::SpatialDomains::Geometry::v_GetDir ( const int  i,
const int  j 
) const
protectedvirtual

Returns the element coordinate direction corresponding to a given face coordinate direction.

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, Nektar::SpatialDomains::QuadGeom, Nektar::SpatialDomains::TetGeom, and Nektar::SpatialDomains::TriGeom.

Definition at line 328 of file Geometry.cpp.

330{
332 "This function has not been defined for this geometry");
333 return 0;
334}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetDir().

◆ v_GetEdge()

Geometry1D * Nektar::SpatialDomains::Geometry::v_GetEdge ( int  i) const
protectedvirtual

Returns edge i of this object.

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, Nektar::SpatialDomains::QuadGeom, Nektar::SpatialDomains::TetGeom, and Nektar::SpatialDomains::TriGeom.

Definition at line 136 of file Geometry.cpp.

137{
139 "This function is only valid for shape type geometries");
140 return nullptr;
141}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetEdge().

◆ v_GetEdgeFaceMap()

int Nektar::SpatialDomains::Geometry::v_GetEdgeFaceMap ( int  i,
int  j 
) const
protectedvirtual

Returns the standard element edge IDs that are connected to a given face.

For example, on a prism, edge 0 is connnected to faces 0 and 1; GetEdgeFaceMap(0,j) would therefore return the values 0 and 1 respectively. We assume that j runs between 0 and 1 inclusive, since every face is connected to precisely two faces for all 3D elements.

This function is used in the construction of the low-energy preconditioner.

Parameters
iThe edge to query connectivity for.
jThe local face index between 0 and 1 connected to this element.
See also
MultiRegions::PreconditionerLowEnergy

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, and Nektar::SpatialDomains::TetGeom.

Definition at line 306 of file Geometry.cpp.

308{
310 "This function has not been defined for this geometry");
311 return 0;
312}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetEdgeFaceMap().

◆ v_GetEdgeNormalToFaceVert()

int Nektar::SpatialDomains::Geometry::v_GetEdgeNormalToFaceVert ( const int  i,
const int  j 
) const
protectedvirtual

Returns the standard lement edge IDs that are normal to a given face vertex.

For example, on a hexahedron, on face 0 at vertices 0,1,2,3 the edges normal to that face are 4,5,6,7, ; so GetEdgeNormalToFaceVert(0,j) would therefore return the values 4, 5, 6 and 7 respectively. We assume that j runs between 0 and 3 inclusive on a quadrilateral face and between 0 and 2 inclusive on a triangular face.

This is used to help set up a length scale normal to an face

Parameters
iThe face to query for the normal edge
jThe local vertex index between 0 and nverts on this face

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, and Nektar::SpatialDomains::TetGeom.

Definition at line 317 of file Geometry.cpp.

319{
321 "This function has not been defined for this geometry");
322 return 0;
323}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetEdgeNormalToFaceVert().

◆ v_GetEorient()

StdRegions::Orientation Nektar::SpatialDomains::Geometry::v_GetEorient ( const int  i) const
protectedvirtual

Returns the orientation of edge i with respect to the ordering of edges in the standard element.

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, Nektar::SpatialDomains::QuadGeom, Nektar::SpatialDomains::TetGeom, and Nektar::SpatialDomains::TriGeom.

Definition at line 166 of file Geometry.cpp.

168{
170 "This function is not valid for this geometry.");
172}

References Nektar::ErrorUtil::efatal, Nektar::StdRegions::eForwards, and NEKERROR.

Referenced by GetEorient().

◆ v_GetFace()

Geometry2D * Nektar::SpatialDomains::Geometry::v_GetFace ( int  i) const
protectedvirtual

Returns face i of this object.

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, and Nektar::SpatialDomains::TetGeom.

Definition at line 146 of file Geometry.cpp.

147{
149 "This function is only valid for shape type geometries");
150 return nullptr;
151}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetFace().

◆ v_GetForient()

StdRegions::Orientation Nektar::SpatialDomains::Geometry::v_GetForient ( const int  i) const
protectedvirtual

Returns the orientation of face i with respect to the ordering of faces in the standard element.

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, and Nektar::SpatialDomains::TetGeom.

Definition at line 177 of file Geometry.cpp.

179{
181 "This function is not valid for this geometry.");
182 return StdRegions::eFwd;
183}

References Nektar::ErrorUtil::efatal, Nektar::StdRegions::eFwd, and NEKERROR.

Referenced by GetForient().

◆ v_GetLocCoords()

NekDouble Nektar::SpatialDomains::Geometry::v_GetLocCoords ( const Array< OneD, const NekDouble > &  coords,
Array< OneD, NekDouble > &  Lcoords 
)
protectedvirtual

Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this geometry object.

For curvilinear and non-affine elements (i.e. where the Jacobian varies as a function of the standard element coordinates), this is a non-linear optimisation problem that requires the use of a Newton iteration. Note therefore that this can be an expensive operation.

Note that, clearly, the provided Cartesian coordinate lie outside the element. The function therefore returns the minimum distance from some position in the element to . Lcoords will also be constrained to fit within the range \([-1,1]^d\) where \( d \) is the dimension of the element.

Parameters
coordsInput Cartesian global coordinates
LcoordsCorresponding local coordinates
Returns
Distance between obtained coordinates and provided ones.

Reimplemented in Nektar::SpatialDomains::Geometry1D, Nektar::SpatialDomains::Geometry2D, and Nektar::SpatialDomains::Geometry3D.

Definition at line 351 of file Geometry.cpp.

354{
356 "This function is only valid for expansion type geometries");
357 return 0.0;
358}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetLocCoords().

◆ v_GetNumEdges()

int Nektar::SpatialDomains::Geometry::v_GetNumEdges ( ) const
protectedvirtual

◆ v_GetNumFaces()

int Nektar::SpatialDomains::Geometry::v_GetNumFaces ( ) const
protectedvirtual

Get the number of faces of this object.

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, and Nektar::SpatialDomains::TetGeom.

Definition at line 196 of file Geometry.cpp.

197{
198 return 0;
199}

Referenced by GetNumFaces().

◆ v_GetNumVerts()

int Nektar::SpatialDomains::Geometry::v_GetNumVerts ( ) const
protectedvirtual

Get the number of vertices of this object.

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, Nektar::SpatialDomains::QuadGeom, Nektar::SpatialDomains::TetGeom, Nektar::SpatialDomains::TriGeom, and Nektar::SpatialDomains::SegGeom.

Definition at line 156 of file Geometry.cpp.

157{
159 "This function is only valid for shape type geometries");
160 return 0;
161}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetNumVerts().

◆ v_GetShapeDim()

int Nektar::SpatialDomains::Geometry::v_GetShapeDim ( ) const
protectedvirtual

Get the object's shape dimension.

For example, a segment is one dimensional and quadrilateral is two dimensional.

Reimplemented in Nektar::SpatialDomains::Geometry0D, Nektar::SpatialDomains::Geometry1D, Nektar::SpatialDomains::Geometry2D, and Nektar::SpatialDomains::Geometry3D.

Definition at line 204 of file Geometry.cpp.

205{
207 "This function is only valid for shape type geometries");
208 return 0;
209}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetShapeDim().

◆ v_GetVertex()

PointGeom * Nektar::SpatialDomains::Geometry::v_GetVertex ( int  i) const
protectedvirtual

Returns vertex i of this object.

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, Nektar::SpatialDomains::PyrGeom, Nektar::SpatialDomains::QuadGeom, Nektar::SpatialDomains::TetGeom, Nektar::SpatialDomains::TriGeom, and Nektar::SpatialDomains::SegGeom.

Definition at line 126 of file Geometry.cpp.

127{
129 "This function is only valid for shape type geometries");
130 return nullptr;
131}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetVertex().

◆ v_GetVertexEdgeMap()

int Nektar::SpatialDomains::Geometry::v_GetVertexEdgeMap ( int  i,
int  j 
) const
protectedvirtual

Returns the standard element edge IDs that are connected to a given vertex.

For example, on a prism, vertex 0 is connnected to edges 0, 3, and 4; GetVertexEdgeMap(0,j) would therefore return the values 0, 1 and 4 respectively. We assume that j runs between 0 and 2 inclusive, which is true for every 3D element asides from the pyramid.

This function is used in the construction of the low-energy preconditioner.

Parameters
iThe vertex to query connectivity for.
jThe local edge index between 0 and 2 connected to this element.
Todo:
Expand to work with pyramid elements.
See also
MultiRegions::PreconditionerLowEnergy

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, and Nektar::SpatialDomains::TetGeom.

Definition at line 284 of file Geometry.cpp.

286{
288 "This function has not been defined for this geometry");
289 return 0;
290}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetVertexEdgeMap().

◆ v_GetVertexFaceMap()

int Nektar::SpatialDomains::Geometry::v_GetVertexFaceMap ( int  i,
int  j 
) const
protectedvirtual

Returns the standard element face IDs that are connected to a given vertex.

For example, on a hexahedron, vertex 0 is connnected to faces 0, 1, and 4; GetVertexFaceMap(0,j) would therefore return the values 0, 1 and 4 respectively. We assume that j runs between 0 and 2 inclusive, which is true for every 3D element asides from the pyramid.

This is used in the construction of the low-energy preconditioner.

Parameters
iThe vertex to query connectivity for.
jThe local face index between 0 and 2 connected to this element.
Todo:
Expand to work with pyramid elements.
See also
MultiRegions::PreconditionerLowEnergy

Reimplemented in Nektar::SpatialDomains::HexGeom, Nektar::SpatialDomains::PrismGeom, and Nektar::SpatialDomains::TetGeom.

Definition at line 295 of file Geometry.cpp.

297{
299 "This function has not been defined for this geometry");
300 return 0;
301}

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetVertexFaceMap().

◆ v_GetVid()

int Nektar::SpatialDomains::Geometry::v_GetVid ( int  i) const
protectedvirtual

Get the ID of vertex i of this object.

Reimplemented in Nektar::SpatialDomains::PointGeom.

Definition at line 102 of file Geometry.cpp.

103{
104 return GetVertex(i)->GetGlobalID();
105}

References GetGlobalID(), and GetVertex().

Referenced by GetVid().

◆ v_GetXmap()

StdRegions::StdExpansionSharedPtr Nektar::SpatialDomains::Geometry::v_GetXmap ( ) const
protectedvirtual

Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standard element to physical element.

Definition at line 226 of file Geometry.cpp.

227{
228 return m_xmap;
229}

References m_xmap.

Referenced by GetXmap().

◆ v_Reset()

void Nektar::SpatialDomains::Geometry::v_Reset ( CurveMap curvedEdges,
CurveMap curvedFaces 
)
protectedvirtual

◆ v_Setup()

void Nektar::SpatialDomains::Geometry::v_Setup ( )
protectedvirtual

◆ ValidateRegGeomFactor()

GeomFactorsSharedPtr Nektar::SpatialDomains::Geometry::ValidateRegGeomFactor ( GeomFactorsSharedPtr  geomFactor)
staticprotected

Check to see if a geometric factor has already been created that contains the same regular information.

The principle behind this is that many regular (i.e. constant Jacobian) elements have identicial geometric factors. Memory may therefore be reduced by storing only the unique factors.

Parameters
geomFactorThe GeomFactor to check.
Returns
Either the cached GeomFactor or geomFactor.
Todo:
Currently this method is disabled since the lookup is very expensive.

Definition at line 81 of file Geometry.cpp.

83{
84 GeomFactorsSharedPtr returnval = geomFactor;
85
86 return returnval;
87}

Referenced by GetGeomFactors().

Member Data Documentation

◆ m_boundingBox

Array<OneD, NekDouble> Nektar::SpatialDomains::Geometry::m_boundingBox
protected

Array containing bounding box.

Definition at line 203 of file Geometry.h.

Referenced by ClearBoundingBox(), and GetBoundingBox().

◆ m_coeffs

Array<OneD, Array<OneD, NekDouble> > Nektar::SpatialDomains::Geometry::m_coeffs
protected

◆ m_coordim

int Nektar::SpatialDomains::Geometry::m_coordim
protected

Coordinate dimension of this geometry object.

Definition at line 183 of file Geometry.h.

Referenced by Nektar::SpatialDomains::PointGeom::Add(), Nektar::SpatialDomains::Geometry2D::Geometry2D(), Nektar::SpatialDomains::Geometry3D::Geometry3D(), GetCoordim(), Nektar::SpatialDomains::PointGeom::GetCoords(), Nektar::SpatialDomains::PointGeom::GetCoords(), MinMaxCheck(), Nektar::SpatialDomains::PointGeom::Mult(), Nektar::SpatialDomains::PointGeom::PointGeom(), Nektar::SpatialDomains::PointGeom::PointGeom(), Nektar::SpatialDomains::PointGeom::PointGeom(), Nektar::SpatialDomains::SegGeom::SegGeom(), SetCoordim(), SetUpCoeffs(), Nektar::SpatialDomains::HexGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::PrismGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::PyrGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::TetGeom::SetUpFaceOrientation(), Nektar::SpatialDomains::PointGeom::Sub(), Nektar::SpatialDomains::TriGeom::TriGeom(), Nektar::SpatialDomains::QuadGeom::v_AllLeftCheck(), Nektar::SpatialDomains::TriGeom::v_AllLeftCheck(), v_ContainsPoint(), Nektar::SpatialDomains::HexGeom::v_FillGeom(), Nektar::SpatialDomains::PrismGeom::v_FillGeom(), Nektar::SpatialDomains::PyrGeom::v_FillGeom(), Nektar::SpatialDomains::QuadGeom::v_FillGeom(), Nektar::SpatialDomains::SegGeom::v_FillGeom(), Nektar::SpatialDomains::TetGeom::v_FillGeom(), Nektar::SpatialDomains::TriGeom::v_FillGeom(), Nektar::SpatialDomains::SegGeom::v_FindDistance(), Nektar::SpatialDomains::HexGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PointGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PrismGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PyrGeom::v_GenGeomFactors(), Nektar::SpatialDomains::QuadGeom::v_GenGeomFactors(), Nektar::SpatialDomains::SegGeom::v_GenGeomFactors(), Nektar::SpatialDomains::TetGeom::v_GenGeomFactors(), Nektar::SpatialDomains::TriGeom::v_GenGeomFactors(), Nektar::SpatialDomains::Geometry1D::v_GetLocCoords(), and Nektar::SpatialDomains::Geometry2D::v_GetLocCoords().

◆ m_geomFactors

GeomFactorsSharedPtr Nektar::SpatialDomains::Geometry::m_geomFactors
protected

◆ m_geomFactorsState

GeomState Nektar::SpatialDomains::Geometry::m_geomFactorsState
protected

◆ m_geomType

GeomType Nektar::SpatialDomains::Geometry::m_geomType
protected

Type of geometry.

Definition at line 195 of file Geometry.h.

◆ m_globalID

int Nektar::SpatialDomains::Geometry::m_globalID
protected

◆ m_invIsoParam

Array<OneD, Array<OneD, NekDouble> > Nektar::SpatialDomains::Geometry::m_invIsoParam
protected

◆ m_isoParameter

Array<OneD, Array<OneD, NekDouble> > Nektar::SpatialDomains::Geometry::m_isoParameter
protected

◆ m_regGeomFactorsManager

GeomFactorsVector Nektar::SpatialDomains::Geometry::m_regGeomFactorsManager
staticprotected

Definition at line 180 of file Geometry.h.

◆ m_setupState

bool Nektar::SpatialDomains::Geometry::m_setupState
protected

◆ m_shapeType

LibUtilities::ShapeType Nektar::SpatialDomains::Geometry::m_shapeType
protected

◆ m_state

GeomState Nektar::SpatialDomains::Geometry::m_state
protected

◆ m_straightEdge

int Nektar::SpatialDomains::Geometry::m_straightEdge
protected

◆ m_xmap

StdRegions::StdExpansionSharedPtr Nektar::SpatialDomains::Geometry::m_xmap
protected

\(\chi\) mapping containing isoparametric transformation.

Definition at line 189 of file Geometry.h.

Referenced by Nektar::SpatialDomains::SegGeom::GenerateOneSpaceDimGeom(), Nektar::SpatialDomains::Geometry3D::NewtonIterationForLocCoord(), Nektar::SpatialDomains::Geometry2D::NewtonIterationForLocCoord(), Nektar::SpatialDomains::Geometry3D::NewtonIterationForLocCoord(), Nektar::SpatialDomains::SegGeom::SegGeom(), Nektar::SpatialDomains::HexGeom::SetUpXmap(), Nektar::SpatialDomains::PrismGeom::SetUpXmap(), Nektar::SpatialDomains::PyrGeom::SetUpXmap(), Nektar::SpatialDomains::QuadGeom::SetUpXmap(), Nektar::SpatialDomains::SegGeom::SetUpXmap(), Nektar::SpatialDomains::TetGeom::SetUpXmap(), Nektar::SpatialDomains::TriGeom::SetUpXmap(), v_ContainsPoint(), Nektar::SpatialDomains::HexGeom::v_FillGeom(), Nektar::SpatialDomains::PrismGeom::v_FillGeom(), Nektar::SpatialDomains::PyrGeom::v_FillGeom(), Nektar::SpatialDomains::QuadGeom::v_FillGeom(), Nektar::SpatialDomains::SegGeom::v_FillGeom(), Nektar::SpatialDomains::TetGeom::v_FillGeom(), Nektar::SpatialDomains::TriGeom::v_FillGeom(), Nektar::SpatialDomains::Geometry2D::v_FindDistance(), Nektar::SpatialDomains::SegGeom::v_FindDistance(), Nektar::SpatialDomains::HexGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PointGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PrismGeom::v_GenGeomFactors(), Nektar::SpatialDomains::PyrGeom::v_GenGeomFactors(), Nektar::SpatialDomains::QuadGeom::v_GenGeomFactors(), Nektar::SpatialDomains::SegGeom::v_GenGeomFactors(), Nektar::SpatialDomains::TetGeom::v_GenGeomFactors(), Nektar::SpatialDomains::TriGeom::v_GenGeomFactors(), Nektar::SpatialDomains::Geometry3D::v_GetCoord(), Nektar::SpatialDomains::QuadGeom::v_GetCoord(), Nektar::SpatialDomains::SegGeom::v_GetCoord(), Nektar::SpatialDomains::TriGeom::v_GetCoord(), Nektar::SpatialDomains::Geometry1D::v_GetLocCoords(), Nektar::SpatialDomains::Geometry2D::v_GetLocCoords(), Nektar::SpatialDomains::Geometry3D::v_GetLocCoords(), v_GetXmap(), Nektar::SpatialDomains::HexGeom::v_Reset(), Nektar::SpatialDomains::PyrGeom::v_Reset(), Nektar::SpatialDomains::QuadGeom::v_Reset(), Nektar::SpatialDomains::SegGeom::v_Reset(), Nektar::SpatialDomains::TriGeom::v_Reset(), Nektar::SpatialDomains::HexGeom::v_Setup(), Nektar::SpatialDomains::PrismGeom::v_Setup(), Nektar::SpatialDomains::PyrGeom::v_Setup(), Nektar::SpatialDomains::QuadGeom::v_Setup(), Nektar::SpatialDomains::SegGeom::v_Setup(), Nektar::SpatialDomains::TetGeom::v_Setup(), and Nektar::SpatialDomains::TriGeom::v_Setup().