v0.9.0
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 (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 ()=default
 
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
 

Static Public Member Functions

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

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 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 std::array< double, 8 > diffShapeFunMBQUADAtCenter
 
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 332 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 192 of file Tools.cpp.

194  {
195  double loc_coord[] = {0, 0, 0};
196  double N[4], diffN[12];
198  CHKERR ShapeDiffMBTET(diffN);
199  CHKERR ShapeMBTET(N, &loc_coord[0], &loc_coord[1], &loc_coord[2], 1);
200  CHKERR ShapeMBTET_inverse(N, diffN, tet_coords, global_coord, loc_coord);
201  CHKERR ShapeMBTET(N, &loc_coord[0], &loc_coord[1], &loc_coord[2], 1);
202  result = true;
203  for (int n = 0; n != 4; ++n) {
204  if (N[n] < -tol || (N[n] - 1) > tol) {
205  result = false;
206  break;
207  }
208  }
210 }
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 212 of file Tools.cpp.

214  {
216  int loc_size;
217  CHKERR VecGetLocalSize(v, &loc_size);
218  int prb_loc_size = 0;
219  boost::shared_ptr<NumeredDofEntity_multiIndex> prb_dofs;
220  switch (row_or_col) {
221  case ROW:
222  prb_loc_size = prb_ptr->getNbLocalDofsRow();
223  prb_dofs = prb_ptr->getNumeredDofsRows();
224  break;
225  case COL:
226  prb_loc_size = prb_ptr->getNbLocalDofsCol();
227  prb_dofs = prb_ptr->getNumeredDofsCols();
228  break;
229  break;
230  default:
231  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
232  "Wrong argument, row_or_col should be row or column");
233  }
234  if (loc_size != prb_loc_size) {
235  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
236  "Inconsistent size of vector and problem %d != %d", loc_size,
237  prb_loc_size);
238  }
239  const double *a;
240  CHKERR VecGetArrayRead(v, &a);
241  MPI_Comm comm = PetscObjectComm((PetscObject)v);
242  for (int ii = 0; ii != loc_size; ++ii) {
243  if (!boost::math::isfinite(a[ii])) {
244  NumeredDofEntityByLocalIdx::iterator dit =
245  prb_dofs->get<PetscLocalIdx_mi_tag>().find(ii);
246  std::ostringstream ss;
247  ss << "Not a number " << a[ii] << " on dof: " << endl
248  << **dit << endl
249  << endl;
250  PetscSynchronizedPrintf(comm, "%s", ss.str().c_str());
251  }
252  }
253  CHKERR VecRestoreArrayRead(v, &a);
254  PetscSynchronizedFlush(comm, PETSC_STDOUT);
256 }
#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 411 of file Tools.cpp.

413  {
414  MoFEM::Interface &m_field = cOre;
415  moab::Interface &moab(m_field.get_moab());
417 
419 
420  auto get_point = [i](auto &t_w, auto &t_delta, auto t) {
422  t = std::max(0., std::min(1., t));
423  t_p(i) = t_w(i) + t * t_delta(i);
424  return t_p;
425  };
426 
427  auto get_distance = [i](auto &t_p, auto &t_n) {
428  FTensor::Tensor1<double, 3> t_dist_vector;
429  t_dist_vector(i) = t_p(i) - t_n(i);
430  return sqrt(t_dist_vector(i) * t_dist_vector(i));
431  };
432 
433  for (auto e : edges) {
434 
435  int num_nodes;
436  const EntityHandle *conn_fixed;
437  CHKERR moab.get_connectivity(e, conn_fixed, num_nodes, true);
438  VectorDouble6 coords_fixed(6);
439  CHKERR moab.get_coords(conn_fixed, num_nodes, &coords_fixed[0]);
441  &coords_fixed[0], &coords_fixed[1], &coords_fixed[2]);
443  &coords_fixed[3], &coords_fixed[4], &coords_fixed[5]);
444 
445  FTensor::Tensor1<double, 3> t_edge_delta;
446  t_edge_delta(i) = t_f1(i) - t_f0(i);
447 
449  v_ptr, v_ptr + 1, v_ptr + 2);
451  o_ptr, o_ptr + 1, o_ptr + 2);
452  FTensor::Tensor0<FTensor::PackPtr<double *, 1>> t_min_dist(min_dist_ptr);
453 
454  EntityHandle *colsest_segment_it = nullptr;
455  if (o_segments)
456  colsest_segment_it = o_segments;
457 
458  for (int n = 0; n != nb; ++n) {
459 
460  double t;
461  if (Tools::minDistancePointFromOnSegment(&t_f0(0), &t_f1(0), &t_n(0),
462  &t) == Tools::SOLUTION_EXIST) {
463  auto t_p = get_point(t_f0, t_edge_delta, t);
464  auto dist_n = get_distance(t_p, t_n);
465  if (dist_n < t_min_dist || t_min_dist < 0) {
466  t_min_dist = dist_n;
467  if (o_ptr)
468  t_min_coords(i) = t_p(i);
469  if (o_segments)
470  *colsest_segment_it = e;
471  }
472  }
473 
474  ++t_n;
475  ++t_min_dist;
476  if (o_ptr)
477  ++t_min_coords;
478  if (o_segments)
479  ++colsest_segment_it;
480  }
481  }
482 
484 }
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:90
#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:322
#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 291 of file Tools.cpp.

291  {
293  edge_coords[2]);
295  edge_coords[5]);
297  t_coords_n0(i) -= t_coords_n1(i);
298  return sqrt(t_coords_n0(i) * t_coords_n0(i));
299 }
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 301 of file Tools.cpp.

301  {
302  MoFEM::Interface &m_field = cOre;
303  moab::Interface &moab(m_field.get_moab());
304  auto get_edge_coords = [edge, &moab](double *const coords) {
306  if (moab.type_from_handle(edge) != MBEDGE) {
307  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "Works only for edge");
308  }
309  const EntityHandle *conn;
310  int num_nodes;
311  CHKERR moab.get_connectivity(edge, conn, num_nodes, true);
312  CHKERR moab.get_coords(conn, 2, coords);
314  };
315  double coords[6];
316  ierr = get_edge_coords(coords);
317  CHKERRABORT(PETSC_COMM_SELF, ierr);
318  return getEdgeLength(coords);
319 }
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:291
#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 97 of file Tools.cpp.

99  {
102 
104  &elem_coords[0], &elem_coords[3], &elem_coords[6], &elem_coords[9]};
105 
108  shapeFunMBTETAt000[3]};
109  FTensor::Tensor1<double, 3> t_coords_at_0;
110 
111  // Build matrix and get coordinates of zero point
112  // ii - global coordinates
113  // jj - local direvatives
114  MatrixDouble3by3 a(3, 3);
115  for (auto ii : {0, 1, 2}) {
118  &diffShapeFunMBTET[9]);
119  for (auto jj : {0, 1, 2}) {
120  a(jj, ii) = t_diff(i) * t_elem_coords(i);
121  ++t_diff;
122  }
123  t_coords_at_0(ii) = t_n(i) * t_elem_coords(i);
124  ++t_elem_coords;
125  }
126 
128  &global_coords[0], &global_coords[1], &global_coords[2]};
130  &local_coords[0], &local_coords[1], &local_coords[2]};
131 
132  // Calculate right hand side
134  for (int ii = 0; ii != nb_nodes; ++ii) {
135  t_local_coords(j) = t_global_coords(j) - t_coords_at_0(j);
136  ++t_local_coords;
137  ++t_global_coords;
138  }
139 
140  // Solve problem
141  int IPIV[3];
142  int info = lapack_dgesv(3, nb_nodes, &a(0, 0), 3, IPIV, local_coords, 3);
143  if (info != 0)
144  SETERRQ1(PETSC_COMM_SELF, 1, "info == %d", info);
145 
147 }
#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:143
#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:202
MatrixBoundedArray< double, 9 > MatrixDouble3by3
Definition: Types.hpp:97
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 149 of file Tools.cpp.

151  {
152  MoFEM::Interface &m_field = cOre;
153  moab::Interface &moab(m_field.get_moab());
155  Range to_write;
156  const EntityHandle *conn;
157  int num_nodes;
158  double coords[12];
159  for (auto tet : tets) {
160  CHKERR m_field.get_moab().get_connectivity(tet, conn, num_nodes, true);
161  if (th) {
162  CHKERR moab.tag_get_data(th, conn, num_nodes, coords);
163  } else {
164  CHKERR moab.get_coords(conn, num_nodes, coords);
165  }
166  double q = Tools::volumeLengthQuality(coords);
167  if (f(q)) {
168  out_tets.insert(tet);
169  }
170  }
172 }
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition: Tools.cpp:34
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 283 of file Tools.cpp.

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

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

258  {
260  double diffN[6];
261  CHKERR ShapeDiffMBTRI(diffN);
262  CHKERR ShapeFaceNormalMBTRI(diffN, coords, normal);
264 }
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 266 of file Tools.cpp.

267  {
268  MoFEM::Interface &m_field = cOre;
269  moab::Interface &moab(m_field.get_moab());
271  if (moab.type_from_handle(tri) != MBTRI) {
272  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA, "Works only for triangle");
273  }
274  const EntityHandle *conn;
275  int num_nodes;
276  double coords[9];
277  CHKERR moab.get_connectivity(tri, conn, num_nodes, true);
278  CHKERR moab.get_coords(conn, num_nodes, coords);
279  CHKERR getTriNormal(coords, normal);
281 }
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:258
#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 345 of file Tools.cpp.

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

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

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

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

69  {
70  MoFEM::Interface &m_field = cOre;
71  moab::Interface &moab(m_field.get_moab());
73  const EntityHandle *conn;
74  int num_nodes;
75  double coords[12];
76  for (auto tet : tets) {
77  CHKERR m_field.get_moab().get_connectivity(tet, conn, num_nodes, true);
78  if (th) {
79  CHKERR moab.tag_get_data(th, conn, num_nodes, coords);
80  } else {
81  CHKERR moab.get_coords(conn, num_nodes, coords);
82  }
83  double q = Tools::volumeLengthQuality(coords);
84  if (!std::isnormal(q))
85  q = -2;
86  min_quality = f(q, min_quality);
87  }
89 }
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition: Tools.cpp:34
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

◆ outerProductOfEdgeIntegrationPtsForQuad()

MoFEMErrorCode MoFEM::Tools::outerProductOfEdgeIntegrationPtsForQuad ( MatrixDouble pts,
const int  edge0,
const int  edge1 
)
static
Examples
gauss_points_on_quad.cpp.

Definition at line 486 of file Tools.cpp.

487  {
489 
490  auto check_rule_edge = [](int rule) {
492  if (rule < 0) {
493  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
494  "Wrong integration rule: %d", rule);
495  }
496  if (rule > QUAD_1D_TABLE_SIZE) {
497  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
498  "rule > quadrature order %d < %d", rule, QUAD_1D_TABLE_SIZE);
499  }
500  if (QUAD_1D_TABLE[rule]->dim != 1) {
501  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
502  }
503  if (QUAD_1D_TABLE[rule]->order < rule) {
504  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
505  "wrong order %d != %d", QUAD_1D_TABLE[rule]->order, rule);
506  }
508  };
509 
510  CHKERR check_rule_edge(rule_ksi);
511  CHKERR check_rule_edge(rule_eta);
512 
513  int nb_gauss_pts_ksi = QUAD_1D_TABLE[rule_ksi]->npoints;
514  int nb_gauss_pts_eta = QUAD_1D_TABLE[rule_eta]->npoints;
515  gauss_pts.resize(3, nb_gauss_pts_ksi * nb_gauss_pts_eta, false);
516  gauss_pts.clear();
517 
518  VectorDouble ones;
519  ones.resize(max(nb_gauss_pts_ksi, nb_gauss_pts_eta), false);
520  fill(ones.begin(), ones.end(), 1.0);
521 
522  cblas_dger(CblasRowMajor, nb_gauss_pts_eta, nb_gauss_pts_ksi, 1, &ones(0), 1,
523  &QUAD_1D_TABLE[rule_ksi]->points[1], 2, &gauss_pts(0, 0),
524  nb_gauss_pts_ksi);
525 
526  cblas_dger(CblasRowMajor, nb_gauss_pts_eta, nb_gauss_pts_ksi, 1,
527  &QUAD_1D_TABLE[rule_eta]->points[1], 2, &ones(0), 1,
528  &gauss_pts(1, 0), nb_gauss_pts_ksi);
529 
530  cblas_dger(CblasRowMajor, nb_gauss_pts_eta, nb_gauss_pts_ksi, 1,
531  QUAD_1D_TABLE[rule_eta]->weights, 1,
532  QUAD_1D_TABLE[rule_ksi]->weights, 1, &gauss_pts(2, 0),
533  nb_gauss_pts_ksi);
534 
536 }
static QUAD *const QUAD_1D_TABLE[]
Definition: quad.h:164
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
int npoints
Definition: quad.h:29
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
void cblas_dger(const enum CBLAS_ORDER order, const int M, const int N, const double alpha, const double *X, const int incX, const double *Y, const int incY, double *A, const int lda)
Definition: cblas_dger.c:12
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
#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
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72
#define QUAD_1D_TABLE_SIZE
Definition: quad.h:163

◆ query_interface()

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

Implements MoFEM::UnknownInterface.

Definition at line 22 of file Tools.cpp.

23  {
25  *iface = NULL;
26  if (uuid == IDD_MOFEMNodeMerger) {
27  *iface = const_cast<Tools *>(this);
29  }
30  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
32 }
#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 448 of file Tools.hpp.

450  {
452  for (int n = 0; n != nb; ++n) {
453  shape[0] = shapeFunMBTET0(*ksi, *eta, *zeta);
454  shape[1] = shapeFunMBTET1(*ksi, *eta, *zeta);
455  shape[2] = shapeFunMBTET2(*ksi, *eta, *zeta);
456  shape[3] = shapeFunMBTET3(*ksi, *eta, *zeta);
457  shape += 4;
458  ksi += LDB;
459  eta += LDB;
460  zeta += LDB;
461  }
463 }
#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:443
static double shapeFunMBTET1(const double x, const double y, const double z)
Definition: Tools.hpp:435
static double shapeFunMBTET2(const double x, const double y, const double z)
Definition: Tools.hpp:439
static double shapeFunMBTET0(const double x, const double y, const double z)
Definition: Tools.hpp:431

◆ shapeFunMBTET0()

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

Definition at line 431 of file Tools.hpp.

431  {
432  return N_MBTET0(x, y, z);
433 }
#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 435 of file Tools.hpp.

435  {
436  return N_MBTET1(x, y, z);
437 }
#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 439 of file Tools.hpp.

439  {
440  return N_MBTET2(x, y, z);
441 }
#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 443 of file Tools.hpp.

443  {
444  return N_MBTET3(x, y, z);
445 };
#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 49 of file Tools.cpp.

49  {
50  double diff_n[12];
51  ShapeDiffMBTET(diff_n);
52  FTensor::Tensor1<double *, 3> t_diff_n(&diff_n[0], &diff_n[1], &diff_n[2], 3);
53  FTensor::Tensor1<const double *, 3> t_coords(&coords[0], &coords[1],
54  &coords[2], 3);
58  jac(i, j) = 0;
59  for (int nn = 0; nn != 4; nn++) {
60  jac(i, j) += t_coords(i) * t_diff_n(j);
61  ++t_coords;
62  ++t_diff_n;
63  }
64  return dEterminant(jac) / 6.;
65 }
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:431

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

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

178  {
179  MoFEM::Interface &m_field = cOre;
180  moab::Interface &moab(m_field.get_moab());
182  Range out_tets;
183  CHKERR getTetsWithQuality(out_tets, tets, th, f);
184  EntityHandle meshset;
185  CHKERR moab.create_meshset(MESHSET_SET, meshset);
186  CHKERR moab.add_entities(meshset, out_tets);
187  CHKERR moab.write_file(file_name, file_type, options, &meshset, 1);
188  CHKERR moab.delete_entities(&meshset, 1);
190 }
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:149
#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_generate_base.cpp.

Definition at line 71 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter

constexpr std::array< double, 8 > MoFEM::Tools::diffShapeFunMBQUADAtCenter
static

◆ diffShapeFunMBQUADAtCenter0x

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

derivative of quad shape function

Definition at line 95 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter0y

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

derivative of quad shape function

Definition at line 97 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter1x

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

derivative of quad shape function

Definition at line 99 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter1y

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

derivative of quad shape function

Definition at line 101 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter2x

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

derivative of quad shape function

Definition at line 103 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter2y

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

derivative of quad shape function

Definition at line 105 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter3x

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

derivative of quad shape function

Definition at line 107 of file Tools.hpp.

◆ diffShapeFunMBQUADAtCenter3y

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

derivative of quad shape function

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

◆ diffShapeFunMBTET0y

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

derivative of tetrahedral shape function

Definition at line 120 of file Tools.hpp.

◆ diffShapeFunMBTET0z

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

derivative of tetrahedral shape function

Definition at line 122 of file Tools.hpp.

◆ diffShapeFunMBTET1x

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

derivative of tetrahedral shape function

Definition at line 124 of file Tools.hpp.

◆ diffShapeFunMBTET1y

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

derivative of tetrahedral shape function

Definition at line 126 of file Tools.hpp.

◆ diffShapeFunMBTET1z

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

derivative of tetrahedral shape function

Definition at line 128 of file Tools.hpp.

◆ diffShapeFunMBTET2x

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

derivative of tetrahedral shape function

Definition at line 130 of file Tools.hpp.

◆ diffShapeFunMBTET2y

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

derivative of tetrahedral shape function

Definition at line 132 of file Tools.hpp.

◆ diffShapeFunMBTET2z

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

derivative of tetrahedral shape function

Definition at line 134 of file Tools.hpp.

◆ diffShapeFunMBTET3x

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

derivative of tetrahedral shape function

Definition at line 136 of file Tools.hpp.

◆ diffShapeFunMBTET3y

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

derivative of tetrahedral shape function

Definition at line 138 of file Tools.hpp.

◆ diffShapeFunMBTET3z

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

derivative of tetrahedral shape function

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

◆ shapeFunMBTET0AtOneThird

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

Definition at line 170 of file Tools.hpp.

◆ shapeFunMBTET1At000

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

Definition at line 166 of file Tools.hpp.

◆ shapeFunMBTET1AtOneThird

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

Definition at line 172 of file Tools.hpp.

◆ shapeFunMBTET2At000

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

Definition at line 167 of file Tools.hpp.

◆ shapeFunMBTET2AtOneThird

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

Definition at line 174 of file Tools.hpp.

◆ shapeFunMBTET3At000

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

Definition at line 168 of file Tools.hpp.

◆ shapeFunMBTET3AtOneThird

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

Definition at line 176 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 202 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 210 of file Tools.hpp.


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