v0.15.0
Loading...
Searching...
No Matches
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#ifndef __TOOLS_HPP__
10#define __TOOLS_HPP__
11
12namespace MoFEM {
13
14/**
15 * \brief Auxiliary tools
16 * \nosubgrouping
17 * \ingroup mofem_tools
18 */
19struct Tools : public UnknownInterface {
20
21 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
22 UnknownInterface **iface) const;
23
25 Tools(const MoFEM::Core &core) : cOre(const_cast<MoFEM::Core &>(core)) {}
26
27 /** \name Computational */
28
29 /**
30 * \brief Calculate tetrahedron volume length quality
31 * @param coords tet coordinates
32 * @return Volume-length quality
33 */
34 static double volumeLengthQuality(const double *coords);
35
36 /**
37 * @brief Calculate volume of tetrahedron
38 *
39 * @param coords
40 * @return double volume
41 */
42 static double tetVolume(const double *coords);
43
44 /**
45 * \brief calculate minimal quality of tetrahedra in range
46 * @param tets range
47 * @param min_quality minimal quality
48 * @return error code
49 */
51 const Range &tets, double &min_quality, Tag th = nullptr,
52 boost::function<double(double, double)> f =
53 [](double a, double b) -> double { return std::min(a, b); });
54
55 static constexpr double shapeFunMBEDGE0At00 = N_MBEDGE0(0);
56 static constexpr double shapeFunMBEDGE1At00 = N_MBEDGE1(0);
57
58 /**
59 * @brief Array of shape function at zero local point on reference element
60 *
61 */
62 static constexpr std::array<double, 2> shapeFunMBEDGEAt00 = {
64
65 static constexpr double diffN_MBEDGE0x = diffN_MBEDGE0;
66 static constexpr double diffN_MBEDGE1x = diffN_MBEDGE1;
67
68 static constexpr std::array<double, 2> diffShapeFunMBEDGE = {diffN_MBEDGE0x,
70
71 static inline double shapeFunMBEDGE0(const double x);
72
73 static inline double shapeFunMBEDGE1(const double x);
74
75 /**
76 * @brief Calculate shape functions on edge
77 *
78 * \note Template parameter is leading dimension of point coordinate arrays,
79 * such that \f$ksi_{n+1} = ksi[n + LDB]\f$
80 *
81 * @tparam 1
82 * @param shape shape functions
83 * @param ksi pointer to first local coordinates
84 * @param nb number of points
85 * @return MoFEMErrorCode
86 */
87 template <int LDB = 1>
88 static MoFEMErrorCode shapeFunMBEDGE(double *shape, const double *ksi,
89 const int nb);
90
91 static constexpr double diffShapeFunMBTRI0x =
92 diffN_MBTRI0x; ///< derivative of triangle shape function
93 static constexpr double diffShapeFunMBTRI0y =
94 diffN_MBTRI0y; ///< derivative of triangle shape function
95 static constexpr double diffShapeFunMBTRI1x =
96 diffN_MBTRI1x; ///< derivative of triangle shape function
97 static constexpr double diffShapeFunMBTRI1y =
98 diffN_MBTRI1y; ///< derivative of triangle shape function
99 static constexpr double diffShapeFunMBTRI2x =
100 diffN_MBTRI2x; ///< derivative of triangle shape function
101 static constexpr double diffShapeFunMBTRI2y =
102 diffN_MBTRI2y; ///< derivative of triangle shape function
103
111
112 static constexpr double shapeFunMBTRI0At00 = N_MBTRI0(0, 0);
113 static constexpr double shapeFunMBTRI1At00 = N_MBTRI1(0, 0);
114 static constexpr double shapeFunMBTRI2At00 = N_MBTRI2(0, 0);
115
116 /**
117 * @brief Array of shape function at zero local point on reference element
118 *
119 */
120 static constexpr std::array<double, 3> shapeFunMBTRIAt00 = {
122
123 static inline double shapeFunMBTRI0(const double x, const double y);
124
125 static inline double shapeFunMBTRI1(const double x, const double y);
126
127 static inline double shapeFunMBTRI2(const double x, const double y);
128
129 /**
130 * @brief Calculate shape functions on triangle
131 *
132 * \note Template parameter is leading dimension of point coordinate arrays,
133 * such that \f$ksi_{n+1} = ksi[n + LDB]\f$
134 *
135 * @tparam 1
136 * @param shape shape functions
137 * @param ksi pointer to first local coordinates
138 * @param eta pointer to second local coordinates
139 * @param nb number of points
140 * @return MoFEMErrorCode
141 */
142 template <int LDB = 1>
143 static MoFEMErrorCode shapeFunMBTRI(double *shape, const double *ksi,
144 const double *eta, const int nb);
145
146 static constexpr double diffShapeFunMBQUADAtCenter0x =
147 diffN_MBQUAD0x(0.5); ///< derivative of quad shape function
148 static constexpr double diffShapeFunMBQUADAtCenter0y =
149 diffN_MBQUAD0y(0.5); ///< derivative of quad shape function
150 static constexpr double diffShapeFunMBQUADAtCenter1x =
151 diffN_MBQUAD1x(0.5); ///< derivative of quad shape function
152 static constexpr double diffShapeFunMBQUADAtCenter1y =
153 diffN_MBQUAD1y(0.5); ///< derivative of quad shape function
154 static constexpr double diffShapeFunMBQUADAtCenter2x =
155 diffN_MBQUAD2x(0.5); ///< derivative of quad shape function
156 static constexpr double diffShapeFunMBQUADAtCenter2y =
157 diffN_MBQUAD2y(0.5); ///< derivative of quad shape function
158 static constexpr double diffShapeFunMBQUADAtCenter3x =
159 diffN_MBQUAD3x(0.5); ///< derivative of quad shape function
160 static constexpr double diffShapeFunMBQUADAtCenter3y =
161 diffN_MBQUAD3y(0.5); ///< derivative of quad shape function
162
163 static constexpr double diffShapeFunMBHEXAtCenter0x =
164 diffN_MBHEX0x(0.5, 0.5); ///< derivative of HEX shape function
165 static constexpr double diffShapeFunMBHEXAtCenter0y =
166 diffN_MBHEX0y(0.5, 0.5); ///< derivative of HEX shape function
167 static constexpr double diffShapeFunMBHEXAtCenter0z =
168 diffN_MBHEX0z(0.5, 0.5); ///< derivative of HEX shape function
169 static constexpr double diffShapeFunMBHEXAtCenter1x =
170 diffN_MBHEX1x(0.5, 0.5); ///< derivative of HEX shape function
171 static constexpr double diffShapeFunMBHEXAtCenter1y =
172 diffN_MBHEX1y(0.5, 0.5); ///< derivative of HEX shape function
173 static constexpr double diffShapeFunMBHEXAtCenter1z =
174 diffN_MBHEX1z(0.5, 0.5); ///< derivative of HEX shape function
175 static constexpr double diffShapeFunMBHEXAtCenter2x =
176 diffN_MBHEX2x(0.5, 0.5); ///< derivative of HEX shape function
177 static constexpr double diffShapeFunMBHEXAtCenter2y =
178 diffN_MBHEX2y(0.5, 0.5); ///< derivative of HEX shape function
179 static constexpr double diffShapeFunMBHEXAtCenter2z =
180 diffN_MBHEX2z(0.5, 0.5); ///< derivative of HEX shape function
181 static constexpr double diffShapeFunMBHEXAtCenter3x =
182 diffN_MBHEX3x(0.5, 0.5); ///< derivative of HEX shape function
183 static constexpr double diffShapeFunMBHEXAtCenter3y =
184 diffN_MBHEX3y(0.5, 0.5); ///< derivative of quad shape function
185 static constexpr double diffShapeFunMBHEXAtCenter3z =
186 diffN_MBHEX3z(0.5, 0.5); ///< derivative of quad shape function
187 static constexpr double diffShapeFunMBHEXAtCenter4x =
188 diffN_MBHEX4x(0.5, 0.5); ///< derivative of HEX shape function
189 static constexpr double diffShapeFunMBHEXAtCenter4y =
190 diffN_MBHEX4y(0.5, 0.5); ///< derivative of quad shape function
191 static constexpr double diffShapeFunMBHEXAtCenter4z =
192 diffN_MBHEX4z(0.5, 0.5); ///< derivative of quad shape function
193 static constexpr double diffShapeFunMBHEXAtCenter5x =
194 diffN_MBHEX5x(0.5, 0.5); ///< derivative of HEX shape function
195 static constexpr double diffShapeFunMBHEXAtCenter5y =
196 diffN_MBHEX5y(0.5, 0.5); ///< derivative of quad shape function
197 static constexpr double diffShapeFunMBHEXAtCenter5z =
198 diffN_MBHEX5z(0.5, 0.5); ///< derivative of quad shape function
199 static constexpr double diffShapeFunMBHEXAtCenter6x =
200 diffN_MBHEX6x(0.5, 0.5); ///< derivative of HEX shape function
201 static constexpr double diffShapeFunMBHEXAtCenter6y =
202 diffN_MBHEX6y(0.5, 0.5); ///< derivative of quad shape function
203 static constexpr double diffShapeFunMBHEXAtCenter6z =
204 diffN_MBHEX6z(0.5, 0.5); ///< derivative of quad shape function
205 static constexpr double diffShapeFunMBHEXAtCenter7x =
206 diffN_MBHEX7x(0.5, 0.5); ///< derivative of HEX shape function
207 static constexpr double diffShapeFunMBHEXAtCenter7y =
208 diffN_MBHEX7y(0.5, 0.5); ///< derivative of quad shape function
209 static constexpr double diffShapeFunMBHEXAtCenter7z =
210 diffN_MBHEX7z(0.5, 0.5); ///< derivative of quad shape function
211
217
245
246 static constexpr double diffShapeFunMBTET0x =
247 diffN_MBTET0x; ///< derivative of tetrahedral shape function
248 static constexpr double diffShapeFunMBTET0y =
249 diffN_MBTET0y; ///< derivative of tetrahedral shape function
250 static constexpr double diffShapeFunMBTET0z =
251 diffN_MBTET0z; ///< derivative of tetrahedral shape function
252 static constexpr double diffShapeFunMBTET1x =
253 diffN_MBTET1x; ///< derivative of tetrahedral shape function
254 static constexpr double diffShapeFunMBTET1y =
255 diffN_MBTET1y; ///< derivative of tetrahedral shape function
256 static constexpr double diffShapeFunMBTET1z =
257 diffN_MBTET1z; ///< derivative of tetrahedral shape function
258 static constexpr double diffShapeFunMBTET2x =
259 diffN_MBTET2x; ///< derivative of tetrahedral shape function
260 static constexpr double diffShapeFunMBTET2y =
261 diffN_MBTET2y; ///< derivative of tetrahedral shape function
262 static constexpr double diffShapeFunMBTET2z =
263 diffN_MBTET2z; ///< derivative of tetrahedral shape function
264 static constexpr double diffShapeFunMBTET3x =
265 diffN_MBTET3x; ///< derivative of tetrahedral shape function
266 static constexpr double diffShapeFunMBTET3y =
267 diffN_MBTET3y; ///< derivative of tetrahedral shape function
268 static constexpr double diffShapeFunMBTET3z =
269 diffN_MBTET3z; ///< derivative of tetrahedral shape function
270
280
281 static inline double shapeFunMBTET0(const double x, const double y,
282 const double z);
283
284 static inline double shapeFunMBTET1(const double x, const double y,
285 const double z);
286
287 static inline double shapeFunMBTET2(const double x, const double y,
288 const double z);
289
290 static inline double shapeFunMBTET3(const double x, const double y,
291 const double z);
292
293 static constexpr double shapeFunMBTET0At000 = N_MBTET0(0, 0, 0);
294 static constexpr double shapeFunMBTET1At000 = N_MBTET1(0, 0, 0);
295 static constexpr double shapeFunMBTET2At000 = N_MBTET2(0, 0, 0);
296 static constexpr double shapeFunMBTET3At000 = N_MBTET3(0, 0, 0);
297
298 static constexpr double shapeFunMBTET0AtOneThird =
299 N_MBTET0(1. / 3., 1. / 3., 1. / 3.);
300 static constexpr double shapeFunMBTET1AtOneThird =
301 N_MBTET1(1. / 3., 1. / 3., 1. / 3.);
302 static constexpr double shapeFunMBTET2AtOneThird =
303 N_MBTET2(1. / 3., 1. / 3., 1. / 3.);
304 static constexpr double shapeFunMBTET3AtOneThird =
305 N_MBTET3(1. / 3., 1. / 3., 1. / 3.);
306
307 /**
308 * @brief Calculate shape functions on tetrahedron
309 *
310 * \note Template parameter is leading dimension of point coordinate arrays,
311 * such that \f$ksi_{n+1} = ksi[n + LDB]\f$
312 *
313 * @tparam 1
314 * @param shape shape functions
315 * @param ksi pointer to first local coordinates
316 * @param eta pointer to second local coordinates
317 * @param zeta pointer to first third coordinates
318 * @param nb number of points
319 * @return MoFEMErrorCode
320 */
321 template <int LDB = 1>
322 static MoFEMErrorCode shapeFunMBTET(double *shape, const double *ksi,
323 const double *eta, const double *zeta,
324 const double nb);
325
326 /**
327 * @brief Array of shape function at zero local point on reference element
328 *
329 */
333
334 /**
335 * @brief Array of shape function at center on reference element
336 *
337 */
341
342 /**
343 * @brief Get the local coordinates on reference four node tet object
344 *
345 * \code
346 * MatrixDouble elem_coords(4, 3);
347 * // Set nodal coordinates
348 * MatrixDouble global_coords(5, 3);
349 * // Set global coordinates
350 * MatrixDouble local_coords(global_coords.size1(), 3);
351 * CHKERR Tools::getLocalCoordinatesOnReferenceFourNodeTet(
352 * &elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
353 * &local_coords(0, 0))
354 * \endcode
355 *
356 * @param elem_coords Global element node coordinates
357 * @param glob_coords Global coordinates
358 * @param nb_nodes Number of points
359 * @param local_coords Result
360 * @return MoFEMErrorCode
361 */
363 const double *elem_coords, const double *glob_coords, const int nb_nodes,
364 double *local_coords);
365
366 /**
367 * @brief Get the local coordinates on reference three node tri object
368 *
369 * \code
370 * MatrixDouble elem_coords(4, 3);
371 * // Set nodal coordinates
372 * MatrixDouble global_coords(5, 3);
373 * // Set global coordinates
374 * MatrixDouble local_coords(global_coords.size1(), 3);
375 * CHKERR Tools::getLocalCoordinatesOnReferenceFourNodeTet(
376 * &elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
377 * &local_coords(0, 0))
378 * \endcode
379 *
380 * @param elem_coords Global element node coordinates
381 * @param glob_coords Global coordinates
382 * @param nb_nodes Number of points
383 * @param local_coords Result
384 * @param d_elem_coords Derivative of local coordinates
385 * @param d_global_coords
386 * @return MoFEMErrorCode
387 */
389 const double *elem_coords, const double *glob_coords, const int nb_nodes,
390 double *local_coords);
391
392 /** @copydoc getLocalCoordinatesOnReferenceThreeNodeTri */
394 const double *elem_coords, const std::complex<double> *glob_coords,
395 const int nb_nodes, std::complex<double> *local_coords);
396
397/** @deprecated use getLocalCoordinatesOnReferenceThreeNodeTri */
399 const double *elem_coords, const double *glob_coords, const int nb_nodes,
400 double *local_coords) {
401 return getLocalCoordinatesOnReferenceThreeNodeTri(elem_coords, glob_coords,
402 nb_nodes, local_coords);
403 }
404
405 /**
406 * @brief Get the local coordinates on reference four node tet object
407 *
408 * \code
409 * MatrixDouble elem_coords(4, 3);
410 * // Set nodal coordinates
411 * MatrixDouble global_coords(5, 3);
412 * // Set global coordinates
413 * MatrixDouble local_coords(global_coords.size1(), 3);
414 * CHKERR Tools::getLocalCoordinatesOnReferenceEdgeNodeEdge(
415 * &elem_coords(0, 0), &global_coords(0, 0), global_coords.size1(),
416 * &local_coords(0, 0))
417 * \endcode
418 *
419 * @param elem_coords Global element node coordinates
420 * @param glob_coords Global coordinates
421 * @param nb_nodes Number of points
422 * @param local_coords Result
423 * @return MoFEMErrorCode
424 */
426 const double *elem_coords, const double *glob_coords, const int nb_nodes,
427 double *local_coords);
428
429 /**
430 * @brief Get the Tets With Quality
431 *
432 * @param out_tets
433 * @param tets
434 * @param th
435 * @param f
436 * @return MoFEMErrorCode
437 */
439 Range &out_tets, const Range &tets, Tag th = nullptr,
440 boost::function<bool(double)> f = [](double q) -> bool {
441 if (q <= 0)
442 return true;
443 else
444 return false;
445 });
446 /**
447 * @brief Write file with tetrahedral of given quality
448 *
449 * @param file_name
450 * @param file_type
451 * @param options
452 * @param tets
453 * @param th
454 * @param f
455 * @return MoFEMErrorCode
456 */
458 const char *file_name, const char *file_type, const char *options,
459 const Range &tets, Tag th = nullptr,
460 boost::function<bool(double)> f = [](double q) -> bool {
461 if (q <= 0)
462 return true;
463 else
464 return false;
465 });
466
467 /**
468 * @brief Check of point is in tetrahedral
469 *
470 * @param tet_coords
471 * @param global_coord
472 * @param tol
473 * @param result
474 * @return MoFEMErrorCode
475 */
476 static MoFEMErrorCode checkIfPointIsInTet(const double tet_coords[],
477 const double global_coord[],
478 const double tol, bool &result);
479
480 /**
481 * @brief Get the Tri Normal objectGet triangle normal
482 *
483 * @param coords
484 * @param normal
485 * @param d_normal derivative, if pointer is null, derivative is not
486 * calculated
487 *
488 * In d_normal, format is 3rd rank tensor, first index is normal direction,
489 * second is node number, third is coordinate.
490 *
491 * This is relation between tensor and pointer
492 * \code {.cpp}
493 * auto t_d_normal = getFTensor3FromPtr<3, 3, 3>(d_normal);
494 * \endcode
495 *
496 * @return MoFEMErrorCode
497 */
498 static MoFEMErrorCode getTriNormal(const double *coords, double *normal,
499 double *d_normal = nullptr);
500
501 /**
502 * @brief Get triangle normal
503 *
504 * @param tri
505 * @param normal
506 * @return MoFEMErrorCode
507 */
508 MoFEMErrorCode getTriNormal(const EntityHandle tri, double *normal) const;
509
510 /**
511 * @brief Get triangle area
512 *
513 * @param tri
514 * @return double
515 */
516 double getTriArea(const EntityHandle tri) const;
517
518 /**
519 * @brief Get edge length
520 *
521 * @param edge_coords
522 * @return double
523 */
524 static double getEdgeLength(const double *edge_coords);
525
526 /**
527 * @brief Get edge length
528 *
529 * @param edge
530 * @return double
531 */
532 double getEdgeLength(const EntityHandle edge);
533
541
542 /**
543 * @brief Find closet point on the segment from the point
544 *
545 * @param w_ptr segment first vertex coordinate
546 * @param v_ptr segment second vertex coordinate
547 * @param p_ptr coordinate of point
548 * @param t_ptr distance on the segment
549 *
550 * \note If t is outside bounds [ 0,-1 ] point is on the line point
551 * beyond segment.
552 *
553 * \code
554 *
555 * double w[] = {-1, 0, 0};
556 * double v[] = {1, 0, 0};
557 * double p[] = {0, 1, 0};
558 * double t;
559 * CHKERR Toolas::minDistancePointFromOnSegment(w, v, p, &t);
560 * double point_on_segment[3];
561 * for (int i = 0; i != 3; ++i)
562 * point_on_segment[i] = w[i] + t * (v[i] - w[i]);
563 *
564 * \endcode
565 *
566 * @return SEGMENT_MIN_DISTANCE
567 */
569 minDistancePointFromOnSegment(const double *w_ptr, const double *v_ptr,
570 const double *p_ptr,
571 double *const t_ptr = nullptr);
572 /**
573 * @brief Find points on two segments in closest distance
574 *
575 * @param w_ptr
576 * @param v_ptr
577 * @param k_ptr
578 * @param l_ptr
579 * @param tvw_ptr
580 * @param tlk_ptr
581 * @return SEGMENT_MIN_DISTANCE
582 *
583 * \note If tvwk or tlk are outside bound [0,-1], it means that points
584 * are on the lines beyond segments, respectively for segment vw and
585 * lk.
586 *
587 */
589 minDistanceFromSegments(const double *w_ptr, const double *v_ptr,
590 const double *k_ptr, const double *l_ptr,
591 double *const tvw_ptr = nullptr,
592 double *const tlk_ptr = nullptr);
593
594 /**
595 * @brief Find minimal distance to edges
596 *
597 * \note Finding only edges with have smaller distance than distance
598 * set on the input by min_dist_ptr
599 *
600 * @param v_ptr point coordinates
601 * @param nb nb points
602 * @param edges range of edges
603 * @param min_dist_ptr on return minimal distance, on input starting
604 * distance
605 * @param o_ptr coordinates of the point on edge
606 * @param o_segments closest segments
607 * @return MoFEMErrorCode
608 */
610 findMinDistanceFromTheEdges(const double *v_ptr, const int nb, Range edges,
611 double *min_dist_ptr, double *o_ptr = nullptr,
612 EntityHandle *o_segments = nullptr) const;
613
614 /** \name Debugging */
615
616 /**@{*/
617
618 /** \brief Print all DOFs for which element of vector is not a number
619 *
620 */
622 const RowColData row_or_col, Vec v);
623
624 /**@}*/
625
626 static MoFEMErrorCode
628 const int edge1);
629
630 static MoFEMErrorCode
632 const int edge1, const int edge2);
633
634 /** \name Mesh refinement */
635
636 /**@{*/
637
638 static constexpr std::array<int, 12> uniformTriangleRefineTriangles = {
639
640 0, 3, 5, // 0
641 3, 1, 4, // 1
642 5, 4, 2, // 2
643 5, 3, 4 // 3
644
645 };
646
648 std::tuple<std::vector<double>, std::vector<int>, std::vector<int>>;
649
650 /**
651 * @brief create uniform triangle mesh of refined elements
652 *
653 * @param nb_levels
654 * @return RefineTrianglesReturn
655 */
657
658 /**
659 * @brief generate integration points for refined triangle mesh for last level
660 *
661 * @param pts
662 * @param refined
663 * @return MatrixDouble
664 */
665 static MatrixDouble
667
668 /**
669 * @brief generate integration points for refined triangle mesh for last level
670 *
671 * @param rule Gauss integration rule
672 * @param refined
673 * @return MatrixDouble
674 */
675 static MatrixDouble
677
678 /**@}*/
679
680};
681
682double Tools::shapeFunMBEDGE0(const double x) {
683 return N_MBEDGE0(x);
684}
685
686double Tools::shapeFunMBEDGE1(const double x) {
687 return N_MBEDGE1(x);
688}
689
690template <int LDB>
691MoFEMErrorCode Tools::shapeFunMBEDGE(double *shape, const double *ksi,
692 const int nb) {
694 for (int n = 0; n != nb; ++n) {
695 shape[0] = shapeFunMBEDGE0(*ksi);
696 shape[1] = shapeFunMBEDGE1(*ksi);
697 shape += 2;
698 ksi += LDB;
699 }
701}
702
703double Tools::shapeFunMBTRI0(const double x, const double y) {
704 return N_MBTRI0(x, y);
705}
706
707double Tools::shapeFunMBTRI1(const double x, const double y) {
708 return N_MBTRI1(x, y);
709}
710
711double Tools::shapeFunMBTRI2(const double x, const double y) {
712 return N_MBTRI2(x, y);
713}
714
715template <int LDB>
716MoFEMErrorCode Tools::shapeFunMBTRI(double *shape, const double *ksi,
717 const double *eta, const int nb) {
719 for (int n = 0; n != nb; ++n) {
720 shape[0] = shapeFunMBTRI0(*ksi, *eta);
721 shape[1] = shapeFunMBTRI1(*ksi, *eta);
722 shape[2] = shapeFunMBTRI2(*ksi, *eta);
723 shape += 3;
724 ksi += LDB;
725 eta += LDB;
726 }
728}
729
730double Tools::shapeFunMBTET0(const double x, const double y, const double z) {
731 return N_MBTET0(x, y, z);
732}
733
734double Tools::shapeFunMBTET1(const double x, const double y, const double z) {
735 return N_MBTET1(x, y, z);
736}
737
738double Tools::shapeFunMBTET2(const double x, const double y, const double z) {
739 return N_MBTET2(x, y, z);
740}
741
742double Tools::shapeFunMBTET3(const double x, const double y, const double z) {
743 return N_MBTET3(x, y, z);
744};
745
746template <int LDB>
747MoFEMErrorCode Tools::shapeFunMBTET(double *shape, const double *ksi,
748 const double *eta, const double *zeta,
749 const double nb) {
751 for (int n = 0; n != nb; ++n) {
752 shape[0] = shapeFunMBTET0(*ksi, *eta, *zeta);
753 shape[1] = shapeFunMBTET1(*ksi, *eta, *zeta);
754 shape[2] = shapeFunMBTET2(*ksi, *eta, *zeta);
755 shape[3] = shapeFunMBTET3(*ksi, *eta, *zeta);
756 shape += 4;
757 ksi += LDB;
758 eta += LDB;
759 zeta += LDB;
760 }
762}
763
764
765} // namespace MoFEM
766
767#endif // __TOOLS_HPP__a
768
769/**
770 * \defgroup mofem_tools Tools interface
771 * \brief Interface for tools
772 *
773 * \ingroup mofem
774 */
constexpr double a
RowColData
RowColData.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define DEPRECATED
Definition definitions.h:17
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define diffN_MBTET0z
derivative of tetrahedral shape function
Definition fem_tools.h:34
#define N_MBTET0(x, y, z)
tetrahedral shape function
Definition fem_tools.h:28
#define diffN_MBTRI2x
derivative of triangle shape function
Definition fem_tools.h:53
#define diffN_MBHEX3z(x, y)
Definition fem_tools.h:98
#define diffN_MBHEX7z(x, y)
Definition fem_tools.h:102
#define N_MBTET1(x, y, z)
tetrahedral shape function
Definition fem_tools.h:29
#define diffN_MBQUAD2y(x)
Definition fem_tools.h:66
#define diffN_MBTRI0y
derivative of triangle shape function
Definition fem_tools.h:50
#define diffN_MBTRI1y
derivative of triangle shape function
Definition fem_tools.h:52
#define diffN_MBQUAD1x(y)
Definition fem_tools.h:63
#define diffN_MBHEX7y(x, z)
Definition fem_tools.h:94
#define N_MBTRI1(x, y)
triangle shape function
Definition fem_tools.h:47
#define diffN_MBHEX6x(y, z)
Definition fem_tools.h:85
#define diffN_MBEDGE0
derivative of edge shape function
Definition fem_tools.h:107
#define diffN_MBQUAD0x(y)
Definition fem_tools.h:61
#define diffN_MBTET2z
derivative of tetrahedral shape function
Definition fem_tools.h:40
#define diffN_MBEDGE1
derivative of edge shape function
Definition fem_tools.h:108
#define diffN_MBQUAD1y(x)
Definition fem_tools.h:64
#define diffN_MBQUAD3y(x)
Definition fem_tools.h:68
#define diffN_MBHEX1y(x, z)
Definition fem_tools.h:88
#define diffN_MBHEX3y(x, z)
Definition fem_tools.h:90
#define N_MBEDGE0(x)
edge shape function
Definition fem_tools.h:105
#define diffN_MBTET3y
derivative of tetrahedral shape function
Definition fem_tools.h:42
#define diffN_MBTET0y
derivative of tetrahedral shape function
Definition fem_tools.h:33
#define diffN_MBHEX4y(x, z)
Definition fem_tools.h:91
#define diffN_MBHEX5x(y, z)
Definition fem_tools.h:84
#define diffN_MBTET1z
derivative of tetrahedral shape function
Definition fem_tools.h:37
#define diffN_MBHEX1z(x, y)
Definition fem_tools.h:96
#define diffN_MBTET3z
derivative of tetrahedral shape function
Definition fem_tools.h:43
#define diffN_MBHEX4z(x, y)
Definition fem_tools.h:99
#define diffN_MBHEX6y(x, z)
Definition fem_tools.h:93
#define diffN_MBHEX6z(x, y)
Definition fem_tools.h:101
#define diffN_MBHEX5z(x, y)
Definition fem_tools.h:100
#define diffN_MBQUAD0y(x)
Definition fem_tools.h:62
#define diffN_MBTET1y
derivative of tetrahedral shape function
Definition fem_tools.h:36
#define diffN_MBHEX0x(y, z)
Definition fem_tools.h:79
#define diffN_MBHEX1x(y, z)
Definition fem_tools.h:80
#define diffN_MBTRI0x
derivative of triangle shape function
Definition fem_tools.h:49
#define diffN_MBHEX2z(x, y)
Definition fem_tools.h:97
#define diffN_MBHEX5y(x, z)
Definition fem_tools.h:92
#define diffN_MBHEX4x(y, z)
Definition fem_tools.h:83
#define diffN_MBTET2y
derivative of tetrahedral shape function
Definition fem_tools.h:39
#define diffN_MBTET2x
derivative of tetrahedral shape function
Definition fem_tools.h:38
#define N_MBTRI0(x, y)
triangle shape function
Definition fem_tools.h:46
#define diffN_MBTRI2y
derivative of triangle shape function
Definition fem_tools.h:54
#define diffN_MBHEX3x(y, z)
Definition fem_tools.h:82
#define diffN_MBQUAD3x(y)
Definition fem_tools.h:67
#define diffN_MBQUAD2x(y)
Definition fem_tools.h:65
#define diffN_MBHEX2y(x, z)
Definition fem_tools.h:89
#define N_MBTRI2(x, y)
triangle shape function
Definition fem_tools.h:48
#define diffN_MBHEX0y(x, z)
Definition fem_tools.h:87
#define diffN_MBHEX0z(x, y)
Definition fem_tools.h:95
#define N_MBTET2(x, y, z)
tetrahedral shape function
Definition fem_tools.h:30
#define N_MBEDGE1(x)
edge shape function
Definition fem_tools.h:106
#define diffN_MBTRI1x
derivative of triangle shape function
Definition fem_tools.h:51
#define diffN_MBTET3x
derivative of tetrahedral shape function
Definition fem_tools.h:41
#define diffN_MBHEX2x(y, z)
Definition fem_tools.h:81
#define diffN_MBTET1x
derivative of tetrahedral shape function
Definition fem_tools.h:35
#define diffN_MBTET0x
derivative of tetrahedral shape function
Definition fem_tools.h:32
#define diffN_MBHEX7x(y, z)
Definition fem_tools.h:86
#define N_MBTET3(x, y, z)
tetrahedral shape function
Definition fem_tools.h:31
static const double edge_coords[6][6]
double eta
const double v
phase velocity of light in medium (cm/ns)
const double n
refractive index of diffusive medium
constexpr int nb_levels
Definition level_set.cpp:58
double tol
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
Definition Types.hpp:77
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
double zeta
Viscous hardening.
Definition plastic.cpp:130
Core (interface) class.
Definition Core.hpp:82
keeps basic data about problem
Auxiliary tools.
Definition Tools.hpp:19
static constexpr double shapeFunMBEDGE1At00
Definition Tools.hpp:56
static constexpr double diffShapeFunMBTET0y
derivative of tetrahedral shape function
Definition Tools.hpp:248
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:472
static RefineTrianglesReturn refineTriangle(int nb_levels)
create uniform triangle mesh of refined elements
Definition Tools.cpp:724
static constexpr double diffShapeFunMBHEXAtCenter6y
derivative of quad shape function
Definition Tools.hpp:201
static constexpr double shapeFunMBTET1At000
Definition Tools.hpp:294
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:88
static constexpr double shapeFunMBTET0At000
Definition Tools.hpp:293
static constexpr double diffShapeFunMBHEXAtCenter0y
derivative of HEX shape function
Definition Tools.hpp:165
static double shapeFunMBTET0(const double x, const double y, const double z)
Definition Tools.hpp:730
static MoFEMErrorCode outerProductOfEdgeIntegrationPtsForQuad(MatrixDouble &pts, const int edge0, const int edge1)
Definition Tools.cpp:613
static constexpr double diffShapeFunMBTRI1x
derivative of triangle shape function
Definition Tools.hpp:95
static constexpr double diffShapeFunMBQUADAtCenter0x
derivative of quad shape function
Definition Tools.hpp:146
static constexpr double diffShapeFunMBQUADAtCenter2y
derivative of quad shape function
Definition Tools.hpp:156
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:449
static constexpr double shapeFunMBTRI1At00
Definition Tools.hpp:113
static double shapeFunMBTET3(const double x, const double y, const double z)
Definition Tools.hpp:742
static constexpr double diffShapeFunMBTET0z
derivative of tetrahedral shape function
Definition Tools.hpp:250
std::tuple< std::vector< double >, std::vector< int >, std::vector< int > > RefineTrianglesReturn
Definition Tools.hpp:647
static constexpr double diffShapeFunMBHEXAtCenter1x
derivative of HEX shape function
Definition Tools.hpp:169
static constexpr double shapeFunMBTRI0At00
Definition Tools.hpp:112
static constexpr double diffShapeFunMBTET3y
derivative of tetrahedral shape function
Definition Tools.hpp:266
static constexpr std::array< double, 8 > diffShapeFunMBQUADAtCenter
Definition Tools.hpp:212
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:538
static constexpr double diffShapeFunMBTET1z
derivative of tetrahedral shape function
Definition Tools.hpp:256
static double shapeFunMBTRI1(const double x, const double y)
Definition Tools.hpp:707
static constexpr double diffShapeFunMBQUADAtCenter2x
derivative of quad shape function
Definition Tools.hpp:154
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition Tools.hpp:104
static constexpr double diffShapeFunMBHEXAtCenter3z
derivative of quad shape function
Definition Tools.hpp:185
static constexpr double diffShapeFunMBHEXAtCenter7y
derivative of quad shape function
Definition Tools.hpp:207
static constexpr double diffShapeFunMBHEXAtCenter4y
derivative of quad shape function
Definition Tools.hpp:189
static constexpr double diffShapeFunMBTRI0x
derivative of triangle shape function
Definition Tools.hpp:91
static double shapeFunMBTET1(const double x, const double y, const double z)
Definition Tools.hpp:734
static double shapeFunMBEDGE1(const double x)
Definition Tools.hpp:686
@ SEGMENT_TWO_AND_TWO_ARE_POINT
Definition Tools.hpp:538
@ SEGMENT_ONE_IS_POINT
Definition Tools.hpp:536
@ SEGMENT_TWO_IS_POINT
Definition Tools.hpp:537
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:307
static constexpr std::array< double, 2 > shapeFunMBEDGEAt00
Array of shape function at zero local point on reference element.
Definition Tools.hpp:62
static constexpr double shapeFunMBEDGE0At00
Definition Tools.hpp:55
static MatrixDouble refineTriangleIntegrationPts(MatrixDouble pts, RefineTrianglesReturn refined)
generate integration points for refined triangle mesh for last level
Definition Tools.cpp:791
static constexpr double diffShapeFunMBHEXAtCenter0z
derivative of HEX shape function
Definition Tools.hpp:167
static constexpr double diffShapeFunMBHEXAtCenter1y
derivative of HEX shape function
Definition Tools.hpp:171
static constexpr double diffShapeFunMBTET2z
derivative of tetrahedral shape function
Definition Tools.hpp:262
static constexpr double diffShapeFunMBTET3x
derivative of tetrahedral shape function
Definition Tools.hpp:264
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:49
static constexpr double diffShapeFunMBHEXAtCenter5y
derivative of quad shape function
Definition Tools.hpp:195
static double shapeFunMBTRI0(const double x, const double y)
Definition Tools.hpp:703
static constexpr double shapeFunMBTET1AtOneThird
Definition Tools.hpp:300
static constexpr std::array< int, 12 > uniformTriangleRefineTriangles
Definition Tools.hpp:638
Tools(const MoFEM::Core &core)
Definition Tools.hpp:25
static constexpr double diffShapeFunMBTET0x
derivative of tetrahedral shape function
Definition Tools.hpp:246
static double shapeFunMBTET2(const double x, const double y, const double z)
Definition Tools.hpp:738
static constexpr double diffShapeFunMBHEXAtCenter2y
derivative of HEX shape function
Definition Tools.hpp:177
static constexpr double diffShapeFunMBHEXAtCenter3x
derivative of HEX shape function
Definition Tools.hpp:181
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:269
static MoFEMErrorCode getLocalCoordinatesOnReferenceThreeNodeTri(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Get the local coordinates on reference three node tri object.
Definition Tools.cpp:188
static constexpr double shapeFunMBTET2AtOneThird
Definition Tools.hpp:302
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:244
static constexpr double diffShapeFunMBQUADAtCenter3y
derivative of quad shape function
Definition Tools.hpp:160
static constexpr double diffShapeFunMBQUADAtCenter1y
derivative of quad shape function
Definition Tools.hpp:152
static MoFEMErrorCode outerProductOfEdgeIntegrationPtsForHex(MatrixDouble &pts, const int edge0, const int edge1, const int edge2)
Definition Tools.cpp:662
static MoFEMErrorCode getTriNormal(const double *coords, double *normal, double *d_normal=nullptr)
Get the Tri Normal objectGet triangle normal.
Definition Tools.cpp:353
static constexpr double shapeFunMBTRI2At00
Definition Tools.hpp:114
static constexpr std::array< double, 24 > diffShapeFunMBHEXAtCenter
Definition Tools.hpp:218
static constexpr double shapeFunMBTET3At000
Definition Tools.hpp:296
static constexpr std::array< double, 4 > shapeFunMBTETAt000
Array of shape function at zero local point on reference element.
Definition Tools.hpp:330
static constexpr double diffShapeFunMBHEXAtCenter0x
derivative of HEX shape function
Definition Tools.hpp:163
static constexpr double diffShapeFunMBQUADAtCenter3x
derivative of quad shape function
Definition Tools.hpp:158
static double shapeFunMBEDGE0(const double x)
Definition Tools.hpp:682
static constexpr double diffN_MBEDGE0x
Definition Tools.hpp:65
static constexpr double diffN_MBEDGE1x
Definition Tools.hpp:66
static constexpr double diffShapeFunMBHEXAtCenter3y
derivative of quad shape function
Definition Tools.hpp:183
static constexpr double diffShapeFunMBTRI2x
derivative of triangle shape function
Definition Tools.hpp:99
static constexpr double diffShapeFunMBHEXAtCenter2x
derivative of HEX shape function
Definition Tools.hpp:175
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:287
static MoFEMErrorCode getLocalCoordinatesOnReferenceEdgeNodeEdge(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:203
static constexpr double diffShapeFunMBHEXAtCenter7z
derivative of quad shape function
Definition Tools.hpp:209
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:747
static constexpr std::array< double, 2 > diffShapeFunMBEDGE
Definition Tools.hpp:68
static constexpr double diffShapeFunMBTET1y
derivative of tetrahedral shape function
Definition Tools.hpp:254
static MoFEMErrorCode shapeFunMBEDGE(double *shape, const double *ksi, const int nb)
Calculate shape functions on edge.
Definition Tools.hpp:691
static double tetVolume(const double *coords)
Calculate volume of tetrahedron.
Definition Tools.cpp:30
static constexpr double diffShapeFunMBTRI2y
derivative of triangle shape function
Definition Tools.hpp:101
static double volumeLengthQuality(const double *coords)
Calculate tetrahedron volume length quality.
Definition Tools.cpp:15
double getTriArea(const EntityHandle tri) const
Get triangle area.
Definition Tools.cpp:408
static constexpr std::array< double, 3 > shapeFunMBTRIAt00
Array of shape function at zero local point on reference element.
Definition Tools.hpp:120
static constexpr double diffShapeFunMBHEXAtCenter7x
derivative of HEX shape function
Definition Tools.hpp:205
static constexpr double diffShapeFunMBTRI0y
derivative of triangle shape function
Definition Tools.hpp:93
static constexpr double shapeFunMBTET3AtOneThird
Definition Tools.hpp:304
static constexpr double diffShapeFunMBTET3z
derivative of tetrahedral shape function
Definition Tools.hpp:268
static constexpr double diffShapeFunMBTET1x
derivative of tetrahedral shape function
Definition Tools.hpp:252
static constexpr double diffShapeFunMBHEXAtCenter6x
derivative of HEX shape function
Definition Tools.hpp:199
static constexpr double diffShapeFunMBHEXAtCenter4x
derivative of HEX shape function
Definition Tools.hpp:187
static constexpr double shapeFunMBTET2At000
Definition Tools.hpp:295
static constexpr double diffShapeFunMBHEXAtCenter5z
derivative of quad shape function
Definition Tools.hpp:197
static constexpr double shapeFunMBTET0AtOneThird
Definition Tools.hpp:298
static double shapeFunMBTRI2(const double x, const double y)
Definition Tools.hpp:711
static constexpr double diffShapeFunMBHEXAtCenter2z
derivative of HEX shape function
Definition Tools.hpp:179
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition Tools.cpp:9
static constexpr double diffShapeFunMBHEXAtCenter4z
derivative of quad shape function
Definition Tools.hpp:191
static constexpr std::array< double, 12 > diffShapeFunMBTET
Definition Tools.hpp:271
static constexpr double diffShapeFunMBTET2x
derivative of tetrahedral shape function
Definition Tools.hpp:258
static double getEdgeLength(const double *edge_coords)
Get edge length.
Definition Tools.cpp:418
static DEPRECATED MoFEMErrorCode getLocalCoordinatesOnReferenceTriNodeTri(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Definition Tools.hpp:398
MoFEM::Core & cOre
Definition Tools.hpp:24
static constexpr double diffShapeFunMBQUADAtCenter1x
derivative of quad shape function
Definition Tools.hpp:150
static constexpr double diffShapeFunMBHEXAtCenter6z
derivative of quad shape function
Definition Tools.hpp:203
static constexpr std::array< double, 4 > shapeFunMBTETAtOneThird
Array of shape function at center on reference element.
Definition Tools.hpp:338
static constexpr double diffShapeFunMBHEXAtCenter1z
derivative of HEX shape function
Definition Tools.hpp:173
static MoFEMErrorCode shapeFunMBTRI(double *shape, const double *ksi, const double *eta, const int nb)
Calculate shape functions on triangle.
Definition Tools.hpp:716
static constexpr double diffShapeFunMBTET2y
derivative of tetrahedral shape function
Definition Tools.hpp:260
static constexpr double diffShapeFunMBHEXAtCenter5x
derivative of HEX shape function
Definition Tools.hpp:193
static constexpr double diffShapeFunMBTRI1y
derivative of triangle shape function
Definition Tools.hpp:97
static constexpr double diffShapeFunMBQUADAtCenter0y
derivative of quad shape function
Definition Tools.hpp:148
base class for all interface classes