Nektar++
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. More...
 
const Geometry1DSharedPtr GetEdge (int i) const
 
Geometry2DSharedPtr GetFace (int i)
 Return face i in this element. More...
 
int GetDir (const int faceidx, const int facedir) const
 Returns the element coordinate direction corresponding to a given face coordinate direction. More...
 
- 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 GetEid (int i) const
 
int GetFid (int i) const
 
int GetTid (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
 
StdRegions::Orientation GetForient (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 More...
 
int GetEdgeFaceMap (int i, int j) const
 
void FillGeom ()
 Put all quadrature information into face/edge structure and backward transform. More...
 
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. More...
 
void SetOwnData ()
 
const LibUtilities::BasisSharedPtr GetBasis (const int i)
 Return the j-th basis of the i-th co-ordinate dimension. More...
 
const LibUtilities::PointsKeyVector GetPointsKeys ()
 
void Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 

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. More...
 
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
 
virtual void v_Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object: unset the current state and remove allocated GeomFactors. More...
 
- 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. More...
 
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 int v_GetShapeDim () const
 Return the dimension of this element. More...
 
virtual int v_GetVid (int i) const
 Return the vertex ID of vertex i. More...
 
virtual PointGeomSharedPtr v_GetVertex (int i) const
 Return vertex i in this element. More...
 
virtual const SegGeomSharedPtr v_GetEdge (int i) const
 Return edge i of this element. More...
 
virtual StdRegions::Orientation v_GetEorient (const int i) const
 Return the orientation of edge i in this element. More...
 
virtual int v_GetEid (int i) const
 Return the ID of edge i in this element. More...
 
virtual const Geometry2DSharedPtr v_GetFace (int i) const
 Return face i in this element. More...
 
virtual StdRegions::Orientation v_GetForient (const int i) const
 Return the orientation of face i in this element. More...
 
virtual int v_GetFid (int i) const
 Return the ID of face i in this element. More...
 
virtual int v_GetEid () const
 Return the ID of this element. More...
 
virtual int v_WhichEdge (SegGeomSharedPtr edge)
 Return the local ID of a given edge. More...
 
virtual int v_WhichFace (Geometry2DSharedPtr face)
 Return the local ID of a given face. More...
 
virtual const LibUtilities::BasisSharedPtr v_GetBasis (const int i)
 Return the j-th basis of the i-th co-ordinate dimension. More...
 
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)
 Initialise the m_coeffs array. More...
 

Private Member Functions

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

Static Private Attributes

static const unsigned int VertexEdgeConnectivity [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::Orientationm_eorient
 
std::vector< StdRegions::Orientationm_forient
 
std::list< CompToElmtm_elmtmap
 
bool m_owndata
 
int m_eid
 
bool m_ownverts
 
- Protected Attributes inherited from Nektar::SpatialDomains::Geometry
int m_coordim
 coordinate dimension More...
 
GeomFactorsSharedPtr m_geomFactors
 
GeomState m_geomFactorsState
 
StdRegions::StdExpansionSharedPtr m_xmap
 
GeomState m_state
 enum identifier to determine if quad points are filled More...
 
GeomType m_geomType
 
LibUtilities::ShapeType m_shapeType
 
int m_globalID
 
Array< OneD, Array< OneD, NekDouble > > m_coeffs
 
- 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 ( )
Nektar::SpatialDomains::HexGeom::HexGeom ( const QuadGeomSharedPtr  faces[])

Copy the face shared pointers

Set up orientation vectors with correct amount of elements.

Definition at line 61 of file HexGeom.cpp.

References Nektar::LibUtilities::eHexahedron, 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(), SetUpLocalVertices(), and SetUpXmap().

61  :
62  Geometry3D(faces[0]->GetEdge(0)->GetVertex(0)->GetCoordim())
63  {
65 
66  /// Copy the face shared pointers
67  m_faces.insert(m_faces.begin(), faces, faces+HexGeom::kNfaces);
68 
69  /// Set up orientation vectors with correct amount of elements.
70  m_eorient.resize(kNedges);
71  m_forient.resize(kNfaces);
72 
77  SetUpXmap();
78  SetUpCoeffs(m_xmap->GetNcoeffs());
79  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
void SetUpXmap()
Set up the m_xmap object by determining the order of each direction from derived faces.
Definition: HexGeom.cpp:873
std::vector< StdRegions::Orientation > m_eorient
Definition: Geometry3D.h:92
static const int kNfaces
Definition: HexGeom.h:64
std::vector< StdRegions::Orientation > m_forient
Definition: Geometry3D.h:93
PointGeomSharedPtr GetVertex(int i) const
Definition: Geometry.h:348
const Geometry1DSharedPtr GetEdge(int i) const
LibUtilities::ShapeType m_shapeType
Definition: Geometry.h:177
void SetUpCoeffs(const int nCoeffs)
Initialise the m_coeffs array.
Definition: Geometry.h:484
static const int kNedges
Definition: HexGeom.h:61
Nektar::SpatialDomains::HexGeom::~HexGeom ( )

Definition at line 81 of file HexGeom.cpp.

82  {
83 
84  }

Member Function Documentation

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

Definition at line 816 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().

817  {
818 
819  // This 2D array holds the local id's of all the vertices
820  // for every edge. For every edge, they are ordered to what we
821  // define as being Forwards
822  const unsigned int edgeVerts[kNedges][2] =
823  { {0,1} ,
824  {1,2} ,
825  {2,3} ,
826  {3,0} ,
827  {0,4} ,
828  {1,5} ,
829  {2,6} ,
830  {3,7} ,
831  {4,5} ,
832  {5,6} ,
833  {6,7} ,
834  {7,4} };
835 
836  int i;
837  for(i = 0; i < kNedges; i++)
838  {
839  if( m_edges[i]->GetVid(0) == m_verts[ edgeVerts[i][0] ]->GetVid() )
840  {
842  }
843  else if( m_edges[i]->GetVid(0) == m_verts[ edgeVerts[i][1] ]->GetVid() )
844  {
846  }
847  else
848  {
849  ASSERTL0(false,"Could not find matching vertex for the edge");
850  }
851  }
852  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::vector< StdRegions::Orientation > m_eorient
Definition: Geometry3D.h:92
int GetVid(int i) const
Definition: Geometry.h:319
static const int kNedges
Definition: HexGeom.h:61
void Nektar::SpatialDomains::HexGeom::SetUpFaceOrientation ( )
private

Definition at line 607 of file HexGeom.cpp.

References ASSERTL0, ASSERTL1, Nektar::SpatialDomains::Geometry::GetVertex(), Nektar::SpatialDomains::Geometry::GetVid(), Nektar::NekConstants::kNekZeroTol, kNfaces, kNqfaces, kNtfaces, Nektar::SpatialDomains::QuadGeom::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().

608  {
609  int f,i;
610 
611  // These arrays represent the vector of the A and B
612  // coordinate of the local elemental coordinate system
613  // where A corresponds with the coordinate direction xi_i
614  // with the lowest index i (for that particular face)
615  // Coordinate 'B' then corresponds to the other local
616  // coordinate (i.e. with the highest index)
617  Array<OneD,NekDouble> elementAaxis(m_coordim);
618  Array<OneD,NekDouble> elementBaxis(m_coordim);
619 
620  // These arrays correspond to the local coordinate
621  // system of the face itself (i.e. the Geometry2D)
622  // faceAaxis correspond to the xi_0 axis
623  // faceBaxis correspond to the xi_1 axis
624  Array<OneD,NekDouble> faceAaxis(m_coordim);
625  Array<OneD,NekDouble> faceBaxis(m_coordim);
626 
627  // This is the base vertex of the face (i.e. the Geometry2D)
628  // This corresponds to thevertex with local ID 0 of the
629  // Geometry2D
630  unsigned int baseVertex;
631 
632  // The lenght of the vectors above
633  NekDouble elementAaxis_length;
634  NekDouble elementBaxis_length;
635  NekDouble faceAaxis_length;
636  NekDouble faceBaxis_length;
637 
638  // This 2D array holds the local id's of all the vertices
639  // for every face. For every face, they are ordered in such
640  // a way that the implementation below allows a unified approach
641  // for all faces.
642  const unsigned int faceVerts[kNfaces][QuadGeom::kNverts] =
643  { {0,1,2,3} ,
644  {0,1,5,4} ,
645  {1,2,6,5} ,
646  {3,2,6,7} ,
647  {0,3,7,4} ,
648  {4,5,6,7} };
649 
650  NekDouble dotproduct1 = 0.0;
651  NekDouble dotproduct2 = 0.0;
652 
653  unsigned int orientation;
654 
655  // Loop over all the faces to set up the orientation
656  for(f = 0; f < kNqfaces + kNtfaces; f++)
657  {
658  // initialisation
659  elementAaxis_length = 0.0;
660  elementBaxis_length = 0.0;
661  faceAaxis_length = 0.0;
662  faceBaxis_length = 0.0;
663 
664  dotproduct1 = 0.0;
665  dotproduct2 = 0.0;
666 
667  baseVertex = m_faces[f]->GetVid(0);
668 
669  // We are going to construct the vectors representing the A and B axis
670  // of every face. These vectors will be constructed as a vector-representation
671  // of the edges of the face. However, for both coordinate directions, we can
672  // represent the vectors by two different edges. That's why we need to make sure that
673  // we pick the edge to which the baseVertex of the Geometry2D-representation of the face
674  // belongs...
675  if( baseVertex == m_verts[ faceVerts[f][0] ]->GetVid() )
676  {
677  for(i = 0; i < m_coordim; i++)
678  {
679  elementAaxis[i] = (*m_verts[ faceVerts[f][1] ])[i] - (*m_verts[ faceVerts[f][0] ])[i];
680  elementBaxis[i] = (*m_verts[ faceVerts[f][3] ])[i] - (*m_verts[ faceVerts[f][0] ])[i];
681  }
682  }
683  else if( baseVertex == m_verts[ faceVerts[f][1] ]->GetVid() )
684  {
685  for(i = 0; i < m_coordim; i++)
686  {
687  elementAaxis[i] = (*m_verts[ faceVerts[f][1] ])[i] - (*m_verts[ faceVerts[f][0] ])[i];
688  elementBaxis[i] = (*m_verts[ faceVerts[f][2] ])[i] - (*m_verts[ faceVerts[f][1] ])[i];
689  }
690  }
691  else if( baseVertex == m_verts[ faceVerts[f][2] ]->GetVid() )
692  {
693  for(i = 0; i < m_coordim; i++)
694  {
695  elementAaxis[i] = (*m_verts[ faceVerts[f][2] ])[i] - (*m_verts[ faceVerts[f][3] ])[i];
696  elementBaxis[i] = (*m_verts[ faceVerts[f][2] ])[i] - (*m_verts[ faceVerts[f][1] ])[i];
697  }
698  }
699  else if( baseVertex == m_verts[ faceVerts[f][3] ]->GetVid() )
700  {
701  for(i = 0; i < m_coordim; i++)
702  {
703  elementAaxis[i] = (*m_verts[ faceVerts[f][2] ])[i] - (*m_verts[ faceVerts[f][3] ])[i];
704  elementBaxis[i] = (*m_verts[ faceVerts[f][3] ])[i] - (*m_verts[ faceVerts[f][0] ])[i];
705  }
706  }
707  else
708  {
709  ASSERTL0(false, "Could not find matching vertex for the face");
710  }
711 
712  // Now, construct the edge-vectors of the local coordinates of
713  // the Geometry2D-representation of the face
714  for(i = 0; i < m_coordim; i++)
715  {
716  faceAaxis[i] = (*m_faces[f]->GetVertex(1))[i] - (*m_faces[f]->GetVertex(0))[i];
717  faceBaxis[i] = (*m_faces[f]->GetVertex(3))[i] - (*m_faces[f]->GetVertex(0))[i];
718 
719  elementAaxis_length += pow(elementAaxis[i],2);
720  elementBaxis_length += pow(elementBaxis[i],2);
721  faceAaxis_length += pow(faceAaxis[i],2);
722  faceBaxis_length += pow(faceBaxis[i],2);
723  }
724 
725  elementAaxis_length = sqrt(elementAaxis_length);
726  elementBaxis_length = sqrt(elementBaxis_length);
727  faceAaxis_length = sqrt(faceAaxis_length);
728  faceBaxis_length = sqrt(faceBaxis_length);
729 
730  // Calculate the inner product of both the A-axis
731  // (i.e. Elemental A axis and face A axis)
732  for(i = 0 ; i < m_coordim; i++)
733  {
734  dotproduct1 += elementAaxis[i]*faceAaxis[i];
735  }
736 
737  orientation = 0;
738  // if the innerproduct is equal to the (absolute value of the ) products of the lengths
739  // of both vectors, then, the coordinate systems will NOT be transposed
740  if( fabs(elementAaxis_length*faceAaxis_length - fabs(dotproduct1)) < NekConstants::kNekZeroTol )
741  {
742  // if the inner product is negative, both A-axis point
743  // in reverse direction
744  if(dotproduct1 < 0.0)
745  {
746  orientation += 2;
747  }
748 
749  // calculate the inner product of both B-axis
750  for(i = 0 ; i < m_coordim; i++)
751  {
752  dotproduct2 += elementBaxis[i]*faceBaxis[i];
753  }
754 
755  // check that both these axis are indeed parallel
756  ASSERTL1(fabs(elementBaxis_length*faceBaxis_length - fabs(dotproduct2)) <
758  "These vectors should be parallel");
759 
760  // if the inner product is negative, both B-axis point
761  // in reverse direction
762  if( dotproduct2 < 0.0 )
763  {
764  orientation++;
765  }
766  }
767  // The coordinate systems are transposed
768  else
769  {
770  orientation = 4;
771 
772  // Calculate the inner product between the elemental A-axis
773  // and the B-axis of the face (which are now the corresponding axis)
774  dotproduct1 = 0.0;
775  for(i = 0 ; i < m_coordim; i++)
776  {
777  dotproduct1 += elementAaxis[i]*faceBaxis[i];
778  }
779 
780  // check that both these axis are indeed parallel
781  ASSERTL1(fabs(elementAaxis_length*faceBaxis_length - fabs(dotproduct1)) <
783  "These vectors should be parallel");
784 
785  // if the result is negative, both axis point in reverse
786  // directions
787  if(dotproduct1 < 0.0)
788  {
789  orientation += 2;
790  }
791 
792  // Do the same for the other two corresponding axis
793  dotproduct2 = 0.0;
794  for(i = 0 ; i < m_coordim; i++)
795  {
796  dotproduct2 += elementBaxis[i]*faceAaxis[i];
797  }
798 
799  // check that both these axis are indeed parallel
800  ASSERTL1(fabs(elementBaxis_length*faceAaxis_length - fabs(dotproduct2)) <
802  "These vectors should be parallel");
803 
804  if( dotproduct2 < 0.0 )
805  {
806  orientation++;
807  }
808  }
809 
810  orientation = orientation + 5;
811  // Fill the m_forient array
812  m_forient[f] = (StdRegions::Orientation) orientation;
813  }
814  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static const int kNtfaces
Definition: HexGeom.h:63
static const int kNfaces
Definition: HexGeom.h:64
std::vector< StdRegions::Orientation > m_forient
Definition: Geometry3D.h:93
static const NekDouble kNekZeroTol
double NekDouble
PointGeomSharedPtr GetVertex(int i) const
Definition: Geometry.h:348
int GetVid(int i) const
Definition: Geometry.h:319
static const int kNqfaces
Definition: HexGeom.h:62
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
int m_coordim
coordinate dimension
Definition: Geometry.h:169
void Nektar::SpatialDomains::HexGeom::SetUpLocalEdges ( )
private

Definition at line 384 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().

385  {
386  // find edge 0
387  int i,j;
388  unsigned int check;
389 
390  SegGeomSharedPtr edge;
391 
392  // First set up the 4 bottom edges
393  int f;
394  for(f = 1; f < 5 ; f++)
395  {
396  check = 0;
397  for(i = 0; i < 4; i++)
398  {
399  for(j = 0; j < 4; j++)
400  {
401  if( (m_faces[0])->GetEid(i) == (m_faces[f])->GetEid(j) )
402  {
403  edge = boost::dynamic_pointer_cast<SegGeom>((m_faces[0])->GetEdge(i));
404  m_edges.push_back(edge);
405  check++;
406  }
407  }
408  }
409 
410  if( check < 1 )
411  {
412  std::ostringstream errstrm;
413  errstrm << "Connected faces do not share an edge. Faces ";
414  errstrm << (m_faces[0])->GetFid() << ", " << (m_faces[f])->GetFid();
415  ASSERTL0(false, errstrm.str());
416  }
417  else if( check > 1)
418  {
419  std::ostringstream errstrm;
420  errstrm << "Connected faces share more than one edge. Faces ";
421  errstrm << (m_faces[0])->GetFid() << ", " << (m_faces[f])->GetFid();
422  ASSERTL0(false, errstrm.str());
423  }
424  }
425 
426  // Then, set up the 4 vertical edges
427  check = 0;
428  for(i = 0; i < 4; i++)
429  {
430  for(j = 0; j < 4; j++)
431  {
432  if( (m_faces[1])->GetEid(i) == (m_faces[4])->GetEid(j) )
433  {
434  edge = boost::dynamic_pointer_cast<SegGeom>((m_faces[1])->GetEdge(i));
435  m_edges.push_back(edge);
436  check++;
437  }
438  }
439  }
440  if( check < 1 )
441  {
442  std::ostringstream errstrm;
443  errstrm << "Connected faces do not share an edge. Faces ";
444  errstrm << (m_faces[1])->GetFid() << ", " << (m_faces[4])->GetFid();
445  ASSERTL0(false, errstrm.str());
446  }
447  else if( check > 1)
448  {
449  std::ostringstream errstrm;
450  errstrm << "Connected faces share more than one edge. Faces ";
451  errstrm << (m_faces[1])->GetFid() << ", " << (m_faces[4])->GetFid();
452  ASSERTL0(false, errstrm.str());
453  }
454  for(f = 1; f < 4 ; f++)
455  {
456  check = 0;
457  for(i = 0; i < 4; i++)
458  {
459  for(j = 0; j < 4; j++)
460  {
461  if( (m_faces[f])->GetEid(i) == (m_faces[f+1])->GetEid(j) )
462  {
463  edge = boost::dynamic_pointer_cast<SegGeom>((m_faces[f])->GetEdge(i));
464  m_edges.push_back(edge);
465  check++;
466  }
467  }
468  }
469 
470  if( check < 1 )
471  {
472  std::ostringstream errstrm;
473  errstrm << "Connected faces do not share an edge. Faces ";
474  errstrm << (m_faces[f])->GetFid() << ", " << (m_faces[f+1])->GetFid();
475  ASSERTL0(false, errstrm.str());
476  }
477  else if( check > 1)
478  {
479  std::ostringstream errstrm;
480  errstrm << "Connected faces share more than one edge. Faces ";
481  errstrm << (m_faces[f])->GetFid() << ", " << (m_faces[f+1])->GetFid();
482  ASSERTL0(false, errstrm.str());
483  }
484  }
485 
486  // Finally, set up the 4 top vertices
487  for(f = 1; f < 5 ; f++)
488  {
489  check = 0;
490  for(i = 0; i < 4; i++)
491  {
492  for(j = 0; j < 4; j++)
493  {
494  if( (m_faces[5])->GetEid(i) == (m_faces[f])->GetEid(j) )
495  {
496  edge = boost::dynamic_pointer_cast<SegGeom>((m_faces[5])->GetEdge(i));
497  m_edges.push_back(edge);
498  check++;
499  }
500  }
501  }
502 
503  if( check < 1 )
504  {
505  std::ostringstream errstrm;
506  errstrm << "Connected faces do not share an edge. Faces ";
507  errstrm << (m_faces[5])->GetFid() << ", " << (m_faces[f])->GetFid();
508  ASSERTL0(false, errstrm.str());
509  }
510  else if( check > 1)
511  {
512  std::ostringstream errstrm;
513  errstrm << "Connected faces share more than one edge. Faces ";
514  errstrm << (m_faces[5])->GetFid() << ", " << (m_faces[f])->GetFid();
515  ASSERTL0(false, errstrm.str());
516  }
517  }
518  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
int GetEid(int i) const
Return the ID of edge i in this element.
Definition: Geometry3D.cpp:71
int GetFid(int i) const
Definition: Geometry.h:329
boost::shared_ptr< SegGeom > SegGeomSharedPtr
Definition: Geometry2D.h:60
const Geometry1DSharedPtr GetEdge(int i) const
void Nektar::SpatialDomains::HexGeom::SetUpLocalVertices ( )
private

Definition at line 520 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().

521  {
522  // Set up the first 2 vertices (i.e. vertex 0,1)
523  if( ( m_edges[0]->GetVid(0) == m_edges[1]->GetVid(0) ) ||
524  ( m_edges[0]->GetVid(0) == m_edges[1]->GetVid(1) ) )
525  {
526  m_verts.push_back(m_edges[0]->GetVertex(1));
527  m_verts.push_back(m_edges[0]->GetVertex(0));
528  }
529  else if( ( m_edges[0]->GetVid(1) == m_edges[1]->GetVid(0) ) ||
530  ( m_edges[0]->GetVid(1) == m_edges[1]->GetVid(1) ) )
531  {
532  m_verts.push_back(m_edges[0]->GetVertex(0));
533  m_verts.push_back(m_edges[0]->GetVertex(1));
534  }
535  else
536  {
537  std::ostringstream errstrm;
538  errstrm << "Connected edges do not share a vertex. Edges ";
539  errstrm << m_edges[0]->GetEid() << ", " << m_edges[1]->GetEid();
540  ASSERTL0(false, errstrm.str());
541  }
542 
543  // set up the other bottom vertices (i.e. vertex 2,3)
544  int i;
545  for(i = 1; i < 3; i++)
546  {
547  if( m_edges[i]->GetVid(0) == m_verts[i]->GetVid() )
548  {
549  m_verts.push_back(m_edges[i]->GetVertex(1));
550  }
551  else if( m_edges[i]->GetVid(1) == m_verts[i]->GetVid() )
552  {
553  m_verts.push_back(m_edges[i]->GetVertex(0));
554  }
555  else
556  {
557  std::ostringstream errstrm;
558  errstrm << "Connected edges do not share a vertex. Edges ";
559  errstrm << m_edges[i]->GetEid() << ", " << m_edges[i-1]->GetEid();
560  ASSERTL0(false, errstrm.str());
561  }
562  }
563 
564  // set up top vertices
565  // First, set up vertices 4,5
566  if( ( m_edges[8]->GetVid(0) == m_edges[9]->GetVid(0) ) ||
567  ( m_edges[8]->GetVid(0) == m_edges[9]->GetVid(1) ) )
568  {
569  m_verts.push_back(m_edges[8]->GetVertex(1));
570  m_verts.push_back(m_edges[8]->GetVertex(0));
571  }
572  else if( ( m_edges[8]->GetVid(1) == m_edges[9]->GetVid(0) ) ||
573  ( m_edges[8]->GetVid(1) == m_edges[9]->GetVid(1) ) )
574  {
575  m_verts.push_back(m_edges[8]->GetVertex(0));
576  m_verts.push_back(m_edges[8]->GetVertex(1));
577  }
578  else
579  {
580  std::ostringstream errstrm;
581  errstrm << "Connected edges do not share a vertex. Edges ";
582  errstrm << m_edges[8]->GetEid() << ", " << m_edges[9]->GetEid();
583  ASSERTL0(false, errstrm.str());
584  }
585 
586  // set up the other top vertices (i.e. vertex 6,7)
587  for(i = 9; i < 11; i++)
588  {
589  if( m_edges[i]->GetVid(0) == m_verts[i-4]->GetVid() )
590  {
591  m_verts.push_back(m_edges[i]->GetVertex(1));
592  }
593  else if( m_edges[i]->GetVid(1) == m_verts[i-4]->GetVid() )
594  {
595  m_verts.push_back(m_edges[i]->GetVertex(0));
596  }
597  else
598  {
599  std::ostringstream errstrm;
600  errstrm << "Connected edges do not share a vertex. Edges ";
601  errstrm << m_edges[i]->GetEid() << ", " << m_edges[i-1]->GetEid();
602  ASSERTL0(false, errstrm.str());
603  }
604  }
605  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
PointGeomSharedPtr GetVertex(int i) const
Definition: Geometry.h:348
int GetVid(int i) const
Definition: Geometry.h:319
void Nektar::SpatialDomains::HexGeom::SetUpXmap ( )
private

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

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

Definition at line 873 of file HexGeom.cpp.

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

Referenced by HexGeom(), and v_Reset().

874  {
875  /// Determine necessary order for standard region. This can almost
876  /// certainly be simplified but works for now!
877  vector<int> tmp1, tmp2;
878 
879  if (m_forient[0] < 9)
880  {
881  tmp1.push_back(m_faces[0]->GetXmap()->GetEdgeNcoeffs (0));
882  tmp1.push_back(m_faces[0]->GetXmap()->GetEdgeNcoeffs (2));
883  tmp2.push_back(m_faces[0]->GetXmap()->GetEdgeNumPoints(0));
884  tmp2.push_back(m_faces[0]->GetXmap()->GetEdgeNumPoints(2));
885  }
886  else
887  {
888  tmp1.push_back(m_faces[0]->GetXmap()->GetEdgeNcoeffs (1));
889  tmp1.push_back(m_faces[0]->GetXmap()->GetEdgeNcoeffs (3));
890  tmp2.push_back(m_faces[0]->GetXmap()->GetEdgeNumPoints(1));
891  tmp2.push_back(m_faces[0]->GetXmap()->GetEdgeNumPoints(3));
892  }
893 
894  if (m_forient[5] < 9)
895  {
896  tmp1.push_back(m_faces[5]->GetXmap()->GetEdgeNcoeffs (0));
897  tmp1.push_back(m_faces[5]->GetXmap()->GetEdgeNcoeffs (2));
898  tmp2.push_back(m_faces[5]->GetXmap()->GetEdgeNumPoints(0));
899  tmp2.push_back(m_faces[5]->GetXmap()->GetEdgeNumPoints(2));
900  }
901  else
902  {
903  tmp1.push_back(m_faces[5]->GetXmap()->GetEdgeNcoeffs (1));
904  tmp1.push_back(m_faces[5]->GetXmap()->GetEdgeNcoeffs (3));
905  tmp2.push_back(m_faces[5]->GetXmap()->GetEdgeNumPoints(1));
906  tmp2.push_back(m_faces[5]->GetXmap()->GetEdgeNumPoints(3));
907  }
908 
909  int order0 = *max_element(tmp1.begin(), tmp1.end());
910  int points0 = *max_element(tmp2.begin(), tmp2.end());
911 
912  tmp1.clear();
913  tmp2.clear();
914 
915  if (m_forient[0] < 9)
916  {
917  tmp1.push_back(m_faces[0]->GetXmap()->GetEdgeNcoeffs (1));
918  tmp1.push_back(m_faces[0]->GetXmap()->GetEdgeNcoeffs (3));
919  tmp2.push_back(m_faces[0]->GetXmap()->GetEdgeNumPoints(1));
920  tmp2.push_back(m_faces[0]->GetXmap()->GetEdgeNumPoints(3));
921  }
922  else
923  {
924  tmp1.push_back(m_faces[0]->GetXmap()->GetEdgeNcoeffs (0));
925  tmp1.push_back(m_faces[0]->GetXmap()->GetEdgeNcoeffs (2));
926  tmp2.push_back(m_faces[0]->GetXmap()->GetEdgeNumPoints(0));
927  tmp2.push_back(m_faces[0]->GetXmap()->GetEdgeNumPoints(2));
928  }
929 
930  if (m_forient[5] < 9)
931  {
932  tmp1.push_back(m_faces[5]->GetXmap()->GetEdgeNcoeffs (1));
933  tmp1.push_back(m_faces[5]->GetXmap()->GetEdgeNcoeffs (3));
934  tmp2.push_back(m_faces[5]->GetXmap()->GetEdgeNumPoints(1));
935  tmp2.push_back(m_faces[5]->GetXmap()->GetEdgeNumPoints(3));
936  }
937  else
938  {
939  tmp1.push_back(m_faces[5]->GetXmap()->GetEdgeNcoeffs (0));
940  tmp1.push_back(m_faces[5]->GetXmap()->GetEdgeNcoeffs (2));
941  tmp2.push_back(m_faces[5]->GetXmap()->GetEdgeNumPoints(0));
942  tmp2.push_back(m_faces[5]->GetXmap()->GetEdgeNumPoints(2));
943  }
944 
945  int order1 = *max_element(tmp1.begin(), tmp1.end());
946  int points1 = *max_element(tmp2.begin(), tmp2.end());
947 
948  tmp1.clear();
949  tmp2.clear();
950 
951  if (m_forient[1] < 9)
952  {
953  tmp1.push_back(m_faces[1]->GetXmap()->GetEdgeNcoeffs (1));
954  tmp1.push_back(m_faces[1]->GetXmap()->GetEdgeNcoeffs (3));
955  tmp2.push_back(m_faces[1]->GetXmap()->GetEdgeNumPoints(1));
956  tmp2.push_back(m_faces[1]->GetXmap()->GetEdgeNumPoints(3));
957  }
958  else
959  {
960  tmp1.push_back(m_faces[1]->GetXmap()->GetEdgeNcoeffs (0));
961  tmp1.push_back(m_faces[1]->GetXmap()->GetEdgeNcoeffs (2));
962  tmp2.push_back(m_faces[1]->GetXmap()->GetEdgeNumPoints(0));
963  tmp2.push_back(m_faces[1]->GetXmap()->GetEdgeNumPoints(2));
964  }
965 
966  if (m_forient[3] < 9)
967  {
968  tmp1.push_back(m_faces[3]->GetXmap()->GetEdgeNcoeffs (1));
969  tmp1.push_back(m_faces[3]->GetXmap()->GetEdgeNcoeffs (3));
970  tmp2.push_back(m_faces[3]->GetXmap()->GetEdgeNumPoints(1));
971  tmp2.push_back(m_faces[3]->GetXmap()->GetEdgeNumPoints(3));
972  }
973  else
974  {
975  tmp1.push_back(m_faces[3]->GetXmap()->GetEdgeNcoeffs (0));
976  tmp1.push_back(m_faces[3]->GetXmap()->GetEdgeNcoeffs (2));
977  tmp2.push_back(m_faces[3]->GetXmap()->GetEdgeNumPoints(0));
978  tmp2.push_back(m_faces[3]->GetXmap()->GetEdgeNumPoints(2));
979  }
980 
981  int order2 = *max_element(tmp1.begin(), tmp1.end());
982  int points2 = *max_element(tmp2.begin(), tmp2.end());
983 
984  const LibUtilities::BasisKey A(
986  LibUtilities::PointsKey(
988  const LibUtilities::BasisKey B(
990  LibUtilities::PointsKey(
992  const LibUtilities::BasisKey C(
994  LibUtilities::PointsKey(
996 
998  A, B, C);
999  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
Principle Modified Functions .
Definition: BasisType.h:49
StdRegions::StdExpansionSharedPtr GetXmap() const
Definition: Geometry.h:383
std::vector< StdRegions::Orientation > m_forient
Definition: Geometry3D.h:93
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:50
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 249 of file HexGeom.cpp.

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

Referenced by v_ContainsPoint().

251  {
252  Array<OneD,NekDouble> locCoord(GetCoordim(),0.0);
253  return v_ContainsPoint(gloCoord,locCoord,tol);
254  }
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...
Definition: HexGeom.cpp:249
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 257 of file HexGeom.cpp.

References v_ContainsPoint().

261  {
262  NekDouble resid;
263  return v_ContainsPoint(gloCoord,locCoord,tol,resid);
264  }
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...
Definition: HexGeom.cpp:249
double NekDouble
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 266 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().

271  {
272  ASSERTL1(gloCoord.num_elements() == 3,
273  "Three dimensional geometry expects three coordinates.");
274 
275  // find min, max point and check if within twice this
276  // distance other false this is advisable since
277  // GetLocCoord is expensive for non regular elements.
278  if(GetMetricInfo()->GetGtype() != eRegular)
279  {
280  int i;
281  Array<OneD, NekDouble> mincoord(3), maxcoord(3);
282  NekDouble diff = 0.0;
283 
284  v_FillGeom();
285 
286  const int npts = m_xmap->GetTotPoints();
287  Array<OneD, NekDouble> pts(npts);
288 
289  for(i = 0; i < 3; ++i)
290  {
291  m_xmap->BwdTrans(m_coeffs[i], pts);
292 
293  mincoord[i] = Vmath::Vmin(pts.num_elements(),pts,1);
294  maxcoord[i] = Vmath::Vmax(pts.num_elements(),pts,1);
295 
296  diff = max(maxcoord[i] - mincoord[i],diff);
297  }
298 
299  for(i = 0; i < 3; ++i)
300  {
301  if((gloCoord[i] < mincoord[i] - 0.2*diff)||
302  (gloCoord[i] > maxcoord[i] + 0.2*diff))
303  {
304  return false;
305  }
306  }
307  }
308 
309  v_GetLocCoords(gloCoord, locCoord);
310 
311  if (locCoord[0] >= -(1+tol) && locCoord[0] <= 1+tol
312  && locCoord[1] >= -(1+tol) && locCoord[1] <= 1+tol
313  && locCoord[2] >= -(1+tol) && locCoord[2] <= 1+tol)
314  {
315  return true;
316  }
317 
318  // If out of range clamp locCoord to be within [-1,1]^3
319  // since any larger value will be very oscillatory if
320  // called by 'returnNearestElmt' option in
321  // ExpList::GetExpIndex
322  for(int i = 0; i < 3; ++i)
323  {
324  if(locCoord[i] <-(1+tol))
325  {
326  locCoord[i] = -(1+tol);
327  }
328 
329  if(locCoord[i] > (1+tol))
330  {
331  locCoord[i] = 1+tol;
332  }
333  }
334 
335  return false;
336  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.cpp:756
T Vmin(int n, const T *x, const int incx)
Return the minimum element in x - called vmin to avoid conflict with min.
Definition: Vmath.cpp:848
virtual NekDouble v_GetLocCoords(const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
Definition: HexGeom.cpp:146
static std::string npts
Definition: InputFld.cpp:43
double NekDouble
virtual void v_FillGeom()
Put all quadrature information into face/edge structure and backward transform.
Definition: Geometry3D.cpp:244
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Definition: Geometry.h:180
Geometry is straight-sided with constant geometric factors.
GeomFactorsSharedPtr GetMetricInfo()
Definition: Geometry.h:299
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
void Nektar::SpatialDomains::HexGeom::v_GenGeomFactors ( )
protectedvirtual

Generate the geometry factors for this element.

Reimplemented from Nektar::SpatialDomains::Geometry3D.

Definition at line 86 of file HexGeom.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::SpatialDomains::eDeformed, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::eRegular, Nektar::NekConstants::kNekZeroTol, kNfaces, Nektar::SpatialDomains::QuadGeom::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().

87  {
89  {
90  int i,f;
91  GeomType Gtype = eRegular;
92 
93  v_FillGeom();
94 
95  // check to see if expansions are linear
96  for(i = 0; i < m_coordim; ++i)
97  {
98  if (m_xmap->GetBasisNumModes(0) != 2 ||
99  m_xmap->GetBasisNumModes(1) != 2 ||
100  m_xmap->GetBasisNumModes(2) != 2 )
101  {
102  Gtype = eDeformed;
103  }
104  }
105 
106  // check to see if all faces are parallelograms
107  if(Gtype == eRegular)
108  {
109  const unsigned int faceVerts[kNfaces][QuadGeom::kNverts] =
110  { {0,1,2,3} ,
111  {0,1,5,4} ,
112  {1,2,6,5} ,
113  {3,2,6,7} ,
114  {0,3,7,4} ,
115  {4,5,6,7} };
116 
117  for(f = 0; f < kNfaces; f++)
118  {
119  // Ensure each face is a parallelogram? Check this.
120  for (i = 0; i < m_coordim; i++)
121  {
122  if( fabs( (*m_verts[ faceVerts[f][0] ])(i) -
123  (*m_verts[ faceVerts[f][1] ])(i) +
124  (*m_verts[ faceVerts[f][2] ])(i) -
125  (*m_verts[ faceVerts[f][3] ])(i) )
127  {
128  Gtype = eDeformed;
129  break;
130  }
131  }
132 
133  if (Gtype == eDeformed)
134  {
135  break;
136  }
137  }
138  }
139 
141  Gtype, m_coordim, m_xmap, m_coeffs);
143  }
144  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
GeomFactorsSharedPtr m_geomFactors
Definition: Geometry.h:170
static const int kNfaces
Definition: HexGeom.h:64
static const NekDouble kNekZeroTol
virtual void v_FillGeom()
Put all quadrature information into face/edge structure and backward transform.
Definition: Geometry3D.cpp:244
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Definition: Geometry.h:180
Geometry is straight-sided with constant geometric factors.
Geometric information has been generated.
GeomType
Indicates the type of element geometry.
Geometry is curved or has non-constant factors.
int m_coordim
coordinate dimension
Definition: Geometry.h:169
int Nektar::SpatialDomains::HexGeom::v_GetDir ( const int  faceidx,
const int  facedir 
) const
protectedvirtual

Implements Nektar::SpatialDomains::Geometry3D.

Definition at line 368 of file HexGeom.cpp.

369  {
370  if (faceidx == 0 || faceidx == 1)
371  {
372  return facedir;
373  }
374  else if (faceidx == 1 || faceidx == 3)
375  {
376  return 2 * facedir;
377  }
378  else
379  {
380  return 1 + facedir;
381  }
382  }
int Nektar::SpatialDomains::HexGeom::v_GetEdgeFaceMap ( const int  i,
const int  j 
) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 363 of file HexGeom.cpp.

References EdgeFaceConnectivity.

364  {
365  return EdgeFaceConnectivity[i][j];
366  }
static const unsigned int EdgeFaceConnectivity[12][2]
Definition: HexGeom.h:107
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 146 of file HexGeom.cpp.

References Nektar::SpatialDomains::PointGeom::dist(), Nektar::SpatialDomains::eRegular, Nektar::SpatialDomains::Geometry::GetMetricInfo(), Vmath::Imin(), Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry3D::m_verts, 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().

149  {
150  NekDouble ptdist = 1e6;
151  int i;
152 
153  v_FillGeom();
154 
155  // calculate local coordinate for coord
156  if(GetMetricInfo()->GetGtype() == eRegular)
157  {
158  NekDouble len0 = 0.0 ;
159  NekDouble len1 = 0.0;
160  NekDouble len2 = 0.0;
161  NekDouble xi0 = 0.0;
162  NekDouble xi1 = 0.0;
163  NekDouble xi2 = 0.0;
164  Array<OneD, NekDouble> pts(m_xmap->GetTotPoints());
165  int nq0, nq1, nq2;
166 
167  // get points;
168  //find end points
169  for(i = 0; i < m_coordim; ++i)
170  {
171  nq0 = m_xmap->GetNumPoints(0);
172  nq1 = m_xmap->GetNumPoints(1);
173  nq2 = m_xmap->GetNumPoints(2);
174 
175  m_xmap->BwdTrans(m_coeffs[i], pts);
176 
177  // use projection to side 1 to determine xi_1 coordinate based on length
178  len0 += (pts[nq0-1]-pts[0])*(pts[nq0-1]-pts[0]);
179  xi0 += (coords[i] -pts[0])*(pts[nq0-1]-pts[0]);
180 
181  // use projection to side 4 to determine xi_2 coordinate based on length
182  len1 += (pts[nq0*(nq1-1)]-pts[0])*(pts[nq0*(nq1-1)]-pts[0]);
183  xi1 += (coords[i] -pts[0])*(pts[nq0*(nq1-1)]-pts[0]);
184 
185  // use projection to side 4 to determine xi_3 coordinate based on length
186  len2 += (pts[nq0*nq1*(nq2-1)]-pts[0])*(pts[nq0*nq1*(nq2-1)]-pts[0]);
187  xi2 += (coords[i] -pts[0])*(pts[nq0*nq1*(nq2-1)]-pts[0]);
188  }
189 
190  Lcoords[0] = 2*xi0/len0-1.0;
191  Lcoords[1] = 2*xi1/len1-1.0;
192  Lcoords[2] = 2*xi2/len2-1.0;
193 
194  // Set ptdist to distance to nearest vertex
195  // Point inside tetrahedron
196  PointGeom r(m_coordim, 0, coords[0], coords[1], coords[2]);
197  for(int i = 0; i < 8; ++i)
198  {
199  ptdist = min(ptdist,r.dist(*m_verts[i]));
200  }
201  }
202  else
203  {
204  // Determine nearest point of coords to values in m_xmap
205  int npts = m_xmap->GetTotPoints();
206  Array<OneD, NekDouble> ptsx(npts), ptsy(npts), ptsz(npts);
207  Array<OneD, NekDouble> tmp1(npts), tmp2(npts);
208 
209  m_xmap->BwdTrans(m_coeffs[0], ptsx);
210  m_xmap->BwdTrans(m_coeffs[1], ptsy);
211  m_xmap->BwdTrans(m_coeffs[2], ptsz);
212 
213  const Array<OneD, const NekDouble> za = m_xmap->GetPoints(0);
214  const Array<OneD, const NekDouble> zb = m_xmap->GetPoints(1);
215  const Array<OneD, const NekDouble> zc = m_xmap->GetPoints(2);
216 
217  //guess the first local coords based on nearest point
218  Vmath::Sadd(npts, -coords[0], ptsx,1,tmp1,1);
219  Vmath::Vmul (npts, tmp1,1,tmp1,1,tmp1,1);
220  Vmath::Sadd(npts, -coords[1], ptsy,1,tmp2,1);
221  Vmath::Vvtvp(npts, tmp2,1,tmp2,1,tmp1,1,tmp1,1);
222  Vmath::Sadd(npts, -coords[2], ptsz,1,tmp2,1);
223  Vmath::Vvtvp(npts, tmp2,1,tmp2,1,tmp1,1,tmp1,1);
224 
225  int min_i = Vmath::Imin(npts,tmp1,1);
226 
227  // distance from coordinate to nearest point for return value.
228  ptdist = sqrt(tmp1[min_i]);
229 
230  // Get Local coordinates
231  int qa = za.num_elements(), qb = zb.num_elements();
232  Lcoords[2] = zc[min_i/(qa*qb)];
233  min_i = min_i%(qa*qb);
234  Lcoords[1] = zb[min_i/qa];
235  Lcoords[0] = za[min_i%qa];
236 
237  // Perform newton iteration to find local coordinates
238  NekDouble resid = 0.0;
239  NewtonIterationForLocCoord(coords, ptsx, ptsy, ptsz, Lcoords,
240  resid);
241  }
242  return ptdist;
243  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
int Imin(int n, const T *x, const int incx)
Return the index of the minimum element in x.
Definition: Vmath.cpp:824
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
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)
Definition: Geometry3D.cpp:96
static std::string npts
Definition: InputFld.cpp:43
double NekDouble
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
Definition: Vmath.cpp:301
virtual void v_FillGeom()
Put all quadrature information into face/edge structure and backward transform.
Definition: Geometry3D.cpp:244
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Definition: Geometry.h:180
Geometry is straight-sided with constant geometric factors.
GeomFactorsSharedPtr GetMetricInfo()
Definition: Geometry.h:299
int m_coordim
coordinate dimension
Definition: Geometry.h:169
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169
int Nektar::SpatialDomains::HexGeom::v_GetNumEdges ( ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 343 of file HexGeom.cpp.

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 348 of file HexGeom.cpp.

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 338 of file HexGeom.cpp.

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

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 353 of file HexGeom.cpp.

References VertexEdgeConnectivity.

354  {
355  return VertexEdgeConnectivity[i][j];
356  }
static const unsigned int VertexEdgeConnectivity[8][3]
Definition: HexGeom.h:105
int Nektar::SpatialDomains::HexGeom::v_GetVertexFaceMap ( const int  i,
const int  j 
) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 358 of file HexGeom.cpp.

References VertexFaceConnectivity.

359  {
360  return VertexFaceConnectivity[i][j];
361  }
static const unsigned int VertexFaceConnectivity[8][3]
Definition: HexGeom.h:106
void Nektar::SpatialDomains::HexGeom::v_Reset ( CurveMap curvedEdges,
CurveMap curvedFaces 
)
protectedvirtual

Reset this geometry object: unset the current state and remove allocated GeomFactors.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 854 of file HexGeom.cpp.

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

857  {
858  Geometry::v_Reset(curvedEdges, curvedFaces);
859 
860  for (int i = 0; i < 6; ++i)
861  {
862  m_faces[i]->Reset(curvedEdges, curvedFaces);
863  }
864 
865  SetUpXmap();
866  SetUpCoeffs(m_xmap->GetNcoeffs());
867  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
void SetUpXmap()
Set up the m_xmap object by determining the order of each direction from derived faces.
Definition: HexGeom.cpp:873
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state and remove allocated GeomFactors.
Definition: Geometry.cpp:307
void SetUpCoeffs(const int nCoeffs)
Initialise the m_coeffs array.
Definition: Geometry.h:484

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 107 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
const int Nektar::SpatialDomains::HexGeom::kNqfaces = 6
static
const int Nektar::SpatialDomains::HexGeom::kNtfaces = 0
static
const int Nektar::SpatialDomains::HexGeom::kNverts = 8
static

Definition at line 60 of file HexGeom.h.

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 105 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 106 of file HexGeom.h.

Referenced by v_GetVertexFaceMap().

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

Definition at line 65 of file HexGeom.h.