v0.14.0
Loading...
Searching...
No Matches
EdgeElementForcesAndSourcesCoreOnParent.cpp
Go to the documentation of this file.
1/** \file EdgeElementForcesAndSourcesCoreOnParent.cpp
2
3\brief Implementation of edge element for parent
4
5*/
6
7namespace MoFEM {
8
10 return -1;
11};
12
16
17 auto ref_fe = refinePtrFE;
18 if (ref_fe == nullptr)
19 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
20 "Pointer to face element is not set");
21
22 const auto ref_ent = ref_fe->getFEEntityHandle();
23 const auto ent = getFEEntityHandle();
24
25 const auto &ref_gauss_pts = ref_fe->gaussPts;
26 const auto nb_integration_points = ref_gauss_pts.size2();
27
28 auto set_integration_pts_for_edge = [&]() {
30
31 auto get_coords = [&](const auto ent) {
32 int num_nodes;
33 const EntityHandle *conn;
34 CHKERR mField.get_moab().get_connectivity(ent, conn, num_nodes, true);
35 std::array<double, 6> node_coords;
36 CHKERR mField.get_moab().get_coords(conn, num_nodes, node_coords.data());
37 return node_coords;
38 };
39
40 auto ref_node_coords = get_coords(ref_ent);
41 auto node_coords = get_coords(ent);
42
43 MatrixDouble ref_shapes(nb_integration_points, 2);
44 CHKERR Tools::shapeFunMBEDGE<1>(&ref_shapes(0, 0), &ref_gauss_pts(0, 0),
45 nb_integration_points);
46
47#ifndef NDEBUG
48 if (ref_shapes.size1() * ref_shapes.size2() !=
49 nb_integration_points * ref_node_coords.size() / 3)
50 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
51 "Wrong number of shape functions");
52#endif
53
54 auto get_glob_coords = [&]() {
55 MatrixDouble glob_coords(nb_integration_points, 3);
56
58 &glob_coords(0, 0), &glob_coords(0, 1), &glob_coords(0, 2)};
59 auto t_shape =
61
62 for (auto gg = 0; gg != nb_integration_points; ++gg) {
64
65 auto t_ref_node_coords =
67 &ref_node_coords[0], &ref_node_coords[1], &ref_node_coords[2]};
68
69 t_glob_coords(i) = 0;
70 for (int nn = 0; nn != ref_node_coords.size() / 3; ++nn) {
71 t_glob_coords(i) += t_shape * t_ref_node_coords(i);
72 ++t_ref_node_coords;
73 ++t_shape;
74 }
75
76 ++t_glob_coords;
77 }
78 return glob_coords;
79 };
80
81 auto glob_coords = get_glob_coords();
82 VectorDouble local_coords(nb_integration_points);
83
85 node_coords.data(), &glob_coords(0, 0), nb_integration_points,
86 &local_coords[0]);
87
88 gaussPts.resize(2, nb_integration_points, false);
90 &gaussPts(0, 0), &gaussPts(1, 0)};
91 auto t_local_coords =
93
94 auto scale_quadarture = [&]() {
95 const double scale = Tools::getEdgeLength(ref_node_coords.data()) /
96 Tools::getEdgeLength(node_coords.data());
97 return scale;
98 };
99
100 const auto sq = scale_quadarture();
101
102 for (auto gg = 0; gg != nb_integration_points; ++gg) {
103 t_gauss_pts(0) = t_local_coords;
104 t_gauss_pts(1) = sq * ref_gauss_pts(1, gg);
105 ++t_gauss_pts;
106 ++t_local_coords;
107 }
109 };
110
111 const auto type = numeredEntFiniteElementPtr->getEntType();
112
113 switch (type) {
114 case MBEDGE:
115 CHKERR set_integration_pts_for_edge();
116 break;
117 default:
118 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
119 "Element type not implemented: %d", type);
120 }
121
123}
124
125} // namespace MoFEM
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
constexpr int order
FTensor::Index< 'i', SPACE_DIM > i
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
double scale
Definition: plastic.cpp:170
virtual moab::Interface & get_moab()=0
EntityHandle getFEEntityHandle() const
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MatrixDouble gaussPts
Matrix of integration points.
ForcesAndSourcesCore * refinePtrFE
Element to integrate parent or child.
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:205
static double getEdgeLength(const double *edge_coords)
Get edge length.
Definition: Tools.cpp:417