Nektar++
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
36
40
42{
43const unsigned int TetGeom::VertexEdgeConnectivity[4][3] = {
44 {0, 2, 3}, {0, 1, 4}, {1, 2, 5}, {3, 4, 5}};
45const unsigned int TetGeom::VertexFaceConnectivity[4][3] = {
46 {0, 1, 3}, {0, 1, 2}, {0, 2, 3}, {1, 2, 3}};
47const unsigned int TetGeom::EdgeFaceConnectivity[6][2] = {
48 {0, 1}, {0, 2}, {0, 3}, {1, 3}, {1, 2}, {2, 3}};
49const unsigned int TetGeom::EdgeNormalToFaceVert[4][3] = {
50 {3, 4, 5}, {1, 2, 5}, {0, 2, 3}, {0, 1, 4}};
51
53{
55}
56
57TetGeom::TetGeom(int id, const TriGeomSharedPtr faces[])
58 : Geometry3D(faces[0]->GetEdge(0)->GetVertex(0)->GetCoordim())
59{
61 m_globalID = id;
62
63 /// Copy the face shared pointers
64 m_faces.insert(m_faces.begin(), faces, faces + TetGeom::kNfaces);
65
66 /// Set up orientation vectors with correct amount of elements.
67 m_eorient.resize(kNedges);
68 m_forient.resize(kNfaces);
69
74}
75
76int TetGeom::v_GetDir(const int faceidx, const int facedir) const
77{
78 if (faceidx == 0)
79 {
80 return facedir;
81 }
82 else if (faceidx == 1)
83 {
84 return 2 * facedir;
85 }
86 else
87 {
88 return 1 + facedir;
89 }
90}
91
92int TetGeom::v_GetVertexEdgeMap(const int i, const int j) const
93{
94 return VertexEdgeConnectivity[i][j];
95}
96
97int TetGeom::v_GetVertexFaceMap(const int i, const int j) const
98{
99 return VertexFaceConnectivity[i][j];
100}
101
102int TetGeom::v_GetEdgeFaceMap(const int i, const int j) const
103{
104 return EdgeFaceConnectivity[i][j];
105}
106
107int TetGeom::v_GetEdgeNormalToFaceVert(const int i, const int j) const
108{
109 return EdgeNormalToFaceVert[i][j];
110}
111
113{
114
115 // find edge 0
116 int i, j;
117 unsigned int check;
118
119 SegGeomSharedPtr edge;
120
121 // First set up the 3 bottom edges
122
123 if (m_faces[0]->GetEid(0) != m_faces[1]->GetEid(0))
124 {
125 std::ostringstream errstrm;
126 errstrm << "Local edge 0 (eid=" << m_faces[0]->GetEid(0);
127 errstrm << ") on face " << m_faces[0]->GetGlobalID();
128 errstrm << " must be the same as local edge 0 (eid="
129 << m_faces[1]->GetEid(0);
130 errstrm << ") on face " << m_faces[1]->GetGlobalID();
131 ASSERTL0(false, errstrm.str());
132 }
133
134 int faceConnected;
135 for (faceConnected = 1; faceConnected < 4; faceConnected++)
136 {
137 check = 0;
138 for (i = 0; i < 3; i++)
139 {
140 if ((m_faces[0])->GetEid(i) == (m_faces[faceConnected])->GetEid(0))
141 {
142 edge = std::dynamic_pointer_cast<SegGeom>(
143 (m_faces[0])->GetEdge(i));
144 m_edges.push_back(edge);
145 check++;
146 }
147 }
148
149 if (check < 1)
150 {
151 std::ostringstream errstrm;
152 errstrm << "Face 0 does not share an edge with first edge of "
153 "adjacent face. Faces ";
154 errstrm << (m_faces[0])->GetGlobalID() << ", "
155 << (m_faces[faceConnected])->GetGlobalID();
156 ASSERTL0(false, errstrm.str());
157 }
158 else if (check > 1)
159 {
160 std::ostringstream errstrm;
161 errstrm << "Connected faces share more than one edge. Faces ";
162 errstrm << (m_faces[0])->GetGlobalID() << ", "
163 << (m_faces[faceConnected])->GetGlobalID();
164 ASSERTL0(false, errstrm.str());
165 }
166 }
167
168 // Then, set up the 3 vertical edges
169 check = 0;
170 for (i = 0; i < 3; i++) // Set up the vertical edge :face(1) and face(3)
171 {
172 for (j = 0; j < 3; j++)
173 {
174 if ((m_faces[1])->GetEid(i) == (m_faces[3])->GetEid(j))
175 {
176 edge = std::dynamic_pointer_cast<SegGeom>(
177 (m_faces[1])->GetEdge(i));
178 m_edges.push_back(edge);
179 check++;
180 }
181 }
182 }
183 if (check < 1)
184 {
185 std::ostringstream errstrm;
186 errstrm << "Connected faces do not share an edge. Faces ";
187 errstrm << (m_faces[1])->GetGlobalID() << ", "
188 << (m_faces[3])->GetGlobalID();
189 ASSERTL0(false, errstrm.str());
190 }
191 else if (check > 1)
192 {
193 std::ostringstream errstrm;
194 errstrm << "Connected faces share more than one edge. Faces ";
195 errstrm << (m_faces[1])->GetGlobalID() << ", "
196 << (m_faces[3])->GetGlobalID();
197 ASSERTL0(false, errstrm.str());
198 }
199 // Set up vertical edges: face(1) through face(3)
200 for (faceConnected = 1; faceConnected < 3; faceConnected++)
201 {
202 check = 0;
203 for (i = 0; i < 3; i++)
204 {
205 for (j = 0; j < 3; j++)
206 {
207 if ((m_faces[faceConnected])->GetEid(i) ==
208 (m_faces[faceConnected + 1])->GetEid(j))
209 {
210 edge = std::dynamic_pointer_cast<SegGeom>(
211 (m_faces[faceConnected])->GetEdge(i));
212 m_edges.push_back(edge);
213 check++;
214 }
215 }
216 }
217
218 if (check < 1)
219 {
220 std::ostringstream errstrm;
221 errstrm << "Connected faces do not share an edge. Faces ";
222 errstrm << (m_faces[faceConnected])->GetGlobalID() << ", "
223 << (m_faces[faceConnected + 1])->GetGlobalID();
224 ASSERTL0(false, errstrm.str());
225 }
226 else if (check > 1)
227 {
228 std::ostringstream errstrm;
229 errstrm << "Connected faces share more than one edge. Faces ";
230 errstrm << (m_faces[faceConnected])->GetGlobalID() << ", "
231 << (m_faces[faceConnected + 1])->GetGlobalID();
232 ASSERTL0(false, errstrm.str());
233 }
234 }
235}
236
238{
239
240 // Set up the first 2 vertices (i.e. vertex 0,1)
241 if ((m_edges[0]->GetVid(0) == m_edges[1]->GetVid(0)) ||
242 (m_edges[0]->GetVid(0) == m_edges[1]->GetVid(1)))
243 {
244 m_verts.push_back(m_edges[0]->GetVertex(1));
245 m_verts.push_back(m_edges[0]->GetVertex(0));
246 }
247 else if ((m_edges[0]->GetVid(1) == m_edges[1]->GetVid(0)) ||
248 (m_edges[0]->GetVid(1) == m_edges[1]->GetVid(1)))
249 {
250 m_verts.push_back(m_edges[0]->GetVertex(0));
251 m_verts.push_back(m_edges[0]->GetVertex(1));
252 }
253 else
254 {
255 std::ostringstream errstrm;
256 errstrm << "Connected edges do not share a vertex. Edges ";
257 errstrm << m_edges[0]->GetGlobalID() << ", "
258 << m_edges[1]->GetGlobalID();
259 ASSERTL0(false, errstrm.str());
260 }
261
262 // set up the other bottom vertices (i.e. vertex 2)
263 for (int i = 1; i < 2; i++)
264 {
265 if (m_edges[i]->GetVid(0) == m_verts[i]->GetGlobalID())
266 {
267 m_verts.push_back(m_edges[i]->GetVertex(1));
268 }
269 else if (m_edges[i]->GetVid(1) == m_verts[i]->GetGlobalID())
270 {
271 m_verts.push_back(m_edges[i]->GetVertex(0));
272 }
273 else
274 {
275 std::ostringstream errstrm;
276 errstrm << "Connected edges do not share a vertex. Edges ";
277 errstrm << m_edges[i]->GetGlobalID() << ", "
278 << m_edges[i - 1]->GetGlobalID();
279 ASSERTL0(false, errstrm.str());
280 }
281 }
282
283 // set up top vertex
284 if (m_edges[3]->GetVid(0) == m_verts[0]->GetGlobalID())
285 {
286 m_verts.push_back(m_edges[3]->GetVertex(1));
287 }
288 else
289 {
290 m_verts.push_back(m_edges[3]->GetVertex(0));
291 }
292
293 // Check the other edges match up.
294 int check = 0;
295 for (int i = 4; i < 6; ++i)
296 {
297 if ((m_edges[i]->GetVid(0) == m_verts[i - 3]->GetGlobalID() &&
298 m_edges[i]->GetVid(1) == m_verts[3]->GetGlobalID()) ||
299 (m_edges[i]->GetVid(1) == m_verts[i - 3]->GetGlobalID() &&
300 m_edges[i]->GetVid(0) == m_verts[3]->GetGlobalID()))
301 {
302 check++;
303 }
304 }
305 if (check != 2)
306 {
307 std::ostringstream errstrm;
308 errstrm << "Connected edges do not share a vertex. Edges ";
309 errstrm << m_edges[3]->GetGlobalID() << ", "
310 << m_edges[2]->GetGlobalID();
311 ASSERTL0(false, errstrm.str());
312 }
313}
314
316{
317
318 // This 2D array holds the local id's of all the vertices
319 // for every edge. For every edge, they are ordered to what we
320 // define as being Forwards
321 const unsigned int edgeVerts[kNedges][2] = {{0, 1}, {1, 2}, {0, 2},
322 {0, 3}, {1, 3}, {2, 3}};
323
324 int i;
325 for (i = 0; i < kNedges; i++)
326 {
327 if (m_edges[i]->GetVid(0) == m_verts[edgeVerts[i][0]]->GetGlobalID())
328 {
330 }
331 else if (m_edges[i]->GetVid(0) ==
332 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 =
476 fabs(dotproduct1) / elementAaxis_length / faceAaxis_length;
477 orientation = 0;
478
479 // if the innerproduct is equal to the (absolute value of the ) products
480 // of the lengths of both vectors, then, the coordinate systems will NOT
481 // be transposed
482 if (fabs(norm - 1.0) < NekConstants::kNekZeroTol)
483 {
484 // if the inner product is negative, both A-axis point
485 // in reverse direction
486 if (dotproduct1 < 0.0)
487 {
488 orientation += 2;
489 }
490
491 // calculate the inner product of both B-axis
492 for (i = 0; i < m_coordim; i++)
493 {
494 dotproduct2 += elementBaxis[i] * faceBaxis[i];
495 }
496
497 norm = fabs(dotproduct2) / elementBaxis_length / faceBaxis_length;
498
499 // check that both these axis are indeed parallel
500 ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
501 "These vectors should be parallel");
502
503 // if the inner product is negative, both B-axis point
504 // in reverse direction
505 if (dotproduct2 < 0.0)
506 {
507 orientation++;
508 }
509 }
510 // The coordinate systems are transposed
511 else
512 {
513 orientation = 4;
514
515 // Calculate the inner product between the elemental A-axis
516 // and the B-axis of the face (which are now the corresponding axis)
517 dotproduct1 = 0.0;
518 for (i = 0; i < m_coordim; i++)
519 {
520 dotproduct1 += elementAaxis[i] * faceBaxis[i];
521 }
522
523 norm = fabs(dotproduct1) / elementAaxis_length / faceBaxis_length;
524
525 // check that both these axis are indeed parallel
526 ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
527 "These vectors should be parallel");
528
529 // if the result is negative, both axis point in reverse
530 // directions
531 if (dotproduct1 < 0.0)
532 {
533 orientation += 2;
534 }
535
536 // Do the same for the other two corresponding axis
537 dotproduct2 = 0.0;
538 for (i = 0; i < m_coordim; i++)
539 {
540 dotproduct2 += elementBaxis[i] * faceAaxis[i];
541 }
542
543 norm = fabs(dotproduct2) / elementBaxis_length / faceAaxis_length;
544
545 // check that both these axis are indeed parallel
546 ASSERTL1(fabs(norm - 1.0) < NekConstants::kNekZeroTol,
547 "These vectors should be parallel");
548
549 if (dotproduct2 < 0.0)
550 {
551 orientation++;
552 }
553 }
554
555 orientation = orientation + 5;
556
558 "Orientation of triangular face (id = " +
559 std::to_string(m_faces[f]->GetGlobalID()) +
560 ") is inconsistent with face " + std::to_string(f) +
561 " of tet element (id = " + std::to_string(m_globalID) +
562 ") since Dir2 is aligned with Dir1. Mesh setup "
563 "needs investigation");
564
565 // Fill the m_forient array
566 m_forient[f] = (StdRegions::Orientation)orientation;
567 }
568}
569
570void TetGeom::v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
571{
572 Geometry::v_Reset(curvedEdges, curvedFaces);
573
574 for (int i = 0; i < 4; ++i)
575 {
576 m_faces[i]->Reset(curvedEdges, curvedFaces);
577 }
578}
579
581{
582 if (!m_setupState)
583 {
584 for (int i = 0; i < 4; ++i)
585 {
586 m_faces[i]->Setup();
587 }
588 SetUpXmap();
589 SetUpCoeffs(m_xmap->GetNcoeffs());
590 m_setupState = true;
591 }
592}
593
594/**
595 * Generate the geometry factors for this element.
596 */
598{
599 if (!m_setupState)
600 {
602 }
603
605 {
606 GeomType Gtype = eRegular;
607
608 v_FillGeom();
609
610 // check to see if expansions are linear
611 m_straightEdge = 1;
612 if (m_xmap->GetBasisNumModes(0) != 2 ||
613 m_xmap->GetBasisNumModes(1) != 2 ||
614 m_xmap->GetBasisNumModes(2) != 2)
615 {
616 Gtype = eDeformed;
617 m_straightEdge = 0;
618 }
619
620 if (Gtype == eRegular)
621 {
623 for (int i = 0; i < 3; ++i)
624 {
626 NekDouble A = (*m_verts[0])(i);
627 NekDouble B = (*m_verts[1])(i);
628 NekDouble C = (*m_verts[2])(i);
629 NekDouble D = (*m_verts[3])(i);
630 m_isoParameter[i][0] = 0.5 * (-A + B + C + D);
631
632 m_isoParameter[i][1] = 0.5 * (-A + B); // xi1
633 m_isoParameter[i][2] = 0.5 * (-A + C); // xi2
634 m_isoParameter[i][3] = 0.5 * (-A + D); // xi3
635 }
637 }
638
640 Gtype, m_coordim, m_xmap, m_coeffs);
642 }
643}
644
645/**
646 * @brief Set up the #m_xmap object by determining the order of each
647 * direction from derived faces.
648 */
650{
651 std::vector<int> tmp;
652 tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(0));
653 int order0 = *std::max_element(tmp.begin(), tmp.end());
654
655 tmp.clear();
656 tmp.push_back(order0);
657 tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(1));
658 tmp.push_back(m_faces[0]->GetXmap()->GetTraceNcoeffs(2));
659 int order1 = *std::max_element(tmp.begin(), tmp.end());
660
661 tmp.clear();
662 tmp.push_back(order0);
663 tmp.push_back(order1);
664 tmp.push_back(m_faces[1]->GetXmap()->GetTraceNcoeffs(1));
665 tmp.push_back(m_faces[1]->GetXmap()->GetTraceNcoeffs(2));
666 tmp.push_back(m_faces[3]->GetXmap()->GetTraceNcoeffs(1));
667 int order2 = *std::max_element(tmp.begin(), tmp.end());
668
671 LibUtilities::PointsKey(order0 + 1,
675 LibUtilities::PointsKey(order1, LibUtilities::eGaussRadauMAlpha1Beta0));
678 LibUtilities::PointsKey(order2, LibUtilities::eGaussRadauMAlpha2Beta0));
679
681}
682
683} // namespace Nektar::SpatialDomains
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
Describes the specification for a Basis.
Definition: Basis.h:45
Defines a specification for a set of points.
Definition: Points.h:50
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
3D geometry information
Definition: Geometry3D.h:65
std::vector< StdRegions::Orientation > m_forient
Definition: Geometry3D.h:81
void v_CalculateInverseIsoParam() override
Definition: Geometry3D.cpp:565
void v_FillGeom() override
Put all quadrature information into face/edge structure and backward transform.
Definition: Geometry3D.cpp:442
std::vector< StdRegions::Orientation > m_eorient
Definition: Geometry3D.h:80
bool m_setupState
Wether or not the setup routines have been run.
Definition: Geometry.h:198
PointGeomSharedPtr GetVertex(int i) const
Returns vertex i of this object.
Definition: Geometry.h:357
void SetUpCoeffs(const int nCoeffs)
Initialise the Geometry::m_coeffs array.
Definition: Geometry.h:700
Geometry1DSharedPtr GetEdge(int i) const
Returns edge i of this object.
Definition: Geometry.h:365
int GetVid(int i) const
Get the ID of vertex i of this object.
Definition: Geometry.cpp:104
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
Definition: Geometry.h:209
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:364
int GetGlobalID(void) const
Get the ID of this object.
Definition: Geometry.h:326
LibUtilities::ShapeType m_shapeType
Type of shape.
Definition: Geometry.h:202
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
Definition: Geometry.h:206
GeomState m_geomFactorsState
State of the geometric factors.
Definition: Geometry.h:192
StdRegions::StdExpansionSharedPtr m_xmap
mapping containing isoparametric transformation.
Definition: Geometry.h:194
StdRegions::StdExpansionSharedPtr GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
Definition: Geometry.h:435
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
Definition: Geometry.h:190
int m_coordim
Coordinate dimension of this geometry object.
Definition: Geometry.h:188
int GetEid(int i) const
Get the ID of edge i of this object.
Definition: Geometry.cpp:112
int v_GetEdgeFaceMap(const int i, const int j) const override
Returns the standard element edge IDs that are connected to a given face.
Definition: TetGeom.cpp:102
int v_GetVertexFaceMap(const int i, const int j) const override
Returns the standard element face IDs that are connected to a given vertex.
Definition: TetGeom.cpp:97
int v_GetVertexEdgeMap(const int i, const int j) const override
Returns the standard element edge IDs that are connected to a given vertex.
Definition: TetGeom.cpp:92
static const unsigned int VertexEdgeConnectivity[4][3]
Definition: TetGeom.h:76
static const unsigned int EdgeNormalToFaceVert[4][3]
Definition: TetGeom.h:79
int v_GetEdgeNormalToFaceVert(const int i, const int j) const override
Returns the standard lement edge IDs that are normal to a given face vertex.
Definition: TetGeom.cpp:107
static const int kNfaces
Definition: TetGeom.h:56
void v_GenGeomFactors() override
Definition: TetGeom.cpp:597
void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces) override
Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated Geo...
Definition: TetGeom.cpp:570
static const int kNqfaces
Definition: TetGeom.h:54
static const int kNtfaces
Definition: TetGeom.h:55
void SetUpXmap()
Set up the m_xmap object by determining the order of each direction from derived faces.
Definition: TetGeom.cpp:649
static const int kNedges
Definition: TetGeom.h:53
int v_GetDir(const int faceidx, const int facedir) const override
Returns the element coordinate direction corresponding to a given face coordinate direction.
Definition: TetGeom.cpp:76
static const unsigned int VertexFaceConnectivity[4][3]
Definition: TetGeom.h:77
static const unsigned int EdgeFaceConnectivity[6][2]
Definition: TetGeom.h:78
static const int kNverts
Definition: TriGeom.h:71
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:51
@ 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:59
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:58
@ ePtsFilled
Geometric information has been generated.
std::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:56
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285