v0.15.0
Loading...
Searching...
No Matches
MoFEM::Tools Struct Reference

Auxiliary tools. More...

#include "src/interfaces/Tools.hpp"

Inheritance diagram for MoFEM::Tools:
[legend]
Collaboration diagram for MoFEM::Tools:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 Tools (const MoFEM::Core &core)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface.
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface reference to pointer of interface.
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface.
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface.
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface.
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface.
 
virtual ~UnknownInterface ()=default
 

Static Public Member Functions

static MoFEMErrorCode outerProductOfEdgeIntegrationPtsForQuad (MatrixDouble &pts, const int edge0, const int edge1)
 
static MoFEMErrorCode outerProductOfEdgeIntegrationPtsForHex (MatrixDouble &pts, const int edge0, const int edge1, const int edge2)
 
- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version.
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version.
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version.
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version.
 

Public Attributes

MoFEM::CorecOre
 

Computational

enum  SEGMENT_MIN_DISTANCE {
  SOLUTION_EXIST , SEGMENT_ONE_IS_POINT , SEGMENT_TWO_IS_POINT , SEGMENT_TWO_AND_TWO_ARE_POINT ,
  NO_SOLUTION
}
 
static constexpr double shapeFunMBEDGE0At00 = N_MBEDGE0(0)
 
static constexpr double shapeFunMBEDGE1At00 = N_MBEDGE1(0)
 
static constexpr std::array< double, 2 > shapeFunMBEDGEAt00
 Array of shape function at zero local point on reference element.
 
static constexpr double diffN_MBEDGE0x = diffN_MBEDGE0
 
static constexpr double diffN_MBEDGE1x = diffN_MBEDGE1
 
static constexpr std::array< double, 2 > diffShapeFunMBEDGE
 
static constexpr double diffShapeFunMBTRI0x
 derivative of triangle shape function
 
static constexpr double diffShapeFunMBTRI0y
 derivative of triangle shape function
 
static constexpr double diffShapeFunMBTRI1x
 derivative of triangle shape function
 
static constexpr double diffShapeFunMBTRI1y
 derivative of triangle shape function
 
static constexpr double diffShapeFunMBTRI2x
 derivative of triangle shape function
 
static constexpr double diffShapeFunMBTRI2y
 derivative of triangle shape function
 
static constexpr std::array< double, 6 > diffShapeFunMBTRI
 
static constexpr double shapeFunMBTRI0At00 = N_MBTRI0(0, 0)
 
static constexpr double shapeFunMBTRI1At00 = N_MBTRI1(0, 0)
 
static constexpr double shapeFunMBTRI2At00 = N_MBTRI2(0, 0)
 
static constexpr std::array< double, 3 > shapeFunMBTRIAt00
 Array of shape function at zero local point on reference element.
 
static constexpr double diffShapeFunMBQUADAtCenter0x
 derivative of quad shape function
 
static constexpr double diffShapeFunMBQUADAtCenter0y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBQUADAtCenter1x
 derivative of quad shape function
 
static constexpr double diffShapeFunMBQUADAtCenter1y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBQUADAtCenter2x
 derivative of quad shape function
 
static constexpr double diffShapeFunMBQUADAtCenter2y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBQUADAtCenter3x
 derivative of quad shape function
 
static constexpr double diffShapeFunMBQUADAtCenter3y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter0x
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter0y
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter0z
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter1x
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter1y
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter1z
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter2x
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter2y
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter2z
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter3x
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter3y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter3z
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter4x
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter4y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter4z
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter5x
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter5y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter5z
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter6x
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter6y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter6z
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter7x
 derivative of HEX shape function
 
static constexpr double diffShapeFunMBHEXAtCenter7y
 derivative of quad shape function
 
static constexpr double diffShapeFunMBHEXAtCenter7z
 derivative of quad shape function
 
static constexpr std::array< double, 8 > diffShapeFunMBQUADAtCenter
 
static constexpr std::array< double, 24 > diffShapeFunMBHEXAtCenter
 
static constexpr double diffShapeFunMBTET0x
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET0y
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET0z
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET1x
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET1y
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET1z
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET2x
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET2y
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET2z
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET3x
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET3y
 derivative of tetrahedral shape function
 
static constexpr double diffShapeFunMBTET3z
 derivative of tetrahedral shape function
 
static constexpr std::array< double, 12 > diffShapeFunMBTET
 
static constexpr double shapeFunMBTET0At000 = N_MBTET0(0, 0, 0)
 
static constexpr double shapeFunMBTET1At000 = N_MBTET1(0, 0, 0)
 
static constexpr double shapeFunMBTET2At000 = N_MBTET2(0, 0, 0)
 
static constexpr double shapeFunMBTET3At000 = N_MBTET3(0, 0, 0)
 
static constexpr double shapeFunMBTET0AtOneThird
 
static constexpr double shapeFunMBTET1AtOneThird
 
static constexpr double shapeFunMBTET2AtOneThird
 
static constexpr double shapeFunMBTET3AtOneThird
 
static constexpr std::array< double, 4 > shapeFunMBTETAt000
 Array of shape function at zero local point on reference element.
 
static constexpr std::array< double, 4 > shapeFunMBTETAtOneThird
 Array of shape function at center on reference element.
 
MoFEMErrorCode minTetsQuality (const Range &tets, double &min_quality, Tag th=nullptr, boost::function< double(double, double)> f=[](double a, double b) -> double { return std::min(a, b);})
 calculate minimal quality of tetrahedra in range
 
MoFEMErrorCode getTetsWithQuality (Range &out_tets, const Range &tets, Tag th=nullptr, boost::function< bool(double)> f=[](double q) -> bool { if(q<=0) return true;else return false;})
 Get the Tets With Quality.
 
MoFEMErrorCode writeTetsWithQuality (const char *file_name, const char *file_type, const char *options, const Range &tets, Tag th=nullptr, boost::function< bool(double)> f=[](double q) -> bool { if(q<=0) return true;else return false;})
 Write file with tetrahedral of given quality.
 
MoFEMErrorCode getTriNormal (const EntityHandle tri, double *normal) const
 Get triangle normal.
 
double getTriArea (const EntityHandle tri) const
 Get triangle area.
 
double getEdgeLength (const EntityHandle edge)
 Get edge length.
 
MoFEMErrorCode findMinDistanceFromTheEdges (const double *v_ptr, const int nb, Range edges, double *min_dist_ptr, double *o_ptr=nullptr, EntityHandle *o_segments=nullptr) const
 Find minimal distance to edges.
 
static double volumeLengthQuality (const double *coords)
 Calculate tetrahedron volume length quality.
 
static double tetVolume (const double *coords)
 Calculate volume of tetrahedron.
 
static double shapeFunMBEDGE0 (const double x)
 
static double shapeFunMBEDGE1 (const double x)
 
template<int LDB = 1>
static MoFEMErrorCode shapeFunMBEDGE (double *shape, const double *ksi, const int nb)
 Calculate shape functions on edge.
 
static double shapeFunMBTRI0 (const double x, const double y)
 
static double shapeFunMBTRI1 (const double x, const double y)
 
static double shapeFunMBTRI2 (const double x, const double y)
 
template<int LDB = 1>
static MoFEMErrorCode shapeFunMBTRI (double *shape, const double *ksi, const double *eta, const int nb)
 Calculate shape functions on triangle.
 
static double shapeFunMBTET0 (const double x, const double y, const double z)
 
static double shapeFunMBTET1 (const double x, const double y, const double z)
 
static double shapeFunMBTET2 (const double x, const double y, const double z)
 
static double shapeFunMBTET3 (const double x, const double y, const double z)
 
template<int LDB = 1>
static MoFEMErrorCode shapeFunMBTET (double *shape, const double *ksi, const double *eta, const double *zeta, const double nb)
 Calculate shape functions on tetrahedron.
 
static MoFEMErrorCode getLocalCoordinatesOnReferenceFourNodeTet (const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
 Get the local coordinates on reference four node tet object.
 
static MoFEMErrorCode getLocalCoordinatesOnReferenceThreeNodeTri (const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
 Get the local coordinates on reference three node tri object.
 
static MoFEMErrorCode getLocalCoordinatesOnReferenceThreeNodeTri (const double *elem_coords, const std::complex< double > *glob_coords, const int nb_nodes, std::complex< double > *local_coords)
 Get the local coordinates on reference three node tri object.
 
static DEPRECATED MoFEMErrorCode getLocalCoordinatesOnReferenceTriNodeTri (const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
 
static MoFEMErrorCode getLocalCoordinatesOnReferenceEdgeNodeEdge (const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
 Get the local coordinates on reference four node tet object.
 
static MoFEMErrorCode checkIfPointIsInTet (const double tet_coords[], const double global_coord[], const double tol, bool &result)
 Check of point is in tetrahedral.
 
static MoFEMErrorCode getTriNormal (const double *coords, double *normal, double *d_normal=nullptr)
 Get the Tri Normal objectGet triangle normal.
 
static double getEdgeLength (const double *edge_coords)
 Get edge length.
 
static SEGMENT_MIN_DISTANCE minDistancePointFromOnSegment (const double *w_ptr, const double *v_ptr, const double *p_ptr, double *const t_ptr=nullptr)
 Find closet point on the segment from the point.
 
static SEGMENT_MIN_DISTANCE minDistanceFromSegments (const double *w_ptr, const double *v_ptr, const double *k_ptr, const double *l_ptr, double *const tvw_ptr=nullptr, double *const tlk_ptr=nullptr)
 Find points on two segments in closest distance.
 

Mesh refinement

using RefineTrianglesReturn
 
static constexpr std::array< int, 12 > uniformTriangleRefineTriangles
 
static RefineTrianglesReturn refineTriangle (int nb_levels)
 create uniform triangle mesh of refined elements
 
static MatrixDouble refineTriangleIntegrationPts (MatrixDouble pts, RefineTrianglesReturn refined)
 generate integration points for refined triangle mesh for last level
 
static MatrixDouble refineTriangleIntegrationPts (int rule, RefineTrianglesReturn refined)
 generate integration points for refined triangle mesh for last level
 

Debugging

MoFEMErrorCode checkVectorForNotANumber (const Problem *prb_ptr, const RowColData row_or_col, Vec v)
 Print all DOFs for which element of vector is not a number.
 

Detailed Description

Auxiliary tools.

Examples
EshelbianPlasticity.cpp.

Definition at line 19 of file Tools.hpp.

Member Typedef Documentation

◆ RefineTrianglesReturn

Initial value:
std::tuple<std::vector<double>, std::vector<int>, std::vector<int>>

Definition at line 647 of file Tools.hpp.

Member Enumeration Documentation

◆ SEGMENT_MIN_DISTANCE

Enumerator
SOLUTION_EXIST 
SEGMENT_ONE_IS_POINT 
SEGMENT_TWO_IS_POINT 
SEGMENT_TWO_AND_TWO_ARE_POINT 
NO_SOLUTION 

Definition at line 534 of file Tools.hpp.

Constructor & Destructor Documentation

◆ Tools()

MoFEM::Tools::Tools ( const MoFEM::Core & core)
inline

Definition at line 25 of file Tools.hpp.

25: cOre(const_cast<MoFEM::Core &>(core)) {}
Core (interface) class.
Definition Core.hpp:82
MoFEM::Core & cOre
Definition Tools.hpp:24

Member Function Documentation

◆ checkIfPointIsInTet()

MoFEMErrorCode MoFEM::Tools::checkIfPointIsInTet ( const double tet_coords[],
const double global_coord[],
const double tol,
bool & result )
static

Check of point is in tetrahedral.

Parameters
tet_coords
global_coord
tol
result
Returns
MoFEMErrorCode

Definition at line 287 of file Tools.cpp.

289 {
290 double loc_coord[] = {0, 0, 0};
291 double N[4], diffN[12];
293 CHKERR ShapeDiffMBTET(diffN);
294 CHKERR ShapeMBTET(N, &loc_coord[0], &loc_coord[1], &loc_coord[2], 1);
295 CHKERR ShapeMBTET_inverse(N, diffN, tet_coords, global_coord, loc_coord);
296 CHKERR ShapeMBTET(N, &loc_coord[0], &loc_coord[1], &loc_coord[2], 1);
297 result = true;
298 for (int n = 0; n != 4; ++n) {
299 if (N[n] < -tol || (N[n] - 1) > tol) {
300 result = false;
301 break;
302 }
303 }
305}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
PetscErrorCode ShapeDiffMBTET(double *diffN)
calculate derivatives of shape functions
Definition fem_tools.c:319
PetscErrorCode ShapeMBTET_inverse(double *N, double *diffN, const double *elem_coords, const double *glob_coords, double *loc_coords)
calculate local coordinates for given global coordinates
Definition fem_tools.c:335
PetscErrorCode ShapeMBTET(double *N, const double *G_X, const double *G_Y, const double *G_Z, int DIM)
calculate shape functions
Definition fem_tools.c:306
const double n
refractive index of diffusive medium
double tol
const int N
Definition speed_test.cpp:3

◆ checkVectorForNotANumber()

MoFEMErrorCode MoFEM::Tools::checkVectorForNotANumber ( const Problem * prb_ptr,
const RowColData row_or_col,
Vec v )

Print all DOFs for which element of vector is not a number.

Definition at line 307 of file Tools.cpp.

309 {
311 int loc_size;
312 CHKERR VecGetLocalSize(v, &loc_size);
313 int prb_loc_size = 0;
314 boost::shared_ptr<NumeredDofEntity_multiIndex> prb_dofs;
315 switch (row_or_col) {
316 case ROW:
317 prb_loc_size = prb_ptr->getNbLocalDofsRow();
318 prb_dofs = prb_ptr->getNumeredRowDofsPtr();
319 break;
320 case COL:
321 prb_loc_size = prb_ptr->getNbLocalDofsCol();
322 prb_dofs = prb_ptr->getNumeredColDofsPtr();
323 break;
324 break;
325 default:
326 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
327 "Wrong argument, row_or_col should be row or column");
328 }
329 if (loc_size != prb_loc_size) {
330 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
331 "Inconsistent size of vector and problem %d != %d", loc_size,
332 prb_loc_size);
333 }
334 const double *a;
335 CHKERR VecGetArrayRead(v, &a);
336 MPI_Comm comm = PetscObjectComm((PetscObject)v);
337 for (int ii = 0; ii != loc_size; ++ii) {
338 if (!boost::math::isfinite(a[ii])) {
339 NumeredDofEntityByLocalIdx::iterator dit =
340 prb_dofs->get<PetscLocalIdx_mi_tag>().find(ii);
341 std::ostringstream ss;
342 ss << "Not a number " << a[ii] << " on dof: " << endl
343 << **dit << endl
344 << endl;
345 PetscSynchronizedPrintf(comm, "%s", ss.str().c_str());
346 }
347 }
348 CHKERR VecRestoreArrayRead(v, &a);
349 PetscSynchronizedFlush(comm, PETSC_STDOUT);
351}
constexpr double a
@ COL
@ ROW
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
const double v
phase velocity of light in medium (cm/ns)

◆ findMinDistanceFromTheEdges()

MoFEMErrorCode MoFEM::Tools::findMinDistanceFromTheEdges ( const double * v_ptr,
const int nb,
Range edges,
double * min_dist_ptr,
double * o_ptr = nullptr,
EntityHandle * o_segments = nullptr ) const

Find minimal distance to edges.

Note
Finding only edges with have smaller distance than distance set on the input by min_dist_ptr
Parameters
v_ptrpoint coordinates
nbnb points
edgesrange of edges
min_dist_ptron return minimal distance, on input starting distance
o_ptrcoordinates of the point on edge
o_segmentsclosest segments
Returns
MoFEMErrorCode

Definition at line 538 of file Tools.cpp.

540 {
541 MoFEM::Interface &m_field = cOre;
542 moab::Interface &moab(m_field.get_moab());
544
545 FTensor::Index<'i', 3> i;
546
547 auto get_point = [i](auto &t_w, auto &t_delta, auto t) {
549 t = std::max(0., std::min(1., t));
550 t_p(i) = t_w(i) + t * t_delta(i);
551 return t_p;
552 };
553
554 auto get_distance = [i](auto &t_p, auto &t_n) {
555 FTensor::Tensor1<double, 3> t_dist_vector;
556 t_dist_vector(i) = t_p(i) - t_n(i);
557 return sqrt(t_dist_vector(i) * t_dist_vector(i));
558 };
559
560 for (auto e : edges) {
561
562 int num_nodes;
563 const EntityHandle *conn_fixed;
564 CHKERR moab.get_connectivity(e, conn_fixed, num_nodes, true);
565 VectorDouble6 coords_fixed(6);
566 CHKERR moab.get_coords(conn_fixed, num_nodes, &coords_fixed[0]);
568 &coords_fixed[0], &coords_fixed[1], &coords_fixed[2]);
570 &coords_fixed[3], &coords_fixed[4], &coords_fixed[5]);
571
572 FTensor::Tensor1<double, 3> t_edge_delta;
573 t_edge_delta(i) = t_f1(i) - t_f0(i);
574
576 v_ptr, v_ptr + 1, v_ptr + 2);
578 o_ptr, o_ptr + 1, o_ptr + 2);
579 FTensor::Tensor0<FTensor::PackPtr<double *, 1>> t_min_dist(min_dist_ptr);
580
581 EntityHandle *colsest_segment_it = nullptr;
582 if (o_segments)
583 colsest_segment_it = o_segments;
584
585 for (int n = 0; n != nb; ++n) {
586
587 double t;
588 if (Tools::minDistancePointFromOnSegment(&t_f0(0), &t_f1(0), &t_n(0),
590 auto t_p = get_point(t_f0, t_edge_delta, t);
591 auto dist_n = get_distance(t_p, t_n);
592 if (dist_n < t_min_dist || t_min_dist < 0) {
593 t_min_dist = dist_n;
594 if (o_ptr)
595 t_min_coords(i) = t_p(i);
596 if (o_segments)
597 *colsest_segment_it = e;
598 }
599 }
600
601 ++t_n;
602 ++t_min_dist;
603 if (o_ptr)
604 ++t_min_coords;
605 if (o_segments)
606 ++colsest_segment_it;
607 }
608 }
609
611}
FTensor::Index< 'i', SPACE_DIM > i
VectorBoundedArray< double, 6 > VectorDouble6
Definition Types.hpp:95
constexpr double t
plate stiffness
Definition plate.cpp:58
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
static SEGMENT_MIN_DISTANCE minDistancePointFromOnSegment(const double *w_ptr, const double *v_ptr, const double *p_ptr, double *const t_ptr=nullptr)
Find closet point on the segment from the point.
Definition Tools.cpp:449

◆ getEdgeLength() [1/2]

double MoFEM::Tools::getEdgeLength ( const double * edge_coords)
static

Get edge length.

Parameters
edge_coords
Returns
double

Definition at line 418 of file Tools.cpp.

418 {
420 edge_coords[2]);
422 edge_coords[5]);
423 FTensor::Index<'i', 3> i;
424 t_coords_n0(i) -= t_coords_n1(i);
425 return sqrt(t_coords_n0(i) * t_coords_n0(i));
426}
static const double edge_coords[6][6]

◆ getEdgeLength() [2/2]

double MoFEM::Tools::getEdgeLength ( const EntityHandle edge)

Get edge length.

Parameters
edge
Returns
double

Definition at line 428 of file Tools.cpp.

428 {
429 MoFEM::Interface &m_field = cOre;
430 moab::Interface &moab(m_field.get_moab());
431 auto get_edge_coords = [edge, &moab](double *const coords) {
433 if (type_from_handle(edge) != MBEDGE) {
434 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "Works only for edge");
435 }
436 const EntityHandle *conn;
437 int num_nodes;
438 CHKERR moab.get_connectivity(edge, conn, num_nodes, true);
439 CHKERR moab.get_coords(conn, 2, coords);
441 };
442 double coords[6];
443 ierr = get_edge_coords(coords);
444 CHKERRABORT(PETSC_COMM_SELF, ierr);
445 return getEdgeLength(coords);
446}
@ MOFEM_INVALID_DATA
Definition definitions.h:36
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
auto type_from_handle(const EntityHandle h)
get type from entity handle
static double getEdgeLength(const double *edge_coords)
Get edge length.
Definition Tools.cpp:418

◆ getLocalCoordinatesOnReferenceEdgeNodeEdge()

MoFEMErrorCode MoFEM::Tools::getLocalCoordinatesOnReferenceEdgeNodeEdge ( const double * elem_coords,
const double * glob_coords,
const int nb_nodes,
double * local_coords )
static

Get the local coordinates on reference four node tet object.

MatrixDouble elem_coords(4, 3);
// Set nodal coordinates
MatrixDouble global_coords(5, 3);
// Set global coordinates
MatrixDouble local_coords(global_coords.size1(), 3);
&elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
&local_coords(0, 0))
static MoFEMErrorCode getLocalCoordinatesOnReferenceEdgeNodeEdge(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Get the local coordinates on reference four node tet object.
Definition Tools.cpp:203
Parameters
elem_coordsGlobal element node coordinates
glob_coordsGlobal coordinates
nb_nodesNumber of points
local_coordsResult
Returns
MoFEMErrorCode

Definition at line 203 of file Tools.cpp.

205 {
206
207 FTensor::Index<'i', 3> i3;
209
211 &elem_coords[0], &elem_coords[3], &elem_coords[6]};
212
213 FTensor::Tensor1<double, 3> t_coords_at_0;
214 // Build matrix and get coordinates of zero point
215 // ii - global coordinates
217 for (auto ii : {0, 1, 2}) {
218 t_a(ii) = diffShapeFunMBEDGE[0] * t_elem_coords(0) +
219 diffShapeFunMBEDGE[1] * t_elem_coords(1);
220 t_coords_at_0(ii) = shapeFunMBEDGEAt00[0] * t_elem_coords(0) +
221 shapeFunMBEDGEAt00[1] * t_elem_coords(1);
222 ++t_elem_coords;
223 }
224
226 &global_coords[0], &global_coords[1], &global_coords[2]};
228 &local_coords[0]};
229
230 const double b = t_a(i3) * t_a(i3);
231 const double inv_b = 1 / b;
232
233 // Calculate right hand side
234 for (int ii = 0; ii != nb_nodes; ++ii) {
235 t_local_coords =
236 inv_b * (t_a(i3) * (t_global_coords(i3) - t_coords_at_0(i3)));
237 ++t_local_coords;
238 ++t_global_coords;
239 }
240
242}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
static constexpr std::array< double, 2 > shapeFunMBEDGEAt00
Array of shape function at zero local point on reference element.
Definition Tools.hpp:62
static constexpr std::array< double, 2 > diffShapeFunMBEDGE
Definition Tools.hpp:68

◆ getLocalCoordinatesOnReferenceFourNodeTet()

MoFEMErrorCode MoFEM::Tools::getLocalCoordinatesOnReferenceFourNodeTet ( const double * elem_coords,
const double * glob_coords,
const int nb_nodes,
double * local_coords )
static

Get the local coordinates on reference four node tet object.

MatrixDouble elem_coords(4, 3);
// Set nodal coordinates
MatrixDouble global_coords(5, 3);
// Set global coordinates
MatrixDouble local_coords(global_coords.size1(), 3);
&elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
&local_coords(0, 0))
static MoFEMErrorCode getLocalCoordinatesOnReferenceFourNodeTet(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Get the local coordinates on reference four node tet object.
Definition Tools.cpp:88
Parameters
elem_coordsGlobal element node coordinates
glob_coordsGlobal coordinates
nb_nodesNumber of points
local_coordsResult
Returns
MoFEMErrorCode

Definition at line 88 of file Tools.cpp.

90 {
91 FTensor::Index<'i', 4> i;
93
95 &elem_coords[0], &elem_coords[3], &elem_coords[6], &elem_coords[9]};
96
100 FTensor::Tensor1<double, 3> t_coords_at_0;
101
102 // Build matrix and get coordinates of zero point
103 // ii - global coordinates
104 // jj - local derivatives
105 MatrixDouble3by3 a(3, 3);
106 for (auto ii : {0, 1, 2}) {
110 for (auto jj : {0, 1, 2}) {
111 a(jj, ii) = t_diff(i) * t_elem_coords(i);
112 ++t_diff;
113 }
114 t_coords_at_0(ii) = t_n(i) * t_elem_coords(i);
115 ++t_elem_coords;
116 }
117
119 &global_coords[0], &global_coords[1], &global_coords[2]};
121 &local_coords[0], &local_coords[1], &local_coords[2]};
122
123 // Calculate right hand side
124 FTensor::Index<'j', 3> j;
125 for (int ii = 0; ii != nb_nodes; ++ii) {
126 t_local_coords(j) = t_global_coords(j) - t_coords_at_0(j);
127 ++t_local_coords;
128 ++t_global_coords;
129 }
130
131 // Solve problem
132 int IPIV[3];
133 int info = lapack_dgesv(3, nb_nodes, &a(0, 0), 3, IPIV, local_coords, 3);
134 if (info != 0)
135 SETERRQ(PETSC_COMM_SELF, 1, "info == %d", info);
136
138}
static __CLPK_integer lapack_dgesv(__CLPK_integer n, __CLPK_integer nrhs, __CLPK_doublereal *a, __CLPK_integer lda, __CLPK_integer *ipiv, __CLPK_doublereal *b, __CLPK_integer ldb)
FTensor::Index< 'j', 3 > j
MatrixBoundedArray< double, 9 > MatrixDouble3by3
Definition Types.hpp:105
static constexpr std::array< double, 4 > shapeFunMBTETAt000
Array of shape function at zero local point on reference element.
Definition Tools.hpp:330
static constexpr std::array< double, 12 > diffShapeFunMBTET
Definition Tools.hpp:271

◆ getLocalCoordinatesOnReferenceThreeNodeTri() [1/2]

MoFEMErrorCode MoFEM::Tools::getLocalCoordinatesOnReferenceThreeNodeTri ( const double * elem_coords,
const double * glob_coords,
const int nb_nodes,
double * local_coords )
static

Get the local coordinates on reference three node tri object.

MatrixDouble elem_coords(4, 3);
// Set nodal coordinates
MatrixDouble global_coords(5, 3);
// Set global coordinates
MatrixDouble local_coords(global_coords.size1(), 3);
&elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
&local_coords(0, 0))
Parameters
elem_coordsGlobal element node coordinates
glob_coordsGlobal coordinates
nb_nodesNumber of points
local_coordsResult
d_elem_coordsDerivative of local coordinates
d_global_coords
Returns
MoFEMErrorCode

Definition at line 188 of file Tools.cpp.

190 {
192 elem_coords, glob_coords, nb_nodes, local_coords);
193}
MoFEMErrorCode getLocalCoordinatesOnReferenceThreeNodeTriImpl(const T1 *elem_coords, const T2 *global_coords, const int nb_nodes, typename FTensor::promote< T1, T2 >::V *local_coords)
Definition Tools.cpp:141

◆ getLocalCoordinatesOnReferenceThreeNodeTri() [2/2]

MoFEMErrorCode MoFEM::Tools::getLocalCoordinatesOnReferenceThreeNodeTri ( const double * elem_coords,
const std::complex< double > * glob_coords,
const int nb_nodes,
std::complex< double > * local_coords )
static

Get the local coordinates on reference three node tri object.

MatrixDouble elem_coords(4, 3);
// Set nodal coordinates
MatrixDouble global_coords(5, 3);
// Set global coordinates
MatrixDouble local_coords(global_coords.size1(), 3);
&elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
&local_coords(0, 0))
Parameters
elem_coordsGlobal element node coordinates
glob_coordsGlobal coordinates
nb_nodesNumber of points
local_coordsResult
d_elem_coordsDerivative of local coordinates
d_global_coords
Returns
MoFEMErrorCode

Definition at line 195 of file Tools.cpp.

197 {
199 std::complex<double>>(
200 elem_coords, glob_coords, nb_nodes, local_coords);
201}

◆ getLocalCoordinatesOnReferenceTriNodeTri()

static DEPRECATED MoFEMErrorCode MoFEM::Tools::getLocalCoordinatesOnReferenceTriNodeTri ( const double * elem_coords,
const double * glob_coords,
const int nb_nodes,
double * local_coords )
inlinestatic
Deprecated
use getLocalCoordinatesOnReferenceThreeNodeTri

Definition at line 398 of file Tools.hpp.

400 {
401 return getLocalCoordinatesOnReferenceThreeNodeTri(elem_coords, glob_coords,
402 nb_nodes, local_coords);
403 }
static MoFEMErrorCode getLocalCoordinatesOnReferenceThreeNodeTri(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Get the local coordinates on reference three node tri object.
Definition Tools.cpp:188

◆ getTetsWithQuality()

MoFEMErrorCode MoFEM::Tools::getTetsWithQuality ( Range & out_tets,
const Range & tets,
Tag th = nullptr,
boost::function< bool(double)> f = [](double q) -> bool { if (q <= 0) return true; else return false; } )

Get the Tets With Quality.

Parameters
out_tets
tets
th
f
Returns
MoFEMErrorCode

Definition at line 244 of file Tools.cpp.

246 {
247 MoFEM::Interface &m_field = cOre;
248 moab::Interface &moab(m_field.get_moab());
250 Range to_write;
251 const EntityHandle *conn;
252 int num_nodes;
253 double coords[12];
254 for (auto tet : tets) {
255 CHKERR m_field.get_moab().get_connectivity(tet, conn, num_nodes, true);
256 if (th) {
257 CHKERR moab.tag_get_data(th, conn, num_nodes, coords);
258 } else {
259 CHKERR moab.get_coords(conn, num_nodes, coords);
260 }
261 double q = Tools::volumeLengthQuality(coords);
262 if (f(q)) {
263 out_tets.insert(tet);
264 }
265 }
267}
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition Tools.cpp:15

◆ getTriArea()

double MoFEM::Tools::getTriArea ( const EntityHandle tri) const

Get triangle area.

Parameters
tri
Returns
double

Definition at line 408 of file Tools.cpp.

408 {
409 if (type_from_handle(tri) != MBTRI) {
410 CHK_THROW_MESSAGE(MOFEM_INVALID_DATA, "Works only for triangle");
411 }
412 FTensor::Index<'i', 3> i;
414 CHK_THROW_MESSAGE(getTriNormal(tri, &t_normal(0)), "calculate area");
415 return sqrt(t_normal(i) * t_normal(i)) * 0.5;
416}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
static MoFEMErrorCode getTriNormal(const double *coords, double *normal, double *d_normal=nullptr)
Get the Tri Normal objectGet triangle normal.
Definition Tools.cpp:353

◆ getTriNormal() [1/2]

MoFEMErrorCode MoFEM::Tools::getTriNormal ( const double * coords,
double * normal,
double * d_normal = nullptr )
static

Get the Tri Normal objectGet triangle normal.

Parameters
coords
normal
d_normalderivative, if pointer is null, derivative is not calculated

In d_normal, format is 3rd rank tensor, first index is normal direction, second is node number, third is coordinate.

This is relation between tensor and pointer

auto t_d_normal = getFTensor3FromPtr<3, 3, 3>(d_normal);
FTensor::Tensor3< FTensor::PackPtr< double *, 27 >, 3, 3, 3 > getFTensor3FromPtr< 3, 3, 3 >(double *ptr)
Returns
MoFEMErrorCode
Examples
EshelbianPlasticity.cpp.

Definition at line 353 of file Tools.cpp.

354 {
356 FTensor::Index<'i', 3> i;
357 FTensor::Index<'j', 3> j;
358 FTensor::Index<'k', 3> k;
359 FTensor::Index<'l', 3> l;
360 FTensor::Index<'n', 3> n;
361 FTensor::Index<'m', 3> m;
362 FTensor::Index<'J', 2> J;
365 auto diff_ptr = Tools::diffShapeFunMBTRI.data();
366 auto t_diff_tensor = getFTensor2FromPtr<3, 2>(const_cast<double *>(diff_ptr));
367 auto t_coords = getFTensor2FromPtr<3, 3>(const_cast<double *>(coords));
369 t_tangent(i, J) = t_coords(n, i) * t_diff_tensor(n, J);
370 auto t_normal = getFTensor1FromPtr<3>(normal);
371 t_normal(j) =
372 (FTensor::levi_civita(i, j, k) * t_tangent(k, N0)) * t_tangent(i, N1);
373 if (d_normal) {
374 constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
376 t_d_coords(i, j, k, n) = t_kd(i, k) * t_kd(j, n);
378 t_d_tangent(i, k, l, J) = t_d_coords(n, i, k, l) * t_diff_tensor(n, J);
379 auto t_d_normal = getFTensor3FromPtr<3, 3, 3>(d_normal);
380 t_d_normal(j, m, n) = (FTensor::levi_civita(i, j, k) * t_tangent(i, N1)) *
381 t_d_tangent(k, m, n, N0)
382
383 +
384
385 (FTensor::levi_civita(i, j, k) * t_tangent(k, N0)) *
386 t_d_tangent(i, m, n, N1);
387 }
389}
Kronecker Delta class.
constexpr auto t_kd
FTensor::Index< 'J', DIM1 > J
Definition level_set.cpp:30
FTensor::Index< 'l', 3 > l
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
auto getFTensor2FromPtr(double *ptr)
Make Tensor2 from pointer.
FTensor::Tensor1< FTensor::PackPtr< double *, S >, DIM > getFTensor1FromPtr(double *ptr)
Make Tensor1 from pointer.
FTensor::Index< 'm', 3 > m
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition Tools.hpp:104

◆ getTriNormal() [2/2]

MoFEMErrorCode MoFEM::Tools::getTriNormal ( const EntityHandle tri,
double * normal ) const

Get triangle normal.

Parameters
tri
normal
Returns
MoFEMErrorCode

Definition at line 391 of file Tools.cpp.

392 {
393 MoFEM::Interface &m_field = cOre;
394 moab::Interface &moab(m_field.get_moab());
396 if (type_from_handle(tri) != MBTRI) {
397 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "Works only for triangle");
398 }
399 const EntityHandle *conn;
400 int num_nodes;
401 double coords[9];
402 CHKERR moab.get_connectivity(tri, conn, num_nodes, true);
403 CHKERR moab.get_coords(conn, num_nodes, coords);
404 CHKERR getTriNormal(coords, normal);
406}

◆ minDistanceFromSegments()

Tools::SEGMENT_MIN_DISTANCE MoFEM::Tools::minDistanceFromSegments ( const double * w_ptr,
const double * v_ptr,
const double * k_ptr,
const double * l_ptr,
double *const tvw_ptr = nullptr,
double *const tlk_ptr = nullptr )
static

Find points on two segments in closest distance.

Parameters
w_ptr
v_ptr
k_ptr
l_ptr
tvw_ptr
tlk_ptr
Returns
SEGMENT_MIN_DISTANCE
Note
If tvwk or tlk are outside bound [0,-1], it means that points are on the lines beyond segments, respectively for segment vw and lk.

Definition at line 472 of file Tools.cpp.

474 {
475
476 FTensor::Index<'i', 3> i;
477 FTensor::Tensor1<const double *, 3> t_w(&w_ptr[0], &w_ptr[1], &w_ptr[2]);
478 FTensor::Tensor1<const double *, 3> t_v(&v_ptr[0], &v_ptr[1], &v_ptr[2]);
479 FTensor::Tensor1<const double *, 3> t_k(&k_ptr[0], &k_ptr[1], &k_ptr[2]);
480 FTensor::Tensor1<const double *, 3> t_l(&l_ptr[0], &l_ptr[1], &l_ptr[2]);
481
482 // First segnent is a point
484 t_vw(i) = t_v(i) - t_w(i);
485 double dot_vw = t_vw(i) * t_vw(i);
486 if (std::abs(dot_vw) < std::numeric_limits<double>::epsilon()) {
487 if (tvw_ptr)
488 *tvw_ptr = 0;
489 if (minDistancePointFromOnSegment(k_ptr, l_ptr, w_ptr, tlk_ptr) ==
492 else
494 }
495
496 // Second segment is a point
498 t_lk(i) = t_l(i) - t_k(i);
499 double dot_lk = t_lk(i) * t_lk(i);
500 if (std::abs(dot_lk) < std::numeric_limits<double>::epsilon()) {
501 if (tlk_ptr)
502 *tlk_ptr = 0;
503 if (minDistancePointFromOnSegment(w_ptr, v_ptr, k_ptr, tvw_ptr) ==
506 else
508 }
509
510 const double a = t_vw(i) * t_vw(i);
511 const double b = -t_vw(i) * t_lk(i);
512 const double c = t_lk(i) * t_lk(i);
513
514 const double det = a * c - b * b;
515 if (std::abs(det) < std::numeric_limits<double>::epsilon()) {
516
517 return NO_SOLUTION;
518
519 } else {
520
522 t_wk(i) = t_w(i) - t_k(i);
523
524 const double ft0 = t_vw(i) * t_wk(i);
525 const double ft1 = -t_lk(i) * t_wk(i);
526 const double t0 = (ft1 * b - ft0 * c) / det;
527 const double t1 = (ft0 * b - ft1 * a) / det;
528
529 if (tvw_ptr)
530 *tvw_ptr = t0;
531 if (tlk_ptr)
532 *tlk_ptr = t1;
533
534 return SOLUTION_EXIST;
535 }
536}
const double c
speed of light (cm/ns)

◆ minDistancePointFromOnSegment()

Tools::SEGMENT_MIN_DISTANCE MoFEM::Tools::minDistancePointFromOnSegment ( const double * w_ptr,
const double * v_ptr,
const double * p_ptr,
double *const t_ptr = nullptr )
static

Find closet point on the segment from the point.

Parameters
w_ptrsegment first vertex coordinate
v_ptrsegment second vertex coordinate
p_ptrcoordinate of point
t_ptrdistance on the segment
Note
If t is outside bounds [ 0,-1 ] point is on the line point beyond segment.
double w[] = {-1, 0, 0};
double v[] = {1, 0, 0};
double p[] = {0, 1, 0};
double t;
CHKERR Toolas::minDistancePointFromOnSegment(w, v, p, &t);
double point_on_segment[3];
for (int i = 0; i != 3; ++i)
point_on_segment[i] = w[i] + t * (v[i] - w[i]);
Returns
SEGMENT_MIN_DISTANCE

Definition at line 449 of file Tools.cpp.

450 {
451 FTensor::Index<'i', 3> i;
452 FTensor::Tensor1<const double *, 3> t_w(&w_ptr[0], &w_ptr[1], &w_ptr[2]);
453 FTensor::Tensor1<const double *, 3> t_v(&v_ptr[0], &v_ptr[1], &v_ptr[2]);
454 FTensor::Tensor1<const double *, 3> t_p(&p_ptr[0], &p_ptr[1], &p_ptr[2]);
456 t_vw(i) = t_v(i) - t_w(i);
457 const double dot_vw = t_vw(i) * t_vw(i);
458 if (std::abs(dot_vw) < std::numeric_limits<double>::epsilon()) {
459 if (t_ptr)
460 *t_ptr = 0;
462 }
464 t_pw(i) = t_p(i) - t_w(i);
465 const double t = t_pw(i) * t_vw(i) / dot_vw;
466 if (t_ptr)
467 *t_ptr = t;
468 return SOLUTION_EXIST;
469}

◆ minTetsQuality()

MoFEMErrorCode MoFEM::Tools::minTetsQuality ( const Range & tets,
double & min_quality,
Tag th = nullptr,
boost::function< double(double, double)> f = [](double adouble b) -> double { return std::min(a, b); } )

calculate minimal quality of tetrahedra in range

Parameters
tetsrange
min_qualityminimal quality
Returns
error code

Definition at line 49 of file Tools.cpp.

50 {
51 MoFEM::Interface &m_field = cOre;
52 moab::Interface &moab(m_field.get_moab());
54 const EntityHandle *conn;
55 int num_nodes;
56 double coords[12];
57 for (auto tet : tets) {
58 CHKERR m_field.get_moab().get_connectivity(tet, conn, num_nodes, true);
59 if (th) {
60 CHKERR moab.tag_get_data(th, conn, num_nodes, coords);
61 } else {
62 CHKERR moab.get_coords(conn, num_nodes, coords);
63 }
64 double q = Tools::volumeLengthQuality(coords);
65 if (!std::isnormal(q))
66 q = -2;
67 min_quality = f(q, min_quality);
68 }
70}

◆ outerProductOfEdgeIntegrationPtsForHex()

MoFEMErrorCode MoFEM::Tools::outerProductOfEdgeIntegrationPtsForHex ( MatrixDouble & pts,
const int edge0,
const int edge1,
const int edge2 )
static

Definition at line 662 of file Tools.cpp.

664 {
666
667 auto check_rule_edge = [](int rule) {
669 if (rule < 0) {
670 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
671 "Wrong integration rule: %d", rule);
672 }
673 if (rule > QUAD_1D_TABLE_SIZE) {
674 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
675 "rule > quadrature order %d < %d", rule, QUAD_1D_TABLE_SIZE);
676 }
677 if (QUAD_1D_TABLE[rule]->dim != 1) {
678 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
679 }
680 if (QUAD_1D_TABLE[rule]->order < rule) {
681 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
682 "wrong order %d != %d", QUAD_1D_TABLE[rule]->order, rule);
683 }
685 };
686
687 CHKERR check_rule_edge(rule_ksi);
688 CHKERR check_rule_edge(rule_eta);
689 CHKERR check_rule_edge(rule_zeta);
690
691 const int nb_gauss_pts_ksi = QUAD_1D_TABLE[rule_ksi]->npoints;
692 const int nb_gauss_pts_eta = QUAD_1D_TABLE[rule_eta]->npoints;
693 const int nb_gauss_pts_zeta = QUAD_1D_TABLE[rule_zeta]->npoints;
694 gauss_pts.resize(4, nb_gauss_pts_ksi * nb_gauss_pts_eta * nb_gauss_pts_zeta,
695 false);
696
697 int gg = 0;
698 for (size_t i = 0; i != nb_gauss_pts_ksi; ++i) {
699 const double wi = QUAD_1D_TABLE[rule_ksi]->weights[i];
700 const double ksi = QUAD_1D_TABLE[rule_ksi]->points[2 * i + 1];
701 for (size_t j = 0; j != nb_gauss_pts_eta; ++j) {
702 const double wj = wi * QUAD_1D_TABLE[rule_eta]->weights[j];
703 const double eta = QUAD_1D_TABLE[rule_eta]->points[2 * j + 1];
704 for (size_t k = 0; k != nb_gauss_pts_zeta; ++k, ++gg) {
705 const double wk = wj * QUAD_1D_TABLE[rule_zeta]->weights[k];
706 const double zeta = QUAD_1D_TABLE[rule_zeta]->points[2 * k + 1];
707 gauss_pts(0, gg) = ksi;
708 gauss_pts(1, gg) = eta;
709 gauss_pts(2, gg) = zeta;
710 gauss_pts(3, gg) = wk;
711 }
712 }
713 }
714
715 if (gg != gauss_pts.size2())
716 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong size of matrix");
717
719}
constexpr int order
double eta
double zeta
Viscous hardening.
Definition plastic.cpp:130
static QUAD *const QUAD_1D_TABLE[]
Definition quad.h:164
#define QUAD_1D_TABLE_SIZE
Definition quad.h:163
double * points
Definition quad.h:30
double * weights
Definition quad.h:31
int npoints
Definition quad.h:29

◆ outerProductOfEdgeIntegrationPtsForQuad()

MoFEMErrorCode MoFEM::Tools::outerProductOfEdgeIntegrationPtsForQuad ( MatrixDouble & pts,
const int edge0,
const int edge1 )
static

Definition at line 613 of file Tools.cpp.

614 {
616
617 auto check_rule_edge = [](int rule) {
619 if (rule < 0) {
620 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
621 "Wrong integration rule: %d", rule);
622 }
623 if (rule > QUAD_1D_TABLE_SIZE) {
624 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
625 "rule > quadrature order %d < %d", rule, QUAD_1D_TABLE_SIZE);
626 }
627 if (QUAD_1D_TABLE[rule]->dim != 1) {
628 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
629 }
630 if (QUAD_1D_TABLE[rule]->order < rule) {
631 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
632 "wrong order %d != %d", QUAD_1D_TABLE[rule]->order, rule);
633 }
635 };
636
637 CHKERR check_rule_edge(rule_ksi);
638 CHKERR check_rule_edge(rule_eta);
639
640 const int nb_gauss_pts_ksi = QUAD_1D_TABLE[rule_ksi]->npoints;
641 const int nb_gauss_pts_eta = QUAD_1D_TABLE[rule_eta]->npoints;
642 gauss_pts.resize(3, nb_gauss_pts_ksi * nb_gauss_pts_eta, false);
643
644 int gg = 0;
645 for (size_t i = 0; i != nb_gauss_pts_ksi; ++i) {
646 const double wi = QUAD_1D_TABLE[rule_ksi]->weights[i];
647 const double ksi = (QUAD_1D_TABLE[rule_ksi]->points[2 * i + 1]);
648 for (size_t j = 0; j != nb_gauss_pts_eta; ++j, ++gg) {
649 const double wk = wi * QUAD_1D_TABLE[rule_eta]->weights[j];
650 const double eta = QUAD_1D_TABLE[rule_eta]->points[2 * j + 1];
651 gauss_pts(0, gg) = ksi;
652 gauss_pts(1, gg) = eta;
653 gauss_pts(2, gg) = wk;
654 }
655 }
656 if (gg != gauss_pts.size2())
657 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong size of matrix");
658
660}

◆ query_interface()

MoFEMErrorCode MoFEM::Tools::query_interface ( boost::typeindex::type_index type_index,
UnknownInterface ** iface ) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 9 of file Tools.cpp.

10 {
11 *iface = const_cast<Tools *>(this);
12 return 0;
13}
Tools(const MoFEM::Core &core)
Definition Tools.hpp:25

◆ refineTriangle()

std::tuple< std::vector< double >, std::vector< int >, std::vector< int > > MoFEM::Tools::refineTriangle ( int nb_levels)
static

create uniform triangle mesh of refined elements

Parameters
nb_levels
Returns
RefineTrianglesReturn

Definition at line 724 of file Tools.cpp.

724 {
725
726 std::vector<int> triangles{0, 1, 2, 3, 4, 5};
727 std::vector<double> nodes{
728
729 0., 0., // 0
730 1., 0., // 1
731 0., 1., // 2
732 0.5, 0., // 3
733 0.5, 0.5, // 4
734 0., 0.5 // 5
735
736 };
737 std::map<std::pair<int, int>, int> edges{
738 {{0, 1}, 3}, {{1, 2}, 4}, {{0, 2}, 5}};
739
740 auto add_edge = [&](auto a, auto b) {
741 if (a > b) {
742 std::swap(a, b);
743 }
744 auto it = edges.find(std::make_pair(a, b));
745 if (it == edges.end()) {
746 int e = 3 + edges.size();
747 edges[std::make_pair(a, b)] = e;
748 for (auto n : {0, 1}) {
749 nodes.push_back((nodes[2 * a + n] + nodes[2 * b + n]) / 2);
750 }
751#ifndef NDEBUG
752 if (e != nodes.size() / 2 - 1)
753 CHK_THROW_MESSAGE(MOFEM_DATA_INCONSISTENCY, "wrong node/edge index");
754#endif
755 return e;
756 }
757 return it->second;
758 };
759
760 auto add_triangle = [&](auto t) {
761 for (auto tt : {0, 1, 2, 3}) {
762 auto last = triangles.size() / 6;
763 for (auto n : {0, 1, 2}) {
764 // add triangle nodes
765 triangles.push_back(
766 triangles[6 * t + uniformTriangleRefineTriangles[3 * tt + n]]);
767 }
768 // add triangle edges
769 auto cycle = std::array<int, 4>{0, 1, 2, 0};
770 for (auto e : {0, 1, 2}) {
771 triangles.push_back(add_edge(triangles[6 * last + cycle[e]],
772 triangles[6 * last + cycle[e + 1]]));
773 }
774 }
775 };
776
777 std::vector<int> level_index{0, 1};
778 auto l = 0;
779 for (; l != nb_levels; ++l) {
780 auto first_tet = level_index[l];
781 auto nb_last_level_test = level_index.back() - level_index[l];
782 for (auto t = first_tet; t != (first_tet + nb_last_level_test); ++t) {
783 add_triangle(t);
784 }
785 level_index.push_back(triangles.size() / 6);
786 }
787
788 return std::make_tuple(nodes, triangles, level_index);
789}
constexpr int nb_levels
Definition level_set.cpp:58
static constexpr std::array< int, 12 > uniformTriangleRefineTriangles
Definition Tools.hpp:638

◆ refineTriangleIntegrationPts() [1/2]

MatrixDouble MoFEM::Tools::refineTriangleIntegrationPts ( int rule,
RefineTrianglesReturn refined )
static

generate integration points for refined triangle mesh for last level

Parameters
ruleGauss integration rule
refined
Returns
MatrixDouble

Definition at line 856 of file Tools.cpp.

856 {
857
858 MatrixDouble gauss_pts;
859
860 const size_t nb_gauss_pts = QUAD_2D_TABLE[rule]->npoints;
861 gauss_pts.resize(3, nb_gauss_pts, false);
862 cblas_dcopy(nb_gauss_pts, &QUAD_2D_TABLE[rule]->points[1], 3,
863 &gauss_pts(0, 0), 1);
864 cblas_dcopy(nb_gauss_pts, &QUAD_2D_TABLE[rule]->points[2], 3,
865 &gauss_pts(1, 0), 1);
866 cblas_dcopy(nb_gauss_pts, QUAD_2D_TABLE[rule]->weights, 1, &gauss_pts(2, 0),
867 1);
868
869 return Tools::refineTriangleIntegrationPts(gauss_pts, refined);
870}
static QUAD *const QUAD_2D_TABLE[]
Definition quad.h:175
static MatrixDouble refineTriangleIntegrationPts(MatrixDouble pts, RefineTrianglesReturn refined)
generate integration points for refined triangle mesh for last level
Definition Tools.cpp:791

◆ refineTriangleIntegrationPts() [2/2]

MatrixDouble MoFEM::Tools::refineTriangleIntegrationPts ( MatrixDouble pts,
RefineTrianglesReturn refined )
static

generate integration points for refined triangle mesh for last level

Parameters
pts
refined
Returns
MatrixDouble

Definition at line 791 of file Tools.cpp.

794 {
795
796 auto [nodes, triangles, level_index] = refined;
797
798 auto get_coords = [&](auto t) {
799 auto [nodes, triangles, level_index] = refined;
800 std::array<double, 9> ele_coords;
801 for (auto n : {0, 1, 2}) {
802 for (auto i : {0, 1}) {
803 ele_coords[3 * n + i] = nodes[2 * triangles[6 * t + n] + i];
804 }
805 ele_coords[3 * n + 2] = 0;
806 }
807 return ele_coords;
808 };
809
810 auto get_normal = [](auto &ele_coords) {
812 Tools::getTriNormal(ele_coords.data(), &t_normal(0));
813 return t_normal;
814 };
815
816 std::vector<double> ele_shape(3 * pts.size2());
817 shapeFunMBTRI<1>(&*ele_shape.begin(), &pts(0, 0), &pts(1, 0), pts.size2());
818
819 int nb_elems = level_index.back() - level_index[level_index.size() - 2];
820 MatrixDouble new_pts(3, pts.size2() * nb_elems);
821 FTensor::Tensor1<FTensor::PackPtr<double *, 1>, 2> t_gauss_pt{&new_pts(0, 0),
822 &new_pts(1, 0)};
824 &new_pts(2, 0)};
825
826 for (auto t = level_index[level_index.size() - 2]; t != level_index.back();
827 ++t) {
828
829 auto ele_coords = get_coords(t);
830 auto t_normal = get_normal(ele_coords);
831 auto area = t_normal.l2();
832
834 &ele_shape[0], &ele_shape[1], &ele_shape[2]};
835 FTensor::Tensor2<double, 3, 2> t_ele_coords{ele_coords[0], ele_coords[1],
836 ele_coords[3], ele_coords[4],
837 ele_coords[6], ele_coords[7]};
838
839 FTensor::Index<'i', 3> i;
840 FTensor::Index<'J', 2> J;
841
842 for (auto gg = 0; gg != pts.size2(); ++gg) {
843 t_gauss_pt(J) = t_ele_shape(i) * t_ele_coords(i, J);
844 t_gauss_weight = area * pts(2, gg);
845 ++t_gauss_pt;
846 ++t_gauss_weight;
847 ++t_ele_shape;
848 }
849
850 }
851
852 return new_pts;
853}
static MoFEMErrorCode shapeFunMBTRI(double *shape, const double *ksi, const double *eta, const int nb)
Calculate shape functions on triangle.
Definition Tools.hpp:716

◆ shapeFunMBEDGE()

template<int LDB>
MoFEMErrorCode MoFEM::Tools::shapeFunMBEDGE ( double * shape,
const double * ksi,
const int nb )
static

Calculate shape functions on edge.

Note
Template parameter is leading dimension of point coordinate arrays, such that \(ksi_{n+1} = ksi[n + LDB]\)
Template Parameters
1
Parameters
shapeshape functions
ksipointer to first local coordinates
nbnumber of points
Returns
MoFEMErrorCode

Definition at line 691 of file Tools.hpp.

692 {
694 for (int n = 0; n != nb; ++n) {
695 shape[0] = shapeFunMBEDGE0(*ksi);
696 shape[1] = shapeFunMBEDGE1(*ksi);
697 shape += 2;
698 ksi += LDB;
699 }
701}
static double shapeFunMBEDGE1(const double x)
Definition Tools.hpp:686
static double shapeFunMBEDGE0(const double x)
Definition Tools.hpp:682

◆ shapeFunMBEDGE0()

double MoFEM::Tools::shapeFunMBEDGE0 ( const double x)
inlinestatic

Definition at line 682 of file Tools.hpp.

682 {
683 return N_MBEDGE0(x);
684}
#define N_MBEDGE0(x)
edge shape function
Definition fem_tools.h:105

◆ shapeFunMBEDGE1()

double MoFEM::Tools::shapeFunMBEDGE1 ( const double x)
inlinestatic

Definition at line 686 of file Tools.hpp.

686 {
687 return N_MBEDGE1(x);
688}
#define N_MBEDGE1(x)
edge shape function
Definition fem_tools.h:106

◆ shapeFunMBTET()

template<int LDB>
MoFEMErrorCode MoFEM::Tools::shapeFunMBTET ( double * shape,
const double * ksi,
const double * eta,
const double * zeta,
const double nb )
static

Calculate shape functions on tetrahedron.

Note
Template parameter is leading dimension of point coordinate arrays, such that \(ksi_{n+1} = ksi[n + LDB]\)
Template Parameters
1
Parameters
shapeshape functions
ksipointer to first local coordinates
etapointer to second local coordinates
zetapointer to first third coordinates
nbnumber of points
Returns
MoFEMErrorCode

Definition at line 747 of file Tools.hpp.

749 {
751 for (int n = 0; n != nb; ++n) {
752 shape[0] = shapeFunMBTET0(*ksi, *eta, *zeta);
753 shape[1] = shapeFunMBTET1(*ksi, *eta, *zeta);
754 shape[2] = shapeFunMBTET2(*ksi, *eta, *zeta);
755 shape[3] = shapeFunMBTET3(*ksi, *eta, *zeta);
756 shape += 4;
757 ksi += LDB;
758 eta += LDB;
759 zeta += LDB;
760 }
762}
static double shapeFunMBTET0(const double x, const double y, const double z)
Definition Tools.hpp:730
static double shapeFunMBTET3(const double x, const double y, const double z)
Definition Tools.hpp:742
static double shapeFunMBTET1(const double x, const double y, const double z)
Definition Tools.hpp:734
static double shapeFunMBTET2(const double x, const double y, const double z)
Definition Tools.hpp:738

◆ shapeFunMBTET0()

double MoFEM::Tools::shapeFunMBTET0 ( const double x,
const double y,
const double z )
inlinestatic

Definition at line 730 of file Tools.hpp.

730 {
731 return N_MBTET0(x, y, z);
732}
#define N_MBTET0(x, y, z)
tetrahedral shape function
Definition fem_tools.h:28

◆ shapeFunMBTET1()

double MoFEM::Tools::shapeFunMBTET1 ( const double x,
const double y,
const double z )
inlinestatic

Definition at line 734 of file Tools.hpp.

734 {
735 return N_MBTET1(x, y, z);
736}
#define N_MBTET1(x, y, z)
tetrahedral shape function
Definition fem_tools.h:29

◆ shapeFunMBTET2()

double MoFEM::Tools::shapeFunMBTET2 ( const double x,
const double y,
const double z )
inlinestatic

Definition at line 738 of file Tools.hpp.

738 {
739 return N_MBTET2(x, y, z);
740}
#define N_MBTET2(x, y, z)
tetrahedral shape function
Definition fem_tools.h:30

◆ shapeFunMBTET3()

double MoFEM::Tools::shapeFunMBTET3 ( const double x,
const double y,
const double z )
inlinestatic

Definition at line 742 of file Tools.hpp.

742 {
743 return N_MBTET3(x, y, z);
744};
#define N_MBTET3(x, y, z)
tetrahedral shape function
Definition fem_tools.h:31

◆ shapeFunMBTRI()

template<int LDB>
MoFEMErrorCode MoFEM::Tools::shapeFunMBTRI ( double * shape,
const double * ksi,
const double * eta,
const int nb )
static

Calculate shape functions on triangle.

Note
Template parameter is leading dimension of point coordinate arrays, such that \(ksi_{n+1} = ksi[n + LDB]\)
Template Parameters
1
Parameters
shapeshape functions
ksipointer to first local coordinates
etapointer to second local coordinates
nbnumber of points
Returns
MoFEMErrorCode

Definition at line 716 of file Tools.hpp.

717 {
719 for (int n = 0; n != nb; ++n) {
720 shape[0] = shapeFunMBTRI0(*ksi, *eta);
721 shape[1] = shapeFunMBTRI1(*ksi, *eta);
722 shape[2] = shapeFunMBTRI2(*ksi, *eta);
723 shape += 3;
724 ksi += LDB;
725 eta += LDB;
726 }
728}
static double shapeFunMBTRI1(const double x, const double y)
Definition Tools.hpp:707
static double shapeFunMBTRI0(const double x, const double y)
Definition Tools.hpp:703
static double shapeFunMBTRI2(const double x, const double y)
Definition Tools.hpp:711

◆ shapeFunMBTRI0()

double MoFEM::Tools::shapeFunMBTRI0 ( const double x,
const double y )
inlinestatic

Definition at line 703 of file Tools.hpp.

703 {
704 return N_MBTRI0(x, y);
705}
#define N_MBTRI0(x, y)
triangle shape function
Definition fem_tools.h:46

◆ shapeFunMBTRI1()

double MoFEM::Tools::shapeFunMBTRI1 ( const double x,
const double y )
inlinestatic

Definition at line 707 of file Tools.hpp.

707 {
708 return N_MBTRI1(x, y);
709}
#define N_MBTRI1(x, y)
triangle shape function
Definition fem_tools.h:47

◆ shapeFunMBTRI2()

double MoFEM::Tools::shapeFunMBTRI2 ( const double x,
const double y )
inlinestatic

Definition at line 711 of file Tools.hpp.

711 {
712 return N_MBTRI2(x, y);
713}
#define N_MBTRI2(x, y)
triangle shape function
Definition fem_tools.h:48

◆ tetVolume()

double MoFEM::Tools::tetVolume ( const double * coords)
static

Calculate volume of tetrahedron.

Parameters
coords
Returns
double volume

Definition at line 30 of file Tools.cpp.

30 {
31 double diff_n[12];
32 ShapeDiffMBTET(diff_n);
33 FTensor::Tensor1<double *, 3> t_diff_n(&diff_n[0], &diff_n[1], &diff_n[2], 3);
34 FTensor::Tensor1<const double *, 3> t_coords(&coords[0], &coords[1],
35 &coords[2], 3);
37 FTensor::Index<'i', 3> i;
38 FTensor::Index<'j', 3> j;
39 jac(i, j) = 0;
40 for (int nn = 0; nn != 4; nn++) {
41 jac(i, j) += t_coords(i) * t_diff_n(j);
42 ++t_coords;
43 ++t_diff_n;
44 }
45 return determinantTensor3by3(jac) / 6.;
46}
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.

◆ volumeLengthQuality()

double MoFEM::Tools::volumeLengthQuality ( const double * coords)
static

Calculate tetrahedron volume length quality.

Parameters
coordstet coordinates
Returns
Volume-length quality

Definition at line 15 of file Tools.cpp.

15 {
16 double lrms = 0;
17 for (int dd = 0; dd != 3; dd++) {
18 lrms += pow(coords[0 * 3 + dd] - coords[1 * 3 + dd], 2) +
19 pow(coords[0 * 3 + dd] - coords[2 * 3 + dd], 2) +
20 pow(coords[0 * 3 + dd] - coords[3 * 3 + dd], 2) +
21 pow(coords[1 * 3 + dd] - coords[2 * 3 + dd], 2) +
22 pow(coords[1 * 3 + dd] - coords[3 * 3 + dd], 2) +
23 pow(coords[2 * 3 + dd] - coords[3 * 3 + dd], 2);
24 }
25 lrms = sqrt((1. / 6.) * lrms);
26 double volume = tetVolume(coords);
27 return 6. * sqrt(2.) * volume / pow(lrms, 3);
28}
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition ddTensor0.hpp:33
static double tetVolume(const double *coords)
Calculate volume of tetrahedron.
Definition Tools.cpp:30

◆ writeTetsWithQuality()

MoFEMErrorCode MoFEM::Tools::writeTetsWithQuality ( const char * file_name,
const char * file_type,
const char * options,
const Range & tets,
Tag th = nullptr,
boost::function< bool(double)> f = [](double q) -> bool { if (q <= 0) return true; else return false; } )

Write file with tetrahedral of given quality.

Parameters
file_name
file_type
options
tets
th
f
Returns
MoFEMErrorCode

Definition at line 269 of file Tools.cpp.

273 {
274 MoFEM::Interface &m_field = cOre;
275 moab::Interface &moab(m_field.get_moab());
277 Range out_tets;
278 CHKERR getTetsWithQuality(out_tets, tets, th, f);
279 EntityHandle meshset;
280 CHKERR moab.create_meshset(MESHSET_SET, meshset);
281 CHKERR moab.add_entities(meshset, out_tets);
282 CHKERR moab.write_file(file_name, file_type, options, &meshset, 1);
283 CHKERR moab.delete_entities(&meshset, 1);
285}
MoFEMErrorCode getTetsWithQuality(Range &out_tets, const Range &tets, Tag th=nullptr, boost::function< bool(double)> f=[](double q) -> bool { if(q<=0) return true;else return false;})
Get the Tets With Quality.
Definition Tools.cpp:244

Member Data Documentation

◆ cOre

MoFEM::Core& MoFEM::Tools::cOre

Definition at line 24 of file Tools.hpp.

◆ diffN_MBEDGE0x

double MoFEM::Tools::diffN_MBEDGE0x = diffN_MBEDGE0
staticconstexpr

Definition at line 65 of file Tools.hpp.

◆ diffN_MBEDGE1x

double MoFEM::Tools::diffN_MBEDGE1x = diffN_MBEDGE1
staticconstexpr

Definition at line 66 of file Tools.hpp.

◆ diffShapeFunMBEDGE

std::array< double, 2 > MoFEM::Tools::diffShapeFunMBEDGE
staticconstexpr
Initial value:
static constexpr double diffN_MBEDGE0x
Definition Tools.hpp:65
static constexpr double diffN_MBEDGE1x
Definition Tools.hpp:66

Definition at line 68 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter

std::array< double, 24 > MoFEM::Tools::diffShapeFunMBHEXAtCenter
staticconstexpr
Initial value:
= {
}
static constexpr double diffShapeFunMBHEXAtCenter6y
derivative of quad shape function
Definition Tools.hpp:201
static constexpr double diffShapeFunMBHEXAtCenter0y
derivative of HEX shape function
Definition Tools.hpp:165
static constexpr double diffShapeFunMBHEXAtCenter1x
derivative of HEX shape function
Definition Tools.hpp:169
static constexpr double diffShapeFunMBHEXAtCenter3z
derivative of quad shape function
Definition Tools.hpp:185
static constexpr double diffShapeFunMBHEXAtCenter7y
derivative of quad shape function
Definition Tools.hpp:207
static constexpr double diffShapeFunMBHEXAtCenter4y
derivative of quad shape function
Definition Tools.hpp:189
static constexpr double diffShapeFunMBHEXAtCenter0z
derivative of HEX shape function
Definition Tools.hpp:167
static constexpr double diffShapeFunMBHEXAtCenter1y
derivative of HEX shape function
Definition Tools.hpp:171
static constexpr double diffShapeFunMBHEXAtCenter5y
derivative of quad shape function
Definition Tools.hpp:195
static constexpr double diffShapeFunMBHEXAtCenter2y
derivative of HEX shape function
Definition Tools.hpp:177
static constexpr double diffShapeFunMBHEXAtCenter3x
derivative of HEX shape function
Definition Tools.hpp:181
static constexpr double diffShapeFunMBHEXAtCenter0x
derivative of HEX shape function
Definition Tools.hpp:163
static constexpr double diffShapeFunMBHEXAtCenter3y
derivative of quad shape function
Definition Tools.hpp:183
static constexpr double diffShapeFunMBHEXAtCenter2x
derivative of HEX shape function
Definition Tools.hpp:175
static constexpr double diffShapeFunMBHEXAtCenter7z
derivative of quad shape function
Definition Tools.hpp:209
static constexpr double diffShapeFunMBHEXAtCenter7x
derivative of HEX shape function
Definition Tools.hpp:205
static constexpr double diffShapeFunMBHEXAtCenter6x
derivative of HEX shape function
Definition Tools.hpp:199
static constexpr double diffShapeFunMBHEXAtCenter4x
derivative of HEX shape function
Definition Tools.hpp:187
static constexpr double diffShapeFunMBHEXAtCenter5z
derivative of quad shape function
Definition Tools.hpp:197
static constexpr double diffShapeFunMBHEXAtCenter2z
derivative of HEX shape function
Definition Tools.hpp:179
static constexpr double diffShapeFunMBHEXAtCenter4z
derivative of quad shape function
Definition Tools.hpp:191
static constexpr double diffShapeFunMBHEXAtCenter6z
derivative of quad shape function
Definition Tools.hpp:203
static constexpr double diffShapeFunMBHEXAtCenter1z
derivative of HEX shape function
Definition Tools.hpp:173
static constexpr double diffShapeFunMBHEXAtCenter5x
derivative of HEX shape function
Definition Tools.hpp:193

Definition at line 218 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter0x

double MoFEM::Tools::diffShapeFunMBHEXAtCenter0x
staticconstexpr
Initial value:
=
diffN_MBHEX0x(0.5, 0.5)
#define diffN_MBHEX0x(y, z)
Definition fem_tools.h:79

derivative of HEX shape function

Definition at line 163 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter0y

double MoFEM::Tools::diffShapeFunMBHEXAtCenter0y
staticconstexpr
Initial value:
=
diffN_MBHEX0y(0.5, 0.5)
#define diffN_MBHEX0y(x, z)
Definition fem_tools.h:87

derivative of HEX shape function

Definition at line 165 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter0z

double MoFEM::Tools::diffShapeFunMBHEXAtCenter0z
staticconstexpr
Initial value:
=
diffN_MBHEX0z(0.5, 0.5)
#define diffN_MBHEX0z(x, y)
Definition fem_tools.h:95

derivative of HEX shape function

Definition at line 167 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter1x

double MoFEM::Tools::diffShapeFunMBHEXAtCenter1x
staticconstexpr
Initial value:
=
diffN_MBHEX1x(0.5, 0.5)
#define diffN_MBHEX1x(y, z)
Definition fem_tools.h:80

derivative of HEX shape function

Definition at line 169 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter1y

double MoFEM::Tools::diffShapeFunMBHEXAtCenter1y
staticconstexpr
Initial value:
=
diffN_MBHEX1y(0.5, 0.5)
#define diffN_MBHEX1y(x, z)
Definition fem_tools.h:88

derivative of HEX shape function

Definition at line 171 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter1z

double MoFEM::Tools::diffShapeFunMBHEXAtCenter1z
staticconstexpr
Initial value:
=
diffN_MBHEX1z(0.5, 0.5)
#define diffN_MBHEX1z(x, y)
Definition fem_tools.h:96

derivative of HEX shape function

Definition at line 173 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter2x

double MoFEM::Tools::diffShapeFunMBHEXAtCenter2x
staticconstexpr
Initial value:
=
diffN_MBHEX2x(0.5, 0.5)
#define diffN_MBHEX2x(y, z)
Definition fem_tools.h:81

derivative of HEX shape function

Definition at line 175 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter2y

double MoFEM::Tools::diffShapeFunMBHEXAtCenter2y
staticconstexpr
Initial value:
=
diffN_MBHEX2y(0.5, 0.5)
#define diffN_MBHEX2y(x, z)
Definition fem_tools.h:89

derivative of HEX shape function

Definition at line 177 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter2z

double MoFEM::Tools::diffShapeFunMBHEXAtCenter2z
staticconstexpr
Initial value:
=
diffN_MBHEX2z(0.5, 0.5)
#define diffN_MBHEX2z(x, y)
Definition fem_tools.h:97

derivative of HEX shape function

Definition at line 179 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter3x

double MoFEM::Tools::diffShapeFunMBHEXAtCenter3x
staticconstexpr
Initial value:
=
diffN_MBHEX3x(0.5, 0.5)
#define diffN_MBHEX3x(y, z)
Definition fem_tools.h:82

derivative of HEX shape function

Definition at line 181 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter3y

double MoFEM::Tools::diffShapeFunMBHEXAtCenter3y
staticconstexpr
Initial value:
=
diffN_MBHEX3y(0.5, 0.5)
#define diffN_MBHEX3y(x, z)
Definition fem_tools.h:90

derivative of quad shape function

Definition at line 183 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter3z

double MoFEM::Tools::diffShapeFunMBHEXAtCenter3z
staticconstexpr
Initial value:
=
diffN_MBHEX3z(0.5, 0.5)
#define diffN_MBHEX3z(x, y)
Definition fem_tools.h:98

derivative of quad shape function

Definition at line 185 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter4x

double MoFEM::Tools::diffShapeFunMBHEXAtCenter4x
staticconstexpr
Initial value:
=
diffN_MBHEX4x(0.5, 0.5)
#define diffN_MBHEX4x(y, z)
Definition fem_tools.h:83

derivative of HEX shape function

Definition at line 187 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter4y

double MoFEM::Tools::diffShapeFunMBHEXAtCenter4y
staticconstexpr
Initial value:
=
diffN_MBHEX4y(0.5, 0.5)
#define diffN_MBHEX4y(x, z)
Definition fem_tools.h:91

derivative of quad shape function

Definition at line 189 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter4z

double MoFEM::Tools::diffShapeFunMBHEXAtCenter4z
staticconstexpr
Initial value:
=
diffN_MBHEX4z(0.5, 0.5)
#define diffN_MBHEX4z(x, y)
Definition fem_tools.h:99

derivative of quad shape function

Definition at line 191 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter5x

double MoFEM::Tools::diffShapeFunMBHEXAtCenter5x
staticconstexpr
Initial value:
=
diffN_MBHEX5x(0.5, 0.5)
#define diffN_MBHEX5x(y, z)
Definition fem_tools.h:84

derivative of HEX shape function

Definition at line 193 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter5y

double MoFEM::Tools::diffShapeFunMBHEXAtCenter5y
staticconstexpr
Initial value:
=
diffN_MBHEX5y(0.5, 0.5)
#define diffN_MBHEX5y(x, z)
Definition fem_tools.h:92

derivative of quad shape function

Definition at line 195 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter5z

double MoFEM::Tools::diffShapeFunMBHEXAtCenter5z
staticconstexpr
Initial value:
=
diffN_MBHEX5z(0.5, 0.5)
#define diffN_MBHEX5z(x, y)
Definition fem_tools.h:100

derivative of quad shape function

Definition at line 197 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter6x

double MoFEM::Tools::diffShapeFunMBHEXAtCenter6x
staticconstexpr
Initial value:
=
diffN_MBHEX6x(0.5, 0.5)
#define diffN_MBHEX6x(y, z)
Definition fem_tools.h:85

derivative of HEX shape function

Definition at line 199 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter6y

double MoFEM::Tools::diffShapeFunMBHEXAtCenter6y
staticconstexpr
Initial value:
=
diffN_MBHEX6y(0.5, 0.5)
#define diffN_MBHEX6y(x, z)
Definition fem_tools.h:93

derivative of quad shape function

Definition at line 201 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter6z

double MoFEM::Tools::diffShapeFunMBHEXAtCenter6z
staticconstexpr
Initial value:
=
diffN_MBHEX6z(0.5, 0.5)
#define diffN_MBHEX6z(x, y)
Definition fem_tools.h:101

derivative of quad shape function

Definition at line 203 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter7x

double MoFEM::Tools::diffShapeFunMBHEXAtCenter7x
staticconstexpr
Initial value:
=
diffN_MBHEX7x(0.5, 0.5)
#define diffN_MBHEX7x(y, z)
Definition fem_tools.h:86

derivative of HEX shape function

Definition at line 205 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter7y

double MoFEM::Tools::diffShapeFunMBHEXAtCenter7y
staticconstexpr
Initial value:
=
diffN_MBHEX7y(0.5, 0.5)
#define diffN_MBHEX7y(x, z)
Definition fem_tools.h:94

derivative of quad shape function

Definition at line 207 of file Tools.hpp.

◆ diffShapeFunMBHEXAtCenter7z

double MoFEM::Tools::diffShapeFunMBHEXAtCenter7z
staticconstexpr
Initial value:
=
diffN_MBHEX7z(0.5, 0.5)
#define diffN_MBHEX7z(x, y)
Definition fem_tools.h:102

derivative of quad shape function

Definition at line 209 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter

std::array< double, 8 > MoFEM::Tools::diffShapeFunMBQUADAtCenter
staticconstexpr
Initial value:
= {
static constexpr double diffShapeFunMBQUADAtCenter0x
derivative of quad shape function
Definition Tools.hpp:146
static constexpr double diffShapeFunMBQUADAtCenter2y
derivative of quad shape function
Definition Tools.hpp:156
static constexpr double diffShapeFunMBQUADAtCenter2x
derivative of quad shape function
Definition Tools.hpp:154
static constexpr double diffShapeFunMBQUADAtCenter3y
derivative of quad shape function
Definition Tools.hpp:160
static constexpr double diffShapeFunMBQUADAtCenter1y
derivative of quad shape function
Definition Tools.hpp:152
static constexpr double diffShapeFunMBQUADAtCenter3x
derivative of quad shape function
Definition Tools.hpp:158
static constexpr double diffShapeFunMBQUADAtCenter1x
derivative of quad shape function
Definition Tools.hpp:150
static constexpr double diffShapeFunMBQUADAtCenter0y
derivative of quad shape function
Definition Tools.hpp:148

Definition at line 212 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter0x

double MoFEM::Tools::diffShapeFunMBQUADAtCenter0x
staticconstexpr
Initial value:
=
#define diffN_MBQUAD0x(y)
Definition fem_tools.h:61

derivative of quad shape function

Definition at line 146 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter0y

double MoFEM::Tools::diffShapeFunMBQUADAtCenter0y
staticconstexpr
Initial value:
=
#define diffN_MBQUAD0y(x)
Definition fem_tools.h:62

derivative of quad shape function

Definition at line 148 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter1x

double MoFEM::Tools::diffShapeFunMBQUADAtCenter1x
staticconstexpr
Initial value:
=
#define diffN_MBQUAD1x(y)
Definition fem_tools.h:63

derivative of quad shape function

Definition at line 150 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter1y

double MoFEM::Tools::diffShapeFunMBQUADAtCenter1y
staticconstexpr
Initial value:
=
#define diffN_MBQUAD1y(x)
Definition fem_tools.h:64

derivative of quad shape function

Definition at line 152 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter2x

double MoFEM::Tools::diffShapeFunMBQUADAtCenter2x
staticconstexpr
Initial value:
=
#define diffN_MBQUAD2x(y)
Definition fem_tools.h:65

derivative of quad shape function

Definition at line 154 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter2y

double MoFEM::Tools::diffShapeFunMBQUADAtCenter2y
staticconstexpr
Initial value:
=
#define diffN_MBQUAD2y(x)
Definition fem_tools.h:66

derivative of quad shape function

Definition at line 156 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter3x

double MoFEM::Tools::diffShapeFunMBQUADAtCenter3x
staticconstexpr
Initial value:
=
#define diffN_MBQUAD3x(y)
Definition fem_tools.h:67

derivative of quad shape function

Definition at line 158 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter3y

double MoFEM::Tools::diffShapeFunMBQUADAtCenter3y
staticconstexpr
Initial value:
=
#define diffN_MBQUAD3y(x)
Definition fem_tools.h:68

derivative of quad shape function

Definition at line 160 of file Tools.hpp.

◆ diffShapeFunMBTET

std::array< double, 12 > MoFEM::Tools::diffShapeFunMBTET
staticconstexpr
Initial value:
= {
static constexpr double diffShapeFunMBTET0y
derivative of tetrahedral shape function
Definition Tools.hpp:248
static constexpr double diffShapeFunMBTET0z
derivative of tetrahedral shape function
Definition Tools.hpp:250
static constexpr double diffShapeFunMBTET3y
derivative of tetrahedral shape function
Definition Tools.hpp:266
static constexpr double diffShapeFunMBTET1z
derivative of tetrahedral shape function
Definition Tools.hpp:256
static constexpr double diffShapeFunMBTET2z
derivative of tetrahedral shape function
Definition Tools.hpp:262
static constexpr double diffShapeFunMBTET3x
derivative of tetrahedral shape function
Definition Tools.hpp:264
static constexpr double diffShapeFunMBTET0x
derivative of tetrahedral shape function
Definition Tools.hpp:246
static constexpr double diffShapeFunMBTET1y
derivative of tetrahedral shape function
Definition Tools.hpp:254
static constexpr double diffShapeFunMBTET3z
derivative of tetrahedral shape function
Definition Tools.hpp:268
static constexpr double diffShapeFunMBTET1x
derivative of tetrahedral shape function
Definition Tools.hpp:252
static constexpr double diffShapeFunMBTET2x
derivative of tetrahedral shape function
Definition Tools.hpp:258
static constexpr double diffShapeFunMBTET2y
derivative of tetrahedral shape function
Definition Tools.hpp:260

Definition at line 271 of file Tools.hpp.

◆ diffShapeFunMBTET0x

double MoFEM::Tools::diffShapeFunMBTET0x
staticconstexpr
Initial value:
=
#define diffN_MBTET0x
derivative of tetrahedral shape function
Definition fem_tools.h:32

derivative of tetrahedral shape function

Definition at line 246 of file Tools.hpp.

◆ diffShapeFunMBTET0y

double MoFEM::Tools::diffShapeFunMBTET0y
staticconstexpr
Initial value:
=
#define diffN_MBTET0y
derivative of tetrahedral shape function
Definition fem_tools.h:33

derivative of tetrahedral shape function

Definition at line 248 of file Tools.hpp.

◆ diffShapeFunMBTET0z

double MoFEM::Tools::diffShapeFunMBTET0z
staticconstexpr
Initial value:
=
#define diffN_MBTET0z
derivative of tetrahedral shape function
Definition fem_tools.h:34

derivative of tetrahedral shape function

Definition at line 250 of file Tools.hpp.

◆ diffShapeFunMBTET1x

double MoFEM::Tools::diffShapeFunMBTET1x
staticconstexpr
Initial value:
=
#define diffN_MBTET1x
derivative of tetrahedral shape function
Definition fem_tools.h:35

derivative of tetrahedral shape function

Definition at line 252 of file Tools.hpp.

◆ diffShapeFunMBTET1y

double MoFEM::Tools::diffShapeFunMBTET1y
staticconstexpr
Initial value:
=
#define diffN_MBTET1y
derivative of tetrahedral shape function
Definition fem_tools.h:36

derivative of tetrahedral shape function

Definition at line 254 of file Tools.hpp.

◆ diffShapeFunMBTET1z

double MoFEM::Tools::diffShapeFunMBTET1z
staticconstexpr
Initial value:
=
#define diffN_MBTET1z
derivative of tetrahedral shape function
Definition fem_tools.h:37

derivative of tetrahedral shape function

Definition at line 256 of file Tools.hpp.

◆ diffShapeFunMBTET2x

double MoFEM::Tools::diffShapeFunMBTET2x
staticconstexpr
Initial value:
=
#define diffN_MBTET2x
derivative of tetrahedral shape function
Definition fem_tools.h:38

derivative of tetrahedral shape function

Definition at line 258 of file Tools.hpp.

◆ diffShapeFunMBTET2y

double MoFEM::Tools::diffShapeFunMBTET2y
staticconstexpr
Initial value:
=
#define diffN_MBTET2y
derivative of tetrahedral shape function
Definition fem_tools.h:39

derivative of tetrahedral shape function

Definition at line 260 of file Tools.hpp.

◆ diffShapeFunMBTET2z

double MoFEM::Tools::diffShapeFunMBTET2z
staticconstexpr
Initial value:
=
#define diffN_MBTET2z
derivative of tetrahedral shape function
Definition fem_tools.h:40

derivative of tetrahedral shape function

Definition at line 262 of file Tools.hpp.

◆ diffShapeFunMBTET3x

double MoFEM::Tools::diffShapeFunMBTET3x
staticconstexpr
Initial value:
=
#define diffN_MBTET3x
derivative of tetrahedral shape function
Definition fem_tools.h:41

derivative of tetrahedral shape function

Definition at line 264 of file Tools.hpp.

◆ diffShapeFunMBTET3y

double MoFEM::Tools::diffShapeFunMBTET3y
staticconstexpr
Initial value:
=
#define diffN_MBTET3y
derivative of tetrahedral shape function
Definition fem_tools.h:42

derivative of tetrahedral shape function

Definition at line 266 of file Tools.hpp.

◆ diffShapeFunMBTET3z

double MoFEM::Tools::diffShapeFunMBTET3z
staticconstexpr
Initial value:
=
#define diffN_MBTET3z
derivative of tetrahedral shape function
Definition fem_tools.h:43

derivative of tetrahedral shape function

Definition at line 268 of file Tools.hpp.

◆ diffShapeFunMBTRI

std::array< double, 6 > MoFEM::Tools::diffShapeFunMBTRI
staticconstexpr
Initial value:
= {
static constexpr double diffShapeFunMBTRI1x
derivative of triangle shape function
Definition Tools.hpp:95
static constexpr double diffShapeFunMBTRI0x
derivative of triangle shape function
Definition Tools.hpp:91
static constexpr double diffShapeFunMBTRI2x
derivative of triangle shape function
Definition Tools.hpp:99
static constexpr double diffShapeFunMBTRI2y
derivative of triangle shape function
Definition Tools.hpp:101
static constexpr double diffShapeFunMBTRI0y
derivative of triangle shape function
Definition Tools.hpp:93
static constexpr double diffShapeFunMBTRI1y
derivative of triangle shape function
Definition Tools.hpp:97

Definition at line 104 of file Tools.hpp.

◆ diffShapeFunMBTRI0x

double MoFEM::Tools::diffShapeFunMBTRI0x
staticconstexpr
Initial value:
=
#define diffN_MBTRI0x
derivative of triangle shape function
Definition fem_tools.h:49

derivative of triangle shape function

Definition at line 91 of file Tools.hpp.

◆ diffShapeFunMBTRI0y

double MoFEM::Tools::diffShapeFunMBTRI0y
staticconstexpr
Initial value:
=
#define diffN_MBTRI0y
derivative of triangle shape function
Definition fem_tools.h:50

derivative of triangle shape function

Definition at line 93 of file Tools.hpp.

◆ diffShapeFunMBTRI1x

double MoFEM::Tools::diffShapeFunMBTRI1x
staticconstexpr
Initial value:
=
#define diffN_MBTRI1x
derivative of triangle shape function
Definition fem_tools.h:51

derivative of triangle shape function

Definition at line 95 of file Tools.hpp.

◆ diffShapeFunMBTRI1y

double MoFEM::Tools::diffShapeFunMBTRI1y
staticconstexpr
Initial value:
=
#define diffN_MBTRI1y
derivative of triangle shape function
Definition fem_tools.h:52

derivative of triangle shape function

Definition at line 97 of file Tools.hpp.

◆ diffShapeFunMBTRI2x

double MoFEM::Tools::diffShapeFunMBTRI2x
staticconstexpr
Initial value:
=
#define diffN_MBTRI2x
derivative of triangle shape function
Definition fem_tools.h:53

derivative of triangle shape function

Definition at line 99 of file Tools.hpp.

◆ diffShapeFunMBTRI2y

double MoFEM::Tools::diffShapeFunMBTRI2y
staticconstexpr
Initial value:
=
#define diffN_MBTRI2y
derivative of triangle shape function
Definition fem_tools.h:54

derivative of triangle shape function

Definition at line 101 of file Tools.hpp.

◆ shapeFunMBEDGE0At00

double MoFEM::Tools::shapeFunMBEDGE0At00 = N_MBEDGE0(0)
staticconstexpr

Definition at line 55 of file Tools.hpp.

◆ shapeFunMBEDGE1At00

double MoFEM::Tools::shapeFunMBEDGE1At00 = N_MBEDGE1(0)
staticconstexpr

Definition at line 56 of file Tools.hpp.

◆ shapeFunMBEDGEAt00

std::array< double, 2 > MoFEM::Tools::shapeFunMBEDGEAt00
staticconstexpr
Initial value:
= {
static constexpr double shapeFunMBEDGE1At00
Definition Tools.hpp:56
static constexpr double shapeFunMBEDGE0At00
Definition Tools.hpp:55

Array of shape function at zero local point on reference element.

Definition at line 62 of file Tools.hpp.

◆ shapeFunMBTET0At000

double MoFEM::Tools::shapeFunMBTET0At000 = N_MBTET0(0, 0, 0)
staticconstexpr

Definition at line 293 of file Tools.hpp.

◆ shapeFunMBTET0AtOneThird

double MoFEM::Tools::shapeFunMBTET0AtOneThird
staticconstexpr
Initial value:
=
N_MBTET0(1. / 3., 1. / 3., 1. / 3.)

Definition at line 298 of file Tools.hpp.

◆ shapeFunMBTET1At000

double MoFEM::Tools::shapeFunMBTET1At000 = N_MBTET1(0, 0, 0)
staticconstexpr

Definition at line 294 of file Tools.hpp.

◆ shapeFunMBTET1AtOneThird

double MoFEM::Tools::shapeFunMBTET1AtOneThird
staticconstexpr
Initial value:
=
N_MBTET1(1. / 3., 1. / 3., 1. / 3.)

Definition at line 300 of file Tools.hpp.

◆ shapeFunMBTET2At000

double MoFEM::Tools::shapeFunMBTET2At000 = N_MBTET2(0, 0, 0)
staticconstexpr

Definition at line 295 of file Tools.hpp.

◆ shapeFunMBTET2AtOneThird

double MoFEM::Tools::shapeFunMBTET2AtOneThird
staticconstexpr
Initial value:
=
N_MBTET2(1. / 3., 1. / 3., 1. / 3.)

Definition at line 302 of file Tools.hpp.

◆ shapeFunMBTET3At000

double MoFEM::Tools::shapeFunMBTET3At000 = N_MBTET3(0, 0, 0)
staticconstexpr

Definition at line 296 of file Tools.hpp.

◆ shapeFunMBTET3AtOneThird

double MoFEM::Tools::shapeFunMBTET3AtOneThird
staticconstexpr
Initial value:
=
N_MBTET3(1. / 3., 1. / 3., 1. / 3.)

Definition at line 304 of file Tools.hpp.

◆ shapeFunMBTETAt000

std::array< double, 4 > MoFEM::Tools::shapeFunMBTETAt000
staticconstexpr
Initial value:
= {
static constexpr double shapeFunMBTET1At000
Definition Tools.hpp:294
static constexpr double shapeFunMBTET0At000
Definition Tools.hpp:293
static constexpr double shapeFunMBTET3At000
Definition Tools.hpp:296
static constexpr double shapeFunMBTET2At000
Definition Tools.hpp:295

Array of shape function at zero local point on reference element.

Definition at line 330 of file Tools.hpp.

◆ shapeFunMBTETAtOneThird

std::array<double, 4> MoFEM::Tools::shapeFunMBTETAtOneThird
staticconstexpr
Initial value:
= {
static constexpr double shapeFunMBTET1AtOneThird
Definition Tools.hpp:300
static constexpr double shapeFunMBTET2AtOneThird
Definition Tools.hpp:302
static constexpr double shapeFunMBTET3AtOneThird
Definition Tools.hpp:304
static constexpr double shapeFunMBTET0AtOneThird
Definition Tools.hpp:298

Array of shape function at center on reference element.

Definition at line 338 of file Tools.hpp.

◆ shapeFunMBTRI0At00

double MoFEM::Tools::shapeFunMBTRI0At00 = N_MBTRI0(0, 0)
staticconstexpr

Definition at line 112 of file Tools.hpp.

◆ shapeFunMBTRI1At00

double MoFEM::Tools::shapeFunMBTRI1At00 = N_MBTRI1(0, 0)
staticconstexpr

Definition at line 113 of file Tools.hpp.

◆ shapeFunMBTRI2At00

double MoFEM::Tools::shapeFunMBTRI2At00 = N_MBTRI2(0, 0)
staticconstexpr

Definition at line 114 of file Tools.hpp.

◆ shapeFunMBTRIAt00

std::array< double, 3 > MoFEM::Tools::shapeFunMBTRIAt00
staticconstexpr
Initial value:
= {
static constexpr double shapeFunMBTRI1At00
Definition Tools.hpp:113
static constexpr double shapeFunMBTRI0At00
Definition Tools.hpp:112
static constexpr double shapeFunMBTRI2At00
Definition Tools.hpp:114

Array of shape function at zero local point on reference element.

Definition at line 120 of file Tools.hpp.

◆ uniformTriangleRefineTriangles

std::array< int, 12 > MoFEM::Tools::uniformTriangleRefineTriangles
staticconstexpr
Initial value:
= {
0, 3, 5,
3, 1, 4,
5, 4, 2,
5, 3, 4
}

Definition at line 638 of file Tools.hpp.

638 {
639
640 0, 3, 5, // 0
641 3, 1, 4, // 1
642 5, 4, 2, // 2
643 5, 3, 4 // 3
644
645 };

The documentation for this struct was generated from the following files: