Nektar++
Public Member Functions | Static Public Attributes | Protected Member Functions | Private Member Functions | Static Private Attributes | List of all members
Nektar::SpatialDomains::TetGeom Class Reference

#include <TetGeom.h>

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

Public Member Functions

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

Static Public Attributes

static const int kNverts = 4
 
static const int kNedges = 6
 
static const int kNqfaces = 0
 
static const int kNtfaces = 4
 
static const int kNfaces = kNqfaces + kNtfaces
 
static const std::string XMLElementType
 
- Static Public Attributes inherited from Nektar::SpatialDomains::Geometry3D
static const int kDim = 3
 

Protected Member Functions

int v_GetVertexEdgeMap (const int i, const int j) const override
 Returns the standard element edge IDs that are connected to a given vertex. More...
 
int v_GetVertexFaceMap (const int i, const int j) const override
 Returns the standard element face IDs that are connected to a given vertex. More...
 
int v_GetEdgeFaceMap (const int i, const int j) const override
 Returns the standard element edge IDs that are connected to a given face. More...
 
int v_GetEdgeNormalToFaceVert (const int i, const int j) const override
 Returns the standard lement edge IDs that are normal to a given face vertex. More...
 
int v_GetDir (const int faceidx, const int facedir) const override
 Returns the element coordinate direction corresponding to a given face coordinate direction. More...
 
void v_Reset (CurveMap &curvedEdges, CurveMap &curvedFaces) override
 Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors. More...
 
void v_Setup () override
 
void v_GenGeomFactors () override
 
- Protected Member Functions inherited from Nektar::SpatialDomains::Geometry3D
NekDouble v_GetLocCoords (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords) override
 Determine the local collapsed coordinates that correspond to a given Cartesian coordinate for this geometry object. More...
 
void NewtonIterationForLocCoord (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &ptsx, const Array< OneD, const NekDouble > &ptsy, const Array< OneD, const NekDouble > &ptsz, Array< OneD, NekDouble > &Lcoords, NekDouble &dist)
 
void NewtonIterationForLocCoord (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
 
void v_FillGeom () override
 Put all quadrature information into face/edge structure and backward transform. More...
 
NekDouble v_GetCoord (const int i, const Array< OneD, const NekDouble > &Lcoord) override
 Given local collapsed coordinate Lcoord return the value of physical coordinate in direction i. More...
 
void v_CalculateInverseIsoParam () override
 
int v_AllLeftCheck (const Array< OneD, const NekDouble > &gloCoord) override
 
int v_GetShapeDim () const override
 Get the object's shape dimension. More...
 
int v_GetNumVerts () const override
 Get the number of vertices of this object. More...
 
int v_GetNumEdges () const override
 Get the number of edges of this object. More...
 
int v_GetNumFaces () const override
 Get the number of faces of this object. More...
 
PointGeomSharedPtr v_GetVertex (int i) const override
 
Geometry1DSharedPtr v_GetEdge (int i) const override
 Returns edge i of this object. More...
 
Geometry2DSharedPtr v_GetFace (int i) const override
 Returns face i of this object. More...
 
StdRegions::Orientation v_GetEorient (const int i) const override
 Returns the orientation of edge i with respect to the ordering of edges in the standard element. More...
 
StdRegions::Orientation v_GetForient (const int i) const override
 Returns the orientation of face i with respect to the ordering of faces in the standard element. More...
 
- Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
virtual PointGeomSharedPtr v_GetVertex (int i) const =0
 
virtual Geometry1DSharedPtr v_GetEdge (int i) const
 Returns edge i of this object. More...
 
virtual Geometry2DSharedPtr v_GetFace (int i) const
 Returns face i of this object. More...
 
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. More...
 
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. More...
 
virtual int v_GetNumVerts () const
 Get the number of vertices of this object. More...
 
virtual int v_GetNumEdges () const
 Get the number of edges of this object. More...
 
virtual int v_GetNumFaces () const
 Get the number of faces of this object. More...
 
virtual int v_GetShapeDim () const
 Get the object's shape dimension. More...
 
virtual StdRegions::StdExpansionSharedPtr v_GetXmap () const
 Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standard element to physical element. More...
 
virtual void v_FillGeom ()
 Populate the coordinate mapping Geometry::m_coeffs information from any children geometry elements. More...
 
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)\). More...
 
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. More...
 
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. More...
 
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. More...
 
virtual int v_GetVertexFaceMap (int i, int j) const
 Returns the standard element face IDs that are connected to a given vertex. More...
 
virtual int v_GetEdgeFaceMap (int i, int j) const
 Returns the standard element edge IDs that are connected to a given face. More...
 
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. More...
 
virtual int v_GetDir (const int faceidx, const int facedir) const
 Returns the element coordinate direction corresponding to a given face coordinate direction. More...
 
virtual void v_Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors. More...
 
virtual void v_Setup ()
 
virtual void v_GenGeomFactors ()=0
 
void SetUpCoeffs (const int nCoeffs)
 Initialise the Geometry::m_coeffs array. More...
 
virtual void v_CalculateInverseIsoParam ()
 

Private Member Functions

void SetUpLocalEdges ()
 
void SetUpLocalVertices ()
 
void SetUpEdgeOrientation ()
 
void SetUpFaceOrientation ()
 
void SetUpXmap ()
 Set up the m_xmap object by determining the order of each direction from derived faces. More...
 

Static Private Attributes

static const unsigned int VertexEdgeConnectivity [4][3]
 
static const unsigned int VertexFaceConnectivity [4][3]
 
static const unsigned int EdgeFaceConnectivity [6][2]
 
static const unsigned int EdgeNormalToFaceVert [4][3]
 

Additional Inherited Members

- Static Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsSharedPtr ValidateRegGeomFactor (GeomFactorsSharedPtr geomFactor)
 Check to see if a geometric factor has already been created that contains the same regular information. More...
 
- Protected Attributes inherited from Nektar::SpatialDomains::Geometry3D
PointGeomVector m_verts
 
SegGeomVector m_edges
 
Geometry2DVector m_faces
 
std::vector< StdRegions::Orientationm_eorient
 
std::vector< StdRegions::Orientationm_forient
 
int m_eid
 
bool m_ownverts
 
- Protected Attributes inherited from Nektar::SpatialDomains::Geometry
int m_coordim
 Coordinate dimension of this geometry object. More...
 
GeomFactorsSharedPtr m_geomFactors
 Geometric factors. More...
 
GeomState m_geomFactorsState
 State of the geometric factors. More...
 
StdRegions::StdExpansionSharedPtr m_xmap
 \(\chi\) mapping containing isoparametric transformation. More...
 
GeomState m_state
 Enumeration to dictate whether coefficients are filled. More...
 
bool m_setupState
 Wether or not the setup routines have been run. More...
 
GeomType m_geomType
 Type of geometry. More...
 
LibUtilities::ShapeType m_shapeType
 Type of shape. More...
 
int m_globalID
 Global ID. More...
 
Array< OneD, Array< OneD, NekDouble > > m_coeffs
 Array containing expansion coefficients of m_xmap. More...
 
Array< OneD, NekDoublem_boundingBox
 Array containing bounding box. More...
 
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
 
Array< OneD, Array< OneD, NekDouble > > m_invIsoParam
 
int m_straightEdge
 
- Static Protected Attributes inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsVector m_regGeomFactorsManager
 

Detailed Description

Definition at line 45 of file TetGeom.h.

Constructor & Destructor Documentation

◆ TetGeom() [1/2]

Nektar::SpatialDomains::TetGeom::TetGeom ( )

◆ TetGeom() [2/2]

Nektar::SpatialDomains::TetGeom::TetGeom ( int  id,
const TriGeomSharedPtr  faces[] 
)

Copy the face shared pointers

Set up orientation vectors with correct amount of elements.

Definition at line 57 of file TetGeom.cpp.

58 : Geometry3D(faces[0]->GetEdge(0)->GetVertex(0)->GetCoordim())
59{
61 m_globalID = id;
62
63 /// Copy the face shared pointers
64 m_faces.insert(m_faces.begin(), faces, faces + TetGeom::kNfaces);
65
66 /// Set up orientation vectors with correct amount of elements.
67 m_eorient.resize(kNedges);
68 m_forient.resize(kNfaces);
69
74}
std::vector< StdRegions::Orientation > m_forient
Definition: Geometry3D.h:81
std::vector< StdRegions::Orientation > m_eorient
Definition: Geometry3D.h:80
PointGeomSharedPtr GetVertex(int i) const
Returns vertex i of this object.
Definition: Geometry.h:357
Geometry1DSharedPtr GetEdge(int i) const
Returns edge i of this object.
Definition: Geometry.h:365
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
static const int kNfaces
Definition: TetGeom.h:56
static const int kNedges
Definition: TetGeom.h:53

References Nektar::LibUtilities::eTetrahedron, kNedges, kNfaces, Nektar::SpatialDomains::Geometry3D::m_eorient, Nektar::SpatialDomains::Geometry3D::m_faces, Nektar::SpatialDomains::Geometry3D::m_forient, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_shapeType, SetUpEdgeOrientation(), SetUpFaceOrientation(), SetUpLocalEdges(), and SetUpLocalVertices().

◆ ~TetGeom()

Nektar::SpatialDomains::TetGeom::~TetGeom ( )
overridedefault

Member Function Documentation

◆ SetUpEdgeOrientation()

void Nektar::SpatialDomains::TetGeom::SetUpEdgeOrientation ( )
private

Definition at line 315 of file TetGeom.cpp.

316{
317
318 // This 2D array holds the local id's of all the vertices
319 // for every edge. For every edge, they are ordered to what we
320 // define as being Forwards
321 const unsigned int edgeVerts[kNedges][2] = {{0, 1}, {1, 2}, {0, 2},
322 {0, 3}, {1, 3}, {2, 3}};
323
324 int i;
325 for (i = 0; i < kNedges; i++)
326 {
327 if (m_edges[i]->GetVid(0) == m_verts[edgeVerts[i][0]]->GetGlobalID())
328 {
330 }
331 else if (m_edges[i]->GetVid(0) ==
332 m_verts[edgeVerts[i][1]]->GetGlobalID())
333 {
335 }
336 else
337 {
338 ASSERTL0(false, "Could not find matching vertex for the edge");
339 }
340 }
341}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
int GetVid(int i) const
Get the ID of vertex i of this object.
Definition: Geometry.cpp:104
int GetGlobalID(void) const
Get the ID of this object.
Definition: Geometry.h:326

References ASSERTL0, Nektar::StdRegions::eBackwards, Nektar::StdRegions::eForwards, Nektar::SpatialDomains::Geometry::GetGlobalID(), Nektar::SpatialDomains::Geometry::GetVid(), kNedges, Nektar::SpatialDomains::Geometry3D::m_edges, Nektar::SpatialDomains::Geometry3D::m_eorient, and Nektar::SpatialDomains::Geometry3D::m_verts.

Referenced by TetGeom().

◆ SetUpFaceOrientation()

void Nektar::SpatialDomains::TetGeom::SetUpFaceOrientation ( )
private

Definition at line 343 of file TetGeom.cpp.

344{
345
346 int f, i;
347
348 // These arrays represent the vector of the A and B
349 // coordinate of the local elemental coordinate system
350 // where A corresponds with the coordinate direction xi_i
351 // with the lowest index i (for that particular face)
352 // Coordinate 'B' then corresponds to the other local
353 // coordinate (i.e. with the highest index)
354 Array<OneD, NekDouble> elementAaxis(m_coordim);
355 Array<OneD, NekDouble> elementBaxis(m_coordim);
356
357 // These arrays correspond to the local coordinate
358 // system of the face itself (i.e. the Geometry2D)
359 // faceAaxis correspond to the xi_0 axis
360 // faceBaxis correspond to the xi_1 axis
361 Array<OneD, NekDouble> faceAaxis(m_coordim);
362 Array<OneD, NekDouble> faceBaxis(m_coordim);
363
364 // This is the base vertex of the face (i.e. the Geometry2D)
365 // This corresponds to thevertex with local ID 0 of the
366 // Geometry2D
367 unsigned int baseVertex;
368
369 // The lenght of the vectors above
370 NekDouble elementAaxis_length;
371 NekDouble elementBaxis_length;
372 NekDouble faceAaxis_length;
373 NekDouble faceBaxis_length;
374
375 // This 2D array holds the local id's of all the vertices
376 // for every face. For every face, they are ordered in such
377 // a way that the implementation below allows a unified approach
378 // for all faces.
379 const unsigned int faceVerts[kNfaces][TriGeom::kNverts] = {
380 {0, 1, 2}, {0, 1, 3}, {1, 2, 3}, {0, 2, 3}};
381
382 NekDouble dotproduct1 = 0.0;
383 NekDouble dotproduct2 = 0.0;
384
385 unsigned int orientation;
386
387 // Loop over all the faces to set up the orientation
388 for (f = 0; f < kNqfaces + kNtfaces; f++)
389 {
390 // initialisation
391 elementAaxis_length = 0.0;
392 elementBaxis_length = 0.0;
393 faceAaxis_length = 0.0;
394 faceBaxis_length = 0.0;
395
396 dotproduct1 = 0.0;
397 dotproduct2 = 0.0;
398
399 baseVertex = m_faces[f]->GetVid(0);
400
401 // We are going to construct the vectors representing the A and B axis
402 // of every face. These vectors will be constructed as a
403 // vector-representation
404 // of the edges of the face. However, for both coordinate directions, we
405 // can
406 // represent the vectors by two different edges. That's why we need to
407 // make sure that
408 // we pick the edge to which the baseVertex of the
409 // Geometry2D-representation of the face
410 // belongs...
411
412 // Compute the length of edges on a base-face
413 if (baseVertex == m_verts[faceVerts[f][0]]->GetGlobalID())
414 {
415 for (i = 0; i < m_coordim; i++)
416 {
417 elementAaxis[i] = (*m_verts[faceVerts[f][1]])[i] -
418 (*m_verts[faceVerts[f][0]])[i];
419 elementBaxis[i] = (*m_verts[faceVerts[f][2]])[i] -
420 (*m_verts[faceVerts[f][0]])[i];
421 }
422 }
423 else if (baseVertex == m_verts[faceVerts[f][1]]->GetGlobalID())
424 {
425 for (i = 0; i < m_coordim; i++)
426 {
427 elementAaxis[i] = (*m_verts[faceVerts[f][1]])[i] -
428 (*m_verts[faceVerts[f][0]])[i];
429 elementBaxis[i] = (*m_verts[faceVerts[f][2]])[i] -
430 (*m_verts[faceVerts[f][1]])[i];
431 }
432 }
433 else if (baseVertex == m_verts[faceVerts[f][2]]->GetGlobalID())
434 {
435 for (i = 0; i < m_coordim; i++)
436 {
437 elementAaxis[i] = (*m_verts[faceVerts[f][1]])[i] -
438 (*m_verts[faceVerts[f][2]])[i];
439 elementBaxis[i] = (*m_verts[faceVerts[f][2]])[i] -
440 (*m_verts[faceVerts[f][0]])[i];
441 }
442 }
443 else
444 {
445 ASSERTL0(false, "Could not find matching vertex for the face");
446 }
447
448 // Now, construct the edge-vectors of the local coordinates of
449 // the Geometry2D-representation of the face
450 for (i = 0; i < m_coordim; i++)
451 {
452 faceAaxis[i] =
453 (*m_faces[f]->GetVertex(1))[i] - (*m_faces[f]->GetVertex(0))[i];
454 faceBaxis[i] =
455 (*m_faces[f]->GetVertex(2))[i] - (*m_faces[f]->GetVertex(0))[i];
456
457 elementAaxis_length += pow(elementAaxis[i], 2);
458 elementBaxis_length += pow(elementBaxis[i], 2);
459 faceAaxis_length += pow(faceAaxis[i], 2);
460 faceBaxis_length += pow(faceBaxis[i], 2);
461 }
462
463 elementAaxis_length = sqrt(elementAaxis_length);
464 elementBaxis_length = sqrt(elementBaxis_length);
465 faceAaxis_length = sqrt(faceAaxis_length);
466 faceBaxis_length = sqrt(faceBaxis_length);
467
468 // Calculate the inner product of both the A-axis
469 // (i.e. Elemental A axis and face A axis)
470 for (i = 0; i < m_coordim; i++)
471 {
472 dotproduct1 += elementAaxis[i] * faceAaxis[i];
473 }
474
475 NekDouble norm =
476 fabs(dotproduct1) / elementAaxis_length / faceAaxis_length;
477 orientation = 0;
478
479 // if the innerproduct is equal to the (absolute value of the ) products
480 // of the lengths of both vectors, then, the coordinate systems will NOT
481 // be transposed
482 if (fabs(norm - 1.0) < NekConstants::kNekZeroTol)
483 {
484 // if the inner product is negative, both A-axis point
485 // in reverse direction
486 if (dotproduct1 < 0.0)
487 {
488 orientation += 2;
489 }
490
491 // calculate the inner product of both B-axis
492 for (i = 0; i < m_coordim; i++)
493 {
494 dotproduct2 += elementBaxis[i] * faceBaxis[i];
495 }
496
497 norm = fabs(dotproduct2) / elementBaxis_length / faceBaxis_length;
498
499 // check that both these axis are indeed parallel
500 ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
501 "These vectors should be parallel");
502
503 // if the inner product is negative, both B-axis point
504 // in reverse direction
505 if (dotproduct2 < 0.0)
506 {
507 orientation++;
508 }
509 }
510 // The coordinate systems are transposed
511 else
512 {
513 orientation = 4;
514
515 // Calculate the inner product between the elemental A-axis
516 // and the B-axis of the face (which are now the corresponding axis)
517 dotproduct1 = 0.0;
518 for (i = 0; i < m_coordim; i++)
519 {
520 dotproduct1 += elementAaxis[i] * faceBaxis[i];
521 }
522
523 norm = fabs(dotproduct1) / elementAaxis_length / faceBaxis_length;
524
525 // check that both these axis are indeed parallel
526 ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
527 "These vectors should be parallel");
528
529 // if the result is negative, both axis point in reverse
530 // directions
531 if (dotproduct1 < 0.0)
532 {
533 orientation += 2;
534 }
535
536 // Do the same for the other two corresponding axis
537 dotproduct2 = 0.0;
538 for (i = 0; i < m_coordim; i++)
539 {
540 dotproduct2 += elementBaxis[i] * faceAaxis[i];
541 }
542
543 norm = fabs(dotproduct2) / elementBaxis_length / faceAaxis_length;
544
545 // check that both these axis are indeed parallel
546 ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
547 "These vectors should be parallel");
548
549 if (dotproduct2 < 0.0)
550 {
551 orientation++;
552 }
553 }
554
555 orientation = orientation + 5;
556
558 "Orientation of triangular face (id = " +
559 std::to_string(m_faces[f]->GetGlobalID()) +
560 ") is inconsistent with face " + std::to_string(f) +
561 " of tet element (id = " + std::to_string(m_globalID) +
562 ") since Dir2 is aligned with Dir1. Mesh setup "
563 "needs investigation");
564
565 // Fill the m_forient array
566 m_forient[f] = (StdRegions::Orientation)orientation;
567 }
568}
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
int m_coordim
Coordinate dimension of this geometry object.
Definition: Geometry.h:188
static const int kNqfaces
Definition: TetGeom.h:54
static const int kNtfaces
Definition: TetGeom.h:55
static const int kNverts
Definition: TriGeom.h:71
static const NekDouble kNekZeroTol
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References ASSERTL0, ASSERTL1, Nektar::StdRegions::eDir1FwdDir2_Dir2FwdDir1, Nektar::SpatialDomains::Geometry::GetGlobalID(), Nektar::SpatialDomains::Geometry::GetVertex(), Nektar::NekConstants::kNekZeroTol, kNfaces, kNqfaces, kNtfaces, Nektar::SpatialDomains::TriGeom::kNverts, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry3D::m_faces, Nektar::SpatialDomains::Geometry3D::m_forient, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry3D::m_verts, and tinysimd::sqrt().

Referenced by TetGeom().

◆ SetUpLocalEdges()

void Nektar::SpatialDomains::TetGeom::SetUpLocalEdges ( )
private

Definition at line 112 of file TetGeom.cpp.

113{
114
115 // find edge 0
116 int i, j;
117 unsigned int check;
118
119 SegGeomSharedPtr edge;
120
121 // First set up the 3 bottom edges
122
123 if (m_faces[0]->GetEid(0) != m_faces[1]->GetEid(0))
124 {
125 std::ostringstream errstrm;
126 errstrm << "Local edge 0 (eid=" << m_faces[0]->GetEid(0);
127 errstrm << ") on face " << m_faces[0]->GetGlobalID();
128 errstrm << " must be the same as local edge 0 (eid="
129 << m_faces[1]->GetEid(0);
130 errstrm << ") on face " << m_faces[1]->GetGlobalID();
131 ASSERTL0(false, errstrm.str());
132 }
133
134 int faceConnected;
135 for (faceConnected = 1; faceConnected < 4; faceConnected++)
136 {
137 check = 0;
138 for (i = 0; i < 3; i++)
139 {
140 if ((m_faces[0])->GetEid(i) == (m_faces[faceConnected])->GetEid(0))
141 {
142 edge = std::dynamic_pointer_cast<SegGeom>(
143 (m_faces[0])->GetEdge(i));
144 m_edges.push_back(edge);
145 check++;
146 }
147 }
148
149 if (check < 1)
150 {
151 std::ostringstream errstrm;
152 errstrm << "Face 0 does not share an edge with first edge of "
153 "adjacent face. Faces ";
154 errstrm << (m_faces[0])->GetGlobalID() << ", "
155 << (m_faces[faceConnected])->GetGlobalID();
156 ASSERTL0(false, errstrm.str());
157 }
158 else if (check > 1)
159 {
160 std::ostringstream errstrm;
161 errstrm << "Connected faces share more than one edge. Faces ";
162 errstrm << (m_faces[0])->GetGlobalID() << ", "
163 << (m_faces[faceConnected])->GetGlobalID();
164 ASSERTL0(false, errstrm.str());
165 }
166 }
167
168 // Then, set up the 3 vertical edges
169 check = 0;
170 for (i = 0; i < 3; i++) // Set up the vertical edge :face(1) and face(3)
171 {
172 for (j = 0; j < 3; j++)
173 {
174 if ((m_faces[1])->GetEid(i) == (m_faces[3])->GetEid(j))
175 {
176 edge = std::dynamic_pointer_cast<SegGeom>(
177 (m_faces[1])->GetEdge(i));
178 m_edges.push_back(edge);
179 check++;
180 }
181 }
182 }
183 if (check < 1)
184 {
185 std::ostringstream errstrm;
186 errstrm << "Connected faces do not share an edge. Faces ";
187 errstrm << (m_faces[1])->GetGlobalID() << ", "
188 << (m_faces[3])->GetGlobalID();
189 ASSERTL0(false, errstrm.str());
190 }
191 else if (check > 1)
192 {
193 std::ostringstream errstrm;
194 errstrm << "Connected faces share more than one edge. Faces ";
195 errstrm << (m_faces[1])->GetGlobalID() << ", "
196 << (m_faces[3])->GetGlobalID();
197 ASSERTL0(false, errstrm.str());
198 }
199 // Set up vertical edges: face(1) through face(3)
200 for (faceConnected = 1; faceConnected < 3; faceConnected++)
201 {
202 check = 0;
203 for (i = 0; i < 3; i++)
204 {
205 for (j = 0; j < 3; j++)
206 {
207 if ((m_faces[faceConnected])->GetEid(i) ==
208 (m_faces[faceConnected + 1])->GetEid(j))
209 {
210 edge = std::dynamic_pointer_cast<SegGeom>(
211 (m_faces[faceConnected])->GetEdge(i));
212 m_edges.push_back(edge);
213 check++;
214 }
215 }
216 }
217
218 if (check < 1)
219 {
220 std::ostringstream errstrm;
221 errstrm << "Connected faces do not share an edge. Faces ";
222 errstrm << (m_faces[faceConnected])->GetGlobalID() << ", "
223 << (m_faces[faceConnected + 1])->GetGlobalID();
224 ASSERTL0(false, errstrm.str());
225 }
226 else if (check > 1)
227 {
228 std::ostringstream errstrm;
229 errstrm << "Connected faces share more than one edge. Faces ";
230 errstrm << (m_faces[faceConnected])->GetGlobalID() << ", "
231 << (m_faces[faceConnected + 1])->GetGlobalID();
232 ASSERTL0(false, errstrm.str());
233 }
234 }
235}
int GetEid(int i) const
Get the ID of edge i of this object.
Definition: Geometry.cpp:112
std::shared_ptr< SegGeom > SegGeomSharedPtr
Definition: Geometry2D.h:58

References ASSERTL0, Nektar::SpatialDomains::Geometry::GetEdge(), Nektar::SpatialDomains::Geometry::GetEid(), Nektar::SpatialDomains::Geometry::GetGlobalID(), Nektar::SpatialDomains::Geometry3D::m_edges, and Nektar::SpatialDomains::Geometry3D::m_faces.

Referenced by TetGeom().

◆ SetUpLocalVertices()

void Nektar::SpatialDomains::TetGeom::SetUpLocalVertices ( )
private

Definition at line 237 of file TetGeom.cpp.

238{
239
240 // Set up the first 2 vertices (i.e. vertex 0,1)
241 if ((m_edges[0]->GetVid(0) == m_edges[1]->GetVid(0)) ||
242 (m_edges[0]->GetVid(0) == m_edges[1]->GetVid(1)))
243 {
244 m_verts.push_back(m_edges[0]->GetVertex(1));
245 m_verts.push_back(m_edges[0]->GetVertex(0));
246 }
247 else if ((m_edges[0]->GetVid(1) == m_edges[1]->GetVid(0)) ||
248 (m_edges[0]->GetVid(1) == m_edges[1]->GetVid(1)))
249 {
250 m_verts.push_back(m_edges[0]->GetVertex(0));
251 m_verts.push_back(m_edges[0]->GetVertex(1));
252 }
253 else
254 {
255 std::ostringstream errstrm;
256 errstrm << "Connected edges do not share a vertex. Edges ";
257 errstrm << m_edges[0]->GetGlobalID() << ", "
258 << m_edges[1]->GetGlobalID();
259 ASSERTL0(false, errstrm.str());
260 }
261
262 // set up the other bottom vertices (i.e. vertex 2)
263 for (int i = 1; i < 2; i++)
264 {
265 if (m_edges[i]->GetVid(0) == m_verts[i]->GetGlobalID())
266 {
267 m_verts.push_back(m_edges[i]->GetVertex(1));
268 }
269 else if (m_edges[i]->GetVid(1) == m_verts[i]->GetGlobalID())
270 {
271 m_verts.push_back(m_edges[i]->GetVertex(0));
272 }
273 else
274 {
275 std::ostringstream errstrm;
276 errstrm << "Connected edges do not share a vertex. Edges ";
277 errstrm << m_edges[i]->GetGlobalID() << ", "
278 << m_edges[i - 1]->GetGlobalID();
279 ASSERTL0(false, errstrm.str());
280 }
281 }
282
283 // set up top vertex
284 if (m_edges[3]->GetVid(0) == m_verts[0]->GetGlobalID())
285 {
286 m_verts.push_back(m_edges[3]->GetVertex(1));
287 }
288 else
289 {
290 m_verts.push_back(m_edges[3]->GetVertex(0));
291 }
292
293 // Check the other edges match up.
294 int check = 0;
295 for (int i = 4; i < 6; ++i)
296 {
297 if ((m_edges[i]->GetVid(0) == m_verts[i - 3]->GetGlobalID() &&
298 m_edges[i]->GetVid(1) == m_verts[3]->GetGlobalID()) ||
299 (m_edges[i]->GetVid(1) == m_verts[i - 3]->GetGlobalID() &&
300 m_edges[i]->GetVid(0) == m_verts[3]->GetGlobalID()))
301 {
302 check++;
303 }
304 }
305 if (check != 2)
306 {
307 std::ostringstream errstrm;
308 errstrm << "Connected edges do not share a vertex. Edges ";
309 errstrm << m_edges[3]->GetGlobalID() << ", "
310 << m_edges[2]->GetGlobalID();
311 ASSERTL0(false, errstrm.str());
312 }
313}

References ASSERTL0, Nektar::SpatialDomains::Geometry::GetGlobalID(), Nektar::SpatialDomains::Geometry::GetVertex(), Nektar::SpatialDomains::Geometry::GetVid(), Nektar::SpatialDomains::Geometry3D::m_edges, and Nektar::SpatialDomains::Geometry3D::m_verts.

Referenced by TetGeom().

◆ SetUpXmap()

void Nektar::SpatialDomains::TetGeom::SetUpXmap ( )
private

Set up the m_xmap object by determining the order of each direction from derived faces.

Definition at line 649 of file TetGeom.cpp.

650{
651 std::vector<int> tmp;
652 tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(0));
653 int order0 = *std::max_element(tmp.begin(), tmp.end());
654
655 tmp.clear();
656 tmp.push_back(order0);
657 tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(1));
658 tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(2));
659 int order1 = *std::max_element(tmp.begin(), tmp.end());
660
661 tmp.clear();
662 tmp.push_back(order0);
663 tmp.push_back(order1);
664 tmp.push_back(m_faces[1]->GetXmap()->GetTraceNcoeffs(1));
665 tmp.push_back(m_faces[1]->GetXmap()->GetTraceNcoeffs(2));
666 tmp.push_back(m_faces[3]->GetXmap()->GetTraceNcoeffs(1));
667 int order2 = *std::max_element(tmp.begin(), tmp.end());
668
669 const LibUtilities::BasisKey A(
671 LibUtilities::PointsKey(order0 + 1,
673 const LibUtilities::BasisKey B(
675 LibUtilities::PointsKey(order1, LibUtilities::eGaussRadauMAlpha1Beta0));
676 const LibUtilities::BasisKey C(
678 LibUtilities::PointsKey(order2, LibUtilities::eGaussRadauMAlpha2Beta0));
679
681}
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
StdRegions::StdExpansionSharedPtr m_xmap
mapping containing isoparametric transformation.
Definition: Geometry.h:194
StdRegions::StdExpansionSharedPtr GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition: Geometry.h:435
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:51
@ eModified_B
Principle Modified Functions .
Definition: BasisType.h:49
@ eModified_C
Principle Modified Functions .
Definition: BasisType.h:50
@ eModified_A
Principle Modified Functions .
Definition: BasisType.h:48

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModified_B, Nektar::LibUtilities::eModified_C, Nektar::SpatialDomains::Geometry::GetXmap(), Nektar::SpatialDomains::Geometry3D::m_faces, and Nektar::SpatialDomains::Geometry::m_xmap.

Referenced by v_Setup().

◆ v_GenGeomFactors()

void Nektar::SpatialDomains::TetGeom::v_GenGeomFactors ( )
overrideprotectedvirtual

Generate the geometry factors for this element.

Implements Nektar::SpatialDomains::Geometry.

Definition at line 597 of file TetGeom.cpp.

598{
599 if (!m_setupState)
600 {
602 }
603
605 {
606 GeomType Gtype = eRegular;
607
608 v_FillGeom();
609
610 // check to see if expansions are linear
611 m_straightEdge = 1;
612 if (m_xmap->GetBasisNumModes(0) != 2 ||
613 m_xmap->GetBasisNumModes(1) != 2 ||
614 m_xmap->GetBasisNumModes(2) != 2)
615 {
616 Gtype = eDeformed;
617 m_straightEdge = 0;
618 }
619
620 if (Gtype == eRegular)
621 {
622 m_isoParameter = Array<OneD, Array<OneD, NekDouble>>(3);
623 for (int i = 0; i < 3; ++i)
624 {
625 m_isoParameter[i] = Array<OneD, NekDouble>(4, 0.);
626 NekDouble A = (*m_verts[0])(i);
627 NekDouble B = (*m_verts[1])(i);
628 NekDouble C = (*m_verts[2])(i);
629 NekDouble D = (*m_verts[3])(i);
630 m_isoParameter[i][0] = 0.5 * (-A + B + C + D);
631
632 m_isoParameter[i][1] = 0.5 * (-A + B); // xi1
633 m_isoParameter[i][2] = 0.5 * (-A + C); // xi2
634 m_isoParameter[i][3] = 0.5 * (-A + D); // xi3
635 }
637 }
638
640 Gtype, m_coordim, m_xmap, m_coeffs);
642 }
643}
void v_CalculateInverseIsoParam() override
Definition: Geometry3D.cpp:565
void v_FillGeom() override
Put all quadrature information into face/edge structure and backward transform.
Definition: Geometry3D.cpp:442
bool m_setupState
Wether or not the setup routines have been run.
Definition: Geometry.h:198
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
Definition: Geometry.h:209
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
Definition: Geometry.h:206
GeomState m_geomFactorsState
State of the geometric factors.
Definition: Geometry.h:192
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
Definition: Geometry.h:190
GeomType
Indicates the type of element geometry.
@ eRegular
Geometry is straight-sided with constant geometric factors.
@ eDeformed
Geometry is curved or has non-constant factors.
@ ePtsFilled
Geometric information has been generated.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::SpatialDomains::eDeformed, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::eRegular, Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry::m_geomFactors, Nektar::SpatialDomains::Geometry::m_geomFactorsState, Nektar::SpatialDomains::Geometry::m_isoParameter, Nektar::SpatialDomains::Geometry::m_setupState, Nektar::SpatialDomains::Geometry::m_straightEdge, Nektar::SpatialDomains::Geometry3D::m_verts, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry3D::v_CalculateInverseIsoParam(), Nektar::SpatialDomains::Geometry3D::v_FillGeom(), and v_Setup().

◆ v_GetDir()

int Nektar::SpatialDomains::TetGeom::v_GetDir ( const int  i,
const int  j 
) const
overrideprotectedvirtual

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 76 of file TetGeom.cpp.

77{
78 if (faceidx == 0)
79 {
80 return facedir;
81 }
82 else if (faceidx == 1)
83 {
84 return 2 * facedir;
85 }
86 else
87 {
88 return 1 + facedir;
89 }
90}

◆ v_GetEdgeFaceMap()

int Nektar::SpatialDomains::TetGeom::v_GetEdgeFaceMap ( const int  i,
const int  j 
) const
overrideprotectedvirtual

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 from Nektar::SpatialDomains::Geometry.

Definition at line 102 of file TetGeom.cpp.

103{
104 return EdgeFaceConnectivity[i][j];
105}
static const unsigned int EdgeFaceConnectivity[6][2]
Definition: TetGeom.h:78

References EdgeFaceConnectivity.

◆ v_GetEdgeNormalToFaceVert()

int Nektar::SpatialDomains::TetGeom::v_GetEdgeNormalToFaceVert ( const int  i,
const int  j 
) const
overrideprotectedvirtual

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 from Nektar::SpatialDomains::Geometry.

Definition at line 107 of file TetGeom.cpp.

108{
109 return EdgeNormalToFaceVert[i][j];
110}
static const unsigned int EdgeNormalToFaceVert[4][3]
Definition: TetGeom.h:79

References EdgeNormalToFaceVert.

◆ v_GetVertexEdgeMap()

int Nektar::SpatialDomains::TetGeom::v_GetVertexEdgeMap ( const int  i,
const int  j 
) const
overrideprotectedvirtual

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 from Nektar::SpatialDomains::Geometry.

Definition at line 92 of file TetGeom.cpp.

93{
94 return VertexEdgeConnectivity[i][j];
95}
static const unsigned int VertexEdgeConnectivity[4][3]
Definition: TetGeom.h:76

References VertexEdgeConnectivity.

◆ v_GetVertexFaceMap()

int Nektar::SpatialDomains::TetGeom::v_GetVertexFaceMap ( const int  i,
const int  j 
) const
overrideprotectedvirtual

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 from Nektar::SpatialDomains::Geometry.

Definition at line 97 of file TetGeom.cpp.

98{
99 return VertexFaceConnectivity[i][j];
100}
static const unsigned int VertexFaceConnectivity[4][3]
Definition: TetGeom.h:77

References VertexFaceConnectivity.

◆ v_Reset()

void Nektar::SpatialDomains::TetGeom::v_Reset ( CurveMap curvedEdges,
CurveMap curvedFaces 
)
overrideprotectedvirtual

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 570 of file TetGeom.cpp.

571{
572 Geometry::v_Reset(curvedEdges, curvedFaces);
573
574 for (int i = 0; i < 4; ++i)
575 {
576 m_faces[i]->Reset(curvedEdges, curvedFaces);
577 }
578}
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

References Nektar::SpatialDomains::Geometry3D::m_faces, and Nektar::SpatialDomains::Geometry::v_Reset().

◆ v_Setup()

void Nektar::SpatialDomains::TetGeom::v_Setup ( )
overrideprotectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 580 of file TetGeom.cpp.

581{
582 if (!m_setupState)
583 {
584 for (int i = 0; i < 4; ++i)
585 {
586 m_faces[i]->Setup();
587 }
588 SetUpXmap();
589 SetUpCoeffs(m_xmap->GetNcoeffs());
590 m_setupState = true;
591 }
592}
void SetUpCoeffs(const int nCoeffs)
Initialise the Geometry::m_coeffs array.
Definition: Geometry.h:700
void SetUpXmap()
Set up the m_xmap object by determining the order of each direction from derived faces.
Definition: TetGeom.cpp:649

References Nektar::SpatialDomains::Geometry3D::m_faces, Nektar::SpatialDomains::Geometry::m_setupState, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry::SetUpCoeffs(), and SetUpXmap().

Referenced by v_GenGeomFactors().

Member Data Documentation

◆ EdgeFaceConnectivity

const unsigned int Nektar::SpatialDomains::TetGeom::EdgeFaceConnectivity
staticprivate
Initial value:
= {
{0, 1}, {0, 2}, {0, 3}, {1, 3}, {1, 2}, {2, 3}}

Definition at line 78 of file TetGeom.h.

Referenced by v_GetEdgeFaceMap().

◆ EdgeNormalToFaceVert

const unsigned int Nektar::SpatialDomains::TetGeom::EdgeNormalToFaceVert
staticprivate
Initial value:
= {
{3, 4, 5}, {1, 2, 5}, {0, 2, 3}, {0, 1, 4}}

Definition at line 79 of file TetGeom.h.

Referenced by v_GetEdgeNormalToFaceVert().

◆ kNedges

const int Nektar::SpatialDomains::TetGeom::kNedges = 6
static

Definition at line 53 of file TetGeom.h.

Referenced by SetUpEdgeOrientation(), and TetGeom().

◆ kNfaces

const int Nektar::SpatialDomains::TetGeom::kNfaces = kNqfaces + kNtfaces
static

◆ kNqfaces

const int Nektar::SpatialDomains::TetGeom::kNqfaces = 0
static

◆ kNtfaces

const int Nektar::SpatialDomains::TetGeom::kNtfaces = 4
static

◆ kNverts

const int Nektar::SpatialDomains::TetGeom::kNverts = 4
static

Definition at line 52 of file TetGeom.h.

◆ VertexEdgeConnectivity

const unsigned int Nektar::SpatialDomains::TetGeom::VertexEdgeConnectivity
staticprivate
Initial value:
= {
{0, 2, 3}, {0, 1, 4}, {1, 2, 5}, {3, 4, 5}}

Definition at line 76 of file TetGeom.h.

Referenced by v_GetVertexEdgeMap().

◆ VertexFaceConnectivity

const unsigned int Nektar::SpatialDomains::TetGeom::VertexFaceConnectivity
staticprivate
Initial value:
= {
{0, 1, 3}, {0, 1, 2}, {0, 2, 3}, {1, 2, 3}}

Definition at line 77 of file TetGeom.h.

Referenced by v_GetVertexFaceMap().

◆ XMLElementType

const std::string Nektar::SpatialDomains::TetGeom::XMLElementType
static

Definition at line 57 of file TetGeom.h.