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

#include <TriGeom.h>

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

Public Member Functions

 TriGeom ()
 
 TriGeom (const TriGeom &in)
 
 TriGeom (const int id, const SegGeomSharedPtr edges[], const CurveSharedPtr curve=CurveSharedPtr())
 
 ~TriGeom () override=default
 
- Public Member Functions inherited from Nektar::SpatialDomains::Geometry2D
 Geometry2D ()
 
 Geometry2D (const int coordim, CurveSharedPtr curve)
 
 ~Geometry2D () override=default
 
CurveSharedPtr GetCurve ()
 
- 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 Member Functions

static StdRegions::Orientation GetFaceOrientation (const TriGeom &face1, const TriGeom &face2, bool doRot, int dir, NekDouble angle, NekDouble tol)
 
static StdRegions::Orientation GetFaceOrientation (const PointGeomVector &face1, const PointGeomVector &face2, bool doRot, int dir, NekDouble angle, NekDouble tol)
 

Static Public Attributes

static const int kNedges = 3
 Get the orientation of face1. More...
 
static const int kNverts = 3
 
- Static Public Attributes inherited from Nektar::SpatialDomains::Geometry2D
static const int kDim = 2
 

Protected Member Functions

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_GenGeomFactors () override
 
void v_FillGeom () override
 
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
 
- Protected Member Functions inherited from Nektar::SpatialDomains::Geometry2D
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, Array< OneD, NekDouble > &Lcoords, NekDouble &dist)
 
void SolveStraightEdgeQuad (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
 
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...
 
PointGeomSharedPtr v_GetVertex (int i) const override
 
Geometry1DSharedPtr v_GetEdge (int i) const override
 Returns edge i of this object. 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...
 
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...
 
NekDouble v_FindDistance (const Array< OneD, const NekDouble > &xs, Array< OneD, NekDouble > &xi) override
 
- 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 SetUpXmap ()
 

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::Geometry2D
PointGeomVector m_verts
 
SegGeomVector m_edges
 
std::vector< StdRegions::Orientationm_eorient
 
CurveSharedPtr m_curve
 
Array< OneD, int > m_manifold
 
Array< OneD, Array< OneD, NekDouble > > m_edgeNormal
 
- 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 59 of file TriGeom.h.

Constructor & Destructor Documentation

◆ TriGeom() [1/3]

Nektar::SpatialDomains::TriGeom::TriGeom ( )

◆ TriGeom() [2/3]

Nektar::SpatialDomains::TriGeom::TriGeom ( const TriGeom in)

Definition at line 81 of file TriGeom.cpp.

81 : Geometry2D(in)
82{
83 // From Geometry
84 m_shapeType = in.m_shapeType;
85
86 // From TriFaceComponent
87 m_globalID = in.m_globalID;
88
89 // From TriGeom
90 m_verts = in.m_verts;
91 m_edges = in.m_edges;
92 for (int i = 0; i < kNedges; i++)
93 {
94 m_eorient[i] = in.m_eorient[i];
95 }
96}
std::vector< StdRegions::Orientation > m_eorient
Definition: Geometry2D.h:81
static const int kNedges
Get the orientation of face1.
Definition: TriGeom.h:70

References kNedges, Nektar::SpatialDomains::Geometry2D::m_edges, Nektar::SpatialDomains::Geometry2D::m_eorient, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_shapeType, and Nektar::SpatialDomains::Geometry2D::m_verts.

◆ TriGeom() [3/3]

Nektar::SpatialDomains::TriGeom::TriGeom ( const int  id,
const SegGeomSharedPtr  edges[],
const CurveSharedPtr  curve = CurveSharedPtr() 
)

Definition at line 52 of file TriGeom.cpp.

54 : Geometry2D(edges[0]->GetVertex(0)->GetCoordim(), curve)
55{
56 int j;
57
59 m_globalID = id;
60
61 // Copy the edge shared pointers.
62 m_edges.insert(m_edges.begin(), edges, edges + TriGeom::kNedges);
63 m_eorient.resize(kNedges);
64
65 for (j = 0; j < kNedges; ++j)
66 {
67 m_eorient[j] =
68 SegGeom::GetEdgeOrientation(*edges[j], *edges[(j + 1) % kNedges]);
69 m_verts.push_back(
70 edges[j]->GetVertex(m_eorient[j] == StdRegions::eForwards ? 0 : 1));
71 }
72
76
77 m_coordim = edges[0]->GetVertex(0)->GetCoordim();
78 ASSERTL0(m_coordim > 1, "Cannot call function with dim == 1");
79}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
PointGeomSharedPtr GetVertex(int i) const
Returns vertex i of this object.
Definition: Geometry.h:357
int GetCoordim() const
Return the coordinate dimension of this object (i.e. the dimension of the space in which this object ...
Definition: Geometry.h:283
int m_coordim
Coordinate dimension of this geometry object.
Definition: Geometry.h:188
static StdRegions::Orientation GetEdgeOrientation(const SegGeom &edge1, const SegGeom &edge2)
Get the orientation of edge1.
Definition: SegGeom.cpp:190

References ASSERTL0, Nektar::StdRegions::eBackwards, Nektar::StdRegions::eForwards, Nektar::LibUtilities::eTriangle, Nektar::SpatialDomains::SegGeom::GetEdgeOrientation(), Nektar::SpatialDomains::Geometry::GetVertex(), kNedges, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry2D::m_edges, Nektar::SpatialDomains::Geometry2D::m_eorient, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_shapeType, and Nektar::SpatialDomains::Geometry2D::m_verts.

◆ ~TriGeom()

Nektar::SpatialDomains::TriGeom::~TriGeom ( )
overridedefault

Member Function Documentation

◆ GetFaceOrientation() [1/2]

StdRegions::Orientation Nektar::SpatialDomains::TriGeom::GetFaceOrientation ( const PointGeomVector face1,
const PointGeomVector face2,
bool  doRot,
int  dir,
NekDouble  angle,
NekDouble  tol 
)
static

Definition at line 119 of file TriGeom.cpp.

122{
123 int i, j, vmap[3] = {-1, -1, -1};
124
125 if (doRot)
126 {
127 PointGeom rotPt;
128
129 for (i = 0; i < 3; ++i)
130 {
131 rotPt.Rotate((*face1[i]), dir, angle);
132 for (j = 0; j < 3; ++j)
133 {
134 if (rotPt.dist(*face2[j]) < tol)
135 {
136 vmap[j] = i;
137 break;
138 }
139 }
140 }
141 }
142 else
143 {
144
145 NekDouble x, y, z, x1, y1, z1, cx = 0.0, cy = 0.0, cz = 0.0;
146
147 // For periodic faces, we calculate the vector between the centre
148 // points of the two faces. (For connected faces this will be
149 // zero). We can then use this to determine alignment later in the
150 // algorithm.
151 for (i = 0; i < 3; ++i)
152 {
153 cx += (*face2[i])(0) - (*face1[i])(0);
154 cy += (*face2[i])(1) - (*face1[i])(1);
155 cz += (*face2[i])(2) - (*face1[i])(2);
156 }
157 cx /= 3;
158 cy /= 3;
159 cz /= 3;
160
161 // Now construct a mapping which takes us from the vertices of one
162 // face to the other. That is, vertex j of face2 corresponds to
163 // vertex vmap[j] of face1.
164 for (i = 0; i < 3; ++i)
165 {
166 x = (*face1[i])(0);
167 y = (*face1[i])(1);
168 z = (*face1[i])(2);
169 for (j = 0; j < 3; ++j)
170 {
171 x1 = (*face2[j])(0) - cx;
172 y1 = (*face2[j])(1) - cy;
173 z1 = (*face2[j])(2) - cz;
174 if (sqrt((x1 - x) * (x1 - x) + (y1 - y) * (y1 - y) +
175 (z1 - z) * (z1 - z)) < 1e-8)
176 {
177 vmap[j] = i;
178 break;
179 }
180 }
181 }
182 }
183
184 if (vmap[1] == (vmap[0] + 1) % 3)
185 {
186 switch (vmap[0])
187 {
188 case 0:
190 break;
191 case 1:
193 break;
194 case 2:
196 break;
197 }
198 }
199 else
200 {
201 switch (vmap[0])
202 {
203 case 0:
205 break;
206 case 1:
208 break;
209 case 2:
211 break;
212 }
213 }
214
215 ASSERTL0(false, "Unable to determine triangle orientation");
217}
std::vector< double > z(NPUPPER)
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References ASSERTL0, Nektar::SpatialDomains::PointGeom::dist(), Nektar::StdRegions::eDir1BwdDir1_Dir2BwdDir2, Nektar::StdRegions::eDir1BwdDir1_Dir2FwdDir2, Nektar::StdRegions::eDir1BwdDir2_Dir2BwdDir1, Nektar::StdRegions::eDir1FwdDir1_Dir2FwdDir2, Nektar::StdRegions::eDir1FwdDir2_Dir2BwdDir1, Nektar::StdRegions::eDir1FwdDir2_Dir2FwdDir1, Nektar::StdRegions::eNoOrientation, Nektar::SpatialDomains::PointGeom::Rotate(), tinysimd::sqrt(), and Nektar::UnitTests::z().

◆ GetFaceOrientation() [2/2]

StdRegions::Orientation Nektar::SpatialDomains::TriGeom::GetFaceOrientation ( const TriGeom face1,
const TriGeom face2,
bool  doRot,
int  dir,
NekDouble  angle,
NekDouble  tol 
)
static

Definition at line 109 of file TriGeom.cpp.

114{
115 return GetFaceOrientation(face1.m_verts, face2.m_verts, doRot, dir, angle,
116 tol);
117}
static StdRegions::Orientation GetFaceOrientation(const TriGeom &face1, const TriGeom &face2, bool doRot, int dir, NekDouble angle, NekDouble tol)
Definition: TriGeom.cpp:109

References GetFaceOrientation(), and Nektar::SpatialDomains::Geometry2D::m_verts.

Referenced by Nektar::MultiRegions::DisContField::FindPeriodicTraces(), and GetFaceOrientation().

◆ SetUpXmap()

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

Definition at line 550 of file TriGeom.cpp.

551{
552 int order0 = m_edges[0]->GetXmap()->GetBasis(0)->GetNumModes();
553 int order1 = std::max(
554 order0, std::max(m_edges[1]->GetXmap()->GetBasis(0)->GetNumModes(),
555 m_edges[2]->GetXmap()->GetBasis(0)->GetNumModes()));
556
557 const LibUtilities::BasisKey B0(
559 LibUtilities::PointsKey(order0 + 1,
561 const LibUtilities::BasisKey B1(
563 LibUtilities::PointsKey(order1, LibUtilities::eGaussRadauMAlpha1Beta0));
564
566}
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_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::SpatialDomains::Geometry::GetXmap(), Nektar::SpatialDomains::Geometry2D::m_edges, and Nektar::SpatialDomains::Geometry::m_xmap.

Referenced by v_Reset(), and v_Setup().

◆ v_FillGeom()

void Nektar::SpatialDomains::TriGeom::v_FillGeom ( )
overrideprotectedvirtual

Note verts and edges are listed according to anticlockwise convention but points in _coeffs have to be in array format from left to right.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 302 of file TriGeom.cpp.

303{
304 // check to see if geometry structure is already filled
305 if (m_state == ePtsFilled)
306 {
307 return;
308 }
309
310 int i, j, k;
311 int nEdgeCoeffs = m_xmap->GetTraceNcoeffs(0);
312
313 if (m_curve)
314 {
315 int pdim = LibUtilities::PointsManager()[LibUtilities::PointsKey(
316 2, m_curve->m_ptype)]
317 ->GetPointsDim();
318
319 // Deal with 2D points type separately
320 // (e.g. electrostatic or Fekete points) to 1D tensor
321 // product.
322 if (pdim == 2)
323 {
324 int N = m_curve->m_points.size();
325 int nEdgePts =
326 (-1 + (int)sqrt(static_cast<NekDouble>(8 * N + 1))) / 2;
327
328 ASSERTL0(nEdgePts * (nEdgePts + 1) / 2 == N,
329 "NUMPOINTS should be a triangle number for"
330 " triangle curved face " +
331 std::to_string(m_globalID));
332
333 // Sanity check 1: are curved vertices consistent with
334 // triangle vertices?
335 for (i = 0; i < 3; ++i)
336 {
337 NekDouble dist = m_verts[i]->dist(*(m_curve->m_points[i]));
339 {
340 std::stringstream ss;
341 ss << "Curved vertex " << i << " of triangle " << m_globalID
342 << " is separated from expansion vertex by"
343 << " more than " << NekConstants::kVertexTheSameDouble
344 << " (dist = " << dist << ")";
345 NEKERROR(ErrorUtil::ewarning, ss.str().c_str());
346 }
347 }
348
349 // Sanity check 2: are curved edges from the face curvature
350 // consistent with curved edges?
351 for (i = 0; i < kNedges; ++i)
352 {
353 CurveSharedPtr edgeCurve = m_edges[i]->GetCurve();
354
355 ASSERTL0(edgeCurve->m_points.size() == nEdgePts,
356 "Number of edge points does not correspond "
357 "to number of face points in triangle " +
358 std::to_string(m_globalID));
359
360 const int offset = 3 + i * (nEdgePts - 2);
361 NekDouble maxDist = 0.0;
362
363 // Account for different ordering of nodal coordinates
364 // vs. Cartesian ordering of element.
366
367 if (i == 2)
368 {
369 orient = orient == StdRegions::eForwards
372 }
373
374 if (orient == StdRegions::eForwards)
375 {
376 for (j = 0; j < nEdgePts - 2; ++j)
377 {
378 NekDouble dist = m_curve->m_points[offset + j]->dist(
379 *(edgeCurve->m_points[j + 1]));
380 maxDist = dist > maxDist ? dist : maxDist;
381 }
382 }
383 else
384 {
385 for (j = 0; j < nEdgePts - 2; ++j)
386 {
387 NekDouble dist = m_curve->m_points[offset + j]->dist(
388 *(edgeCurve->m_points[nEdgePts - 2 - j]));
389 maxDist = dist > maxDist ? dist : maxDist;
390 }
391 }
392
394 {
395 std::stringstream ss;
396 ss << "Curved edge " << i << " of triangle " << m_globalID
397 << " has a point separated from edge interior"
398 << " points by more than "
400 << " (maxdist = " << maxDist << ")";
401 NEKERROR(ErrorUtil::ewarning, ss.str().c_str());
402 }
403 }
404
405 const LibUtilities::PointsKey P0(
407 const LibUtilities::PointsKey P1(
408 nEdgePts, LibUtilities::eGaussRadauMAlpha1Beta0);
409 const LibUtilities::BasisKey T0(LibUtilities::eOrtho_A, nEdgePts,
410 P0);
411 const LibUtilities::BasisKey T1(LibUtilities::eOrtho_B, nEdgePts,
412 P1);
413 Array<OneD, NekDouble> phys(
414 std::max(nEdgePts * nEdgePts, m_xmap->GetTotPoints()));
415 Array<OneD, NekDouble> tmp(nEdgePts * nEdgePts);
416
417 for (i = 0; i < m_coordim; ++i)
418 {
419 // Create a StdNodalTriExp.
422 AllocateSharedPtr(T0, T1, m_curve->m_ptype);
423
424 for (j = 0; j < N; ++j)
425 {
426 phys[j] = (m_curve->m_points[j]->GetPtr())[i];
427 }
428
429 t->BwdTrans(phys, tmp);
430
431 // Interpolate points to standard region.
433 P0, P1, tmp, m_xmap->GetBasis(0)->GetPointsKey(),
434 m_xmap->GetBasis(1)->GetPointsKey(), phys);
435
436 // Forwards transform to get coefficient space.
437 m_xmap->FwdTrans(phys, m_coeffs[i]);
438 }
439 }
440 else if (pdim == 1)
441 {
442 int npts = m_curve->m_points.size();
443 int nEdgePts = (int)sqrt(static_cast<NekDouble>(npts));
444 Array<OneD, NekDouble> tmp(npts);
445 Array<OneD, NekDouble> phys(m_xmap->GetTotPoints());
446 LibUtilities::PointsKey curveKey(nEdgePts, m_curve->m_ptype);
447
448 // Sanity checks:
449 // - Curved faces should have square number of points;
450 // - Each edge should have sqrt(npts) points.
451 ASSERTL0(nEdgePts * nEdgePts == npts,
452 "NUMPOINTS should be a square number for"
453 " triangle " +
454 std::to_string(m_globalID));
455
456 for (i = 0; i < kNedges; ++i)
457 {
458 ASSERTL0(m_edges[i]->GetXmap()->GetNcoeffs() == nEdgePts,
459 "Number of edge points does not correspond to "
460 "number of face points in triangle " +
461 std::to_string(m_globalID));
462 }
463
464 for (i = 0; i < m_coordim; ++i)
465 {
466 for (j = 0; j < npts; ++j)
467 {
468 tmp[j] = (m_curve->m_points[j]->GetPtr())[i];
469 }
470
471 // Interpolate curve points to standard triangle
472 // points.
473 LibUtilities::Interp2D(curveKey, curveKey, tmp,
474 m_xmap->GetBasis(0)->GetPointsKey(),
475 m_xmap->GetBasis(1)->GetPointsKey(),
476 phys);
477
478 // Forwards transform to get coefficient space.
479 m_xmap->FwdTrans(phys, m_coeffs[i]);
480 }
481 }
482 else
483 {
484 ASSERTL0(false, "Only 1D/2D points distributions "
485 "supported.");
486 }
487 }
488
489 Array<OneD, unsigned int> mapArray(nEdgeCoeffs);
490 Array<OneD, int> signArray(nEdgeCoeffs);
491
492 for (i = 0; i < kNedges; i++)
493 {
494 m_edges[i]->FillGeom();
495 m_xmap->GetTraceToElementMap(i, mapArray, signArray, m_eorient[i]);
496
497 nEdgeCoeffs = m_edges[i]->GetXmap()->GetNcoeffs();
498
499 for (j = 0; j < m_coordim; j++)
500 {
501 for (k = 0; k < nEdgeCoeffs; k++)
502 {
503 m_coeffs[j][mapArray[k]] =
504 signArray[k] * m_edges[i]->GetCoeffs(j)[k];
505 }
506 }
507 }
508
510}
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:202
GeomState m_state
Enumeration to dictate whether coefficients are filled.
Definition: Geometry.h:196
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
Definition: Geometry.h:206
PointsManagerT & PointsManager(void)
void Interp2D(const BasisKey &fbasis0, const BasisKey &fbasis1, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, Array< OneD, NekDouble > &to)
this function interpolates a 2D function evaluated at the quadrature points of the 2D basis,...
Definition: Interp.cpp:101
@ eOrtho_A
Principle Orthogonal Functions .
Definition: BasisType.h:42
@ eOrtho_B
Principle Orthogonal Functions .
Definition: BasisType.h:44
static const NekDouble kVertexTheSameDouble
std::shared_ptr< Curve > CurveSharedPtr
Definition: Curve.hpp:58
@ ePtsFilled
Geometric information has been generated.
std::shared_ptr< StdNodalTriExp > StdNodalTriExpSharedPtr

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::StdRegions::eBackwards, Nektar::StdRegions::eForwards, Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eOrtho_A, Nektar::LibUtilities::eOrtho_B, Nektar::SpatialDomains::ePtsFilled, Nektar::ErrorUtil::ewarning, Nektar::SpatialDomains::Geometry::GetXmap(), Nektar::LibUtilities::Interp2D(), kNedges, Nektar::NekConstants::kVertexTheSameDouble, Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry2D::m_curve, Nektar::SpatialDomains::Geometry2D::m_edges, Nektar::SpatialDomains::Geometry2D::m_eorient, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_state, Nektar::SpatialDomains::Geometry2D::m_verts, Nektar::SpatialDomains::Geometry::m_xmap, NEKERROR, Nektar::LibUtilities::PointsManager(), and tinysimd::sqrt().

Referenced by v_GenGeomFactors().

◆ v_GenGeomFactors()

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

Implements Nektar::SpatialDomains::Geometry.

Definition at line 219 of file TriGeom.cpp.

220{
221 if (!m_setupState)
222 {
224 }
225
227 {
228 GeomType Gtype = eRegular;
229
231
232 // check to see if expansions are linear
233 m_straightEdge = 1;
234 if (m_xmap->GetBasisNumModes(0) != 2 ||
235 m_xmap->GetBasisNumModes(1) != 2)
236 {
237 Gtype = eDeformed;
238 m_straightEdge = 0;
239 }
240
241 m_manifold = Array<OneD, int>(m_coordim);
242 m_manifold[0] = 0;
243 m_manifold[1] = 1;
244 if (m_coordim == 3)
245 {
246 PointGeom e01, e21, norm;
247 e01.Sub(*m_verts[0], *m_verts[1]);
248 e21.Sub(*m_verts[2], *m_verts[1]);
249 norm.Mult(e01, e21);
250 int tmpi = 0;
251 double tmp = std::fabs(norm[0]);
252 if (tmp < fabs(norm[1]))
253 {
254 tmp = fabs(norm[1]);
255 tmpi = 1;
256 }
257 if (tmp < fabs(norm[2]))
258 {
259 tmpi = 2;
260 }
261 m_manifold[0] = (tmpi + 1) % 3;
262 m_manifold[1] = (tmpi + 2) % 3;
263 m_manifold[2] = (tmpi + 3) % 3;
264 }
265 if (Gtype == eRegular)
266 {
267 Array<OneD, Array<OneD, NekDouble>> verts(m_verts.size());
268 for (int i = 0; i < m_verts.size(); ++i)
269 {
270 verts[i] = Array<OneD, NekDouble>(3);
271 m_verts[i]->GetCoords(verts[i]);
272 }
273 // a00 + a01 xi1 + a02 xi2
274 // a10 + a11 xi1 + a12 xi2
275 m_isoParameter = Array<OneD, Array<OneD, NekDouble>>(2);
276 for (int i = 0; i < 2; ++i)
277 {
278 unsigned int d = m_manifold[i];
279 m_isoParameter[i] = Array<OneD, NekDouble>(3, 0.);
280 NekDouble A = verts[0][d];
281 NekDouble B = verts[1][d];
282 NekDouble C = verts[2][d];
283 m_isoParameter[i][0] = 0.5 * (B + C); // 1
284 m_isoParameter[i][1] = 0.5 * (-A + B); // xi1
285 m_isoParameter[i][2] = 0.5 * (-A + C); // xi2
286 }
288 }
289
291 Gtype, m_coordim, m_xmap, m_coeffs);
292
294 }
295}
Array< OneD, int > m_manifold
Definition: Geometry2D.h:83
void v_CalculateInverseIsoParam() override
Definition: Geometry2D.cpp:629
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
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.
std::vector< double > d(NPUPPER *NPUPPER)

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::UnitTests::d(), 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::Geometry2D::m_manifold, Nektar::SpatialDomains::Geometry::m_setupState, Nektar::SpatialDomains::Geometry::m_straightEdge, Nektar::SpatialDomains::Geometry2D::m_verts, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::PointGeom::Mult(), Nektar::SpatialDomains::PointGeom::Sub(), Nektar::SpatialDomains::Geometry2D::v_CalculateInverseIsoParam(), v_FillGeom(), and v_Setup().

◆ v_GetCoord()

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

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 98 of file TriGeom.cpp.

100{
101 ASSERTL1(m_state == ePtsFilled, "Geometry is not in physical space");
102
103 Array<OneD, NekDouble> tmp(m_xmap->GetTotPoints());
104 m_xmap->BwdTrans(m_coeffs[i], tmp);
105
106 return m_xmap->PhysEvaluate(Lcoord, tmp);
107}
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242

References ASSERTL1, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_state, and Nektar::SpatialDomains::Geometry::m_xmap.

◆ v_GetDir()

int Nektar::SpatialDomains::TriGeom::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 512 of file TriGeom.cpp.

513{
514 return i == 0 ? 0 : 1;
515}

◆ v_Reset()

void Nektar::SpatialDomains::TriGeom::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 517 of file TriGeom.cpp.

518{
519 Geometry::v_Reset(curvedEdges, curvedFaces);
520 CurveMap::iterator it = curvedFaces.find(m_globalID);
521
522 if (it != curvedFaces.end())
523 {
524 m_curve = it->second;
525 }
526
527 for (int i = 0; i < 3; ++i)
528 {
529 m_edges[i]->Reset(curvedEdges, curvedFaces);
530 }
531
532 SetUpXmap();
533 SetUpCoeffs(m_xmap->GetNcoeffs());
534}
void SetUpCoeffs(const int nCoeffs)
Initialise the Geometry::m_coeffs array.
Definition: Geometry.h:700
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::Geometry2D::m_curve, Nektar::SpatialDomains::Geometry2D::m_edges, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry::SetUpCoeffs(), SetUpXmap(), and Nektar::SpatialDomains::Geometry::v_Reset().

◆ v_Setup()

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 536 of file TriGeom.cpp.

537{
538 if (!m_setupState)
539 {
540 for (int i = 0; i < 3; ++i)
541 {
542 m_edges[i]->Setup();
543 }
544 SetUpXmap();
545 SetUpCoeffs(m_xmap->GetNcoeffs());
546 m_setupState = true;
547 }
548}

References Nektar::SpatialDomains::Geometry2D::m_edges, Nektar::SpatialDomains::Geometry::m_setupState, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry::SetUpCoeffs(), and SetUpXmap().

Referenced by v_GenGeomFactors().

Member Data Documentation

◆ kNedges

const int Nektar::SpatialDomains::TriGeom::kNedges = 3
static

◆ kNverts

const int Nektar::SpatialDomains::TriGeom::kNverts = 3
static

Definition at line 71 of file TriGeom.h.

Referenced by Nektar::SpatialDomains::TetGeom::SetUpFaceOrientation().