149 {
151
154 const int nb_gauss_pts = pts.size2();
155
157 (unsigned int)nb_gauss_pts)
159 "Base functions or nodes has wrong number of integration points "
160 "for base %s",
163
165 auto &ptr = data.getBBAlphaIndicesSharedPtr(
field_name);
166 if (!ptr)
168 return *ptr;
169 };
170
173 if (!ptr)
175 return *ptr;
176 };
177
179 auto &ptr = data.getBBDiffNSharedPtr(
field_name);
180 if (!ptr)
182 return *ptr;
183 };
184
185 auto get_alpha_by_name_ptr =
186 [](auto &data,
187 const std::string &
field_name) -> boost::shared_ptr<MatrixInt> & {
188 return data.getBBAlphaIndicesSharedPtr(
field_name);
189 };
190
191 auto get_base_by_name_ptr =
192 [](auto &data,
193 const std::string &
field_name) -> boost::shared_ptr<MatrixDouble> & {
195 };
196
197 auto get_diff_base_by_name_ptr =
198 [](auto &data,
199 const std::string &
field_name) -> boost::shared_ptr<MatrixDouble> & {
201 };
202
203 auto get_alpha_by_order_ptr =
204 [](
auto &data,
const size_t o) -> boost::shared_ptr<MatrixInt> & {
205 return data.getBBAlphaIndicesByOrderSharedPtr(
o);
206 };
207
208 auto get_base_by_order_ptr =
209 [](
auto &data,
const size_t o) -> boost::shared_ptr<MatrixDouble> & {
210 return data.getBBNByOrderSharedPtr(
o);
211 };
212
213 auto get_diff_base_by_order_ptr =
214 [](
auto &data,
const size_t o) -> boost::shared_ptr<MatrixDouble> & {
215 return data.getBBDiffNByOrderSharedPtr(
o);
216 };
217
219 auto &vertex_alpha = get_alpha(vert_ent_data);
220 vertex_alpha.resize(4, 4, false);
221 vertex_alpha.clear();
222 for (
int n = 0;
n != 4; ++
n)
224
225 auto &vert_get_n = get_base(vert_ent_data);
226 auto &vert_get_diff_n = get_diff_base(vert_ent_data);
227 vert_get_n.resize(nb_gauss_pts, 4, false);
228 vert_get_diff_n.resize(nb_gauss_pts, 12, false);
230 1,
lambda.size1(), vertex_alpha.size1(), &vertex_alpha(0, 0),
232 &vert_get_diff_n(0, 0));
233 for (
int n = 0;
n != 4; ++
n) {
234 const double f = boost::math::factorial<double>(
236 for (
int g = 0;
g != nb_gauss_pts; ++
g) {
237 vert_get_n(
g,
n) *=
f;
238 for (
int d = 0;
d != 3; ++
d)
239 vert_get_diff_n(
g, 3 *
n + d) *=
f;
240 }
241 }
242
243
247 "Wrong size of ent data");
248
249 constexpr int edges_nodes[6][2] = {{0, 1}, {1, 2}, {2, 0},
250 {0, 3}, {1, 3}, {2, 3}};
251 for (int ee = 0; ee != 6; ++ee) {
253 const int sense = ent_data.getSense();
254 if (sense == 0)
256 "Sense of the edge unknown");
257 const int order = ent_data.getOrder();
259
260 if (nb_dofs) {
261 if (get_alpha_by_order_ptr(ent_data,
order)) {
263 get_alpha_by_order_ptr(ent_data,
order);
265 get_base_by_order_ptr(ent_data,
order);
266 get_diff_base_by_name_ptr(ent_data,
field_name) =
267 get_diff_base_by_order_ptr(ent_data,
order);
268 } else {
269 auto &get_n = get_base(ent_data);
270 auto &get_diff_n = get_diff_base(ent_data);
271 get_n.resize(nb_gauss_pts, nb_dofs, false);
272 get_diff_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
273
275 edge_alpha.resize(nb_dofs, 4, false);
277 &edge_alpha(0, 0));
278 if (sense == -1) {
279 for (
int i = 0;
i != edge_alpha.size1(); ++
i) {
280 int a = edge_alpha(
i, edges_nodes[ee][0]);
281 edge_alpha(
i, edges_nodes[ee][0]) =
282 edge_alpha(
i, edges_nodes[ee][1]);
283 edge_alpha(
i, edges_nodes[ee][1]) =
a;
284 }
285 }
287 order,
lambda.size1(), edge_alpha.size1(), &edge_alpha(0, 0),
289 &get_diff_n(0, 0));
290
291 get_alpha_by_order_ptr(ent_data,
order) =
293 get_base_by_order_ptr(ent_data,
order) =
295 get_diff_base_by_order_ptr(ent_data,
order) =
296 get_diff_base_by_name_ptr(ent_data,
field_name);
297 }
298 }
299 }
300 } else {
301 for (int ee = 0; ee != 6; ++ee) {
303 ent_data.getBBAlphaIndicesSharedPtr(
field_name).reset();
304 auto &get_n = get_base(ent_data);
305 auto &get_diff_n = get_diff_base(ent_data);
306 get_n.resize(nb_gauss_pts, 0, false);
307 get_diff_n.resize(nb_gauss_pts, 0, false);
308 }
309 }
310
311
315 "Wrong size of ent data");
320
321 for (int ff = 0; ff != 4; ++ff) {
323 const int order = ent_data.getOrder();
325
326 if (nb_dofs) {
327 if (get_alpha_by_order_ptr(ent_data,
order)) {
329 get_alpha_by_order_ptr(ent_data,
order);
331 get_base_by_order_ptr(ent_data,
order);
332 get_diff_base_by_name_ptr(ent_data,
field_name) =
333 get_diff_base_by_order_ptr(ent_data,
order);
334 } else {
335
336 auto &get_n = get_base(ent_data);
337 auto &get_diff_n = get_diff_base(ent_data);
338 get_n.resize(nb_gauss_pts, nb_dofs, false);
339 get_diff_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
340
341 auto &face_alpha = get_alpha(ent_data);
342 face_alpha.resize(nb_dofs, 4, false);
343
345 &face_alpha(0, 0));
346 senseFaceAlpha.resize(face_alpha.size1(), face_alpha.size2(),
false);
348 constexpr int tri_nodes[4][3] = {
349 {0, 1, 3}, {1, 2, 3}, {0, 2, 3}, {0, 1, 2}};
350 for (
int d = 0;
d != nb_dofs; ++
d)
351 for (
int n = 0;
n != 3; ++
n)
353 face_alpha(d, tri_nodes[ff][
n]);
356 order,
lambda.size1(), face_alpha.size1(), &face_alpha(0, 0),
358 &get_diff_n(0, 0));
359
360 get_alpha_by_order_ptr(ent_data,
order) =
362 get_base_by_order_ptr(ent_data,
order) =
364 get_diff_base_by_order_ptr(ent_data,
order) =
365 get_diff_base_by_name_ptr(ent_data,
field_name);
366 }
367 }
368 }
369 } else {
370 for (int ff = 0; ff != 4; ++ff) {
372 ent_data.getBBAlphaIndicesSharedPtr(
field_name).reset();
373 auto &get_n = get_base(ent_data);
374 auto &get_diff_n = get_diff_base(ent_data);
375 get_n.resize(nb_gauss_pts, 0, false);
376 get_diff_n.resize(nb_gauss_pts, 0, false);
377 }
378 }
379
383 "Wrong size ent of ent data");
384
386 const int order = ent_data.getOrder();
388 if (get_alpha_by_order_ptr(ent_data,
order)) {
390 get_alpha_by_order_ptr(ent_data,
order);
392 get_base_by_order_ptr(ent_data,
order);
393 get_diff_base_by_name_ptr(ent_data,
field_name) =
394 get_diff_base_by_order_ptr(ent_data,
order);
395 } else {
396
397 auto &get_n = get_base(ent_data);
398 auto &get_diff_n = get_diff_base(ent_data);
399 get_n.resize(nb_gauss_pts, nb_dofs, false);
400 get_diff_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
401 if (nb_dofs) {
402 auto &tet_alpha = get_alpha(ent_data);
403 tet_alpha.resize(nb_dofs, 4, false);
404
407 order,
lambda.size1(), tet_alpha.size1(), &tet_alpha(0, 0),
409 &get_diff_n(0, 0));
410
411 get_alpha_by_order_ptr(ent_data,
order) =
413 get_base_by_order_ptr(ent_data,
order) =
415 get_diff_base_by_order_ptr(ent_data,
order) =
416 get_diff_base_by_name_ptr(ent_data,
field_name);
417 }
418 }
419 } else {
421 ent_data.getBBAlphaIndicesSharedPtr(
field_name).reset();
422 auto &get_n = get_base(ent_data);
423 auto &get_diff_n = get_diff_base(ent_data);
424 get_n.resize(nb_gauss_pts, 0, false);
425 get_diff_n.resize(nb_gauss_pts, 0, false);
426 }
427
429}
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'i', SPACE_DIM > i
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
UBlasMatrix< double > MatrixDouble
UBlasMatrix< int > MatrixInt
constexpr auto field_name
static MoFEMErrorCode generateIndicesTriTet(const int N[], int *alpha[])
static MoFEMErrorCode baseFunctionsTet(const int N, const int gdim, const int n_alpha, const int *alpha, const double *lambda, const double *grad_lambda, double *base, double *grad_base)
static MoFEMErrorCode generateIndicesEdgeTet(const int N[], int *alpha[])
static MoFEMErrorCode generateIndicesTetTet(const int N, int *alpha)
const std::string fieldName
Public Member Functions inherited from MoFEM::BaseFunction