v0.8.23
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 
42  /**
43  * \brief Calculate tetrahedron volume length quality
44  * @param coords tet coordinates
45  * @return Volume-length quality
46  */
47  static double volumeLengthQuality(const double *coords);
48 
49  /**
50  * @brief Calculate volume of tetrahedron
51  *
52  * @param coords
53  * @return double volume
54  */
55  static double tetVolume(const double *coords);
56 
57  /**
58  * \brief calculate minimal quality of tetrahedra in range
59  * @param tets range
60  * @param min_quality mimimal quality
61  * @return error code
62  */
63  MoFEMErrorCode minTetsQuality(const Range &tets, double &min_quality,
64  Tag th = nullptr,
65  boost::function<double(double, double)> f =
66  [](double a, double b) -> double {
67  return std::min(a, b);
68  });
69 
70  static constexpr double diffShapeFunMBTRI0x =
71  diffN_MBTRI0x; ///< derivative of triangle shape function
72  static constexpr double diffShapeFunMBTRI0y =
73  diffN_MBTRI0y; ///< derivative of triangle shape function
74  static constexpr double diffShapeFunMBTRI1x =
75  diffN_MBTRI1x; ///< derivative of triangle shape function
76  static constexpr double diffShapeFunMBTRI1y =
77  diffN_MBTRI1y; ///< derivative of triangle shape function
78  static constexpr double diffShapeFunMBTRI2x =
79  diffN_MBTRI2x; ///< derivative of triangle shape function
80  static constexpr double diffShapeFunMBTRI2y =
81  diffN_MBTRI2y; ///< derivative of triangle shape function
82 
83  static constexpr std::array<double, 6> diffShapeFunMBTRI = {
84 
86 
88 
90 
91  static constexpr double diffShapeFunMBTET0x =
92  diffN_MBTET0x; ///< derivative of tetrahedral shape function
93  static constexpr double diffShapeFunMBTET0y =
94  diffN_MBTET0y; ///< derivative of tetrahedral shape function
95  static constexpr double diffShapeFunMBTET0z =
96  diffN_MBTET0z; ///< derivative of tetrahedral shape function
97  static constexpr double diffShapeFunMBTET1x =
98  diffN_MBTET1x; ///< derivative of tetrahedral shape function
99  static constexpr double diffShapeFunMBTET1y =
100  diffN_MBTET1y; ///< derivative of tetrahedral shape function
101  static constexpr double diffShapeFunMBTET1z =
102  diffN_MBTET1z; ///< derivative of tetrahedral shape function
103  static constexpr double diffShapeFunMBTET2x =
104  diffN_MBTET2x; ///< derivative of tetrahedral shape function
105  static constexpr double diffShapeFunMBTET2y =
106  diffN_MBTET2y; ///< derivative of tetrahedral shape function
107  static constexpr double diffShapeFunMBTET2z =
108  diffN_MBTET2z; ///< derivative of tetrahedral shape function
109  static constexpr double diffShapeFunMBTET3x =
110  diffN_MBTET3x; ///< derivative of tetrahedral shape function
111  static constexpr double diffShapeFunMBTET3y =
112  diffN_MBTET3y; ///< derivative of tetrahedral shape function
113  static constexpr double diffShapeFunMBTET3z =
114  diffN_MBTET3z; ///< derivative of tetrahedral shape function
115 
116  static constexpr std::array<double, 12> diffShapeFunMBTET = {
117 
119 
121 
123 
125 
126  static inline double shapeFunMBTET0(const double x, const double y,
127  const double z) {
128  return N_MBTET0(x, y, z);
129  }
130 
131  static inline double shapeFunMBTET1(const double x, const double y,
132  const double z) {
133  return N_MBTET1(x, y, z);
134  }
135 
136  static inline double shapeFunMBTET2(const double x, const double y,
137  const double z) {
138  return N_MBTET2(x, y, z);
139  }
140 
141  static inline double shapeFunMBTET3(const double x, const double y,
142  const double z) {
143  return N_MBTET3(x, y, z);
144  ;
145  };
146 
147  static constexpr double shapeFunMBTET0At000 = N_MBTET0(0, 0, 0);
148  static constexpr double shapeFunMBTET1At000 = N_MBTET1(0, 0, 0);
149  static constexpr double shapeFunMBTET2At000 = N_MBTET2(0, 0, 0);
150  static constexpr double shapeFunMBTET3At000 = N_MBTET3(0, 0, 0);
151 
152  static constexpr double shapeFunMBTET0AtOneThird =
153  N_MBTET0(1. / 3., 1. / 3., 1. / 3.);
154  static constexpr double shapeFunMBTET1AtOneThird =
155  N_MBTET1(1. / 3., 1. / 3., 1. / 3.);
156  static constexpr double shapeFunMBTET2AtOneThird =
157  N_MBTET2(1. / 3., 1. / 3., 1. / 3.);
158  static constexpr double shapeFunMBTET3AtOneThird =
159  N_MBTET3(1. / 3., 1. / 3., 1. / 3.);
160 
161  /**
162  * @brief Calculate shape functions on tetrahedron
163  *
164  * \note Template parameter is leading dimension of point coordinate arrays,
165  * such that \f$ksi_{n+1} = ksi[n + LDB]\f$
166  *
167  * @tparam 1
168  * @param shape shape functions
169  * @param ksi pointer to first local coordinates
170  * @param eta pointer to second local coordinates
171  * @param zeta pointer to first third coordinates
172  * @param nb number of points
173  * @return MoFEMErrorCode
174  */
175  template <int LDB = 1>
176  static MoFEMErrorCode shapeFunMBTET(double *shape, const double *ksi,
177  const double *eta, const double *zeta,
178  const double nb) {
180  for (int n = 0; n != nb; ++n) {
181  shape[0] = shapeFunMBTET0(*ksi, *eta, *zeta);
182  shape[1] = shapeFunMBTET1(*ksi, *eta, *zeta);
183  shape[2] = shapeFunMBTET2(*ksi, *eta, *zeta);
184  shape[3] = shapeFunMBTET3(*ksi, *eta, *zeta);
185  shape += 4;
186  ksi += LDB;
187  eta += LDB;
188  zeta += LDB;
189  }
191  }
192 
193  /**
194  * @brief Array of shape function at zero local point on reference element
195  *
196  */
197  static constexpr std::array<double, 4> shapeFunMBTETAt000 = {
200 
201  /**
202  * @brief Array of shape function at center on reference element
203  *
204  */
205  static constexpr std::array<double, 4> shapeFunMBTETAtOneThird = {
208 
209  /**
210  * @brief Get the Local Coordinates On Reference Four Node Tet object
211  *
212  * \code
213  * MatrixDouble elem_coords(4, 3);
214  * // Set nodal coordinates
215  * MatrixDouble global_coords(5, 3);
216  * // Set global coordinates
217  * MatrixDouble local_coords(global_coords.size1(), 3);
218  * CHKERR Tools::getLocalCoordinatesOnReferenceFourNodeTet(
219  * &elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
220  * &local_coords(0, 0))
221  * \endcode
222  *
223  * @param elem_coords Global element node coordinates
224  * @param glob_coords Globale coordinates
225  * @param nb_nodes Number of points
226  * @param local_coords Result
227  * @return MoFEMErrorCode
228  */
230  const double *elem_coords, const double *glob_coords, const int nb_nodes,
231  double *local_coords);
232 
233  /**
234  * @brief Get the Tets With Quality
235  *
236  * @param out_tets
237  * @param tets
238  * @param th
239  * @param f
240  * @return MoFEMErrorCode
241  */
243  getTetsWithQuality(Range &out_tets, const Range &tets, Tag th = nullptr,
244  boost::function<bool(double)> f = [](double q) -> bool {
245  if (q <= 0)
246  return true;
247  else
248  return false;
249  });
250  /**
251  * @brief Write file with tetrahedral of given quality
252  *
253  * @param file_name
254  * @param file_type
255  * @param options
256  * @param tets
257  * @param th
258  * @param f
259  * @return MoFEMErrorCode
260  */
262  writeTetsWithQuality(const char *file_name, const char *file_type,
263  const char *options, const Range &tets, Tag th = nullptr,
264  boost::function<bool(double)> f = [](double q) -> bool {
265  if (q <= 0)
266  return true;
267  else
268  return false;
269  });
270 
271  /**
272  * @brief Check of point is in tetrahedral
273  *
274  * @param tet_coords
275  * @param global_coord
276  * @param tol
277  * @param result
278  * @return MoFEMErrorCode
279  */
280  static MoFEMErrorCode checkIfPointIsInTet(const double tet_coords[],
281  const double global_coord[],
282  const double tol, bool &result);
283 
284  /**
285  * @brief Get the Tri Normal objectGet triangle normal
286  *
287  * @param coords
288  * @param normal
289  * @return MoFEMErrorCode
290  */
291  static MoFEMErrorCode getTriNormal(const double *coords, double *normal);
292 
293  /**
294  * @brief Get triangle normal
295  *
296  * @param tri
297  * @param normal
298  * @return MoFEMErrorCode
299  */
300  MoFEMErrorCode getTriNormal(const EntityHandle tri, double *normal) const;
301 
302  /**
303  * @brief Get triangle area
304  *
305  * @param tri
306  * @return double
307  */
308  double getTriArea(const EntityHandle tri) const;
309 
310  /**
311  * @brief Get edge length
312  *
313  * @param edge_coords
314  * @return double
315  */
316  static double getEdgeLength(const double *edge_coords);
317 
318  /**
319  * @brief Get edge length
320  *
321  * @param edge
322  * @return double
323  */
324  double getEdgeLength(const EntityHandle edge);
325 
332  };
333 
334  /**
335  * @brief Find closet point on the segment from the point
336  *
337  * @param w_ptr segment first vertex coordinate
338  * @param v_ptr segment second vertex coordinate
339  * @param p_ptr coordinate of point
340  * @param t_ptr distance on the segment
341  *
342  * \note If t is outside bounds [ 0,-1 ] point is on the line point
343  * beyond segment.
344  *
345  * \code
346  *
347  * double w[] = {-1, 0, 0};
348  * double v[] = {1, 0, 0};
349  * double p[] = {0, 1, 0};
350  * double t;
351  * CHKERR Toolas::minDistancePointFromOnSegment(w, v, p, &t);
352  * double point_on_segment[3];
353  * for (int i = 0; i != 3; ++i)
354  * point_on_segment[i] = w[i] + t * (v[i] - w[i]);
355  *
356  * \endcode
357  *
358  * @return SEGMENT_MIN_DISTANCE
359  */
360  static SEGMENT_MIN_DISTANCE
361  minDistancePointFromOnSegment(const double *w_ptr, const double *v_ptr,
362  const double *p_ptr,
363  double *const t_ptr = nullptr);
364  /**
365  * @brief Find points on two segments in closest distance
366  *
367  * @param w_ptr
368  * @param v_ptr
369  * @param k_ptr
370  * @param l_ptr
371  * @param tvw_ptr
372  * @param tlk_ptr
373  * @return SEGMENT_MIN_DISTANCE
374  *
375  * \note If tvwk or tlk are outside bound [0,-1], it means that points
376  * are on the lines beyond segments, respectively for segment vw and
377  * lk.
378  *
379  */
380  static SEGMENT_MIN_DISTANCE
381  minDistanceFromSegments(const double *w_ptr, const double *v_ptr,
382  const double *k_ptr, const double *l_ptr,
383  double *const tvw_ptr = nullptr,
384  double *const tlk_ptr = nullptr);
385 
386  /**
387  * @brief Find minimal distance to edges
388  *
389  * \note Finding only edges with have smaller distance than distance
390  * set on the input by min_dist_ptr
391  *
392  * @param v_ptr point coordinates
393  * @param nb nb points
394  * @param edges range of edges
395  * @param min_dist_ptr on return minimal distance, on input starting
396  * distance
397  * @param o_ptr coordinates of the point on edge
398  * @param o_segments closest segments
399  * @return MoFEMErrorCode
400  */
402  findMinDistanceFromTheEdges(const double *v_ptr, const int nb, Range edges,
403  double *min_dist_ptr, double *o_ptr,
404  EntityHandle *o_segments = nullptr) const;
405 
406  /**@}*/
407 
408  /** \name Debugging */
409 
410  /**@{*/
411 
412  /** \brief Print all DOFs for which element of vector is not a number
413  *
414  */
416  const RowColData row_or_col, Vec v);
417 
418  /**@}*/
419 };
420 
421 } // namespace MoFEM
422 
423 #endif // __TOOLS_HPP__a
424 
425 /**
426  * \defgroup mofem_tools Tools interface
427  * \brief Interface for tools
428  *
429  * \ingroup mofem
430  */
#define diffN_MBTRI1y
derivative of triangle shape function
Definition: fem_tools.h:60
static double shapeFunMBTET2(const double x, const double y, const double z)
Definition: Tools.hpp:136
#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:80
MoFEMErrorCode findMinDistanceFromTheEdges(const double *v_ptr, const int nb, Range edges, double *min_dist_ptr, double *o_ptr, EntityHandle *o_segments=nullptr) const
Find minimal distance to edges.
Definition: Tools.cpp:407
MoFEM interface unique ID.
static constexpr double shapeFunMBTET3AtOneThird
Definition: Tools.hpp:158
static constexpr double diffShapeFunMBTET1y
derivative of tetrahedral shape function
Definition: Tools.hpp:99
static constexpr double shapeFunMBTET1AtOneThird
Definition: Tools.hpp:154
#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:32
static constexpr double diffShapeFunMBTET1x
derivative of tetrahedral shape function
Definition: Tools.hpp:97
#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:72
static constexpr double diffShapeFunMBTET1z
derivative of tetrahedral shape function
Definition: Tools.hpp:101
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:145
static constexpr double shapeFunMBTET1At000
Definition: Tools.hpp:148
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:500
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:341
static constexpr double diffShapeFunMBTRI1x
derivative of triangle shape function
Definition: Tools.hpp:74
#define diffN_MBTET2x
derivative of tetrahedral shape function
Definition: fem_tools.h:46
static double shapeFunMBTET3(const double x, const double y, const double z)
Definition: Tools.hpp:141
base class for all interface classes
static constexpr double diffShapeFunMBTET0y
derivative of tetrahedral shape function
Definition: Tools.hpp:93
static constexpr std::array< double, 12 > diffShapeFunMBTET
Definition: Tools.hpp:116
double getTriArea(const EntityHandle tri) const
Get triangle area.
Definition: Tools.cpp:279
#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:93
Core (interface) class.
Definition: Core.hpp:50
static constexpr double diffShapeFunMBTET2y
derivative of tetrahedral shape function
Definition: Tools.hpp:105
static constexpr double diffShapeFunMBTET2x
derivative of tetrahedral shape function
Definition: Tools.hpp:103
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:507
#define diffN_MBTET1x
derivative of tetrahedral shape function
Definition: fem_tools.h:43
keeps basic data about problemThis is low level structure with information about problem,...
static double shapeFunMBTET0(const double x, const double y, const double z)
Definition: Tools.hpp:126
static constexpr double shapeFunMBTET0AtOneThird
Definition: Tools.hpp:152
RowColData
RowColData.
Definition: definitions.h:185
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:170
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:91
static constexpr double diffShapeFunMBTET2z
derivative of tetrahedral shape function
Definition: Tools.hpp:107
static constexpr double diffShapeFunMBTET3z
derivative of tetrahedral shape function
Definition: Tools.hpp:113
static constexpr double shapeFunMBTET3At000
Definition: Tools.hpp:150
static constexpr std::array< double, 4 > shapeFunMBTETAt000
Array of shape function at zero local point on reference element.
Definition: Tools.hpp:197
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:188
static double tetVolume(const double *coords)
Calculate volume of tetrahedron.
Definition: Tools.cpp:47
#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
static constexpr double shapeFunMBTET2AtOneThird
Definition: Tools.hpp:156
#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:66
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
static constexpr double shapeFunMBTET0At000
Definition: Tools.hpp:147
static constexpr double diffShapeFunMBTRI1y
derivative of triangle shape function
Definition: Tools.hpp:76
#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:208
static MoFEMErrorCode getTriNormal(const double *coords, double *normal)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:254
#define N_MBTET3(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:39
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: Tools.cpp:20
static const double edge_coords[6][6]
static const MOFEMuuid IDD_MOFEMTools
Definition: Tools.hpp:23
static constexpr double diffShapeFunMBTRI0x
derivative of triangle shape function
Definition: Tools.hpp:70
static constexpr double diffShapeFunMBTET3x
derivative of tetrahedral shape function
Definition: Tools.hpp:109
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition: Tools.hpp:83
#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:318
static constexpr double diffShapeFunMBTET0z
derivative of tetrahedral shape function
Definition: Tools.hpp:95
static constexpr double shapeFunMBTET2At000
Definition: Tools.hpp:149
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:176
#define N_MBTET2(x, y, z)
tetrahedral shape function
Definition: fem_tools.h:38
#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 getEdgeLength(const double *edge_coords)
Get edge length.
Definition: Tools.cpp:287
static constexpr double diffShapeFunMBTRI2x
derivative of triangle shape function
Definition: Tools.hpp:78
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:111
#define diffN_MBTET1y
derivative of tetrahedral shape function
Definition: fem_tools.h:44
static double shapeFunMBTET1(const double x, const double y, const double z)
Definition: Tools.hpp:131
#define diffN_MBTET2z
derivative of tetrahedral shape function
Definition: fem_tools.h:48
#define diffN_MBTET0x
derivative of tetrahedral shape function
Definition: fem_tools.h:40
MoFEM::Core & cOre
Definition: Tools.hpp:35
static constexpr std::array< double, 4 > shapeFunMBTETAtOneThird
Array of shape function at center on reference element.
Definition: Tools.hpp:205
SEGMENT_MIN_DISTANCE
Definition: Tools.hpp:326