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

#include <src/post_proc/PostProcBrokenMeshInMoabBase.hpp>

Inheritance diagram for MoFEM::PostProcGenerateRefMesh< MBEDGE >:
[legend]
Collaboration diagram for MoFEM::PostProcGenerateRefMesh< MBEDGE >:
[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 89 of file PostProcBrokenMeshInMoabBase.hpp.

Member Function Documentation

◆ generateReferenceElementMesh()

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

Implements MoFEM::PostProcGenerateRefMeshBase.

Definition at line 528 of file PostProcBrokenMeshInMoabBase.cpp.

528 {
530
531#ifndef NDEBUG
532 if (defMaxLevel > 0)
533 MOFEM_LOG("WORLD", Sev::warning)
534 << "Refinement for edges is not implemented";
535#endif
536
537 auto set_gauss_pts = [&](std::map<EntityHandle, int> &little_map) {
539
540 int nb_nodes = 2;
541 if (hoNodes)
542 nb_nodes = 3;
543
544 auto &gauss_pts = levelGaussPtsOnRefMesh[0];
545 gauss_pts.resize(2, nb_nodes, false);
546 gauss_pts.clear();
547
548 int nn = 0;
549 for (; nn != 2; ++nn) {
550 gauss_pts(0, nn) = static_cast<double>(nn);
551 little_map[nn] = nn;
552 }
553
554 if (nn < nb_nodes) {
555 gauss_pts(0, nn) = 0.5;
556 little_map[nn] = 2;
557 }
558
560 };
561
562 auto set_ref_edges = [&](std::map<EntityHandle, int> &little_map) {
564
565 int level = 0;
566 int nb_edges = level + 1;
567
568 int nb_nodes = 2;
569 if (hoNodes)
570 nb_nodes = 3;
571
572 auto &ref_edges = levelRef[level];
573 ref_edges.resize(nb_edges, nb_nodes, false);
574
575 for (int ee = 0; ee != nb_edges; ++ee) {
576 int nn = 0;
577 for (; nn != 2; ++nn) {
578 ref_edges(ee, nn) = nb_nodes * ee + nn;
579 }
580 if (nn < nb_nodes) {
581 ref_edges(ee, nn) = nb_nodes * ee + 2;
582 }
583 }
584
586 };
587
588 auto set_shape_functions = [&]() {
590 auto &gauss_pts = levelGaussPtsOnRefMesh[0];
591 auto &shape_functions = levelShapeFunctions[0];
592 const auto nb_gauss_pts = gauss_pts.size2();
593 shape_functions.resize(nb_gauss_pts, 2);
594 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
595 const double ksi = gauss_pts(0, gg);
596 shape_functions(gg, 0) = N_MBEDGE0(ksi);
597 shape_functions(gg, 1) = N_MBEDGE1(ksi);
598 }
600 };
601
602 levelRef.resize(1);
603 levelGaussPtsOnRefMesh.resize(1);
604 levelShapeFunctions.resize(1);
605
606 std::map<EntityHandle, int> little_map;
607 CHKERR set_gauss_pts(little_map);
608 CHKERR set_ref_edges(little_map);
609 CHKERR set_shape_functions();
610
612}
#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
#define N_MBEDGE0(x)
edge shape function
Definition: fem_tools.h:105
#define N_MBEDGE1(x)
edge shape function
Definition: fem_tools.h:106
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
std::vector< ublas::matrix< int > > levelRef
PetscBool hoNodes
int defMaxLevel
std::vector< MatrixDouble > levelShapeFunctions
std::vector< MatrixDouble > levelGaussPtsOnRefMesh

◆ PostProcGenerateRefMeshBase()

MoFEM::PostProcGenerateRefMeshBase::PostProcGenerateRefMeshBase ( )

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