#include <TriGeom.h>

Inheritance diagram for Nektar::SpatialDomains::TriGeom:

Public Member Functions
	TriGeom ()

	TriGeom (const TriGeom &in)

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

	~TriGeom ()

Public Member Functions inherited from Nektar::SpatialDomains::Geometry2D
	Geometry2D ()

	Geometry2D (const int coordim, CurveSharedPtr curve)

virtual	~Geometry2D ()

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 destructor. More...

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...

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 &resid)
	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...

bool	MinMaxCheck (const Array< OneD, const NekDouble > &gloCoord)
	Check if given global coord is within twice the min/max distance of the element. More...

void	ClampLocCoords (Array< OneD, NekDouble > &locCoord, NekDouble tol)
	Clamp local coords to be within standard regions [-1, 1]^dim. More...

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...

void	Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
	Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated GeomFactors. More...

void	Setup ()

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
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 void	v_GenGeomFactors ()

virtual void	v_FillGeom ()

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 bool	v_ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &resid)

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 ()

Protected Member Functions inherited from Nektar::SpatialDomains::Geometry2D
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 &resid)

Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
void	GenGeomFactors ()
	Handles generation of geometry factors. More...

virtual Geometry2DSharedPtr	v_GetFace (int i) const
	Returns face `i` of this object. 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_GetNumFaces () const
	Get the number of faces of this object. 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 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...

void	SetUpCoeffs (const int nCoeffs)
	Initialise the Geometry::m_coeffs array. More...

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::Orientation >	m_eorient

CurveSharedPtr	m_curve

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...

Static Protected Attributes inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsVector	m_regGeomFactorsManager

Detailed Description

Definition at line 61 of file TriGeom.h.

Constructor & Destructor Documentation

◆ TriGeom() [1/3]

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

Definition at line 51 of file TriGeom.cpp.

References Nektar::LibUtilities::eTriangle.

 {
     m_shapeType = LibUtilities::eTriangle;
 }

◆ TriGeom() [2/3]

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

Definition at line 85 of file TriGeom.cpp.

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.

     : Geometry2D(in)
 {
     // From Geometry
     m_shapeType = in.m_shapeType;
 
     // From TriFaceComponent
     m_globalID = in.m_globalID;
 
     // From TriGeom
     m_verts = in.m_verts;
     m_edges = in.m_edges;
     for (int i = 0; i < kNedges; i++)
     {
         m_eorient[i] = in.m_eorient[i];
     }
 }

◆ TriGeom() [3/3]

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

Definition at line 56 of file TriGeom.cpp.

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.

     : Geometry2D(edges[0]->GetVertex(0)->GetCoordim(), curve)
 {
     int j;
 
     m_shapeType = LibUtilities::eTriangle;
     m_globalID = id;
 
     // Copy the edge shared pointers.
     m_edges.insert(m_edges.begin(), edges, edges + TriGeom::kNedges);
     m_eorient.resize(kNedges);
 
     for (j = 0; j < kNedges; ++j)
     {
         m_eorient[j] =
             SegGeom::GetEdgeOrientation(*edges[j], *edges[(j + 1) % kNedges]);
         m_verts.push_back(
             edges[j]->GetVertex(m_eorient[j] == StdRegions::eForwards ? 0 : 1));
     }
 
     m_eorient[2] = m_eorient[2] == StdRegions::eBackwards ?
         StdRegions::eForwards : StdRegions::eBackwards;
 
     m_coordim = edges[0]->GetVertex(0)->GetCoordim();
     ASSERTL0(m_coordim > 1, "Cannot call function with dim == 1");
 }

◆ ~TriGeom()

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

Definition at line 103 of file TriGeom.cpp.

104 {

105 }

Member Function Documentation

◆ GetFaceOrientation() [1/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 118 of file TriGeom.cpp.

References Nektar::SpatialDomains::Geometry2D::m_verts.

Referenced by Nektar::MultiRegions::DisContField3D::FindPeriodicFaces().

 {
     return GetFaceOrientation(face1.m_verts, face2.m_verts,
                               doRot, dir, angle, tol);
 }

◆ GetFaceOrientation() [2/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 126 of file TriGeom.cpp.

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, and Nektar::SpatialDomains::PointGeom::Rotate().

 {
     int i, j, vmap[3] = {-1, -1, -1};
 
     if(doRot)
     {
         PointGeom rotPt;
         
         for (i = 0; i < 3; ++i)
         {
             rotPt.Rotate((*face1[i]), dir, angle);
             for (j = 0; j < 3; ++j)
             {
                 if (rotPt.dist(*face2[j]) < tol)
                 {
                     vmap[j] = i;
                     break;
                 }
             }
         }
     }
     else
     {
         
         NekDouble x, y, z, x1, y1, z1, cx = 0.0, cy = 0.0, cz = 0.0;
 
         // For periodic faces, we calculate the vector between the centre
         // points of the two faces. (For connected faces this will be
         // zero). We can then use this to determine alignment later in the
         // algorithm.
         for (i = 0; i < 3; ++i)
         {
             cx += (*face2[i])(0) - (*face1[i])(0);
             cy += (*face2[i])(1) - (*face1[i])(1);
             cz += (*face2[i])(2) - (*face1[i])(2);
         }
         cx /= 3;
         cy /= 3;
         cz /= 3;
         
         // Now construct a mapping which takes us from the vertices of one
         // face to the other. That is, vertex j of face2 corresponds to
         // vertex vmap[j] of face1.
         for (i = 0; i < 3; ++i)
         {
             x = (*face1[i])(0);
             y = (*face1[i])(1);
             z = (*face1[i])(2);
             for (j = 0; j < 3; ++j)
             {
                 x1 = (*face2[j])(0) - cx;
                 y1 = (*face2[j])(1) - cy;
                 z1 = (*face2[j])(2) - cz;
                 if (sqrt((x1 - x) * (x1 - x) + (y1 - y) * (y1 - y) +
                          (z1 - z) * (z1 - z)) < 1e-8)
                 {
                     vmap[j] = i;
                 break;
                 }
             }
         }
     }
     
     if (vmap[1] == (vmap[0] + 1) % 3)
     {
         switch (vmap[0])
         {
         case 0:
             return StdRegions::eDir1FwdDir1_Dir2FwdDir2;
             break;
         case 1:
             return StdRegions::eDir1FwdDir2_Dir2BwdDir1;
             break;
         case 2:
             return StdRegions::eDir1BwdDir1_Dir2BwdDir2;
             break;
         }
     }
     else
     {
         switch (vmap[0])
         {
         case 0:
             return StdRegions::eDir1FwdDir2_Dir2FwdDir1;
             break;
         case 1:
             return StdRegions::eDir1BwdDir1_Dir2FwdDir2;
             break;
         case 2:
             return StdRegions::eDir1BwdDir2_Dir2BwdDir1;
             break;
         }
     }
     
     ASSERTL0(false, "Unable to determine triangle orientation");
     return StdRegions::eNoOrientation;
 }

◆ SetUpXmap()

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

private

Definition at line 608 of file TriGeom.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eGaussRadauMAlpha1Beta0, 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().

 {
     int order0 = m_edges[0]->GetXmap()->GetBasis(0)->GetNumModes();
     int order1 = max(order0,
                      max(m_edges[1]->GetXmap()->GetBasis(0)->GetNumModes(),
                          m_edges[2]->GetXmap()->GetBasis(0)->GetNumModes()));
 
     const LibUtilities::BasisKey B0(
         LibUtilities::eModified_A,
         order0,
         LibUtilities::PointsKey(order0+1, LibUtilities::eGaussLobattoLegendre));
     const LibUtilities::BasisKey B1(
         LibUtilities::eModified_B,
         order1,
         LibUtilities::PointsKey(order1,
                                 LibUtilities::eGaussRadauMAlpha1Beta0));
 
     m_xmap = MemoryManager<StdRegions::StdTriExp>::AllocateSharedPtr(B0, B1);
 }

◆ v_ContainsPoint()

bool Nektar::SpatialDomains::TriGeom::v_ContainsPoint	(	const Array< OneD, const NekDouble > &	gloCoord,
		Array< OneD, NekDouble > &	locCoord,
		NekDouble	tol,
		NekDouble &	resid
	)

protectedvirtual

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

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

Parameters

gloCoord	The coordinate \( (x,y,z) \).
locCoord	On exit, this is the local coordinate \(\vec{\xi}\) such that \(\chi(\vec{\xi}) = \vec{x}\).
tol	The tolerance used to dictate the bounding box of the standard coordinates \(\vec{\xi}\).
resid	On exit, returns the minimum distance between `gloCoord` and the quadrature points inside the element.

Returns: true if the coordinate gloCoord is contained in the element; false otherwise.

See also: Geometry::GetLocCoords.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 546 of file TriGeom.cpp.

References Nektar::SpatialDomains::Geometry::ClampLocCoords(), Nektar::SpatialDomains::eRegular, Nektar::SpatialDomains::Geometry::GetLocCoords(), Nektar::SpatialDomains::Geometry::GetMetricInfo(), and Nektar::SpatialDomains::Geometry::MinMaxCheck().

 {
     //Rough check if within twice min/max point
     if (GetMetricInfo()->GetGtype() != eRegular)
     {
         if (!MinMaxCheck(gloCoord))
         {
             return false;
         }
     }
 
     // Convert to the local (eta) coordinates.
     resid = GetLocCoords(gloCoord, stdCoord);
 
     if (stdCoord[0] >= -(1 + tol) && stdCoord[1] >= -(1 + tol) &&
         stdCoord[0] + stdCoord[1] <= tol)
     {
         return true;
     }
 
     //Clamp local coords
     ClampLocCoords(stdCoord, tol);
 
     return false;
 }

◆ v_FillGeom()

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

protectedvirtual

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 261 of file TriGeom.cpp.

Referenced by v_GenGeomFactors(), and v_GetLocCoords().

 {
     // check to see if geometry structure is already filled
     if (m_state == ePtsFilled)
     {
         return;
     }
 
     int i, j, k;
     int nEdgeCoeffs = m_xmap->GetEdgeNcoeffs(0);
 
     if (m_curve)
     {
         int pdim = LibUtilities::PointsManager()[LibUtilities::PointsKey(
                                                      2, m_curve->m_ptype)]
                        ->GetPointsDim();
 
         // Deal with 2D points type separately
         // (e.g. electrostatic or Fekete points) to 1D tensor
         // product.
         if (pdim == 2)
         {
             int N = m_curve->m_points.size();
             int nEdgePts =
                 (-1 + (int)sqrt(static_cast<NekDouble>(8 * N + 1))) / 2;
 
             ASSERTL0(nEdgePts * (nEdgePts + 1) / 2 == N,
                      "NUMPOINTS should be a triangle number for"
                      " triangle curved face " +
                          boost::lexical_cast<string>(m_globalID));
 
             // Sanity check 1: are curved vertices consistent with
             // triangle vertices?
             for (i = 0; i < 3; ++i)
             {
                 NekDouble dist = m_verts[i]->dist(*(m_curve->m_points[i]));
                 if (dist > NekConstants::kVertexTheSameDouble)
                 {
                     std::stringstream ss;
                     ss << "Curved vertex " << i << " of triangle " << m_globalID
                        << " is separated from expansion vertex by"
                        << " more than " << NekConstants::kVertexTheSameDouble
                        << " (dist = " << dist << ")";
                     NEKERROR(ErrorUtil::ewarning, ss.str().c_str());
                 }
             }
 
             // Sanity check 2: are curved edges from the face curvature
             // consistent with curved edges?
             for (i = 0; i < kNedges; ++i)
             {
                 CurveSharedPtr edgeCurve = m_edges[i]->GetCurve();
 
                 ASSERTL0(edgeCurve->m_points.size() == nEdgePts,
                          "Number of edge points does not correspond "
                          "to number of face points in triangle " +
                              boost::lexical_cast<string>(m_globalID));
 
                 const int offset = 3 + i * (nEdgePts - 2);
                 NekDouble maxDist = 0.0;
 
                 // Account for different ordering of nodal coordinates
                 // vs. Cartesian ordering of element.
                 StdRegions::Orientation orient = m_eorient[i];
 
                 if (i == 2)
                 {
                     orient = orient == StdRegions::eForwards ?
                         StdRegions::eBackwards : StdRegions::eForwards;
                 }
 
                 if (orient == StdRegions::eForwards)
                 {
                     for (j = 0; j < nEdgePts - 2; ++j)
                     {
                         NekDouble dist = m_curve->m_points[offset + j]->dist(
                             *(edgeCurve->m_points[j + 1]));
                         maxDist = dist > maxDist ? dist : maxDist;
                     }
                 }
                 else
                 {
                     for (j = 0; j < nEdgePts - 2; ++j)
                     {
                         NekDouble dist = m_curve->m_points[offset + j]->dist(
                             *(edgeCurve->m_points[nEdgePts - 2 - j]));
                         maxDist = dist > maxDist ? dist : maxDist;
                     }
                 }
 
                 if (maxDist > NekConstants::kVertexTheSameDouble)
                 {
                     std::stringstream ss;
                     ss << "Curved edge " << i << " of triangle " << m_globalID
                        << " has a point separated from edge interior"
                        << " points by more than "
                        << NekConstants::kVertexTheSameDouble
                        << " (maxdist = " << maxDist << ")";
                     NEKERROR(ErrorUtil::ewarning, ss.str().c_str());
                 }
             }
 
             const LibUtilities::PointsKey P0(
                 nEdgePts, LibUtilities::eGaussLobattoLegendre);
             const LibUtilities::PointsKey P1(
                 nEdgePts, LibUtilities::eGaussRadauMAlpha1Beta0);
             const LibUtilities::BasisKey T0(
                 LibUtilities::eOrtho_A, nEdgePts, P0);
             const LibUtilities::BasisKey T1(
                 LibUtilities::eOrtho_B, nEdgePts, P1);
             Array<OneD, NekDouble> phys(
                 max(nEdgePts * nEdgePts, m_xmap->GetTotPoints()));
             Array<OneD, NekDouble> tmp(nEdgePts * nEdgePts);
 
             for (i = 0; i < m_coordim; ++i)
             {
                 // Create a StdNodalTriExp.
                 StdRegions::StdNodalTriExpSharedPtr t =
                     MemoryManager<StdRegions::StdNodalTriExp>::
                         AllocateSharedPtr(T0, T1, m_curve->m_ptype);
 
                 for (j = 0; j < N; ++j)
                 {
                     phys[j] = (m_curve->m_points[j]->GetPtr())[i];
                 }
 
                 t->BwdTrans(phys, tmp);
 
                 // Interpolate points to standard region.
                 LibUtilities::Interp2D(P0,
                                        P1,
                                        tmp,
                                        m_xmap->GetBasis(0)->GetPointsKey(),
                                        m_xmap->GetBasis(1)->GetPointsKey(),
                                        phys);
 
                 // Forwards transform to get coefficient space.
                 m_xmap->FwdTrans(phys, m_coeffs[i]);
             }
         }
         else if (pdim == 1)
         {
             int npts = m_curve->m_points.size();
             int nEdgePts = (int)sqrt(static_cast<NekDouble>(npts));
             Array<OneD, NekDouble> tmp(npts);
             Array<OneD, NekDouble> phys(m_xmap->GetTotPoints());
             LibUtilities::PointsKey curveKey(nEdgePts, m_curve->m_ptype);
 
             // Sanity checks:
             // - Curved faces should have square number of points;
             // - Each edge should have sqrt(npts) points.
             ASSERTL0(nEdgePts * nEdgePts == npts,
                      "NUMPOINTS should be a square number for"
                      " triangle " +
                          boost::lexical_cast<string>(m_globalID));
 
             for (i = 0; i < kNedges; ++i)
             {
                 ASSERTL0(m_edges[i]->GetXmap()->GetNcoeffs() == nEdgePts,
                          "Number of edge points does not correspond to "
                          "number of face points in triangle " +
                              boost::lexical_cast<string>(m_globalID));
             }
 
             for (i = 0; i < m_coordim; ++i)
             {
                 for (j = 0; j < npts; ++j)
                 {
                     tmp[j] = (m_curve->m_points[j]->GetPtr())[i];
                 }
 
                 // Interpolate curve points to standard triangle
                 // points.
                 LibUtilities::Interp2D(curveKey,
                                        curveKey,
                                        tmp,
                                        m_xmap->GetBasis(0)->GetPointsKey(),
                                        m_xmap->GetBasis(1)->GetPointsKey(),
                                        phys);
 
                 // Forwards transform to get coefficient space.
                 m_xmap->FwdTrans(phys, m_coeffs[i]);
             }
         }
         else
         {
             ASSERTL0(false,
                      "Only 1D/2D points distributions "
                      "supported.");
         }
     }
 
     Array<OneD, unsigned int> mapArray(nEdgeCoeffs);
     Array<OneD, int> signArray(nEdgeCoeffs);
 
     for (i = 0; i < kNedges; i++)
     {
         m_edges[i]->FillGeom();
         m_xmap->GetEdgeToElementMap(i, m_eorient[i], mapArray, signArray);
 
         nEdgeCoeffs = m_edges[i]->GetXmap()->GetNcoeffs();
 
         for (j = 0; j < m_coordim; j++)
         {
             for (k = 0; k < nEdgeCoeffs; k++)
             {
                 m_coeffs[j][mapArray[k]] =
                     signArray[k] * m_edges[i]->GetCoeffs(j)[k];
             }
         }
     }
 
     m_state = ePtsFilled;
 }

◆ v_GenGeomFactors()

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

protectedvirtual

Implements Nektar::SpatialDomains::Geometry.

Definition at line 226 of file TriGeom.cpp.

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_setupState, Nektar::SpatialDomains::Geometry::m_xmap, v_FillGeom(), and v_Setup().

 {
     if(!m_setupState)
     {
         TriGeom::v_Setup();
     }
 
     if (m_geomFactorsState != ePtsFilled)
     {
         GeomType Gtype = eRegular;
 
         TriGeom::v_FillGeom();
 
         // check to see if expansions are linear
         for (int i = 0; i < m_coordim; ++i)
         {
             if (m_xmap->GetBasisNumModes(0) != 2 ||
                 m_xmap->GetBasisNumModes(1) != 2)
             {
                 Gtype = eDeformed;
             }
         }
 
         m_geomFactors = MemoryManager<GeomFactors>::AllocateSharedPtr(
             Gtype, m_coordim, m_xmap, m_coeffs);
 
         m_geomFactorsState = ePtsFilled;
     }
 }

◆ v_GetCoord()

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

protectedvirtual

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 107 of file TriGeom.cpp.

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

 {
     ASSERTL1(m_state == ePtsFilled, "Geometry is not in physical space");
 
     Array<OneD, NekDouble> tmp(m_xmap->GetTotPoints());
     m_xmap->BwdTrans(m_coeffs[i], tmp);
 
     return m_xmap->PhysEvaluate(Lcoord, tmp);
 }

◆ v_GetLocCoords()

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

protectedvirtual

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

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

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

Parameters

coords	Input Cartesian global coordinates
Lcoords	Corresponding local coordinates

Returns: Distance between obtained coordinates and provided ones.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 479 of file TriGeom.cpp.

 {
     NekDouble resid = 0.0;
     TriGeom::v_FillGeom();
 
     // calculate local coordinate for coord
     if (GetMetricInfo()->GetGtype() == eRegular)
     {
         NekDouble coords2 = (m_coordim == 3) ? coords[2] : 0.0;
         PointGeom dv1, dv2, norm, orth1, orth2;
         PointGeom xin(m_coordim, 0, coords[0], coords[1], coords2);
 
         // Calculate edge vectors from 0-1 and 0-2 edges.
         dv1.Sub(*m_verts[1], *m_verts[0]);
         dv2.Sub(*m_verts[2], *m_verts[0]);
 
         // Obtain normal to plane in which dv1 and dv2 lie
         norm.Mult(dv1, dv2);
 
         // Obtain vector which are proportional to normal of dv1 and dv2.
         orth1.Mult(norm, dv1);
         orth2.Mult(norm, dv2);
 
         // Start with vector of desired points minus vertex_0
         xin -= *m_verts[0];
 
         // Calculate length using L/|dv1| = (x-v0).n1/(dv1.n1) for coordiante 1
         // Then rescale to [-1,1].
         Lcoords[0] = xin.dot(orth2) / dv1.dot(orth2);
         Lcoords[0] = 2 * Lcoords[0] - 1;
         Lcoords[1] = xin.dot(orth1) / dv2.dot(orth1);
         Lcoords[1] = 2 * Lcoords[1] - 1;
     }
     else
     {
         // Determine nearest point of coords  to values in m_xmap
         int npts = m_xmap->GetTotPoints();
         Array<OneD, NekDouble> ptsx(npts), ptsy(npts);
         Array<OneD, NekDouble> tmpx(npts), tmpy(npts);
 
         m_xmap->BwdTrans(m_coeffs[0], ptsx);
         m_xmap->BwdTrans(m_coeffs[1], ptsy);
 
         const Array<OneD, const NekDouble> za = m_xmap->GetPoints(0);
         const Array<OneD, const NekDouble> zb = m_xmap->GetPoints(1);
 
         // guess the first local coords based on nearest point
         Vmath::Sadd(npts, -coords[0], ptsx, 1, tmpx, 1);
         Vmath::Sadd(npts, -coords[1], ptsy, 1, tmpy, 1);
         Vmath::Vmul(npts, tmpx, 1, tmpx, 1, tmpx, 1);
         Vmath::Vvtvp(npts, tmpy, 1, tmpy, 1, tmpx, 1, tmpx, 1);
 
         int min_i = Vmath::Imin(npts, tmpx, 1);
 
         Lcoords[0] = za[min_i % za.num_elements()];
         Lcoords[1] = zb[min_i / za.num_elements()];
 
         // recover cartesian coordinate from collapsed coordinate.
         Lcoords[0] = (1.0 + Lcoords[0]) * (1.0 - Lcoords[1]) / 2 - 1.0;
 
         // Perform newton iteration to find local coordinates
         NewtonIterationForLocCoord(coords, ptsx, ptsy, Lcoords, resid);
     }
     return resid;
 }

◆ v_Reset()

void Nektar::SpatialDomains::TriGeom::v_Reset	(	CurveMap &	curvedEdges,
		CurveMap &	curvedFaces
	)

protectedvirtual

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 575 of file TriGeom.cpp.

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().

 {
     Geometry::v_Reset(curvedEdges, curvedFaces);
     CurveMap::iterator it = curvedFaces.find(m_globalID);
 
     if (it != curvedFaces.end())
     {
         m_curve = it->second;
     }
 
     for (int i = 0; i < 3; ++i)
     {
         m_edges[i]->Reset(curvedEdges, curvedFaces);
     }
 
     SetUpXmap();
     SetUpCoeffs(m_xmap->GetNcoeffs());
 }

◆ v_Setup()

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

protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 594 of file TriGeom.cpp.

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().

 {
     if(!m_setupState)
     {
         for (int i = 0; i < 3; ++i)
         {
             m_edges[i]->Setup();
         }
         SetUpXmap();
         SetUpCoeffs(m_xmap->GetNcoeffs());
         m_setupState = true;
     }
 }

Member Data Documentation

◆ kNedges

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

static

Get the orientation of face1.

Definition at line 73 of file TriGeom.h.

Referenced by TriGeom(), and v_FillGeom().

◆ kNverts

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

static

Definition at line 74 of file TriGeom.h.

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

Public Member Functions

Static Public Member Functions

Static Public Attributes

Protected Member Functions

Private Member Functions

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ TriGeom() [1/3]

◆ TriGeom() [2/3]

◆ TriGeom() [3/3]

◆ ~TriGeom()

Member Function Documentation

◆ GetFaceOrientation() [1/2]

◆ GetFaceOrientation() [2/2]

◆ SetUpXmap()

◆ v_ContainsPoint()

◆ v_FillGeom()

◆ v_GenGeomFactors()

◆ v_GetCoord()

◆ v_GetLocCoords()

◆ v_Reset()

◆ v_Setup()

Member Data Documentation

◆ kNedges

◆ kNverts