v0.10.0
Tools.hpp
Go to the documentation of this file.
1 /** \file Tools.hpp
2  * \brief Tools interface
3  *
4  * Implementation of some useful and very often used methods * in MoFEM.
5  *
6  * \ingroup mofem_tools
7  */
8 
9 /* MoFEM is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
12  * License for more details.
13  *
14  * You should have received a copy of the GNU Lesser General Public
15  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
16  */
17 
18 #ifndef __TOOLS_HPP__
19 #define __TOOLS_HPP__
20 
21 namespace MoFEM {
22 
24 
25 /**
26  * \brief Auxiliary tools
27  * \nosubgrouping
28  * \ingroup mofem_tools
29  */
30 struct Tools : public UnknownInterface {
31 
33  UnknownInterface **iface) const;
34 
36  Tools(const MoFEM::Core &core) : cOre(const_cast<MoFEM::Core &>(core)) {}
37 
38  /** \name Computational */
39 
40  /**
41  * \brief Calculate tetrahedron volume length quality
42  * @param coords tet coordinates
43  * @return Volume-length quality
44  */
45  static double volumeLengthQuality(const double *coords);
46 
47  /**
48  * @brief Calculate volume of tetrahedron
49  *
50  * @param coords
51  * @return double volume
52  */
53  static double tetVolume(const double *coords);
54 
55  /**
56  * \brief calculate minimal quality of tetrahedra in range
57  * @param tets range
58  * @param min_quality mimimal quality
59  * @return error code
60  */
61  MoFEMErrorCode minTetsQuality(const Range &tets, double &min_quality,
62  Tag th = nullptr,
63  boost::function<double(double, double)> f =
64  [](double a, double b) -> double {
65  return std::min(a, b);
66  });
67 
68  static constexpr double diffN_MBEDGE0x = diffN_MBEDGE0;
69  static constexpr double diffN_MBEDGE1x = diffN_MBEDGE1;
70 
71  static constexpr std::array<double, 2> diffShapeFunMBEDGE = {diffN_MBEDGE0x,
73 
74  static constexpr double diffShapeFunMBTRI0x =
75  diffN_MBTRI0x; ///< derivative of triangle shape function
76  static constexpr double diffShapeFunMBTRI0y =
77  diffN_MBTRI0y; ///< derivative of triangle shape function
78  static constexpr double diffShapeFunMBTRI1x =
79  diffN_MBTRI1x; ///< derivative of triangle shape function
80  static constexpr double diffShapeFunMBTRI1y =
81  diffN_MBTRI1y; ///< derivative of triangle shape function
82  static constexpr double diffShapeFunMBTRI2x =
83  diffN_MBTRI2x; ///< derivative of triangle shape function
84  static constexpr double diffShapeFunMBTRI2y =
85  diffN_MBTRI2y; ///< derivative of triangle shape function
86 
87  static constexpr std::array<double, 6> diffShapeFunMBTRI = {
88 
90 
92 
94 
95  static inline double shapeFunMBTRI0(const double x, const double y);
96 
97  static inline double shapeFunMBTRI1(const double x, const double y);
98 
99  static inline double shapeFunMBTRI2(const double x, const double y);
100 
101  /**
102  * @brief Calculate shape functions on triangle
103  *
104  * \note Template parameter is leading dimension of point coordinate arrays,
105  * such that \f$ksi_{n+1} = ksi[n + LDB]\f$
106  *
107  * @tparam 1
108  * @param shape shape functions
109  * @param ksi pointer to first local coordinates
110  * @param eta pointer to second local coordinates
111  * @param nb number of points
112  * @return MoFEMErrorCode
113  */
114  template <int LDB = 1>
115  static MoFEMErrorCode shapeFunMBTRI(double *shape, const double *ksi,
116  const double *eta, const double nb);
117 
118  static constexpr double diffShapeFunMBQUADAtCenter0x =
119  diffN_MBQUAD0x(0.5); ///< derivative of quad shape function
120  static constexpr double diffShapeFunMBQUADAtCenter0y =
121  diffN_MBQUAD0y(0.5); ///< derivative of quad shape function
122  static constexpr double diffShapeFunMBQUADAtCenter1x =
123  diffN_MBQUAD1x(0.5); ///< derivative of quad shape function
124  static constexpr double diffShapeFunMBQUADAtCenter1y =
125  diffN_MBQUAD1y(0.5); ///< derivative of quad shape function
126  static constexpr double diffShapeFunMBQUADAtCenter2x =
127  diffN_MBQUAD2x(0.5); ///< derivative of quad shape function
128  static constexpr double diffShapeFunMBQUADAtCenter2y =
129  diffN_MBQUAD2y(0.5); ///< derivative of quad shape function
130  static constexpr double diffShapeFunMBQUADAtCenter3x =
131  diffN_MBQUAD3x(0.5); ///< derivative of quad shape function
132  static constexpr double diffShapeFunMBQUADAtCenter3y =
133  diffN_MBQUAD3y(0.5); ///< derivative of quad shape function
134 
135  static constexpr std::array<double, 8> diffShapeFunMBQUADAtCenter = {
140 
141  static constexpr double diffShapeFunMBTET0x =
142  diffN_MBTET0x; ///< derivative of tetrahedral shape function
143  static constexpr double diffShapeFunMBTET0y =
144  diffN_MBTET0y; ///< derivative of tetrahedral shape function
145  static constexpr double diffShapeFunMBTET0z =
146  diffN_MBTET0z; ///< derivative of tetrahedral shape function
147  static constexpr double diffShapeFunMBTET1x =
148  diffN_MBTET1x; ///< derivative of tetrahedral shape function
149  static constexpr double diffShapeFunMBTET1y =
150  diffN_MBTET1y; ///< derivative of tetrahedral shape function
151  static constexpr double diffShapeFunMBTET1z =
152  diffN_MBTET1z; ///< derivative of tetrahedral shape function
153  static constexpr double diffShapeFunMBTET2x =
154  diffN_MBTET2x; ///< derivative of tetrahedral shape function
155  static constexpr double diffShapeFunMBTET2y =
156  diffN_MBTET2y; ///< derivative of tetrahedral shape function
157  static constexpr double diffShapeFunMBTET2z =
158  diffN_MBTET2z; ///< derivative of tetrahedral shape function
159  static constexpr double diffShapeFunMBTET3x =
160  diffN_MBTET3x; ///< derivative of tetrahedral shape function
161  static constexpr double diffShapeFunMBTET3y =
162  diffN_MBTET3y; ///< derivative of tetrahedral shape function
163  static constexpr double diffShapeFunMBTET3z =
164  diffN_MBTET3z; ///< derivative of tetrahedral shape function
165 
166  static constexpr std::array<double, 12> diffShapeFunMBTET = {
167 
169 
171 
173 
175 
176  static inline double shapeFunMBTET0(const double x, const double y,
177  const double z);
178 
179  static inline double shapeFunMBTET1(const double x, const double y,
180  const double z);
181 
182  static inline double shapeFunMBTET2(const double x, const double y,
183  const double z);
184 
185  static inline double shapeFunMBTET3(const double x, const double y,
186  const double z);
187 
188  static constexpr double shapeFunMBTET0At000 = N_MBTET0(0, 0, 0);
189  static constexpr double shapeFunMBTET1At000 = N_MBTET1(0, 0, 0);
190  static constexpr double shapeFunMBTET2At000 = N_MBTET2(0, 0, 0);
191  static constexpr double shapeFunMBTET3At000 = N_MBTET3(0, 0, 0);
192 
193  static constexpr double shapeFunMBTET0AtOneThird =
194  N_MBTET0(1. / 3., 1. / 3., 1. / 3.);
195  static constexpr double shapeFunMBTET1AtOneThird =
196  N_MBTET1(1. / 3., 1. / 3., 1. / 3.);
197  static constexpr double shapeFunMBTET2AtOneThird =
198  N_MBTET2(1. / 3., 1. / 3., 1. / 3.);
199  static constexpr double shapeFunMBTET3AtOneThird =
200  N_MBTET3(1. / 3., 1. / 3., 1. / 3.);
201 
202  /**
203  * @brief Calculate shape functions on tetrahedron
204  *
205  * \note Template parameter is leading dimension of point coordinate arrays,
206  * such that \f$ksi_{n+1} = ksi[n + LDB]\f$
207  *
208  * @tparam 1
209  * @param shape shape functions
210  * @param ksi pointer to first local coordinates
211  * @param eta pointer to second local coordinates
212  * @param zeta pointer to first third coordinates
213  * @param nb number of points
214  * @return MoFEMErrorCode
215  */
216  template <int LDB = 1>
217  static MoFEMErrorCode shapeFunMBTET(double *shape, const double *ksi,
218  const double *eta, const double *zeta,
219  const double nb);
220 
221  /**
222  * @brief Array of shape function at zero local point on reference element
223  *
224  */
225  static constexpr std::array<double, 4> shapeFunMBTETAt000 = {
228 
229  /**
230  * @brief Array of shape function at center on reference element
231  *
232  */
233  static constexpr std::array<double, 4> shapeFunMBTETAtOneThird = {
236 
237 
238  /**
239  * @brief Get the Local Coordinates On Reference Four Node Tet object
240  *
241  * \code
242  * MatrixDouble elem_coords(4, 3);
243  * // Set nodal coordinates
244  * MatrixDouble global_coords(5, 3);
245  * // Set global coordinates
246  * MatrixDouble local_coords(global_coords.size1(), 3);
247  * CHKERR Tools::getLocalCoordinatesOnReferenceFourNodeTet(
248  * &elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
249  * &local_coords(0, 0))
250  * \endcode
251  *
252  * @param elem_coords Global element node coordinates
253  * @param glob_coords Globale coordinates
254  * @param nb_nodes Number of points
255  * @param local_coords Result
256  * @return MoFEMErrorCode
257  */
259  const double *elem_coords, const double *glob_coords, const int nb_nodes,
260  double *local_coords);
261 
262  /**
263  * @brief Get the Tets With Quality
264  *
265  * @param out_tets
266  * @param tets
267  * @param th
268  * @param f
269  * @return MoFEMErrorCode
270  */
272  getTetsWithQuality(Range &out_tets, const Range &tets, Tag th = nullptr,
273  boost::function<bool(double)> f = [](double q) -> bool {
274  if (q <= 0)
275  return true;
276  else
277  return false;
278  });
279  /**
280  * @brief Write file with tetrahedral of given quality
281  *
282  * @param file_name
283  * @param file_type
284  * @param options
285  * @param tets
286  * @param th
287  * @param f
288  * @return MoFEMErrorCode
289  */
291  writeTetsWithQuality(const char *file_name, const char *file_type,
292  const char *options, const Range &tets, Tag th = nullptr,
293  boost::function<bool(double)> f = [](double q) -> bool {
294  if (q <= 0)
295  return true;
296  else
297  return false;
298  });
299 
300  /**
301  * @brief Check of point is in tetrahedral
302  *
303  * @param tet_coords
304  * @param global_coord
305  * @param tol
306  * @param result
307  * @return MoFEMErrorCode
308  */
309  static MoFEMErrorCode checkIfPointIsInTet(const double tet_coords[],
310  const double global_coord[],
311  const double tol, bool &result);
312 
313  /**
314  * @brief Get the Tri Normal objectGet triangle normal
315  *
316  * @param coords
317  * @param normal
318  * @return MoFEMErrorCode
319  */
320  static MoFEMErrorCode getTriNormal(const double *coords, double *normal);
321 
322  /**
323  * @brief Get triangle normal
324  *
325  * @param tri
326  * @param normal
327  * @return MoFEMErrorCode
328  */
329  MoFEMErrorCode getTriNormal(const EntityHandle tri, double *normal) const;
330 
331  /**
332  * @brief Get triangle area
333  *
334  * @param tri
335  * @return double
336  */
337  double getTriArea(const EntityHandle tri) const;
338 
339  /**
340  * @brief Get edge length
341  *
342  * @param edge_coords
343  * @return double
344  */
345  static double getEdgeLength(const double *edge_coords);
346 
347  /**
348  * @brief Get edge length
349  *
350  * @param edge
351  * @return double
352  */
353  double getEdgeLength(const EntityHandle edge);
354 
361  };
362 
363  /**
364  * @brief Find closet point on the segment from the point
365  *
366  * @param w_ptr segment first vertex coordinate
367  * @param v_ptr segment second vertex coordinate
368  * @param p_ptr coordinate of point
369  * @param t_ptr distance on the segment
370  *
371  * \note If t is outside bounds [ 0,-1 ] point is on the line point
372  * beyond segment.
373  *
374  * \code
375  *
376  * double w[] = {-1, 0, 0};
377  * double v[] = {1, 0, 0};
378  * double p[] = {0, 1, 0};
379  * double t;
380  * CHKERR Toolas::minDistancePointFromOnSegment(w, v, p, &t);
381  * double point_on_segment[3];
382  * for (int i = 0; i != 3; ++i)
383  * point_on_segment[i] = w[i] + t * (v[i] - w[i]);
384  *
385  * \endcode
386  *
387  * @return SEGMENT_MIN_DISTANCE
388  */
389  static SEGMENT_MIN_DISTANCE
390  minDistancePointFromOnSegment(const double *w_ptr, const double *v_ptr,
391  const double *p_ptr,
392  double *const t_ptr = nullptr);
393  /**
394  * @brief Find points on two segments in closest distance
395  *
396  * @param w_ptr
397  * @param v_ptr
398  * @param k_ptr
399  * @param l_ptr
400  * @param tvw_ptr
401  * @param tlk_ptr
402  * @return SEGMENT_MIN_DISTANCE
403  *
404  * \note If tvwk or tlk are outside bound [0,-1], it means that points
405  * are on the lines beyond segments, respectively for segment vw and
406  * lk.
407  *
408  */
409  static SEGMENT_MIN_DISTANCE
410  minDistanceFromSegments(const double *w_ptr, const double *v_ptr,
411  const double *k_ptr, const double *l_ptr,
412  double *const tvw_ptr = nullptr,
413  double *const tlk_ptr = nullptr);
414 
415  /**
416  * @brief Find minimal distance to edges
417  *
418  * \note Finding only edges with have smaller distance than distance
419  * set on the input by min_dist_ptr
420  *
421  * @param v_ptr point coordinates
422  * @param nb nb points
423  * @param edges range of edges
424  * @param min_dist_ptr on return minimal distance, on input starting
425  * distance
426  * @param o_ptr coordinates of the point on edge
427  * @param o_segments closest segments
428  * @return MoFEMErrorCode
429  */
431  findMinDistanceFromTheEdges(const double *v_ptr, const int nb, Range edges,
432  double *min_dist_ptr, double *o_ptr = nullptr,
433  EntityHandle *o_segments = nullptr) const;
434 
435  /**@}*/
436 
437  /** \name Debugging */
438 
439  /**@{*/
440 
441  /** \brief Print all DOFs for which element of vector is not a number
442  *
443  */
445  const RowColData row_or_col, Vec v);
446 
447  /**@}*/
448 
450  const int edge0,
451  const int edge1);
452 };
453 
454 double Tools::shapeFunMBTRI0(const double x, const double y) {
455  return N_MBTRI0(x, y);
456 }
457 
458 double Tools::shapeFunMBTRI1(const double x, const double y) {
459  return N_MBTRI1(x, y);
460 }
461 
462 double Tools::shapeFunMBTRI2(const double x, const double y) {
463  return N_MBTRI2(x, y);
464 }
465 
466 template <int LDB>
467 MoFEMErrorCode Tools::shapeFunMBTRI(double *shape, const double *ksi,
468  const double *eta, const double nb) {
470  for (int n = 0; n != nb; ++n) {
471  shape[0] = shapeFunMBTRI0(*ksi, *eta);
472  shape[1] = shapeFunMBTRI1(*ksi, *eta);
473  shape[2] = shapeFunMBTRI2(*ksi, *eta);
474  shape += 3;
475  ksi += LDB;
476  eta += LDB;
477  }
479 }
480 
481 double Tools::shapeFunMBTET0(const double x, const double y, const double z) {
482  return N_MBTET0(x, y, z);
483 }
484 
485 double Tools::shapeFunMBTET1(const double x, const double y, const double z) {
486  return N_MBTET1(x, y, z);
487 }
488 
489 double Tools::shapeFunMBTET2(const double x, const double y, const double z) {
490  return N_MBTET2(x, y, z);
491 }
492 
493 double Tools::shapeFunMBTET3(const double x, const double y, const double z) {
494  return N_MBTET3(x, y, z);
495 };
496 
497 template <int LDB>
498 MoFEMErrorCode Tools::shapeFunMBTET(double *shape, const double *ksi,
499  const double *eta, const double *zeta,
500  const double nb) {
502  for (int n = 0; n != nb; ++n) {
503  shape[0] = shapeFunMBTET0(*ksi, *eta, *zeta);
504  shape[1] = shapeFunMBTET1(*ksi, *eta, *zeta);
505  shape[2] = shapeFunMBTET2(*ksi, *eta, *zeta);
506  shape[3] = shapeFunMBTET3(*ksi, *eta, *zeta);
507  shape += 4;
508  ksi += LDB;
509  eta += LDB;
510  zeta += LDB;
511  }
513 }
514 
515 } // namespace MoFEM
516 
517 #endif // __TOOLS_HPP__a
518 
519 /**
520  * \defgroup mofem_tools Tools interface
521  * \brief Interface for tools
522  *
523  * \ingroup mofem
524  */
static constexpr double diffShapeFunMBQUADAtCenter2y
derivative of quad shape function
Definition: Tools.hpp:128
#define diffN_MBTRI1y
derivative of triangle shape function
Definition: fem_tools.h:60
#define diffN_MBTET1z
derivative of tetrahedral shape function
Definition: fem_tools.h:45
#define diffN_MBTET0z
derivative of tetrahedral shape function
Definition: fem_tools.h:42
static constexpr double diffShapeFunMBTRI2y
derivative of triangle shape function
Definition: Tools.hpp:84
static constexpr double diffShapeFunMBQUADAtCenter3x
derivative of quad shape function
Definition: Tools.hpp:130
MoFEM interface unique ID.
static constexpr double shapeFunMBTET3AtOneThird
Definition: Tools.hpp:199
static constexpr double diffShapeFunMBTET1y
derivative of tetrahedral shape function
Definition: Tools.hpp:149
static constexpr double shapeFunMBTET1AtOneThird
Definition: Tools.hpp:195
#define diffN_MBTRI0x
derivative of triangle shape function
Definition: fem_tools.h:57
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition: Tools.cpp:34
static constexpr double diffShapeFunMBTET1x
derivative of tetrahedral shape function
Definition: Tools.hpp:147
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.
Definition: Tools.cpp:411
#define diffN_MBQUAD3x(y)
Definition: fem_tools.h:75
#define diffN_MBTRI1x
derivative of triangle shape function
Definition: fem_tools.h:59
static double shapeFunMBTRI1(const double x, const double y)
Definition: Tools.hpp:458
static constexpr double diffShapeFunMBTRI0y
derivative of triangle shape function
Definition: Tools.hpp:76
static constexpr double diffShapeFunMBTET1z
derivative of tetrahedral shape function
Definition: Tools.hpp:151
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
static constexpr double shapeFunMBTET1At000
Definition: Tools.hpp:189
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
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.
Definition: Tools.cpp:345
static constexpr double diffShapeFunMBTRI1x
derivative of triangle shape function
Definition: Tools.hpp:78
static constexpr std::array< double, 2 > diffShapeFunMBEDGE
Definition: Tools.hpp:71
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
#define diffN_MBTET2x
derivative of tetrahedral shape function
Definition: fem_tools.h:46
base class for all interface classes
static constexpr double diffShapeFunMBTET0y
derivative of tetrahedral shape function
Definition: Tools.hpp:143
FTensor::Tensor1< double, 2 > normal(FTensor::Tensor1< T1, 3 > &t_coords, FTensor::Tensor1< T2, 2 > &t_disp)
Definition: ContactOps.hpp:137
static constexpr std::array< double, 12 > diffShapeFunMBTET
Definition: Tools.hpp:166
double getTriArea(const EntityHandle tri) const
Get triangle area.
Definition: Tools.cpp:283
#define diffN_MBTET3y
derivative of tetrahedral shape function
Definition: fem_tools.h:50
Auxiliary tools.
Definition: Tools.hpp:30
#define diffN_MBTRI2x
derivative of triangle shape function
Definition: fem_tools.h:61
static MoFEMErrorCode shapeFunMBTRI(double *shape, const double *ksi, const double *eta, const double nb)
Calculate shape functions on triangle.
Definition: Tools.hpp:467
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:97
static constexpr double diffShapeFunMBTET2y
derivative of tetrahedral shape function
Definition: Tools.hpp:155
static constexpr double diffShapeFunMBTET2x
derivative of tetrahedral shape function
Definition: Tools.hpp:153
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
#define diffN_MBTET1x
derivative of tetrahedral shape function
Definition: fem_tools.h:43
static double shapeFunMBTRI2(const double x, const double y)
Definition: Tools.hpp:462
#define diffN_MBQUAD0x(y)
Definition: fem_tools.h:69
keeps basic data about problemThis is low level structure with information about problem,...
static constexpr double diffShapeFunMBQUADAtCenter3y
derivative of quad shape function
Definition: Tools.hpp:132
static constexpr double shapeFunMBTET0AtOneThird
Definition: Tools.hpp:193
static constexpr double diffN_MBEDGE0x
Definition: Tools.hpp:68
#define diffN_MBQUAD0y(x)
Definition: fem_tools.h:70
RowColData
RowColData.
Definition: definitions.h:192
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.
Definition: Tools.cpp:174
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
static constexpr double diffShapeFunMBTET0x
derivative of tetrahedral shape function
Definition: Tools.hpp:141
static constexpr double diffShapeFunMBTET2z
derivative of tetrahedral shape function
Definition: Tools.hpp:157
static constexpr double diffShapeFunMBTET3z
derivative of tetrahedral shape function
Definition: Tools.hpp:163
static constexpr double shapeFunMBTET3At000
Definition: Tools.hpp:191
#define diffN_MBEDGE1
derivative of edge shape function
Definition: fem_tools.h:82
static constexpr std::array< double, 4 > shapeFunMBTETAt000
Array of shape function at zero local point on reference element.
Definition: Tools.hpp:225
static MoFEMErrorCode checkIfPointIsInTet(const double tet_coords[], const double global_coord[], const double tol, bool &result)
Check of point is in tetrahedral.
Definition: Tools.cpp:192
#define diffN_MBQUAD2y(x)
Definition: fem_tools.h:74
static Index< 'n', 3 > n
static double tetVolume(const double *coords)
Calculate volume of tetrahedron.
Definition: Tools.cpp:49
#define diffN_MBTET3z
derivative of tetrahedral shape function
Definition: fem_tools.h:51
#define diffN_MBTET0y
derivative of tetrahedral shape function
Definition: fem_tools.h:41
#define diffN_MBQUAD2x(y)
Definition: fem_tools.h:73
static constexpr double diffShapeFunMBQUADAtCenter2x
derivative of quad shape function
Definition: Tools.hpp:126
static constexpr double shapeFunMBTET2AtOneThird
Definition: Tools.hpp:197
#define diffN_MBTET2y
derivative of tetrahedral shape function
Definition: fem_tools.h:47
double tol
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
Definition: Tools.cpp:68
static double shapeFunMBTET3(const double x, const double y, const double z)
Definition: Tools.hpp:493
static constexpr double diffN_MBEDGE1x
Definition: Tools.hpp:69
static double shapeFunMBTRI0(const double x, const double y)
Definition: Tools.hpp:454
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
static constexpr double shapeFunMBTET0At000
Definition: Tools.hpp:188
static constexpr double diffShapeFunMBTRI1y
derivative of triangle shape function
Definition: Tools.hpp:80
#define diffN_MBTET3x
derivative of tetrahedral shape function
Definition: fem_tools.h:49
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.
Definition: Tools.cpp:212
static MoFEMErrorCode getTriNormal(const double *coords, double *normal)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:258
static double shapeFunMBTET1(const double x, const double y, const double z)
Definition: Tools.hpp:485
#define N_MBTET3(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:39
#define diffN_MBEDGE0
derivative of edge shape function
Definition: fem_tools.h:81
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: Tools.cpp:22
static const double edge_coords[6][6]
static const MOFEMuuid IDD_MOFEMTools
Definition: Tools.hpp:23
static constexpr std::array< double, 8 > diffShapeFunMBQUADAtCenter
Definition: Tools.hpp:135
static constexpr double diffShapeFunMBTRI0x
derivative of triangle shape function
Definition: Tools.hpp:74
static constexpr double diffShapeFunMBTET3x
derivative of tetrahedral shape function
Definition: Tools.hpp:159
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition: Tools.hpp:87
#define N_MBTET0(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:36
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
static constexpr double diffShapeFunMBTET0z
derivative of tetrahedral shape function
Definition: Tools.hpp:145
static constexpr double shapeFunMBTET2At000
Definition: Tools.hpp:190
#define N_MBTRI2(x, y)
triangle shape function
Definition: fem_tools.h:56
static constexpr double diffShapeFunMBQUADAtCenter1y
derivative of quad shape function
Definition: Tools.hpp:124
static MoFEMErrorCode outerProductOfEdgeIntegrationPtsForQuad(MatrixDouble &pts, const int edge0, const int edge1)
Definition: Tools.cpp:486
#define N_MBTET2(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:38
#define diffN_MBQUAD1x(y)
Definition: fem_tools.h:71
#define diffN_MBTRI2y
derivative of triangle shape function
Definition: fem_tools.h:62
#define diffN_MBTRI0y
derivative of triangle shape function
Definition: fem_tools.h:58
static double shapeFunMBTET2(const double x, const double y, const double z)
Definition: Tools.hpp:489
static double shapeFunMBTET0(const double x, const double y, const double z)
Definition: Tools.hpp:481
static double getEdgeLength(const double *edge_coords)
Get edge length.
Definition: Tools.cpp:291
static constexpr double diffShapeFunMBQUADAtCenter0x
derivative of quad shape function
Definition: Tools.hpp:118
Core (interface) class.
Definition: Core.hpp:77
static constexpr double diffShapeFunMBTRI2x
derivative of triangle shape function
Definition: Tools.hpp:82
static constexpr double diffShapeFunMBQUADAtCenter0y
derivative of quad shape function
Definition: Tools.hpp:120
Tools(const MoFEM::Core &core)
Definition: Tools.hpp:36
#define N_MBTET1(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:37
static constexpr double diffShapeFunMBTET3y
derivative of tetrahedral shape function
Definition: Tools.hpp:161
#define diffN_MBTET1y
derivative of tetrahedral shape function
Definition: fem_tools.h:44
#define diffN_MBTET2z
derivative of tetrahedral shape function
Definition: fem_tools.h:48
#define diffN_MBQUAD3y(x)
Definition: fem_tools.h:76
#define N_MBTRI1(x, y)
triangle shape function
Definition: fem_tools.h:55
#define diffN_MBTET0x
derivative of tetrahedral shape function
Definition: fem_tools.h:40
MoFEM::Core & cOre
Definition: Tools.hpp:35
static MoFEMErrorCode shapeFunMBTET(double *shape, const double *ksi, const double *eta, const double *zeta, const double nb)
Calculate shape functions on tetrahedron.
Definition: Tools.hpp:498
static constexpr std::array< double, 4 > shapeFunMBTETAtOneThird
Array of shape function at center on reference element.
Definition: Tools.hpp:233
#define N_MBTRI0(x, y)
triangle shape function
Definition: fem_tools.h:54
static constexpr double diffShapeFunMBQUADAtCenter1x
derivative of quad shape function
Definition: Tools.hpp:122
SEGMENT_MIN_DISTANCE
Definition: Tools.hpp:355
#define diffN_MBQUAD1y(x)
Definition: fem_tools.h:72