44 {0, 2, 3}, {0, 1, 4}, {1, 2, 5}, {3, 4, 5}};
46 {0, 1, 3}, {0, 1, 2}, {0, 2, 3}, {1, 2, 3}};
48 {0, 1}, {0, 2}, {0, 3}, {1, 3}, {1, 2}, {2, 3}};
50 {3, 4, 5}, {1, 2, 5}, {0, 2, 3}, {0, 1, 4}};
58 :
Geometry3D(faces[0]->GetEdge(0)->GetVertex(0)->GetCoordim())
82 else if (faceidx == 1)
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="
130 errstrm <<
") on face " <<
m_faces[1]->GetGlobalID();
135 for (faceConnected = 1; faceConnected < 4; faceConnected++)
138 for (i = 0; i < 3; i++)
142 edge = std::dynamic_pointer_cast<SegGeom>(
151 std::ostringstream errstrm;
152 errstrm <<
"Face 0 does not share an edge with first edge of "
153 "adjacent face. Faces ";
160 std::ostringstream errstrm;
161 errstrm <<
"Connected faces share more than one edge. Faces ";
170 for (i = 0; i < 3; i++)
172 for (j = 0; j < 3; j++)
176 edge = std::dynamic_pointer_cast<SegGeom>(
185 std::ostringstream errstrm;
186 errstrm <<
"Connected faces do not share an edge. Faces ";
193 std::ostringstream errstrm;
194 errstrm <<
"Connected faces share more than one edge. Faces ";
200 for (faceConnected = 1; faceConnected < 3; faceConnected++)
203 for (i = 0; i < 3; i++)
205 for (j = 0; j < 3; j++)
210 edge = std::dynamic_pointer_cast<SegGeom>(
220 std::ostringstream errstrm;
221 errstrm <<
"Connected faces do not share an edge. Faces ";
228 std::ostringstream errstrm;
229 errstrm <<
"Connected faces share more than one edge. Faces ";
255 std::ostringstream errstrm;
256 errstrm <<
"Connected edges do not share a vertex. Edges ";
257 errstrm <<
m_edges[0]->GetGlobalID() <<
", "
263 for (
int i = 1; i < 2; i++)
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();
295 for (
int i = 4; i < 6; ++i)
307 std::ostringstream errstrm;
308 errstrm <<
"Connected edges do not share a vertex. Edges ";
309 errstrm <<
m_edges[3]->GetGlobalID() <<
", "
321 const unsigned int edgeVerts[
kNedges][2] = {{0, 1}, {1, 2}, {0, 2},
322 {0, 3}, {1, 3}, {2, 3}};
338 ASSERTL0(
false,
"Could not find matching vertex for the edge");
367 unsigned int baseVertex;
380 {0, 1, 2}, {0, 1, 3}, {1, 2, 3}, {0, 2, 3}};
385 unsigned int orientation;
391 elementAaxis_length = 0.0;
392 elementBaxis_length = 0.0;
393 faceAaxis_length = 0.0;
394 faceBaxis_length = 0.0;
399 baseVertex =
m_faces[f]->GetVid(0);
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];
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];
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];
445 ASSERTL0(
false,
"Could not find matching vertex for the face");
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);
463 elementAaxis_length =
sqrt(elementAaxis_length);
464 elementBaxis_length =
sqrt(elementBaxis_length);
465 faceAaxis_length =
sqrt(faceAaxis_length);
466 faceBaxis_length =
sqrt(faceBaxis_length);
472 dotproduct1 += elementAaxis[i] * faceAaxis[i];
476 fabs(dotproduct1) / elementAaxis_length / faceAaxis_length;
486 if (dotproduct1 < 0.0)
494 dotproduct2 += elementBaxis[i] * faceBaxis[i];
497 norm = fabs(dotproduct2) / elementBaxis_length / faceBaxis_length;
501 "These vectors should be parallel");
505 if (dotproduct2 < 0.0)
520 dotproduct1 += elementAaxis[i] * faceBaxis[i];
523 norm = fabs(dotproduct1) / elementAaxis_length / faceBaxis_length;
527 "These vectors should be parallel");
531 if (dotproduct1 < 0.0)
540 dotproduct2 += elementBaxis[i] * faceAaxis[i];
543 norm = fabs(dotproduct2) / elementBaxis_length / faceAaxis_length;
547 "These vectors should be parallel");
549 if (dotproduct2 < 0.0)
555 orientation = orientation + 5;
558 "Orientation of triangular face (id = " +
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");
574 for (
int i = 0; i < 4; ++i)
576 m_faces[i]->Reset(curvedEdges, curvedFaces);
584 for (
int i = 0; i < 4; ++i)
612 if (
m_xmap->GetBasisNumModes(0) != 2 ||
613 m_xmap->GetBasisNumModes(1) != 2 ||
614 m_xmap->GetBasisNumModes(2) != 2)
623 for (
int i = 0; i < 3; ++i)
651 std::vector<int> tmp;
653 int order0 = *std::max_element(tmp.begin(), tmp.end());
656 tmp.push_back(order0);
659 int order1 = *std::max_element(tmp.begin(), tmp.end());
662 tmp.push_back(order0);
663 tmp.push_back(order1);
667 int order2 = *std::max_element(tmp.begin(), tmp.end());
#define ASSERTL0(condition, msg)
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Describes the specification for a Basis.
Defines a specification for a set of points.
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::vector< StdRegions::Orientation > m_forient
void v_CalculateInverseIsoParam() override
void v_FillGeom() override
Put all quadrature information into face/edge structure and backward transform.
std::vector< StdRegions::Orientation > m_eorient
bool m_setupState
Wether or not the setup routines have been run.
PointGeomSharedPtr GetVertex(int i) const
Returns vertex i of this object.
void SetUpCoeffs(const int nCoeffs)
Initialise the Geometry::m_coeffs array.
Geometry1DSharedPtr GetEdge(int i) const
Returns edge i of this object.
int GetVid(int i) const
Get the ID of vertex i of this object.
Array< OneD, Array< OneD, NekDouble > > m_isoParameter
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated Geo...
int GetGlobalID(void) const
Get the ID of this object.
LibUtilities::ShapeType m_shapeType
Type of shape.
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Array containing expansion coefficients of m_xmap.
GeomState m_geomFactorsState
State of the geometric factors.
StdRegions::StdExpansionSharedPtr m_xmap
mapping containing isoparametric transformation.
StdRegions::StdExpansionSharedPtr GetXmap() const
Return the mapping object Geometry::m_xmap that represents the coordinate transformation from standar...
GeomFactorsSharedPtr m_geomFactors
Geometric factors.
int m_coordim
Coordinate dimension of this geometry object.
int GetEid(int i) const
Get the ID of edge i of this object.
int v_GetEdgeFaceMap(const int i, const int j) const override
Returns the standard element edge IDs that are connected to a given face.
int v_GetVertexFaceMap(const int i, const int j) const override
Returns the standard element face IDs that are connected to a given vertex.
int v_GetVertexEdgeMap(const int i, const int j) const override
Returns the standard element edge IDs that are connected to a given vertex.
static const unsigned int VertexEdgeConnectivity[4][3]
static const unsigned int EdgeNormalToFaceVert[4][3]
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.
void SetUpFaceOrientation()
void SetUpEdgeOrientation()
void SetUpLocalVertices()
void v_GenGeomFactors() override
void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces) override
Reset this geometry object: unset the current state, zero Geometry::m_coeffs and remove allocated Geo...
static const int kNqfaces
static const int kNtfaces
void SetUpXmap()
Set up the m_xmap object by determining the order of each direction from derived faces.
int v_GetDir(const int faceidx, const int facedir) const override
Returns the element coordinate direction corresponding to a given face coordinate direction.
static const unsigned int VertexFaceConnectivity[4][3]
static const unsigned int EdgeFaceConnectivity[6][2]
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
@ eModified_B
Principle Modified Functions .
@ eModified_C
Principle Modified Functions .
@ eModified_A
Principle Modified Functions .
static const NekDouble kNekZeroTol
std::unordered_map< int, CurveSharedPtr > CurveMap
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
@ ePtsFilled
Geometric information has been generated.
std::shared_ptr< TriGeom > TriGeomSharedPtr
@ eDir1FwdDir2_Dir2FwdDir1
scalarT< T > sqrt(scalarT< T > in)