v0.15.0
Loading...
Searching...
No Matches
PostProcBrokenMeshInMoabBase.cpp
Go to the documentation of this file.
1/**
2 * @file PostProc.cpp
3 * @brief Post processing elements and operators
4 *
5 * @copyright Copyright (c) 2022
6 *
7 */
8
9#include <MoFEM.hpp>
10
11namespace MoFEM {
12
13PostProcGenerateRefMeshBase::PostProcGenerateRefMeshBase()
14 : hoNodes(PETSC_TRUE), defMaxLevel(0), countEle(0), countVertEle(0),
15 nbVertices(0) {}
16
17MoFEMErrorCode PostProcGenerateRefMeshBase::getOptions(std::string prefix) {
18 MoFEMFunctionBegin;
19
20 std::string opt1 = prefix.size() ? "-" + prefix + "_max_post_proc_ref_level"
21 : "-max_post_proc_ref_level";
22 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, opt1.c_str(), &defMaxLevel, PETSC_NULLPTR);
23
24 std::string opt2 = prefix.size() ? "-" + prefix + "_max_post_ho_nodes"
25 : "-max_post_ho_nodes";
26 CHKERR PetscOptionsGetBool(PETSC_NULLPTR, opt2.c_str(), &hoNodes, PETSC_NULLPTR);
27
28 if (defMaxLevel < 0)
29 SETERRQ(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
30 "Wrong parameter -max_post_proc_ref_level "
31 "should be positive number");
32
33 MoFEMFunctionReturn(0);
34};
35
36MoFEMErrorCode PostProcGenerateRefMesh<MBTET>::generateReferenceElementMesh() {
37 MoFEMFunctionBegin;
38
39 const int max_level = defMaxLevel;
40
41 moab::Core core_ref;
42 moab::Interface &moab_ref = core_ref;
43
44 auto create_reference_element = [&moab_ref]() {
45 MoFEMFunctionBegin;
46 constexpr double base_coords[] = {0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1};
47 EntityHandle nodes[4];
48 for (int nn = 0; nn < 4; nn++) {
49 CHKERR
50 moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
51 }
52 EntityHandle tet;
53 CHKERR moab_ref.create_element(MBTET, nodes, 4, tet);
54 MoFEMFunctionReturn(0);
55 };
56
57 MoFEM::CoreTmp<-1> m_core_ref(moab_ref, PETSC_COMM_SELF, -2);
58 MoFEM::Interface &m_field_ref = m_core_ref;
59
60 auto refine_ref_tetrahedron = [this, &m_field_ref, max_level]() {
61 MoFEMFunctionBegin;
62 // seed ref mofem database by setting bit ref level to reference
63 // tetrahedron
64 CHKERR
65 m_field_ref.getInterface<BitRefManager>()->setBitRefLevelByDim(
66 0, 3, BitRefLevel().set(0));
67 for (int ll = 0; ll != max_level; ++ll) {
68 MOFEM_TAG_AND_LOG_C("WORLD", Sev::noisy, "PostProc", "Refine Level %d",
69 ll);
70 Range edges;
71 CHKERR m_field_ref.getInterface<BitRefManager>()
72 ->getEntitiesByTypeAndRefLevel(BitRefLevel().set(ll),
73 BitRefLevel().set(), MBEDGE, edges);
74 Range tets;
75 CHKERR m_field_ref.getInterface<BitRefManager>()
76 ->getEntitiesByTypeAndRefLevel(BitRefLevel().set(ll),
77 BitRefLevel(ll).set(), MBTET, tets);
78 // refine mesh
79 MeshRefinement *m_ref;
80 CHKERR m_field_ref.getInterface(m_ref);
81 CHKERR m_ref->addVerticesInTheMiddleOfEdges(edges,
82 BitRefLevel().set(ll + 1));
83 CHKERR m_ref->refineTets(tets, BitRefLevel().set(ll + 1));
84 }
85 MoFEMFunctionReturn(0);
86 };
87
88 auto get_ref_gauss_pts_and_shape_functions = [this, max_level, &moab_ref,
89 &m_field_ref]() {
90 MoFEMFunctionBegin;
91 for (int ll = 0; ll != max_level + 1; ++ll) {
92 Range tets;
93 CHKERR
94 m_field_ref.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
95 BitRefLevel().set(ll), BitRefLevel().set(ll), MBTET, tets);
96 if (hoNodes) {
97 EntityHandle meshset;
98 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
99 CHKERR moab_ref.add_entities(meshset, tets);
100 CHKERR moab_ref.convert_entities(meshset, true, false, false);
101 CHKERR moab_ref.delete_entities(&meshset, 1);
102 }
103 Range elem_nodes;
104 CHKERR moab_ref.get_connectivity(tets, elem_nodes, false);
105
106 auto &gauss_pts = levelGaussPtsOnRefMesh[ll];
107 gauss_pts.resize(elem_nodes.size(), 4, false);
108 std::map<EntityHandle, int> little_map;
109 Range::iterator nit = elem_nodes.begin();
110 for (int gg = 0; nit != elem_nodes.end(); nit++, gg++) {
111 CHKERR moab_ref.get_coords(&*nit, 1, &gauss_pts(gg, 0));
112 little_map[*nit] = gg;
113 }
114 gauss_pts = trans(gauss_pts);
115
116 auto &ref_tets = levelRef[ll];
117 Range::iterator tit = tets.begin();
118 for (int tt = 0; tit != tets.end(); ++tit, ++tt) {
119 const EntityHandle *conn;
120 int num_nodes;
121 CHKERR moab_ref.get_connectivity(*tit, conn, num_nodes, false);
122 if (tt == 0) {
123 ref_tets.resize(tets.size(), num_nodes);
124 }
125 for (int nn = 0; nn != num_nodes; ++nn) {
126 ref_tets(tt, nn) = little_map[conn[nn]];
127 }
128 }
129
130 auto &shape_functions = levelShapeFunctions[ll];
131 shape_functions.resize(elem_nodes.size(), 4);
132 CHKERR ShapeMBTET(&*shape_functions.data().begin(), &gauss_pts(0, 0),
133 &gauss_pts(1, 0), &gauss_pts(2, 0), elem_nodes.size());
134 }
135 MoFEMFunctionReturn(0);
136 };
137
138 levelRef.resize(max_level + 1);
139 levelGaussPtsOnRefMesh.resize(max_level + 1);
140 levelShapeFunctions.resize(max_level + 1);
141
142 CHKERR create_reference_element();
143 CHKERR refine_ref_tetrahedron();
144 CHKERR get_ref_gauss_pts_and_shape_functions();
145
146 MoFEMFunctionReturn(0);
147}
148
149MoFEMErrorCode PostProcGenerateRefMesh<MBHEX>::generateReferenceElementMesh() {
150 MoFEMFunctionBegin;
151
152#ifndef NDEBUG
153 if (defMaxLevel > 0)
154 MOFEM_LOG("WORLD", Sev::warning)
155 << "Refinement for hexes is not implemented";
156#endif
157
158 moab::Core core_ref;
159 moab::Interface &moab_ref = core_ref;
160
161 auto create_reference_element = [&moab_ref]() {
162 MoFEMFunctionBegin;
163 constexpr double base_coords[] = {0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0,
164 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1};
165 EntityHandle nodes[8];
166 for (int nn = 0; nn < 8; nn++)
167 CHKERR moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
168 EntityHandle hex;
169 CHKERR moab_ref.create_element(MBHEX, nodes, 8, hex);
170 MoFEMFunctionReturn(0);
171 };
172
173 auto add_ho_nodes = [&]() {
174 MoFEMFunctionBegin;
175 Range hexes;
176 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes, true);
177 EntityHandle meshset;
178 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
179 CHKERR moab_ref.add_entities(meshset, hexes);
180 CHKERR moab_ref.convert_entities(meshset, true, true, true);
181 CHKERR moab_ref.delete_entities(&meshset, 1);
182 MoFEMFunctionReturn(0);
183 };
184
185 auto set_gauss_pts = [&](std::map<EntityHandle, int> &little_map) {
186 MoFEMFunctionBegin;
187 Range hexes;
188 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes, true);
189 Range hexes_nodes;
190 CHKERR moab_ref.get_connectivity(hexes, hexes_nodes, false);
191 auto &gauss_pts = levelGaussPtsOnRefMesh[0];
192 gauss_pts.resize(hexes_nodes.size(), 4, false);
193 size_t gg = 0;
194 for (auto node : hexes_nodes) {
195 CHKERR moab_ref.get_coords(&node, 1, &gauss_pts(gg, 0));
196 little_map[node] = gg;
197 ++gg;
198 }
199 gauss_pts = trans(gauss_pts);
200 MoFEMFunctionReturn(0);
201 };
202
203 auto set_ref_hexes = [&](std::map<EntityHandle, int> &little_map) {
204 MoFEMFunctionBegin;
205 Range hexes;
206 CHKERR moab_ref.get_entities_by_type(0, MBHEX, hexes, true);
207 size_t hh = 0;
208 auto &ref_hexes = levelRef[0];
209 for (auto hex : hexes) {
210 const EntityHandle *conn;
211 int num_nodes;
212 CHKERR moab_ref.get_connectivity(hex, conn, num_nodes, false);
213 if (ref_hexes.size2() != num_nodes) {
214 ref_hexes.resize(hexes.size(), num_nodes);
215 }
216 for (int nn = 0; nn != num_nodes; ++nn) {
217 ref_hexes(hh, nn) = little_map[conn[nn]];
218 }
219 ++hh;
220 }
221 MoFEMFunctionReturn(0);
222 };
223
224 auto set_shape_functions = [&]() {
225 MoFEMFunctionBegin;
226 auto &gauss_pts = levelGaussPtsOnRefMesh[0];
227 auto &shape_functions = levelShapeFunctions[0];
228 const auto nb_gauss_pts = gauss_pts.size2();
229 shape_functions.resize(nb_gauss_pts, 8);
230 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
231 const double ksi = gauss_pts(0, gg);
232 const double zeta = gauss_pts(1, gg);
233 const double eta = gauss_pts(2, gg);
234 shape_functions(gg, 0) = N_MBHEX0(ksi, zeta, eta);
235 shape_functions(gg, 1) = N_MBHEX1(ksi, zeta, eta);
236 shape_functions(gg, 2) = N_MBHEX2(ksi, zeta, eta);
237 shape_functions(gg, 3) = N_MBHEX3(ksi, zeta, eta);
238 shape_functions(gg, 4) = N_MBHEX4(ksi, zeta, eta);
239 shape_functions(gg, 5) = N_MBHEX5(ksi, zeta, eta);
240 shape_functions(gg, 6) = N_MBHEX6(ksi, zeta, eta);
241 shape_functions(gg, 7) = N_MBHEX7(ksi, zeta, eta);
242 }
243 MoFEMFunctionReturn(0);
244 };
245
246 levelRef.resize(1);
247 levelGaussPtsOnRefMesh.resize(1);
248 levelShapeFunctions.resize(1);
249
250 CHKERR create_reference_element();
251 if (hoNodes)
252 CHKERR add_ho_nodes();
253 std::map<EntityHandle, int> little_map;
254 CHKERR set_gauss_pts(little_map);
255 CHKERR set_ref_hexes(little_map);
256 CHKERR set_shape_functions();
257
258 MoFEMFunctionReturn(0);
259}
260
261MoFEMErrorCode PostProcGenerateRefMesh<MBTRI>::generateReferenceElementMesh() {
262 MoFEMFunctionBegin;
263
264 const int max_level = defMaxLevel;
265
266 moab::Core core_ref;
267 moab::Interface &moab_ref = core_ref;
268
269 auto create_reference_element = [&moab_ref]() {
270 MoFEMFunctionBegin;
271 constexpr double base_coords[] = {
272
273 0, 0,
274 0,
275
276 1, 0,
277 0,
278
279 0, 1,
280 0
281
282 };
283 EntityHandle nodes[3];
284 for (int nn = 0; nn != 3; ++nn)
285 CHKERR moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
286 EntityHandle tri;
287 CHKERR moab_ref.create_element(MBTRI, nodes, 3, tri);
288
289 Range edges;
290 CHKERR moab_ref.get_adjacencies(&tri, 1, 1, true, edges,
291 moab::Interface::UNION);
292
293 MoFEMFunctionReturn(0);
294 };
295
296 CHKERR create_reference_element();
297
298 MoFEM::CoreTmp<-1> m_core_ref(moab_ref, PETSC_COMM_SELF, -2);
299 MoFEM::Interface &m_field_ref = m_core_ref;
300
301 auto refine_ref_triangles = [this, &m_field_ref, max_level]() {
302 MoFEMFunctionBegin;
303 // seed ref mofem database by setting bit ref level to reference
304 // tetrahedron
305 CHKERR
306 m_field_ref.getInterface<BitRefManager>()->setBitRefLevelByDim(
307 0, 2, BitRefLevel().set(0));
308
309 for (int ll = 0; ll != max_level; ++ll) {
310 MOFEM_TAG_AND_LOG_C("WORLD", Sev::noisy, "PostProc", "Refine Level %d",
311 ll);
312 Range edges;
313 CHKERR m_field_ref.getInterface<BitRefManager>()
314 ->getEntitiesByTypeAndRefLevel(BitRefLevel().set(ll),
315 BitRefLevel().set(), MBEDGE, edges);
316 Range tris;
317 CHKERR m_field_ref.getInterface<BitRefManager>()
318 ->getEntitiesByTypeAndRefLevel(BitRefLevel().set(ll),
319 BitRefLevel(ll).set(), MBTRI, tris);
320 // refine mesh
321 auto m_ref = m_field_ref.getInterface<MeshRefinement>();
322 CHKERR m_ref->addVerticesInTheMiddleOfEdges(edges,
323 BitRefLevel().set(ll + 1));
324 CHKERR m_ref->refineTris(tris, BitRefLevel().set(ll + 1));
325 }
326 MoFEMFunctionReturn(0);
327 };
328
329 // auto set_gauss_pts = [&](std::map<EntityHandle, int> &little_map) {
330 // MoFEMFunctionBegin;
331 // Range faces;
332 // CHKERR moab_ref.get_entities_by_type(0, MBTRI, faces, true);
333 // Range faces_nodes;
334 // CHKERR moab_ref.get_connectivity(faces, faces_nodes, false);
335 // auto &gauss_pts = levelGaussPtsOnRefMesh[0];
336 // gauss_pts.resize(faces_nodes.size(), 4, false);
337 // size_t gg = 0;
338 // for (auto node : faces_nodes) {
339 // CHKERR moab_ref.get_coords(&node, 1, &gauss_pts(gg, 0));
340 // little_map[node] = gg;
341 // ++gg;
342 // }
343 // gauss_pts = trans(gauss_pts);
344 // MoFEMFunctionReturn(0);
345 // };
346
347 auto get_ref_gauss_pts_and_shape_functions = [this, max_level, &moab_ref,
348 &m_field_ref]() {
349 MoFEMFunctionBegin;
350 for (int ll = 0; ll != max_level + 1; ++ll) {
351
352 Range tris;
353 CHKERR
354 m_field_ref.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
355 BitRefLevel().set(ll), BitRefLevel().set(ll), MBTRI, tris);
356
357 if (hoNodes) {
358 EntityHandle meshset;
359 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
360 CHKERR moab_ref.add_entities(meshset, tris);
361 CHKERR moab_ref.convert_entities(meshset, true, false, false);
362 CHKERR moab_ref.delete_entities(&meshset, 1);
363 }
364
365 Range elem_nodes;
366 CHKERR moab_ref.get_connectivity(tris, elem_nodes, false);
367
368 auto &gauss_pts = levelGaussPtsOnRefMesh[ll];
369 gauss_pts.resize(elem_nodes.size(), 3, false);
370 std::map<EntityHandle, int> little_map;
371 Range::iterator nit = elem_nodes.begin();
372 for (int gg = 0; nit != elem_nodes.end(); nit++, gg++) {
373 CHKERR moab_ref.get_coords(&*nit, 1, &gauss_pts(gg, 0));
374 little_map[*nit] = gg;
375 }
376 gauss_pts = trans(gauss_pts);
377
378 auto &ref_tris = levelRef[ll];
379 Range::iterator tit = tris.begin();
380 for (int tt = 0; tit != tris.end(); ++tit, ++tt) {
381 const EntityHandle *conn;
382 int num_nodes;
383 CHKERR moab_ref.get_connectivity(*tit, conn, num_nodes, false);
384 if (tt == 0) {
385 ref_tris.resize(tris.size(), num_nodes);
386 }
387 for (int nn = 0; nn != num_nodes; ++nn) {
388 ref_tris(tt, nn) = little_map[conn[nn]];
389 }
390 }
391
392 auto &shape_functions = levelShapeFunctions[ll];
393 shape_functions.resize(elem_nodes.size(), 3);
394 CHKERR ShapeMBTRI(&*shape_functions.data().begin(), &gauss_pts(0, 0),
395 &gauss_pts(1, 0), elem_nodes.size());
396 }
397 MoFEMFunctionReturn(0);
398 };
399
400 levelRef.resize(max_level + 1);
401 levelGaussPtsOnRefMesh.resize(max_level + 1);
402 levelShapeFunctions.resize(max_level + 1);
403
404 CHKERR refine_ref_triangles();
405 CHKERR get_ref_gauss_pts_and_shape_functions();
406
407 MoFEMFunctionReturn(0);
408}
409
410MoFEMErrorCode PostProcGenerateRefMesh<MBQUAD>::generateReferenceElementMesh() {
411 MoFEMFunctionBegin;
412
413#ifndef NDEBUG
414 if (defMaxLevel > 0)
415 MOFEM_LOG("WORLD", Sev::warning)
416 << "Refinement for quad is not implemented";
417#endif
418
419 moab::Core core_ref;
420 moab::Interface &moab_ref = core_ref;
421
422 auto create_reference_element = [&moab_ref]() {
423 MoFEMFunctionBegin;
424 constexpr double base_coords[] = {
425
426 0, 0,
427 0,
428
429 1, 0,
430 0,
431
432 1, 1,
433 0,
434
435 0, 1,
436 0
437
438 };
439 EntityHandle nodes[4];
440 for (int nn = 0; nn < 4; nn++)
441 CHKERR moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
442 EntityHandle quad;
443 CHKERR moab_ref.create_element(MBQUAD, nodes, 4, quad);
444 MoFEMFunctionReturn(0);
445 };
446
447 auto add_ho_nodes = [&]() {
448 MoFEMFunctionBegin;
449 Range quads;
450 CHKERR moab_ref.get_entities_by_type(0, MBQUAD, quads, true);
451 EntityHandle meshset;
452 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
453 CHKERR moab_ref.add_entities(meshset, quads);
454 CHKERR moab_ref.convert_entities(meshset, true, true, true);
455 CHKERR moab_ref.delete_entities(&meshset, 1);
456 MoFEMFunctionReturn(0);
457 };
458
459 auto set_gauss_pts = [&](std::map<EntityHandle, int> &little_map) {
460 MoFEMFunctionBegin;
461 Range quads;
462 CHKERR moab_ref.get_entities_by_type(0, MBQUAD, quads, true);
463 Range quads_nodes;
464 CHKERR moab_ref.get_connectivity(quads, quads_nodes, false);
465 auto &gauss_pts = levelGaussPtsOnRefMesh[0];
466 gauss_pts.resize(quads_nodes.size(), 4, false);
467 size_t gg = 0;
468 for (auto node : quads_nodes) {
469 CHKERR moab_ref.get_coords(&node, 1, &gauss_pts(gg, 0));
470 little_map[node] = gg;
471 ++gg;
472 }
473 gauss_pts = trans(gauss_pts);
474 MoFEMFunctionReturn(0);
475 };
476
477 auto set_ref_quads = [&](std::map<EntityHandle, int> &little_map) {
478 MoFEMFunctionBegin;
479 Range quads;
480 CHKERR moab_ref.get_entities_by_type(0, MBQUAD, quads, true);
481 size_t hh = 0;
482 auto &ref_quads = levelRef[0];
483 for (auto quad : quads) {
484 const EntityHandle *conn;
485 int num_nodes;
486 CHKERR moab_ref.get_connectivity(quad, conn, num_nodes, false);
487 if (ref_quads.size2() != num_nodes) {
488 ref_quads.resize(quads.size(), num_nodes);
489 }
490 for (int nn = 0; nn != num_nodes; ++nn) {
491 ref_quads(hh, nn) = little_map[conn[nn]];
492 }
493 ++hh;
494 }
495 MoFEMFunctionReturn(0);
496 };
497
498 auto set_shape_functions = [&]() {
499 MoFEMFunctionBegin;
500 auto &gauss_pts = levelGaussPtsOnRefMesh[0];
501 auto &shape_functions = levelShapeFunctions[0];
502 const auto nb_gauss_pts = gauss_pts.size2();
503 shape_functions.resize(nb_gauss_pts, 4);
504 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
505 const double ksi = gauss_pts(0, gg);
506 const double zeta = gauss_pts(1, gg);
507 shape_functions(gg, 0) = N_MBQUAD0(ksi, zeta);
508 shape_functions(gg, 1) = N_MBQUAD1(ksi, zeta);
509 shape_functions(gg, 2) = N_MBQUAD2(ksi, zeta);
510 shape_functions(gg, 3) = N_MBQUAD3(ksi, zeta);
511 }
512 MoFEMFunctionReturn(0);
513 };
514
515 levelRef.resize(1);
516 levelGaussPtsOnRefMesh.resize(1);
517 levelShapeFunctions.resize(1);
518
519 CHKERR create_reference_element();
520 if (hoNodes)
521 CHKERR add_ho_nodes();
522 std::map<EntityHandle, int> little_map;
523 CHKERR set_gauss_pts(little_map);
524 CHKERR set_ref_quads(little_map);
525 CHKERR set_shape_functions();
526
527 MoFEMFunctionReturn(0);
528}
529
530MoFEMErrorCode PostProcGenerateRefMesh<MBEDGE>::generateReferenceElementMesh() {
531 MoFEMFunctionBegin;
532
533#ifndef NDEBUG
534 if (defMaxLevel > 0)
535 MOFEM_LOG("WORLD", Sev::warning)
536 << "Refinement for edges is not implemented";
537#endif
538
539 auto set_gauss_pts = [&](std::map<EntityHandle, int> &little_map) {
540 MoFEMFunctionBegin;
541
542 int nb_nodes = 2;
543 if (hoNodes)
544 nb_nodes = 3;
545
546 auto &gauss_pts = levelGaussPtsOnRefMesh[0];
547 gauss_pts.resize(2, nb_nodes, false);
548 gauss_pts.clear();
549
550 int nn = 0;
551 for (; nn != 2; ++nn) {
552 gauss_pts(0, nn) = static_cast<double>(nn);
553 little_map[nn] = nn;
554 }
555
556 if (nn < nb_nodes) {
557 gauss_pts(0, nn) = 0.5;
558 little_map[nn] = 2;
559 }
560
561 MoFEMFunctionReturn(0);
562 };
563
564 auto set_ref_edges = [&](std::map<EntityHandle, int> &little_map) {
565 MoFEMFunctionBegin;
566
567 int level = 0;
568 int nb_edges = level + 1;
569
570 int nb_nodes = 2;
571 if (hoNodes)
572 nb_nodes = 3;
573
574 auto &ref_edges = levelRef[level];
575 ref_edges.resize(nb_edges, nb_nodes, false);
576
577 for (int ee = 0; ee != nb_edges; ++ee) {
578 int nn = 0;
579 for (; nn != 2; ++nn) {
580 ref_edges(ee, nn) = nb_nodes * ee + nn;
581 }
582 if (nn < nb_nodes) {
583 ref_edges(ee, nn) = nb_nodes * ee + 2;
584 }
585 }
586
587 MoFEMFunctionReturn(0);
588 };
589
590 auto set_shape_functions = [&]() {
591 MoFEMFunctionBegin;
592 auto &gauss_pts = levelGaussPtsOnRefMesh[0];
593 auto &shape_functions = levelShapeFunctions[0];
594 const auto nb_gauss_pts = gauss_pts.size2();
595 shape_functions.resize(nb_gauss_pts, 2);
596 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
597 const double ksi = gauss_pts(0, gg);
598 shape_functions(gg, 0) = N_MBEDGE0(ksi);
599 shape_functions(gg, 1) = N_MBEDGE1(ksi);
600 }
601 MoFEMFunctionReturn(0);
602 };
603
604 levelRef.resize(1);
605 levelGaussPtsOnRefMesh.resize(1);
606 levelShapeFunctions.resize(1);
607
608 std::map<EntityHandle, int> little_map;
609 CHKERR set_gauss_pts(little_map);
610 CHKERR set_ref_edges(little_map);
611 CHKERR set_shape_functions();
612
613 MoFEMFunctionReturn(0);
614}
615
616} // namespace MoFEM
MOFEM_TAG_AND_LOG_C
#define MOFEM_TAG_AND_LOG_C(channel, severity, tag, format,...)
Tag and log in channel.
Definition LogManager.hpp:384
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition definitions.h:359
MOFEM_INVALID_DATA
@ MOFEM_INVALID_DATA
Definition definitions.h:36
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition definitions.h:432
CHKERR
#define CHKERR
Inline error check.
Definition definitions.h:551
N_MBQUAD3
#define N_MBQUAD3(x, y)
quad shape function
Definition fem_tools.h:60
N_MBHEX7
#define N_MBHEX7(x, y, z)
Definition fem_tools.h:78
N_MBHEX3
#define N_MBHEX3(x, y, z)
Definition fem_tools.h:74
N_MBHEX5
#define N_MBHEX5(x, y, z)
Definition fem_tools.h:76
N_MBEDGE0
#define N_MBEDGE0(x)
edge shape function
Definition fem_tools.h:105
N_MBHEX4
#define N_MBHEX4(x, y, z)
Definition fem_tools.h:75
ShapeMBTET
PetscErrorCode ShapeMBTET(double *N, const double *G_X, const double *G_Y, const double *G_Z, int DIM)
calculate shape functions
Definition fem_tools.c:306
N_MBHEX0
#define N_MBHEX0(x, y, z)
Definition fem_tools.h:71
N_MBHEX6
#define N_MBHEX6(x, y, z)
Definition fem_tools.h:77
N_MBHEX2
#define N_MBHEX2(x, y, z)
Definition fem_tools.h:73
N_MBQUAD0
#define N_MBQUAD0(x, y)
quad shape function
Definition fem_tools.h:57
N_MBHEX1
#define N_MBHEX1(x, y, z)
Definition fem_tools.h:72
N_MBQUAD2
#define N_MBQUAD2(x, y)
quad shape function
Definition fem_tools.h:59
N_MBQUAD1
#define N_MBQUAD1(x, y)
quad shape function
Definition fem_tools.h:58
N_MBEDGE1
#define N_MBEDGE1(x)
edge shape function
Definition fem_tools.h:106
ShapeMBTRI
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition fem_tools.c:182
eta
double eta
Definition free_surface.cpp:196
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition LogManager.hpp:316
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition Exceptions.hpp:56
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition Types.hpp:40
MoFEM
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
MoFEM::PetscOptionsGetInt
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
Definition DeprecatedPetsc.hpp:142
MoFEM::PetscOptionsGetBool
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
Definition DeprecatedPetsc.hpp:182
zeta
double zeta
Viscous hardening.
Definition plastic.cpp:130
MoFEM::BitRefManager
Managing BitRefLevels.
Definition BitRefManager.hpp:21
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition DeprecatedCoreInterface.hpp:16
MoFEM::MeshRefinement
Mesh refinement interface.
Definition MeshRefinement.hpp:26
MoFEM::PostProcGenerateRefMeshBase::getOptions
virtual MoFEMErrorCode getOptions(std::string prefix)
Definition PostProcBrokenMeshInMoabBase.cpp:17
MoFEM::PostProcGenerateRefMesh
Element for postprocessing. Uses MoAB to generate post-processing mesh.
Definition PostProcBrokenMeshInMoabBase.hpp:62
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition UnknownInterface.hpp:93
Generated by Doxygen 1.9.8 and hosted at University of Glasgow