Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
TetGeom.cpp
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // File: TetGeom.cpp
4 //
5 // For more information, please see: http://www.nektar.info/
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // Permission is hereby granted, free of charge, to any person obtaining a
14 // copy of this software and associated documentation files (the "Software"),
15 // to deal in the Software without restriction, including without limitation
16 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
17 // and/or sell copies of the Software, and to permit persons to whom the
18 // Software is furnished to do so, subject to the following conditions:
19 //
20 // The above copyright notice and this permission notice shall be included
21 // in all copies or substantial portions of the Software.
22 //
23 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29 // DEALINGS IN THE SOFTWARE.
30 //
31 // Description: Tetrahedral geometry information.
32 //
33 ////////////////////////////////////////////////////////////////////////////////
34 
35 #include <SpatialDomains/TetGeom.h>
36 
38 #include <StdRegions/StdTetExp.h>
39 #include <SpatialDomains/SegGeom.h>
40 
41 using namespace std;
42 
43 namespace Nektar
44 {
45 namespace SpatialDomains
46 {
47 const unsigned int TetGeom::VertexEdgeConnectivity[4][3] = {
48  {0, 2, 3}, {0, 1, 4}, {1, 2, 5}, {3, 4, 5}};
49 const unsigned int TetGeom::VertexFaceConnectivity[4][3] = {
50  {0, 1, 3}, {0, 1, 2}, {0, 2, 3}, {1, 2, 3}};
51 const unsigned int TetGeom::EdgeFaceConnectivity[6][2] = {
52  {0, 1}, {0, 2}, {0, 3}, {1, 3}, {1, 2}, {2, 3}};
53 
54 TetGeom::TetGeom()
55 {
56  m_shapeType = LibUtilities::eTetrahedron;
57 }
58 
59 TetGeom::TetGeom(int id, const TriGeomSharedPtr faces[])
60  : Geometry3D(faces[0]->GetEdge(0)->GetVertex(0)->GetCoordim())
61 {
63  m_globalID = id;
64 
65  /// Copy the face shared pointers
66  m_faces.insert(m_faces.begin(), faces, faces + TetGeom::kNfaces);
67 
68  /// Set up orientation vectors with correct amount of elements.
69  m_eorient.resize(kNedges);
70  m_forient.resize(kNfaces);
71 
76 }
77 
79 {
80 }
81 
82 int TetGeom::v_GetDir(const int faceidx, const int facedir) const
83 {
84  if (faceidx == 0)
85  {
86  return facedir;
87  }
88  else if (faceidx == 1)
89  {
90  return 2 * facedir;
91  }
92  else
93  {
94  return 1 + facedir;
95  }
96 }
97 
98 int TetGeom::v_GetVertexEdgeMap(const int i, const int j) const
99 {
100  return VertexEdgeConnectivity[i][j];
101 }
102 
103 int TetGeom::v_GetVertexFaceMap(const int i, const int j) const
104 {
105  return VertexFaceConnectivity[i][j];
106 }
107 
108 int TetGeom::v_GetEdgeFaceMap(const int i, const int j) const
109 {
110  return EdgeFaceConnectivity[i][j];
111 }
112 
114 {
115 
116  // find edge 0
117  int i, j;
118  unsigned int check;
119 
120  SegGeomSharedPtr edge;
121 
122  // First set up the 3 bottom edges
123 
124  if (m_faces[0]->GetEid(0) != m_faces[1]->GetEid(0))
125  {
126  std::ostringstream errstrm;
127  errstrm << "Local edge 0 (eid=" << m_faces[0]->GetEid(0);
128  errstrm << ") on face " << m_faces[0]->GetGlobalID();
129  errstrm << " must be the same as local edge 0 (eid="
130  << m_faces[1]->GetEid(0);
131  errstrm << ") on face " << m_faces[1]->GetGlobalID();
132  ASSERTL0(false, errstrm.str());
133  }
134 
135  int faceConnected;
136  for (faceConnected = 1; faceConnected < 4; faceConnected++)
137  {
138  check = 0;
139  for (i = 0; i < 3; i++)
140  {
141  if ((m_faces[0])->GetEid(i) == (m_faces[faceConnected])->GetEid(0))
142  {
143  edge = dynamic_pointer_cast<SegGeom>(
144  (m_faces[0])->GetEdge(i));
145  m_edges.push_back(edge);
146  check++;
147  }
148  }
149 
150  if (check < 1)
151  {
152  std::ostringstream errstrm;
153  errstrm << "Face 0 does not share an edge with first edge of "
154  "adjacent face. Faces ";
155  errstrm << (m_faces[0])->GetGlobalID() << ", "
156  << (m_faces[faceConnected])->GetGlobalID();
157  ASSERTL0(false, errstrm.str());
158  }
159  else if (check > 1)
160  {
161  std::ostringstream errstrm;
162  errstrm << "Connected faces share more than one edge. Faces ";
163  errstrm << (m_faces[0])->GetGlobalID() << ", "
164  << (m_faces[faceConnected])->GetGlobalID();
165  ASSERTL0(false, errstrm.str());
166  }
167  }
168 
169  // Then, set up the 3 vertical edges
170  check = 0;
171  for (i = 0; i < 3; i++) // Set up the vertical edge :face(1) and face(3)
172  {
173  for (j = 0; j < 3; j++)
174  {
175  if ((m_faces[1])->GetEid(i) == (m_faces[3])->GetEid(j))
176  {
177  edge = dynamic_pointer_cast<SegGeom>(
178  (m_faces[1])->GetEdge(i));
179  m_edges.push_back(edge);
180  check++;
181  }
182  }
183  }
184  if (check < 1)
185  {
186  std::ostringstream errstrm;
187  errstrm << "Connected faces do not share an edge. Faces ";
188  errstrm << (m_faces[1])->GetGlobalID() << ", "
189  << (m_faces[3])->GetGlobalID();
190  ASSERTL0(false, errstrm.str());
191  }
192  else if (check > 1)
193  {
194  std::ostringstream errstrm;
195  errstrm << "Connected faces share more than one edge. Faces ";
196  errstrm << (m_faces[1])->GetGlobalID() << ", "
197  << (m_faces[3])->GetGlobalID();
198  ASSERTL0(false, errstrm.str());
199  }
200  // Set up vertical edges: face(1) through face(3)
201  for (faceConnected = 1; faceConnected < 3; faceConnected++)
202  {
203  check = 0;
204  for (i = 0; i < 3; i++)
205  {
206  for (j = 0; j < 3; j++)
207  {
208  if ((m_faces[faceConnected])->GetEid(i) ==
209  (m_faces[faceConnected + 1])->GetEid(j))
210  {
211  edge = dynamic_pointer_cast<SegGeom>(
212  (m_faces[faceConnected])->GetEdge(i));
213  m_edges.push_back(edge);
214  check++;
215  }
216  }
217  }
218 
219  if (check < 1)
220  {
221  std::ostringstream errstrm;
222  errstrm << "Connected faces do not share an edge. Faces ";
223  errstrm << (m_faces[faceConnected])->GetGlobalID() << ", "
224  << (m_faces[faceConnected + 1])->GetGlobalID();
225  ASSERTL0(false, errstrm.str());
226  }
227  else if (check > 1)
228  {
229  std::ostringstream errstrm;
230  errstrm << "Connected faces share more than one edge. Faces ";
231  errstrm << (m_faces[faceConnected])->GetGlobalID() << ", "
232  << (m_faces[faceConnected + 1])->GetGlobalID();
233  ASSERTL0(false, errstrm.str());
234  }
235  }
236 }
237 
239 {
240 
241  // Set up the first 2 vertices (i.e. vertex 0,1)
242  if ((m_edges[0]->GetVid(0) == m_edges[1]->GetVid(0)) ||
243  (m_edges[0]->GetVid(0) == m_edges[1]->GetVid(1)))
244  {
245  m_verts.push_back(m_edges[0]->GetVertex(1));
246  m_verts.push_back(m_edges[0]->GetVertex(0));
247  }
248  else if ((m_edges[0]->GetVid(1) == m_edges[1]->GetVid(0)) ||
249  (m_edges[0]->GetVid(1) == m_edges[1]->GetVid(1)))
250  {
251  m_verts.push_back(m_edges[0]->GetVertex(0));
252  m_verts.push_back(m_edges[0]->GetVertex(1));
253  }
254  else
255  {
256  std::ostringstream errstrm;
257  errstrm << "Connected edges do not share a vertex. Edges ";
258  errstrm << m_edges[0]->GetGlobalID() << ", "
259  << m_edges[1]->GetGlobalID();
260  ASSERTL0(false, errstrm.str());
261  }
262 
263  // set up the other bottom vertices (i.e. vertex 2)
264  for (int i = 1; i < 2; i++)
265  {
266  if (m_edges[i]->GetVid(0) == m_verts[i]->GetGlobalID())
267  {
268  m_verts.push_back(m_edges[i]->GetVertex(1));
269  }
270  else if (m_edges[i]->GetVid(1) == m_verts[i]->GetGlobalID())
271  {
272  m_verts.push_back(m_edges[i]->GetVertex(0));
273  }
274  else
275  {
276  std::ostringstream errstrm;
277  errstrm << "Connected edges do not share a vertex. Edges ";
278  errstrm << m_edges[i]->GetGlobalID() << ", "
279  << m_edges[i - 1]->GetGlobalID();
280  ASSERTL0(false, errstrm.str());
281  }
282  }
283 
284  // set up top vertex
285  if (m_edges[3]->GetVid(0) == m_verts[0]->GetGlobalID())
286  {
287  m_verts.push_back(m_edges[3]->GetVertex(1));
288  }
289  else
290  {
291  m_verts.push_back(m_edges[3]->GetVertex(0));
292  }
293 
294  // Check the other edges match up.
295  int check = 0;
296  for (int i = 4; i < 6; ++i)
297  {
298  if ((m_edges[i]->GetVid(0) == m_verts[i - 3]->GetGlobalID() &&
299  m_edges[i]->GetVid(1) == m_verts[3]->GetGlobalID()) ||
300  (m_edges[i]->GetVid(1) == m_verts[i - 3]->GetGlobalID() &&
301  m_edges[i]->GetVid(0) == m_verts[3]->GetGlobalID()))
302  {
303  check++;
304  }
305  }
306  if (check != 2)
307  {
308  std::ostringstream errstrm;
309  errstrm << "Connected edges do not share a vertex. Edges ";
310  errstrm << m_edges[3]->GetGlobalID() << ", "
311  << m_edges[2]->GetGlobalID();
312  ASSERTL0(false, errstrm.str());
313  }
314 }
315 
317 {
318 
319  // This 2D array holds the local id's of all the vertices
320  // for every edge. For every edge, they are ordered to what we
321  // define as being Forwards
322  const unsigned int edgeVerts[kNedges][2] = {
323  {0, 1}, {1, 2}, {0, 2}, {0, 3}, {1, 3}, {2, 3}};
324 
325  int i;
326  for (i = 0; i < kNedges; i++)
327  {
328  if (m_edges[i]->GetVid(0) == m_verts[edgeVerts[i][0]]->GetGlobalID())
329  {
331  }
332  else if (m_edges[i]->GetVid(0) == m_verts[edgeVerts[i][1]]->GetGlobalID())
333  {
335  }
336  else
337  {
338  ASSERTL0(false, "Could not find matching vertex for the edge");
339  }
340  }
341 }
342 
344 {
345 
346  int f, i;
347 
348  // These arrays represent the vector of the A and B
349  // coordinate of the local elemental coordinate system
350  // where A corresponds with the coordinate direction xi_i
351  // with the lowest index i (for that particular face)
352  // Coordinate 'B' then corresponds to the other local
353  // coordinate (i.e. with the highest index)
354  Array<OneD, NekDouble> elementAaxis(m_coordim);
355  Array<OneD, NekDouble> elementBaxis(m_coordim);
356 
357  // These arrays correspond to the local coordinate
358  // system of the face itself (i.e. the Geometry2D)
359  // faceAaxis correspond to the xi_0 axis
360  // faceBaxis correspond to the xi_1 axis
363 
364  // This is the base vertex of the face (i.e. the Geometry2D)
365  // This corresponds to thevertex with local ID 0 of the
366  // Geometry2D
367  unsigned int baseVertex;
368 
369  // The lenght of the vectors above
370  NekDouble elementAaxis_length;
371  NekDouble elementBaxis_length;
372  NekDouble faceAaxis_length;
373  NekDouble faceBaxis_length;
374 
375  // This 2D array holds the local id's of all the vertices
376  // for every face. For every face, they are ordered in such
377  // a way that the implementation below allows a unified approach
378  // for all faces.
379  const unsigned int faceVerts[kNfaces][TriGeom::kNverts] = {
380  {0, 1, 2}, {0, 1, 3}, {1, 2, 3}, {0, 2, 3}};
381 
382  NekDouble dotproduct1 = 0.0;
383  NekDouble dotproduct2 = 0.0;
384 
385  unsigned int orientation;
386 
387  // Loop over all the faces to set up the orientation
388  for (f = 0; f < kNqfaces + kNtfaces; f++)
389  {
390  // initialisation
391  elementAaxis_length = 0.0;
392  elementBaxis_length = 0.0;
393  faceAaxis_length = 0.0;
394  faceBaxis_length = 0.0;
395 
396  dotproduct1 = 0.0;
397  dotproduct2 = 0.0;
398 
399  baseVertex = m_faces[f]->GetVid(0);
400 
401  // We are going to construct the vectors representing the A and B axis
402  // of every face. These vectors will be constructed as a
403  // vector-representation
404  // of the edges of the face. However, for both coordinate directions, we
405  // can
406  // represent the vectors by two different edges. That's why we need to
407  // make sure that
408  // we pick the edge to which the baseVertex of the
409  // Geometry2D-representation of the face
410  // belongs...
411 
412  // Compute the length of edges on a base-face
413  if (baseVertex == m_verts[faceVerts[f][0]]->GetGlobalID())
414  {
415  for (i = 0; i < m_coordim; i++)
416  {
417  elementAaxis[i] = (*m_verts[faceVerts[f][1]])[i] -
418  (*m_verts[faceVerts[f][0]])[i];
419  elementBaxis[i] = (*m_verts[faceVerts[f][2]])[i] -
420  (*m_verts[faceVerts[f][0]])[i];
421  }
422  }
423  else if (baseVertex == m_verts[faceVerts[f][1]]->GetGlobalID())
424  {
425  for (i = 0; i < m_coordim; i++)
426  {
427  elementAaxis[i] = (*m_verts[faceVerts[f][1]])[i] -
428  (*m_verts[faceVerts[f][0]])[i];
429  elementBaxis[i] = (*m_verts[faceVerts[f][2]])[i] -
430  (*m_verts[faceVerts[f][1]])[i];
431  }
432  }
433  else if (baseVertex == m_verts[faceVerts[f][2]]->GetGlobalID())
434  {
435  for (i = 0; i < m_coordim; i++)
436  {
437  elementAaxis[i] = (*m_verts[faceVerts[f][1]])[i] -
438  (*m_verts[faceVerts[f][2]])[i];
439  elementBaxis[i] = (*m_verts[faceVerts[f][2]])[i] -
440  (*m_verts[faceVerts[f][0]])[i];
441  }
442  }
443  else
444  {
445  ASSERTL0(false, "Could not find matching vertex for the face");
446  }
447 
448  // Now, construct the edge-vectors of the local coordinates of
449  // the Geometry2D-representation of the face
450  for (i = 0; i < m_coordim; i++)
451  {
452  faceAaxis[i] =
453  (*m_faces[f]->GetVertex(1))[i] - (*m_faces[f]->GetVertex(0))[i];
454  faceBaxis[i] =
455  (*m_faces[f]->GetVertex(2))[i] - (*m_faces[f]->GetVertex(0))[i];
456 
457  elementAaxis_length += pow(elementAaxis[i], 2);
458  elementBaxis_length += pow(elementBaxis[i], 2);
459  faceAaxis_length += pow(faceAaxis[i], 2);
460  faceBaxis_length += pow(faceBaxis[i], 2);
461  }
462 
463  elementAaxis_length = sqrt(elementAaxis_length);
464  elementBaxis_length = sqrt(elementBaxis_length);
465  faceAaxis_length = sqrt(faceAaxis_length);
466  faceBaxis_length = sqrt(faceBaxis_length);
467 
468  // Calculate the inner product of both the A-axis
469  // (i.e. Elemental A axis and face A axis)
470  for (i = 0; i < m_coordim; i++)
471  {
472  dotproduct1 += elementAaxis[i] * faceAaxis[i];
473  }
474 
475  NekDouble norm = fabs(dotproduct1) / elementAaxis_length / faceAaxis_length;
476  orientation = 0;
477 
478  // if the innerproduct is equal to the (absolute value of the ) products
479  // of the lengths of both vectors, then, the coordinate systems will NOT
480  // be transposed
481  if (fabs(norm - 1.0) < NekConstants::kNekZeroTol)
482  {
483  // if the inner product is negative, both A-axis point
484  // in reverse direction
485  if (dotproduct1 < 0.0)
486  {
487  orientation += 2;
488  }
489 
490  // calculate the inner product of both B-axis
491  for (i = 0; i < m_coordim; i++)
492  {
493  dotproduct2 += elementBaxis[i] * faceBaxis[i];
494  }
495 
496  norm = fabs(dotproduct2) / elementBaxis_length / faceBaxis_length;
497 
498  // check that both these axis are indeed parallel
499  ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
500  "These vectors should be parallel");
501 
502  // if the inner product is negative, both B-axis point
503  // in reverse direction
504  if (dotproduct2 < 0.0)
505  {
506  orientation++;
507  }
508  }
509  // The coordinate systems are transposed
510  else
511  {
512  orientation = 4;
513 
514  // Calculate the inner product between the elemental A-axis
515  // and the B-axis of the face (which are now the corresponding axis)
516  dotproduct1 = 0.0;
517  for (i = 0; i < m_coordim; i++)
518  {
519  dotproduct1 += elementAaxis[i] * faceBaxis[i];
520  }
521 
522  norm = fabs(dotproduct1) / elementAaxis_length / faceBaxis_length;
523 
524  // check that both these axis are indeed parallel
525  ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
526  "These vectors should be parallel");
527 
528  // if the result is negative, both axis point in reverse
529  // directions
530  if (dotproduct1 < 0.0)
531  {
532  orientation += 2;
533  }
534 
535  // Do the same for the other two corresponding axis
536  dotproduct2 = 0.0;
537  for (i = 0; i < m_coordim; i++)
538  {
539  dotproduct2 += elementBaxis[i] * faceAaxis[i];
540  }
541 
542  norm = fabs(dotproduct2) / elementBaxis_length / faceAaxis_length;
543 
544  // check that both these axis are indeed parallel
545  ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
546  "These vectors should be parallel");
547 
548  if (dotproduct2 < 0.0)
549  {
550  orientation++;
551  }
552  }
553 
554  orientation = orientation + 5;
555 
557  "Orientation of triangular face (id = " +
558  boost::lexical_cast<string>(m_faces[f]->GetGlobalID()) +
559  ") is inconsistent with face "+
560  boost::lexical_cast<string>(f) +
561  " of tet element (id = "+
562  boost::lexical_cast<string>(m_globalID) +
563  ") since Dir2 is aligned with Dir1. Mesh setup "
564  "needs investigation");
565 
566  // Fill the m_forient array
567  m_forient[f] = (StdRegions::Orientation)orientation;
568  }
569 }
570 
571 void TetGeom::v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
572 {
573  Geometry::v_Reset(curvedEdges, curvedFaces);
574 
575  for (int i = 0; i < 4; ++i)
576  {
577  m_faces[i]->Reset(curvedEdges, curvedFaces);
578  }
579 }
580 
582 {
583  if(!m_setupState)
584  {
585  for (int i = 0; i < 4; ++i)
586  {
587  m_faces[i]->Setup();
588  }
589  SetUpXmap();
590  SetUpCoeffs(m_xmap->GetNcoeffs());
591  m_setupState = true;
592  }
593 }
594 
595 /**
596  * Generate the geometry factors for this element.
597  */
599 {
600  if(!m_setupState)
601  {
603  }
604 
606  {
607  GeomType Gtype = eRegular;
608 
609  v_FillGeom();
610 
611  // check to see if expansions are linear
612  for (int i = 0; i < m_coordim; ++i)
613  {
614  if (m_xmap->GetBasisNumModes(0) != 2 ||
615  m_xmap->GetBasisNumModes(1) != 2 ||
616  m_xmap->GetBasisNumModes(2) != 2)
617  {
618  Gtype = eDeformed;
619  }
620  }
621 
623  Gtype, m_coordim, m_xmap, m_coeffs);
625  }
626 }
627 
628 /**
629  * @brief Set up the #m_xmap object by determining the order of each
630  * direction from derived faces.
631  */
633 {
634  vector<int> tmp;
635  tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(0));
636  int order0 = *max_element(tmp.begin(), tmp.end());
637 
638  tmp.clear();
639  tmp.push_back(order0);
640  tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(1));
641  tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(2));
642  int order1 = *max_element(tmp.begin(), tmp.end());
643 
644  tmp.clear();
645  tmp.push_back(order0);
646  tmp.push_back(order1);
647  tmp.push_back(m_faces[1]->GetXmap()->GetTraceNcoeffs(1));
648  tmp.push_back(m_faces[1]->GetXmap()->GetTraceNcoeffs(2));
649  tmp.push_back(m_faces[3]->GetXmap()->GetTraceNcoeffs(1));
650  int order2 = *max_element(tmp.begin(), tmp.end());
651 
654  order0,
656  const LibUtilities::BasisKey B(
658  order1,
661  const LibUtilities::BasisKey C(
663  order2,
666 
668 }
669 
670 }
671 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:216
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:250
Describes the specification for a Basis.
Definition: Basis.h:50
Defines a specification for a set of points.
Definition: Points.h:60
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
3D geometry information
Definition: Geometry3D.h:68
virtual void v_FillGeom()
Put all quadrature information into face/edge structure and backward transform.
Definition: Geometry3D.cpp:365
std::vector< StdRegions::Orientation > m_forient
Definition: Geometry3D.h:84
std::vector< StdRegions::Orientation > m_eorient
Definition: Geometry3D.h:83
bool m_setupState
Wether or not the setup routines have been run.
Definition: Geometry.h:195
PointGeomSharedPtr GetVertex(int i) const
Returns vertex i of this object.
Definition: Geometry.h:348
void SetUpCoeffs(const int nCoeffs)
Initialise the Geometry::m_coeffs array.
Definition: Geometry.h:658
Geometry1DSharedPtr GetEdge(int i) const
Returns edge i of this object.
Definition: Geometry.h:356
int GetVid(int i) const
Get the ID of vertex i of this object.
Definition: Geometry.cpp:137
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated Geo...
Definition: Geometry.cpp:355
int GetGlobalID(void) const
Get the ID of this object.
Definition: Geometry.h:319
LibUtilities::ShapeType m_shapeType
Type of shape.
Definition: Geometry.h:199
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
Definition: Geometry.h:203
GeomState m_geomFactorsState
State of the geometric factors.
Definition: Geometry.h:189
StdRegions::StdExpansionSharedPtr m_xmap
mapping containing isoparametric transformation.
Definition: Geometry.h:191
StdRegions::StdExpansionSharedPtr GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition: Geometry.h:426
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
Definition: Geometry.h:187
int m_coordim
Coordinate dimension of this geometry object.
Definition: Geometry.h:185
int GetEid(int i) const
Get the ID of edge i of this object.
Definition: Geometry.cpp:145
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated Geo...
Definition: TetGeom.cpp:571
static const unsigned int VertexEdgeConnectivity[4][3]
Definition: TetGeom.h:77
static const int kNfaces
Definition: TetGeom.h:58
virtual void v_GenGeomFactors()
Definition: TetGeom.cpp:598
virtual int v_GetVertexEdgeMap(const int i, const int j) const
Returns the standard element edge IDs that are connected to a given vertex.
Definition: TetGeom.cpp:98
static const int kNqfaces
Definition: TetGeom.h:56
static const int kNtfaces
Definition: TetGeom.h:57
virtual int v_GetDir(const int faceidx, const int facedir) const
Returns the element coordinate direction corresponding to a given face coordinate direction.
Definition: TetGeom.cpp:82
void SetUpXmap()
Set up the m_xmap object by determining the order of each direction from derived faces.
Definition: TetGeom.cpp:632
static const int kNedges
Definition: TetGeom.h:55
virtual int v_GetVertexFaceMap(const int i, const int j) const
Returns the standard element face IDs that are connected to a given vertex.
Definition: TetGeom.cpp:103
static const unsigned int VertexFaceConnectivity[4][3]
Definition: TetGeom.h:78
virtual int v_GetEdgeFaceMap(const int i, const int j) const
Returns the standard element edge IDs that are connected to a given face.
Definition: TetGeom.cpp:108
static const unsigned int EdgeFaceConnectivity[6][2]
Definition: TetGeom.h:79
static const int kNverts
Definition: TriGeom.h:74
@ eGaussRadauMAlpha2Beta0
Gauss Radau pinned at x=-1, .
Definition: PointsType.h:59
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:51
@ eGaussRadauMAlpha1Beta0
Gauss Radau pinned at x=-1, .
Definition: PointsType.h:58
@ eModified_B
Principle Modified Functions .
Definition: BasisType.h:49
@ eModified_C
Principle Modified Functions .
Definition: BasisType.h:50
@ eModified_A
Principle Modified Functions .
Definition: BasisType.h:48
static const NekDouble kNekZeroTol
std::unordered_map< int, CurveSharedPtr > CurveMap
Definition: Curve.hpp:62
GeomType
Indicates the type of element geometry.
@ eRegular
Geometry is straight-sided with constant geometric factors.
@ eDeformed
Geometry is curved or has non-constant factors.
std::shared_ptr< SegGeom > SegGeomSharedPtr
Definition: Geometry2D.h:62
@ ePtsFilled
Geometric information has been generated.
std::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:58
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:1
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:267