efem.py

Define the classes, methods and functions for Edge Finite Element Method (EFEM) of lowest order in tetrahedral meshes, namely, Nedelec elements.

petgem.efem.efem.NiTetrahedralSecondOrder(a1, a2, a3, a4, b1, b2, b3, b4, c1, c2, c3, c4, D, E, F, G, H, II, JJ, K, x, y, z, r)[source]

Computation of Ni (Nedelec basis functions of second order) in a tetrahedral element with vertices (x,y,z) for point r.

Parameters:
  • coefficients (float) – a1,a2,a3,a4,b1,b2,b3,b4,c1,c2,c3,c4,D,E, F,G,H,II,JJ,K.
  • x (float-array) – x-coordinates of reference element.
  • y (float-array) – y-coordinates of reference element.
  • z (float-array) – z-coordinates of reference element.
  • r (float-array) – xyz coordinates of the evaluation point.
Returns:

basis nedelec funcions of second order.

Return type:

ndarray.

petgem.efem.efem.NiTetrahedralThirdOrder(a1, a2, a3, a4, a5, a6, a7, a8, a9, a10, b1, b2, b3, b4, b5, b6, b7, b8, b9, b10, c1, c2, c3, c4, c5, c6, c7, c8, c9, c10, D, E, F, G, H, II, JJ, K, L, M, N, OO, P, Q, R, x, y, z, r)[source]

Computation of Ni (Nedelec basis functions of second order) in a tetrahedral element with vertices (x,y,z) for point r.

Parameters:
  • coefficients (float) – a1,a2,a3,a4,b1,b2,b3,b4,c1,c2,c3,c4,D,E,F, G,H,II,JJ,K.
  • x (float-array) – x-coordinates of reference element.
  • y (float-array) – y-coordinates of reference element.
  • z (float-array) – z-coordinates of reference element.
  • r (float-array) – xyz coordinates of the evaluation point.
Returns:

basis nedelec funcions of second order.

Return type:

ndarray.

petgem.efem.efem.analyticIntegral(input_functions)[source]

Compute the integrals of the product of monomials in the reference triangle.

Parameters:input_functions (float-array) – matrix of functions with coefficients of Ni before integration process.
Returns:integrals of the product of monomials in the reference triangle.
Return type:ndarray
petgem.efem.efem.aristaMappingTetrahedral(edge_vertices_old, edge_vertices_new)[source]

Mapping edges of a real element defined by edge_vertices to a reference element por edge_vertices_ref. Furthermore, compute relations between “tau” tangential vectors to edges.

Parameters:
  • edge_vertices_old (int-array) – edges definition based on vertices (old).
  • edge_vertices_new (int-array) – edges definition based on vertices (new).
Returns:

signs of unitary tangential vectors for each edge and edges mapping.

petgem.efem.efem.componentProductPolynomial(nxNi_term, function_n, component_q, face)[source]

Compute the product for a given component of nxNi.

Parameters:
  • nxNi_term (float-array) – matrix of coefficients.
  • function_n (float-array) – matrix function.
  • component_q (float-array) – id of the q component.
  • face (float-array) – face of the aforementioned parameters.
Returns:

product of nxNi, face and function_q.

Return type:

ndarray

petgem.efem.efem.computeBoundaries(elemsN, nElems, edgesNodes, nedelec_order)[source]

Compute boundaries of the domain for edge finite element computations.

Parameters:
  • elemsN (ndarray) – elements-nodes connectivity.
  • nElems (int) – number of elements in the mesh.
  • edgesN (ndarray) – edges-nodes connectivity.
  • nedelec_order (int) – nedelec element order.
Returns:

boundaries, ndofs

Return type:

ndarray

petgem.efem.efem.computeBoundaryEdges(edgesN, bfacesN)[source]

Compute boundary edges of a tetrahedral mesh.

Parameters:
  • edgesN (ndarray) – edges-nodes connectivity.
  • bfacesN (ndarray) – boundary-faces-nodes connectivity.
Returns:

boundary-edges connectivity.

Return type:

ndarray

petgem.efem.efem.computeBoundaryFaces(elemsF, facesN)[source]

Compute boundary faces of a tetrahedral mesh.

Parameters:
  • elemsF (ndarray) – elements-face connectivity.
  • facesN (ndarray) – faces-nodes connectivity.
Returns:

nodal-connectivity and indexes of boundary-faces.

Return type:

ndarray

petgem.efem.efem.computeCoefficientsSecondOrder(qface1, qface2, qface3, qface4, r_vertices_ref, edge_vertices, face_vertices)[source]

Compute the coefficients for edge basis functions of second order on reference element.

Parameters:
  • qface1 (float-array) – vector q on face 1.
  • qface2 (float-array) – vector q on face 2.
  • qface3 (float-array) – vector q on face 3.
  • qface4 (float-array) – vector q on face 4.
  • r_vertices (float-array) – vertices of reference element.
  • edge_vertices (int-array) – initialization of edges connectivity.
  • face_vertices (int-array) – initialization of faces connectivity.
Returns:

coefficients.

Return type:

float.

petgem.efem.efem.computeCoefficientsThirdOrder(qface1, qface2, qface3, qface4)[source]

Compute the coefficients for edge basis functions of third order on reference element.

Parameters:
  • qface1 (float-array) – vector q on face 1.
  • qface2 (float-array) – vector q on face 2.
  • qface3 (float-array) – vector q on face 3.
  • qface4 (float-array) – vector q on face 4.
Returns:

coefficients.

Return type:

float.

petgem.efem.efem.computeDerivativesSecondOrder(a2, a3, a4, b2, b3, b4, c2, c3, c4, D, E, F, G, H, II, JJ, K, point)[source]

Compute partial derivatives of basis functions for tetrahedral edge elements of second order (reference element)

Parameters:
  • coefficients (float) – a2,a3,a4,b2,b3,b4,c2,c3,c4,D,E,F,G,H,II,JJ,K.
  • point (float-array) – coordinates of the gaussian point.
Returns:

partial derivatives for edge element of second order.

Return type:

ndarray.

petgem.efem.efem.computeDerivativesThirdOrder(a2, a3, a4, a5, a6, a7, a8, a9, a10, b2, b3, b4, b5, b6, b7, b8, b9, b10, c2, c3, c4, c5, c6, c7, c8, c9, c10, D, E, F, G, H, II, JJ, K, L, M, N, O, P, Q, R, point)[source]

Compute partial derivatives of basis functions for tetrahedral edge elements of third order (reference element)

Parameters:
  • coefficients (float) – a2,a3,a4,a5,a6,a7,a8,a9,a10,b2,b3,b4,b5,b6, b7,b8,b9,b10,c2,c3,c4,c5,c6,c7,c8,c9,c10,D,E,F,G,H,II,JJ,K,L,M,N,O,P,Q,R.
  • point (float-array) – coordinates of the gaussian point.
Returns:

partial derivatives for edge element of third order.

Return type:

ndarray.

petgem.efem.efem.computeEdges(elemsN, nElems, nedelec_order)[source]

Compute edges of a 3D tetrahedral mesh. For edge finite element of linear order the edges are the dofs. For edge finite element of second order the dofs are computed in runtime based on edges and faces on each tetrahedral element.

Parameters:
  • elemsN (ndarray) – elements-nodes connectivity.
  • nElems (int) – number of tetrahedral elements in the mesh.
  • nedelec_order (int) – nedelec element order.
Returns:

edges connectivity and edgesNodes connectivity.

Return type:

ndarray

petgem.efem.efem.computeElementDOFs(iEle, nodesEle, edgesEle, facesEle, nodesFace, nEdges, nFaces, nedelec_order)[source]

Compute global DOFs for tetrahedral element.

Parameters:
  • iEle (int) – id of the tetrahedral element-th.
  • nodesEle (ndarray) – array with nodal indexes of element.
  • edgesEle (ndarray) – array with edge indexes of element.
  • facesEle (ndarray) – array with face indexes of element.
  • nodesFace (ndarray) – array with node indexes of faces.
  • nEdges (ndarray) – total number of edges in the mesh.
  • nFaces (ndarray) – total number of faces in the mesh.
  • nedelec_order (int) – nedelec element order.
Returns:

DOFs indexes.

Return type:

ndarray

petgem.efem.efem.computeFaces(elemsN, nElems, nedelec_order)[source]

Compute the element’s faces of a 3D tetrahedral mesh.

Parameters:
  • matrix (ndarray) – elements-nodes connectivity.
  • nElems (int) – number of elements in the mesh.
  • nedelec_order (int) – nedelec element order.
Returns:

element/faces connectivity.

Return type:

ndarray

Note

References:

Rognes, Marie E., Robert Cndarray. Kirby, and Anders Logg. “Efficient assembly of H(div) and H(curl) conforming finite elements.” SIAM Journal on Scientific Computing 31.6 (2009): 4130-4151.

petgem.efem.efem.computeMassMatrixSecondOrder(a1, a2, a3, a4, b1, b2, b3, b4, c1, c2, c3, c4, D, E, F, G, H, II, JJ, K, ref_ele, GR, signs, DetJacob, Wi, rx, ry, rz, ngaussP)[source]

Compute mass matrix for tetrahedral edge elements of second order.

Parameters:
  • coefficients (float) – a1,a2,a3,a4,b1,b2,b3,b4,c1,c2,c3,c4,D,E,F,G, H,II,JJ,K.
  • ref_ele (float-array) – nodal coordinates of reference element.
  • GR (float-array) – tensor.
  • signs (int-array) – dofs signs.
  • DetJacob (float) – determinant of the jacobian.
  • Wi (float) – gauss weigths.
  • rx (float-array) – x-coordinates of gauss points.
  • ry (float-array) – y-coordinates of gauss points.
  • rz (float-array) – z-coordinates of gauss points.
  • ngaussP (int) – number of gauss points.
Returns:

mass matrix for edge element of second order.

Return type:

ndarray.

petgem.efem.efem.computeMassMatrixThirdOrder(a1, a2, a3, a4, a5, a6, a7, a8, a9, a10, b1, b2, b3, b4, b5, b6, b7, b8, b9, b10, c1, c2, c3, c4, c5, c6, c7, c8, c9, c10, D, E, F, G, H, II, JJ, K, L, M, N, O, P, Q, R, ref_ele, GR, signs, DetJacob)[source]

Compute mass matrix for tetrahedral edge elements of third order.

Parameters:
  • coefficients (float) – a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,b1,b2,b3,b4,b5, b6,b7,b8,b9,b10,c1,c2,c3,c4,c5,c6,c7,c8,c9,c10,D,E,F,G,H,II,JJ,K,L,M,N, O,P,Q,R
  • ref_ele (float-array) – nodal coordinates of reference element.
  • GR (float-array) – tensor.
  • signs (int-array) – dofs signs.
  • DetJacob (float) – determinant of the jacobian.
Returns:

mass matrix for edge element of third order.

Return type:

ndarray.

petgem.efem.efem.computeSignsJacobLegth(eleNodes, edge_vertices, face_vertices, nedelec_order)[source]

Compute signs, edges length and jacobian for mapping.

Parameters:
  • nodesEle (ndarray) – array with nodal indexes of element.
  • edge_vertices (ndarray) – edges connectivity based on vertices.
  • face_vertices (ndarray) –

    .

  • nedelec_order (int) – nedelec element order.
Returns:

tangential unitary vectors for each edge, edges length, dofs signs, normal face vectors, jacobiano, determinant of jacobian, faces area.

Return type:

ndarray

petgem.efem.efem.computeStiffnessMatrixSecondOrder(a2, a3, a4, b2, b3, b4, c2, c3, c4, D, E, F, G, H, II, JJ, K, invjj, signs, DetJacob, Wi, rx, ry, rz, ngaussP)[source]

Compute stiffness matrix for tetrahedral edge elements of second order.

Parameters:
  • coefficients (float) – a2,a3,a4,b2,b3,b4,c2,c3,c4,D,E,F,G,H,II,JJ,K.
  • invjj (float-array) – inverse jacobian.
  • signs (int-array) – dofs signs.
  • DetJacob (float) – determinant of the jacobian.
  • Wi (float) – gauss weigths.
  • rx (float-array) – x-coordinates of gauss points.
  • ry (float-array) – y-coordinates of gauss points.
  • rz (float-array) – z-coordinates of gauss points.
  • ngaussP (int) – number of gauss points.
Returns:

stiffness matrix for edge element of second order.

Return type:

ndarray.

petgem.efem.efem.computeStiffnessMatrixThirdOrder(a2, a3, a4, a5, a6, a7, a8, a9, a10, b2, b3, b4, b5, b6, b7, b8, b9, b10, c2, c3, c4, c5, c6, c7, c8, c9, c10, D, E, F, G, H, II, JJ, K, L, M, N, O, P, Q, R, invjj, signs, DetJacob)[source]

Compute stiffness matrix for tetrahedral edge elements of third order.

Parameters:
  • coefficients (float) – a2,a3,a4,a5,a6,a7,a8,a9,a10,b2,b3,b4,b5,b6,b7, b8,b9,b10,c2,c3,c4,c5,c6,c7,c8,c9,c10,D,E,F,G,H,II,JJ,K,L,M,N,O,P,Q,R.
  • invjj (float-array) – inverse jacobian.
  • signs (int-array) – dofs signs.
  • DetJacob (float) – determinant of the jacobian.
Returns:

stiffness matrix for edge element of third order.

Return type:

ndarray.

petgem.efem.efem.defineEFEMConstants(nedelec_order)[source]

Set constants for edge finite element computations, namely, order of edge finite elements (edgeOrder), order of nodal finite elements (nodalOrder) and number of dimensions (numDimensions).

Parameters:nedelec_order (int) – nedelec element order.
Returns:edgeOrder, nodalOrder and numDimensions.
Return type:integer.
petgem.efem.efem.definitionHighOrderTet(nreal, signs, jacob, nedelec_order)[source]

Compute “q” vectors on faces -dof definition of edge tetrahedral element of second and third order.

Parameters:
  • nreal (float-array) – normal face vectors.
  • signs (float-array) – dofs signs.
  • jacob (float) – jacobian.
  • nedelec_order (int) – nedelec element order.
Returns:

q vectors on faces, inverse jacobian, 3x3 tensor.

Return type:

ndarray.

petgem.efem.efem.edgeFaceVerticesInit()[source]

Initialization of arrays that define edges and faces of tetrahedral elements. Vector basis functions are based on these arrays.

Param:None.
Returns:edge_vertices and face_vertices.
Return type:ndarray
petgem.efem.efem.faceMappingTetrahedral(face_vertices_old, face_vertices_new, stage)[source]

Compute mapping between faces of a tetrahedron with face_vertices to a tetrahedron with faces_vertices_ref.

Parameters:
  • face_vertices_old (int-array) – faces definition based on vertices (old).
  • face_vertices_new (int-array) – faces definition based on vertices (new).
  • face (int) – number of stage of mapping (1 or 2).
Returns:

signs of unitary tangential vectors for each face and faces mapping.

petgem.efem.efem.isigno(xx, yy, zz)[source]

Compute signs for edges in a given tetrahedral element.

Parameters:
  • xx (float) – x-coordinates of the nodes.
  • yy (float) – x-coordinates of the nodes.
  • zz (float) – x-coordinates of the nodes.
Returns:

signs of edges.

Return type:

ndarray

petgem.efem.efem.nedelecBasisIterative(eleNodes, points, eleVol, lengthEdges, edgeOrder)[source]

Compute the basis Nedelec functions in an iterative way for a set of points in a given element.

Parameters:
  • eleNodes (ndarray) – nodal spatial coordinates of the element
  • points (ndarray) – spatial coordinates of the evaluation points
  • eleVol (float) – element’s volume
  • lengthEdges (ndarray) – element’s edges defined by their length
  • edgeOrder (int) – order of tetrahedral edge element
Returns:

values of Nedelec functions.

Return type:

ndarray.

petgem.efem.efem.nedelecBasisSecondOrder(a1, a2, a3, a4, b1, b2, b3, b4, c1, c2, c3, c4, D, E, F, G, H, II, JJ, K, ref_ele, points)[source]

This function computes the basis nedelec functions of second order for a set of points in a given tetrahedral element.

Parameters:
  • coefficients (float) – a1,a2,a3,a4,b1,b2,b3,b4,c1,c2,c3,c4, D,E,F,G,H,II,JJ,K.
  • ref_ele (float-array) – nodal coordinates of reference element.
  • points (float-array) – spatial coordinates of the evaluation points.
Returns:

basis nedelec functions of second order.

Return type:

ndarray.

petgem.efem.efem.nedelecBasisThirdOrder(a1, a2, a3, a4, a5, a6, a7, a8, a9, a10, b1, b2, b3, b4, b5, b6, b7, b8, b9, b10, c1, c2, c3, c4, c5, c6, c7, c8, c9, c10, D, E, F, G, H, II, JJ, K, L, M, N, O, P, Q, R, ref_ele, points)[source]

This function computes the basis nedelec functions of third order for a set of points in a given tetrahedral element.

Parameters:
  • coefficients (float) – a1,a2,a3,a4,a5,a6,a7,a8,a9,a10,b1,b2,b3,b4,b5, b6,b7,b8,b9,b10,c1,c2,c3,c4,c5,c6,c7,c8,c9,c10,D,E,F,G,H,II,JJ,K,L,M,N, O,P,Q,R.
  • ref_ele (float-array) – nodal coordinates of reference element.
  • points (float-array) – spatial coordinates of the evaluation points.
Returns:

basis nedelec functions of second order.

Return type:

ndarray.

petgem.efem.efem.polynomialProduct(poly1, poly2)[source]

Compute the product of two polynomials according to convention for tetrahedral edge element of third order.

Parameters:
  • poly1 (float-array) – first polynomial to be computed.
  • poly2 (float-array) – second polynomial to be computed.
Returns:

polynomial product.

Return type:

ndarray

petgem.efem.efem.unitary_test()[source]

Unitary test for efem.py script.