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::Prism Class Reference

A 3-dimensional five-faced element (2 triangles, 3 quadrilaterals). More...

#include <Prism.h>

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

Public Member Functions

NEKMESHUTILS_EXPORT Prism (ElmtConfig pConf, std::vector< NodeSharedPtr > pNodeList, std::vector< int > pTagList)
 Create a prism element. More...
 
NEKMESHUTILS_EXPORT Prism (const Prism &pSrc)
 
virtual NEKMESHUTILS_EXPORT ~Prism ()
 
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 prism. More...
 

Public Attributes

unsigned int m_orientation
 
- 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 OrientPrism ()
 Orient prism to align degenerate vertices. More...
 

Static Private Attributes

static int m_faceIds [5][4]
 Vertex IDs that make up prism 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 five-faced element (2 triangles, 3 quadrilaterals).

Definition at line 50 of file Prism.h.

Constructor & Destructor Documentation

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

Create a prism element.

Definition at line 64 of file Prism.cpp.

References ASSERTL1, Nektar::LibUtilities::eNodalTriEvenlySpaced, Nektar::LibUtilities::ePolyEvenlySpaced, 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_orientation, Nektar::NekMeshUtils::ElmtConfig::m_reorient, Nektar::NekMeshUtils::Element::m_tag, Nektar::NekMeshUtils::Element::m_taglist, Nektar::NekMeshUtils::Element::m_vertex, and OrientPrism().

Referenced by create().

67  : Element(pConf, GetNumNodes(pConf), pNodeList.size())
68 {
69  m_tag = "R";
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. This is based on the ordering produced by
76  // gmsh.
77  map<pair<int, int>, int> edgeNodeMap;
78  map<pair<int, int>, int>::iterator it;
79 
80  // This edge-node map is based on Nektar++ ordering.
81  edgeNodeMap[pair<int, int>(1, 2)] = 7;
82  edgeNodeMap[pair<int, int>(2, 3)] = 7 + n;
83  edgeNodeMap[pair<int, int>(4, 3)] = 7 + 2 * n;
84  edgeNodeMap[pair<int, int>(1, 4)] = 7 + 3 * n;
85  edgeNodeMap[pair<int, int>(1, 5)] = 7 + 4 * n;
86  edgeNodeMap[pair<int, int>(2, 5)] = 7 + 5 * n;
87  edgeNodeMap[pair<int, int>(3, 6)] = 7 + 6 * n;
88  edgeNodeMap[pair<int, int>(4, 6)] = 7 + 7 * n;
89  edgeNodeMap[pair<int, int>(5, 6)] = 7 + 8 * n;
90 
91  // Add vertices
92  for (int i = 0; i < 6; ++i)
93  {
94  m_vertex.push_back(pNodeList[i]);
95  }
96 
97  int eid = 0;
98  // Create edges (with corresponding set of edge points)
99  for (it = edgeNodeMap.begin(); it != edgeNodeMap.end(); ++it)
100  {
101  vector<NodeSharedPtr> edgeNodes;
102  if (m_conf.m_order > 1)
103  {
104  for (int j = it->second; j < it->second + n; ++j)
105  {
106  edgeNodes.push_back(pNodeList[j - 1]);
107  }
108  }
109  m_edge.push_back(EdgeSharedPtr(new Edge(pNodeList[it->first.first - 1],
110  pNodeList[it->first.second - 1],
111  edgeNodes,
113  m_edge.back()->m_id = eid++;
114  }
115 
116  if (m_conf.m_reorient)
117  {
118  OrientPrism();
119  }
120  else
121  {
122  m_orientation = 0;
123  }
124 
125  // Create faces
126  int face_edges[5][4];
127 
128  int face_offset[5];
129  face_offset[0] = 6 + 9 * n;
130  for (int j = 0; j < 4; ++j)
131  {
132  int facenodes = j % 2 == 0 ? n * n : n * (n - 1) / 2;
133  face_offset[j + 1] = face_offset[j] + facenodes;
134  }
135 
136  for (int j = 0; j < 5; ++j)
137  {
138  vector<NodeSharedPtr> faceVertices;
139  vector<EdgeSharedPtr> faceEdges;
140  vector<NodeSharedPtr> faceNodes;
141  int nEdge = 3 - (j % 2 - 1);
142 
143  for (int k = 0; k < nEdge; ++k)
144  {
145  faceVertices.push_back(m_vertex[m_faceIds[j][k]]);
146  NodeSharedPtr a = m_vertex[m_faceIds[j][k]];
147  NodeSharedPtr b = m_vertex[m_faceIds[j][(k + 1) % nEdge]];
148  unsigned int i;
149  for (i = 0; i < m_edge.size(); ++i)
150  {
151  if ((m_edge[i]->m_n1->m_id == a->m_id &&
152  m_edge[i]->m_n2->m_id == b->m_id) ||
153  (m_edge[i]->m_n1->m_id == b->m_id &&
154  m_edge[i]->m_n2->m_id == a->m_id))
155  {
156  faceEdges.push_back(m_edge[i]);
157  face_edges[j][k] = i;
158  break;
159  }
160  }
161 
162  if (i == m_edge.size())
163  {
164  face_edges[j][k] = -1;
165  }
166  }
167 
168  if (m_conf.m_faceNodes)
169  {
170  int face = j, facenodes;
171 
172  if (j % 2 == 0)
173  {
174  facenodes = n * n;
175  if (m_orientation == 1)
176  {
177  face = (face + 4) % 6;
178  }
179  else if (m_orientation == 2)
180  {
181  face = (face + 2) % 6;
182  }
183  }
184  else
185  {
186  // TODO: need to rotate these too.
187  facenodes = n * (n - 1) / 2;
188  }
189 
190  for (int i = 0; i < facenodes; ++i)
191  {
192  faceNodes.push_back(pNodeList[face_offset[face] + i]);
193  }
194  }
195 
196  // Try to translate between common face curve types
198 
199  if (pType == LibUtilities::ePolyEvenlySpaced && (j == 1 || j == 3))
200  {
202  }
203 
204  m_face.push_back(
205  FaceSharedPtr(new Face(faceVertices, faceNodes, faceEdges, pType)));
206  }
207 
208  // Re-order edge array to be consistent with Nektar++ ordering.
209  vector<EdgeSharedPtr> tmp(9);
210  ASSERTL1(face_edges[0][0] != -1, "face_edges[0][0] == -1");
211  tmp[0] = m_edge[face_edges[0][0]];
212  ASSERTL1(face_edges[0][1] != -1, "face_edges[0][1] == -1");
213  tmp[1] = m_edge[face_edges[0][1]];
214  ASSERTL1(face_edges[0][2] != -1, "face_edges[0][2] == -1");
215  tmp[2] = m_edge[face_edges[0][2]];
216  ASSERTL1(face_edges[0][3] != -1, "face_edges[0][3] == -1");
217  tmp[3] = m_edge[face_edges[0][3]];
218  ASSERTL1(face_edges[1][2] != -1, "face_edges[1][2] == -1");
219  tmp[4] = m_edge[face_edges[1][2]];
220  ASSERTL1(face_edges[1][1] != -1, "face_edges[1][1] == -1");
221  tmp[5] = m_edge[face_edges[1][1]];
222  ASSERTL1(face_edges[2][1] != -1, "face_edges[2][1] == -1");
223  tmp[6] = m_edge[face_edges[2][1]];
224  ASSERTL1(face_edges[3][2] != -1, "face_edges[3][2] == -1");
225  tmp[7] = m_edge[face_edges[3][2]];
226  ASSERTL1(face_edges[4][2] != -1, "face_edges[4][2] == -1");
227  tmp[8] = m_edge[face_edges[4][2]];
228  m_edge = tmp;
229 }
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
LibUtilities::PointsType m_faceCurveType
Distribution of points in faces.
Definition: ElementConfig.h:94
unsigned int m_orientation
Definition: Prism.h:96
static int m_faceIds[5][4]
Vertex IDs that make up prism faces.
Definition: Prism.h:102
ElmtConfig m_conf
Contains configuration of the element.
Definition: Element.h:380
void OrientPrism()
Orient prism to align degenerate vertices.
Definition: Prism.cpp:522
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
1D Evenly-spaced points using Lagrange polynomial
Definition: PointsType.h:65
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 NEKMESHUTILS_EXPORT unsigned int GetNumNodes(ElmtConfig pConf)
Return the number of nodes defining a prism.
Definition: Prism.cpp:234
std::vector< EdgeSharedPtr > m_edge
List of element edges.
Definition: Element.h:388
boost::shared_ptr< Node > NodeSharedPtr
Definition: Node.h:50
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
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
2D Evenly-spaced points on a Triangle
Definition: PointsType.h:72
boost::shared_ptr< Face > FaceSharedPtr
Definition: Face.h:153
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
NEKMESHUTILS_EXPORT Nektar::NekMeshUtils::Prism::Prism ( const Prism pSrc)
virtual NEKMESHUTILS_EXPORT Nektar::NekMeshUtils::Prism::~Prism ( )
inlinevirtual

Definition at line 68 of file Prism.h.

69  {
70  }

Member Function Documentation

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

Creates an instance of this class.

Definition at line 54 of file Prism.h.

References Prism().

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

get list of volume interior nodes

Reimplemented from Nektar::NekMeshUtils::Element.

Definition at line 353 of file Prism.cpp.

References Nektar::NekMeshUtils::HOTriangle< T >::Align(), Nektar::NekMeshUtils::HOQuadrilateral< 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, Nektar::NekMeshUtils::HOTriangle< T >::surfVerts, and Nektar::NekMeshUtils::HOQuadrilateral< T >::surfVerts.

354 {
355  int n = m_edge[0]->GetNodeCount();
356  nodeList.resize(n*n*(n+1)/2);
357 
358  nodeList[0] = m_vertex[0];
359  nodeList[1] = m_vertex[1];
360  nodeList[2] = m_vertex[2];
361  nodeList[3] = m_vertex[3];
362  nodeList[4] = m_vertex[4];
363  nodeList[5] = m_vertex[5];
364  int k = 6;
365 
366  for(int i = 0; i < 4; i++)
367  {
368  bool reverseEdge = m_edge[i]->m_n1 == m_vertex[i];
369  if (reverseEdge)
370  {
371  for(int j = 0; j < n-2; j++)
372  {
373  nodeList[k++] = m_edge[i]->m_edgeNodes[j];
374  }
375  }
376  else
377  {
378  for(int j = n-3; j >= 0; j--)
379  {
380  nodeList[k++] = m_edge[i]->m_edgeNodes[j];
381  }
382  }
383  }
384 
385  for(int i = 4; i < 8; i++)
386  {
387  bool reverseEdge = m_edge[i]->m_n1 == m_vertex[i-4];
388  if (reverseEdge)
389  {
390  for(int j = 0; j < n-2; j++)
391  {
392  nodeList[k++] = m_edge[i]->m_edgeNodes[j];
393  }
394  }
395  else
396  {
397  for(int j = n-3; j >= 0; j--)
398  {
399  nodeList[k++] = m_edge[i]->m_edgeNodes[j];
400  }
401  }
402  }
403  bool reverseEdge = m_edge[8]->m_n1 == m_vertex[4];
404  if (reverseEdge)
405  {
406  for(int j = 0; j < n-2; j++)
407  {
408  nodeList[k++] = m_edge[8]->m_edgeNodes[j];
409  }
410  }
411  else
412  {
413  for(int j = n-3; j >= 0; j--)
414  {
415  nodeList[k++] = m_edge[8]->m_edgeNodes[j];
416  }
417  }
418 
419  vector<vector<int> > ts;
420  {
421  vector<int> t(4);
422  t[0] = m_vertex[0]->m_id;
423  t[1] = m_vertex[1]->m_id;
424  t[2] = m_vertex[2]->m_id;
425  t[3] = m_vertex[3]->m_id;
426  ts.push_back(t);
427  }
428  {
429  vector<int> t(3);
430  t[0] = m_vertex[0]->m_id;
431  t[1] = m_vertex[1]->m_id;
432  t[2] = m_vertex[4]->m_id;
433  ts.push_back(t);
434  }
435  {
436  vector<int> t(4);
437  t[0] = m_vertex[1]->m_id;
438  t[1] = m_vertex[2]->m_id;
439  t[2] = m_vertex[5]->m_id;
440  t[3] = m_vertex[4]->m_id;
441  ts.push_back(t);
442  }
443  {
444  vector<int> t(3);
445  t[0] = m_vertex[3]->m_id;
446  t[1] = m_vertex[2]->m_id;
447  t[2] = m_vertex[5]->m_id;
448  ts.push_back(t);
449  }
450  {
451  vector<int> t(4);
452  t[0] = m_vertex[0]->m_id;
453  t[1] = m_vertex[3]->m_id;
454  t[2] = m_vertex[5]->m_id;
455  t[3] = m_vertex[4]->m_id;
456  ts.push_back(t);
457  }
458 
459  for(int i = 0; i < ts.size(); i++)
460  {
461  if(ts[i].size() == 3)
462  {
463  vector<int> fcid;
464  fcid.push_back(m_face[i]->m_vertexList[0]->m_id);
465  fcid.push_back(m_face[i]->m_vertexList[1]->m_id);
466  fcid.push_back(m_face[i]->m_vertexList[2]->m_id);
467 
468  HOTriangle<NodeSharedPtr> hot(fcid, m_face[i]->m_faceNodes);
469 
470  hot.Align(ts[i]);
471 
472  std::copy(hot.surfVerts.begin(),
473  hot.surfVerts.end(),
474  nodeList.begin() + k);
475  k+= hot.surfVerts.size();
476  }
477  else
478  {
479  vector<int> fcid;
480  fcid.push_back(m_face[i]->m_vertexList[0]->m_id);
481  fcid.push_back(m_face[i]->m_vertexList[1]->m_id);
482  fcid.push_back(m_face[i]->m_vertexList[2]->m_id);
483  fcid.push_back(m_face[i]->m_vertexList[3]->m_id);
484 
485  HOQuadrilateral<NodeSharedPtr> hoq(fcid, m_face[i]->m_faceNodes);
486 
487  hoq.Align(ts[i]);
488 
489  std::copy(hoq.surfVerts.begin(),
490  hoq.surfVerts.end(),
491  nodeList.begin() + k);
492  k+= hoq.surfVerts.size();
493  }
494  }
495 
496  std::copy(m_volumeNodes.begin(),
497  m_volumeNodes.end(),
498  nodeList.begin() + k);
499 }
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::Prism::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 261 of file Prism.cpp.

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

263 {
264  static int edgeVerts[9][2] = {
265  {0,1}, {1,2}, {3,2}, {0,3}, {0,4}, {1,4}, {2,5}, {3,5}, {4,5}
266  };
267 
268  if (edge->m_n1 == m_vertex[edgeVerts[edgeId][0]])
269  {
270  return StdRegions::eForwards;
271  }
272  else if (edge->m_n1 == m_vertex[edgeVerts[edgeId][1]])
273  {
274  return StdRegions::eBackwards;
275  }
276  else
277  {
278  ASSERTL1(false, "Edge is not connected to this quadrilateral.");
279  }
280 
282 }
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::Prism::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 Prism.h.

References m_faceIds.

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

Generate a Nektar++ geometry object for this element.

Reimplemented from Nektar::NekMeshUtils::Element.

Definition at line 246 of file Prism.cpp.

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

247 {
250 
251  for (int i = 0; i < 5; ++i)
252  {
253  faces[i] = m_face[i]->GetGeom(coordDim);
254  }
255 
257 
258  return ret;
259 }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
boost::shared_ptr< Geometry2D > Geometry2DSharedPtr
Definition: Geometry2D.h:59
boost::shared_ptr< PrismGeom > PrismGeomSharedPtr
Definition: PrismGeom.h:109
std::vector< FaceSharedPtr > m_face
List of element faces.
Definition: Element.h:390
unsigned int Nektar::NekMeshUtils::Prism::GetNumNodes ( ElmtConfig  pConf)
static

Return the number of nodes defining a prism.

Definition at line 234 of file Prism.cpp.

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

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

235 {
236  int n = pConf.m_order;
237  if (pConf.m_faceNodes && pConf.m_volumeNodes)
238  return (n + 1) * (n + 1) * (n + 2) / 2;
239  else if (pConf.m_faceNodes && !pConf.m_volumeNodes)
240  return 3 * (n + 1) * (n + 1) + 2 * (n + 1) * (n + 2) / 2 - 9 * (n + 1) +
241  6;
242  else
243  return 9 * (n + 1) - 12;
244 }
void Nektar::NekMeshUtils::Prism::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 284 of file Prism.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().

290 {
291  m_conf.m_order = order;
292  m_curveType = pType;
293  m_volumeNodes.clear();
294 
295  if (order == 1)
296  {
298  return;
299  }
300  else if (order == 2)
301  {
302  m_conf.m_faceNodes = true;
303  m_conf.m_volumeNodes = false;
304  return;
305  }
306 
307  m_conf.m_faceNodes = true;
308  m_conf.m_volumeNodes = true;
309 
310  if (justConfig)
311  {
312  return;
313  }
314 
315  int nPoints = order + 1;
316  StdRegions::StdExpansionSharedPtr xmap = geom->GetXmap();
317 
318  Array<OneD, NekDouble> px, py, pz;
319  LibUtilities::PointsKey pKey(nPoints, pType);
320  ASSERTL1(pKey.GetPointsDim() == 3, "Points distribution must be 3D");
321  LibUtilities::PointsManager()[pKey]->GetPoints(px, py, pz);
322 
323  Array<OneD, Array<OneD, NekDouble> > phys(coordDim);
324 
325  for (int i = 0; i < coordDim; ++i)
326  {
327  phys[i] = Array<OneD, NekDouble>(xmap->GetTotPoints());
328  xmap->BwdTrans(geom->GetCoeffs(i), phys[i]);
329  }
330 
331  const int nPrismPts = nPoints * nPoints * (nPoints + 1) / 2;
332  const int nPrismIntPts = (nPoints - 2) * (nPoints - 3) * (nPoints - 2) / 2;
333  m_volumeNodes.resize(nPrismIntPts);
334 
335  for (int i = nPrismPts - nPrismIntPts, cnt = 0; i < nPrismPts; ++i, ++cnt)
336  {
337  Array<OneD, NekDouble> xp(3);
338  xp[0] = px[i];
339  xp[1] = py[i];
340  xp[2] = pz[i];
341 
342  Array<OneD, NekDouble> x(3, 0.0);
343  for (int j = 0; j < coordDim; ++j)
344  {
345  x[j] = xmap->PhysEvaluate(xp, phys[j]);
346  }
347 
348  m_volumeNodes[cnt] = boost::shared_ptr<Node>(
349  new Node(id++, x[0], x[1], x[2]));
350  }
351 }
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::Prism::OrientPrism ( )
protected

Orient prism to align degenerate vertices.

Orientation of prismatric 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 points are re-ordered so that the highest global IDs represent the two singular points of the prism. Then, if necessary, the nodes are rotated either clockwise or counter-clockwise (w.r.t to the p-r plane) to correctly align the prism. The #orientation variable is set to:

  • 0 if the prism is not rotated;
  • 1 if the prism is rotated clockwise;
  • 2 if the prism is rotated counter-clockwise.

This is necessary for some input modules (e.g. #InputNek) which add high-order information or bounary conditions to faces.

Definition at line 522 of file Prism.cpp.

References m_orientation, and Nektar::NekMeshUtils::Element::m_vertex.

Referenced by Prism().

523 {
524  int lid[6], gid[6];
525 
526  // Re-order vertices.
527  for (int i = 0; i < 6; ++i)
528  {
529  lid[i] = i;
530  gid[i] = m_vertex[i]->m_id;
531  }
532 
533  gid[0] = gid[3] = max(gid[0], gid[3]);
534  gid[1] = gid[2] = max(gid[1], gid[2]);
535  gid[4] = gid[5] = max(gid[4], gid[5]);
536 
537  for (int i = 1; i < 6; ++i)
538  {
539  if (gid[0] < gid[i])
540  {
541  swap(gid[i], gid[0]);
542  swap(lid[i], lid[0]);
543  }
544  }
545 
546  if (lid[0] == 4 || lid[0] == 5)
547  {
548  m_orientation = 0;
549  }
550  else if (lid[0] == 1 || lid[0] == 2)
551  {
552  // Rotate prism clockwise in p-r plane
553  vector<NodeSharedPtr> vertexmap(6);
554  vertexmap[0] = m_vertex[4];
555  vertexmap[1] = m_vertex[0];
556  vertexmap[2] = m_vertex[3];
557  vertexmap[3] = m_vertex[5];
558  vertexmap[4] = m_vertex[1];
559  vertexmap[5] = m_vertex[2];
560  m_vertex = vertexmap;
561  m_orientation = 1;
562  }
563  else if (lid[0] == 0 || lid[0] == 3)
564  {
565  // Rotate prism counter-clockwise in p-r plane
566  vector<NodeSharedPtr> vertexmap(6);
567  vertexmap[0] = m_vertex[1];
568  vertexmap[1] = m_vertex[4];
569  vertexmap[2] = m_vertex[5];
570  vertexmap[3] = m_vertex[2];
571  vertexmap[4] = m_vertex[0];
572  vertexmap[5] = m_vertex[3];
573  m_vertex = vertexmap;
574  m_orientation = 2;
575  }
576  else
577  {
578  cerr << "Warning: possible prism orientation problem." << endl;
579  }
580 }
unsigned int m_orientation
Definition: Prism.h:96
std::vector< NodeSharedPtr > m_vertex
List of element vertex nodes.
Definition: Element.h:386

Member Data Documentation

int Nektar::NekMeshUtils::Prism::m_faceIds
staticprivate
Initial value:
= {
{0, 1, 2, 3}, {0, 1, 4, -1}, {1, 2, 5, 4}, {3, 2, 5, -1}, {0, 3, 5, 4}
}

Vertex IDs that make up prism faces.

Definition at line 102 of file Prism.h.

Referenced by GetFaceVertex(), and Prism().

unsigned int Nektar::NekMeshUtils::Prism::m_orientation

Orientation of prism; unchanged = 0; clockwise = 1; counter-clockwise = 2. This is set by OrientPrism.

Definition at line 96 of file Prism.h.

Referenced by OrientPrism(), and Prism().

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

Element type.

Definition at line 62 of file Prism.h.