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

#include <QuadGeom.h>

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

Public Member Functions

 QuadGeom ()
 
 QuadGeom (const QuadGeom &in)
 
 QuadGeom (const int id, const SegGeomSharedPtr edges[], const CurveSharedPtr curve=CurveSharedPtr())
 
 ~QuadGeom () 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 QuadGeom &face1, const QuadGeom &face2, bool doRot=false, int dir=0, NekDouble angle=0.0, NekDouble tol=1e-8)
 Get the orientation of face1. More...
 
static StdRegions::Orientation GetFaceOrientation (const PointGeomVector &face1, const PointGeomVector &face2, bool doRot=false, int dir=0, NekDouble angle=0.0, NekDouble tol=1e-8)
 

Static Public Attributes

static const int kNverts = 4
 
static const int kNedges = 4
 
static const std::string XMLElementType
 
- 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
 
void PreSolveStraightEdge ()
 
- 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 55 of file QuadGeom.h.

Constructor & Destructor Documentation

◆ QuadGeom() [1/3]

Nektar::SpatialDomains::QuadGeom::QuadGeom ( )

◆ QuadGeom() [2/3]

Nektar::SpatialDomains::QuadGeom::QuadGeom ( const QuadGeom in)

Definition at line 84 of file QuadGeom.cpp.

84 : Geometry2D(in)
85{
86 // From Geometry
87 m_shapeType = in.m_shapeType;
88 m_globalID = in.m_globalID;
89
90 // From QuadGeom
91 m_verts = in.m_verts;
92 m_edges = in.m_edges;
93 for (int i = 0; i < kNedges; i++)
94 {
95 m_eorient[i] = in.m_eorient[i];
96 }
97}
std::vector< StdRegions::Orientation > m_eorient
Definition: Geometry2D.h:81
static const int kNedges
Definition: QuadGeom.h:75

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.

◆ QuadGeom() [3/3]

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

Copy the edge shared pointers.

Definition at line 52 of file QuadGeom.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 + QuadGeom::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
73 for (j = 2; j < kNedges; ++j)
74 {
78 }
79
80 m_coordim = edges[0]->GetVertex(0)->GetCoordim();
81 ASSERTL0(m_coordim > 1, "Cannot call function with dim == 1");
82}
#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::eQuadrilateral, 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.

◆ ~QuadGeom()

Nektar::SpatialDomains::QuadGeom::~QuadGeom ( )
overridedefault

Member Function Documentation

◆ GetFaceOrientation() [1/2]

StdRegions::Orientation Nektar::SpatialDomains::QuadGeom::GetFaceOrientation ( const PointGeomVector face1,
const PointGeomVector face2,
bool  doRot = false,
int  dir = 0,
NekDouble  angle = 0.0,
NekDouble  tol = 1e-8 
)
static

Calculate the orientation of face2 to face1 (note this is not face1 to face2!).

Definition at line 143 of file QuadGeom.cpp.

146{
147 int i, j, vmap[4] = {-1, -1, -1, -1};
148
149 if (doRot)
150 {
151 PointGeom rotPt;
152
153 for (i = 0; i < 4; ++i)
154 {
155 rotPt.Rotate((*face1[i]), dir, angle);
156 for (j = 0; j < 4; ++j)
157 {
158 if (rotPt.dist(*face2[j]) < tol)
159 {
160 vmap[j] = i;
161 break;
162 }
163 }
164 }
165 }
166 else
167 {
168
169 NekDouble x, y, z, x1, y1, z1, cx = 0.0, cy = 0.0, cz = 0.0;
170
171 // For periodic faces, we calculate the vector between the centre
172 // points of the two faces. (For connected faces this will be
173 // zero). We can then use this to determine alignment later in the
174 // algorithm.
175 for (i = 0; i < 4; ++i)
176 {
177 cx += (*face2[i])(0) - (*face1[i])(0);
178 cy += (*face2[i])(1) - (*face1[i])(1);
179 cz += (*face2[i])(2) - (*face1[i])(2);
180 }
181 cx /= 4;
182 cy /= 4;
183 cz /= 4;
184
185 // Now construct a mapping which takes us from the vertices of one
186 // face to the other. That is, vertex j of face2 corresponds to
187 // vertex vmap[j] of face1.
188 for (i = 0; i < 4; ++i)
189 {
190 x = (*face1[i])(0);
191 y = (*face1[i])(1);
192 z = (*face1[i])(2);
193 for (j = 0; j < 4; ++j)
194 {
195 x1 = (*face2[j])(0) - cx;
196 y1 = (*face2[j])(1) - cy;
197 z1 = (*face2[j])(2) - cz;
198 if (sqrt((x1 - x) * (x1 - x) + (y1 - y) * (y1 - y) +
199 (z1 - z) * (z1 - z)) < 1e-8)
200 {
201 vmap[j] = i;
202 break;
203 }
204 }
205 }
206 }
207
208 // Use the mapping to determine the eight alignment options between
209 // faces.
210 if (vmap[1] == (vmap[0] + 1) % 4)
211 {
212 switch (vmap[0])
213 {
214 case 0:
216 break;
217 case 1:
219 break;
220 case 2:
222 break;
223 case 3:
225 break;
226 }
227 }
228 else
229 {
230 switch (vmap[0])
231 {
232 case 0:
234 break;
235 case 1:
237 break;
238 case 2:
240 break;
241 case 3:
243 break;
244 }
245 }
246
247 ASSERTL0(false, "unable to determine face orientation");
249}
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::eDir1BwdDir2_Dir2FwdDir1, Nektar::StdRegions::eDir1FwdDir1_Dir2BwdDir2, Nektar::StdRegions::eDir1FwdDir1_Dir2FwdDir2, Nektar::StdRegions::eDir1FwdDir2_Dir2BwdDir1, Nektar::StdRegions::eDir1FwdDir2_Dir2FwdDir1, Nektar::SpatialDomains::PointGeom::Rotate(), tinysimd::sqrt(), and Nektar::UnitTests::z().

◆ GetFaceOrientation() [2/2]

StdRegions::Orientation Nektar::SpatialDomains::QuadGeom::GetFaceOrientation ( const QuadGeom face1,
const QuadGeom face2,
bool  doRot = false,
int  dir = 0,
NekDouble  angle = 0.0,
NekDouble  tol = 1e-8 
)
static

Get the orientation of face1.

Definition at line 129 of file QuadGeom.cpp.

134{
135 return GetFaceOrientation(face1.m_verts, face2.m_verts, doRot, dir, angle,
136 tol);
137}
static StdRegions::Orientation GetFaceOrientation(const QuadGeom &face1, const QuadGeom &face2, bool doRot=false, int dir=0, NekDouble angle=0.0, NekDouble tol=1e-8)
Get the orientation of face1.
Definition: QuadGeom.cpp:129

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

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

◆ PreSolveStraightEdge()

void Nektar::SpatialDomains::QuadGeom::PreSolveStraightEdge ( )
protected

Definition at line 458 of file QuadGeom.cpp.

459{
460 int i0, i1, j1, j2;
461 if (fabs(m_isoParameter[0][3]) >= fabs(m_isoParameter[1][3]))
462 {
463 i0 = 0;
464 i1 = 1;
465 }
466 else
467 {
468 i1 = 0;
469 i0 = 1;
470 m_straightEdge |= 2;
471 }
472 NekDouble gamma = m_isoParameter[i1][3] / m_isoParameter[i0][3];
473 std::vector<NekDouble> c(3);
474 for (int i = 0; i < 3; ++i)
475 {
476 c[i] = m_isoParameter[i1][i] - gamma * m_isoParameter[i0][i];
477 }
478 if (fabs(c[1]) >= fabs(c[2]))
479 {
480 j1 = 1;
481 j2 = 2;
482 }
483 else
484 {
485 j1 = 2;
486 j2 = 1;
487 m_straightEdge |= 4;
488 }
489 NekDouble beta = c[j2] / c[j1];
490 if (i0 == 1)
491 {
493 }
494 if (j1 == 2)
495 {
496 NekDouble temp = m_isoParameter[0][j1];
497 m_isoParameter[0][j1] = m_isoParameter[0][j2];
498 m_isoParameter[0][j2] = temp;
499 }
500 m_isoParameter[0][2] -= m_isoParameter[0][1] * beta;
501 m_isoParameter[1][0] = c[0];
502 m_isoParameter[1][1] = 1. / c[j1];
503 m_isoParameter[1][2] = beta;
504 m_isoParameter[1][3] = gamma;
505}
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
Definition: Geometry.h:209
@ beta
Gauss Radau pinned at x=-1,.
Definition: PointsType.h:59
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

References Nektar::LibUtilities::beta, Nektar::SpatialDomains::Geometry::m_isoParameter, Nektar::SpatialDomains::Geometry::m_straightEdge, and Vmath::Vcopy().

Referenced by v_GenGeomFactors().

◆ SetUpXmap()

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

Definition at line 99 of file QuadGeom.cpp.

100{
101 int order0 = std::max(m_edges[0]->GetXmap()->GetBasis(0)->GetNumModes(),
102 m_edges[2]->GetXmap()->GetBasis(0)->GetNumModes());
103 int order1 = std::max(m_edges[1]->GetXmap()->GetBasis(0)->GetNumModes(),
104 m_edges[3]->GetXmap()->GetBasis(0)->GetNumModes());
105
106 const LibUtilities::BasisKey B0(
108 LibUtilities::PointsKey(order0 + 1,
110 const LibUtilities::BasisKey B1(
112 LibUtilities::PointsKey(order1 + 1,
114
116}
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_A
Principle Modified Functions .
Definition: BasisType.h:48

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eModified_A, 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::QuadGeom::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 383 of file QuadGeom.cpp.

384{
385 // check to see if geometry structure is already filled
386 if (m_state != ePtsFilled)
387 {
388 int i, j, k;
389 int nEdgeCoeffs;
390
391 if (m_curve)
392 {
393 int npts = m_curve->m_points.size();
394 int nEdgePts = (int)sqrt(static_cast<NekDouble>(npts));
395 Array<OneD, NekDouble> tmp(npts);
396 Array<OneD, NekDouble> tmp2(m_xmap->GetTotPoints());
397 LibUtilities::PointsKey curveKey(nEdgePts, m_curve->m_ptype);
398
399 // Sanity checks:
400 // - Curved faces should have square number of points;
401 // - Each edge should have sqrt(npts) points.
402 ASSERTL0(nEdgePts * nEdgePts == npts,
403 "NUMPOINTS should be a square number in"
404 " quadrilteral " +
405 std::to_string(m_globalID));
406
407 for (i = 0; i < kNedges; ++i)
408 {
409 ASSERTL0(m_edges[i]->GetXmap()->GetNcoeffs() == nEdgePts,
410 "Number of edge points does not correspond to "
411 "number of face points in quadrilateral " +
412 std::to_string(m_globalID));
413 }
414
415 for (i = 0; i < m_coordim; ++i)
416 {
417 for (j = 0; j < npts; ++j)
418 {
419 tmp[j] = (m_curve->m_points[j]->GetPtr())[i];
420 }
421
422 // Interpolate m_curve points to GLL points
423 LibUtilities::Interp2D(curveKey, curveKey, tmp,
424 m_xmap->GetBasis(0)->GetPointsKey(),
425 m_xmap->GetBasis(1)->GetPointsKey(),
426 tmp2);
427
428 // Forwards transform to get coefficient space.
429 m_xmap->FwdTrans(tmp2, m_coeffs[i]);
430 }
431 }
432
433 // Now fill in edges.
434 Array<OneD, unsigned int> mapArray;
435 Array<OneD, int> signArray;
436
437 for (i = 0; i < kNedges; i++)
438 {
439 m_edges[i]->FillGeom();
440 m_xmap->GetTraceToElementMap(i, mapArray, signArray, m_eorient[i]);
441
442 nEdgeCoeffs = m_edges[i]->GetXmap()->GetNcoeffs();
443
444 for (j = 0; j < m_coordim; j++)
445 {
446 for (k = 0; k < nEdgeCoeffs; k++)
447 {
448 m_coeffs[j][mapArray[k]] =
449 signArray[k] * (m_edges[i]->GetCoeffs(j))[k];
450 }
451 }
452 }
453
455 }
456}
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
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
@ ePtsFilled
Geometric information has been generated.

References ASSERTL0, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::Geometry::GetXmap(), Nektar::LibUtilities::Interp2D(), kNedges, 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::Geometry::m_xmap, and tinysimd::sqrt().

Referenced by v_GenGeomFactors().

◆ v_GenGeomFactors()

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

Set up GeoFac for this geometry using Coord quadrature distribution

Implements Nektar::SpatialDomains::Geometry.

Definition at line 254 of file QuadGeom.cpp.

255{
256 if (!m_setupState)
257 {
259 }
260
262 {
263 GeomType Gtype = eRegular;
264
266
267 // We will first check whether we have a regular or deformed
268 // geometry. We will define regular as those cases where the
269 // Jacobian and the metric terms of the derivative are constants
270 // (i.e. not coordinate dependent)
271
272 // Check to see if expansions are linear
273 // If not linear => deformed geometry
274 m_straightEdge = 1;
275 if ((m_xmap->GetBasisNumModes(0) != 2) ||
276 (m_xmap->GetBasisNumModes(1) != 2))
277 {
278 Gtype = eDeformed;
279 m_straightEdge = 0;
280 }
281
282 // For linear expansions, the mapping from standard to local
283 // element is given by the relation:
284 // x_i = 0.25 * [ ( x_i^A + x_i^B + x_i^C + x_i^D) +
285 // (-x_i^A + x_i^B + x_i^C - x_i^D)*xi_1 +
286 // (-x_i^A - x_i^B + x_i^C + x_i^D)*xi_2 +
287 // ( x_i^A - x_i^B + x_i^C - x_i^D)*xi_1*xi_2 ]
288 //
289 // The jacobian of the transformation and the metric terms
290 // dxi_i/dx_j, involve only terms of the form dx_i/dxi_j (both
291 // for coordim == 2 or 3). Inspecting the formula above, it can
292 // be appreciated that the derivatives dx_i/dxi_j will be
293 // constant, if the coefficient of the non-linear term is zero.
294 //
295 // That is why for regular geometry, we require
296 //
297 // x_i^A - x_i^B + x_i^C - x_i^D = 0
298 //
299 // or equivalently
300 //
301 // x_i^A - x_i^B = x_i^D - x_i^C
302 //
303 // This corresponds to quadrilaterals which are paralellograms.
304 m_manifold = Array<OneD, int>(m_coordim);
305 m_manifold[0] = 0;
306 m_manifold[1] = 1;
307 if (m_coordim == 3)
308 {
309 PointGeom e01, e21, norm;
310 e01.Sub(*m_verts[0], *m_verts[1]);
311 e21.Sub(*m_verts[3], *m_verts[1]);
312 norm.Mult(e01, e21);
313 int tmpi = 0;
314 double tmp = std::fabs(norm[0]);
315 if (tmp < fabs(norm[1]))
316 {
317 tmp = fabs(norm[1]);
318 tmpi = 1;
319 }
320 if (tmp < fabs(norm[2]))
321 {
322 tmpi = 2;
323 }
324 m_manifold[0] = (tmpi + 1) % 3;
325 m_manifold[1] = (tmpi + 2) % 3;
326 m_manifold[2] = (tmpi + 3) % 3;
327 }
328
329 if (Gtype == eRegular)
330 {
331 Array<OneD, Array<OneD, NekDouble>> verts(m_verts.size());
332 for (int i = 0; i < m_verts.size(); ++i)
333 {
334 verts[i] = Array<OneD, NekDouble>(3);
335 m_verts[i]->GetCoords(verts[i]);
336 }
337 // a00 + a01 xi1 + a02 xi2 + a03 xi1 xi2
338 // a10 + a11 xi1 + a12 xi2 + a03 xi1 xi2
339 m_isoParameter = Array<OneD, Array<OneD, NekDouble>>(2);
340 for (int i = 0; i < 2; i++)
341 {
342 unsigned int d = m_manifold[i];
343 m_isoParameter[i] = Array<OneD, NekDouble>(4, 0.);
344 // Karniadakis, Sherwin 2005, Appendix D
345 NekDouble A = verts[0][d];
346 NekDouble B = verts[1][d];
347 NekDouble D = verts[2][d];
348 NekDouble C = verts[3][d];
349 m_isoParameter[i][0] = 0.25 * (A + B + C + D); // 1
350 m_isoParameter[i][1] = 0.25 * (-A + B - C + D); // xi1
351 m_isoParameter[i][2] = 0.25 * (-A - B + C + D); // xi2
352 m_isoParameter[i][3] = 0.25 * (A - B - C + D); // xi1*xi2
353 NekDouble tmp =
354 fabs(m_isoParameter[i][1]) + fabs(m_isoParameter[i][2]);
355 if (fabs(m_isoParameter[i][3]) >
357 {
358 Gtype = eDeformed;
359 }
360 }
361 }
362
363 if (Gtype == eRegular)
364 {
366 }
367 else if (m_straightEdge)
368 {
370 }
371
373 Gtype, m_coordim, m_xmap, m_coeffs);
375 }
376}
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
GeomState m_geomFactorsState
State of the geometric factors.
Definition: Geometry.h:192
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
Definition: Geometry.h:190
static const NekDouble kNekZeroTol
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::NekConstants::kNekZeroTol, 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(), PreSolveStraightEdge(), Nektar::SpatialDomains::PointGeom::Sub(), Nektar::SpatialDomains::Geometry2D::v_CalculateInverseIsoParam(), v_FillGeom(), and v_Setup().

◆ v_GetCoord()

NekDouble Nektar::SpatialDomains::QuadGeom::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 118 of file QuadGeom.cpp.

120{
121 ASSERTL1(m_state == ePtsFilled, "Geometry is not in physical space");
122
123 Array<OneD, NekDouble> tmp(m_xmap->GetTotPoints());
124 m_xmap->BwdTrans(m_coeffs[i], tmp);
125
126 return m_xmap->PhysEvaluate(Lcoord, tmp);
127}
#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::QuadGeom::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 507 of file QuadGeom.cpp.

508{
509 return i % 2;
510}

◆ v_Reset()

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

513{
514 Geometry::v_Reset(curvedEdges, curvedFaces);
515 CurveMap::iterator it = curvedFaces.find(m_globalID);
516
517 if (it != curvedFaces.end())
518 {
519 m_curve = it->second;
520 }
521
522 for (int i = 0; i < 4; ++i)
523 {
524 m_edges[i]->Reset(curvedEdges, curvedFaces);
525 }
526
527 SetUpXmap();
528 SetUpCoeffs(m_xmap->GetNcoeffs());
529}
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::QuadGeom::v_Setup ( )
overrideprotectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 531 of file QuadGeom.cpp.

532{
533 if (!m_setupState)
534 {
535 for (int i = 0; i < 4; ++i)
536 {
537 m_edges[i]->Setup();
538 }
539 SetUpXmap();
540 SetUpCoeffs(m_xmap->GetNcoeffs());
541 m_setupState = true;
542 }
543}

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::QuadGeom::kNedges = 4
static

◆ kNverts

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

◆ XMLElementType

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

Definition at line 76 of file QuadGeom.h.