Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Static Public Attributes | Protected Member Functions | Private Member Functions | Static Private Attributes | List of all members
Nektar::SpatialDomains::HexGeom Class Reference

#include <HexGeom.h>

Inheritance diagram for Nektar::SpatialDomains::HexGeom:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SpatialDomains::HexGeom:
Collaboration graph
[legend]

Public Member Functions

 HexGeom ()
 HexGeom (const QuadGeomSharedPtr faces[])
 ~HexGeom ()
- Public Member Functions inherited from Nektar::SpatialDomains::Geometry3D
 Geometry3D ()
 Geometry3D (const int coordim)
virtual ~Geometry3D ()
int GetEid (int i) const
 Return the ID of edge i in this element.
const Geometry1DSharedPtr GetEdge (int i) const
Geometry2DSharedPtr GetFace (int i)
 Return face i in this element.
StdRegions::Orientation GetFaceOrient (const int i) const
 Return the orientation of face i in this element.
int GetDir (const int faceidx, const int facedir) const
 Returns the element coordinate direction corresponding to a given face coordinate direction.
- Public Member Functions inherited from Nektar::SpatialDomains::Geometry
 Geometry ()
 Geometry (int coordim)
virtual ~Geometry ()
bool IsElmtConnected (int gvo_id, int locid) const
void AddElmtConnected (int gvo_id, int locid)
int NumElmtConnected () const
int GetCoordim () const
void SetCoordim (int coordim)
GeomFactorsSharedPtr GetGeomFactors ()
GeomFactorsSharedPtr GetRefGeomFactors (const Array< OneD, const LibUtilities::BasisSharedPtr > &tbasis)
GeomFactorsSharedPtr GetMetricInfo ()
LibUtilities::ShapeType GetShapeType (void)
int GetGlobalID (void)
void SetGlobalID (int globalid)
int GetVid (int i) const
int GetFid (int i) const
int GetNumVerts () const
PointGeomSharedPtr GetVertex (int i) const
StdRegions::Orientation GetEorient (const int i) const
StdRegions::Orientation GetPorient (const int i) const
int GetNumEdges () const
int GetNumFaces () const
int GetShapeDim () const
StdRegions::StdExpansionSharedPtr GetXmap () const
const Array< OneD, const
NekDouble > & 		GetCoeffs (const int i) const
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol)
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &resid)
int GetVertexEdgeMap (int i, int j) const
int GetVertexFaceMap (int i, int j) const
 return the id of the $j^{th}$ face attached to the $ i^{th}$ vertex
int GetEdgeFaceMap (int i, int j) const
void FillGeom ()
 Put all quadrature information into face/edge structure and backward transform.
NekDouble GetLocCoords (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
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.
void SetOwnData ()
const LibUtilities::BasisSharedPtr GetBasis (const int i)
 Return the j-th basis of the i-th co-ordinate dimension.
const LibUtilities::PointsKeyVector GetPointsKeys ()

Static Public Attributes

static const int kNverts = 8
static const int kNedges = 12
static const int kNqfaces = 6
static const int kNtfaces = 0
static const int kNfaces = kNqfaces + kNtfaces
static const std::string XMLElementType

Protected Member Functions

virtual void v_GenGeomFactors ()
virtual NekDouble v_GetLocCoords (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
virtual bool v_ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
 Determines if a point specified in global coordinates is located within this tetrahedral geometry.
virtual bool v_ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol)
virtual bool v_ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &resid)
virtual int v_GetNumVerts () const
virtual int v_GetNumEdges () const
virtual int v_GetNumFaces () const
virtual int v_GetVertexEdgeMap (const int i, const int j) const
virtual int v_GetVertexFaceMap (const int i, const int j) const
virtual int v_GetEdgeFaceMap (const int i, const int j) const
virtual int v_GetDir (const int faceidx, const int facedir) const
- Protected Member Functions inherited from Nektar::SpatialDomains::Geometry3D
void NewtonIterationForLocCoord (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &ptsx, const Array< OneD, const NekDouble > &ptsy, const Array< OneD, const NekDouble > &ptsz, Array< OneD, NekDouble > &Lcoords, NekDouble &resid)
virtual void v_FillGeom ()
 Put all quadrature information into face/edge structure and backward transform.
virtual NekDouble v_GetCoord (const int i, const Array< OneD, const NekDouble > &Lcoord)
 Given local collapsed coordinate Lcoord return the value of physical coordinate in direction i.
virtual int v_GetShapeDim () const
 Return the dimension of this element.
virtual int v_GetVid (int i) const
 Return the vertex ID of vertex i.
virtual PointGeomSharedPtr v_GetVertex (int i) const
 Return vertex i in this element.
virtual const SegGeomSharedPtr v_GetEdge (int i) const
 Return edge i of this element.
virtual StdRegions::Orientation v_GetEorient (const int i) const
 Return the orientation of edge i in this element.
virtual int v_GetEid (int i) const
 Return the ID of edge i in this element.
virtual const Geometry2DSharedPtr v_GetFace (int i) const
 Return face i in this element.
virtual StdRegions::Orientation v_GetFaceOrient (const int i) const
 Return the orientation of face i in this element.
virtual int v_GetFid (int i) const
 Return the ID of face i in this element.
virtual int v_GetEid () const
 Return the ID of this element.
virtual int v_WhichEdge (SegGeomSharedPtr edge)
 Return the local ID of a given edge.
virtual int v_WhichFace (Geometry2DSharedPtr face)
 Return the local ID of a given face.
virtual const
LibUtilities::BasisSharedPtr 		v_GetBasis (const int i)
 Return the j-th basis of the i-th co-ordinate dimension.
virtual void v_AddElmtConnected (int gvo_id, int locid)
virtual bool v_IsElmtConnected (int gvo_id, int locid) const
virtual int v_NumElmtConnected () const
virtual void v_SetOwnData ()
- Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
void GenGeomFactors ()
virtual StdRegions::Orientation v_GetPorient (const int i) const
virtual
StdRegions::StdExpansionSharedPtr 		v_GetXmap () const
virtual int v_GetCoordim () const
void SetUpCoeffs (const int nCoeffs)

Private Member Functions

void SetUpLocalEdges ()
void SetUpLocalVertices ()
void SetUpEdgeOrientation ()
void SetUpFaceOrientation ()

Static Private Attributes

static const unsigned int VertexEdgeConnectivity [8][3]
static const unsigned int VertexFaceConnectivity [8][3]
static const unsigned int EdgeFaceConnectivity [12][2]

Additional Inherited Members

- Static Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsSharedPtr ValidateRegGeomFactor (GeomFactorsSharedPtr geomFactor)
- Protected Attributes inherited from Nektar::SpatialDomains::Geometry3D
PointGeomVector m_verts
SegGeomVector m_edges
Geometry2DVector m_faces
std::vector
< StdRegions::Orientation > 		m_eorient
std::vector
< StdRegions::Orientation > 		m_forient
std::list< CompToElmt > m_elmtmap
bool m_owndata
int m_eid
bool m_ownverts
- Static Protected Attributes inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsVector m_regGeomFactorsManager

Detailed Description

Definition at line 47 of file HexGeom.h.

Constructor & Destructor Documentation

Nektar::SpatialDomains::HexGeom::HexGeom ( )

Definition at line 56 of file HexGeom.cpp.

References Nektar::LibUtilities::eHexahedron, and Nektar::SpatialDomains::Geometry::m_shapeType.

{
m_shapeType = LibUtilities::eHexahedron;
}
Nektar::SpatialDomains::HexGeom::HexGeom ( const QuadGeomSharedPtr  faces[])

Copy the face shared pointers

Set up orientation vectors with correct amount of elements.

Determine necessary order for standard region. This can almost certainly be simplified but works for now!

Definition at line 61 of file HexGeom.cpp.

References Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eHexahedron, Nektar::LibUtilities::eModified_A, Nektar::SpatialDomains::Geometry::GetXmap(), kNedges, kNfaces, Nektar::SpatialDomains::Geometry3D::m_eorient, Nektar::SpatialDomains::Geometry3D::m_faces, Nektar::SpatialDomains::Geometry3D::m_forient, Nektar::SpatialDomains::Geometry::m_shapeType, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry::SetUpCoeffs(), SetUpEdgeOrientation(), SetUpFaceOrientation(), SetUpLocalEdges(), and SetUpLocalVertices().

:
Geometry3D(faces[0]->GetEdge(0)->GetVertex(0)->GetCoordim())
{
m_shapeType = LibUtilities::eHexahedron;
/// Copy the face shared pointers
m_faces.insert(m_faces.begin(), faces, faces+HexGeom::kNfaces);
/// Set up orientation vectors with correct amount of elements.
m_eorient.resize(kNedges);
m_forient.resize(kNfaces);
SetUpLocalEdges();
SetUpLocalVertices();
SetUpEdgeOrientation();
SetUpFaceOrientation();
/// Determine necessary order for standard region. This can almost
/// certainly be simplified but works for now!
vector<int> tmp1, tmp2;
if (m_forient[0] < 9)
{
tmp1.push_back(faces[0]->GetXmap()->GetEdgeNcoeffs (0));
tmp1.push_back(faces[0]->GetXmap()->GetEdgeNcoeffs (2));
tmp2.push_back(faces[0]->GetXmap()->GetEdgeNumPoints(0));
tmp2.push_back(faces[0]->GetXmap()->GetEdgeNumPoints(2));
}
else
{
tmp1.push_back(faces[0]->GetXmap()->GetEdgeNcoeffs (1));
tmp1.push_back(faces[0]->GetXmap()->GetEdgeNcoeffs (3));
tmp2.push_back(faces[0]->GetXmap()->GetEdgeNumPoints(1));
tmp2.push_back(faces[0]->GetXmap()->GetEdgeNumPoints(3));
}
if (m_forient[5] < 9)
{
tmp1.push_back(faces[5]->GetXmap()->GetEdgeNcoeffs (0));
tmp1.push_back(faces[5]->GetXmap()->GetEdgeNcoeffs (2));
tmp2.push_back(faces[5]->GetXmap()->GetEdgeNumPoints(0));
tmp2.push_back(faces[5]->GetXmap()->GetEdgeNumPoints(2));
}
else
{
tmp1.push_back(faces[5]->GetXmap()->GetEdgeNcoeffs (1));
tmp1.push_back(faces[5]->GetXmap()->GetEdgeNcoeffs (3));
tmp2.push_back(faces[5]->GetXmap()->GetEdgeNumPoints(1));
tmp2.push_back(faces[5]->GetXmap()->GetEdgeNumPoints(3));
}
int order0 = *max_element(tmp1.begin(), tmp1.end());
int points0 = *max_element(tmp2.begin(), tmp2.end());
tmp1.clear();
tmp2.clear();
if (m_forient[0] < 9)
{
tmp1.push_back(faces[0]->GetXmap()->GetEdgeNcoeffs (1));
tmp1.push_back(faces[0]->GetXmap()->GetEdgeNcoeffs (3));
tmp2.push_back(faces[0]->GetXmap()->GetEdgeNumPoints(1));
tmp2.push_back(faces[0]->GetXmap()->GetEdgeNumPoints(3));
}
else
{
tmp1.push_back(faces[0]->GetXmap()->GetEdgeNcoeffs (0));
tmp1.push_back(faces[0]->GetXmap()->GetEdgeNcoeffs (2));
tmp2.push_back(faces[0]->GetXmap()->GetEdgeNumPoints(0));
tmp2.push_back(faces[0]->GetXmap()->GetEdgeNumPoints(2));
}
if (m_forient[5] < 9)
{
tmp1.push_back(faces[5]->GetXmap()->GetEdgeNcoeffs (1));
tmp1.push_back(faces[5]->GetXmap()->GetEdgeNcoeffs (3));
tmp2.push_back(faces[5]->GetXmap()->GetEdgeNumPoints(1));
tmp2.push_back(faces[5]->GetXmap()->GetEdgeNumPoints(3));
}
else
{
tmp1.push_back(faces[5]->GetXmap()->GetEdgeNcoeffs (0));
tmp1.push_back(faces[5]->GetXmap()->GetEdgeNcoeffs (2));
tmp2.push_back(faces[5]->GetXmap()->GetEdgeNumPoints(0));
tmp2.push_back(faces[5]->GetXmap()->GetEdgeNumPoints(2));
}
int order1 = *max_element(tmp1.begin(), tmp1.end());
int points1 = *max_element(tmp2.begin(), tmp2.end());
tmp1.clear();
tmp2.clear();
if (m_forient[1] < 9)
{
tmp1.push_back(faces[1]->GetXmap()->GetEdgeNcoeffs (1));
tmp1.push_back(faces[1]->GetXmap()->GetEdgeNcoeffs (3));
tmp2.push_back(faces[1]->GetXmap()->GetEdgeNumPoints(1));
tmp2.push_back(faces[1]->GetXmap()->GetEdgeNumPoints(3));
}
else
{
tmp1.push_back(faces[1]->GetXmap()->GetEdgeNcoeffs (0));
tmp1.push_back(faces[1]->GetXmap()->GetEdgeNcoeffs (2));
tmp2.push_back(faces[1]->GetXmap()->GetEdgeNumPoints(0));
tmp2.push_back(faces[1]->GetXmap()->GetEdgeNumPoints(2));
}
if (m_forient[3] < 9)
{
tmp1.push_back(faces[3]->GetXmap()->GetEdgeNcoeffs (1));
tmp1.push_back(faces[3]->GetXmap()->GetEdgeNcoeffs (3));
tmp2.push_back(faces[3]->GetXmap()->GetEdgeNumPoints(1));
tmp2.push_back(faces[3]->GetXmap()->GetEdgeNumPoints(3));
}
else
{
tmp1.push_back(faces[3]->GetXmap()->GetEdgeNcoeffs (0));
tmp1.push_back(faces[3]->GetXmap()->GetEdgeNcoeffs (2));
tmp2.push_back(faces[3]->GetXmap()->GetEdgeNumPoints(0));
tmp2.push_back(faces[3]->GetXmap()->GetEdgeNumPoints(2));
}
int order2 = *max_element(tmp1.begin(), tmp1.end());
int points2 = *max_element(tmp2.begin(), tmp2.end());
const LibUtilities::BasisKey A(
LibUtilities::eModified_A, order0,
LibUtilities::PointsKey(points0,LibUtilities::eGaussLobattoLegendre));
const LibUtilities::BasisKey B(
LibUtilities::eModified_A, order1,
LibUtilities::PointsKey(points1,LibUtilities::eGaussLobattoLegendre));
const LibUtilities::BasisKey C(
LibUtilities::eModified_A, order2,
LibUtilities::PointsKey(points2,LibUtilities::eGaussLobattoLegendre));
m_xmap = MemoryManager<StdRegions::StdHexExp>::AllocateSharedPtr(A,B,C);
SetUpCoeffs(m_xmap->GetNcoeffs());
}
Nektar::SpatialDomains::HexGeom::~HexGeom ( )

Definition at line 201 of file HexGeom.cpp.

{
}

Member Function Documentation

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

Definition at line 906 of file HexGeom.cpp.

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

Referenced by HexGeom().

{
// This 2D array holds the local id's of all the vertices
// for every edge. For every edge, they are ordered to what we
// define as being Forwards
const unsigned int edgeVerts[kNedges][2] =
{ {0,1} ,
{1,2} ,
{2,3} ,
{3,0} ,
{0,4} ,
{1,5} ,
{2,6} ,
{3,7} ,
{4,5} ,
{5,6} ,
{6,7} ,
{7,4} };
int i;
for(i = 0; i < kNedges; i++)
{
if( m_edges[i]->GetVid(0) == m_verts[ edgeVerts[i][0] ]->GetVid() )
{
m_eorient[i] = StdRegions::eForwards;
}
else if( m_edges[i]->GetVid(0) == m_verts[ edgeVerts[i][1] ]->GetVid() )
{
m_eorient[i] = StdRegions::eBackwards;
}
else
{
ASSERTL0(false,"Could not find matching vertex for the edge");
}
}
}
void Nektar::SpatialDomains::HexGeom::SetUpFaceOrientation ( )
private

Definition at line 697 of file HexGeom.cpp.

References ASSERTL0, ASSERTL1, Nektar::SpatialDomains::Geometry::GetVertex(), Nektar::SpatialDomains::Geometry::GetVid(), Nektar::NekConstants::kNekZeroTol, kNfaces, kNqfaces, kNtfaces, kNverts, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry3D::m_faces, Nektar::SpatialDomains::Geometry3D::m_forient, and Nektar::SpatialDomains::Geometry3D::m_verts.

Referenced by HexGeom().

{
int f,i;
// These arrays represent the vector of the A and B
// coordinate of the local elemental coordinate system
// where A corresponds with the coordinate direction xi_i
// with the lowest index i (for that particular face)
// Coordinate 'B' then corresponds to the other local
// coordinate (i.e. with the highest index)
Array<OneD,NekDouble> elementAaxis(m_coordim);
Array<OneD,NekDouble> elementBaxis(m_coordim);
// These arrays correspond to the local coordinate
// system of the face itself (i.e. the Geometry2D)
// faceAaxis correspond to the xi_0 axis
// faceBaxis correspond to the xi_1 axis
Array<OneD,NekDouble> faceAaxis(m_coordim);
Array<OneD,NekDouble> faceBaxis(m_coordim);
// This is the base vertex of the face (i.e. the Geometry2D)
// This corresponds to thevertex with local ID 0 of the
// Geometry2D
unsigned int baseVertex;
// The lenght of the vectors above
NekDouble elementAaxis_length;
NekDouble elementBaxis_length;
NekDouble faceAaxis_length;
NekDouble faceBaxis_length;
// This 2D array holds the local id's of all the vertices
// for every face. For every face, they are ordered in such
// a way that the implementation below allows a unified approach
// for all faces.
const unsigned int faceVerts[kNfaces][QuadGeom::kNverts] =
{ {0,1,2,3} ,
{0,1,5,4} ,
{1,2,6,5} ,
{3,2,6,7} ,
{0,3,7,4} ,
{4,5,6,7} };
NekDouble dotproduct1 = 0.0;
NekDouble dotproduct2 = 0.0;
unsigned int orientation;
// Loop over all the faces to set up the orientation
for(f = 0; f < kNqfaces + kNtfaces; f++)
{
// initialisation
elementAaxis_length = 0.0;
elementBaxis_length = 0.0;
faceAaxis_length = 0.0;
faceBaxis_length = 0.0;
dotproduct1 = 0.0;
dotproduct2 = 0.0;
baseVertex = m_faces[f]->GetVid(0);
// We are going to construct the vectors representing the A and B axis
// of every face. These vectors will be constructed as a vector-representation
// of the edges of the face. However, for both coordinate directions, we can
// represent the vectors by two different edges. That's why we need to make sure that
// we pick the edge to which the baseVertex of the Geometry2D-representation of the face
// belongs...
if( baseVertex == m_verts[ faceVerts[f][0] ]->GetVid() )
{
for(i = 0; i < m_coordim; i++)
{
elementAaxis[i] = (*m_verts[ faceVerts[f][1] ])[i] - (*m_verts[ faceVerts[f][0] ])[i];
elementBaxis[i] = (*m_verts[ faceVerts[f][3] ])[i] - (*m_verts[ faceVerts[f][0] ])[i];
}
}
else if( baseVertex == m_verts[ faceVerts[f][1] ]->GetVid() )
{
for(i = 0; i < m_coordim; i++)
{
elementAaxis[i] = (*m_verts[ faceVerts[f][1] ])[i] - (*m_verts[ faceVerts[f][0] ])[i];
elementBaxis[i] = (*m_verts[ faceVerts[f][2] ])[i] - (*m_verts[ faceVerts[f][1] ])[i];
}
}
else if( baseVertex == m_verts[ faceVerts[f][2] ]->GetVid() )
{
for(i = 0; i < m_coordim; i++)
{
elementAaxis[i] = (*m_verts[ faceVerts[f][2] ])[i] - (*m_verts[ faceVerts[f][3] ])[i];
elementBaxis[i] = (*m_verts[ faceVerts[f][2] ])[i] - (*m_verts[ faceVerts[f][1] ])[i];
}
}
else if( baseVertex == m_verts[ faceVerts[f][3] ]->GetVid() )
{
for(i = 0; i < m_coordim; i++)
{
elementAaxis[i] = (*m_verts[ faceVerts[f][2] ])[i] - (*m_verts[ faceVerts[f][3] ])[i];
elementBaxis[i] = (*m_verts[ faceVerts[f][3] ])[i] - (*m_verts[ faceVerts[f][0] ])[i];
}
}
else
{
ASSERTL0(false, "Could not find matching vertex for the face");
}
// Now, construct the edge-vectors of the local coordinates of
// the Geometry2D-representation of the face
for(i = 0; i < m_coordim; i++)
{
faceAaxis[i] = (*m_faces[f]->GetVertex(1))[i] - (*m_faces[f]->GetVertex(0))[i];
faceBaxis[i] = (*m_faces[f]->GetVertex(3))[i] - (*m_faces[f]->GetVertex(0))[i];
elementAaxis_length += pow(elementAaxis[i],2);
elementBaxis_length += pow(elementBaxis[i],2);
faceAaxis_length += pow(faceAaxis[i],2);
faceBaxis_length += pow(faceBaxis[i],2);
}
elementAaxis_length = sqrt(elementAaxis_length);
elementBaxis_length = sqrt(elementBaxis_length);
faceAaxis_length = sqrt(faceAaxis_length);
faceBaxis_length = sqrt(faceBaxis_length);
// Calculate the inner product of both the A-axis
// (i.e. Elemental A axis and face A axis)
for(i = 0 ; i < m_coordim; i++)
{
dotproduct1 += elementAaxis[i]*faceAaxis[i];
}
orientation = 0;
// if the innerproduct is equal to the (absolute value of the ) products of the lengths
// of both vectors, then, the coordinate systems will NOT be transposed
if( fabs(elementAaxis_length*faceAaxis_length - fabs(dotproduct1)) < NekConstants::kNekZeroTol )
{
// if the inner product is negative, both A-axis point
// in reverse direction
if(dotproduct1 < 0.0)
{
orientation += 2;
}
// calculate the inner product of both B-axis
for(i = 0 ; i < m_coordim; i++)
{
dotproduct2 += elementBaxis[i]*faceBaxis[i];
}
// check that both these axis are indeed parallel
ASSERTL1(fabs(elementBaxis_length*faceBaxis_length - fabs(dotproduct2)) <
NekConstants::kNekZeroTol,
"These vectors should be parallel");
// if the inner product is negative, both B-axis point
// in reverse direction
if( dotproduct2 < 0.0 )
{
orientation++;
}
}
// The coordinate systems are transposed
else
{
orientation = 4;
// Calculate the inner product between the elemental A-axis
// and the B-axis of the face (which are now the corresponding axis)
dotproduct1 = 0.0;
for(i = 0 ; i < m_coordim; i++)
{
dotproduct1 += elementAaxis[i]*faceBaxis[i];
}
// check that both these axis are indeed parallel
ASSERTL1(fabs(elementAaxis_length*faceBaxis_length - fabs(dotproduct1)) <
NekConstants::kNekZeroTol,
"These vectors should be parallel");
// if the result is negative, both axis point in reverse
// directions
if(dotproduct1 < 0.0)
{
orientation += 2;
}
// Do the same for the other two corresponding axis
dotproduct2 = 0.0;
for(i = 0 ; i < m_coordim; i++)
{
dotproduct2 += elementBaxis[i]*faceAaxis[i];
}
// check that both these axis are indeed parallel
ASSERTL1(fabs(elementBaxis_length*faceAaxis_length - fabs(dotproduct2)) <
NekConstants::kNekZeroTol,
"These vectors should be parallel");
if( dotproduct2 < 0.0 )
{
orientation++;
}
}
orientation = orientation + 5;
// Fill the m_forient array
m_forient[f] = (StdRegions::Orientation) orientation;
}
}
void Nektar::SpatialDomains::HexGeom::SetUpLocalEdges ( )
private

Definition at line 474 of file HexGeom.cpp.

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

Referenced by HexGeom().

{
// find edge 0
int i,j;
unsigned int check;
SegGeomSharedPtr edge;
// First set up the 4 bottom edges
int f;
for(f = 1; f < 5 ; f++)
{
check = 0;
for(i = 0; i < 4; i++)
{
for(j = 0; j < 4; j++)
{
if( (m_faces[0])->GetEid(i) == (m_faces[f])->GetEid(j) )
{
edge = boost::dynamic_pointer_cast<SegGeom>((m_faces[0])->GetEdge(i));
m_edges.push_back(edge);
check++;
}
}
}
if( check < 1 )
{
std::ostringstream errstrm;
errstrm << "Connected faces do not share an edge. Faces ";
errstrm << (m_faces[0])->GetFid() << ", " << (m_faces[f])->GetFid();
ASSERTL0(false, errstrm.str());
}
else if( check > 1)
{
std::ostringstream errstrm;
errstrm << "Connected faces share more than one edge. Faces ";
errstrm << (m_faces[0])->GetFid() << ", " << (m_faces[f])->GetFid();
ASSERTL0(false, errstrm.str());
}
}
// Then, set up the 4 vertical edges
check = 0;
for(i = 0; i < 4; i++)
{
for(j = 0; j < 4; j++)
{
if( (m_faces[1])->GetEid(i) == (m_faces[4])->GetEid(j) )
{
edge = boost::dynamic_pointer_cast<SegGeom>((m_faces[1])->GetEdge(i));
m_edges.push_back(edge);
check++;
}
}
}
if( check < 1 )
{
std::ostringstream errstrm;
errstrm << "Connected faces do not share an edge. Faces ";
errstrm << (m_faces[1])->GetFid() << ", " << (m_faces[4])->GetFid();
ASSERTL0(false, errstrm.str());
}
else if( check > 1)
{
std::ostringstream errstrm;
errstrm << "Connected faces share more than one edge. Faces ";
errstrm << (m_faces[1])->GetFid() << ", " << (m_faces[4])->GetFid();
ASSERTL0(false, errstrm.str());
}
for(f = 1; f < 4 ; f++)
{
check = 0;
for(i = 0; i < 4; i++)
{
for(j = 0; j < 4; j++)
{
if( (m_faces[f])->GetEid(i) == (m_faces[f+1])->GetEid(j) )
{
edge = boost::dynamic_pointer_cast<SegGeom>((m_faces[f])->GetEdge(i));
m_edges.push_back(edge);
check++;
}
}
}
if( check < 1 )
{
std::ostringstream errstrm;
errstrm << "Connected faces do not share an edge. Faces ";
errstrm << (m_faces[f])->GetFid() << ", " << (m_faces[f+1])->GetFid();
ASSERTL0(false, errstrm.str());
}
else if( check > 1)
{
std::ostringstream errstrm;
errstrm << "Connected faces share more than one edge. Faces ";
errstrm << (m_faces[f])->GetFid() << ", " << (m_faces[f+1])->GetFid();
ASSERTL0(false, errstrm.str());
}
}
// Finally, set up the 4 top vertices
for(f = 1; f < 5 ; f++)
{
check = 0;
for(i = 0; i < 4; i++)
{
for(j = 0; j < 4; j++)
{
if( (m_faces[5])->GetEid(i) == (m_faces[f])->GetEid(j) )
{
edge = boost::dynamic_pointer_cast<SegGeom>((m_faces[5])->GetEdge(i));
m_edges.push_back(edge);
check++;
}
}
}
if( check < 1 )
{
std::ostringstream errstrm;
errstrm << "Connected faces do not share an edge. Faces ";
errstrm << (m_faces[5])->GetFid() << ", " << (m_faces[f])->GetFid();
ASSERTL0(false, errstrm.str());
}
else if( check > 1)
{
std::ostringstream errstrm;
errstrm << "Connected faces share more than one edge. Faces ";
errstrm << (m_faces[5])->GetFid() << ", " << (m_faces[f])->GetFid();
ASSERTL0(false, errstrm.str());
}
}
}
void Nektar::SpatialDomains::HexGeom::SetUpLocalVertices ( )
private

Definition at line 610 of file HexGeom.cpp.

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

Referenced by HexGeom().

{
// Set up the first 2 vertices (i.e. vertex 0,1)
if( ( m_edges[0]->GetVid(0) == m_edges[1]->GetVid(0) ) ||
( m_edges[0]->GetVid(0) == m_edges[1]->GetVid(1) ) )
{
m_verts.push_back(m_edges[0]->GetVertex(1));
m_verts.push_back(m_edges[0]->GetVertex(0));
}
else if( ( m_edges[0]->GetVid(1) == m_edges[1]->GetVid(0) ) ||
( m_edges[0]->GetVid(1) == m_edges[1]->GetVid(1) ) )
{
m_verts.push_back(m_edges[0]->GetVertex(0));
m_verts.push_back(m_edges[0]->GetVertex(1));
}
else
{
std::ostringstream errstrm;
errstrm << "Connected edges do not share a vertex. Edges ";
errstrm << m_edges[0]->GetEid() << ", " << m_edges[1]->GetEid();
ASSERTL0(false, errstrm.str());
}
// set up the other bottom vertices (i.e. vertex 2,3)
int i;
for(i = 1; i < 3; i++)
{
if( m_edges[i]->GetVid(0) == m_verts[i]->GetVid() )
{
m_verts.push_back(m_edges[i]->GetVertex(1));
}
else if( m_edges[i]->GetVid(1) == m_verts[i]->GetVid() )
{
m_verts.push_back(m_edges[i]->GetVertex(0));
}
else
{
std::ostringstream errstrm;
errstrm << "Connected edges do not share a vertex. Edges ";
errstrm << m_edges[i]->GetEid() << ", " << m_edges[i-1]->GetEid();
ASSERTL0(false, errstrm.str());
}
}
// set up top vertices
// First, set up vertices 4,5
if( ( m_edges[8]->GetVid(0) == m_edges[9]->GetVid(0) ) ||
( m_edges[8]->GetVid(0) == m_edges[9]->GetVid(1) ) )
{
m_verts.push_back(m_edges[8]->GetVertex(1));
m_verts.push_back(m_edges[8]->GetVertex(0));
}
else if( ( m_edges[8]->GetVid(1) == m_edges[9]->GetVid(0) ) ||
( m_edges[8]->GetVid(1) == m_edges[9]->GetVid(1) ) )
{
m_verts.push_back(m_edges[8]->GetVertex(0));
m_verts.push_back(m_edges[8]->GetVertex(1));
}
else
{
std::ostringstream errstrm;
errstrm << "Connected edges do not share a vertex. Edges ";
errstrm << m_edges[8]->GetEid() << ", " << m_edges[9]->GetEid();
ASSERTL0(false, errstrm.str());
}
// set up the other top vertices (i.e. vertex 6,7)
for(i = 9; i < 11; i++)
{
if( m_edges[i]->GetVid(0) == m_verts[i-4]->GetVid() )
{
m_verts.push_back(m_edges[i]->GetVertex(1));
}
else if( m_edges[i]->GetVid(1) == m_verts[i-4]->GetVid() )
{
m_verts.push_back(m_edges[i]->GetVertex(0));
}
else
{
std::ostringstream errstrm;
errstrm << "Connected edges do not share a vertex. Edges ";
errstrm << m_edges[i]->GetEid() << ", " << m_edges[i-1]->GetEid();
ASSERTL0(false, errstrm.str());
}
}
}
bool Nektar::SpatialDomains::HexGeom::v_ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
NekDouble  tol = 0.0 
)
protectedvirtual

Determines if a point specified in global coordinates is located within this tetrahedral geometry.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 357 of file HexGeom.cpp.

References Nektar::SpatialDomains::Geometry::GetCoordim().

Referenced by v_ContainsPoint().

{
Array<OneD,NekDouble> locCoord(GetCoordim(),0.0);
return v_ContainsPoint(gloCoord,locCoord,tol);
}
bool Nektar::SpatialDomains::HexGeom::v_ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
Array< OneD, NekDouble > &  locCoord,
NekDouble  tol 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 365 of file HexGeom.cpp.

References v_ContainsPoint().

{
NekDouble resid;
return v_ContainsPoint(gloCoord,locCoord,tol,resid);
}
bool Nektar::SpatialDomains::HexGeom::v_ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
Array< OneD, NekDouble > &  locCoord,
NekDouble  tol,
NekDouble &  resid 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 374 of file HexGeom.cpp.

References ASSERTL1, Nektar::SpatialDomains::eRegular, Nektar::SpatialDomains::Geometry::GetMetricInfo(), Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_xmap, npts, Nektar::SpatialDomains::Geometry3D::v_FillGeom(), v_GetLocCoords(), Vmath::Vmax(), and Vmath::Vmin().

{
ASSERTL1(gloCoord.num_elements() == 3,
"Three dimensional geometry expects three coordinates.");
// find min, max point and check if within twice this
// distance other false this is advisable since
// GetLocCoord is expensive for non regular elements.
if(GetMetricInfo()->GetGtype() != eRegular)
{
int i;
Array<OneD, NekDouble> mincoord(3), maxcoord(3);
NekDouble diff = 0.0;
v_FillGeom();
const int npts = m_xmap->GetTotPoints();
Array<OneD, NekDouble> pts(npts);
for(i = 0; i < 3; ++i)
{
m_xmap->BwdTrans(m_coeffs[i], pts);
mincoord[i] = Vmath::Vmin(pts.num_elements(),pts,1);
maxcoord[i] = Vmath::Vmax(pts.num_elements(),pts,1);
diff = max(maxcoord[i] - mincoord[i],diff);
}
for(i = 0; i < 3; ++i)
{
if((gloCoord[i] < mincoord[i] - 0.2*diff)||
(gloCoord[i] > maxcoord[i] + 0.2*diff))
{
return false;
}
}
}
v_GetLocCoords(gloCoord, locCoord);
if (locCoord[0] >= -(1+tol) && locCoord[0] <= 1+tol
&& locCoord[1] >= -(1+tol) && locCoord[1] <= 1+tol
&& locCoord[2] >= -(1+tol) && locCoord[2] <= 1+tol)
{
return true;
}
return false;
}
void Nektar::SpatialDomains::HexGeom::v_GenGeomFactors ( )
protectedvirtual

Generate the geometry factors for this element.

Reimplemented from Nektar::SpatialDomains::Geometry3D.

Definition at line 206 of file HexGeom.cpp.

References Nektar::SpatialDomains::eDeformed, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::eRegular, Nektar::NekConstants::kNekZeroTol, kNfaces, kNverts, Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry::m_geomFactors, Nektar::SpatialDomains::Geometry::m_geomFactorsState, Nektar::SpatialDomains::Geometry3D::m_verts, Nektar::SpatialDomains::Geometry::m_xmap, and Nektar::SpatialDomains::Geometry3D::v_FillGeom().

{
if (m_geomFactorsState != ePtsFilled)
{
int i,f;
GeomType Gtype = eRegular;
v_FillGeom();
// check to see if expansions are linear
for(i = 0; i < m_coordim; ++i)
{
if (m_xmap->GetBasisNumModes(0) != 2 ||
m_xmap->GetBasisNumModes(1) != 2 ||
m_xmap->GetBasisNumModes(2) != 2 )
{
Gtype = eDeformed;
}
}
// check to see if all faces are parallelograms
if(Gtype == eRegular)
{
const unsigned int faceVerts[kNfaces][QuadGeom::kNverts] =
{ {0,1,2,3} ,
{0,1,5,4} ,
{1,2,6,5} ,
{3,2,6,7} ,
{0,3,7,4} ,
{4,5,6,7} };
for(f = 0; f < kNfaces; f++)
{
// Ensure each face is a parallelogram? Check this.
for (i = 0; i < m_coordim; i++)
{
if( fabs( (*m_verts[ faceVerts[f][0] ])(i) -
(*m_verts[ faceVerts[f][1] ])(i) +
(*m_verts[ faceVerts[f][2] ])(i) -
(*m_verts[ faceVerts[f][3] ])(i) )
> NekConstants::kNekZeroTol )
{
Gtype = eDeformed;
break;
}
}
if (Gtype == eDeformed)
{
break;
}
}
}
m_geomFactors = MemoryManager<GeomFactors>::AllocateSharedPtr(
Gtype, m_coordim, m_xmap, m_coeffs);
m_geomFactorsState = ePtsFilled;
}
}
int Nektar::SpatialDomains::HexGeom::v_GetDir ( const int  faceidx,
const int  facedir 
) const
protectedvirtual

Implements Nektar::SpatialDomains::Geometry3D.

Definition at line 458 of file HexGeom.cpp.

{
if (faceidx == 0 || faceidx == 1)
{
return facedir;
}
else if (faceidx == 1 || faceidx == 3)
{
return 2 * facedir;
}
else
{
return 1 + facedir;
}
}
int Nektar::SpatialDomains::HexGeom::v_GetEdgeFaceMap ( const int  i,
const int  j 
) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 453 of file HexGeom.cpp.

References EdgeFaceConnectivity.

{
return EdgeFaceConnectivity[i][j];
}
NekDouble Nektar::SpatialDomains::HexGeom::v_GetLocCoords ( const Array< OneD, const NekDouble > &  coords,
Array< OneD, NekDouble > &  Lcoords 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 266 of file HexGeom.cpp.

References Nektar::SpatialDomains::eRegular, Nektar::SpatialDomains::Geometry::GetMetricInfo(), Vmath::Imin(), Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry3D::NewtonIterationForLocCoord(), npts, Vmath::Sadd(), Nektar::SpatialDomains::Geometry3D::v_FillGeom(), Vmath::Vmul(), and Vmath::Vvtvp().

Referenced by v_ContainsPoint().

{
NekDouble resid = 0.0;
int i;
v_FillGeom();
// calculate local coordinate for coord
if(GetMetricInfo()->GetGtype() == eRegular)
{
NekDouble len0 = 0.0 ;
NekDouble len1 = 0.0;
NekDouble len2 = 0.0;
NekDouble xi0 = 0.0;
NekDouble xi1 = 0.0;
NekDouble xi2 = 0.0;
Array<OneD, NekDouble> pts(m_xmap->GetTotPoints());
int nq0, nq1, nq2;
// get points;
//find end points
for(i = 0; i < m_coordim; ++i)
{
nq0 = m_xmap->GetNumPoints(0);
nq1 = m_xmap->GetNumPoints(1);
nq2 = m_xmap->GetNumPoints(2);
m_xmap->BwdTrans(m_coeffs[i], pts);
// use projection to side 1 to determine xi_1 coordinate based on length
len0 += (pts[nq0-1]-pts[0])*(pts[nq0-1]-pts[0]);
xi0 += (coords[i] -pts[0])*(pts[nq0-1]-pts[0]);
// use projection to side 4 to determine xi_2 coordinate based on length
len1 += (pts[nq0*(nq1-1)]-pts[0])*(pts[nq0*(nq1-1)]-pts[0]);
xi1 += (coords[i] -pts[0])*(pts[nq0*(nq1-1)]-pts[0]);
// use projection to side 4 to determine xi_3 coordinate based on length
len2 += (pts[nq0*nq1*(nq2-1)]-pts[0])*(pts[nq0*nq1*(nq2-1)]-pts[0]);
xi2 += (coords[i] -pts[0])*(pts[nq0*nq1*(nq2-1)]-pts[0]);
}
Lcoords[0] = 2*xi0/len0-1.0;
Lcoords[1] = 2*xi1/len1-1.0;
Lcoords[2] = 2*xi2/len2-1.0;
}
else
{
// Determine nearest point of coords to values in m_xmap
int npts = m_xmap->GetTotPoints();
Array<OneD, NekDouble> ptsx(npts), ptsy(npts), ptsz(npts);
Array<OneD, NekDouble> tmp1(npts), tmp2(npts);
m_xmap->BwdTrans(m_coeffs[0], ptsx);
m_xmap->BwdTrans(m_coeffs[1], ptsy);
m_xmap->BwdTrans(m_coeffs[2], ptsz);
const Array<OneD, const NekDouble> za = m_xmap->GetPoints(0);
const Array<OneD, const NekDouble> zb = m_xmap->GetPoints(1);
const Array<OneD, const NekDouble> zc = m_xmap->GetPoints(2);
//guess the first local coords based on nearest point
Vmath::Sadd(npts, -coords[0], ptsx,1,tmp1,1);
Vmath::Vmul (npts, tmp1,1,tmp1,1,tmp1,1);
Vmath::Sadd(npts, -coords[1], ptsy,1,tmp2,1);
Vmath::Vvtvp(npts, tmp2,1,tmp2,1,tmp1,1,tmp1,1);
Vmath::Sadd(npts, -coords[2], ptsz,1,tmp2,1);
Vmath::Vvtvp(npts, tmp2,1,tmp2,1,tmp1,1,tmp1,1);
int min_i = Vmath::Imin(npts,tmp1,1);
// Get Local coordinates
int qa = za.num_elements(), qb = zb.num_elements();
Lcoords[2] = zc[min_i/(qa*qb)];
min_i = min_i%(qa*qb);
Lcoords[1] = zb[min_i/qa];
Lcoords[0] = za[min_i%qa];
// Perform newton iteration to find local coordinates
NewtonIterationForLocCoord(coords, ptsx, ptsy, ptsz, Lcoords,
resid);
}
return resid;
}
int Nektar::SpatialDomains::HexGeom::v_GetNumEdges ( ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 433 of file HexGeom.cpp.

{
return 12;
}
int Nektar::SpatialDomains::HexGeom::v_GetNumFaces ( ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 438 of file HexGeom.cpp.

{
return 6;
}
int Nektar::SpatialDomains::HexGeom::v_GetNumVerts ( ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 428 of file HexGeom.cpp.

{
return 8;
}
int Nektar::SpatialDomains::HexGeom::v_GetVertexEdgeMap ( const int  i,
const int  j 
) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 443 of file HexGeom.cpp.

References VertexEdgeConnectivity.

{
return VertexEdgeConnectivity[i][j];
}
int Nektar::SpatialDomains::HexGeom::v_GetVertexFaceMap ( const int  i,
const int  j 
) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 448 of file HexGeom.cpp.

References VertexFaceConnectivity.

{
return VertexFaceConnectivity[i][j];
}

Member Data Documentation

const unsigned int Nektar::SpatialDomains::HexGeom::EdgeFaceConnectivity
staticprivate
Initial value:
{
{0,1},{0,2},{0,3},{0,4},{1,4},{1,2},{2,3},{3,4},
{1,5},{2,5},{3,5},{4,5}}

Definition at line 103 of file HexGeom.h.

Referenced by v_GetEdgeFaceMap().

const int Nektar::SpatialDomains::HexGeom::kNedges = 12
static

Definition at line 61 of file HexGeom.h.

Referenced by HexGeom(), and SetUpEdgeOrientation().

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

Definition at line 64 of file HexGeom.h.

Referenced by HexGeom(), Nektar::SpatialDomains::MeshGraph3D::ReadElements(), SetUpFaceOrientation(), and v_GenGeomFactors().

const int Nektar::SpatialDomains::HexGeom::kNqfaces = 6
static

Definition at line 62 of file HexGeom.h.

Referenced by Nektar::SpatialDomains::MeshGraph3D::ReadElements(), and SetUpFaceOrientation().

const int Nektar::SpatialDomains::HexGeom::kNtfaces = 0
static

Definition at line 63 of file HexGeom.h.

Referenced by Nektar::SpatialDomains::MeshGraph3D::ReadElements(), and SetUpFaceOrientation().

const int Nektar::SpatialDomains::HexGeom::kNverts = 8
static

Definition at line 60 of file HexGeom.h.

Referenced by SetUpFaceOrientation(), and v_GenGeomFactors().

const unsigned int Nektar::SpatialDomains::HexGeom::VertexEdgeConnectivity
staticprivate
Initial value:
{
{0,3,4},{0,1,5},{1,2,6},{2,3,7},
{4,8,11},{5,8,9},{6,9,10},{7,10,11}}

Definition at line 101 of file HexGeom.h.

Referenced by v_GetVertexEdgeMap().

const unsigned int Nektar::SpatialDomains::HexGeom::VertexFaceConnectivity
staticprivate
Initial value:
{
{0,1,4},{0,1,2},{0,2,3},{0,3,4},
{1,4,5},{1,2,5},{2,3,5},{3,4,5}}

Definition at line 102 of file HexGeom.h.

Referenced by v_GetVertexFaceMap().

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

Definition at line 65 of file HexGeom.h.

Generated on Sun Mar 15 2015 00:12:03 for Nektar++ by   doxygen 1.8.1.2