v0.8.23
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 (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 Tools (const MoFEM::Core &core)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of 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 ()
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

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

Additional Inherited Members

- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Auxiliary tools.

Definition at line 30 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 308 of file Tools.hpp.

Constructor & Destructor Documentation

◆ Tools()

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

Definition at line 36 of file Tools.hpp.

36 : cOre(const_cast<MoFEM::Core &>(core)) {}
MoFEM::Core & cOre
Definition: Tools.hpp:35

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

191  {
192  double loc_coord[] = {0, 0, 0};
193  double N[4], diffN[12];
195  CHKERR ShapeDiffMBTET(diffN);
196  CHKERR ShapeMBTET(N, &loc_coord[0], &loc_coord[1], &loc_coord[2], 1);
197  CHKERR ShapeMBTET_inverse(N, diffN, tet_coords, global_coord, loc_coord);
198  CHKERR ShapeMBTET(N, &loc_coord[0], &loc_coord[1], &loc_coord[2], 1);
199  result = true;
200  for (int n = 0; n != 4; ++n) {
201  if (N[n] < -tol || (N[n] - 1) > tol) {
202  result = false;
203  break;
204  }
205  }
207 }
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:318
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
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:347
PetscErrorCode ShapeDiffMBTET(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:331
double tol
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
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 209 of file Tools.cpp.

211  {
213  int loc_size;
214  CHKERR VecGetLocalSize(v, &loc_size);
215  int prb_loc_size = 0;
216  boost::shared_ptr<NumeredDofEntity_multiIndex> prb_dofs;
217  switch (row_or_col) {
218  case ROW:
219  prb_loc_size = prb_ptr->getNbLocalDofsRow();
220  prb_dofs = prb_ptr->getNumeredDofsRows();
221  break;
222  case COL:
223  prb_loc_size = prb_ptr->getNbLocalDofsCol();
224  prb_dofs = prb_ptr->getNumeredDofsCols();
225  break;
226  break;
227  default:
228  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
229  "Wrong argument, row_or_col should be row or column");
230  }
231  if (loc_size != prb_loc_size) {
232  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
233  "Inconsistent size of vector and problem %d != %d", loc_size,
234  prb_loc_size);
235  }
236  const double *a;
237  CHKERR VecGetArrayRead(v, &a);
238  MPI_Comm comm = PetscObjectComm((PetscObject)v);
239  for (int ii = 0; ii != loc_size; ++ii) {
240  if (!boost::math::isfinite(a[ii])) {
241  NumeredDofEntityByLocalIdx::iterator dit =
242  prb_dofs->get<PetscLocalIdx_mi_tag>().find(ii);
243  std::ostringstream ss;
244  ss << "Not a number " << a[ii] << " on dof: " << endl
245  << **dit << endl
246  << endl;
247  PetscSynchronizedPrintf(comm, "%s", ss.str().c_str());
248  }
249  }
250  CHKERR VecRestoreArrayRead(v, &a);
251  PetscSynchronizedFlush(comm, PETSC_STDOUT);
253 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

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

410  {
411  MoFEM::Interface &m_field = cOre;
412  moab::Interface &moab(m_field.get_moab());
414 
416 
417  auto get_point = [i](auto &t_w, auto &t_delta, auto t) {
419  t = std::max(0., std::min(1., t));
420  t_p(i) = t_w(i) + t * t_delta(i);
421  return t_p;
422  };
423 
424  auto get_distance = [i](auto &t_p, auto &t_n) {
425  FTensor::Tensor1<double, 3> t_dist_vector;
426  t_dist_vector(i) = t_p(i) - t_n(i);
427  return sqrt(t_dist_vector(i) * t_dist_vector(i));
428  };
429 
430  for (auto e : edges) {
431 
432  int num_nodes;
433  const EntityHandle *conn_fixed;
434  CHKERR moab.get_connectivity(e, conn_fixed, num_nodes, true);
435  VectorDouble6 coords_fixed(6);
436  CHKERR moab.get_coords(conn_fixed, num_nodes, &coords_fixed[0]);
438  &coords_fixed[0], &coords_fixed[1], &coords_fixed[2]);
440  &coords_fixed[3], &coords_fixed[4], &coords_fixed[5]);
441 
442  FTensor::Tensor1<double, 3> t_edge_delta;
443  t_edge_delta(i) = t_f1(i) - t_f0(i);
444 
446  v_ptr, v_ptr + 1, v_ptr + 2);
448  o_ptr, o_ptr + 1, o_ptr + 2);
449  FTensor::Tensor0<FTensor::PackPtr<double *, 1>> t_min_dist(min_dist_ptr);
450 
451  EntityHandle *colsest_segment_it = nullptr;
452  if (o_segments)
453  colsest_segment_it = o_segments;
454 
455  for (int n = 0; n != nb; ++n) {
456 
457  double t;
458  if (Tools::minDistancePointFromOnSegment(&t_f0(0), &t_f1(0), &t_n(0),
459  &t) == Tools::SOLUTION_EXIST) {
460  auto t_p = get_point(t_f0, t_edge_delta, t);
461  auto dist_n = get_distance(t_p, t_n);
462  if (dist_n < t_min_dist || t_min_dist < 0) {
463  t_min_dist = dist_n;
464  if (o_ptr)
465  t_min_coords(i) = t_p(i);
466  if (o_segments)
467  *colsest_segment_it = e;
468  }
469  }
470 
471  ++t_n;
472  ++t_min_dist;
473  if (o_ptr)
474  ++t_min_coords;
475  if (o_segments)
476  ++colsest_segment_it;
477  }
478  }
479 
481 }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
VectorBoundedArray< double, 6 > VectorDouble6
Definition: Types.hpp:88
#define CHKERR
Inline error check.
Definition: definitions.h:596
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:319
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MoFEM::Core & cOre
Definition: Tools.hpp:35

◆ getEdgeLength() [1/2]

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

Get edge length.

Parameters
edge_coords
Returns
double

Definition at line 288 of file Tools.cpp.

288  {
290  edge_coords[2]);
292  edge_coords[5]);
294  t_coords_n0(i) -= t_coords_n1(i);
295  return sqrt(t_coords_n0(i) * t_coords_n0(i));
296 }
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 298 of file Tools.cpp.

298  {
299  MoFEM::Interface &m_field = cOre;
300  moab::Interface &moab(m_field.get_moab());
301  auto get_edge_coords = [edge, &moab](double *const coords) {
303  if (moab.type_from_handle(edge) != MBEDGE) {
304  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "Works only for edge");
305  }
306  const EntityHandle *conn;
307  int num_nodes;
308  CHKERR moab.get_connectivity(edge, conn, num_nodes, true);
309  CHKERR moab.get_coords(conn, 2, coords);
311  };
312  double coords[6];
313  ierr = get_edge_coords(coords);
314  CHKERRABORT(PETSC_COMM_SELF, ierr);
315  return getEdgeLength(coords);
316 }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
static double getEdgeLength(const double *edge_coords)
Get edge length.
Definition: Tools.cpp:288
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MoFEM::Core & cOre
Definition: Tools.hpp:35

◆ 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
Examples
find_local_coordinates_on_tet.cpp.

Definition at line 94 of file Tools.cpp.

96  {
99 
101  &elem_coords[0], &elem_coords[3], &elem_coords[6], &elem_coords[9]};
102 
105  shapeFunMBTETAt000[3]};
106  FTensor::Tensor1<double, 3> t_coords_at_0;
107 
108  // Build matrix and get coordinates of zero point
109  // ii - global coordinates
110  // jj - local direvatives
111  MatrixDouble3by3 a(3, 3);
112  for (auto ii : {0, 1, 2}) {
115  &diffShapeFunMBTET[9]);
116  for (auto jj : {0, 1, 2}) {
117  a(jj, ii) = t_diff(i) * t_elem_coords(i);
118  ++t_diff;
119  }
120  t_coords_at_0(ii) = t_n(i) * t_elem_coords(i);
121  ++t_elem_coords;
122  }
123 
125  &global_coords[0], &global_coords[1], &global_coords[2]};
127  &local_coords[0], &local_coords[1], &local_coords[2]};
128 
129  // Calculate right hand side
131  for (int ii = 0; ii != nb_nodes; ++ii) {
132  t_local_coords(j) = t_global_coords(j) - t_coords_at_0(j);
133  ++t_local_coords;
134  ++t_global_coords;
135  }
136 
137  // Solve problem
138  int IPIV[3];
139  int info = lapack_dgesv(3, nb_nodes, &a(0, 0), 3, IPIV, local_coords, 3);
140  if (info != 0)
141  SETERRQ1(PETSC_COMM_SELF, 1, "info == %d", info);
142 
144 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
static constexpr std::array< double, 12 > diffShapeFunMBTET
Definition: Tools.hpp:120
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
static constexpr std::array< double, 4 > shapeFunMBTETAt000
Array of shape function at zero local point on reference element.
Definition: Tools.hpp:179
MatrixBoundedArray< double, 9 > MatrixDouble3by3
Definition: Types.hpp:95
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: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 146 of file Tools.cpp.

148  {
149  MoFEM::Interface &m_field = cOre;
150  moab::Interface &moab(m_field.get_moab());
152  Range to_write;
153  const EntityHandle *conn;
154  int num_nodes;
155  double coords[12];
156  for (auto tet : tets) {
157  CHKERR m_field.get_moab().get_connectivity(tet, conn, num_nodes, true);
158  if (th) {
159  CHKERR moab.tag_get_data(th, conn, num_nodes, coords);
160  } else {
161  CHKERR moab.get_coords(conn, num_nodes, coords);
162  }
163  double q = Tools::volumeLengthQuality(coords);
164  if (f(q)) {
165  out_tets.insert(tet);
166  }
167  }
169 }
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition: Tools.cpp:32
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MoFEM::Core & cOre
Definition: Tools.hpp:35

◆ getTriArea()

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

Get triangle area.

Parameters
tri
Returns
double

Definition at line 280 of file Tools.cpp.

280  {
282  ierr = getTriNormal(tri, &t_normal(0));
283  CHKERRABORT(PETSC_COMM_SELF, ierr);
285  return sqrt(t_normal(i) * t_normal(i)) * 0.5;
286 }
static MoFEMErrorCode getTriNormal(const double *coords, double *normal)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:255

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

255  {
257  double diffN[6];
258  CHKERR ShapeDiffMBTRI(diffN);
259  CHKERR ShapeFaceNormalMBTRI(diffN, coords, normal);
261 }
PetscErrorCode ShapeFaceNormalMBTRI(double *diffN, const double *coords, double *normal)
Definition: fem_tools.c:241
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscErrorCode ShapeDiffMBTRI(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:206
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ getTriNormal() [2/2]

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

Get triangle normal.

Parameters
tri
normal
Returns
MoFEMErrorCode

Definition at line 263 of file Tools.cpp.

264  {
265  MoFEM::Interface &m_field = cOre;
266  moab::Interface &moab(m_field.get_moab());
268  if (moab.type_from_handle(tri) != MBTRI) {
269  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "Works only for triangle");
270  }
271  const EntityHandle *conn;
272  int num_nodes;
273  double coords[9];
274  CHKERR moab.get_connectivity(tri, conn, num_nodes, true);
275  CHKERR moab.get_coords(conn, num_nodes, coords);
276  CHKERR getTriNormal(coords, normal);
278 }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
static MoFEMErrorCode getTriNormal(const double *coords, double *normal)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:255
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MoFEM::Core & cOre
Definition: Tools.hpp:35

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

344  {
345 
347  FTensor::Tensor1<const double *, 3> t_w(&w_ptr[0], &w_ptr[1], &w_ptr[2]);
348  FTensor::Tensor1<const double *, 3> t_v(&v_ptr[0], &v_ptr[1], &v_ptr[2]);
349  FTensor::Tensor1<const double *, 3> t_k(&k_ptr[0], &k_ptr[1], &k_ptr[2]);
350  FTensor::Tensor1<const double *, 3> t_l(&l_ptr[0], &l_ptr[1], &l_ptr[2]);
351 
352  // First segnent is a point
354  t_vw(i) = t_v(i) - t_w(i);
355  double dot_vw = t_vw(i) * t_vw(i);
356  if (dot_vw == 0) {
357  if (tvw_ptr)
358  *tvw_ptr = 0;
359  if (minDistancePointFromOnSegment(k_ptr, l_ptr, w_ptr, tlk_ptr) ==
362  else
363  return SEGMENT_ONE_IS_POINT;
364  }
365 
366  // Second segment is a point
368  t_lk(i) = t_l(i) - t_k(i);
369  double dot_lk = t_lk(i) * t_lk(i);
370  if (dot_lk == 0) {
371  if (tlk_ptr)
372  *tlk_ptr = 0;
373  if (minDistancePointFromOnSegment(w_ptr, v_ptr, k_ptr, tvw_ptr) ==
376  else
377  return SEGMENT_TWO_IS_POINT;
378  }
379 
380  const double a = t_vw(i) * t_vw(i);
381  const double b = -t_vw(i) * t_lk(i);
382  const double c = t_lk(i) * t_lk(i);
383 
384  const double det = a * c - b * b;
385  if (det == 0) {
386 
387  return NO_SOLUTION;
388 
389  } else {
390 
392  t_wk(i) = t_w(i) - t_k(i);
393 
394  const double ft0 = t_vw(i) * t_wk(i);
395  const double ft1 = -t_lk(i) * t_wk(i);
396  const double t0 = (ft1 * b - ft0 * c) / det;
397  const double t1 = (ft0 * b - ft1 * a) / det;
398 
399  if (tvw_ptr)
400  *tvw_ptr = t0;
401  if (tlk_ptr)
402  *tlk_ptr = t1;
403 
404  return SOLUTION_EXIST;
405  }
406 }
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:319

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

320  {
322  FTensor::Tensor1<const double *, 3> t_w(&w_ptr[0], &w_ptr[1], &w_ptr[2]);
323  FTensor::Tensor1<const double *, 3> t_v(&v_ptr[0], &v_ptr[1], &v_ptr[2]);
324  FTensor::Tensor1<const double *, 3> t_p(&p_ptr[0], &p_ptr[1], &p_ptr[2]);
326  t_vw(i) = t_v(i) - t_w(i);
327  const double dot_vw = t_vw(i) * t_vw(i);
328  if (dot_vw == 0) {
329  if (t_ptr)
330  *t_ptr = 0;
331  return SEGMENT_ONE_IS_POINT;
332  }
334  t_pw(i) = t_p(i) - t_w(i);
335  const double t = t_pw(i) * t_vw(i) / dot_vw;
336  if (t_ptr)
337  *t_ptr = t;
338  return SOLUTION_EXIST;
339 }

◆ minTetsQuality()

MoFEMErrorCode MoFEM::Tools::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

Parameters
tetsrange
min_qualitymimimal quality
Returns
error code

Definition at line 66 of file Tools.cpp.

67  {
68  MoFEM::Interface &m_field = cOre;
69  moab::Interface &moab(m_field.get_moab());
71  const EntityHandle *conn;
72  int num_nodes;
73  double coords[12];
74  for (auto tet : tets) {
75  CHKERR m_field.get_moab().get_connectivity(tet, conn, num_nodes, true);
76  if (th) {
77  CHKERR moab.tag_get_data(th, conn, num_nodes, coords);
78  } else {
79  CHKERR moab.get_coords(conn, num_nodes, coords);
80  }
81  double q = Tools::volumeLengthQuality(coords);
82  if (!std::isnormal(q))
83  q = -2;
84  min_quality = f(q, min_quality);
85  }
87 }
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition: Tools.cpp:32
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MoFEM::Core & cOre
Definition: Tools.hpp:35

◆ query_interface()

MoFEMErrorCode MoFEM::Tools::query_interface ( const MOFEMuuid uuid,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 20 of file Tools.cpp.

21  {
23  *iface = NULL;
24  if (uuid == IDD_MOFEMNodeMerger) {
25  *iface = const_cast<Tools *>(this);
27  }
28  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
30 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
static const MOFEMuuid IDD_MOFEMNodeMerger
Definition: NodeMerger.hpp:23
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

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

422  {
424  for (int n = 0; n != nb; ++n) {
425  shape[0] = shapeFunMBTET0(*ksi, *eta, *zeta);
426  shape[1] = shapeFunMBTET1(*ksi, *eta, *zeta);
427  shape[2] = shapeFunMBTET2(*ksi, *eta, *zeta);
428  shape[3] = shapeFunMBTET3(*ksi, *eta, *zeta);
429  shape += 4;
430  ksi += LDB;
431  eta += LDB;
432  zeta += LDB;
433  }
435 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
static double shapeFunMBTET3(const double x, const double y, const double z)
Definition: Tools.hpp:415
static double shapeFunMBTET1(const double x, const double y, const double z)
Definition: Tools.hpp:407
static double shapeFunMBTET2(const double x, const double y, const double z)
Definition: Tools.hpp:411
static double shapeFunMBTET0(const double x, const double y, const double z)
Definition: Tools.hpp:403

◆ shapeFunMBTET0()

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

Definition at line 403 of file Tools.hpp.

403  {
404  return N_MBTET0(x, y, z);
405 }
#define N_MBTET0(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:36

◆ shapeFunMBTET1()

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

Definition at line 407 of file Tools.hpp.

407  {
408  return N_MBTET1(x, y, z);
409 }
#define N_MBTET1(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:37

◆ shapeFunMBTET2()

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

Definition at line 411 of file Tools.hpp.

411  {
412  return N_MBTET2(x, y, z);
413 }
#define N_MBTET2(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:38

◆ shapeFunMBTET3()

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

Definition at line 415 of file Tools.hpp.

415  {
416  return N_MBTET3(x, y, z);
417 };
#define N_MBTET3(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:39

◆ tetVolume()

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

Calculate volume of tetrahedron.

Parameters
coords
Returns
double volume

Definition at line 47 of file Tools.cpp.

47  {
48  double diff_n[12];
49  ShapeDiffMBTET(diff_n);
50  FTensor::Tensor1<double *, 3> t_diff_n(&diff_n[0], &diff_n[1], &diff_n[2], 3);
51  FTensor::Tensor1<const double *, 3> t_coords(&coords[0], &coords[1],
52  &coords[2], 3);
56  jac(i, j) = 0;
57  for (int nn = 0; nn != 4; nn++) {
58  jac(i, j) += t_coords(i) * t_diff_n(j);
59  ++t_coords;
60  ++t_diff_n;
61  }
62  return dEterminant(jac) / 6.;
63 }
PetscErrorCode ShapeDiffMBTET(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:331
static double dEterminant(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:385

◆ volumeLengthQuality()

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

Calculate tetrahedron volume length quality.

Parameters
coordstet coordinates
Returns
Volume-length quality
Examples
mesh_smoothing.cpp.

Definition at line 32 of file Tools.cpp.

32  {
33  double lrms = 0;
34  for (int dd = 0; dd != 3; dd++) {
35  lrms += pow(coords[0 * 3 + dd] - coords[1 * 3 + dd], 2) +
36  pow(coords[0 * 3 + dd] - coords[2 * 3 + dd], 2) +
37  pow(coords[0 * 3 + dd] - coords[3 * 3 + dd], 2) +
38  pow(coords[1 * 3 + dd] - coords[2 * 3 + dd], 2) +
39  pow(coords[1 * 3 + dd] - coords[3 * 3 + dd], 2) +
40  pow(coords[2 * 3 + dd] - coords[3 * 3 + dd], 2);
41  }
42  lrms = sqrt((1. / 6.) * lrms);
43  double volume = tetVolume(coords);
44  return 6. * sqrt(2.) * volume / pow(lrms, 3);
45 }
static double tetVolume(const double *coords)
Calculate volume of tetrahedron.
Definition: Tools.cpp:47
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

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

175  {
176  MoFEM::Interface &m_field = cOre;
177  moab::Interface &moab(m_field.get_moab());
179  Range out_tets;
180  CHKERR getTetsWithQuality(out_tets, tets, th, f);
181  EntityHandle meshset;
182  CHKERR moab.create_meshset(MESHSET_SET, meshset);
183  CHKERR moab.add_entities(meshset, out_tets);
184  CHKERR moab.write_file(file_name, file_type, options, &meshset, 1);
185  CHKERR moab.delete_entities(&meshset, 1);
187 }
virtual moab::Interface & get_moab()=0
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:146
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MoFEM::Core & cOre
Definition: Tools.hpp:35

Member Data Documentation

◆ cOre

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

Definition at line 35 of file Tools.hpp.

◆ diffN_MBEDGE0x

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

Definition at line 68 of file Tools.hpp.

◆ diffN_MBEDGE1x

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

Definition at line 69 of file Tools.hpp.

◆ diffShapeFunMBEDGE

constexpr std::array< double, 2 > MoFEM::Tools::diffShapeFunMBEDGE
static
Initial value:
Examples
bernstein_bezier_genrate_base.cpp.

Definition at line 71 of file Tools.hpp.

◆ diffShapeFunMBTET

constexpr std::array< double, 12 > MoFEM::Tools::diffShapeFunMBTET
static

◆ diffShapeFunMBTET0x

constexpr double MoFEM::Tools::diffShapeFunMBTET0x
static
Initial value:

derivative of tetrahedral shape function

Definition at line 95 of file Tools.hpp.

◆ diffShapeFunMBTET0y

constexpr double MoFEM::Tools::diffShapeFunMBTET0y
static
Initial value:

derivative of tetrahedral shape function

Definition at line 97 of file Tools.hpp.

◆ diffShapeFunMBTET0z

constexpr double MoFEM::Tools::diffShapeFunMBTET0z
static
Initial value:

derivative of tetrahedral shape function

Definition at line 99 of file Tools.hpp.

◆ diffShapeFunMBTET1x

constexpr double MoFEM::Tools::diffShapeFunMBTET1x
static
Initial value:

derivative of tetrahedral shape function

Definition at line 101 of file Tools.hpp.

◆ diffShapeFunMBTET1y

constexpr double MoFEM::Tools::diffShapeFunMBTET1y
static
Initial value:

derivative of tetrahedral shape function

Definition at line 103 of file Tools.hpp.

◆ diffShapeFunMBTET1z

constexpr double MoFEM::Tools::diffShapeFunMBTET1z
static
Initial value:

derivative of tetrahedral shape function

Definition at line 105 of file Tools.hpp.

◆ diffShapeFunMBTET2x

constexpr double MoFEM::Tools::diffShapeFunMBTET2x
static
Initial value:

derivative of tetrahedral shape function

Definition at line 107 of file Tools.hpp.

◆ diffShapeFunMBTET2y

constexpr double MoFEM::Tools::diffShapeFunMBTET2y
static
Initial value:

derivative of tetrahedral shape function

Definition at line 109 of file Tools.hpp.

◆ diffShapeFunMBTET2z

constexpr double MoFEM::Tools::diffShapeFunMBTET2z
static
Initial value:

derivative of tetrahedral shape function

Definition at line 111 of file Tools.hpp.

◆ diffShapeFunMBTET3x

constexpr double MoFEM::Tools::diffShapeFunMBTET3x
static
Initial value:

derivative of tetrahedral shape function

Definition at line 113 of file Tools.hpp.

◆ diffShapeFunMBTET3y

constexpr double MoFEM::Tools::diffShapeFunMBTET3y
static
Initial value:

derivative of tetrahedral shape function

Definition at line 115 of file Tools.hpp.

◆ diffShapeFunMBTET3z

constexpr double MoFEM::Tools::diffShapeFunMBTET3z
static
Initial value:

derivative of tetrahedral shape function

Definition at line 117 of file Tools.hpp.

◆ diffShapeFunMBTRI

constexpr std::array< double, 6 > MoFEM::Tools::diffShapeFunMBTRI
static

◆ diffShapeFunMBTRI0x

constexpr double MoFEM::Tools::diffShapeFunMBTRI0x
static
Initial value:

derivative of triangle shape function

Definition at line 74 of file Tools.hpp.

◆ diffShapeFunMBTRI0y

constexpr double MoFEM::Tools::diffShapeFunMBTRI0y
static
Initial value:

derivative of triangle shape function

Definition at line 76 of file Tools.hpp.

◆ diffShapeFunMBTRI1x

constexpr double MoFEM::Tools::diffShapeFunMBTRI1x
static
Initial value:

derivative of triangle shape function

Definition at line 78 of file Tools.hpp.

◆ diffShapeFunMBTRI1y

constexpr double MoFEM::Tools::diffShapeFunMBTRI1y
static
Initial value:

derivative of triangle shape function

Definition at line 80 of file Tools.hpp.

◆ diffShapeFunMBTRI2x

constexpr double MoFEM::Tools::diffShapeFunMBTRI2x
static
Initial value:

derivative of triangle shape function

Definition at line 82 of file Tools.hpp.

◆ diffShapeFunMBTRI2y

constexpr double MoFEM::Tools::diffShapeFunMBTRI2y
static
Initial value:

derivative of triangle shape function

Definition at line 84 of file Tools.hpp.

◆ shapeFunMBTET0At000

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

Definition at line 142 of file Tools.hpp.

◆ shapeFunMBTET0AtOneThird

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

Definition at line 147 of file Tools.hpp.

◆ shapeFunMBTET1At000

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

Definition at line 143 of file Tools.hpp.

◆ shapeFunMBTET1AtOneThird

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

Definition at line 149 of file Tools.hpp.

◆ shapeFunMBTET2At000

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

Definition at line 144 of file Tools.hpp.

◆ shapeFunMBTET2AtOneThird

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

Definition at line 151 of file Tools.hpp.

◆ shapeFunMBTET3At000

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

Definition at line 145 of file Tools.hpp.

◆ shapeFunMBTET3AtOneThird

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

Definition at line 153 of file Tools.hpp.

◆ shapeFunMBTETAt000

constexpr std::array< double, 4 > MoFEM::Tools::shapeFunMBTETAt000
static
Initial value:

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

Definition at line 179 of file Tools.hpp.

◆ shapeFunMBTETAtOneThird

constexpr std::array<double, 4> MoFEM::Tools::shapeFunMBTETAtOneThird
static
Initial value:

Array of shape function at center on reference element.

Definition at line 187 of file Tools.hpp.


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