v0.13.1
Public Member Functions | Private Member Functions | Private Attributes | List of all members
MoFEM::Simple Struct Reference

Simple interface for fast problem set-up. More...

#include <src/interfaces/Simple.hpp>

Inheritance diagram for MoFEM::Simple:
[legend]
Collaboration diagram for MoFEM::Simple:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 Simple (const MoFEM::Core &core)
 
virtual ~Simple ()=default
 
MoFEMErrorCode getOptions ()
 get options More...
 
MoFEMErrorCode loadFile (const std::string options, const std::string mesh_file_name)
 Load mesh file. More...
 
MoFEMErrorCode loadFile (const std::string mesh_file_name="")
 Load mesh file with parallel options if number of cores > 1. More...
 
MoFEMErrorCode addDomainField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
 Add field on domain. More...
 
MoFEMErrorCode addBoundaryField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
 Add field on boundary. More...
 
MoFEMErrorCode addSkeletonField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
 Add field on skeleton. More...
 
MoFEMErrorCode addDataField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
 Add field on domain. More...
 
MoFEMErrorCode removeDomainField (const std::string &name)
 Remove field form domain. More...
 
MoFEMErrorCode removeBoundaryField (const std::string &name)
 Remove field form boundary. More...
 
MoFEMErrorCode removeSkeletonField (const std::string &name)
 Remove field form skeleton. More...
 
MoFEMErrorCode defineFiniteElements ()
 Define finite elements. More...
 
MoFEMErrorCode defineProblem (const PetscBool is_partitioned=PETSC_TRUE)
 define problem More...
 
MoFEMErrorCode setFieldOrder (const std::string field_name, const int order, const Range *ents=NULL)
 Set field order. More...
 
MoFEMErrorCode buildFields ()
 Build fields. More...
 
MoFEMErrorCode buildFiniteElements ()
 Build finite elements. More...
 
MoFEMErrorCode setSkeletonAdjacency (int dim=-1)
 Set the skeleton adjacency object. More...
 
MoFEMErrorCode buildProblem ()
 Build problem. More...
 
MoFEMErrorCode setUp (const PetscBool is_partitioned=PETSC_TRUE)
 Setup problem. More...
 
MoFEMErrorCode reSetUp (bool only_dm=false)
 Rebuild internal MoFEM data structures. More...
 
MoFEMErrorCode getDM (DM *dm)
 Get DM. More...
 
SmartPetscObj< DM > getDM ()
 Return smart DM object. More...
 
int getDim () const
 Get the problem dimension. More...
 
void setDim (int dim)
 Set the problem dimension. More...
 
EntityHandle & getMeshSet ()
 Get the MeshSet object. More...
 
EntityHandle & getBoundaryMeshSet ()
 Get the BoundaryMeshSet object. More...
 
EntityHandle & getSkeletonMeshSet ()
 Get the SkeletonMeshSet object. More...
 
BitRefLevelgetBitRefLevel ()
 Get the BitRefLevel. More...
 
BitRefLevelgetBitRefLevelMask ()
 Get the BitRefLevel. More...
 
const std::string getDomainFEName () const
 Get the Domain FE Name. More...
 
const std::string getBoundaryFEName () const
 Get the Boundary FE Name. More...
 
const std::string getSkeletonFEName () const
 Get the Skeleton FE Name. More...
 
const std::string getProblemName () const
 Get the Problem Name. More...
 
std::string & getDomainFEName ()
 Get the Domain FE Name. More...
 
std::string & getBoundaryFEName ()
 Get the Boundary FE Name. More...
 
std::string & getSkeletonFEName ()
 Get the Skeleton FE Name. More...
 
std::string & getProblemName ()
 Get the Problem Name. More...
 
std::vector< std::string > & getOtherFiniteElements ()
 Get the Other Finite Elements. More...
 
MoFEMErrorCode deleteDM ()
 Delete dm. More...
 
MoFEMErrorCode deleteFiniteElements ()
 Delete finite elements. More...
 
boolgetAddSkeletonFE ()
 Get the addSkeletonFE. More...
 
boolgetAddBoundaryFE ()
 Get the addSkeletonFE. More...
 
boolgetParentAdjacencies ()
 Get the addParentAdjacencies. More...
 
auto & getBitAdjParent ()
 bit ref level for parent More...
 
auto & getBitAdjParentMask ()
 bit ref level for parent marent More...
 
auto & getBitAdjEnt ()
 bit ref level for parent More...
 
auto & getBitAdjEntMask ()
 bit ref level for parent marent More...
 
MoFEMErrorCode addFieldToEmptyFieldBlocks (const std::string row_field, const std::string col_field) const
 add empty block to problem More...
 
- Public Member Functions inherited from MoFEM::UnknownInterface
virtual MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 

Private Member Functions

MoFEMErrorCode createBoundaryMeshset ()
 
MoFEMErrorCode createSkeletonMeshset ()
 
MoFEMErrorCode exchangeGhostCells ()
 
template<int DIM = -1>
MoFEMErrorCode setSkeletonAdjacency ()
 
template<int DIM = -1>
MoFEMErrorCode setParentAdjacency ()
 
template<>
MoFEMErrorCode setSkeletonAdjacency ()
 
template<>
MoFEMErrorCode setSkeletonAdjacency ()
 
template<>
MoFEMErrorCode setSkeletonAdjacency ()
 
template<>
MoFEMErrorCode setParentAdjacency ()
 
template<>
MoFEMErrorCode setParentAdjacency ()
 
template<>
MoFEMErrorCode setParentAdjacency ()
 

Private Attributes

MoFEM::CorecOre
 
BitRefLevel bitLevel
 BitRefLevel of the probelm. More...
 
BitRefLevel bitLevelMask
 BitRefLevel of the probelm. More...
 
PetscLogEvent MOFEM_EVENT_SimpleSetUP
 
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
 
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
 
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
 
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
 
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
 
EntityHandle meshSet
 domain meshset More...
 
EntityHandle boundaryMeshset
 meshset with boundary More...
 
EntityHandle skeletonMeshset
 skeleton meshset with boundary More...
 
bool addSkeletonFE
 Add skeleton FE. More...
 
bool addBoundaryFE
 Add boundary FE. More...
 
bool addParentAdjacencies
 If set to true parent adjacencies are build. More...
 
BitRefLevel bitAdjParent
 bit ref level for parent More...
 
BitRefLevel bitAdjParentMask
 bit ref level for parent marent More...
 
BitRefLevel bitAdjEnt
 bit ref level for parent More...
 
BitRefLevel bitAdjEntMask
 bit ref level for parent marent More...
 
std::vector< std::string > domainFields
 domain fields More...
 
std::vector< std::string > boundaryFields
 boundary fields More...
 
std::vector< std::string > skeletonFields
 fields on the skeleton More...
 
std::vector< std::string > dataFields
 Data fields. More...
 
std::vector< std::string > noFieldFields
 NOFIELD field name. More...
 
std::vector< std::string > noFieldDataFields
 NOFIELD field name. More...
 
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
 fields order More...
 
std::string nameOfProblem
 problem name More...
 
std::string domainFE
 domain finite element More...
 
std::string boundaryFE
 boundary finite element More...
 
std::string skeletonFE
 skeleton finite element More...
 
std::vector< std::string > otherFEs
 Other finite elements. More...
 
char meshFileName [255]
 mesh file name More...
 
int dIm
 dimension of problem More...
 
SmartPetscObj< DM > dM
 Discrete manager (interface to PETSc/MoFEM functions) More...
 
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
 
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 2 > > parentAdjFunctionDim2
 
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 1 > > parentAdjFunctionDim1
 

Additional Inherited Members

- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version. More...
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version. More...
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version. More...
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version. More...
 

Detailed Description

Simple interface for fast problem set-up.

Examples
approx_sphere.cpp, boundary_marker.cpp, child_and_parent.cpp, contact.cpp, continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, eigen_elastic.cpp, field_evaluator.cpp, free_surface.cpp, hanging_node_approx.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, heat_method.cpp, hello_world.cpp, helmholtz.cpp, higher_derivatives.cpp, loop_entities.cpp, mesh_smoothing.cpp, mixed_poisson.cpp, nonlinear_elastic.cpp, partition_mesh.cpp, phase.cpp, photon_diffusion.cpp, plastic.cpp, plate.cpp, plot_base.cpp, poisson_2d_dis_galerkin.cpp, poisson_2d_homogeneous.cpp, reaction_diffusion.cpp, scalar_check_approximation.cpp, shallow_wave.cpp, simple_elasticity.cpp, simple_interface.cpp, simple_l2_only.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, thermo_plastic.cpp, and wave_equation.cpp.

Definition at line 35 of file Simple.hpp.

Constructor & Destructor Documentation

◆ Simple()

MoFEM::Simple::Simple ( const MoFEM::Core core)

Definition at line 276 of file Simple.cpp.

277 : cOre(const_cast<Core &>(core)), bitLevel(BitRefLevel().set(0)),
279 skeletonMeshset(0), nameOfProblem("SimpleProblem"), domainFE("dFE"),
280 boundaryFE("bFE"), skeletonFE("sFE"), dIm(-1), addSkeletonFE(false),
281 addBoundaryFE(false), addParentAdjacencies(false),
284 PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleSetUP);
285 PetscLogEventRegister("SimpleLoadMesh", 0, &MOFEM_EVENT_SimpleLoadMesh);
286 PetscLogEventRegister("SimpleBuildFields", 0, &MOFEM_EVENT_SimpleBuildFields);
287 PetscLogEventRegister("SimpleBuildFiniteElements", 0,
289 PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleBuildProblem);
290 PetscLogEventRegister("SimpleKSPSolve", 0, &MOFEM_EVENT_SimpleKSPSolve);
291 strcpy(meshFileName, "mesh.h5m");
292}
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
CoreTmp< 0 > Core
Definition: Core.hpp:1096
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:485
MoFEM::Core & cOre
Definition: Simple.hpp:472
BitRefLevel bitAdjEnt
bit ref level for parent
Definition: Simple.hpp:494
char meshFileName[255]
mesh file name
Definition: Simple.hpp:514
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:484
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:479
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:478
bool addSkeletonFE
Add skeleton FE.
Definition: Simple.hpp:488
BitRefLevel bitAdjParent
bit ref level for parent
Definition: Simple.hpp:492
bool addBoundaryFE
Add boundary FE.
Definition: Simple.hpp:489
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:510
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:509
std::string nameOfProblem
problem name
Definition: Simple.hpp:507
std::string domainFE
domain finite element
Definition: Simple.hpp:508
BitRefLevel bitAdjEntMask
bit ref level for parent marent
Definition: Simple.hpp:495
BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:474
int dIm
dimension of problem
Definition: Simple.hpp:515
PetscLogEvent MOFEM_EVENT_SimpleSetUP
Definition: Simple.hpp:477
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:482
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:480
BitRefLevel bitLevelMask
BitRefLevel of the probelm.
Definition: Simple.hpp:475
BitRefLevel bitAdjParentMask
bit ref level for parent marent
Definition: Simple.hpp:493
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:481
EntityHandle skeletonMeshset
skeleton meshset with boundary
Definition: Simple.hpp:486
bool addParentAdjacencies
If set to true parent adjacencies are build.
Definition: Simple.hpp:490

◆ ~Simple()

virtual MoFEM::Simple::~Simple ( )
virtualdefault

Member Function Documentation

◆ addBoundaryField()

MoFEMErrorCode MoFEM::Simple::addBoundaryField ( const std::string &  name,
const FieldSpace  space,
const FieldApproximationBase  base,
const FieldCoefficientsNumber  nb_of_coefficients,
const TagType  tag_type = MB_TAG_SPARSE,
const enum MoFEMTypes  bh = MF_ZERO,
int  verb = -1 
)

Add field on boundary.

Parameters
namename of the filed
spacespace (L2,H1,Hdiv,Hcurl)
baseapproximation base, see FieldApproximationBase
nb_of_coefficientsnumber of field coefficients
tag_typetype of the tag MB_TAG_DENSE or MB_TAG_SPARSE (DENSE is faster and uses less memory, SPARSE is more flexible if you define field on subdomains)
bhif MF_EXCL throws error if field exits, MF_ZERO no error if field exist
verbverbosity level
Returns
error code
Examples
hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, mesh_smoothing.cpp, and simple_interface.cpp.

Definition at line 396 of file Simple.cpp.

400 {
401 Interface &m_field = cOre;
403 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
404 verb);
405 boundaryFields.push_back(name);
406 if (space == NOFIELD)
407 SETERRQ(
408 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
409 "NOFIELD space for boundary filed not implemented in Simple interface");
411}
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:97
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:498

◆ addDataField()

MoFEMErrorCode MoFEM::Simple::addDataField ( const std::string &  name,
const FieldSpace  space,
const FieldApproximationBase  base,
const FieldCoefficientsNumber  nb_of_coefficients,
const TagType  tag_type = MB_TAG_SPARSE,
const enum MoFEMTypes  bh = MF_ZERO,
int  verb = -1 
)

Add field on domain.

Parameters
namename of the filed
spacespace (L2,H1,Hdiv,Hcurl)
baseapproximation base, see FieldApproximationBase
nb_of_coefficientsnumber of field coefficients
tag_typetype of the tag MB_TAG_DENSE or MB_TAG_SPARSE (DENSE is faster and uses less memory, SPARSE is more flexible if you define field on subdomains)
bhif MF_EXCL throws error if field exits, MF_ZERO no error if field exist
verbverbosity level
Returns
error code
Examples
hcurl_curl_operator.cpp, and hdiv_divergence_operator.cpp.

Definition at line 434 of file Simple.cpp.

438 {
439
440 Interface &m_field = cOre;
442 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
443 verb);
444 if (space == NOFIELD)
445 noFieldDataFields.push_back(name);
446 else
447 dataFields.push_back(name);
449}
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:500
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:502

◆ addDomainField()

MoFEMErrorCode MoFEM::Simple::addDomainField ( const std::string &  name,
const FieldSpace  space,
const FieldApproximationBase  base,
const FieldCoefficientsNumber  nb_of_coefficients,
const TagType  tag_type = MB_TAG_SPARSE,
const enum MoFEMTypes  bh = MF_ZERO,
int  verb = -1 
)

Add field on domain.

Parameters
namename of the filed
spacespace (L2,H1,Hdiv,Hcurl)
baseapproximation base, see FieldApproximationBase
nb_of_coefficientsnumber of field coefficients
tag_typetype of the tag MB_TAG_DENSE or MB_TAG_SPARSE (DENSE is faster and uses less memory, SPARSE is more flexible if you define field on subdomains)
bhif MF_EXCL throws error if field exits, MF_ZERO no error if field exist
verbverbosity level
Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation.cpp, simple_elasticity.cpp, simple_interface.cpp, simple_l2_only.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 378 of file Simple.cpp.

382 {
383
384 Interface &m_field = cOre;
386 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
387 verb);
388 if (space == NOFIELD)
389 noFieldFields.push_back(name);
390 else
391 domainFields.push_back(name);
393}
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:497
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:501

◆ addFieldToEmptyFieldBlocks()

MoFEMErrorCode MoFEM::Simple::addFieldToEmptyFieldBlocks ( const std::string  row_field,
const std::string  col_field 
) const

add empty block to problem

MatrixManager assumes that all blocks, i.e. all fields combinations are non zero. This is not always the case, to optimise code and reduce memory usage user can specifi which blocks are empty.

Parameters
row_fieldrow filed name
col_fieldcol field name
Returns
MoFEMErrorCode

Definition at line 845 of file Simple.cpp.

846 {
847 Interface &m_field = cOre;
849 CHKERR m_field.getInterface<ProblemsManager>()->addFieldToEmptyFieldBlocks(
850 getProblemName(), row_field, col_field);
852}
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
add empty block to problem
Definition: Simple.cpp:845
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:342

◆ addSkeletonField()

MoFEMErrorCode MoFEM::Simple::addSkeletonField ( const std::string &  name,
const FieldSpace  space,
const FieldApproximationBase  base,
const FieldCoefficientsNumber  nb_of_coefficients,
const TagType  tag_type = MB_TAG_SPARSE,
const enum MoFEMTypes  bh = MF_ZERO,
int  verb = -1 
)

Add field on skeleton.

Parameters
namename of the filed
spacespace (L2,H1,Hdiv,Hcurl)
baseapproximation base, see FieldApproximationBase
nb_of_coefficientsnumber of field coefficients
tag_typetype of the tag MB_TAG_DENSE or MB_TAG_SPARSE (DENSE is faster and uses less memory, SPARSE is more flexible if you define field on subdomains)
bhif MF_EXCL throws error if field exits, MF_ZERO no error if field exist
verbverbosity level
Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, and simple_interface.cpp.

Definition at line 414 of file Simple.cpp.

418 {
419
420 Interface &m_field = cOre;
422 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
423 verb);
424 skeletonFields.push_back(name);
425 if (space == NOFIELD)
426 SETERRQ(
427 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
428 "NOFIELD space for boundary filed not implemented in Simple interface");
429
431}
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:499

◆ buildFields()

MoFEMErrorCode MoFEM::Simple::buildFields ( )

Build fields.

Returns
error code
Examples
mesh_smoothing.cpp, and simple_elasticity.cpp.

Definition at line 601 of file Simple.cpp.

601 {
602 Interface &m_field = cOre;
604 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
605
606 auto comm_interface_ptr = m_field.getInterface<CommInterface>();
607 auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
608
609 // Add entities to the fields
610 auto add_ents_to_field = [&](auto meshset, auto dim, auto &fields) {
612 for (auto &field : fields) {
613 CHKERR m_field.add_ents_to_field_by_dim(meshset, dim, field);
614 CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
615 }
617 };
618
619 auto make_no_field_ents = [&](auto &fields) {
621 for (auto &field : fields) {
622 std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
623 CHKERR m_field.create_vertices_and_add_to_field(field, coords.data(), 1);
624 CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
625 CHKERR bit_ref_ptr->setFieldEntitiesBitRefLevel(field, bitLevel);
626 }
628 };
629
630 CHKERR add_ents_to_field(meshSet, dIm, domainFields);
631 CHKERR add_ents_to_field(meshSet, dIm, dataFields);
632 CHKERR add_ents_to_field(boundaryMeshset, dIm - 1, boundaryFields);
633 CHKERR add_ents_to_field(skeletonMeshset, dIm - 1, skeletonFields);
634
635 std::set<std::string> nofield_fields;
636 for (auto &f : noFieldFields)
637 nofield_fields.insert(f);
638 for (auto &f : noFieldDataFields)
639 nofield_fields.insert(f);
640
641 CHKERR make_no_field_ents(nofield_fields);
642
643 // Set order
644 auto set_order = [&](auto meshset, auto dim, auto &fields) {
646 for (auto &field : fields) {
647 if (fieldsOrder.find(field) == fieldsOrder.end()) {
648 SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
649 "Order for field not set %s", field.c_str());
650 }
651 int dds = 0;
652 const Field *field_ptr = m_field.get_field_structure(field);
653 switch (field_ptr->getSpace()) {
654 case L2:
655 dds = dim;
656 break;
657 case HDIV:
658 dds = 2;
659 break;
660 case HCURL:
661 dds = 1;
662 break;
663 case H1:
664 dds = 1;
665 break;
666 default:
667 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
668 "Glasgow we have a problem");
669 }
670
671 auto set_order = [&](auto field, auto &ents) {
673
674 auto range = fieldsOrder.equal_range(field);
675 for (auto o = range.first; o != range.second; ++o) {
676 if (!o->second.second.empty())
677 ents = intersect(ents, o->second.second);
678 CHKERR m_field.set_field_order(ents, field, o->second.first);
679 }
680
682 };
683
684 if (field_ptr->getSpace() == H1) {
685 if (field_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
686 Range ents;
687 CHKERR m_field.get_field_entities_by_dimension(field, 0, ents);
688 CHKERR set_order(field, ents);
689 } else {
690 CHKERR m_field.set_field_order(meshSet, MBVERTEX, field, 1);
691 }
692 }
693 for (int dd = dds; dd <= dim; dd++) {
694 Range ents;
695 CHKERR m_field.get_field_entities_by_dimension(field, dd, ents);
696 CHKERR set_order(field, ents);
697 }
698 }
700 };
701
702 CHKERR set_order(meshSet, dIm, domainFields);
703 CHKERR set_order(meshSet, dIm, dataFields);
704 CHKERR set_order(meshSet, dIm - 1, boundaryFields);
705 CHKERR set_order(meshSet, dIm - 1, skeletonFields);
706
707 // Build fields
708 CHKERR m_field.build_fields();
709 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
711}
static Index< 'o', 3 > o
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:77
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ L2
field with C-1 continuity
Definition: definitions.h:101
@ H1
continuous field
Definition: definitions.h:98
@ HCURL
field with continuous tangents
Definition: definitions.h:99
@ HDIV
field with continuous normal traction
Definition: definitions.h:100
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
@ MOFEM_INVALID_DATA
Definition: definitions.h:49
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
const int dim
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:505

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

Returns
error code
Examples
mesh_smoothing.cpp, and simple_elasticity.cpp.

Definition at line 713 of file Simple.cpp.

713 {
714 Interface &m_field = cOre;
716 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
717 // Add finite elements
718 CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
719 true);
720 CHKERR m_field.build_finite_elements(domainFE);
721 if (addBoundaryFE || !boundaryFields.empty()) {
722 CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
723 boundaryFE, true);
724 CHKERR m_field.build_finite_elements(boundaryFE);
725 }
726 if (addSkeletonFE || !skeletonFields.empty()) {
727 CHKERR m_field.add_ents_to_finite_element_by_dim(skeletonMeshset, dIm - 1,
728 skeletonFE, true);
729 CHKERR m_field.build_finite_elements(skeletonFE);
730 }
731 for (std::vector<std::string>::iterator fit = otherFEs.begin();
732 fit != otherFEs.end(); ++fit) {
733 CHKERR m_field.build_finite_elements(*fit);
734 }
735 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
737}
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:512

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

Returns
error code
Examples
mesh_smoothing.cpp, and simple_elasticity.cpp.

Definition at line 739 of file Simple.cpp.

739 {
740 Interface &m_field = cOre;
742 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
743 CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
744 // Set problem by the DOFs on the fields rather that by adding DOFs on the
745 // elements
746 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
748 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
749 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
751}
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1224
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:518

◆ createBoundaryMeshset()

MoFEMErrorCode MoFEM::Simple::createBoundaryMeshset ( )
private

Definition at line 854 of file Simple.cpp.

854 {
855 Interface &m_field = cOre;
857 ParallelComm *pcomm =
858 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
859
860 auto get_skin = [&](auto meshset) {
861 // filter not owned entities, those are not on boundary
862
863 Range domain_ents;
864 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
865 domain_ents, true);
866 CHKERR pcomm->filter_pstatus(domain_ents,
867 PSTATUS_SHARED | PSTATUS_MULTISHARED,
868 PSTATUS_NOT, -1, nullptr);
869
870 Skinner skin(&m_field.get_moab());
871 Range domain_skin;
872 CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
873 CHKERR pcomm->filter_pstatus(domain_skin,
874 PSTATUS_SHARED | PSTATUS_MULTISHARED,
875 PSTATUS_NOT, -1, nullptr);
876 return domain_skin;
877 };
878
879 auto create_boundary_meshset = [&](auto &&domain_skin) {
881 // create boundary meshset
882 if (boundaryMeshset != 0) {
884 }
885 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
886 CHKERR m_field.get_moab().add_entities(boundaryMeshset, domain_skin);
887 for (int dd = 0; dd != dIm - 1; dd++) {
888 Range adj;
889 CHKERR m_field.get_moab().get_adjacencies(domain_skin, dd, false, adj,
890 moab::Interface::UNION);
891 CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
892 }
894 };
895
896 CHKERR create_boundary_meshset(get_skin(meshSet));
897
899}
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:228

◆ createSkeletonMeshset()

MoFEMErrorCode MoFEM::Simple::createSkeletonMeshset ( )
private

Definition at line 901 of file Simple.cpp.

901 {
902 Interface &m_field = cOre;
904
905 ParallelComm *pcomm =
906 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
907
908 auto create_skeleton_meshset = [&](auto meshset) {
910 // create boundary meshset
911 if (skeletonMeshset != 0) {
913 }
914 Range boundary_ents, skeleton_ents;
915 CHKERR m_field.get_moab().get_entities_by_dimension(boundaryMeshset,
916 dIm - 1, boundary_ents);
917 Range domain_ents;
918 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
919 domain_ents, true);
920 CHKERR m_field.get_moab().get_adjacencies(
921 domain_ents, dIm - 1, false, skeleton_ents, moab::Interface::UNION);
922 skeleton_ents = subtract(skeleton_ents, boundary_ents);
923 CHKERR pcomm->filter_pstatus(skeleton_ents, PSTATUS_NOT_OWNED, PSTATUS_NOT,
924 -1, nullptr);
925 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, skeletonMeshset);
926 CHKERR m_field.get_moab().add_entities(skeletonMeshset, skeleton_ents);
928 };
929
930 CHKERR create_skeleton_meshset(meshSet);
931
933}

◆ defineFiniteElements()

MoFEMErrorCode MoFEM::Simple::defineFiniteElements ( )

Define finite elements.

Returns
Error code
Examples
mesh_smoothing.cpp, and simple_elasticity.cpp.

Definition at line 497 of file Simple.cpp.

497 {
498 Interface &m_field = cOre;
500
501 auto clear_rows_and_cols = [&](auto &fe_name) {
503 auto fe_ptr = m_field.get_finite_elements();
504 auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
505 ->get<FiniteElement_name_mi_tag>();
506 auto it_fe = fe_by_name.find(fe_name);
507 if (it_fe != fe_by_name.end()) {
508
509 if (!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
510 SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
511 "modification unsuccessful");
512
513 if (!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
514 SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
515 "modification unsuccessful");
516 }
518 };
519 CHKERR clear_rows_and_cols(domainFE);
520 CHKERR clear_rows_and_cols(boundaryFE);
521 CHKERR clear_rows_and_cols(skeletonFE);
522
523 // Define finite elements
524 CHKERR m_field.add_finite_element(domainFE, MF_ZERO);
525
526 auto add_fields = [&](auto &fe_name, auto &fields) {
528 for (auto &field : fields) {
529 CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
530 CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
531 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
532 }
534 };
535
536 auto add_data_fields = [&](auto &fe_name, auto &fields) {
538 for (auto &field : fields)
539 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
541 };
542
543 CHKERR add_fields(domainFE, domainFields);
544 CHKERR add_data_fields(domainFE, dataFields);
545 CHKERR add_fields(domainFE, noFieldFields);
546 CHKERR add_data_fields(domainFE, noFieldDataFields);
547
548 if (addBoundaryFE || !boundaryFields.empty()) {
549 CHKERR m_field.add_finite_element(boundaryFE, MF_ZERO);
550 CHKERR add_fields(boundaryFE, domainFields);
551 if (!boundaryFields.empty())
552 CHKERR add_fields(boundaryFE, boundaryFields);
553 CHKERR add_data_fields(boundaryFE, dataFields);
554 CHKERR add_data_fields(boundaryFE, noFieldDataFields);
555 CHKERR add_fields(boundaryFE, noFieldFields);
556 }
557 if (addSkeletonFE || !skeletonFields.empty()) {
558 CHKERR m_field.add_finite_element(skeletonFE, MF_ZERO);
559 CHKERR add_fields(skeletonFE, domainFields);
560 if (!skeletonFields.empty())
561 CHKERR add_fields(skeletonFE, skeletonFields);
562 CHKERR add_data_fields(skeletonFE, dataFields);
563 CHKERR add_data_fields(skeletonFE, noFieldDataFields);
564 CHKERR add_fields(skeletonFE, noFieldFields);
565 }
567}
@ MF_ZERO
Definition: definitions.h:111
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition: definitions.h:47
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.

◆ defineProblem()

MoFEMErrorCode MoFEM::Simple::defineProblem ( const PetscBool  is_partitioned = PETSC_TRUE)

define problem

Returns
error code
Examples
mesh_smoothing.cpp, and simple_elasticity.cpp.

Definition at line 569 of file Simple.cpp.

569 {
570 Interface &m_field = cOre;
572 // Create dm instance
573 dM = createSmartDM(m_field.get_comm(), "DMMOFEM");
574 // set dm data structure which created mofem data structures
577 CHKERR DMSetFromOptions(dM);
579 if (addBoundaryFE || !boundaryFields.empty()) {
581 }
582 if (addSkeletonFE || !skeletonFields.empty()) {
584 }
585 for (std::vector<std::string>::iterator fit = otherFEs.begin();
586 fit != otherFEs.end(); ++fit) {
587 CHKERR DMMoFEMAddElement(dM, fit->c_str());
588 }
589 CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
591}
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:1082
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:130
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:462
auto createSmartDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.

◆ deleteDM()

MoFEMErrorCode MoFEM::Simple::deleteDM ( )

Delete dm.

Returns
MoFEMErrorCode

Definition at line 822 of file Simple.cpp.

822 {
824 dM.reset();
826}

◆ deleteFiniteElements()

MoFEMErrorCode MoFEM::Simple::deleteFiniteElements ( )

Delete finite elements.

Returns
MoFEMErrorCode

Definition at line 833 of file Simple.cpp.

833 {
834 Interface &m_field = cOre;
836 for (auto fe : {domainFE, boundaryFE, skeletonFE}) {
837 if (m_field.check_finite_element(fe)) {
838 CHKERR m_field.delete_finite_element(fe);
839 }
840 }
842}

◆ exchangeGhostCells()

MoFEMErrorCode MoFEM::Simple::exchangeGhostCells ( )
private

Definition at line 935 of file Simple.cpp.

935 {
936 Interface &m_field = cOre;
938 MOFEM_LOG("WORLD", Sev::verbose) << "Exchange ghost cells";
939
940 ParallelComm *pcomm =
941 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
942 if (pcomm == NULL)
943 pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
944
945 Range verts;
946 CHKERR m_field.get_moab().get_entities_by_type(0, MBVERTEX, verts);
947 CHKERR pcomm->exchange_ghost_cells(getDim(), getDim() - 1, 1,
948 3 /**get all adjacent ghosted entities */,
949 true, false, meshSet ? &meshSet : nullptr);
950
951 Range shared;
952 CHKERR m_field.get_moab().get_entities_by_dimension(0, dIm, shared);
953 for (auto d = dIm - 1; d >= 1; --d) {
954 CHKERR m_field.get_moab().get_adjacencies(shared, d, false, shared,
955 moab::Interface::UNION);
956 }
957 CHKERR pcomm->filter_pstatus(shared, PSTATUS_SHARED | PSTATUS_MULTISHARED,
958 PSTATUS_OR, -1, &shared);
959 Tag part_tag = pcomm->part_tag();
960 CHKERR pcomm->exchange_tags(part_tag, shared);
961
963}
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
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)
Definition: dTensor0.hpp:27
int getDim() const
Get the problem dimension.
Definition: Simple.hpp:272

◆ getAddBoundaryFE()

bool & MoFEM::Simple::getAddBoundaryFE ( )

Get the addSkeletonFE.

If variable set to true, boundary element is created regardless field on skelton is added or not.

Returns
true
false
Examples
simple_l2_only.cpp.

Definition at line 416 of file Simple.hpp.

416{ return addBoundaryFE; }

◆ getAddSkeletonFE()

bool & MoFEM::Simple::getAddSkeletonFE ( )

Get the addSkeletonFE.

If variable set to true, skeleton element is created regardless field on skelton is added or not.

Returns
true
false
Examples
continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, and simple_l2_only.cpp.

Definition at line 405 of file Simple.hpp.

405{ return addSkeletonFE; }

◆ getBitAdjEnt()

auto & MoFEM::Simple::getBitAdjEnt ( )

bit ref level for parent

Returns
auto&

Definition at line 447 of file Simple.hpp.

447{ return bitAdjEnt; }

◆ getBitAdjEntMask()

auto & MoFEM::Simple::getBitAdjEntMask ( )

bit ref level for parent marent

Returns
auto&

Definition at line 454 of file Simple.hpp.

454{ return bitAdjEntMask; }

◆ getBitAdjParent()

auto & MoFEM::Simple::getBitAdjParent ( )

bit ref level for parent

Returns
auto&

Definition at line 433 of file Simple.hpp.

433{ return bitAdjParent; }

◆ getBitAdjParentMask()

auto & MoFEM::Simple::getBitAdjParentMask ( )

bit ref level for parent marent

Returns
auto&

Definition at line 440 of file Simple.hpp.

440{ return bitAdjParentMask; }

◆ getBitRefLevel()

BitRefLevel & MoFEM::Simple::getBitRefLevel ( )

Get the BitRefLevel.

Returns
BitRefLevel

Definition at line 307 of file Simple.hpp.

307{ return bitLevel; }

◆ getBitRefLevelMask()

BitRefLevel & MoFEM::Simple::getBitRefLevelMask ( )

Get the BitRefLevel.

Returns
BitRefLevel

Definition at line 314 of file Simple.hpp.

314{ return bitLevelMask; }

◆ getBoundaryFEName() [1/2]

std::string & MoFEM::Simple::getBoundaryFEName ( )

Get the Boundary FE Name.

Returns
std::string&

Definition at line 356 of file Simple.hpp.

356{ return boundaryFE; }

◆ getBoundaryFEName() [2/2]

const std::string MoFEM::Simple::getBoundaryFEName ( ) const

Get the Boundary FE Name.

Returns
const std::string
Examples
mesh_smoothing.cpp, and simple_interface.cpp.

Definition at line 328 of file Simple.hpp.

328{ return boundaryFE; }

◆ getBoundaryMeshSet()

EntityHandle & MoFEM::Simple::getBoundaryMeshSet ( )

Get the BoundaryMeshSet object.

Returns
EntityHandle&

Definition at line 293 of file Simple.hpp.

293{ return boundaryMeshset; }

◆ getDim()

int MoFEM::Simple::getDim ( ) const

Get the problem dimension.

Problem dimension is determined by highest dimension of entities on the mesh.

Returns
int
Examples
field_evaluator.cpp.

Definition at line 272 of file Simple.hpp.

272{ return dIm; }

◆ getDM() [1/2]

SmartPetscObj< DM > MoFEM::Simple::getDM ( )

Return smart DM object.

{
auto dm = simple_interface->getDM();
// ...
// dm is automatically destroyed when out of the scope
}
Returns
SmartPetscObj<DM>

Definition at line 262 of file Simple.hpp.

262{ return dM; }

◆ getDM() [2/2]

MoFEMErrorCode MoFEM::Simple::getDM ( DM *  dm)

◆ getDomainFEName() [1/2]

std::string & MoFEM::Simple::getDomainFEName ( )

Get the Domain FE Name.

Returns
std::string&

Definition at line 349 of file Simple.hpp.

349{ return domainFE; }

◆ getDomainFEName() [2/2]

const std::string MoFEM::Simple::getDomainFEName ( ) const

◆ getMeshSet()

EntityHandle & MoFEM::Simple::getMeshSet ( )

Get the MeshSet object.

Returns
EntityHandle&

Definition at line 286 of file Simple.hpp.

286{ return meshSet; }

◆ getOptions()

MoFEMErrorCode MoFEM::Simple::getOptions ( )

get options

Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation.cpp, simple_elasticity.cpp, simple_interface.cpp, simple_l2_only.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 294 of file Simple.cpp.

294 {
295 PetscBool flg = PETSC_TRUE;
297 ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "Simple interface options",
298 "none");
299 CHKERRG(ierr);
300 ierr = PetscOptionsString("-file_name", "file name", "", "mesh.h5m",
301 meshFileName, 255, &flg);
302 CHKERRG(ierr);
303 ierr = PetscOptionsEnd();
304 CHKERRG(ierr);
306}
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87

◆ getOtherFiniteElements()

std::vector< std::string > & MoFEM::Simple::getOtherFiniteElements ( )

Get the Other Finite Elements.

User can create finite elements using directly core interface and and them to simple problem by this function

Returns
std::vector<std::string>&
Examples
mesh_smoothing.cpp.

Definition at line 380 of file Simple.hpp.

380{ return otherFEs; }

◆ getParentAdjacencies()

bool & MoFEM::Simple::getParentAdjacencies ( )

Get the addParentAdjacencies.

If set true add parent adjacencies

Returns
true
false

Definition at line 426 of file Simple.hpp.

426{ return addParentAdjacencies; }

◆ getProblemName() [1/2]

std::string & MoFEM::Simple::getProblemName ( )

Get the Problem Name.

Returns
std::string&

Definition at line 370 of file Simple.hpp.

370{ return nameOfProblem; }

◆ getProblemName() [2/2]

const std::string MoFEM::Simple::getProblemName ( ) const

Get the Problem Name.

Returns
const std::string
Examples
reaction_diffusion.cpp, and test_cache_on_entities.cpp.

Definition at line 342 of file Simple.hpp.

342{ return nameOfProblem; }

◆ getSkeletonFEName() [1/2]

std::string & MoFEM::Simple::getSkeletonFEName ( )

Get the Skeleton FE Name.

Returns
std::string&

Definition at line 363 of file Simple.hpp.

363{ return skeletonFE; }

◆ getSkeletonFEName() [2/2]

const std::string MoFEM::Simple::getSkeletonFEName ( ) const

◆ getSkeletonMeshSet()

EntityHandle & MoFEM::Simple::getSkeletonMeshSet ( )

Get the SkeletonMeshSet object.

Returns
EntityHandle&

Definition at line 300 of file Simple.hpp.

300{ return skeletonMeshset; }

◆ loadFile() [1/2]

MoFEMErrorCode MoFEM::Simple::loadFile ( const std::string  mesh_file_name = "")

Load mesh file with parallel options if number of cores > 1.

Parameters
mesh_file_namefile name if not set default or set by command line is used.
Returns
error code

Definition at line 364 of file Simple.cpp.

364 {
366 Interface &m_field = cOre;
367 if (m_field.get_comm_size() == 1)
369 else
370 CHKERR loadFile("PARALLEL=READ_PART;"
371 "PARALLEL_RESOLVE_SHARED_ENTS;"
372 "PARTITION=PARALLEL_PARTITION;",
375}
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:308

◆ loadFile() [2/2]

MoFEMErrorCode MoFEM::Simple::loadFile ( const std::string  options,
const std::string  mesh_file_name 
)

Load mesh file.

Parameters
optionsfile load options
mesh_file_namefile name if not set default or set by command line is used.
Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion.cpp, scalar_check_approximation.cpp, simple_elasticity.cpp, simple_interface.cpp, simple_l2_only.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 308 of file Simple.cpp.

309 {
310 Interface &m_field = cOre;
312 PetscLogEventBegin(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
313
314 if (!mesh_file_name.empty())
315 strcpy(meshFileName, mesh_file_name.c_str());
316
317 // This is a case of distributed mesh and algebra. In that case each processor
318 // keep only part of the problem.
319 CHKERR m_field.get_moab().load_file(meshFileName, 0, options.c_str());
320 CHKERR m_field.rebuild_database();
321
322 // determine problem dimension
323 if (dIm == -1) {
324 int nb_ents_3d;
325 CHKERR m_field.get_moab().get_number_entities_by_dimension(
326 meshSet, 3, nb_ents_3d, true);
327 if (nb_ents_3d > 0) {
328 dIm = 3;
329 } else {
330 int nb_ents_2d;
331 CHKERR m_field.get_moab().get_number_entities_by_dimension(
332 meshSet, 2, nb_ents_2d, true);
333 if (nb_ents_2d > 0) {
334 dIm = 2;
335 } else {
336 dIm = 1;
337 }
338 }
339 }
340
341 if (!boundaryMeshset)
343 if (!skeletonMeshset)
345 if (addSkeletonFE)
347
348 if (bitLevel.any()) {
349 Range ents;
350 CHKERR m_field.get_moab().get_entities_by_dimension(meshSet, dIm, ents,
351 true);
352 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(ents, bitLevel,
353 false);
354 } else {
355 MOFEM_LOG("WORLD", Sev::warning) << "BitRefLevel is none and not set";
356 CHKERR m_field.getInterface<BitRefManager>()->addToDatabaseBitRefLevelByDim(
357 dIm, BitRefLevel().set(), BitRefLevel().set());
358 }
359
360 PetscLogEventEnd(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
362}
char mesh_file_name[255]
MoFEMErrorCode exchangeGhostCells()
Definition: Simple.cpp:935
MoFEMErrorCode createBoundaryMeshset()
Definition: Simple.cpp:854
MoFEMErrorCode createSkeletonMeshset()
Definition: Simple.cpp:901

◆ query_interface()

MoFEMErrorCode MoFEM::Simple::query_interface ( boost::typeindex::type_index  type_index,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 22 of file Simple.cpp.

23 {
24 *iface = const_cast<Simple *>(this);
25 return 0;
26}
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:276

◆ removeBoundaryField()

MoFEMErrorCode MoFEM::Simple::removeBoundaryField ( const std::string &  name)

Remove field form boundary.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 467 of file Simple.cpp.

467 {
468 Interface &m_field = cOre;
470
471 auto remove_field_from_list = [&](auto &vec) {
472 auto it = std::find(vec.begin(), vec.end(), name);
473 if (it != vec.end())
474 vec.erase(it);
475 };
476
477 remove_field_from_list(boundaryFields);
478
480}

◆ removeDomainField()

MoFEMErrorCode MoFEM::Simple::removeDomainField ( const std::string &  name)

Remove field form domain.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 451 of file Simple.cpp.

451 {
452 Interface &m_field = cOre;
454
455 auto remove_field_from_list = [&](auto &vec) {
456 auto it = std::find(vec.begin(), vec.end(), name);
457 if (it != vec.end())
458 vec.erase(it);
459 };
460
461 remove_field_from_list(noFieldFields);
462 remove_field_from_list(domainFields);
463
465}

◆ removeSkeletonField()

MoFEMErrorCode MoFEM::Simple::removeSkeletonField ( const std::string &  name)

Remove field form skeleton.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 482 of file Simple.cpp.

482 {
483 Interface &m_field = cOre;
485
486 auto remove_field_from_list = [&](auto &vec) {
487 auto it = std::find(vec.begin(), vec.end(), name);
488 if (it != vec.end())
489 vec.erase(it);
490 };
491
492 remove_field_from_list(skeletonFields);
493
495}

◆ reSetUp()

MoFEMErrorCode MoFEM::Simple::reSetUp ( bool  only_dm = false)

Rebuild internal MoFEM data structures.

Call this function after you add field or remove it.

Note
If you add field, or remove it, finite element and problem needs to be rebuild. However DM can remain the same.
Returns
MoFEMErrorCode

Definition at line 776 of file Simple.cpp.

776 {
777 Interface &m_field = cOre;
779 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
780
781 if (!only_dm) {
784 if (addSkeletonFE || !skeletonFields.empty())
789 }
790
791 CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
792
793 const Problem *problem_ptr;
794 CHKERR DMMoFEMGetProblemPtr(dM, &problem_ptr);
795 const auto problem_name = problem_ptr->getName();
796 CHKERR m_field.modify_problem_ref_level_set_bit(problem_name, bitLevel);
797 CHKERR m_field.modify_problem_mask_ref_level_set_bit(problem_name,
799
800 // Set problem by the DOFs on the fields rather that by adding DOFs on the
801 // elements
802 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
804 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
805
806 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
808}
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:385
MoFEMErrorCode setParentAdjacency()
Definition: Simple.cpp:185
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:497
MoFEMErrorCode setSkeletonAdjacency()
Definition: Simple.cpp:28
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:713
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:601

◆ setDim()

void MoFEM::Simple::setDim ( int  dim)

Set the problem dimension.

Returns
int

Definition at line 279 of file Simple.hpp.

279{ dIm = dim; };

◆ setFieldOrder()

MoFEMErrorCode MoFEM::Simple::setFieldOrder ( const std::string  field_name,
const int  order,
const Range *  ents = NULL 
)

◆ setParentAdjacency() [1/4]

template<>
MoFEMErrorCode MoFEM::Simple::setParentAdjacency ( )
private

Definition at line 190 of file Simple.cpp.

190 {
191 Interface &m_field = cOre;
193
195 boost::make_shared<ParentFiniteElementAdjacencyFunction<3>>(
198 boost::make_shared<ParentFiniteElementAdjacencyFunction<2>>(
200
201 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBTET,
203 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBHEX,
205 if (addBoundaryFE || !boundaryFields.empty()) {
206 CHKERR m_field.modify_finite_element_adjacency_table(
208 CHKERR m_field.modify_finite_element_adjacency_table(
210 }
211 if (addSkeletonFE || !skeletonFields.empty()) {
212 CHKERR m_field.modify_finite_element_adjacency_table(
214 CHKERR m_field.modify_finite_element_adjacency_table(
216 }
217
219}
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 2 > > parentAdjFunctionDim2
Definition: Simple.hpp:531
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
Definition: Simple.hpp:529

◆ setParentAdjacency() [2/4]

template<>
MoFEMErrorCode MoFEM::Simple::setParentAdjacency ( )
private

Definition at line 221 of file Simple.cpp.

221 {
222 Interface &m_field = cOre;
224
226 boost::make_shared<ParentFiniteElementAdjacencyFunction<2>>(
229 boost::make_shared<ParentFiniteElementAdjacencyFunction<1>>(
231
232 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBTRI,
234 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBQUAD,
236 if (addBoundaryFE || !boundaryFields.empty())
237 CHKERR m_field.modify_finite_element_adjacency_table(
239 if (addSkeletonFE || !skeletonFields.empty())
240 CHKERR m_field.modify_finite_element_adjacency_table(
242
244}
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 1 > > parentAdjFunctionDim1
Definition: Simple.hpp:533

◆ setParentAdjacency() [3/4]

template<>
MoFEMErrorCode MoFEM::Simple::setParentAdjacency ( )
private

Definition at line 246 of file Simple.cpp.

246 {
248 switch (getDim()) {
249 case 1:
250 THROW_MESSAGE("Not implemented");
251 case 2:
252 return setParentAdjacency<2>();
253 case 3:
254 return setParentAdjacency<3>();
255 default:
256 THROW_MESSAGE("Not implemented");
257 }
259}
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:574

◆ setParentAdjacency() [4/4]

template<int DIM>
MoFEMErrorCode MoFEM::Simple::setParentAdjacency
private

Definition at line 185 of file Simple.cpp.

185 {
186 static_assert(DIM == 2 || DIM == 3, "not implemented");
188}

◆ setSkeletonAdjacency() [1/5]

template<>
MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( )
private

Definition at line 33 of file Simple.cpp.

33 {
34 Interface &m_field = cOre;
36
37 auto defaultSkeletonEdge =
38 [&](moab::Interface &moab, const Field &field, const EntFiniteElement &fe,
39 std::vector<EntityHandle> &adjacency) -> MoFEMErrorCode {
41
42 CHKERR DefaultElementAdjacency::defaultEdge(moab, field, fe, adjacency);
43
44 if (std::find(domainFields.begin(), domainFields.end(), field.getName()) !=
45 domainFields.end()) {
46
47 const EntityHandle fe_ent = fe.getEnt();
48 std::vector<EntityHandle> bride_adjacency_edge;
49 CHKERR moab.get_adjacencies(&fe_ent, 1, 2, false, bride_adjacency_edge);
50
51 switch (field.getSpace()) {
52 case H1:
53 CHKERR moab.get_connectivity(&*bride_adjacency_edge.begin(),
54 bride_adjacency_edge.size(), adjacency,
55 true);
56 case HCURL:
57 case HDIV:
58 CHKERR moab.get_adjacencies(&*bride_adjacency_edge.begin(),
59 bride_adjacency_edge.size(), 1, false,
60 adjacency, moab::Interface::UNION);
61 case L2:
62 adjacency.insert(adjacency.end(), bride_adjacency_edge.begin(),
63 bride_adjacency_edge.end());
64 break;
65 case NOFIELD:
66 break;
67 default:
68 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
69 "this field is not implemented for TRI finite element");
70 }
71
72 std::sort(adjacency.begin(), adjacency.end());
73 auto it = std::unique(adjacency.begin(), adjacency.end());
74
75 std::vector<EntityHandle> new_adjacency(
76 std::distance(adjacency.begin(), it));
77 std::copy(adjacency.begin(), it, new_adjacency.begin());
78
79 for (auto e : new_adjacency) {
80 auto side_table = fe.getSideNumberTable();
81 if (side_table.find(e) == side_table.end())
82 const_cast<SideNumber_multiIndex &>(fe.getSideNumberTable())
83 .insert(
84 boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
85 }
86
87 adjacency.swap(new_adjacency);
88 }
89
91 };
92
93 CHKERR m_field.modify_finite_element_adjacency_table(skeletonFE, MBEDGE,
94 defaultSkeletonEdge);
95
97}
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, signed char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
static MoFEMErrorCode defaultEdge(Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)

◆ setSkeletonAdjacency() [2/5]

template<>
MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( )
private

Definition at line 99 of file Simple.cpp.

99 {
100 Interface &m_field = cOre;
102
103 auto defaultSkeletonEdge =
104 [&](moab::Interface &moab, const Field &field, const EntFiniteElement &fe,
105 std::vector<EntityHandle> &adjacency) -> MoFEMErrorCode {
107
108 CHKERR DefaultElementAdjacency::defaultFace(moab, field, fe, adjacency);
109
110 if (std::find(domainFields.begin(), domainFields.end(), field.getName()) !=
111 domainFields.end()) {
112
113 const EntityHandle fe_ent = fe.getEnt();
114 std::vector<EntityHandle> bride_adjacency_edge;
115 CHKERR moab.get_adjacencies(&fe_ent, 1, 2, false, bride_adjacency_edge);
116
117 switch (field.getSpace()) {
118 case H1:
119 CHKERR moab.get_connectivity(&*bride_adjacency_edge.begin(),
120 bride_adjacency_edge.size(), adjacency,
121 true);
122 case HCURL:
123 CHKERR moab.get_adjacencies(&*bride_adjacency_edge.begin(),
124 bride_adjacency_edge.size(), 1, false,
125 adjacency, moab::Interface::UNION);
126 case HDIV:
127 CHKERR moab.get_adjacencies(&*bride_adjacency_edge.begin(),
128 bride_adjacency_edge.size(), 2, false,
129 adjacency, moab::Interface::UNION);
130 case L2:
131 adjacency.insert(adjacency.end(), bride_adjacency_edge.begin(),
132 bride_adjacency_edge.end());
133 break;
134 case NOFIELD:
135 break;
136 default:
137 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
138 "this field is not implemented for TRI finite element");
139 }
140
141 std::sort(adjacency.begin(), adjacency.end());
142 auto it = std::unique(adjacency.begin(), adjacency.end());
143
144 std::vector<EntityHandle> new_adjacency(
145 std::distance(adjacency.begin(), it));
146 std::copy(adjacency.begin(), it, new_adjacency.begin());
147
148 for (auto e : new_adjacency) {
149 auto side_table = fe.getSideNumberTable();
150 if (side_table.find(e) == side_table.end())
151 const_cast<SideNumber_multiIndex &>(fe.getSideNumberTable())
152 .insert(
153 boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
154 }
155
156 adjacency.swap(new_adjacency);
157 }
158
160 };
161
162 CHKERR m_field.modify_finite_element_adjacency_table(skeletonFE, MBTRI,
163 defaultSkeletonEdge);
164 CHKERR m_field.modify_finite_element_adjacency_table(skeletonFE, MBQUAD,
165 defaultSkeletonEdge);
166
168}
static MoFEMErrorCode defaultFace(Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)

◆ setSkeletonAdjacency() [3/5]

template<>
MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( )
private

Definition at line 170 of file Simple.cpp.

170 {
172 switch (getDim()) {
173 case 1:
174 THROW_MESSAGE("Not implemented");
175 case 2:
176 return setSkeletonAdjacency<2>();
177 case 3:
178 return setSkeletonAdjacency<3>();
179 default:
180 THROW_MESSAGE("Not implemented");
181 }
183}

◆ setSkeletonAdjacency() [4/5]

template<int DIM>
MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency
private

Definition at line 28 of file Simple.cpp.

28 {
29 static_assert(DIM == 2 || DIM == 3, "not implemented");
31}

◆ setSkeletonAdjacency() [5/5]

MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( int  dim = -1)

Set the skeleton adjacency object.

Parameters
dim
Returns
MoFEMErrorCode

Definition at line 261 of file Simple.cpp.

261 {
263 if (dim == -1)
264 dim = getDim();
265 switch (dim) {
266 case 2:
267 return setSkeletonAdjacency<2>();
268 case 3:
269 return setSkeletonAdjacency<3>();
270 default:
271 SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_IMPLEMENTED, "Not implemented");
272 }
274}

◆ setUp()

MoFEMErrorCode MoFEM::Simple::setUp ( const PetscBool  is_partitioned = PETSC_TRUE)

Member Data Documentation

◆ addBoundaryFE

bool MoFEM::Simple::addBoundaryFE
private

Add boundary FE.

Definition at line 489 of file Simple.hpp.

◆ addParentAdjacencies

bool MoFEM::Simple::addParentAdjacencies
private

If set to true parent adjacencies are build.

Definition at line 490 of file Simple.hpp.

◆ addSkeletonFE

bool MoFEM::Simple::addSkeletonFE
private

Add skeleton FE.

Definition at line 488 of file Simple.hpp.

◆ bitAdjEnt

BitRefLevel MoFEM::Simple::bitAdjEnt
private

bit ref level for parent

Definition at line 494 of file Simple.hpp.

◆ bitAdjEntMask

BitRefLevel MoFEM::Simple::bitAdjEntMask
private

bit ref level for parent marent

Definition at line 495 of file Simple.hpp.

◆ bitAdjParent

BitRefLevel MoFEM::Simple::bitAdjParent
private

bit ref level for parent

Definition at line 492 of file Simple.hpp.

◆ bitAdjParentMask

BitRefLevel MoFEM::Simple::bitAdjParentMask
private

bit ref level for parent marent

Definition at line 493 of file Simple.hpp.

◆ bitLevel

BitRefLevel MoFEM::Simple::bitLevel
private

BitRefLevel of the probelm.

Definition at line 474 of file Simple.hpp.

◆ bitLevelMask

BitRefLevel MoFEM::Simple::bitLevelMask
private

BitRefLevel of the probelm.

Definition at line 475 of file Simple.hpp.

◆ boundaryFE

std::string MoFEM::Simple::boundaryFE
private

boundary finite element

Definition at line 509 of file Simple.hpp.

◆ boundaryFields

std::vector<std::string> MoFEM::Simple::boundaryFields
private

boundary fields

Definition at line 498 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 485 of file Simple.hpp.

◆ cOre

MoFEM::Core& MoFEM::Simple::cOre
private

Definition at line 472 of file Simple.hpp.

◆ dataFields

std::vector<std::string> MoFEM::Simple::dataFields
private

Data fields.

Definition at line 500 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 515 of file Simple.hpp.

◆ dM

SmartPetscObj<DM> MoFEM::Simple::dM
private

Discrete manager (interface to PETSc/MoFEM functions)

Definition at line 518 of file Simple.hpp.

◆ domainFE

std::string MoFEM::Simple::domainFE
private

domain finite element

Definition at line 508 of file Simple.hpp.

◆ domainFields

std::vector<std::string> MoFEM::Simple::domainFields
private

domain fields

Definition at line 497 of file Simple.hpp.

◆ fieldsOrder

std::multimap<std::string, std::pair<int, Range> > MoFEM::Simple::fieldsOrder
private

fields order

Definition at line 505 of file Simple.hpp.

◆ meshFileName

char MoFEM::Simple::meshFileName[255]
private

mesh file name

Definition at line 514 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 484 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 479 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 480 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 481 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 482 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 478 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleSetUP

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleSetUP
private

Definition at line 477 of file Simple.hpp.

◆ nameOfProblem

std::string MoFEM::Simple::nameOfProblem
private

problem name

Definition at line 507 of file Simple.hpp.

◆ noFieldDataFields

std::vector<std::string> MoFEM::Simple::noFieldDataFields
private

NOFIELD field name.

Definition at line 502 of file Simple.hpp.

◆ noFieldFields

std::vector<std::string> MoFEM::Simple::noFieldFields
private

NOFIELD field name.

Definition at line 501 of file Simple.hpp.

◆ otherFEs

std::vector<std::string> MoFEM::Simple::otherFEs
private

Other finite elements.

Definition at line 512 of file Simple.hpp.

◆ parentAdjFunctionDim1

boost::shared_ptr<ParentFiniteElementAdjacencyFunction<1> > MoFEM::Simple::parentAdjFunctionDim1
private

Definition at line 533 of file Simple.hpp.

◆ parentAdjFunctionDim2

boost::shared_ptr<ParentFiniteElementAdjacencyFunction<2> > MoFEM::Simple::parentAdjFunctionDim2
private

Definition at line 531 of file Simple.hpp.

◆ parentAdjFunctionDim3

boost::shared_ptr<ParentFiniteElementAdjacencyFunction<3> > MoFEM::Simple::parentAdjFunctionDim3
private

Definition at line 529 of file Simple.hpp.

◆ skeletonFE

std::string MoFEM::Simple::skeletonFE
private

skeleton finite element

Definition at line 510 of file Simple.hpp.

◆ skeletonFields

std::vector<std::string> MoFEM::Simple::skeletonFields
private

fields on the skeleton

Definition at line 499 of file Simple.hpp.

◆ skeletonMeshset

EntityHandle MoFEM::Simple::skeletonMeshset
private

skeleton meshset with boundary

Definition at line 486 of file Simple.hpp.


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