v0.13.1
Loading...
Searching...
No Matches
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, std::string fe_name="")
 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...
 
EntityHandlegetMeshSet ()
 Get the MeshSet object. More...
 
EntityHandlegetBoundaryMeshSet ()
 Get the BoundaryMeshSet object. More...
 
EntityHandlegetSkeletonMeshSet ()
 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 (std::string fe_name)
 
template<int DIM = -1>
MoFEMErrorCode setParentAdjacency ()
 
template<>
MoFEMErrorCode setSkeletonAdjacency (std::string fe_name)
 
template<>
MoFEMErrorCode setSkeletonAdjacency (std::string fe_name)
 
template<>
MoFEMErrorCode setSkeletonAdjacency (std::string fe_name)
 
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::list< std::tuple< std::string, 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_elastic.cpp, and wave_equation.cpp.

Definition at line 26 of file Simple.hpp.

Constructor & Destructor Documentation

◆ Simple()

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

Definition at line 271 of file Simple.cpp.

272 : cOre(const_cast<Core &>(core)), bitLevel(BitRefLevel().set(0)),
274 skeletonMeshset(0), nameOfProblem("SimpleProblem"), domainFE("dFE"),
275 boundaryFE("bFE"), skeletonFE("sFE"), dIm(-1), addSkeletonFE(false),
276 addBoundaryFE(false), addParentAdjacencies(false),
279 PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleSetUP);
280 PetscLogEventRegister("SimpleLoadMesh", 0, &MOFEM_EVENT_SimpleLoadMesh);
281 PetscLogEventRegister("SimpleBuildFields", 0, &MOFEM_EVENT_SimpleBuildFields);
282 PetscLogEventRegister("SimpleBuildFiniteElements", 0,
284 PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleBuildProblem);
285 PetscLogEventRegister("SimpleKSPSolve", 0, &MOFEM_EVENT_SimpleKSPSolve);
286 strcpy(meshFileName, "mesh.h5m");
287}
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
CoreTmp< 0 > Core
Definition: Core.hpp:1086
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:477
MoFEM::Core & cOre
Definition: Simple.hpp:464
BitRefLevel bitAdjEnt
bit ref level for parent
Definition: Simple.hpp:486
char meshFileName[255]
mesh file name
Definition: Simple.hpp:506
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:476
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:471
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:470
bool addSkeletonFE
Add skeleton FE.
Definition: Simple.hpp:480
BitRefLevel bitAdjParent
bit ref level for parent
Definition: Simple.hpp:484
bool addBoundaryFE
Add boundary FE.
Definition: Simple.hpp:481
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:502
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:501
std::string nameOfProblem
problem name
Definition: Simple.hpp:499
std::string domainFE
domain finite element
Definition: Simple.hpp:500
BitRefLevel bitAdjEntMask
bit ref level for parent marent
Definition: Simple.hpp:487
BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:466
int dIm
dimension of problem
Definition: Simple.hpp:507
PetscLogEvent MOFEM_EVENT_SimpleSetUP
Definition: Simple.hpp:469
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:474
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:472
BitRefLevel bitLevelMask
BitRefLevel of the probelm.
Definition: Simple.hpp:467
BitRefLevel bitAdjParentMask
bit ref level for parent marent
Definition: Simple.hpp:485
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:473
EntityHandle skeletonMeshset
skeleton meshset with boundary
Definition: Simple.hpp:478
bool addParentAdjacencies
If set to true parent adjacencies are build.
Definition: Simple.hpp:482

◆ ~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 391 of file Simple.cpp.

395 {
396 Interface &m_field = cOre;
398 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
399 verb);
400 boundaryFields.push_back(name);
401 if (space == NOFIELD)
402 SETERRQ(
403 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
404 "NOFIELD space for boundary filed not implemented in Simple interface");
406}
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:84
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#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
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1955
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:490

◆ 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 429 of file Simple.cpp.

433 {
434
435 Interface &m_field = cOre;
437 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
438 verb);
439 if (space == NOFIELD)
440 noFieldDataFields.push_back(name);
441 else
442 dataFields.push_back(name);
444}
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:492
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:494

◆ 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, 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 373 of file Simple.cpp.

377 {
378
379 Interface &m_field = cOre;
381 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
382 verb);
383 if (space == NOFIELD)
384 noFieldFields.push_back(name);
385 else
386 domainFields.push_back(name);
388}
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:489
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:493

◆ 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 815 of file Simple.cpp.

816 {
817 Interface &m_field = cOre;
819 CHKERR m_field.getInterface<ProblemsManager>()->addFieldToEmptyFieldBlocks(
820 getProblemName(), row_field, col_field);
822}
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
add empty block to problem
Definition: Simple.cpp:815
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:334

◆ 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 409 of file Simple.cpp.

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

◆ buildFields()

MoFEMErrorCode MoFEM::Simple::buildFields ( )

Build fields.

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

Definition at line 596 of file Simple.cpp.

596 {
597 Interface &m_field = cOre;
599 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
600
601 auto comm_interface_ptr = m_field.getInterface<CommInterface>();
602 auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
603
604 // Add entities to the fields
605 auto add_ents_to_field = [&](auto meshset, auto dim, auto &fields) {
607 for (auto &field : fields) {
608 CHKERR m_field.add_ents_to_field_by_dim(meshset, dim, field);
609 CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
610 }
612 };
613
614 auto make_no_field_ents = [&](auto &fields) {
616 for (auto &field : fields) {
617 std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
618 CHKERR m_field.create_vertices_and_add_to_field(field, coords.data(), 1);
619 CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
620 CHKERR bit_ref_ptr->setFieldEntitiesBitRefLevel(field, bitLevel);
621 }
623 };
624
625 CHKERR add_ents_to_field(meshSet, dIm, domainFields);
626 CHKERR add_ents_to_field(meshSet, dIm, dataFields);
627 CHKERR add_ents_to_field(boundaryMeshset, dIm - 1, boundaryFields);
628 CHKERR add_ents_to_field(skeletonMeshset, dIm - 1, skeletonFields);
629
630 std::set<std::string> nofield_fields;
631 for (auto &f : noFieldFields)
632 nofield_fields.insert(f);
633 for (auto &f : noFieldDataFields)
634 nofield_fields.insert(f);
635
636 CHKERR make_no_field_ents(nofield_fields);
637
638 // Set order
639
640 auto get_field_ptr = [&](auto &f) { return m_field.get_field_structure(f); };
641 for (auto &t : fieldsOrder) {
642 const auto f = std::get<0>(t);
643 const auto order = std::get<1>(t);
644
645 MOFEM_TAG_AND_LOG("WORLD", Sev::inform, "Simple")
646 << "Set order to field " << f << " order " << order;
647 if (!std::get<2>(t).empty()) {
648 MOFEM_LOG_CHANNEL("SYNC");
649 MOFEM_TAG_AND_LOG("SYNC", Sev::verbose, "Simple")
650 << "To ents: " << std::endl
651 << std::get<2>(t) << std::endl;
652 MOFEM_LOG_SYNCHRONISE(m_field.get_comm());
653 }
654
655 if (std::get<2>(t).empty()) {
656 auto f_ptr = get_field_ptr(f);
657
658 if (f_ptr->getSpace() == H1) {
659 if (f_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
660 CHKERR m_field.set_field_order(meshSet, MBVERTEX, f, order);
661 } else {
662 CHKERR m_field.set_field_order(meshSet, MBVERTEX, f, 1);
663 }
664 }
665
666 for (auto d = 1; d <= dIm; ++d) {
667 for (EntityType t = CN::TypeDimensionMap[d].first;
668 t <= CN::TypeDimensionMap[d].second; ++t) {
669 CHKERR m_field.set_field_order(meshSet, t, f, order);
670 }
671 }
672 } else {
673 CHKERR m_field.set_field_order(std::get<2>(t), f, order);
674 }
675 }
676 MOFEM_LOG_CHANNEL("WORLD");
677 // Build fields
678 CHKERR m_field.build_fields();
679 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
681}
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:338
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:345
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:64
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ H1
continuous field
Definition: definitions.h:85
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
const int dim
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:277
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
auto f
Definition: HenckyOps.hpp:5
constexpr double t
plate stiffness
Definition: plate.cpp:59
std::list< std::tuple< std::string, int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:497

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

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

Definition at line 683 of file Simple.cpp.

683 {
684 Interface &m_field = cOre;
686 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
687 // Add finite elements
688 CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
689 true);
690 CHKERR m_field.build_finite_elements(domainFE);
691 if (addBoundaryFE || !boundaryFields.empty()) {
692 CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
693 boundaryFE, true);
694 CHKERR m_field.build_finite_elements(boundaryFE);
695 }
696 if (addSkeletonFE || !skeletonFields.empty()) {
697 CHKERR m_field.add_ents_to_finite_element_by_dim(skeletonMeshset, dIm - 1,
698 skeletonFE, true);
699 CHKERR m_field.build_finite_elements(skeletonFE);
700 }
701 for (std::vector<std::string>::iterator fit = otherFEs.begin();
702 fit != otherFEs.end(); ++fit) {
703 CHKERR m_field.build_finite_elements(*fit);
704 }
705 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
707}
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:504

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

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

Definition at line 709 of file Simple.cpp.

709 {
710 Interface &m_field = cOre;
712 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
713 CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
714 // Set problem by the DOFs on the fields rather that by adding DOFs on the
715 // elements
716 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
718 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
719 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
721}
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1228
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:510

◆ createBoundaryMeshset()

MoFEMErrorCode MoFEM::Simple::createBoundaryMeshset ( )
private

Definition at line 824 of file Simple.cpp.

824 {
825 Interface &m_field = cOre;
827 ParallelComm *pcomm =
828 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
829
830 auto get_skin = [&](auto meshset) {
831 // filter not owned entities, those are not on boundary
832
833 Range domain_ents;
834 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
835 domain_ents, true);
836 CHKERR pcomm->filter_pstatus(domain_ents,
837 PSTATUS_SHARED | PSTATUS_MULTISHARED,
838 PSTATUS_NOT, -1, nullptr);
839
840 Skinner skin(&m_field.get_moab());
841 Range domain_skin;
842 CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
843 CHKERR pcomm->filter_pstatus(domain_skin,
844 PSTATUS_SHARED | PSTATUS_MULTISHARED,
845 PSTATUS_NOT, -1, nullptr);
846 return domain_skin;
847 };
848
849 auto create_boundary_meshset = [&](auto &&domain_skin) {
851 // create boundary meshset
852 if (boundaryMeshset != 0) {
854 }
855 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
856 CHKERR m_field.get_moab().add_entities(boundaryMeshset, domain_skin);
857 for (int dd = 0; dd != dIm - 1; dd++) {
858 Range adj;
859 CHKERR m_field.get_moab().get_adjacencies(domain_skin, dd, false, adj,
860 moab::Interface::UNION);
861 CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
862 }
864 };
865
866 CHKERR create_boundary_meshset(get_skin(meshSet));
867
869}
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:215
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

◆ createSkeletonMeshset()

MoFEMErrorCode MoFEM::Simple::createSkeletonMeshset ( )
private

Definition at line 871 of file Simple.cpp.

871 {
872 Interface &m_field = cOre;
874
875 ParallelComm *pcomm =
876 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
877
878 auto create_skeleton_meshset = [&](auto meshset) {
880 // create boundary meshset
881 if (skeletonMeshset != 0) {
883 }
884 Range boundary_ents, skeleton_ents;
885 CHKERR m_field.get_moab().get_entities_by_dimension(boundaryMeshset,
886 dIm - 1, boundary_ents);
887 Range domain_ents;
888 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
889 domain_ents, true);
890 CHKERR m_field.get_moab().get_adjacencies(
891 domain_ents, dIm - 1, false, skeleton_ents, moab::Interface::UNION);
892 skeleton_ents = subtract(skeleton_ents, boundary_ents);
893 CHKERR pcomm->filter_pstatus(skeleton_ents, PSTATUS_NOT_OWNED, PSTATUS_NOT,
894 -1, nullptr);
895 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, skeletonMeshset);
896 CHKERR m_field.get_moab().add_entities(skeletonMeshset, skeleton_ents);
898 };
899
900 CHKERR create_skeleton_meshset(meshSet);
901
903}

◆ defineFiniteElements()

MoFEMErrorCode MoFEM::Simple::defineFiniteElements ( )

Define finite elements.

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

Definition at line 492 of file Simple.cpp.

492 {
493 Interface &m_field = cOre;
495
496 auto clear_rows_and_cols = [&](auto &fe_name) {
498 auto fe_ptr = m_field.get_finite_elements();
499 auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
500 ->get<FiniteElement_name_mi_tag>();
501 auto it_fe = fe_by_name.find(fe_name);
502 if (it_fe != fe_by_name.end()) {
503
504 if (!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
505 SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
506 "modification unsuccessful");
507
508 if (!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
509 SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
510 "modification unsuccessful");
511 }
513 };
514 CHKERR clear_rows_and_cols(domainFE);
515 CHKERR clear_rows_and_cols(boundaryFE);
516 CHKERR clear_rows_and_cols(skeletonFE);
517
518 // Define finite elements
519 CHKERR m_field.add_finite_element(domainFE, MF_ZERO);
520
521 auto add_fields = [&](auto &fe_name, auto &fields) {
523 for (auto &field : fields) {
524 CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
525 CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
526 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
527 }
529 };
530
531 auto add_data_fields = [&](auto &fe_name, auto &fields) {
533 for (auto &field : fields)
534 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
536 };
537
538 CHKERR add_fields(domainFE, domainFields);
539 CHKERR add_data_fields(domainFE, dataFields);
540 CHKERR add_fields(domainFE, noFieldFields);
541 CHKERR add_data_fields(domainFE, noFieldDataFields);
542
543 if (addBoundaryFE || !boundaryFields.empty()) {
544 CHKERR m_field.add_finite_element(boundaryFE, MF_ZERO);
545 CHKERR add_fields(boundaryFE, domainFields);
546 if (!boundaryFields.empty())
547 CHKERR add_fields(boundaryFE, boundaryFields);
548 CHKERR add_data_fields(boundaryFE, dataFields);
549 CHKERR add_data_fields(boundaryFE, noFieldDataFields);
550 CHKERR add_fields(boundaryFE, noFieldFields);
551 }
552 if (addSkeletonFE || !skeletonFields.empty()) {
553 CHKERR m_field.add_finite_element(skeletonFE, MF_ZERO);
554 CHKERR add_fields(skeletonFE, domainFields);
555 if (!skeletonFields.empty())
556 CHKERR add_fields(skeletonFE, skeletonFields);
557 CHKERR add_data_fields(skeletonFE, dataFields);
558 CHKERR add_data_fields(skeletonFE, noFieldDataFields);
559 CHKERR add_fields(skeletonFE, noFieldFields);
560 }
562}
@ MF_ZERO
Definition: definitions.h:98
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition: definitions.h:34
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 564 of file Simple.cpp.

564 {
565 Interface &m_field = cOre;
567 // Create dm instance
568 dM = createSmartDM(m_field.get_comm(), "DMMOFEM");
569 // set dm data structure which created mofem data structures
572 CHKERR DMSetFromOptions(dM);
574 if (addBoundaryFE || !boundaryFields.empty()) {
576 }
577 if (addSkeletonFE || !skeletonFields.empty()) {
579 }
580 for (std::vector<std::string>::iterator fit = otherFEs.begin();
581 fit != otherFEs.end(); ++fit) {
582 CHKERR DMMoFEMAddElement(dM, fit->c_str());
583 }
584 CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
586}
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:1070
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:118
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:450
auto createSmartDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.

◆ deleteDM()

MoFEMErrorCode MoFEM::Simple::deleteDM ( )

Delete dm.

Returns
MoFEMErrorCode
Returns
MoFEMErrorCode

Definition at line 792 of file Simple.cpp.

792 {
794 dM.reset();
796}

◆ deleteFiniteElements()

MoFEMErrorCode MoFEM::Simple::deleteFiniteElements ( )

Delete finite elements.

Returns
MoFEMErrorCode
Returns
MoFEMErrorCode

Definition at line 803 of file Simple.cpp.

803 {
804 Interface &m_field = cOre;
806 for (auto fe : {domainFE, boundaryFE, skeletonFE}) {
807 if (m_field.check_finite_element(fe)) {
808 CHKERR m_field.delete_finite_element(fe);
809 }
810 }
812}

◆ exchangeGhostCells()

MoFEMErrorCode MoFEM::Simple::exchangeGhostCells ( )
private

get all adjacent ghosted entities

Definition at line 905 of file Simple.cpp.

905 {
906 Interface &m_field = cOre;
908 MOFEM_LOG("WORLD", Sev::verbose) << "Exchange ghost cells";
909
910 ParallelComm *pcomm =
911 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
912 if (pcomm == NULL)
913 pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
914
915 Range verts;
916 CHKERR m_field.get_moab().get_entities_by_type(0, MBVERTEX, verts);
917 CHKERR pcomm->exchange_ghost_cells(getDim(), getDim() - 1, 1,
918 3 /**get all adjacent ghosted entities */,
919 true, false, meshSet ? &meshSet : nullptr);
920
921 Range shared;
922 CHKERR m_field.get_moab().get_entities_by_dimension(0, dIm, shared);
923 for (auto d = dIm - 1; d >= 1; --d) {
924 CHKERR m_field.get_moab().get_adjacencies(shared, d, false, shared,
925 moab::Interface::UNION);
926 }
927 CHKERR pcomm->filter_pstatus(shared, PSTATUS_SHARED | PSTATUS_MULTISHARED,
928 PSTATUS_OR, -1, &shared);
929 Tag part_tag = pcomm->part_tag();
930 CHKERR pcomm->exchange_tags(part_tag, shared);
931
933}
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
int getDim() const
Get the problem dimension.
Definition: Simple.hpp:264

◆ getAddBoundaryFE()

bool & MoFEM::Simple::getAddBoundaryFE ( )
inline

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 408 of file Simple.hpp.

408{ return addBoundaryFE; }

◆ getAddSkeletonFE()

bool & MoFEM::Simple::getAddSkeletonFE ( )
inline

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 397 of file Simple.hpp.

397{ return addSkeletonFE; }

◆ getBitAdjEnt()

auto & MoFEM::Simple::getBitAdjEnt ( )
inline

bit ref level for parent

Returns
auto&

Definition at line 439 of file Simple.hpp.

439{ return bitAdjEnt; }

◆ getBitAdjEntMask()

auto & MoFEM::Simple::getBitAdjEntMask ( )
inline

bit ref level for parent marent

Returns
auto&

Definition at line 446 of file Simple.hpp.

446{ return bitAdjEntMask; }

◆ getBitAdjParent()

auto & MoFEM::Simple::getBitAdjParent ( )
inline

bit ref level for parent

Returns
auto&

Definition at line 425 of file Simple.hpp.

425{ return bitAdjParent; }

◆ getBitAdjParentMask()

auto & MoFEM::Simple::getBitAdjParentMask ( )
inline

bit ref level for parent marent

Returns
auto&

Definition at line 432 of file Simple.hpp.

432{ return bitAdjParentMask; }

◆ getBitRefLevel()

BitRefLevel & MoFEM::Simple::getBitRefLevel ( )
inline

Get the BitRefLevel.

Returns
BitRefLevel

Definition at line 299 of file Simple.hpp.

299{ return bitLevel; }

◆ getBitRefLevelMask()

BitRefLevel & MoFEM::Simple::getBitRefLevelMask ( )
inline

Get the BitRefLevel.

Returns
BitRefLevel

Definition at line 306 of file Simple.hpp.

306{ return bitLevelMask; }

◆ getBoundaryFEName() [1/2]

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

Get the Boundary FE Name.

Returns
std::string&

Definition at line 348 of file Simple.hpp.

348{ return boundaryFE; }

◆ getBoundaryFEName() [2/2]

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

Get the Boundary FE Name.

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

Definition at line 320 of file Simple.hpp.

320{ return boundaryFE; }

◆ getBoundaryMeshSet()

EntityHandle & MoFEM::Simple::getBoundaryMeshSet ( )
inline

Get the BoundaryMeshSet object.

Returns
EntityHandle&

Definition at line 285 of file Simple.hpp.

285{ return boundaryMeshset; }

◆ getDim()

int MoFEM::Simple::getDim ( ) const
inline

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 264 of file Simple.hpp.

264{ return dIm; }

◆ getDM() [1/2]

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

Return smart DM object.

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

Definition at line 254 of file Simple.hpp.

254{ return dM; }

◆ getDM() [2/2]

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

◆ getDomainFEName() [1/2]

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

Get the Domain FE Name.

Returns
std::string&

Definition at line 341 of file Simple.hpp.

341{ return domainFE; }

◆ getDomainFEName() [2/2]

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

◆ getMeshSet()

EntityHandle & MoFEM::Simple::getMeshSet ( )
inline

Get the MeshSet object.

Returns
EntityHandle&

Definition at line 278 of file Simple.hpp.

278{ 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, 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 289 of file Simple.cpp.

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

◆ getOtherFiniteElements()

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

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 372 of file Simple.hpp.

372{ return otherFEs; }

◆ getParentAdjacencies()

bool & MoFEM::Simple::getParentAdjacencies ( )
inline

Get the addParentAdjacencies.

If set true add parent adjacencies

Returns
true
false

Definition at line 418 of file Simple.hpp.

418{ return addParentAdjacencies; }

◆ getProblemName() [1/2]

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

Get the Problem Name.

Returns
std::string&

Definition at line 362 of file Simple.hpp.

362{ return nameOfProblem; }

◆ getProblemName() [2/2]

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

Get the Problem Name.

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

Definition at line 334 of file Simple.hpp.

334{ return nameOfProblem; }

◆ getSkeletonFEName() [1/2]

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

Get the Skeleton FE Name.

Returns
std::string&

Definition at line 355 of file Simple.hpp.

355{ return skeletonFE; }

◆ getSkeletonFEName() [2/2]

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

◆ getSkeletonMeshSet()

EntityHandle & MoFEM::Simple::getSkeletonMeshSet ( )
inline

Get the SkeletonMeshSet object.

Returns
EntityHandle&

Definition at line 292 of file Simple.hpp.

292{ 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 359 of file Simple.cpp.

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

◆ 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, 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 303 of file Simple.cpp.

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

◆ query_interface()

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

Implements MoFEM::UnknownInterface.

Definition at line 8 of file Simple.cpp.

9 {
10 *iface = const_cast<Simple *>(this);
11 return 0;
12}
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:271

◆ removeBoundaryField()

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

Remove field form boundary.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 462 of file Simple.cpp.

462 {
463 Interface &m_field = cOre;
465
466 auto remove_field_from_list = [&](auto &vec) {
467 auto it = std::find(vec.begin(), vec.end(), name);
468 if (it != vec.end())
469 vec.erase(it);
470 };
471
472 remove_field_from_list(boundaryFields);
473
475}

◆ removeDomainField()

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

Remove field form domain.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 446 of file Simple.cpp.

446 {
447 Interface &m_field = cOre;
449
450 auto remove_field_from_list = [&](auto &vec) {
451 auto it = std::find(vec.begin(), vec.end(), name);
452 if (it != vec.end())
453 vec.erase(it);
454 };
455
456 remove_field_from_list(noFieldFields);
457 remove_field_from_list(domainFields);
458
460}

◆ removeSkeletonField()

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

Remove field form skeleton.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 477 of file Simple.cpp.

477 {
478 Interface &m_field = cOre;
480
481 auto remove_field_from_list = [&](auto &vec) {
482 auto it = std::find(vec.begin(), vec.end(), name);
483 if (it != vec.end())
484 vec.erase(it);
485 };
486
487 remove_field_from_list(skeletonFields);
488
490}

◆ 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 746 of file Simple.cpp.

746 {
747 Interface &m_field = cOre;
749 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
750
751 if (!only_dm) {
754 if (addSkeletonFE || !skeletonFields.empty())
759 }
760
761 CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
762
763 const Problem *problem_ptr;
764 CHKERR DMMoFEMGetProblemPtr(dM, &problem_ptr);
765 const auto problem_name = problem_ptr->getName();
766 CHKERR m_field.modify_problem_ref_level_set_bit(problem_name, bitLevel);
767 CHKERR m_field.modify_problem_mask_ref_level_set_bit(problem_name,
769
770 // Set problem by the DOFs on the fields rather that by adding DOFs on the
771 // elements
772 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
774 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
775
776 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
778}
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:373
MoFEMErrorCode setParentAdjacency()
Definition: Simple.cpp:176
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:492
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:683
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:596
MoFEMErrorCode setSkeletonAdjacency(int dim=-1, std::string fe_name="")
Set the skeleton adjacency object.
Definition: Simple.cpp:252

◆ setDim()

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

Set the problem dimension.

Returns
int

Definition at line 271 of file Simple.hpp.

271{ 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 181 of file Simple.cpp.

181 {
182 Interface &m_field = cOre;
184
186 boost::make_shared<ParentFiniteElementAdjacencyFunction<3>>(
189 boost::make_shared<ParentFiniteElementAdjacencyFunction<2>>(
191
192 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBTET,
194 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBHEX,
196 if (addBoundaryFE || !boundaryFields.empty()) {
197 CHKERR m_field.modify_finite_element_adjacency_table(
199 CHKERR m_field.modify_finite_element_adjacency_table(
201 }
202 if (addSkeletonFE || !skeletonFields.empty()) {
203 CHKERR m_field.modify_finite_element_adjacency_table(
205 CHKERR m_field.modify_finite_element_adjacency_table(
207 }
208
210}
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 2 > > parentAdjFunctionDim2
Definition: Simple.hpp:524
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
Definition: Simple.hpp:522

◆ setParentAdjacency() [2/4]

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

Definition at line 212 of file Simple.cpp.

212 {
213 Interface &m_field = cOre;
215
217 boost::make_shared<ParentFiniteElementAdjacencyFunction<2>>(
220 boost::make_shared<ParentFiniteElementAdjacencyFunction<1>>(
222
223 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBTRI,
225 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBQUAD,
227 if (addBoundaryFE || !boundaryFields.empty())
228 CHKERR m_field.modify_finite_element_adjacency_table(
230 if (addSkeletonFE || !skeletonFields.empty())
231 CHKERR m_field.modify_finite_element_adjacency_table(
233
235}
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 1 > > parentAdjFunctionDim1
Definition: Simple.hpp:526

◆ setParentAdjacency() [3/4]

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

Definition at line 237 of file Simple.cpp.

237 {
239 switch (getDim()) {
240 case 1:
241 THROW_MESSAGE("Not implemented");
242 case 2:
243 return setParentAdjacency<2>();
244 case 3:
245 return setParentAdjacency<3>();
246 default:
247 THROW_MESSAGE("Not implemented");
248 }
250}
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561

◆ setParentAdjacency() [4/4]

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

Definition at line 176 of file Simple.cpp.

176 {
177 static_assert(DIM == 2 || DIM == 3, "not implemented");
179}

◆ setSkeletonAdjacency() [1/5]

MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( int  dim = -1,
std::string  fe_name = "" 
)

Set the skeleton adjacency object.

Parameters
dim
fe_name
Returns
MoFEMErrorCode

Definition at line 252 of file Simple.cpp.

252 {
254 if (dim == -1)
255 dim = getDim();
256
257 if (fe_name.empty())
258 fe_name = skeletonFE;
259
260 switch (dim) {
261 case 2:
262 return setSkeletonAdjacency<2>(fe_name);
263 case 3:
264 return setSkeletonAdjacency<3>(fe_name);
265 default:
266 SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_IMPLEMENTED, "Not implemented");
267 }
269}

◆ setSkeletonAdjacency() [2/5]

template<>
MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( std::string  fe_name)
private

Definition at line 21 of file Simple.cpp.

21 {
22 Interface &m_field = cOre;
24
25 auto defaultSkeletonEdge =
26 [&](moab::Interface &moab, const Field &field, const EntFiniteElement &fe,
27 std::vector<EntityHandle> &adjacency) -> MoFEMErrorCode {
29
30 CHKERR DefaultElementAdjacency::defaultEdge(moab, field, fe, adjacency);
31
32 if (std::find(domainFields.begin(), domainFields.end(), field.getName()) !=
33 domainFields.end()) {
34
35 const EntityHandle fe_ent = fe.getEnt();
36 std::vector<EntityHandle> bride_adjacency_edge;
37 CHKERR moab.get_adjacencies(&fe_ent, 1, 2, false, bride_adjacency_edge);
38
39 switch (field.getSpace()) {
40 case H1:
41 CHKERR moab.get_connectivity(&*bride_adjacency_edge.begin(),
42 bride_adjacency_edge.size(), adjacency,
43 true);
44 case HCURL:
45 case HDIV:
46 CHKERR moab.get_adjacencies(&*bride_adjacency_edge.begin(),
47 bride_adjacency_edge.size(), 1, false,
48 adjacency, moab::Interface::UNION);
49 case L2:
50 adjacency.insert(adjacency.end(), bride_adjacency_edge.begin(),
51 bride_adjacency_edge.end());
52 break;
53 case NOFIELD:
54 break;
55 default:
56 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
57 "this field is not implemented for TRI finite element");
58 }
59
60 std::sort(adjacency.begin(), adjacency.end());
61 auto it = std::unique(adjacency.begin(), adjacency.end());
62
63 std::vector<EntityHandle> new_adjacency(
64 std::distance(adjacency.begin(), it));
65 std::copy(adjacency.begin(), it, new_adjacency.begin());
66
67 for (auto e : new_adjacency) {
68 auto side_table = fe.getSideNumberTable();
69 if (side_table.find(e) == side_table.end())
70 const_cast<SideNumber_multiIndex &>(fe.getSideNumberTable())
71 .insert(
72 boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
73 }
74
75 adjacency.swap(new_adjacency);
76 }
77
79 };
80
81 CHKERR m_field.modify_finite_element_adjacency_table(fe_name, MBEDGE,
82 defaultSkeletonEdge);
83
85}
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ HCURL
field with continuous tangents
Definition: definitions.h:86
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
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:56
static MoFEMErrorCode defaultEdge(Interface &moab, const Field &field, const EntFiniteElement &fe, std::vector< EntityHandle > &adjacency)

◆ setSkeletonAdjacency() [3/5]

template<>
MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( std::string  fe_name)
private

Definition at line 88 of file Simple.cpp.

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

◆ setSkeletonAdjacency() [4/5]

template<>
MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( std::string  fe_name)
private

Definition at line 160 of file Simple.cpp.

160 {
162
163 switch (getDim()) {
164 case 1:
165 THROW_MESSAGE("Not implemented");
166 case 2:
167 return setSkeletonAdjacency<2>(fe_name);
168 case 3:
169 return setSkeletonAdjacency<3>(fe_name);
170 default:
171 THROW_MESSAGE("Not implemented");
172 }
174}

◆ setSkeletonAdjacency() [5/5]

template<int DIM>
MoFEMErrorCode MoFEM::Simple::setSkeletonAdjacency ( std::string  fe_name)
private

Definition at line 15 of file Simple.cpp.

15 {
16 static_assert(DIM == 2 || DIM == 3, "not implemented");
18}

◆ 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 481 of file Simple.hpp.

◆ addParentAdjacencies

bool MoFEM::Simple::addParentAdjacencies
private

If set to true parent adjacencies are build.

Definition at line 482 of file Simple.hpp.

◆ addSkeletonFE

bool MoFEM::Simple::addSkeletonFE
private

Add skeleton FE.

Definition at line 480 of file Simple.hpp.

◆ bitAdjEnt

BitRefLevel MoFEM::Simple::bitAdjEnt
private

bit ref level for parent

Definition at line 486 of file Simple.hpp.

◆ bitAdjEntMask

BitRefLevel MoFEM::Simple::bitAdjEntMask
private

bit ref level for parent marent

Definition at line 487 of file Simple.hpp.

◆ bitAdjParent

BitRefLevel MoFEM::Simple::bitAdjParent
private

bit ref level for parent

Definition at line 484 of file Simple.hpp.

◆ bitAdjParentMask

BitRefLevel MoFEM::Simple::bitAdjParentMask
private

bit ref level for parent marent

Definition at line 485 of file Simple.hpp.

◆ bitLevel

BitRefLevel MoFEM::Simple::bitLevel
private

BitRefLevel of the probelm.

Definition at line 466 of file Simple.hpp.

◆ bitLevelMask

BitRefLevel MoFEM::Simple::bitLevelMask
private

BitRefLevel of the probelm.

Definition at line 467 of file Simple.hpp.

◆ boundaryFE

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

boundary finite element

Definition at line 501 of file Simple.hpp.

◆ boundaryFields

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

boundary fields

Definition at line 490 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 477 of file Simple.hpp.

◆ cOre

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

Definition at line 464 of file Simple.hpp.

◆ dataFields

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

Data fields.

Definition at line 492 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 507 of file Simple.hpp.

◆ dM

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

Discrete manager (interface to PETSc/MoFEM functions)

Definition at line 510 of file Simple.hpp.

◆ domainFE

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

domain finite element

Definition at line 500 of file Simple.hpp.

◆ domainFields

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

domain fields

Definition at line 489 of file Simple.hpp.

◆ fieldsOrder

std::list<std::tuple<std::string, int, Range> > MoFEM::Simple::fieldsOrder
private

fields order

Definition at line 497 of file Simple.hpp.

◆ meshFileName

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

mesh file name

Definition at line 506 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 476 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 471 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 472 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 473 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 474 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 470 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleSetUP

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleSetUP
private

Definition at line 469 of file Simple.hpp.

◆ nameOfProblem

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

problem name

Definition at line 499 of file Simple.hpp.

◆ noFieldDataFields

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

NOFIELD field name.

Definition at line 494 of file Simple.hpp.

◆ noFieldFields

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

NOFIELD field name.

Definition at line 493 of file Simple.hpp.

◆ otherFEs

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

Other finite elements.

Definition at line 504 of file Simple.hpp.

◆ parentAdjFunctionDim1

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

Definition at line 526 of file Simple.hpp.

◆ parentAdjFunctionDim2

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

Definition at line 524 of file Simple.hpp.

◆ parentAdjFunctionDim3

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

Definition at line 522 of file Simple.hpp.

◆ skeletonFE

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

skeleton finite element

Definition at line 502 of file Simple.hpp.

◆ skeletonFields

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

fields on the skeleton

Definition at line 491 of file Simple.hpp.

◆ skeletonMeshset

EntityHandle MoFEM::Simple::skeletonMeshset
private

skeleton meshset with boundary

Definition at line 478 of file Simple.hpp.


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