v0.13.2
Loading...
Searching...
No Matches
Public Member Functions | Static Public Member Functions | Public Attributes | List of all members
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
virtual MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
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. More...
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version. More...
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version. More...
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version. More...
 

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. More...
 
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 More...
 
static constexpr double diffShapeFunMBTRI0y
 derivative of triangle shape function More...
 
static constexpr double diffShapeFunMBTRI1x
 derivative of triangle shape function More...
 
static constexpr double diffShapeFunMBTRI1y
 derivative of triangle shape function More...
 
static constexpr double diffShapeFunMBTRI2x
 derivative of triangle shape function More...
 
static constexpr double diffShapeFunMBTRI2y
 derivative of triangle shape function More...
 
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. More...
 
static constexpr double diffShapeFunMBQUADAtCenter0x
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBQUADAtCenter0y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBQUADAtCenter1x
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBQUADAtCenter1y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBQUADAtCenter2x
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBQUADAtCenter2y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBQUADAtCenter3x
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBQUADAtCenter3y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter0x
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter0y
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter0z
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter1x
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter1y
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter1z
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter2x
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter2y
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter2z
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter3x
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter3y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter3z
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter4x
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter4y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter4z
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter5x
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter5y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter5z
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter6x
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter6y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter6z
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter7x
 derivative of HEX shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter7y
 derivative of quad shape function More...
 
static constexpr double diffShapeFunMBHEXAtCenter7z
 derivative of quad shape function More...
 
static constexpr std::array< double, 8 > diffShapeFunMBQUADAtCenter
 
static constexpr std::array< double, 24 > diffShapeFunMBHEXAtCenter
 
static constexpr double diffShapeFunMBTET0x
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET0y
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET0z
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET1x
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET1y
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET1z
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET2x
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET2y
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET2z
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET3x
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET3y
 derivative of tetrahedral shape function More...
 
static constexpr double diffShapeFunMBTET3z
 derivative of tetrahedral shape function More...
 
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. More...
 
static constexpr std::array< double, 4 > shapeFunMBTETAtOneThird
 Array of shape function at center on reference element. More...
 
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 More...
 
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. More...
 
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. More...
 
MoFEMErrorCode getTriNormal (const EntityHandle tri, double *normal) const
 Get triangle normal. More...
 
double getTriArea (const EntityHandle tri) const
 Get triangle area. More...
 
double getEdgeLength (const EntityHandle edge)
 Get edge length. More...
 
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. More...
 
static double volumeLengthQuality (const double *coords)
 Calculate tetrahedron volume length quality. More...
 
static double tetVolume (const double *coords)
 Calculate volume of tetrahedron. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
static MoFEMErrorCode getLocalCoordinatesOnReferenceTriNodeTri (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. More...
 
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. More...
 
static MoFEMErrorCode checkIfPointIsInTet (const double tet_coords[], const double global_coord[], const double tol, bool &result)
 Check of point is in tetrahedral. More...
 
static MoFEMErrorCode getTriNormal (const double *coords, double *normal)
 Get the Tri Normal objectGet triangle normal. More...
 
static double getEdgeLength (const double *edge_coords)
 Get edge length. More...
 
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. More...
 
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. More...
 

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. More...
 

Detailed Description

Auxiliary tools.

Definition at line 19 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 507 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 286 of file Tools.cpp.

288 {
289 double loc_coord[] = {0, 0, 0};
290 double N[4], diffN[12];
292 CHKERR ShapeDiffMBTET(diffN);
293 CHKERR ShapeMBTET(N, &loc_coord[0], &loc_coord[1], &loc_coord[2], 1);
294 CHKERR ShapeMBTET_inverse(N, diffN, tet_coords, global_coord, loc_coord);
295 CHKERR ShapeMBTET(N, &loc_coord[0], &loc_coord[1], &loc_coord[2], 1);
296 result = true;
297 for (int n = 0; n != 4; ++n) {
298 if (N[n] < -tol || (N[n] - 1) > tol) {
299 result = false;
300 break;
301 }
302 }
304}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
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
FTensor::Index< 'n', SPACE_DIM > n
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 306 of file Tools.cpp.

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

506 {
507 MoFEM::Interface &m_field = cOre;
508 moab::Interface &moab(m_field.get_moab());
510
512
513 auto get_point = [i](auto &t_w, auto &t_delta, auto t) {
515 t = std::max(0., std::min(1., t));
516 t_p(i) = t_w(i) + t * t_delta(i);
517 return t_p;
518 };
519
520 auto get_distance = [i](auto &t_p, auto &t_n) {
521 FTensor::Tensor1<double, 3> t_dist_vector;
522 t_dist_vector(i) = t_p(i) - t_n(i);
523 return sqrt(t_dist_vector(i) * t_dist_vector(i));
524 };
525
526 for (auto e : edges) {
527
528 int num_nodes;
529 const EntityHandle *conn_fixed;
530 CHKERR moab.get_connectivity(e, conn_fixed, num_nodes, true);
531 VectorDouble6 coords_fixed(6);
532 CHKERR moab.get_coords(conn_fixed, num_nodes, &coords_fixed[0]);
534 &coords_fixed[0], &coords_fixed[1], &coords_fixed[2]);
536 &coords_fixed[3], &coords_fixed[4], &coords_fixed[5]);
537
538 FTensor::Tensor1<double, 3> t_edge_delta;
539 t_edge_delta(i) = t_f1(i) - t_f0(i);
540
542 v_ptr, v_ptr + 1, v_ptr + 2);
544 o_ptr, o_ptr + 1, o_ptr + 2);
545 FTensor::Tensor0<FTensor::PackPtr<double *, 1>> t_min_dist(min_dist_ptr);
546
547 EntityHandle *colsest_segment_it = nullptr;
548 if (o_segments)
549 colsest_segment_it = o_segments;
550
551 for (int n = 0; n != nb; ++n) {
552
553 double t;
554 if (Tools::minDistancePointFromOnSegment(&t_f0(0), &t_f1(0), &t_n(0),
556 auto t_p = get_point(t_f0, t_edge_delta, t);
557 auto dist_n = get_distance(t_p, t_n);
558 if (dist_n < t_min_dist || t_min_dist < 0) {
559 t_min_dist = dist_n;
560 if (o_ptr)
561 t_min_coords(i) = t_p(i);
562 if (o_segments)
563 *colsest_segment_it = e;
564 }
565 }
566
567 ++t_n;
568 ++t_min_dist;
569 if (o_ptr)
570 ++t_min_coords;
571 if (o_segments)
572 ++colsest_segment_it;
573 }
574 }
575
577}
FTensor::Index< 'i', SPACE_DIM > i
VectorBoundedArray< double, 6 > VectorDouble6
Definition: Types.hpp:95
constexpr double t
plate stiffness
Definition: plate.cpp:59
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:415

◆ getEdgeLength() [1/2]

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

Get edge length.

Parameters
edge_coords
Returns
double

Definition at line 384 of file Tools.cpp.

384 {
386 edge_coords[2]);
388 edge_coords[5]);
390 t_coords_n0(i) -= t_coords_n1(i);
391 return sqrt(t_coords_n0(i) * t_coords_n0(i));
392}
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 394 of file Tools.cpp.

394 {
395 MoFEM::Interface &m_field = cOre;
396 moab::Interface &moab(m_field.get_moab());
397 auto get_edge_coords = [edge, &moab](double *const coords) {
399 if (type_from_handle(edge) != MBEDGE) {
400 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "Works only for edge");
401 }
402 const EntityHandle *conn;
403 int num_nodes;
404 CHKERR moab.get_connectivity(edge, conn, num_nodes, true);
405 CHKERR moab.get_coords(conn, 2, coords);
407 };
408 double coords[6];
409 ierr = get_edge_coords(coords);
410 CHKERRABORT(PETSC_COMM_SELF, ierr);
411 return getEdgeLength(coords);
412}
@ MOFEM_INVALID_DATA
Definition: definitions.h:36
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1634
static double getEdgeLength(const double *edge_coords)
Get edge length.
Definition: Tools.cpp:384

◆ 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 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:90
Parameters
elem_coordsGlobal element node coordinates
glob_coordsGlobale coordinates
nb_nodesNumber of points
local_coordsResult
Returns
MoFEMErrorCode

Definition at line 200 of file Tools.cpp.

202 {
203
206
208 &elem_coords[0], &elem_coords[3], &elem_coords[6]};
209
212 FTensor::Tensor1<double, 3> t_coords_at_0;
213 // Build matrix and get coordinates of zero point
214 // ii - global coordinates
216 for (auto ii : {0, 1, 2}) {
217 t_a(ii) = diffShapeFunMBEDGE[0] * t_elem_coords(0) +
218 diffShapeFunMBEDGE[1] * t_elem_coords(1);
219 t_coords_at_0(ii) = shapeFunMBEDGEAt00[0] * t_elem_coords(0) +
220 shapeFunMBEDGEAt00[1] * t_elem_coords(1);
221 ++t_elem_coords;
222 }
223
225 &global_coords[0], &global_coords[1], &global_coords[2]};
227 &local_coords[0]};
228
229 const double b = t_a(i3) * t_a(i3);
230 const double inv_b = 1 / b;
231
232 // Calculate right hand side
233 for (int ii = 0; ii != nb_nodes; ++ii) {
234 t_local_coords =
235 inv_b * (t_a(i3) * (t_global_coords(i3) - t_coords_at_0(i3)));
236 ++t_local_coords;
237 ++t_global_coords;
238 }
239
241}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
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))
Parameters
elem_coordsGlobal element node coordinates
glob_coordsGlobale coordinates
nb_nodesNumber of points
local_coordsResult
Returns
MoFEMErrorCode

Definition at line 90 of file Tools.cpp.

92 {
95
97 &elem_coords[0], &elem_coords[3], &elem_coords[6], &elem_coords[9]};
98
102 FTensor::Tensor1<double, 3> t_coords_at_0;
103
104 // Build matrix and get coordinates of zero point
105 // ii - global coordinates
106 // jj - local derivatives
107 MatrixDouble3by3 a(3, 3);
108 for (auto ii : {0, 1, 2}) {
112 for (auto jj : {0, 1, 2}) {
113 a(jj, ii) = t_diff(i) * t_elem_coords(i);
114 ++t_diff;
115 }
116 t_coords_at_0(ii) = t_n(i) * t_elem_coords(i);
117 ++t_elem_coords;
118 }
119
121 &global_coords[0], &global_coords[1], &global_coords[2]};
123 &local_coords[0], &local_coords[1], &local_coords[2]};
124
125 // Calculate right hand side
127 for (int ii = 0; ii != nb_nodes; ++ii) {
128 t_local_coords(j) = t_global_coords(j) - t_coords_at_0(j);
129 ++t_local_coords;
130 ++t_global_coords;
131 }
132
133 // Solve problem
134 int IPIV[3];
135 int info = lapack_dgesv(3, nb_nodes, &a(0, 0), 3, IPIV, local_coords, 3);
136 if (info != 0)
137 SETERRQ1(PETSC_COMM_SELF, 1, "info == %d", info);
138
140}
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)
Definition: lapack_wrap.h:176
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

◆ getLocalCoordinatesOnReferenceTriNodeTri()

MoFEMErrorCode MoFEM::Tools::getLocalCoordinatesOnReferenceTriNodeTri ( 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))
Parameters
elem_coordsGlobal element node coordinates
glob_coordsGlobale coordinates
nb_nodesNumber of points
local_coordsResult
Returns
MoFEMErrorCode

Definition at line 142 of file Tools.cpp.

144 {
145
151
153 &elem_coords[0], &elem_coords[3], &elem_coords[6]};
154
157 FTensor::Tensor1<double, 3> t_coords_at_0;
158
159 // Build matrix and get coordinates of zero point
160 // ii - global coordinates
161 // jj - local derivatives
163 for (auto ii : {0, 1, 2}) {
166 for (auto jj : {0, 1}) {
167 t_a(jj, ii) = t_diff(i3) * t_elem_coords(i3);
168 ++t_diff;
169 }
170 t_coords_at_0(ii) = t_n(i3) * t_elem_coords(i3);
171 ++t_elem_coords;
172 }
173
175 &global_coords[0], &global_coords[1], &global_coords[2]};
177 &local_coords[0], &local_coords[1]};
178
180 t_b(K2, L2) = t_a(K2, j3) ^ t_a(L2, j3);
181
182 // Solve problem
183 const auto inv_det = 1. / (t_b(0, 0) * t_b(1, 1) - t_b(0, 1) * t_b(1, 0));
184 t_inv_b(0, 0) = t_b(1, 1) * inv_det;
185 t_inv_b(0, 1) = -t_b(0, 1) * inv_det;
186 t_inv_b(1, 1) = t_b(0, 0) * inv_det;
187
188 // Calculate right hand side
189 for (int ii = 0; ii != nb_nodes; ++ii) {
190 t_local_coords(L2) =
191 t_inv_b(L2, K2) *
192 (t_a(K2, j3) * (t_global_coords(j3) - t_coords_at_0(j3)));
193 ++t_local_coords;
194 ++t_global_coords;
195 }
196
198}
@ L2
field with C-1 continuity
Definition: definitions.h:88
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition: Tools.hpp:104
static constexpr std::array< double, 3 > shapeFunMBTRIAt00
Array of shape function at zero local point on reference element.
Definition: Tools.hpp:120

◆ 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 243 of file Tools.cpp.

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

◆ getTriArea()

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

Get triangle area.

Parameters
tri
Returns
double

Definition at line 377 of file Tools.cpp.

377 {
380 CHK_THROW_MESSAGE(getTriNormal(tri, &t_normal(0)), "calculate area");
381 return sqrt(t_normal(i) * t_normal(i)) * 0.5;
382}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
static MoFEMErrorCode getTriNormal(const double *coords, double *normal)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:352

◆ getTriNormal() [1/2]

MoFEMErrorCode MoFEM::Tools::getTriNormal ( const double coords,
double normal 
)
static

Get the Tri Normal objectGet triangle normal.

Parameters
coords
normal
Returns
MoFEMErrorCode

Definition at line 352 of file Tools.cpp.

352 {
354 double diffN[6];
355 CHKERR ShapeDiffMBTRI(diffN);
356 CHKERR ShapeFaceNormalMBTRI(diffN, coords, normal);
358}
PetscErrorCode ShapeFaceNormalMBTRI(double *diffN, const double *coords, double *normal)
Definition: fem_tools.c:229
PetscErrorCode ShapeDiffMBTRI(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:194

◆ getTriNormal() [2/2]

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

Get triangle normal.

Parameters
tri
normal
Returns
MoFEMErrorCode

Definition at line 360 of file Tools.cpp.

361 {
362 MoFEM::Interface &m_field = cOre;
363 moab::Interface &moab(m_field.get_moab());
365 if (type_from_handle(tri) != MBTRI) {
366 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "Works only for triangle");
367 }
368 const EntityHandle *conn;
369 int num_nodes;
370 double coords[9];
371 CHKERR moab.get_connectivity(tri, conn, num_nodes, true);
372 CHKERR moab.get_coords(conn, num_nodes, coords);
373 CHKERR getTriNormal(coords, normal);
375}

◆ 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 438 of file Tools.cpp.

440 {
441
443 FTensor::Tensor1<const double *, 3> t_w(&w_ptr[0], &w_ptr[1], &w_ptr[2]);
444 FTensor::Tensor1<const double *, 3> t_v(&v_ptr[0], &v_ptr[1], &v_ptr[2]);
445 FTensor::Tensor1<const double *, 3> t_k(&k_ptr[0], &k_ptr[1], &k_ptr[2]);
446 FTensor::Tensor1<const double *, 3> t_l(&l_ptr[0], &l_ptr[1], &l_ptr[2]);
447
448 // First segnent is a point
450 t_vw(i) = t_v(i) - t_w(i);
451 double dot_vw = t_vw(i) * t_vw(i);
452 if (std::abs(dot_vw) < std::numeric_limits<double>::epsilon()) {
453 if (tvw_ptr)
454 *tvw_ptr = 0;
455 if (minDistancePointFromOnSegment(k_ptr, l_ptr, w_ptr, tlk_ptr) ==
458 else
460 }
461
462 // Second segment is a point
464 t_lk(i) = t_l(i) - t_k(i);
465 double dot_lk = t_lk(i) * t_lk(i);
466 if (std::abs(dot_lk) < std::numeric_limits<double>::epsilon()) {
467 if (tlk_ptr)
468 *tlk_ptr = 0;
469 if (minDistancePointFromOnSegment(w_ptr, v_ptr, k_ptr, tvw_ptr) ==
472 else
474 }
475
476 const double a = t_vw(i) * t_vw(i);
477 const double b = -t_vw(i) * t_lk(i);
478 const double c = t_lk(i) * t_lk(i);
479
480 const double det = a * c - b * b;
481 if (std::abs(det) < std::numeric_limits<double>::epsilon()) {
482
483 return NO_SOLUTION;
484
485 } else {
486
488 t_wk(i) = t_w(i) - t_k(i);
489
490 const double ft0 = t_vw(i) * t_wk(i);
491 const double ft1 = -t_lk(i) * t_wk(i);
492 const double t0 = (ft1 * b - ft0 * c) / det;
493 const double t1 = (ft0 * b - ft1 * a) / det;
494
495 if (tvw_ptr)
496 *tvw_ptr = t0;
497 if (tlk_ptr)
498 *tlk_ptr = t1;
499
500 return SOLUTION_EXIST;
501 }
502}
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]);
static Index< 'p', 3 > p
Returns
SEGMENT_MIN_DISTANCE

Definition at line 415 of file Tools.cpp.

416 {
418 FTensor::Tensor1<const double *, 3> t_w(&w_ptr[0], &w_ptr[1], &w_ptr[2]);
419 FTensor::Tensor1<const double *, 3> t_v(&v_ptr[0], &v_ptr[1], &v_ptr[2]);
420 FTensor::Tensor1<const double *, 3> t_p(&p_ptr[0], &p_ptr[1], &p_ptr[2]);
422 t_vw(i) = t_v(i) - t_w(i);
423 const double dot_vw = t_vw(i) * t_vw(i);
424 if (std::abs(dot_vw) < std::numeric_limits<double>::epsilon()) {
425 if (t_ptr)
426 *t_ptr = 0;
428 }
430 t_pw(i) = t_p(i) - t_w(i);
431 const double t = t_pw(i) * t_vw(i) / dot_vw;
432 if (t_ptr)
433 *t_ptr = t;
434 return SOLUTION_EXIST;
435}

◆ 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_qualitymimimal quality
Returns
error code

Definition at line 51 of file Tools.cpp.

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

◆ outerProductOfEdgeIntegrationPtsForHex()

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

Definition at line 628 of file Tools.cpp.

630 {
632
633 auto check_rule_edge = [](int rule) {
635 if (rule < 0) {
636 SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
637 "Wrong integration rule: %d", rule);
638 }
639 if (rule > QUAD_1D_TABLE_SIZE) {
640 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
641 "rule > quadrature order %d < %d", rule, QUAD_1D_TABLE_SIZE);
642 }
643 if (QUAD_1D_TABLE[rule]->dim != 1) {
644 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
645 }
646 if (QUAD_1D_TABLE[rule]->order < rule) {
647 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
648 "wrong order %d != %d", QUAD_1D_TABLE[rule]->order, rule);
649 }
651 };
652
653 CHKERR check_rule_edge(rule_ksi);
654 CHKERR check_rule_edge(rule_eta);
655 CHKERR check_rule_edge(rule_zeta);
656
657 const int nb_gauss_pts_ksi = QUAD_1D_TABLE[rule_ksi]->npoints;
658 const int nb_gauss_pts_eta = QUAD_1D_TABLE[rule_eta]->npoints;
659 const int nb_gauss_pts_zeta = QUAD_1D_TABLE[rule_zeta]->npoints;
660 gauss_pts.resize(4, nb_gauss_pts_ksi * nb_gauss_pts_eta * nb_gauss_pts_zeta,
661 false);
662
663 int gg = 0;
664 for (size_t i = 0; i != nb_gauss_pts_ksi; ++i) {
665 const double wi = QUAD_1D_TABLE[rule_ksi]->weights[i];
666 const double ksi = QUAD_1D_TABLE[rule_ksi]->points[2 * i + 1];
667 for (size_t j = 0; j != nb_gauss_pts_eta; ++j) {
668 const double wj = wi * QUAD_1D_TABLE[rule_eta]->weights[j];
669 const double eta = QUAD_1D_TABLE[rule_eta]->points[2 * j + 1];
670 for (size_t k = 0; k != nb_gauss_pts_zeta; ++k, ++gg) {
671 const double wk = wj * QUAD_1D_TABLE[rule_zeta]->weights[k];
672 const double zeta = QUAD_1D_TABLE[rule_zeta]->points[2 * k + 1];
673 gauss_pts(0, gg) = ksi;
674 gauss_pts(1, gg) = eta;
675 gauss_pts(2, gg) = zeta;
676 gauss_pts(3, gg) = wk;
677 }
678 }
679 }
680
681 if (gg != gauss_pts.size2())
682 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong size of matrix");
683
685}
const int dim
constexpr double eta
FTensor::Index< 'k', 3 > k
double zeta
Definition: plastic.cpp:107
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 579 of file Tools.cpp.

580 {
582
583 auto check_rule_edge = [](int rule) {
585 if (rule < 0) {
586 SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
587 "Wrong integration rule: %d", rule);
588 }
589 if (rule > QUAD_1D_TABLE_SIZE) {
590 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
591 "rule > quadrature order %d < %d", rule, QUAD_1D_TABLE_SIZE);
592 }
593 if (QUAD_1D_TABLE[rule]->dim != 1) {
594 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
595 }
596 if (QUAD_1D_TABLE[rule]->order < rule) {
597 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
598 "wrong order %d != %d", QUAD_1D_TABLE[rule]->order, rule);
599 }
601 };
602
603 CHKERR check_rule_edge(rule_ksi);
604 CHKERR check_rule_edge(rule_eta);
605
606 const int nb_gauss_pts_ksi = QUAD_1D_TABLE[rule_ksi]->npoints;
607 const int nb_gauss_pts_eta = QUAD_1D_TABLE[rule_eta]->npoints;
608 gauss_pts.resize(3, nb_gauss_pts_ksi * nb_gauss_pts_eta, false);
609
610 int gg = 0;
611 for (size_t i = 0; i != nb_gauss_pts_ksi; ++i) {
612 const double wi = QUAD_1D_TABLE[rule_ksi]->weights[i];
613 const double ksi = (QUAD_1D_TABLE[rule_ksi]->points[2 * i + 1]);
614 for (size_t j = 0; j != nb_gauss_pts_eta; ++j, ++gg) {
615 const double wk = wi * QUAD_1D_TABLE[rule_eta]->weights[j];
616 const double eta = QUAD_1D_TABLE[rule_eta]->points[2 * j + 1];
617 gauss_pts(0, gg) = ksi;
618 gauss_pts(1, gg) = eta;
619 gauss_pts(2, gg) = wk;
620 }
621 }
622 if (gg != gauss_pts.size2())
623 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong size of matrix");
624
626}

◆ query_interface()

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

Implements MoFEM::UnknownInterface.

Definition at line 11 of file Tools.cpp.

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

◆ 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 617 of file Tools.hpp.

618 {
620 for (int n = 0; n != nb; ++n) {
621 shape[0] = shapeFunMBEDGE0(*ksi);
622 shape[1] = shapeFunMBEDGE1(*ksi);
623 shape += 2;
624 ksi += LDB;
625 }
627}
static double shapeFunMBEDGE1(const double x)
Definition: Tools.hpp:612
static double shapeFunMBEDGE0(const double x)
Definition: Tools.hpp:608

◆ shapeFunMBEDGE0()

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

Definition at line 608 of file Tools.hpp.

608 {
609 return N_MBEDGE0(x);
610}
#define N_MBEDGE0(x)
edge shape function
Definition: fem_tools.h:105

◆ shapeFunMBEDGE1()

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

Definition at line 612 of file Tools.hpp.

612 {
613 return N_MBEDGE1(x);
614}
#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 673 of file Tools.hpp.

675 {
677 for (int n = 0; n != nb; ++n) {
678 shape[0] = shapeFunMBTET0(*ksi, *eta, *zeta);
679 shape[1] = shapeFunMBTET1(*ksi, *eta, *zeta);
680 shape[2] = shapeFunMBTET2(*ksi, *eta, *zeta);
681 shape[3] = shapeFunMBTET3(*ksi, *eta, *zeta);
682 shape += 4;
683 ksi += LDB;
684 eta += LDB;
685 zeta += LDB;
686 }
688}
static double shapeFunMBTET0(const double x, const double y, const double z)
Definition: Tools.hpp:656
static double shapeFunMBTET3(const double x, const double y, const double z)
Definition: Tools.hpp:668
static double shapeFunMBTET1(const double x, const double y, const double z)
Definition: Tools.hpp:660
static double shapeFunMBTET2(const double x, const double y, const double z)
Definition: Tools.hpp:664

◆ shapeFunMBTET0()

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

Definition at line 656 of file Tools.hpp.

656 {
657 return N_MBTET0(x, y, z);
658}
#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 660 of file Tools.hpp.

660 {
661 return N_MBTET1(x, y, z);
662}
#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 664 of file Tools.hpp.

664 {
665 return N_MBTET2(x, y, z);
666}
#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 668 of file Tools.hpp.

668 {
669 return N_MBTET3(x, y, z);
670};
#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 642 of file Tools.hpp.

643 {
645 for (int n = 0; n != nb; ++n) {
646 shape[0] = shapeFunMBTRI0(*ksi, *eta);
647 shape[1] = shapeFunMBTRI1(*ksi, *eta);
648 shape[2] = shapeFunMBTRI2(*ksi, *eta);
649 shape += 3;
650 ksi += LDB;
651 eta += LDB;
652 }
654}
static double shapeFunMBTRI1(const double x, const double y)
Definition: Tools.hpp:633
static double shapeFunMBTRI0(const double x, const double y)
Definition: Tools.hpp:629
static double shapeFunMBTRI2(const double x, const double y)
Definition: Tools.hpp:637

◆ shapeFunMBTRI0()

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

Definition at line 629 of file Tools.hpp.

629 {
630 return N_MBTRI0(x, y);
631}
#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 633 of file Tools.hpp.

633 {
634 return N_MBTRI1(x, y);
635}
#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 637 of file Tools.hpp.

637 {
638 return N_MBTRI2(x, y);
639}
#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 32 of file Tools.cpp.

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

◆ volumeLengthQuality()

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

Calculate tetrahedron volume length quality.

Parameters
coordstet coordinates
Returns
Volume-length quality

Definition at line 17 of file Tools.cpp.

17 {
18 double lrms = 0;
19 for (int dd = 0; dd != 3; dd++) {
20 lrms += pow(coords[0 * 3 + dd] - coords[1 * 3 + dd], 2) +
21 pow(coords[0 * 3 + dd] - coords[2 * 3 + dd], 2) +
22 pow(coords[0 * 3 + dd] - coords[3 * 3 + dd], 2) +
23 pow(coords[1 * 3 + dd] - coords[2 * 3 + dd], 2) +
24 pow(coords[1 * 3 + dd] - coords[3 * 3 + dd], 2) +
25 pow(coords[2 * 3 + dd] - coords[3 * 3 + dd], 2);
26 }
27 lrms = sqrt((1. / 6.) * lrms);
28 double volume = tetVolume(coords);
29 return 6. * sqrt(2.) * volume / pow(lrms, 3);
30}
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:32

◆ 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 268 of file Tools.cpp.

272 {
273 MoFEM::Interface &m_field = cOre;
274 moab::Interface &moab(m_field.get_moab());
276 Range out_tets;
277 CHKERR getTetsWithQuality(out_tets, tets, th, f);
278 EntityHandle meshset;
279 CHKERR moab.create_meshset(MESHSET_SET, meshset);
280 CHKERR moab.add_entities(meshset, out_tets);
281 CHKERR moab.write_file(file_name, file_type, options, &meshset, 1);
282 CHKERR moab.delete_entities(&meshset, 1);
284}
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:243

Member Data Documentation

◆ cOre

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

Definition at line 24 of file Tools.hpp.

◆ diffN_MBEDGE0x

constexpr double MoFEM::Tools::diffN_MBEDGE0x = diffN_MBEDGE0
staticconstexpr

Definition at line 65 of file Tools.hpp.

◆ diffN_MBEDGE1x

constexpr double MoFEM::Tools::diffN_MBEDGE1x = diffN_MBEDGE1
staticconstexpr

Definition at line 66 of file Tools.hpp.

◆ diffShapeFunMBEDGE

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

constexpr 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

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

Definition at line 55 of file Tools.hpp.

◆ shapeFunMBEDGE1At00

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

Definition at line 56 of file Tools.hpp.

◆ shapeFunMBEDGEAt00

constexpr 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

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

Definition at line 293 of file Tools.hpp.

◆ shapeFunMBTET0AtOneThird

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

Definition at line 298 of file Tools.hpp.

◆ shapeFunMBTET1At000

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

Definition at line 294 of file Tools.hpp.

◆ shapeFunMBTET1AtOneThird

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

Definition at line 300 of file Tools.hpp.

◆ shapeFunMBTET2At000

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

Definition at line 295 of file Tools.hpp.

◆ shapeFunMBTET2AtOneThird

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

Definition at line 302 of file Tools.hpp.

◆ shapeFunMBTET3At000

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

Definition at line 296 of file Tools.hpp.

◆ shapeFunMBTET3AtOneThird

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

Definition at line 304 of file Tools.hpp.

◆ shapeFunMBTETAt000

constexpr 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

constexpr 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

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

Definition at line 112 of file Tools.hpp.

◆ shapeFunMBTRI1At00

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

Definition at line 113 of file Tools.hpp.

◆ shapeFunMBTRI2At00

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

Definition at line 114 of file Tools.hpp.

◆ shapeFunMBTRIAt00

constexpr 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.


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