v0.9.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 constexpr double diffShapeFunMBQUADAtCenter0x =
96  diffN_MBQUAD0x(0.5); ///< derivative of quad shape function
97  static constexpr double diffShapeFunMBQUADAtCenter0y =
98  diffN_MBQUAD0y(0.5); ///< derivative of quad shape function
99  static constexpr double diffShapeFunMBQUADAtCenter1x =
100  diffN_MBQUAD1x(0.5); ///< derivative of quad shape function
101  static constexpr double diffShapeFunMBQUADAtCenter1y =
102  diffN_MBQUAD1y(0.5); ///< derivative of quad shape function
103  static constexpr double diffShapeFunMBQUADAtCenter2x =
104  diffN_MBQUAD2x(0.5); ///< derivative of quad shape function
105  static constexpr double diffShapeFunMBQUADAtCenter2y =
106  diffN_MBQUAD2y(0.5); ///< derivative of quad shape function
107  static constexpr double diffShapeFunMBQUADAtCenter3x =
108  diffN_MBQUAD3x(0.5); ///< derivative of quad shape function
109  static constexpr double diffShapeFunMBQUADAtCenter3y =
110  diffN_MBQUAD3y(0.5); ///< derivative of quad shape function
111 
112  static constexpr std::array<double, 8> diffShapeFunMBQUADAtCenter = {
117 
118  static constexpr double diffShapeFunMBTET0x =
119  diffN_MBTET0x; ///< derivative of tetrahedral shape function
120  static constexpr double diffShapeFunMBTET0y =
121  diffN_MBTET0y; ///< derivative of tetrahedral shape function
122  static constexpr double diffShapeFunMBTET0z =
123  diffN_MBTET0z; ///< derivative of tetrahedral shape function
124  static constexpr double diffShapeFunMBTET1x =
125  diffN_MBTET1x; ///< derivative of tetrahedral shape function
126  static constexpr double diffShapeFunMBTET1y =
127  diffN_MBTET1y; ///< derivative of tetrahedral shape function
128  static constexpr double diffShapeFunMBTET1z =
129  diffN_MBTET1z; ///< derivative of tetrahedral shape function
130  static constexpr double diffShapeFunMBTET2x =
131  diffN_MBTET2x; ///< derivative of tetrahedral shape function
132  static constexpr double diffShapeFunMBTET2y =
133  diffN_MBTET2y; ///< derivative of tetrahedral shape function
134  static constexpr double diffShapeFunMBTET2z =
135  diffN_MBTET2z; ///< derivative of tetrahedral shape function
136  static constexpr double diffShapeFunMBTET3x =
137  diffN_MBTET3x; ///< derivative of tetrahedral shape function
138  static constexpr double diffShapeFunMBTET3y =
139  diffN_MBTET3y; ///< derivative of tetrahedral shape function
140  static constexpr double diffShapeFunMBTET3z =
141  diffN_MBTET3z; ///< derivative of tetrahedral shape function
142 
143  static constexpr std::array<double, 12> diffShapeFunMBTET = {
144 
146 
148 
150 
152 
153  static inline double shapeFunMBTET0(const double x, const double y,
154  const double z);
155 
156  static inline double shapeFunMBTET1(const double x, const double y,
157  const double z);
158 
159  static inline double shapeFunMBTET2(const double x, const double y,
160  const double z);
161 
162  static inline double shapeFunMBTET3(const double x, const double y,
163  const double z);
164 
165  static constexpr double shapeFunMBTET0At000 = N_MBTET0(0, 0, 0);
166  static constexpr double shapeFunMBTET1At000 = N_MBTET1(0, 0, 0);
167  static constexpr double shapeFunMBTET2At000 = N_MBTET2(0, 0, 0);
168  static constexpr double shapeFunMBTET3At000 = N_MBTET3(0, 0, 0);
169 
170  static constexpr double shapeFunMBTET0AtOneThird =
171  N_MBTET0(1. / 3., 1. / 3., 1. / 3.);
172  static constexpr double shapeFunMBTET1AtOneThird =
173  N_MBTET1(1. / 3., 1. / 3., 1. / 3.);
174  static constexpr double shapeFunMBTET2AtOneThird =
175  N_MBTET2(1. / 3., 1. / 3., 1. / 3.);
176  static constexpr double shapeFunMBTET3AtOneThird =
177  N_MBTET3(1. / 3., 1. / 3., 1. / 3.);
178 
179  /**
180  * @brief Calculate shape functions on tetrahedron
181  *
182  * \note Template parameter is leading dimension of point coordinate arrays,
183  * such that \f$ksi_{n+1} = ksi[n + LDB]\f$
184  *
185  * @tparam 1
186  * @param shape shape functions
187  * @param ksi pointer to first local coordinates
188  * @param eta pointer to second local coordinates
189  * @param zeta pointer to first third coordinates
190  * @param nb number of points
191  * @return MoFEMErrorCode
192  */
193  template <int LDB = 1>
194  static MoFEMErrorCode shapeFunMBTET(double *shape, const double *ksi,
195  const double *eta, const double *zeta,
196  const double nb);
197 
198  /**
199  * @brief Array of shape function at zero local point on reference element
200  *
201  */
202  static constexpr std::array<double, 4> shapeFunMBTETAt000 = {
205 
206  /**
207  * @brief Array of shape function at center on reference element
208  *
209  */
210  static constexpr std::array<double, 4> shapeFunMBTETAtOneThird = {
213 
214 
215  /**
216  * @brief Get the Local Coordinates On Reference Four Node Tet object
217  *
218  * \code
219  * MatrixDouble elem_coords(4, 3);
220  * // Set nodal coordinates
221  * MatrixDouble global_coords(5, 3);
222  * // Set global coordinates
223  * MatrixDouble local_coords(global_coords.size1(), 3);
224  * CHKERR Tools::getLocalCoordinatesOnReferenceFourNodeTet(
225  * &elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
226  * &local_coords(0, 0))
227  * \endcode
228  *
229  * @param elem_coords Global element node coordinates
230  * @param glob_coords Globale coordinates
231  * @param nb_nodes Number of points
232  * @param local_coords Result
233  * @return MoFEMErrorCode
234  */
236  const double *elem_coords, const double *glob_coords, const int nb_nodes,
237  double *local_coords);
238 
239  /**
240  * @brief Get the Tets With Quality
241  *
242  * @param out_tets
243  * @param tets
244  * @param th
245  * @param f
246  * @return MoFEMErrorCode
247  */
249  getTetsWithQuality(Range &out_tets, const Range &tets, Tag th = nullptr,
250  boost::function<bool(double)> f = [](double q) -> bool {
251  if (q <= 0)
252  return true;
253  else
254  return false;
255  });
256  /**
257  * @brief Write file with tetrahedral of given quality
258  *
259  * @param file_name
260  * @param file_type
261  * @param options
262  * @param tets
263  * @param th
264  * @param f
265  * @return MoFEMErrorCode
266  */
268  writeTetsWithQuality(const char *file_name, const char *file_type,
269  const char *options, const Range &tets, Tag th = nullptr,
270  boost::function<bool(double)> f = [](double q) -> bool {
271  if (q <= 0)
272  return true;
273  else
274  return false;
275  });
276 
277  /**
278  * @brief Check of point is in tetrahedral
279  *
280  * @param tet_coords
281  * @param global_coord
282  * @param tol
283  * @param result
284  * @return MoFEMErrorCode
285  */
286  static MoFEMErrorCode checkIfPointIsInTet(const double tet_coords[],
287  const double global_coord[],
288  const double tol, bool &result);
289 
290  /**
291  * @brief Get the Tri Normal objectGet triangle normal
292  *
293  * @param coords
294  * @param normal
295  * @return MoFEMErrorCode
296  */
297  static MoFEMErrorCode getTriNormal(const double *coords, double *normal);
298 
299  /**
300  * @brief Get triangle normal
301  *
302  * @param tri
303  * @param normal
304  * @return MoFEMErrorCode
305  */
306  MoFEMErrorCode getTriNormal(const EntityHandle tri, double *normal) const;
307 
308  /**
309  * @brief Get triangle area
310  *
311  * @param tri
312  * @return double
313  */
314  double getTriArea(const EntityHandle tri) const;
315 
316  /**
317  * @brief Get edge length
318  *
319  * @param edge_coords
320  * @return double
321  */
322  static double getEdgeLength(const double *edge_coords);
323 
324  /**
325  * @brief Get edge length
326  *
327  * @param edge
328  * @return double
329  */
330  double getEdgeLength(const EntityHandle edge);
331 
338  };
339 
340  /**
341  * @brief Find closet point on the segment from the point
342  *
343  * @param w_ptr segment first vertex coordinate
344  * @param v_ptr segment second vertex coordinate
345  * @param p_ptr coordinate of point
346  * @param t_ptr distance on the segment
347  *
348  * \note If t is outside bounds [ 0,-1 ] point is on the line point
349  * beyond segment.
350  *
351  * \code
352  *
353  * double w[] = {-1, 0, 0};
354  * double v[] = {1, 0, 0};
355  * double p[] = {0, 1, 0};
356  * double t;
357  * CHKERR Toolas::minDistancePointFromOnSegment(w, v, p, &t);
358  * double point_on_segment[3];
359  * for (int i = 0; i != 3; ++i)
360  * point_on_segment[i] = w[i] + t * (v[i] - w[i]);
361  *
362  * \endcode
363  *
364  * @return SEGMENT_MIN_DISTANCE
365  */
366  static SEGMENT_MIN_DISTANCE
367  minDistancePointFromOnSegment(const double *w_ptr, const double *v_ptr,
368  const double *p_ptr,
369  double *const t_ptr = nullptr);
370  /**
371  * @brief Find points on two segments in closest distance
372  *
373  * @param w_ptr
374  * @param v_ptr
375  * @param k_ptr
376  * @param l_ptr
377  * @param tvw_ptr
378  * @param tlk_ptr
379  * @return SEGMENT_MIN_DISTANCE
380  *
381  * \note If tvwk or tlk are outside bound [0,-1], it means that points
382  * are on the lines beyond segments, respectively for segment vw and
383  * lk.
384  *
385  */
386  static SEGMENT_MIN_DISTANCE
387  minDistanceFromSegments(const double *w_ptr, const double *v_ptr,
388  const double *k_ptr, const double *l_ptr,
389  double *const tvw_ptr = nullptr,
390  double *const tlk_ptr = nullptr);
391 
392  /**
393  * @brief Find minimal distance to edges
394  *
395  * \note Finding only edges with have smaller distance than distance
396  * set on the input by min_dist_ptr
397  *
398  * @param v_ptr point coordinates
399  * @param nb nb points
400  * @param edges range of edges
401  * @param min_dist_ptr on return minimal distance, on input starting
402  * distance
403  * @param o_ptr coordinates of the point on edge
404  * @param o_segments closest segments
405  * @return MoFEMErrorCode
406  */
408  findMinDistanceFromTheEdges(const double *v_ptr, const int nb, Range edges,
409  double *min_dist_ptr, double *o_ptr = nullptr,
410  EntityHandle *o_segments = nullptr) const;
411 
412  /**@}*/
413 
414  /** \name Debugging */
415 
416  /**@{*/
417 
418  /** \brief Print all DOFs for which element of vector is not a number
419  *
420  */
422  const RowColData row_or_col, Vec v);
423 
424  /**@}*/
425 
427  const int edge0,
428  const int edge1);
429 };
430 
431 double Tools::shapeFunMBTET0(const double x, const double y, const double z) {
432  return N_MBTET0(x, y, z);
433 }
434 
435 double Tools::shapeFunMBTET1(const double x, const double y, const double z) {
436  return N_MBTET1(x, y, z);
437 }
438 
439 double Tools::shapeFunMBTET2(const double x, const double y, const double z) {
440  return N_MBTET2(x, y, z);
441 }
442 
443 double Tools::shapeFunMBTET3(const double x, const double y, const double z) {
444  return N_MBTET3(x, y, z);
445 };
446 
447 template <int LDB>
448 MoFEMErrorCode Tools::shapeFunMBTET(double *shape, const double *ksi,
449  const double *eta, const double *zeta,
450  const double nb) {
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 }
464 
465 } // namespace MoFEM
466 
467 #endif // __TOOLS_HPP__a
468 
469 /**
470  * \defgroup mofem_tools Tools interface
471  * \brief Interface for tools
472  *
473  * \ingroup mofem
474  */
static constexpr double diffShapeFunMBQUADAtCenter2y
derivative of quad shape function
Definition: Tools.hpp:105
#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:107
MoFEM interface unique ID.
static constexpr double shapeFunMBTET3AtOneThird
Definition: Tools.hpp:176
static constexpr double diffShapeFunMBTET1y
derivative of tetrahedral shape function
Definition: Tools.hpp:126
static constexpr double shapeFunMBTET1AtOneThird
Definition: Tools.hpp:172
#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:124
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 constexpr double diffShapeFunMBTRI0y
derivative of triangle shape function
Definition: Tools.hpp:76
static constexpr double diffShapeFunMBTET1z
derivative of tetrahedral shape function
Definition: Tools.hpp:128
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:166
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
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:74
#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:120
static constexpr std::array< double, 12 > diffShapeFunMBTET
Definition: Tools.hpp:143
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 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
Core (interface) class.
Definition: Core.hpp:50
static constexpr double diffShapeFunMBTET2y
derivative of tetrahedral shape function
Definition: Tools.hpp:132
static constexpr double diffShapeFunMBTET2x
derivative of tetrahedral shape function
Definition: Tools.hpp:130
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
#define diffN_MBTET1x
derivative of tetrahedral shape function
Definition: fem_tools.h:43
#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:109
static constexpr double shapeFunMBTET0AtOneThird
Definition: Tools.hpp:170
static constexpr double diffN_MBEDGE0x
Definition: Tools.hpp:68
#define diffN_MBQUAD0y(x)
Definition: fem_tools.h:70
RowColData
RowColData.
Definition: definitions.h:186
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:118
static constexpr double diffShapeFunMBTET2z
derivative of tetrahedral shape function
Definition: Tools.hpp:134
static constexpr double diffShapeFunMBTET3z
derivative of tetrahedral shape function
Definition: Tools.hpp:140
static constexpr double shapeFunMBTET3At000
Definition: Tools.hpp:168
#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:202
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 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:103
static constexpr double shapeFunMBTET2AtOneThird
Definition: Tools.hpp:174
#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:443
static constexpr double diffN_MBEDGE1x
Definition: Tools.hpp:69
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
static constexpr double shapeFunMBTET0At000
Definition: Tools.hpp:165
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:435
#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:112
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:136
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:122
static constexpr double shapeFunMBTET2At000
Definition: Tools.hpp:167
static constexpr double diffShapeFunMBQUADAtCenter1y
derivative of quad shape function
Definition: Tools.hpp:101
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:439
static double shapeFunMBTET0(const double x, const double y, const double z)
Definition: Tools.hpp:431
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:95
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:97
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:138
#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 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:448
static constexpr std::array< double, 4 > shapeFunMBTETAtOneThird
Array of shape function at center on reference element.
Definition: Tools.hpp:210
static constexpr double diffShapeFunMBQUADAtCenter1x
derivative of quad shape function
Definition: Tools.hpp:99
SEGMENT_MIN_DISTANCE
Definition: Tools.hpp:332
#define diffN_MBQUAD1y(x)
Definition: fem_tools.h:72