v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | List of all members
MoFEM::PostProcGenerateRefMesh< MBTRI > Struct Reference

#include <src/post_proc/PostProcBrokenMeshInMoabBase.hpp>

Inheritance diagram for MoFEM::PostProcGenerateRefMesh< MBTRI >:
[legend]
Collaboration diagram for MoFEM::PostProcGenerateRefMesh< MBTRI >:
[legend]

Public Member Functions

MoFEMErrorCode generateReferenceElementMesh ()
 
 PostProcGenerateRefMeshBase ()
 
- Public Member Functions inherited from MoFEM::PostProcGenerateRefMeshBase
 PostProcGenerateRefMeshBase ()
 
virtual ~PostProcGenerateRefMeshBase ()=default
 
virtual MoFEMErrorCode getOptions (std::string prefix)
 
virtual MoFEMErrorCode generateReferenceElementMesh ()=0
 

Additional Inherited Members

- Public Attributes inherited from MoFEM::PostProcGenerateRefMeshBase
std::vector< MatrixDoublelevelShapeFunctions
 
std::vector< MatrixDoublelevelGaussPtsOnRefMesh
 
std::vector< ublas::matrix< int > > levelRef
 
EntityHandle startingVertEleHandle
 
std::vector< double * > verticesOnEleArrays
 
EntityHandle startingEleHandle
 
EntityHandleeleConn
 
int countEle
 
int countVertEle
 
int nbVertices
 
int nbEles
 
PetscBool hoNodes
 
int defMaxLevel
 

Detailed Description

Definition at line 77 of file PostProcBrokenMeshInMoabBase.hpp.

Member Function Documentation

◆ generateReferenceElementMesh()

MoFEMErrorCode MoFEM::PostProcGenerateRefMesh< MBTRI >::generateReferenceElementMesh ( )
virtual

Implements MoFEM::PostProcGenerateRefMeshBase.

Definition at line 259 of file PostProcBrokenMeshInMoabBase.cpp.

259 {
261
262 const int max_level = defMaxLevel;
263
264 moab::Core core_ref;
265 moab::Interface &moab_ref = core_ref;
266
267 auto create_reference_element = [&moab_ref]() {
269 constexpr double base_coords[] = {
270
271 0, 0,
272 0,
273
274 1, 0,
275 0,
276
277 0, 1,
278 0
279
280 };
281 EntityHandle nodes[3];
282 for (int nn = 0; nn != 3; ++nn)
283 CHKERR moab_ref.create_vertex(&base_coords[3 * nn], nodes[nn]);
284 EntityHandle tri;
285 CHKERR moab_ref.create_element(MBTRI, nodes, 3, tri);
286
287 Range edges;
288 CHKERR moab_ref.get_adjacencies(&tri, 1, 1, true, edges,
289 moab::Interface::UNION);
290
292 };
293
294 CHKERR create_reference_element();
295
296 MoFEM::CoreTmp<-1> m_core_ref(moab_ref, PETSC_COMM_SELF, -2);
297 MoFEM::Interface &m_field_ref = m_core_ref;
298
299 auto refine_ref_triangles = [this, &m_field_ref, max_level]() {
301 // seed ref mofem database by setting bit ref level to reference
302 // tetrahedron
303 CHKERR
304 m_field_ref.getInterface<BitRefManager>()->setBitRefLevelByDim(
305 0, 2, BitRefLevel().set(0));
306
307 for (int ll = 0; ll != max_level; ++ll) {
308 MOFEM_TAG_AND_LOG_C("WORLD", Sev::noisy, "PostProc", "Refine Level %d",
309 ll);
310 Range edges;
311 CHKERR m_field_ref.getInterface<BitRefManager>()
312 ->getEntitiesByTypeAndRefLevel(BitRefLevel().set(ll),
313 BitRefLevel().set(), MBEDGE, edges);
314 Range tris;
315 CHKERR m_field_ref.getInterface<BitRefManager>()
316 ->getEntitiesByTypeAndRefLevel(BitRefLevel().set(ll),
317 BitRefLevel(ll).set(), MBTRI, tris);
318 // refine mesh
319 auto m_ref = m_field_ref.getInterface<MeshRefinement>();
320 CHKERR m_ref->addVerticesInTheMiddleOfEdges(edges,
321 BitRefLevel().set(ll + 1));
322 CHKERR m_ref->refineTris(tris, BitRefLevel().set(ll + 1));
323 }
325 };
326
327 // auto set_gauss_pts = [&](std::map<EntityHandle, int> &little_map) {
328 // MoFEMFunctionBegin;
329 // Range faces;
330 // CHKERR moab_ref.get_entities_by_type(0, MBTRI, faces, true);
331 // Range faces_nodes;
332 // CHKERR moab_ref.get_connectivity(faces, faces_nodes, false);
333 // auto &gauss_pts = levelGaussPtsOnRefMesh[0];
334 // gauss_pts.resize(faces_nodes.size(), 4, false);
335 // size_t gg = 0;
336 // for (auto node : faces_nodes) {
337 // CHKERR moab_ref.get_coords(&node, 1, &gauss_pts(gg, 0));
338 // little_map[node] = gg;
339 // ++gg;
340 // }
341 // gauss_pts = trans(gauss_pts);
342 // MoFEMFunctionReturn(0);
343 // };
344
345 auto get_ref_gauss_pts_and_shape_functions = [this, max_level, &moab_ref,
346 &m_field_ref]() {
348 for (int ll = 0; ll != max_level + 1; ++ll) {
349
350 Range tris;
351 CHKERR
352 m_field_ref.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
353 BitRefLevel().set(ll), BitRefLevel().set(ll), MBTRI, tris);
354
355 if (hoNodes) {
356 EntityHandle meshset;
357 CHKERR moab_ref.create_meshset(MESHSET_SET, meshset);
358 CHKERR moab_ref.add_entities(meshset, tris);
359 CHKERR moab_ref.convert_entities(meshset, true, false, false);
360 CHKERR moab_ref.delete_entities(&meshset, 1);
361 }
362
363 Range elem_nodes;
364 CHKERR moab_ref.get_connectivity(tris, elem_nodes, false);
365
366 auto &gauss_pts = levelGaussPtsOnRefMesh[ll];
367 gauss_pts.resize(elem_nodes.size(), 3, false);
368 std::map<EntityHandle, int> little_map;
369 Range::iterator nit = elem_nodes.begin();
370 for (int gg = 0; nit != elem_nodes.end(); nit++, gg++) {
371 CHKERR moab_ref.get_coords(&*nit, 1, &gauss_pts(gg, 0));
372 little_map[*nit] = gg;
373 }
374 gauss_pts = trans(gauss_pts);
375
376 auto &ref_tris = levelRef[ll];
377 Range::iterator tit = tris.begin();
378 for (int tt = 0; tit != tris.end(); ++tit, ++tt) {
379 const EntityHandle *conn;
380 int num_nodes;
381 CHKERR moab_ref.get_connectivity(*tit, conn, num_nodes, false);
382 if (tt == 0) {
383 ref_tris.resize(tris.size(), num_nodes);
384 }
385 for (int nn = 0; nn != num_nodes; ++nn) {
386 ref_tris(tt, nn) = little_map[conn[nn]];
387 }
388 }
389
390 auto &shape_functions = levelShapeFunctions[ll];
391 shape_functions.resize(elem_nodes.size(), 3);
392 CHKERR ShapeMBTRI(&*shape_functions.data().begin(), &gauss_pts(0, 0),
393 &gauss_pts(1, 0), elem_nodes.size());
394 }
396 };
397
398 levelRef.resize(max_level + 1);
399 levelGaussPtsOnRefMesh.resize(max_level + 1);
400 levelShapeFunctions.resize(max_level + 1);
401
402 CHKERR refine_ref_triangles();
403 CHKERR get_ref_gauss_pts_and_shape_functions();
404
406}
#define MOFEM_TAG_AND_LOG_C(channel, severity, tag, format,...)
Tag and log in channel.
Definition: LogManager.hpp:370
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#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
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:182
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
Deprecated interface functions.
std::vector< ublas::matrix< int > > levelRef
PetscBool hoNodes
int defMaxLevel
std::vector< MatrixDouble > levelShapeFunctions
std::vector< MatrixDouble > levelGaussPtsOnRefMesh
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ PostProcGenerateRefMeshBase()

MoFEM::PostProcGenerateRefMeshBase::PostProcGenerateRefMeshBase ( )

The documentation for this struct was generated from the following files: