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

A 3-dimensional four-faced element. More...

#include <Tetrahedron.h>

Inheritance diagram for Nektar::NekMeshUtils::Tetrahedron:
Inheritance graph
[legend]
Collaboration diagram for Nektar::NekMeshUtils::Tetrahedron:
Collaboration graph
[legend]

Public Member Functions

NEKMESHUTILS_EXPORT Tetrahedron (ElmtConfig pConf, std::vector< NodeSharedPtr > pNodeList, std::vector< int > pTagList)
 Create a tetrahedron element. More...
 
NEKMESHUTILS_EXPORT Tetrahedron (const Tetrahedron &pSrc)
 
virtual NEKMESHUTILS_EXPORT ~Tetrahedron ()
 
virtual NEKMESHUTILS_EXPORT
SpatialDomains::GeometrySharedPtr 
GetGeom (int coordDim)
 Generate a Nektar++ geometry object for this element. More...
 
virtual NEKMESHUTILS_EXPORT void GetCurvedNodes (std::vector< NodeSharedPtr > &nodeList) const
 get list of volume interior nodes More...
 
virtual NEKMESHUTILS_EXPORT
StdRegions::Orientation 
GetEdgeOrient (int edgeId, EdgeSharedPtr edge)
 Get the edge orientation of edge with respect to the local element, which lies at edge index edgeId. More...
 
virtual NEKMESHUTILS_EXPORT void MakeOrder (int order, SpatialDomains::GeometrySharedPtr geom, LibUtilities::PointsType pType, int coordDim, int &id, bool justConfig=false)
 Insert interior (i.e. volume) points into this element to make the geometry an order order representation. More...
 
virtual NEKMESHUTILS_EXPORT int GetFaceVertex (int i, int j)
 Returns the local index of vertex j of face i. More...
 
- Public Member Functions inherited from Nektar::NekMeshUtils::Element
NEKMESHUTILS_EXPORT Element (ElmtConfig pConf, unsigned int pNumNodes, unsigned int pGotNodes)
 
NEKMESHUTILS_EXPORT unsigned int GetId () const
 Returns the ID of the element (or associated edge or face for boundary elements). More...
 
NEKMESHUTILS_EXPORT unsigned int GetDim () const
 Returns the expansion dimension of the element. More...
 
NEKMESHUTILS_EXPORT ElmtConfig GetConf () const
 Returns the configuration of the element. More...
 
NEKMESHUTILS_EXPORT
LibUtilities::ShapeType 
GetShapeType () const
 returns the shapetype More...
 
NEKMESHUTILS_EXPORT std::string GetTag () const
 Returns the tag which defines the element shape. More...
 
NEKMESHUTILS_EXPORT NodeSharedPtr GetVertex (unsigned int i) const
 Access a vertex node. More...
 
NEKMESHUTILS_EXPORT EdgeSharedPtr GetEdge (unsigned int i) const
 Access an edge. More...
 
NEKMESHUTILS_EXPORT FaceSharedPtr GetFace (unsigned int i) const
 Access a face. More...
 
NEKMESHUTILS_EXPORT
std::vector< NodeSharedPtr
GetVertexList () const
 Access the list of vertex nodes. More...
 
NEKMESHUTILS_EXPORT
std::vector< EdgeSharedPtr
GetEdgeList () const
 Access the list of edges. More...
 
NEKMESHUTILS_EXPORT
std::vector< FaceSharedPtr
GetFaceList () const
 Access the list of faces. More...
 
NEKMESHUTILS_EXPORT
std::vector< NodeSharedPtr
GetVolumeNodes () const
 Access the list of volume nodes. More...
 
NEKMESHUTILS_EXPORT void SetVolumeNodes (std::vector< NodeSharedPtr > &nodes)
 
NEKMESHUTILS_EXPORT
LibUtilities::PointsType 
GetCurveType () const
 
NEKMESHUTILS_EXPORT void SetCurveType (LibUtilities::PointsType cT)
 
NEKMESHUTILS_EXPORT unsigned int GetNodeCount ()
 Returns the total number of nodes (vertices, edge nodes and face nodes and volume nodes). More...
 
NEKMESHUTILS_EXPORT
std::vector< int > 
GetTagList () const
 Access the list of tags associated with this element. More...
 
NEKMESHUTILS_EXPORT unsigned int GetVertexCount () const
 Returns the number of vertices. More...
 
NEKMESHUTILS_EXPORT unsigned int GetEdgeCount () const
 Returns the number of edges. More...
 
NEKMESHUTILS_EXPORT unsigned int GetFaceCount () const
 Returns the number of faces. More...
 
NEKMESHUTILS_EXPORT void SetId (unsigned int p)
 Change the ID of the element. More...
 
NEKMESHUTILS_EXPORT void SetVertex (unsigned int p, NodeSharedPtr pNew, bool descend=true)
 Replace a vertex in the element. More...
 
NEKMESHUTILS_EXPORT void SetEdge (unsigned int p, EdgeSharedPtr pNew, bool descend=true)
 Replace an edge in the element. More...
 
NEKMESHUTILS_EXPORT void SetFace (unsigned int p, FaceSharedPtr pNew)
 Replace a face in the element. More...
 
NEKMESHUTILS_EXPORT void SetEdgeLink (EdgeSharedPtr pLink)
 Set a correspondence between this element and an edge (2D boundary element). More...
 
NEKMESHUTILS_EXPORT EdgeSharedPtr GetEdgeLink ()
 Get correspondence between this element and an edge. More...
 
NEKMESHUTILS_EXPORT void SetFaceLink (FaceSharedPtr pLink)
 Set a correspondence between this element and a face (3D boundary element). More...
 
NEKMESHUTILS_EXPORT FaceSharedPtr GetFaceLink ()
 Get correspondence between this element and a face. More...
 
NEKMESHUTILS_EXPORT void SetBoundaryLink (int i, int j)
 Set a correspondence between edge or face i and its representative boundary element m->element[expDim-1][j]. More...
 
NEKMESHUTILS_EXPORT int GetBoundaryLink (int i)
 Get the location of the boundary face/edge i for this element. More...
 
NEKMESHUTILS_EXPORT void SetTagList (const std::vector< int > &tags)
 Set the list of tags associated with this element. More...
 
virtual NEKMESHUTILS_EXPORT
std::string 
GetXmlString ()
 Generate a list of vertices (1D), edges (2D), or faces (3D). More...
 
NEKMESHUTILS_EXPORT std::string GetXmlCurveString ()
 Generates a string listing the coordinates of all nodes associated with this element. More...
 
NEKMESHUTILS_EXPORT int GetMaxOrder ()
 Obtain the order of an element by looking at edges. More...
 
NEKMESHUTILS_EXPORT void Print ()
 
virtual NEKMESHUTILS_EXPORT
Array< OneD, NekDouble
Normal (bool inward=false)
 returns the normal to the element More...
 

Static Public Member Functions

static ElementSharedPtr create (ElmtConfig pConf, std::vector< NodeSharedPtr > pNodeList, std::vector< int > pTagList)
 Creates an instance of this class. More...
 
static NEKMESHUTILS_EXPORT
unsigned int 
GetNumNodes (ElmtConfig pConf)
 Return the number of nodes defining a tetrahedron. More...
 

Public Attributes

int m_orientationMap [4]
 
int m_origVertMap [4]
 
- Public Attributes inherited from Nektar::NekMeshUtils::Element
CADObjectSharedPtr m_parentCAD
 

Static Public Attributes

static LibUtilities::ShapeType m_type
 Element type. More...
 

Protected Member Functions

void OrientTet ()
 Orient tetrahedron to align degenerate vertices. More...
 

Static Private Attributes

static int m_faceIds [4][3]
 Vertex IDs that make up tetrahedron faces. More...
 

Additional Inherited Members

- Protected Attributes inherited from Nektar::NekMeshUtils::Element
unsigned int m_id
 ID of the element. More...
 
unsigned int m_dim
 Dimension of the element. More...
 
ElmtConfig m_conf
 Contains configuration of the element. More...
 
std::string m_tag
 Tag character describing the element. More...
 
std::vector< int > m_taglist
 List of integers specifying properties of the element. More...
 
std::vector< NodeSharedPtrm_vertex
 List of element vertex nodes. More...
 
std::vector< EdgeSharedPtrm_edge
 List of element edges. More...
 
std::vector< FaceSharedPtrm_face
 List of element faces. More...
 
std::vector< NodeSharedPtrm_volumeNodes
 List of element volume nodes. More...
 
LibUtilities::PointsType m_curveType
 Volume curve type. More...
 
EdgeSharedPtr m_edgeLink
 Pointer to the corresponding edge if element is a 2D boundary. More...
 
FaceSharedPtr m_faceLink
 Pointer to the corresponding face if element is a 3D boundary. More...
 
std::map< int, int > m_boundaryLinks
 Array mapping faces/edges to the location of the appropriate boundary elements in m->element. More...
 
SpatialDomains::GeometrySharedPtr m_geom
 Nektar++ geometry object for this element. More...
 

Detailed Description

A 3-dimensional four-faced element.

Definition at line 50 of file Tetrahedron.h.

Constructor & Destructor Documentation

Nektar::NekMeshUtils::Tetrahedron::Tetrahedron ( ElmtConfig  pConf,
std::vector< NodeSharedPtr pNodeList,
std::vector< int >  pTagList 
)

Create a tetrahedron element.

Definition at line 64 of file Tetrahedron.cpp.

References Nektar::NekMeshUtils::HOTriangle< T >::Align(), ASSERTL0, Nektar::iterator, Nektar::NekMeshUtils::Element::m_conf, Nektar::NekMeshUtils::Element::m_dim, Nektar::NekMeshUtils::Element::m_edge, Nektar::NekMeshUtils::ElmtConfig::m_edgeCurveType, Nektar::NekMeshUtils::Element::m_face, Nektar::NekMeshUtils::ElmtConfig::m_faceCurveType, m_faceIds, Nektar::NekMeshUtils::ElmtConfig::m_faceNodes, Nektar::NekMeshUtils::ElmtConfig::m_order, m_orientationMap, Nektar::NekMeshUtils::ElmtConfig::m_reorient, Nektar::NekMeshUtils::Element::m_tag, Nektar::NekMeshUtils::Element::m_taglist, Nektar::NekMeshUtils::Element::m_vertex, Nektar::NekMeshUtils::ElmtConfig::m_volumeNodes, Nektar::NekMeshUtils::Element::m_volumeNodes, OrientTet(), and Nektar::NekMeshUtils::HOTriangle< T >::surfVerts.

Referenced by create().

67  : Element(pConf, GetNumNodes(pConf), pNodeList.size())
68 {
69  m_tag = "A";
70  m_dim = 3;
71  m_taglist = pTagList;
72  int n = m_conf.m_order - 1;
73 
74  // Create a map to relate edge nodes to a pair of vertices
75  // defining an edge.
76  map<pair<int, int>, int> edgeNodeMap;
77  map<pair<int, int>, int>::iterator it;
78  edgeNodeMap[pair<int, int>(1, 2)] = 5;
79  edgeNodeMap[pair<int, int>(2, 3)] = 5 + n;
80  edgeNodeMap[pair<int, int>(1, 3)] = 5 + 2 * n;
81  edgeNodeMap[pair<int, int>(1, 4)] = 5 + 3 * n;
82  edgeNodeMap[pair<int, int>(2, 4)] = 5 + 4 * n;
83  edgeNodeMap[pair<int, int>(3, 4)] = 5 + 5 * n;
84 
85  // Add vertices
86  for (int i = 0; i < 4; ++i)
87  {
88  m_vertex.push_back(pNodeList[i]);
89  }
90 
91  // Create edges (with corresponding set of edge points)
92  int eid = 0;
93  for (it = edgeNodeMap.begin(); it != edgeNodeMap.end(); ++it)
94  {
95  vector<NodeSharedPtr> edgeNodes;
96  if (m_conf.m_order > 1)
97  {
98  for (int j = it->second; j < it->second + n; ++j)
99  {
100  edgeNodes.push_back(pNodeList[j - 1]);
101  }
102  }
103  m_edge.push_back(EdgeSharedPtr(new Edge(pNodeList[it->first.first - 1],
104  pNodeList[it->first.second - 1],
105  edgeNodes,
107  m_edge.back()->m_id = eid++;
108  }
109 
110  // Reorient the tet to ensure collapsed coordinates align between
111  // adjacent elements.
112  if (m_conf.m_reorient)
113  {
114  OrientTet();
115  }
116  else
117  {
118  // If we didn't need to orient the tet then set up the
119  // orientation map as the identity mapping.
120  for (int i = 0; i < 4; ++i)
121  {
122  m_orientationMap[i] = i;
123  }
124  }
125 
126  // Face-edge IDs
127  int face_edges[4][3];
128 
129  // Create faces
130  for (int j = 0; j < 4; ++j)
131  {
132  vector<NodeSharedPtr> faceVertices;
133  vector<EdgeSharedPtr> faceEdges;
134  vector<NodeSharedPtr> faceNodes;
135 
136  // Extract the edges for this face. We need to do this because
137  // of the reorientation which might have been applied (see the
138  // additional note below).
139  for (int k = 0; k < 3; ++k)
140  {
141  faceVertices.push_back(m_vertex[m_faceIds[j][k]]);
142  NodeSharedPtr a = m_vertex[m_faceIds[j][k]];
143  NodeSharedPtr b = m_vertex[m_faceIds[j][(k + 1) % 3]];
144  for (unsigned int i = 0; i < m_edge.size(); ++i)
145  {
146  if (((*(m_edge[i]->m_n1) == *a) &&
147  (*(m_edge[i]->m_n2) == *b)) ||
148  ((*(m_edge[i]->m_n1) == *b) && (*(m_edge[i]->m_n2) == *a)))
149  {
150  face_edges[j][k] = i;
151  faceEdges.push_back(m_edge[i]);
152  break;
153  }
154  }
155  }
156 
157  // When face curvature is supplied, it may have been the case
158  // that our tetrahedron was reoriented. In this case, faces have
159  // different vertex IDs and so we have to rotate the face
160  // curvature so that the two align appropriately.
161  if (m_conf.m_faceNodes)
162  {
163  const int nFaceNodes = n * (n - 1) / 2;
164 
165  // Get the vertex IDs of whatever face we are processing.
166  vector<int> faceVertIds(3);
167  faceVertIds[0] = faceVertices[0]->m_id;
168  faceVertIds[1] = faceVertices[1]->m_id;
169  faceVertIds[2] = faceVertices[2]->m_id;
170 
171  // Find out the original face number as we were given it in
172  // the constructor using the orientation map.
173  int origFace = -1;
174  for (int i = 0; i < 4; ++i)
175  {
176  if (m_orientationMap[i] == j)
177  {
178  origFace = i;
179  break;
180  }
181  }
182 
183  ASSERTL0(origFace >= 0, "Couldn't find face");
184 
185  // Now get the face nodes for the original face.
186  int N = 4 + 6 * n + origFace * nFaceNodes;
187  for (int i = 0; i < nFaceNodes; ++i)
188  {
189  faceNodes.push_back(pNodeList[N + i]);
190  }
191 
192  // Find the original face vertex IDs.
193  vector<int> origFaceIds(3);
194  origFaceIds[0] = pNodeList[m_faceIds[origFace][0]]->m_id;
195  origFaceIds[1] = pNodeList[m_faceIds[origFace][1]]->m_id;
196  origFaceIds[2] = pNodeList[m_faceIds[origFace][2]]->m_id;
197 
198  // Construct a HOTriangle object which performs the
199  // orientation magically for us.
200  HOTriangle<NodeSharedPtr> hoTri(origFaceIds, faceNodes);
201  hoTri.Align(faceVertIds);
202 
203  // Copy the face nodes back again.
204  faceNodes = hoTri.surfVerts;
205  }
206 
207  m_face.push_back(FaceSharedPtr(new Face(
208  faceVertices, faceNodes, faceEdges, m_conf.m_faceCurveType)));
209  }
210 
211  if (m_conf.m_volumeNodes)
212  {
213  const int nFaceNodes = n * (n - 1) / 2;
214  for (int i = 4 + 6 * n + 4 * nFaceNodes; i < pNodeList.size(); ++i)
215  {
216  m_volumeNodes.push_back(pNodeList[i]);
217  }
218  }
219 
220  vector<EdgeSharedPtr> tmp(6);
221  tmp[0] = m_edge[face_edges[0][0]];
222  tmp[1] = m_edge[face_edges[0][1]];
223  tmp[2] = m_edge[face_edges[0][2]];
224  tmp[3] = m_edge[face_edges[1][2]];
225  tmp[4] = m_edge[face_edges[1][1]];
226  tmp[5] = m_edge[face_edges[2][1]];
227  m_edge = tmp;
228 }
bool m_faceNodes
Denotes whether the element contains face nodes. For 2D elements, if this is true then the element co...
Definition: ElementConfig.h:80
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
LibUtilities::PointsType m_faceCurveType
Distribution of points in faces.
Definition: ElementConfig.h:94
ElmtConfig m_conf
Contains configuration of the element.
Definition: Element.h:380
std::vector< int > m_taglist
List of integers specifying properties of the element.
Definition: Element.h:384
LibUtilities::PointsType m_edgeCurveType
Distribution of points in edges.
Definition: ElementConfig.h:92
unsigned int m_order
Order of the element.
Definition: ElementConfig.h:87
std::vector< NodeSharedPtr > m_vertex
List of element vertex nodes.
Definition: Element.h:386
unsigned int m_dim
Dimension of the element.
Definition: Element.h:378
static int m_faceIds[4][3]
Vertex IDs that make up tetrahedron faces.
Definition: Tetrahedron.h:99
bool m_volumeNodes
Denotes whether the element contains volume (i.e. interior) nodes. These are not supported by either ...
Definition: ElementConfig.h:85
std::vector< EdgeSharedPtr > m_edge
List of element edges.
Definition: Element.h:388
boost::shared_ptr< Node > NodeSharedPtr
Definition: Node.h:50
std::vector< NodeSharedPtr > m_volumeNodes
List of element volume nodes.
Definition: Element.h:392
void OrientTet()
Orient tetrahedron to align degenerate vertices.
std::string m_tag
Tag character describing the element.
Definition: Element.h:382
boost::shared_ptr< Edge > EdgeSharedPtr
Shared pointer to an edge.
Definition: Edge.h:136
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
NEKMESHUTILS_EXPORT Element(ElmtConfig pConf, unsigned int pNumNodes, unsigned int pGotNodes)
Definition: Element.cpp:54
static NEKMESHUTILS_EXPORT unsigned int GetNumNodes(ElmtConfig pConf)
Return the number of nodes defining a tetrahedron.
std::vector< FaceSharedPtr > m_face
List of element faces.
Definition: Element.h:390
bool m_reorient
Denotes whether the element needs to be re-orientated for a spectral element framework.
Definition: ElementConfig.h:90
boost::shared_ptr< Face > FaceSharedPtr
Definition: Face.h:153
NEKMESHUTILS_EXPORT Nektar::NekMeshUtils::Tetrahedron::Tetrahedron ( const Tetrahedron pSrc)
virtual NEKMESHUTILS_EXPORT Nektar::NekMeshUtils::Tetrahedron::~Tetrahedron ( )
inlinevirtual

Definition at line 68 of file Tetrahedron.h.

69  {
70  }

Member Function Documentation

static ElementSharedPtr Nektar::NekMeshUtils::Tetrahedron::create ( ElmtConfig  pConf,
std::vector< NodeSharedPtr pNodeList,
std::vector< int >  pTagList 
)
inlinestatic

Creates an instance of this class.

Definition at line 54 of file Tetrahedron.h.

References Tetrahedron().

57  {
58  return boost::shared_ptr<Element>(
59  new Tetrahedron(pConf, pNodeList, pTagList));
60  }
NEKMESHUTILS_EXPORT Tetrahedron(ElmtConfig pConf, std::vector< NodeSharedPtr > pNodeList, std::vector< int > pTagList)
Create a tetrahedron element.
Definition: Tetrahedron.cpp:64
void Nektar::NekMeshUtils::Tetrahedron::GetCurvedNodes ( std::vector< NodeSharedPtr > &  nodeList) const
virtual

get list of volume interior nodes

Reimplemented from Nektar::NekMeshUtils::Element.

Definition at line 391 of file Tetrahedron.cpp.

References Nektar::NekMeshUtils::HOTriangle< T >::Align(), CellMLToNektar.pycml::copy(), Nektar::NekMeshUtils::Element::m_edge, Nektar::NekMeshUtils::Element::m_face, Nektar::NekMeshUtils::Element::m_id, Nektar::NekMeshUtils::Element::m_vertex, Nektar::NekMeshUtils::Element::m_volumeNodes, and Nektar::NekMeshUtils::HOTriangle< T >::surfVerts.

392 {
393  int n = m_edge[0]->GetNodeCount();
394  nodeList.resize(n*(n+1)*(n+2)/6);
395 
396  nodeList[0] = m_vertex[0];
397  nodeList[1] = m_vertex[1];
398  nodeList[2] = m_vertex[2];
399  nodeList[3] = m_vertex[3];
400  int k = 4;
401 
402  for(int i = 0; i < 6; i++)
403  {
404  bool reverseEdge = false;
405  if(i < 3)
406  {
407  reverseEdge = m_edge[i]->m_n1 == m_vertex[i];
408  }
409  else
410  {
411  reverseEdge = m_edge[i]->m_n1 == m_vertex[i-3];
412  }
413 
414  if (reverseEdge)
415  {
416  for(int j = 0; j < n-2; j++)
417  {
418  nodeList[k++] = m_edge[i]->m_edgeNodes[j];
419  }
420  }
421  else
422  {
423  for(int j = n-3; j >= 0; j--)
424  {
425  nodeList[k++] = m_edge[i]->m_edgeNodes[j];
426  }
427  }
428  }
429 
430  vector<vector<int> > ts;
431  vector<int> t(3);
432  t[0] = m_vertex[0]->m_id;
433  t[1] = m_vertex[1]->m_id;
434  t[2] = m_vertex[2]->m_id;
435  ts.push_back(t);
436  t[0] = m_vertex[0]->m_id;
437  t[1] = m_vertex[1]->m_id;
438  t[2] = m_vertex[3]->m_id;
439  ts.push_back(t);
440  t[0] = m_vertex[1]->m_id;
441  t[1] = m_vertex[2]->m_id;
442  t[2] = m_vertex[3]->m_id;
443  ts.push_back(t);
444  t[0] = m_vertex[0]->m_id;
445  t[1] = m_vertex[2]->m_id;
446  t[2] = m_vertex[3]->m_id;
447  ts.push_back(t);
448 
449  for(int i = 0; i < 4; i++)
450  {
451  vector<int> fcid;
452  fcid.push_back(m_face[i]->m_vertexList[0]->m_id);
453  fcid.push_back(m_face[i]->m_vertexList[1]->m_id);
454  fcid.push_back(m_face[i]->m_vertexList[2]->m_id);
455 
456  HOTriangle<NodeSharedPtr> hot(fcid, m_face[i]->m_faceNodes);
457 
458  hot.Align(ts[i]);
459 
460  std::copy(hot.surfVerts.begin(),
461  hot.surfVerts.end(),
462  nodeList.begin() + k);
463  k+= hot.surfVerts.size();
464  }
465 
466  std::copy(m_volumeNodes.begin(),
467  m_volumeNodes.end(),
468  nodeList.begin() + k);
469 }
std::vector< NodeSharedPtr > m_vertex
List of element vertex nodes.
Definition: Element.h:386
std::vector< EdgeSharedPtr > m_edge
List of element edges.
Definition: Element.h:388
std::vector< NodeSharedPtr > m_volumeNodes
List of element volume nodes.
Definition: Element.h:392
unsigned int m_id
ID of the element.
Definition: Element.h:376
std::vector< FaceSharedPtr > m_face
List of element faces.
Definition: Element.h:390
StdRegions::Orientation Nektar::NekMeshUtils::Tetrahedron::GetEdgeOrient ( int  edgeId,
EdgeSharedPtr  edge 
)
virtual

Get the edge orientation of edge with respect to the local element, which lies at edge index edgeId.

Reimplemented from Nektar::NekMeshUtils::Element.

Definition at line 246 of file Tetrahedron.cpp.

References ASSERTL1, Nektar::StdRegions::eBackwards, Nektar::StdRegions::eForwards, Nektar::StdRegions::eNoOrientation, and Nektar::NekMeshUtils::Element::m_vertex.

248 {
249  static int edgeVerts[6][2] = { {0,1}, {1,2}, {0,2}, {0,3}, {1,3}, {2,3} };
250 
251  if (edge->m_n1 == m_vertex[edgeVerts[edgeId][0]])
252  {
253  return StdRegions::eForwards;
254  }
255  else if (edge->m_n1 == m_vertex[edgeVerts[edgeId][1]])
256  {
257  return StdRegions::eBackwards;
258  }
259  else
260  {
261  ASSERTL1(false, "Edge is not connected to this quadrilateral.");
262  }
263 
265 }
std::vector< NodeSharedPtr > m_vertex
List of element vertex nodes.
Definition: Element.h:386
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
virtual NEKMESHUTILS_EXPORT int Nektar::NekMeshUtils::Tetrahedron::GetFaceVertex ( int  i,
int  j 
)
inlinevirtual

Returns the local index of vertex j of face i.

Reimplemented from Nektar::NekMeshUtils::Element.

Definition at line 87 of file Tetrahedron.h.

References m_faceIds.

88  {
89  return m_faceIds[i][j];
90  }
static int m_faceIds[4][3]
Vertex IDs that make up tetrahedron faces.
Definition: Tetrahedron.h:99
SpatialDomains::GeometrySharedPtr Nektar::NekMeshUtils::Tetrahedron::GetGeom ( int  coordDim)
virtual

Generate a Nektar++ geometry object for this element.

Reimplemented from Nektar::NekMeshUtils::Element.

Definition at line 230 of file Tetrahedron.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and Nektar::NekMeshUtils::Element::m_face.

231 {
234 
235  for (int i = 0; i < 4; ++i)
236  {
237  tfaces[i] = boost::dynamic_pointer_cast<SpatialDomains::TriGeom>(
238  m_face[i]->GetGeom(coordDim));
239  }
240 
242 
243  return ret;
244 }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
boost::shared_ptr< TetGeom > TetGeomSharedPtr
Definition: TetGeom.h:106
boost::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:58
std::vector< FaceSharedPtr > m_face
List of element faces.
Definition: Element.h:390
unsigned int Nektar::NekMeshUtils::Tetrahedron::GetNumNodes ( ElmtConfig  pConf)
static

Return the number of nodes defining a tetrahedron.

Definition at line 344 of file Tetrahedron.cpp.

References Nektar::NekMeshUtils::ElmtConfig::m_faceNodes, Nektar::NekMeshUtils::ElmtConfig::m_order, and Nektar::NekMeshUtils::ElmtConfig::m_volumeNodes.

Referenced by Nektar::Utilities::InputGmsh::GetNnodes().

345 {
346  int n = pConf.m_order;
347  if (pConf.m_volumeNodes && pConf.m_faceNodes)
348  return (n + 1) * (n + 2) * (n + 3) / 6;
349  else if (!pConf.m_volumeNodes && pConf.m_faceNodes)
350  return 4 * (n + 1) * (n + 2) / 2 - 6 * (n + 1) + 4;
351  else
352  return 6 * (n + 1) - 8;
353 }
void Nektar::NekMeshUtils::Tetrahedron::MakeOrder ( int  order,
SpatialDomains::GeometrySharedPtr  geom,
LibUtilities::PointsType  edgeType,
int  coordDim,
int &  id,
bool  justConfig = false 
)
virtual

Insert interior (i.e. volume) points into this element to make the geometry an order order representation.

Parameters
orderThe desired polynomial order.
geomThe geometry object used to describe the curvature mapping.
edgeTypeThe points distribution to use on the volume.
coordDimThe coordinate (i.e. space) dimension.
idCounter which should be incremented to supply consistent vertex IDs.
justConfigIf true, then the configuration Element::m_conf will be updated but no nodes will be generated. This is used when considering boundary elements, which just require copying of face or edge interior nodes.

Reimplemented from Nektar::NekMeshUtils::Element.

Definition at line 267 of file Tetrahedron.cpp.

References ASSERTL1, Nektar::LibUtilities::PointsKey::GetPointsDim(), Nektar::NekMeshUtils::Element::m_conf, Nektar::NekMeshUtils::Element::m_curveType, Nektar::NekMeshUtils::ElmtConfig::m_faceNodes, Nektar::NekMeshUtils::ElmtConfig::m_order, Nektar::NekMeshUtils::ElmtConfig::m_volumeNodes, Nektar::NekMeshUtils::Element::m_volumeNodes, class_topology::Node, and Nektar::LibUtilities::PointsManager().

273 {
274  m_conf.m_order = order;
275  m_curveType = pType;
276  m_volumeNodes.clear();
277 
278  if (order == 1 || order == 2)
279  {
281  return;
282  }
283  else if (order == 2)
284  {
285  m_conf.m_faceNodes = true;
286  m_conf.m_volumeNodes = false;
287  return;
288  }
289  else if (order == 3)
290  {
291  m_conf.m_volumeNodes = false;
292  return;
293  }
294 
295  m_conf.m_faceNodes = true;
296  m_conf.m_volumeNodes = true;
297 
298  if (justConfig)
299  {
300  return;
301  }
302 
303  int nPoints = order + 1;
304  StdRegions::StdExpansionSharedPtr xmap = geom->GetXmap();
305 
306  Array<OneD, NekDouble> px, py, pz;
307  LibUtilities::PointsKey pKey(nPoints, pType);
308  ASSERTL1(pKey.GetPointsDim() == 3, "Points distribution must be 3D");
309  LibUtilities::PointsManager()[pKey]->GetPoints(px, py, pz);
310 
311  Array<OneD, Array<OneD, NekDouble> > phys(coordDim);
312 
313  for (int i = 0; i < coordDim; ++i)
314  {
315  phys[i] = Array<OneD, NekDouble>(xmap->GetTotPoints());
316  xmap->BwdTrans(geom->GetCoeffs(i), phys[i]);
317  }
318 
319  const int nTetPts = nPoints * (nPoints + 1) * (nPoints + 2) / 6;
320  const int nTetIntPts = (nPoints - 4) * (nPoints - 3) * (nPoints - 2) / 6;
321  m_volumeNodes.resize(nTetIntPts);
322 
323  for (int i = nTetPts - nTetIntPts, cnt = 0; i < nTetPts; ++i, ++cnt)
324  {
325  Array<OneD, NekDouble> xp(3);
326  xp[0] = px[i];
327  xp[1] = py[i];
328  xp[2] = pz[i];
329 
330  Array<OneD, NekDouble> x(3, 0.0);
331  for (int j = 0; j < coordDim; ++j)
332  {
333  x[j] = xmap->PhysEvaluate(xp, phys[j]);
334  }
335 
336  m_volumeNodes[cnt] = boost::shared_ptr<Node>(
337  new Node(id++, x[0], x[1], x[2]));
338  }
339 }
bool m_faceNodes
Denotes whether the element contains face nodes. For 2D elements, if this is true then the element co...
Definition: ElementConfig.h:80
ElmtConfig m_conf
Contains configuration of the element.
Definition: Element.h:380
unsigned int m_order
Order of the element.
Definition: ElementConfig.h:87
bool m_volumeNodes
Denotes whether the element contains volume (i.e. interior) nodes. These are not supported by either ...
Definition: ElementConfig.h:85
PointsManagerT & PointsManager(void)
std::vector< NodeSharedPtr > m_volumeNodes
List of element volume nodes.
Definition: Element.h:392
LibUtilities::PointsType m_curveType
Volume curve type.
Definition: Element.h:394
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Nektar::NekMeshUtils::Tetrahedron::OrientTet ( )
protected

Orient tetrahedron to align degenerate vertices.

Orientation of tetrahedral elements is required so that the singular vertices of triangular faces (which occur as a part of the collapsed co-ordinate system) align. The algorithm is based on that used in T. Warburton's thesis and in the original Nektar source.

First the vertices are ordered with the highest global vertex at the top degenerate point, and the base degenerate point has second lowest ID. These vertices are swapped if the element is incorrectly oriented.

Definition at line 484 of file Tetrahedron.cpp.

References Nektar::iterator, m_faceIds, Nektar::NekMeshUtils::Element::m_id, m_orientationMap, m_origVertMap, and Nektar::NekMeshUtils::Element::m_vertex.

Referenced by Tetrahedron().

485 {
486  TetOrientSet faces;
487 
488  // Create a copy of the original vertex ordering. This is used to
489  // construct a mapping, #orientationMap, which maps the original
490  // face ordering to the new face ordering.
491  for (int i = 0; i < 4; ++i)
492  {
493  vector<int> nodes(3);
494 
495  nodes[0] = m_vertex[m_faceIds[i][0]]->m_id;
496  nodes[1] = m_vertex[m_faceIds[i][1]]->m_id;
497  nodes[2] = m_vertex[m_faceIds[i][2]]->m_id;
498 
499  sort(nodes.begin(), nodes.end());
500  struct TetOrient faceNodes(nodes, i);
501  faces.insert(faceNodes);
502  }
503 
504  // Store a copy of the original vertex ordering so we can create a
505  // permutation map later.
506  vector<NodeSharedPtr> origVert = m_vertex;
507 
508  // Order vertices with highest global vertex at top degenerate
509  // point. Place second highest global vertex at base degenerate
510  // point.
511  sort(m_vertex.begin(), m_vertex.end());
512 
513  // Calculate a.(b x c) if negative, reverse order of
514  // non-degenerate points to correctly orientate the tet.
515 
516  NekDouble ax = m_vertex[1]->m_x - m_vertex[0]->m_x;
517  NekDouble ay = m_vertex[1]->m_y - m_vertex[0]->m_y;
518  NekDouble az = m_vertex[1]->m_z - m_vertex[0]->m_z;
519  NekDouble bx = m_vertex[2]->m_x - m_vertex[0]->m_x;
520  NekDouble by = m_vertex[2]->m_y - m_vertex[0]->m_y;
521  NekDouble bz = m_vertex[2]->m_z - m_vertex[0]->m_z;
522  NekDouble cx = m_vertex[3]->m_x - m_vertex[0]->m_x;
523  NekDouble cy = m_vertex[3]->m_y - m_vertex[0]->m_y;
524  NekDouble cz = m_vertex[3]->m_z - m_vertex[0]->m_z;
525 
526  NekDouble nx = (ay * bz - az * by);
527  NekDouble ny = (az * bx - ax * bz);
528  NekDouble nz = (ax * by - ay * bx);
529  NekDouble nmag = sqrt(nx * nx + ny * ny + nz * nz);
530  nx /= nmag;
531  ny /= nmag;
532  nz /= nmag;
533 
534  // distance of top vertex from base
535  NekDouble dist = cx * nx + cy * ny + cz * nz;
536 
537  if (fabs(dist) <= 1e-4)
538  {
539  cerr << "Warning: degenerate tetrahedron, 3rd vertex is = " << dist
540  << " from face" << endl;
541  }
542 
543  if (dist < 0)
544  {
545  swap(m_vertex[0], m_vertex[1]);
546  }
547 
548  nx = (ay * cz - az * cy);
549  ny = (az * cx - ax * cz);
550  nz = (ax * cy - ay * cx);
551  nmag = sqrt(nx * nx + ny * ny + nz * nz);
552  nx /= nmag;
553  ny /= nmag;
554  nz /= nmag;
555 
556  // distance of top vertex from base
557  dist = bx * nx + by * ny + bz * nz;
558 
559  if (fabs(dist) <= 1e-4)
560  {
561  cerr << "Warning: degenerate tetrahedron, 2nd vertex is = " << dist
562  << " from face" << endl;
563  }
564 
565  nx = (by * cz - bz * cy);
566  ny = (bz * cx - bx * cz);
567  nz = (bx * cy - by * cx);
568  nmag = sqrt(nx * nx + ny * ny + nz * nz);
569  nx /= nmag;
570  ny /= nmag;
571  nz /= nmag;
572 
573  // distance of top vertex from base
574  dist = ax * nx + ay * ny + az * nz;
575 
576  if (fabs(dist) <= 1e-4)
577  {
578  cerr << "Warning: degenerate tetrahedron, 1st vertex is = " << dist
579  << " from face" << endl;
580  }
581 
583 
584  // Search for the face in the original set of face nodes. Then use
585  // this to construct the #orientationMap.
586  for (int i = 0; i < 4; ++i)
587  {
588  vector<int> nodes(3);
589 
590  nodes[0] = m_vertex[m_faceIds[i][0]]->m_id;
591  nodes[1] = m_vertex[m_faceIds[i][1]]->m_id;
592  nodes[2] = m_vertex[m_faceIds[i][2]]->m_id;
593  sort(nodes.begin(), nodes.end());
594 
595  struct TetOrient faceNodes(nodes, 0);
596 
597  it = faces.find(faceNodes);
598  m_orientationMap[it->fid] = i;
599 
600  for (int j = 0; j < 4; ++j)
601  {
602  if (m_vertex[i]->m_id == origVert[j]->m_id)
603  {
604  m_origVertMap[j] = i;
605  break;
606  }
607  }
608  }
609 }
boost::unordered_set< struct TetOrient, TetOrientHash > TetOrientSet
std::vector< NodeSharedPtr > m_vertex
List of element vertex nodes.
Definition: Element.h:386
static int m_faceIds[4][3]
Vertex IDs that make up tetrahedron faces.
Definition: Tetrahedron.h:99
double NekDouble
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
unsigned int m_id
ID of the element.
Definition: Element.h:376

Member Data Documentation

int Nektar::NekMeshUtils::Tetrahedron::m_faceIds
staticprivate
Initial value:
= {
{0, 1, 2}, {0, 1, 3}, {1, 2, 3}, {0, 2, 3}
}

Vertex IDs that make up tetrahedron faces.

Definition at line 99 of file Tetrahedron.h.

Referenced by GetFaceVertex(), OrientTet(), and Tetrahedron().

int Nektar::NekMeshUtils::Tetrahedron::m_orientationMap[4]

Definition at line 92 of file Tetrahedron.h.

Referenced by OrientTet(), Nektar::Utilities::InputNek::Process(), and Tetrahedron().

int Nektar::NekMeshUtils::Tetrahedron::m_origVertMap[4]

Definition at line 93 of file Tetrahedron.h.

Referenced by OrientTet().

LibUtilities::ShapeType Nektar::NekMeshUtils::Tetrahedron::m_type
static
Initial value:

Element type.

Definition at line 62 of file Tetrahedron.h.