v0.13.2
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...
 
DEPRECATED EntityHandlegetMeshSet ()
 
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 parent More...
 
auto & getBitAdjEnt ()
 bit ref level for parent More...
 
auto & getBitAdjEntMask ()
 bit ref level for parent parent 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 problem. More...
 
BitRefLevel bitLevelMask
 BitRefLevel of the problem. 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 parent More...
 
BitRefLevel bitAdjEnt
 bit ref level for parent More...
 
BitRefLevel bitAdjEntMask
 bit ref level for parent parent 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, bool > > fieldsOrder
 
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
 
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 2 > > parentAdjSkeletonFunctionDim2
 
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 1 > > parentAdjSkeletonFunctionDim1
 

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, 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, level_set.cpp, loop_entities.cpp, mesh_smoothing.cpp, mixed_poisson.cpp, nonlinear_elastic.cpp, operators_tests.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, poisson_3d_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 27 of file Simple.hpp.

Constructor & Destructor Documentation

◆ Simple()

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

Definition at line 162 of file Simple.cpp.

163 : cOre(const_cast<Core &>(core)), bitLevel(BitRefLevel().set(0)),
165 skeletonMeshset(0), nameOfProblem("SimpleProblem"), domainFE("dFE"),
166 boundaryFE("bFE"), skeletonFE("sFE"), dIm(-1), addSkeletonFE(false),
167 addBoundaryFE(false), addParentAdjacencies(false),
170 PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleSetUP);
171 PetscLogEventRegister("SimpleLoadMesh", 0, &MOFEM_EVENT_SimpleLoadMesh);
172 PetscLogEventRegister("SimpleBuildFields", 0, &MOFEM_EVENT_SimpleBuildFields);
173 PetscLogEventRegister("SimpleBuildFiniteElements", 0,
175 PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleBuildProblem);
176 PetscLogEventRegister("SimpleKSPSolve", 0, &MOFEM_EVENT_SimpleKSPSolve);
177 strcpy(meshFileName, "mesh.h5m");
178}
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:1094
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:491
MoFEM::Core & cOre
Definition: Simple.hpp:478
BitRefLevel bitAdjEnt
bit ref level for parent
Definition: Simple.hpp:500
char meshFileName[255]
mesh file name
Definition: Simple.hpp:521
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:490
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:485
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:484
bool addSkeletonFE
Add skeleton FE.
Definition: Simple.hpp:494
BitRefLevel bitAdjParent
bit ref level for parent
Definition: Simple.hpp:498
bool addBoundaryFE
Add boundary FE.
Definition: Simple.hpp:495
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:517
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:516
std::string nameOfProblem
problem name
Definition: Simple.hpp:514
std::string domainFE
domain finite element
Definition: Simple.hpp:515
BitRefLevel bitAdjEntMask
bit ref level for parent parent
Definition: Simple.hpp:501
BitRefLevel bitLevel
BitRefLevel of the problem.
Definition: Simple.hpp:480
int dIm
dimension of problem
Definition: Simple.hpp:522
PetscLogEvent MOFEM_EVENT_SimpleSetUP
Definition: Simple.hpp:483
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:488
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:486
BitRefLevel bitLevelMask
BitRefLevel of the problem.
Definition: Simple.hpp:481
BitRefLevel bitAdjParentMask
bit ref level for parent parent
Definition: Simple.hpp:499
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:487
EntityHandle skeletonMeshset
skeleton meshset with boundary
Definition: Simple.hpp:492
bool addParentAdjacencies
If set to true parent adjacencies are build.
Definition: Simple.hpp:496

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

286 {
287 Interface &m_field = cOre;
289 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
290 verb);
291 boundaryFields.push_back(name);
292 if (space == NOFIELD)
293 SETERRQ(
294 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
295 "NOFIELD space for boundary filed not implemented in Simple interface");
297}
@ 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:1975
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:504

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

324 {
325
326 Interface &m_field = cOre;
328 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
329 verb);
330 if (space == NOFIELD)
331 noFieldDataFields.push_back(name);
332 else
333 dataFields.push_back(name);
335}
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:506
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:508

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

268 {
269
270 Interface &m_field = cOre;
272 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
273 verb);
274 if (space == NOFIELD)
275 noFieldFields.push_back(name);
276 else
277 domainFields.push_back(name);
279}
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:503
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:507

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

707 {
708 Interface &m_field = cOre;
710 CHKERR m_field.getInterface<ProblemsManager>()->addFieldToEmptyFieldBlocks(
711 getProblemName(), row_field, col_field);
713}
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
add empty block to problem
Definition: Simple.cpp:706
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:348

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

304 {
305
306 Interface &m_field = cOre;
308 CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
309 verb);
310 skeletonFields.push_back(name);
311 if (space == NOFIELD)
312 SETERRQ(
313 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
314 "NOFIELD space for boundary filed not implemented in Simple interface");
315
317}
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:505

◆ buildFields()

MoFEMErrorCode MoFEM::Simple::buildFields ( )

Build fields.

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

Definition at line 485 of file Simple.cpp.

485 {
486 Interface &m_field = cOre;
488 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
489
490 auto comm_interface_ptr = m_field.getInterface<CommInterface>();
491 auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
492
493 // Add entities to the fields
494 auto add_ents_to_field = [&](auto meshset, auto dim, auto &fields) {
496 for (auto &field : fields) {
497 CHKERR m_field.add_ents_to_field_by_dim(meshset, dim, field);
498 CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
499 }
501 };
502
503 auto make_no_field_ents = [&](auto &fields) {
505 for (auto &field : fields) {
506 std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
507 CHKERR m_field.create_vertices_and_add_to_field(field, coords.data(), 1);
508 CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
509 CHKERR bit_ref_ptr->setFieldEntitiesBitRefLevel(field, bitLevel);
510 }
512 };
513
514 CHKERR add_ents_to_field(meshSet, dIm, domainFields);
515 CHKERR add_ents_to_field(meshSet, dIm, dataFields);
516 CHKERR add_ents_to_field(boundaryMeshset, dIm - 1, boundaryFields);
517 CHKERR add_ents_to_field(skeletonMeshset, dIm - 1, skeletonFields);
518
519 std::set<std::string> nofield_fields;
520 for (auto &f : noFieldFields)
521 nofield_fields.insert(f);
522 for (auto &f : noFieldDataFields)
523 nofield_fields.insert(f);
524
525 CHKERR make_no_field_ents(nofield_fields);
526
527 // Set order
528
529 auto get_field_ptr = [&](auto &f) { return m_field.get_field_structure(f); };
530 for (auto &t : fieldsOrder) {
531 const auto f = std::get<0>(t);
532 const auto order = std::get<1>(t);
533
534 MOFEM_TAG_AND_LOG("WORLD", Sev::inform, "Simple")
535 << "Set order to field " << f << " order " << order;
536 MOFEM_LOG_CHANNEL("SYNC");
537 if (std::get<3>(t)) {
538 MOFEM_TAG_AND_LOG("SYNC", Sev::verbose, "Simple")
539 << "To ents: " << std::endl
540 << std::get<2>(t) << std::endl;
541 }
542 MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::verbose);
543
544 if (std::get<3>(t)) {
545
546 CHKERR m_field.set_field_order(std::get<2>(t), f, order);
547
548 } else {
549 auto f_ptr = get_field_ptr(f);
550
551 if (f_ptr->getSpace() == H1) {
552 if (f_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
553 CHKERR m_field.set_field_order(meshSet, MBVERTEX, f, order);
554 } else {
555 CHKERR m_field.set_field_order(meshSet, MBVERTEX, f, 1);
556 }
557 }
558
559 for (auto d = 1; d <= dIm; ++d) {
560 for (EntityType t = CN::TypeDimensionMap[d].first;
561 t <= CN::TypeDimensionMap[d].second; ++t) {
562 CHKERR m_field.set_field_order(meshSet, t, f, order);
563 }
564 }
565 }
566 }
567 MOFEM_LOG_CHANNEL("WORLD");
568 // Build fields
569 CHKERR m_field.build_fields();
570 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
572}
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
Definition: LogManager.hpp:352
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:359
@ 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:284
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
constexpr double t
plate stiffness
Definition: plate.cpp:59
std::list< std::tuple< std::string, int, Range, bool > > fieldsOrder
Definition: Simple.hpp:511

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

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

Definition at line 574 of file Simple.cpp.

574 {
575 Interface &m_field = cOre;
577 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
578 // Add finite elements
579 CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
580 true);
581 CHKERR m_field.build_finite_elements(domainFE);
582 if (addBoundaryFE || boundaryFields.size()) {
583 CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
584 boundaryFE, true);
585 CHKERR m_field.build_finite_elements(boundaryFE);
586 }
587 if (addSkeletonFE || skeletonFields.size()) {
588 CHKERR m_field.add_ents_to_finite_element_by_dim(skeletonMeshset, dIm - 1,
589 skeletonFE, true);
590 CHKERR m_field.build_finite_elements(skeletonFE);
591 }
592 for (std::vector<std::string>::iterator fit = otherFEs.begin();
593 fit != otherFEs.end(); ++fit) {
594 CHKERR m_field.build_finite_elements(*fit);
595 }
596 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
598}
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:519

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

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

Definition at line 600 of file Simple.cpp.

600 {
601 Interface &m_field = cOre;
603 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
604 CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
605 // Set problem by the DOFs on the fields rather that by adding DOFs on the
606 // elements
607 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
609 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
610 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
612}
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMoFEM.cpp:1269
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:525

◆ createBoundaryMeshset()

MoFEMErrorCode MoFEM::Simple::createBoundaryMeshset ( )
private

Definition at line 715 of file Simple.cpp.

715 {
716 Interface &m_field = cOre;
718 ParallelComm *pcomm =
719 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
720
721 auto get_skin = [&](auto meshset) {
722 // filter not owned entities, those are not on boundary
723
724 Range domain_ents;
725 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
726 domain_ents, true);
727 CHKERR pcomm->filter_pstatus(domain_ents,
728 PSTATUS_SHARED | PSTATUS_MULTISHARED,
729 PSTATUS_NOT, -1, nullptr);
730
731 Skinner skin(&m_field.get_moab());
732 Range domain_skin;
733 CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
734 CHKERR pcomm->filter_pstatus(domain_skin,
735 PSTATUS_SHARED | PSTATUS_MULTISHARED,
736 PSTATUS_NOT, -1, nullptr);
737 return domain_skin;
738 };
739
740 auto create_boundary_meshset = [&](auto &&domain_skin) {
742 // create boundary meshset
743 if (boundaryMeshset != 0) {
745 }
746 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
747 CHKERR m_field.get_moab().add_entities(boundaryMeshset, domain_skin);
748 for (int dd = 0; dd != dIm - 1; dd++) {
749 Range adj;
750 CHKERR m_field.get_moab().get_adjacencies(domain_skin, dd, false, adj,
751 moab::Interface::UNION);
752 CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
753 }
755 };
756
757 CHKERR create_boundary_meshset(get_skin(meshSet));
758
760}
#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 762 of file Simple.cpp.

762 {
763 Interface &m_field = cOre;
765
766 ParallelComm *pcomm =
767 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
768
769 auto create_skeleton_meshset = [&](auto meshset) {
771 // create boundary meshset
772 if (skeletonMeshset != 0) {
774 }
775 Range boundary_ents, skeleton_ents;
776 CHKERR m_field.get_moab().get_entities_by_dimension(boundaryMeshset,
777 dIm - 1, boundary_ents);
778 Range domain_ents;
779 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
780 domain_ents, true);
781 CHKERR m_field.get_moab().get_adjacencies(
782 domain_ents, dIm - 1, false, skeleton_ents, moab::Interface::UNION);
783 skeleton_ents = subtract(skeleton_ents, boundary_ents);
784 CHKERR pcomm->filter_pstatus(skeleton_ents, PSTATUS_NOT_OWNED, PSTATUS_NOT,
785 -1, nullptr);
786 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, skeletonMeshset);
787 CHKERR m_field.get_moab().add_entities(skeletonMeshset, skeleton_ents);
789 };
790
791 CHKERR create_skeleton_meshset(meshSet);
792
794}

◆ defineFiniteElements()

MoFEMErrorCode MoFEM::Simple::defineFiniteElements ( )

Define finite elements.

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

Definition at line 383 of file Simple.cpp.

383 {
384 Interface &m_field = cOre;
386
387 auto clear_rows_and_cols = [&](auto &fe_name) {
389 auto fe_ptr = m_field.get_finite_elements();
390 auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
391 ->get<FiniteElement_name_mi_tag>();
392 auto it_fe = fe_by_name.find(fe_name);
393 if (it_fe != fe_by_name.end()) {
394
395 if (!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
396 SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
397 "modification unsuccessful");
398
399 if (!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
400 SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
401 "modification unsuccessful");
402 }
404 };
405 CHKERR clear_rows_and_cols(domainFE);
406 CHKERR clear_rows_and_cols(boundaryFE);
407 CHKERR clear_rows_and_cols(skeletonFE);
408
409 // Define finite elements
410 CHKERR m_field.add_finite_element(domainFE, MF_ZERO);
411
412 auto add_fields = [&](auto &fe_name, auto &fields) {
414 for (auto &field : fields) {
415 CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
416 CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
417 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
418 }
420 };
421
422 auto add_data_fields = [&](auto &fe_name, auto &fields) {
424 for (auto &field : fields)
425 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
427 };
428
429 CHKERR add_fields(domainFE, domainFields);
430 CHKERR add_data_fields(domainFE, dataFields);
431 CHKERR add_fields(domainFE, noFieldFields);
432 CHKERR add_data_fields(domainFE, noFieldDataFields);
433
434 if (addBoundaryFE || !boundaryFields.empty()) {
435 CHKERR m_field.add_finite_element(boundaryFE, MF_ZERO);
436 CHKERR add_fields(boundaryFE, domainFields);
437 if (!boundaryFields.empty())
438 CHKERR add_fields(boundaryFE, boundaryFields);
439 CHKERR add_data_fields(boundaryFE, dataFields);
440 CHKERR add_data_fields(boundaryFE, noFieldDataFields);
441 CHKERR add_fields(boundaryFE, noFieldFields);
442 }
443 if (addSkeletonFE || !skeletonFields.empty()) {
444 CHKERR m_field.add_finite_element(skeletonFE, MF_ZERO);
445 CHKERR add_fields(skeletonFE, domainFields);
446 if (!skeletonFields.empty())
447 CHKERR add_fields(skeletonFE, skeletonFields);
448 CHKERR add_data_fields(skeletonFE, dataFields);
449 CHKERR add_data_fields(skeletonFE, noFieldDataFields);
450 CHKERR add_fields(skeletonFE, noFieldFields);
451 }
453}
@ 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 455 of file Simple.cpp.

455 {
456 Interface &m_field = cOre;
458 // Create dm instance
459 dM = createSmartDM(m_field.get_comm(), "DMMOFEM");
460 // set dm data structure which created mofem data structures
463 CHKERR DMSetFromOptions(dM);
465 if (addBoundaryFE || !boundaryFields.empty()) {
467 }
468 if (addSkeletonFE || !skeletonFields.empty()) {
470 }
472 CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
474}
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMoFEM.cpp:1111
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:485
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: DMMoFEM.cpp:118
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 683 of file Simple.cpp.

683 {
685 dM.reset();
687}

◆ deleteFiniteElements()

MoFEMErrorCode MoFEM::Simple::deleteFiniteElements ( )

Delete finite elements.

Returns
MoFEMErrorCode
Returns
MoFEMErrorCode

Definition at line 694 of file Simple.cpp.

694 {
695 Interface &m_field = cOre;
697 for (auto fe : {domainFE, boundaryFE, skeletonFE}) {
698 if (m_field.check_finite_element(fe)) {
699 CHKERR m_field.delete_finite_element(fe);
700 }
701 }
703}

◆ exchangeGhostCells()

MoFEMErrorCode MoFEM::Simple::exchangeGhostCells ( )
private

get all adjacent ghosted entities

Definition at line 796 of file Simple.cpp.

796 {
797 Interface &m_field = cOre;
799 MOFEM_LOG("WORLD", Sev::verbose) << "Exchange ghost cells";
800
801 ParallelComm *pcomm =
802 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
803 if (pcomm == NULL)
804 pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
805
806 CHKERR pcomm->exchange_ghost_cells(getDim(), getDim() - 1, 1,
807 3 /**get all adjacent ghosted entities */,
808 true, false, meshSet ? &meshSet : nullptr);
809
810 Range shared;
811 CHKERR m_field.get_moab().get_entities_by_dimension(0, dIm, shared);
812 for (auto d = dIm - 1; d >= 0; --d) {
813 CHKERR m_field.get_moab().get_adjacencies(shared, d, false, shared,
814 moab::Interface::UNION);
815 }
816 CHKERR pcomm->filter_pstatus(shared, PSTATUS_SHARED | PSTATUS_MULTISHARED,
817 PSTATUS_OR, -1, &shared);
818 Tag part_tag = pcomm->part_tag();
819 CHKERR pcomm->exchange_tags(part_tag, shared);
820 CHKERR m_field.getInterface<CommInterface>()->resolveParentEntities(shared,
821 VERBOSE);
822
824}
@ VERBOSE
Definition: definitions.h:209
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
int getDim() const
Get the problem dimension.
Definition: Simple.hpp:271

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

422{ 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 411 of file Simple.hpp.

411{ return addSkeletonFE; }

◆ getBitAdjEnt()

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

bit ref level for parent

Returns
auto&

Definition at line 453 of file Simple.hpp.

453{ return bitAdjEnt; }

◆ getBitAdjEntMask()

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

bit ref level for parent parent

Returns
auto&

Definition at line 460 of file Simple.hpp.

460{ return bitAdjEntMask; }

◆ getBitAdjParent()

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

bit ref level for parent

Returns
auto&

Definition at line 439 of file Simple.hpp.

439{ return bitAdjParent; }

◆ getBitAdjParentMask()

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

bit ref level for parent parent

Returns
auto&

Definition at line 446 of file Simple.hpp.

446{ return bitAdjParentMask; }

◆ getBitRefLevel()

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

Get the BitRefLevel.

Returns
BitRefLevel
Examples
level_set.cpp.

Definition at line 313 of file Simple.hpp.

313{ return bitLevel; }

◆ getBitRefLevelMask()

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

Get the BitRefLevel.

Returns
BitRefLevel

Definition at line 320 of file Simple.hpp.

320{ return bitLevelMask; }

◆ getBoundaryFEName() [1/2]

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

Get the Boundary FE Name.

Returns
std::string&

Definition at line 362 of file Simple.hpp.

362{ 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 334 of file Simple.hpp.

334{ return boundaryFE; }

◆ getBoundaryMeshSet()

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

Get the BoundaryMeshSet object.

Returns
EntityHandle&

Definition at line 299 of file Simple.hpp.

299{ 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 271 of file Simple.hpp.

271{ 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 261 of file Simple.hpp.

261{ 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 355 of file Simple.hpp.

355{ return domainFE; }

◆ getDomainFEName() [2/2]

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

◆ getMeshSet()

DEPRECATED EntityHandle & MoFEM::Simple::getMeshSet ( )
inline
Deprecated:
Use getMeshset
Returns
EntityHandle&

Definition at line 285 of file Simple.hpp.

285{ return meshSet; }

◆ getMeshset()

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

Get the MeshSet object.

Returns
EntityHandle&

Definition at line 292 of file Simple.hpp.

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

180 {
181 PetscBool flg = PETSC_TRUE;
183 ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "Simple interface options",
184 "none");
185 CHKERRG(ierr);
186 ierr = PetscOptionsString("-file_name", "file name", "", "mesh.h5m",
187 meshFileName, 255, &flg);
188 CHKERRG(ierr);
189 ierr = PetscOptionsEnd();
190 CHKERRG(ierr);
192}
#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 386 of file Simple.hpp.

386{ return otherFEs; }

◆ getParentAdjacencies()

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

Get the addParentAdjacencies.

If set true add parent adjacencies

Returns
true
false

Definition at line 432 of file Simple.hpp.

432{ return addParentAdjacencies; }

◆ getProblemName() [1/2]

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

Get the Problem Name.

Returns
std::string&

Definition at line 376 of file Simple.hpp.

376{ 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 348 of file Simple.hpp.

348{ return nameOfProblem; }

◆ getSkeletonFEName() [1/2]

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

Get the Skeleton FE Name.

Returns
std::string&

Definition at line 369 of file Simple.hpp.

369{ 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 306 of file Simple.hpp.

306{ 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 250 of file Simple.cpp.

250 {
252 Interface &m_field = cOre;
253 if (m_field.get_comm_size() == 1)
255 else
256 CHKERR loadFile("PARALLEL=READ_PART;"
257 "PARALLEL_RESOLVE_SHARED_ENTS;"
258 "PARTITION=PARALLEL_PARTITION;",
261}
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:194

◆ 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.
Note
If bitRefLevel is set to any, bit ref level of loaded entities is not changed. After mesh is load, bit ref level should be set to create problem. Default setting of bit ref level is on first bit, and if is set all mesh entities on load are set to set 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 194 of file Simple.cpp.

195 {
196 Interface &m_field = cOre;
198 PetscLogEventBegin(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
199
200 if (!mesh_file_name.empty())
201 strcpy(meshFileName, mesh_file_name.c_str());
202
203 // This is a case of distributed mesh and algebra. In that case each processor
204 // keep only part of the problem.
205 CHKERR m_field.get_moab().load_file(meshFileName, 0, options.c_str());
206 CHKERR m_field.rebuild_database();
207
208 // determine problem dimension
209 if (dIm == -1) {
210 int nb_ents_3d;
211 CHKERR m_field.get_moab().get_number_entities_by_dimension(
212 meshSet, 3, nb_ents_3d, true);
213 if (nb_ents_3d > 0) {
214 dIm = 3;
215 } else {
216 int nb_ents_2d;
217 CHKERR m_field.get_moab().get_number_entities_by_dimension(
218 meshSet, 2, nb_ents_2d, true);
219 if (nb_ents_2d > 0) {
220 dIm = 2;
221 } else {
222 dIm = 1;
223 }
224 }
225 }
226
227 if (!boundaryMeshset)
229 if (!skeletonMeshset)
231 if (addSkeletonFE)
233
234 if (bitLevel.any()) {
235 Range ents;
236 CHKERR m_field.get_moab().get_entities_by_dimension(meshSet, dIm, ents,
237 true);
238 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(ents, bitLevel,
239 false);
240 } else {
241 MOFEM_LOG("WORLD", Sev::warning) << "BitRefLevel is none and not set";
242 CHKERR m_field.getInterface<BitRefManager>()->addToDatabaseBitRefLevelByDim(
243 dIm, BitRefLevel().set(), BitRefLevel().set());
244 }
245
246 PetscLogEventEnd(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
248}
char mesh_file_name[255]
MoFEMErrorCode exchangeGhostCells()
Definition: Simple.cpp:796
MoFEMErrorCode createBoundaryMeshset()
Definition: Simple.cpp:715
MoFEMErrorCode createSkeletonMeshset()
Definition: Simple.cpp:762

◆ 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:162

◆ removeBoundaryField()

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

Remove field form boundary.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 353 of file Simple.cpp.

353 {
354 Interface &m_field = cOre;
356
357 auto remove_field_from_list = [&](auto &vec) {
358 auto it = std::find(vec.begin(), vec.end(), name);
359 if (it != vec.end())
360 vec.erase(it);
361 };
362
363 remove_field_from_list(boundaryFields);
364
366}

◆ removeDomainField()

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

Remove field form domain.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 337 of file Simple.cpp.

337 {
338 Interface &m_field = cOre;
340
341 auto remove_field_from_list = [&](auto &vec) {
342 auto it = std::find(vec.begin(), vec.end(), name);
343 if (it != vec.end())
344 vec.erase(it);
345 };
346
347 remove_field_from_list(noFieldFields);
348 remove_field_from_list(domainFields);
349
351}

◆ removeSkeletonField()

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

Remove field form skeleton.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 368 of file Simple.cpp.

368 {
369 Interface &m_field = cOre;
371
372 auto remove_field_from_list = [&](auto &vec) {
373 auto it = std::find(vec.begin(), vec.end(), name);
374 if (it != vec.end())
375 vec.erase(it);
376 };
377
378 remove_field_from_list(skeletonFields);
379
381}

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

637 {
638 Interface &m_field = cOre;
640 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
641
642 if (!only_dm) {
645 if (addSkeletonFE || !skeletonFields.empty())
650 }
651
652 CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
653
654 const Problem *problem_ptr;
655 CHKERR DMMoFEMGetProblemPtr(dM, &problem_ptr);
656 const auto problem_name = problem_ptr->getName();
657 CHKERR m_field.modify_problem_ref_level_set_bit(problem_name, bitLevel);
658 CHKERR m_field.modify_problem_mask_ref_level_set_bit(problem_name,
660
661 // Set problem by the DOFs on the fields rather that by adding DOFs on the
662 // elements
663 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
665 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
666
667 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
669}
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMoFEM.cpp:414
MoFEMErrorCode setParentAdjacency()
Definition: Simple.cpp:67
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:383
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:574
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:485
MoFEMErrorCode setSkeletonAdjacency(int dim=-1, std::string fe_name="")
Set the skeleton adjacency object.
Definition: Simple.cpp:143

◆ setDim()

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

Set the problem dimension.

Returns
int

Definition at line 278 of file Simple.hpp.

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

72 {
73 Interface &m_field = cOre;
75
77 boost::make_shared<ParentFiniteElementAdjacencyFunction<3>>(
80 boost::make_shared<ParentFiniteElementAdjacencyFunction<2>>(
82
83 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBTET,
85 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBHEX,
87 if (addBoundaryFE || !boundaryFields.empty()) {
88 CHKERR m_field.modify_finite_element_adjacency_table(
90 CHKERR m_field.modify_finite_element_adjacency_table(
92 }
93 if (addSkeletonFE || !skeletonFields.empty()) {
94 CHKERR m_field.modify_finite_element_adjacency_table(
96 CHKERR m_field.modify_finite_element_adjacency_table(
98 }
99
101}
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 2 > > parentAdjFunctionDim2
Definition: Simple.hpp:539
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
Definition: Simple.hpp:537

◆ setParentAdjacency() [2/4]

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

Definition at line 103 of file Simple.cpp.

103 {
104 Interface &m_field = cOre;
106
108 boost::make_shared<ParentFiniteElementAdjacencyFunction<2>>(
111 boost::make_shared<ParentFiniteElementAdjacencyFunction<1>>(
113
114 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBTRI,
116 CHKERR m_field.modify_finite_element_adjacency_table(domainFE, MBQUAD,
118 if (addBoundaryFE || !boundaryFields.empty())
119 CHKERR m_field.modify_finite_element_adjacency_table(
121 if (addSkeletonFE || !skeletonFields.empty())
122 CHKERR m_field.modify_finite_element_adjacency_table(
124
126}
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 1 > > parentAdjFunctionDim1
Definition: Simple.hpp:541

◆ setParentAdjacency() [3/4]

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

Definition at line 128 of file Simple.cpp.

128 {
130 switch (getDim()) {
131 case 1:
132 THROW_MESSAGE("Not implemented");
133 case 2:
134 return setParentAdjacency<2>();
135 case 3:
136 return setParentAdjacency<3>();
137 default:
138 THROW_MESSAGE("Not implemented");
139 }
141}
#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 67 of file Simple.cpp.

67 {
68 static_assert(DIM == 2 || DIM == 3, "not implemented");
70}

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

143 {
145 if (dim == -1)
146 dim = getDim();
147
148 if (fe_name.empty())
149 fe_name = skeletonFE;
150
151 switch (dim) {
152 case 2:
153 return setSkeletonAdjacency<2>(fe_name);
154 case 3:
155 return setSkeletonAdjacency<3>(fe_name);
156 default:
157 SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_IMPLEMENTED, "Not implemented");
158 }
160}

◆ 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
26 boost::make_shared<ParentFiniteElementAdjacencyFunctionSkeleton<1>>(
28 CHKERR m_field.modify_finite_element_adjacency_table(
29 fe_name, MBEDGE, *parentAdjSkeletonFunctionDim1);
30
32}
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 1 > > parentAdjSkeletonFunctionDim1
Definition: Simple.hpp:546

◆ setSkeletonAdjacency() [3/5]

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

Definition at line 35 of file Simple.cpp.

35 {
36 Interface &m_field = cOre;
38
40 boost::make_shared<ParentFiniteElementAdjacencyFunctionSkeleton<2>>(
42
43 CHKERR m_field.modify_finite_element_adjacency_table(
44 fe_name, MBTRI, *parentAdjSkeletonFunctionDim2);
45 CHKERR m_field.modify_finite_element_adjacency_table(
46 fe_name, MBQUAD, *parentAdjSkeletonFunctionDim2);
47
49}
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 2 > > parentAdjSkeletonFunctionDim2
Definition: Simple.hpp:544

◆ setSkeletonAdjacency() [4/5]

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

Definition at line 52 of file Simple.cpp.

52 {
54
55 switch (getDim()) {
56 case 1:
57 THROW_MESSAGE("Not implemented");
58 case 2:
59 return setSkeletonAdjacency<2>(fe_name);
60 case 3:
61 return setSkeletonAdjacency<3>(fe_name);
62 default:
63 THROW_MESSAGE("Not implemented");
64 }
66}

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

◆ addParentAdjacencies

bool MoFEM::Simple::addParentAdjacencies
private

If set to true parent adjacencies are build.

Definition at line 496 of file Simple.hpp.

◆ addSkeletonFE

bool MoFEM::Simple::addSkeletonFE
private

Add skeleton FE.

Definition at line 494 of file Simple.hpp.

◆ bitAdjEnt

BitRefLevel MoFEM::Simple::bitAdjEnt
private

bit ref level for parent

Definition at line 500 of file Simple.hpp.

◆ bitAdjEntMask

BitRefLevel MoFEM::Simple::bitAdjEntMask
private

bit ref level for parent parent

Definition at line 501 of file Simple.hpp.

◆ bitAdjParent

BitRefLevel MoFEM::Simple::bitAdjParent
private

bit ref level for parent

Definition at line 498 of file Simple.hpp.

◆ bitAdjParentMask

BitRefLevel MoFEM::Simple::bitAdjParentMask
private

bit ref level for parent parent

Definition at line 499 of file Simple.hpp.

◆ bitLevel

BitRefLevel MoFEM::Simple::bitLevel
private

BitRefLevel of the problem.

Definition at line 480 of file Simple.hpp.

◆ bitLevelMask

BitRefLevel MoFEM::Simple::bitLevelMask
private

BitRefLevel of the problem.

Definition at line 481 of file Simple.hpp.

◆ boundaryFE

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

boundary finite element

Definition at line 516 of file Simple.hpp.

◆ boundaryFields

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

boundary fields

Definition at line 504 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 491 of file Simple.hpp.

◆ cOre

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

Definition at line 478 of file Simple.hpp.

◆ dataFields

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

Data fields.

Definition at line 506 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 522 of file Simple.hpp.

◆ dM

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

Discrete manager (interface to PETSc/MoFEM functions)

Definition at line 525 of file Simple.hpp.

◆ domainFE

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

domain finite element

Definition at line 515 of file Simple.hpp.

◆ domainFields

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

domain fields

Definition at line 503 of file Simple.hpp.

◆ fieldsOrder

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

fields order. 1: field name, order, range, set by range if true

Definition at line 511 of file Simple.hpp.

◆ meshFileName

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

mesh file name

Definition at line 521 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 490 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 485 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 486 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 487 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 488 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 484 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleSetUP

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleSetUP
private

Definition at line 483 of file Simple.hpp.

◆ nameOfProblem

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

problem name

Definition at line 514 of file Simple.hpp.

◆ noFieldDataFields

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

NOFIELD field name.

Definition at line 508 of file Simple.hpp.

◆ noFieldFields

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

NOFIELD field name.

Definition at line 507 of file Simple.hpp.

◆ otherFEs

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

Other finite elements.

Definition at line 519 of file Simple.hpp.

◆ parentAdjFunctionDim1

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

Definition at line 541 of file Simple.hpp.

◆ parentAdjFunctionDim2

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

Definition at line 539 of file Simple.hpp.

◆ parentAdjFunctionDim3

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

Definition at line 537 of file Simple.hpp.

◆ parentAdjSkeletonFunctionDim1

boost::shared_ptr<ParentFiniteElementAdjacencyFunctionSkeleton<1> > MoFEM::Simple::parentAdjSkeletonFunctionDim1
private

Definition at line 546 of file Simple.hpp.

◆ parentAdjSkeletonFunctionDim2

boost::shared_ptr<ParentFiniteElementAdjacencyFunctionSkeleton<2> > MoFEM::Simple::parentAdjSkeletonFunctionDim2
private

Definition at line 544 of file Simple.hpp.

◆ skeletonFE

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

skeleton finite element

Definition at line 517 of file Simple.hpp.

◆ skeletonFields

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

fields on the skeleton

Definition at line 505 of file Simple.hpp.

◆ skeletonMeshset

EntityHandle MoFEM::Simple::skeletonMeshset
private

skeleton meshset with boundary

Definition at line 492 of file Simple.hpp.


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