v0.15.0
Loading...
Searching...
No Matches
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 Types

typedef boost::function< MoFEMErrorCode(Interface &, const char *, const char *)> LoadFileFunc
 

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
 
MoFEMErrorCode loadFile (const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
 Load mesh file.
 
MoFEMErrorCode loadFile (const std::string mesh_file_name="")
 Load mesh file with parallel options if number of cores > 1.
 
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.
 
MoFEMErrorCode addDomainBrokenField (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 broken field on domain.
 
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.
 
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.
 
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.
 
MoFEMErrorCode addMeshsetField (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 meshset field.
 
MoFEMErrorCode removeDomainField (const std::string name)
 Remove field form domain.
 
MoFEMErrorCode removeBoundaryField (const std::string name)
 Remove field form boundary.
 
MoFEMErrorCode removeSkeletonField (const std::string name)
 Remove field form skeleton.
 
MoFEMErrorCode defineFiniteElements ()
 Define finite elements.
 
MoFEMErrorCode defineProblem (const PetscBool is_partitioned=PETSC_TRUE)
 define problem
 
MoFEMErrorCode setFieldOrder (const std::string field_name, const int order, const Range *ents=NULL)
 Set field order.
 
MoFEMErrorCode buildFields ()
 Build fields.
 
MoFEMErrorCode buildFiniteElements ()
 Build finite elements.
 
MoFEMErrorCode setSkeletonAdjacency (int dim=-1, std::string fe_name="")
 Set the skeleton adjacency object.
 
MoFEMErrorCode buildProblem ()
 Build problem.
 
MoFEMErrorCode setUp (const PetscBool is_partitioned=PETSC_TRUE)
 Setup problem.
 
MoFEMErrorCode reSetUp (bool only_dm=false)
 Rebuild internal MoFEM data structures.
 
MoFEMErrorCode getDM (DM *dm)
 Get DM.
 
SmartPetscObj< DM > getDM ()
 Return smart DM object.
 
int getDim () const
 Get the problem dimension.
 
void setDim (int dim)
 Set the problem dimension.
 
DEPRECATED EntityHandlegetMeshSet ()
 
EntityHandlegetMeshset ()
 Get the MeshSet object.
 
EntityHandlegetBoundaryMeshSet ()
 Get the BoundaryMeshSet object.
 
EntityHandlegetSkeletonMeshSet ()
 Get the SkeletonMeshSet object.
 
BitRefLevelgetBitRefLevel ()
 Get the BitRefLevel.
 
BitRefLevelgetBitRefLevelMask ()
 Get the BitRefLevel.
 
const std::string getDomainFEName () const
 Get the Domain FE Name.
 
const std::string getBoundaryFEName () const
 Get the Boundary FE Name.
 
const std::string getSkeletonFEName () const
 Get the Skeleton FE Name.
 
const std::string getProblemName () const
 Get the Problem Name.
 
std::string & getDomainFEName ()
 Get the Domain FE Name.
 
std::string & getBoundaryFEName ()
 Get the Boundary FE Name.
 
std::string & getSkeletonFEName ()
 Get the Skeleton FE Name.
 
std::string & getMeshsetFEName ()
 Get the Skeleton FE Name.
 
auto & getMeshsetFiniteElementEntities ()
 Range of entities added to meshset finite element.
 
std::string & getProblemName ()
 Get the Problem Name.
 
std::vector< std::string > & getOtherFiniteElements ()
 Get the Other Finite Elements.
 
MoFEMErrorCode deleteDM ()
 Delete dm.
 
MoFEMErrorCode deleteFiniteElements ()
 Delete finite elements.
 
boolgetAddSkeletonFE ()
 Get the addSkeletonFE.
 
boolgetAddBoundaryFE ()
 Get the addSkeletonFE.
 
boolgetParentAdjacencies ()
 Get the addParentAdjacencies.
 
auto & getBitAdjParent ()
 bit ref level for parent
 
auto & getBitAdjParentMask ()
 bit ref level for parent parent
 
auto & getBitAdjEnt ()
 bit ref level for parent
 
auto & getBitAdjEntMask ()
 bit ref level for parent parent
 
MoFEMErrorCode addFieldToEmptyFieldBlocks (const std::string row_field, const std::string col_field) const
 add empty block to problem
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface.
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface reference to pointer of interface.
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface.
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface.
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface.
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface.
 
virtual ~UnknownInterface ()=default
 

Static Public Member Functions

static MoFEMErrorCode defaultLoadFileFunc (Interface &m_field, const char *file_name, const char *options)
 Set function for loading mesh file.
 
- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version.
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version.
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version.
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version.
 

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.
 
BitRefLevel bitLevelMask
 BitRefLevel of the problem.
 
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
 
EntityHandle boundaryMeshset
 meshset with boundary
 
EntityHandle skeletonMeshset
 skeleton meshset with boundary
 
bool addSkeletonFE
 Add skeleton FE.
 
bool addBoundaryFE
 Add boundary FE.
 
bool addParentAdjacencies
 If set to true parent adjacencies are build.
 
BitRefLevel bitAdjParent
 bit ref level for parent
 
BitRefLevel bitAdjParentMask
 bit ref level for parent parent
 
BitRefLevel bitAdjEnt
 bit ref level for parent
 
BitRefLevel bitAdjEntMask
 bit ref level for parent parent
 
std::vector< std::string > domainFields
 domain fields
 
std::vector< std::string > boundaryFields
 boundary fields
 
std::vector< std::string > skeletonFields
 fields on the skeleton
 
std::vector< std::string > dataFields
 Data fields.
 
std::vector< std::string > meshsetFields
 meshset fields
 
std::vector< std::string > noFieldFields
 NOFIELD field name.
 
std::vector< std::string > noFieldDataFields
 NOFIELD field name.
 
std::list< std::tuple< std::string, int, Range, bool > > fieldsOrder
 
std::string nameOfProblem = "SimpleProblem"
 problem name
 
std::string domainFE = "dFE"
 domain finite element
 
std::string boundaryFE = "bFE"
 boundary finite element
 
std::string skeletonFE = "sFE"
 skeleton finite element
 
std::string meshsetFE = "mFE"
 meshset finite element
 
std::vector< std::string > otherFEs
 Other finite elements.
 
std::vector< RangemeshsetFiniteElementEntities
 Meshset element entities.
 
char meshFileName [255]
 mesh file name
 
int dIm
 dimension of problem
 
SmartPetscObj< DM > dM
 Discrete manager (interface to PETSc/MoFEM functions)
 
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
 

Detailed Description

Member Typedef Documentation

◆ LoadFileFunc

boost::function<MoFEMErrorCode(Interface &, const char *, const char *)> MoFEM::Simple::LoadFileFunc

Definition at line 43 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("SimpSetUp", 0, &MOFEM_EVENT_SimpleSetUP);
171 PetscLogEventRegister("SimpLoadMesh", 0, &MOFEM_EVENT_SimpleLoadMesh);
172 PetscLogEventRegister("SimpBuildFields", 0, &MOFEM_EVENT_SimpleBuildFields);
173 PetscLogEventRegister("SimpBuildFEs", 0,
175 PetscLogEventRegister("SimpSetUp", 0, &MOFEM_EVENT_SimpleBuildProblem);
176 PetscLogEventRegister("SimpKSPSolve", 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
EntityHandle boundaryMeshset
meshset with boundary
Definition Simple.hpp:572
MoFEM::Core & cOre
Definition Simple.hpp:559
BitRefLevel bitAdjEnt
bit ref level for parent
Definition Simple.hpp:581
char meshFileName[255]
mesh file name
Definition Simple.hpp:606
EntityHandle meshSet
domain meshset
Definition Simple.hpp:571
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition Simple.hpp:566
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition Simple.hpp:565
bool addSkeletonFE
Add skeleton FE.
Definition Simple.hpp:575
BitRefLevel bitAdjParent
bit ref level for parent
Definition Simple.hpp:579
bool addBoundaryFE
Add boundary FE.
Definition Simple.hpp:576
std::string skeletonFE
skeleton finite element
Definition Simple.hpp:599
std::string boundaryFE
boundary finite element
Definition Simple.hpp:598
std::string nameOfProblem
problem name
Definition Simple.hpp:596
std::string domainFE
domain finite element
Definition Simple.hpp:597
BitRefLevel bitAdjEntMask
bit ref level for parent parent
Definition Simple.hpp:582
BitRefLevel bitLevel
BitRefLevel of the problem.
Definition Simple.hpp:561
int dIm
dimension of problem
Definition Simple.hpp:607
PetscLogEvent MOFEM_EVENT_SimpleSetUP
Definition Simple.hpp:564
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition Simple.hpp:569
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition Simple.hpp:567
BitRefLevel bitLevelMask
BitRefLevel of the problem.
Definition Simple.hpp:562
BitRefLevel bitAdjParentMask
bit ref level for parent parent
Definition Simple.hpp:580
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition Simple.hpp:568
EntityHandle skeletonMeshset
skeleton meshset with boundary
Definition Simple.hpp:573
bool addParentAdjacencies
If set to true parent adjacencies are build.
Definition Simple.hpp:577

◆ ~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 field
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
child_and_parent.cpp, hanging_node_approx.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, helmholtz.cpp, mesh_smoothing.cpp, phase.cpp, and simple_interface.cpp.

Definition at line 355 of file Simple.cpp.

359 {
360 Interface &m_field = cOre;
362 CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
363 verb);
364 boundaryFields.push_back(name);
365 if (space == NOFIELD)
366 SETERRQ(
367 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
368 "NOFIELD space for boundary field not implemented in Simple interface");
370}
@ 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 ...
@ MOFEM_NOT_IMPLEMENTED
Definition definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
DeprecatedCoreInterface Interface
virtual MoFEMErrorCode add_field(const std::string name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)
Add field.
std::vector< std::string > boundaryFields
boundary fields
Definition Simple.hpp:585

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

397 {
398
399 Interface &m_field = cOre;
401 CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
402 verb);
403 if (space == NOFIELD)
404 noFieldDataFields.push_back(name);
405 else
406 dataFields.push_back(name);
408}
std::vector< std::string > dataFields
Data fields.
Definition Simple.hpp:587
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition Simple.hpp:590

◆ addDomainBrokenField()

MoFEMErrorCode MoFEM::Simple::addDomainBrokenField ( 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 broken field on domain.

Parameters
namename of the field
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
plot_base.cpp.

Definition at line 282 of file Simple.cpp.

286 {
287
288 Interface &m_field = cOre;
290
291 auto get_side_map_hcurl = [&]() {
292 return std::vector<
293
294 std::pair<EntityType,
296
297 >>{
298
299 {MBTRI,
300 [&](BaseFunction::DofsSideMap &dofs_side_map) -> MoFEMErrorCode {
302 dofs_side_map);
303 }}
304
305 };
306 };
307
308 auto get_side_map_hdiv = [&]() {
309 return std::vector<
310
311 std::pair<EntityType,
313
314 >>{
315
316 {MBTET,
317 [&](BaseFunction::DofsSideMap &dofs_side_map) -> MoFEMErrorCode {
319 dofs_side_map);
320 }}
321
322 };
323 };
324
325 auto get_side_map = [&](auto space) {
326 switch (space) {
327 case HCURL:
328 return get_side_map_hcurl();
329 case HDIV:
330 return get_side_map_hdiv();
331 default:
332 MOFEM_LOG_CHANNEL("WORLD");
333 MOFEM_LOG("WORLD", Sev::warning)
334 << "Side dof map for broken space <" << space << "> not implemented";
335 }
336 return std::vector<
337
338 std::pair<EntityType,
340
341 >>{};
342 };
343
344 CHKERR m_field.add_broken_field(name, space, base, nb_of_coefficients,
345 get_side_map(space), tag_type, bh, verb);
346 domainFields.push_back(name);
347 if (space == NOFIELD)
348 SETERRQ(
349 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
350 "NOFIELD space for boundary field not implemented in Simple interface");
352}
@ HCURL
field with continuous tangents
Definition definitions.h:86
@ HDIV
field with continuous normal traction
Definition definitions.h:87
@ DISCONTINUOUS
Broken continuity (No effect on L2 space)
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
multi_index_container< DofsSideMapData, indexed_by< ordered_non_unique< tag< TypeSide_mi_tag >, composite_key< DofsSideMapData, member< DofsSideMapData, EntityType, &DofsSideMapData::type >, member< DofsSideMapData, int, &DofsSideMapData::side > > >, ordered_unique< tag< EntDofIdx_mi_tag >, member< DofsSideMapData, int, &DofsSideMapData::dof > > > > DofsSideMap
Map entity stype and side to element/entity dof index.
virtual MoFEMErrorCode add_broken_field(const std::string name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_coefficients, std::vector< std::pair< EntityType, std::function< MoFEMErrorCode(BaseFunction::DofsSideMap &)> > > list_dof_side_map, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)
Add field.
std::vector< std::string > domainFields
domain fields
Definition Simple.hpp:584
static MoFEMErrorCode setDofsSideMap(const FieldSpace space, const FieldContinuity continuity, const FieldApproximationBase base, DofsSideMap &)
Set map of dof to side number.
static MoFEMErrorCode setDofsSideMap(const FieldSpace space, const FieldContinuity continuity, const FieldApproximationBase base, DofsSideMap &)
Set map of dof to side number.

◆ 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 field
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
approx_sphere.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, dg_projection.cpp, eigen_elastic.cpp, field_evaluator.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, helmholtz.cpp, higher_derivatives.cpp, level_set.cpp, loop_entities.cpp, mesh_smoothing.cpp, mixed_poisson.cpp, phase.cpp, photon_diffusion.cpp, plot_base.cpp, poisson_2d_dis_galerkin.cpp, poisson_2d_homogeneous.cpp, reaction_diffusion.cpp, seepage.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 261 of file Simple.cpp.

265 {
266
267 Interface &m_field = cOre;
269 CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
270 verb);
271 domainFields.push_back(name);
272 if (space == NOFIELD)
273 SETERRQ(
274 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
275 "NOFIELD space for boundary field not implemented in Simple interface");
276
277
279}

◆ 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 field name
col_fieldcol field name
Returns
MoFEMErrorCode
Examples
phase.cpp.

Definition at line 835 of file Simple.cpp.

836 {
837 Interface &m_field = cOre;
839 CHKERR m_field.getInterface<ProblemsManager>()->addFieldToEmptyFieldBlocks(
840 getProblemName(), row_field, col_field);
842}
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
add empty block to problem
Definition Simple.cpp:835
const std::string getProblemName() const
Get the Problem Name.
Definition Simple.hpp:410
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.

◆ addMeshsetField()

MoFEMErrorCode MoFEM::Simple::addMeshsetField ( 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 meshset field.

Add fields managed by meshset finite element

Parameters
name
space
base
nb_of_coefficients
tag_type
bh
verb
Returns
MoFEMErrorCode
Examples
simple_interface.cpp, and simple_l2_only.cpp.

Definition at line 411 of file Simple.cpp.

415 {
416
417 Interface &m_field = cOre;
419 CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
420 verb);
421 if (space != NOFIELD) {
422 meshsetFields.push_back(name);
423 } else {
424 noFieldFields.push_back(name);
425 }
427}
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition Simple.hpp:589
std::vector< std::string > meshsetFields
meshset fields
Definition Simple.hpp:588

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

377 {
378
379 Interface &m_field = cOre;
381 CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
382 verb);
383 skeletonFields.push_back(name);
384 if (space == NOFIELD)
385 SETERRQ(
386 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
387 "NOFIELD space for boundary field not implemented in Simple interface");
388
390}
std::vector< std::string > skeletonFields
fields on the skeleton
Definition Simple.hpp:586

◆ buildFields()

MoFEMErrorCode MoFEM::Simple::buildFields ( )

Build fields.

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

Definition at line 584 of file Simple.cpp.

584 {
585 Interface &m_field = cOre;
587 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
588
589 auto comm_interface_ptr = m_field.getInterface<CommInterface>();
590 auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
591
592 // Add entities to the fields
593 auto add_ents_to_field = [&](auto meshset, auto dim, auto &fields) {
595 for (auto &field : fields) {
596 CHKERR m_field.add_ents_to_field_by_dim(meshset, dim, field);
597 CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
598 }
600 };
601
602 CHKERR add_ents_to_field(meshSet, dIm, domainFields);
603 CHKERR add_ents_to_field(meshSet, dIm, dataFields);
604 CHKERR add_ents_to_field(boundaryMeshset, dIm - 1, boundaryFields);
605 CHKERR add_ents_to_field(skeletonMeshset, dIm - 1, skeletonFields);
606
607 std::set<std::string> nofield_fields;
608 for (auto &f : noFieldFields)
609 nofield_fields.insert(f);
610 for (auto &f : noFieldDataFields)
611 nofield_fields.insert(f);
612
613 // Set order
614
615 auto get_field_ptr = [&](auto &f) { return m_field.get_field_structure(f); };
616 for (auto &t : fieldsOrder) {
617 const auto f = std::get<0>(t);
618 const auto order = std::get<1>(t);
619
620 MOFEM_TAG_AND_LOG("WORLD", Sev::inform, "Simple")
621 << "Set order to field " << f << " order " << order;
622 MOFEM_LOG_CHANNEL("SYNC");
623 if (std::get<3>(t)) {
624 MOFEM_TAG_AND_LOG("SYNC", Sev::verbose, "Simple")
625 << "To ents: " << std::endl
626 << std::get<2>(t) << std::endl;
627 }
628 MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::verbose);
629
630 if (std::get<3>(t)) {
631
632 CHKERR m_field.set_field_order(std::get<2>(t), f, order);
633
634 } else {
635 auto f_ptr = get_field_ptr(f);
636
637 if (f_ptr->getSpace() == H1) {
638 if (f_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
639 CHKERR m_field.set_field_order(meshSet, MBVERTEX, f, order);
640 } else {
641 CHKERR m_field.set_field_order(meshSet, MBVERTEX, f, 1);
642 }
643 }
644
645 for (auto d = 1; d <= dIm; ++d) {
646 for (EntityType t = CN::TypeDimensionMap[d].first;
647 t <= CN::TypeDimensionMap[d].second; ++t) {
648 CHKERR m_field.set_field_order(meshSet, t, f, order);
649 }
650 }
651 }
652 }
653 MOFEM_LOG_CHANNEL("WORLD");
654 // Build fields
655 CHKERR m_field.build_fields();
656 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
658}
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition definitions.h:64
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ 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 ...
constexpr int order
constexpr double t
plate stiffness
Definition plate.cpp:58
std::list< std::tuple< std::string, int, Range, bool > > fieldsOrder
Definition Simple.hpp:593

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

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

Definition at line 660 of file Simple.cpp.

660 {
661 Interface &m_field = cOre;
663 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
664 // Add finite elements
665 CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
666 true);
667 CHKERR m_field.build_finite_elements(domainFE);
668 if (addBoundaryFE || boundaryFields.size()) {
669 CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
670 boundaryFE, true);
671 CHKERR m_field.build_finite_elements(boundaryFE);
672 }
673 if (addSkeletonFE || skeletonFields.size()) {
674 CHKERR m_field.add_ents_to_finite_element_by_dim(skeletonMeshset, dIm - 1,
675 skeletonFE, true);
676 CHKERR m_field.build_finite_elements(skeletonFE);
677 }
678
679 if (noFieldFields.size() || meshsetFields.size()) {
680
681 auto create_meshset = [&]() {
682 EntityHandle meshset;
683 auto create = [&]() {
685 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, meshset);
686 {
687 ParallelComm *pcomm =
688 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
689 Tag part_tag = pcomm->part_tag();
690 int rank = pcomm->rank();
691 CHKERR m_field.get_moab().tag_set_data(part_tag, &meshset, 1,
692 &rank);
693 CHKERR m_field.getInterface<BitRefManager>()->setBitLevelToMeshset(
694 meshset, BitRefLevel().set());
695 }
697 };
698 CHK_THROW_MESSAGE(create(), "Failed to create meshset");
699 return meshset;
700 };
701
702 auto fe_meshset = create_meshset();
703 CHKERR m_field.add_ents_to_finite_element_by_MESHSET(fe_meshset,
704 meshsetFE, false);
705 // Add entities to the meshset, that will make DOFs on those entities
706 // accessible on meshset finite element
707 for (auto ents : meshsetFiniteElementEntities) {
708 CHKERR m_field.get_moab().add_entities(fe_meshset, ents);
709 }
710 CHKERR m_field.build_finite_elements(meshsetFE);
711
712 }
713
714 for (std::vector<std::string>::iterator fit = otherFEs.begin();
715 fit != otherFEs.end(); ++fit) {
716 CHKERR m_field.build_finite_elements(*fit);
717 }
718 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
720}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MYPCOMM_INDEX
default communicator number PCOMM
std::vector< std::string > otherFEs
Other finite elements.
Definition Simple.hpp:602
std::string meshsetFE
meshset finite element
Definition Simple.hpp:600
std::vector< Range > meshsetFiniteElementEntities
Meshset element entities.
Definition Simple.hpp:604

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

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

Definition at line 722 of file Simple.cpp.

722 {
723 Interface &m_field = cOre;
725 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
726 CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
727 // Set problem by the DOFs on the fields rather that by adding DOFs on the
728 // elements
729 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
731 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
732 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
734}
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition DMMoFEM.cpp:1297
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition Simple.hpp:610

◆ createBoundaryMeshset()

MoFEMErrorCode MoFEM::Simple::createBoundaryMeshset ( )
private

Definition at line 844 of file Simple.cpp.

844 {
845 Interface &m_field = cOre;
847 ParallelComm *pcomm =
848 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
849
850 auto get_skin = [&](auto meshset) {
851 // filter not owned entities, those are not on boundary
852
853 Range domain_ents;
854 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
855 domain_ents, true);
856 CHKERR pcomm->filter_pstatus(domain_ents,
857 PSTATUS_SHARED | PSTATUS_MULTISHARED,
858 PSTATUS_NOT, -1, nullptr);
859
860 Skinner skin(&m_field.get_moab());
861 Range domain_skin;
862 CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
863 CHKERR pcomm->filter_pstatus(domain_skin,
864 PSTATUS_SHARED | PSTATUS_MULTISHARED,
865 PSTATUS_NOT, -1, nullptr);
866 return domain_skin;
867 };
868
869 auto create_boundary_meshset = [&](auto &&domain_skin) {
871 // create boundary meshset
872 if (boundaryMeshset != 0) {
874 }
875 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
876 CHKERR m_field.get_moab().add_entities(boundaryMeshset, domain_skin);
877 for (int dd = 0; dd != dIm - 1; dd++) {
878 Range adj;
879 CHKERR m_field.get_moab().get_adjacencies(domain_skin, dd, false, adj,
880 moab::Interface::UNION);
881 CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
882 }
884 };
885
886 CHKERR create_boundary_meshset(get_skin(meshSet));
887
889}
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 891 of file Simple.cpp.

891 {
892 Interface &m_field = cOre;
894
895 ParallelComm *pcomm =
896 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
897
898 auto create_skeleton_meshset = [&](auto meshset) {
900 // create boundary meshset
901 if (skeletonMeshset != 0) {
903 }
904 Range boundary_ents, skeleton_ents;
905 CHKERR m_field.get_moab().get_entities_by_dimension(boundaryMeshset,
906 dIm - 1, boundary_ents);
907 Range domain_ents;
908 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
909 domain_ents, true);
910 CHKERR m_field.get_moab().get_adjacencies(
911 domain_ents, dIm - 1, true, skeleton_ents, moab::Interface::UNION);
912 skeleton_ents = subtract(skeleton_ents, boundary_ents);
913 CHKERR pcomm->filter_pstatus(skeleton_ents, PSTATUS_NOT_OWNED, PSTATUS_NOT,
914 -1, nullptr);
915 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, skeletonMeshset);
916 CHKERR m_field.get_moab().add_entities(skeletonMeshset, skeleton_ents);
918 };
919
920 CHKERR create_skeleton_meshset(meshSet);
921
923}

◆ defaultLoadFileFunc()

static MoFEMErrorCode MoFEM::Simple::defaultLoadFileFunc ( Interface & m_field,
const char * file_name,
const char * options )
inlinestatic

Set function for loading mesh file.

Parameters
m_fieldinterface
file_namefile name
optionsloader options
Returns
error code

Definition at line 51 of file Simple.hpp.

53 {
54 // Default behavior using the member moab.load_file
55 return m_field.get_moab().load_file(file_name, 0, options);
56 }

◆ defineFiniteElements()

MoFEMErrorCode MoFEM::Simple::defineFiniteElements ( )

Define finite elements.

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

Definition at line 471 of file Simple.cpp.

471 {
472 Interface &m_field = cOre;
474
475 auto clear_rows_and_cols = [&](auto &fe_name) {
477 auto fe_ptr = m_field.get_finite_elements();
478 auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
479 ->get<FiniteElement_name_mi_tag>();
480 auto it_fe = fe_by_name.find(fe_name);
481 if (it_fe != fe_by_name.end()) {
482
483 if (!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
484 SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
485 "modification unsuccessful");
486
487 if (!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
488 SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
489 "modification unsuccessful");
490 }
492 };
493 CHKERR clear_rows_and_cols(domainFE);
494 CHKERR clear_rows_and_cols(boundaryFE);
495 CHKERR clear_rows_and_cols(skeletonFE);
496
497 // Define finite elements
498 CHKERR m_field.add_finite_element(domainFE, MF_ZERO);
499
500 auto add_fields = [&](auto &fe_name, auto &fields) {
502 for (auto &field : fields) {
503 CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
504 CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
505 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
506 }
508 };
509
510 auto add_data_fields = [&](auto &fe_name, auto &fields) {
512 for (auto &field : fields)
513 CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
515 };
516
517 CHKERR add_fields(domainFE, domainFields);
518 CHKERR add_data_fields(domainFE, dataFields);
519 CHKERR add_fields(domainFE, noFieldFields);
520 CHKERR add_data_fields(domainFE, noFieldDataFields);
521
522 if (addBoundaryFE || !boundaryFields.empty()) {
523 CHKERR m_field.add_finite_element(boundaryFE, MF_ZERO);
524 CHKERR add_fields(boundaryFE, domainFields);
525 if (!boundaryFields.empty())
526 CHKERR add_fields(boundaryFE, boundaryFields);
527 CHKERR add_data_fields(boundaryFE, dataFields);
528 CHKERR add_data_fields(boundaryFE, noFieldDataFields);
529 CHKERR add_fields(boundaryFE, noFieldFields);
530 }
531 if (addSkeletonFE || !skeletonFields.empty()) {
532 CHKERR m_field.add_finite_element(skeletonFE, MF_ZERO);
533 CHKERR add_fields(skeletonFE, domainFields);
534 if (!skeletonFields.empty())
535 CHKERR add_fields(skeletonFE, skeletonFields);
536 CHKERR add_data_fields(skeletonFE, dataFields);
537 CHKERR add_data_fields(skeletonFE, noFieldDataFields);
538 CHKERR add_fields(skeletonFE, noFieldFields);
539 }
540
541 if (noFieldFields.size() || meshsetFields.size()) {
542 CHKERR m_field.add_finite_element(meshsetFE, MF_ZERO);
543 CHKERR add_fields(meshsetFE, meshsetFields);
544 CHKERR add_fields(meshsetFE, noFieldFields);
545 CHKERR add_data_fields(meshsetFE, noFieldDataFields);
546 }
547
549}
@ MF_ZERO
@ 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 551 of file Simple.cpp.

551 {
552 Interface &m_field = cOre;
554 // Create dm instance
555 dM = createDM(m_field.get_comm(), "DMMOFEM");
556 // set dm data structure which created mofem data structures
559 CHKERR DMSetFromOptions(dM);
561 if (addBoundaryFE || !boundaryFields.empty()) {
563 }
564 if (addSkeletonFE || !skeletonFields.empty()) {
566 }
567 if (noFieldFields.size() || meshsetFields.size()) {
569 }
571 CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
573}
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition DMMoFEM.cpp:1113
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition DMMoFEM.cpp:488
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:114
auto createDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.

◆ deleteDM()

MoFEMErrorCode MoFEM::Simple::deleteDM ( )

Delete dm.

Returns
MoFEMErrorCode

Definition at line 812 of file Simple.cpp.

812 {
814 dM.reset();
816}

◆ deleteFiniteElements()

MoFEMErrorCode MoFEM::Simple::deleteFiniteElements ( )

Delete finite elements.

Returns
MoFEMErrorCode

Definition at line 823 of file Simple.cpp.

823 {
824 Interface &m_field = cOre;
826 for (auto fe : {domainFE, boundaryFE, skeletonFE}) {
827 if (m_field.check_finite_element(fe)) {
828 CHKERR m_field.delete_finite_element(fe);
829 }
830 }
832}

◆ exchangeGhostCells()

MoFEMErrorCode MoFEM::Simple::exchangeGhostCells ( )
private

get all adjacent ghosted entities

Definition at line 925 of file Simple.cpp.

925 {
926 Interface &m_field = cOre;
928 MOFEM_LOG("WORLD", Sev::verbose) << "Exchange ghost cells";
929
930 ParallelComm *pcomm =
931 ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
932 if (pcomm == NULL)
933 pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
934
935 CHKERR pcomm->exchange_ghost_cells(getDim(), getDim() - 1, 1,
936 3 /**get all adjacent ghosted entities */,
937 true, true, meshSet ? &meshSet : nullptr);
938
939 Range shared;
940 CHKERR m_field.get_moab().get_entities_by_dimension(0, dIm, shared);
941 for (auto d = dIm - 1; d >= 0; --d) {
942 CHKERR m_field.get_moab().get_adjacencies(shared, d, false, shared,
943 moab::Interface::UNION);
944 }
945 CHKERR pcomm->filter_pstatus(shared, PSTATUS_SHARED | PSTATUS_MULTISHARED,
946 PSTATUS_OR, -1, &shared);
947 Tag part_tag = pcomm->part_tag();
948 CHKERR pcomm->exchange_tags(part_tag, shared);
949 CHKERR m_field.getInterface<CommInterface>()->resolveParentEntities(shared,
950 VERBOSE);
951
953}
@ VERBOSE
int getDim() const
Get the problem dimension.
Definition Simple.hpp:333

◆ 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
poisson_2d_dis_galerkin.cpp, and simple_l2_only.cpp.

Definition at line 504 of file Simple.hpp.

504{ 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, poisson_2d_dis_galerkin.cpp, and simple_l2_only.cpp.

Definition at line 493 of file Simple.hpp.

493{ return addSkeletonFE; }

◆ getBitAdjEnt()

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

bit ref level for parent

Returns
auto&

Definition at line 535 of file Simple.hpp.

535{ return bitAdjEnt; }

◆ getBitAdjEntMask()

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

bit ref level for parent parent

Returns
auto&

Definition at line 542 of file Simple.hpp.

542{ return bitAdjEntMask; }

◆ getBitAdjParent()

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

bit ref level for parent

Returns
auto&

Definition at line 521 of file Simple.hpp.

521{ return bitAdjParent; }

◆ getBitAdjParentMask()

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

bit ref level for parent parent

Returns
auto&

Definition at line 528 of file Simple.hpp.

528{ return bitAdjParentMask; }

◆ getBitRefLevel()

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

Get the BitRefLevel.

Returns
BitRefLevel
Examples
child_and_parent.cpp, hanging_node_approx.cpp, level_set.cpp, mixed_poisson.cpp, and plot_base.cpp.

Definition at line 375 of file Simple.hpp.

375{ return bitLevel; }

◆ getBitRefLevelMask()

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

Get the BitRefLevel.

Returns
BitRefLevel
Examples
child_and_parent.cpp, and hanging_node_approx.cpp.

Definition at line 382 of file Simple.hpp.

382{ return bitLevelMask; }

◆ getBoundaryFEName() [1/2]

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

Get the Boundary FE Name.

Returns
std::string&

Definition at line 424 of file Simple.hpp.

424{ 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 396 of file Simple.hpp.

396{ return boundaryFE; }

◆ getBoundaryMeshSet()

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

Get the BoundaryMeshSet object.

Returns
EntityHandle&
Examples
hanging_node_approx.cpp.

Definition at line 361 of file Simple.hpp.

361{ 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
child_and_parent.cpp, field_evaluator.cpp, and hanging_node_approx.cpp.

Definition at line 333 of file Simple.hpp.

333{ 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 323 of file Simple.hpp.

323{ return dM; }

◆ getDM() [2/2]

◆ getDomainFEName() [1/2]

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

Get the Domain FE Name.

Returns
std::string&

Definition at line 417 of file Simple.hpp.

417{ 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 347 of file Simple.hpp.

347{ return meshSet; }

◆ getMeshset()

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

Get the MeshSet object.

Returns
EntityHandle&

Definition at line 354 of file Simple.hpp.

354{ return meshSet; }

◆ getMeshsetFEName()

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

Get the Skeleton FE Name.

Returns
std::string&
Examples
simple_interface.cpp.

Definition at line 438 of file Simple.hpp.

438{ return meshsetFE; }

◆ getMeshsetFiniteElementEntities()

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

Range of entities added to meshset finite element.

Meshset element have not defined integration rule, if you like to integrate you have to add mesh entity to meshset finite element. In principle meshset element become standard FE.

Returns
std::vector<std::string>&
Examples
simple_interface.cpp, and simple_l2_only.cpp.

Definition at line 449 of file Simple.hpp.

449 {
451 }

◆ getOptions()

MoFEMErrorCode MoFEM::Simple::getOptions ( )

get options

Returns
error code
Examples
adolc_plasticity.cpp, 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, dg_projection.cpp, dynamic_first_order_con_law.cpp, eigen_elastic.cpp, field_evaluator.cpp, hanging_node_approx.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_method.cpp, helmholtz.cpp, higher_derivatives.cpp, level_set.cpp, loop_entities.cpp, mesh_smoothing.cpp, mixed_poisson.cpp, mofem/tutorials/vec-7/adjoint.cpp, operators_tests.cpp, phase.cpp, plastic.cpp, plate.cpp, plot_base.cpp, poisson_2d_dis_galerkin.cpp, poisson_2d_homogeneous.cpp, reaction_diffusion.cpp, schur_test_diag_mat.cpp, seepage.cpp, shallow_wave.cpp, simple_elasticity.cpp, simple_interface.cpp, simple_l2_only.cpp, tensor_divergence_operator.cpp, test_broken_space.cpp, test_cache_on_entities.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and thermo_elastic.cpp.

Definition at line 180 of file Simple.cpp.

180 {
181 PetscBool flg = PETSC_TRUE;
183 PetscOptionsBegin(PETSC_COMM_WORLD, "", "Simple interface options",
184 "none");
185 ierr = PetscOptionsString("-file_name", "file name", "", "mesh.h5m",
186 meshFileName, 255, &flg);
187 PetscOptionsEnd();
189}
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr

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

468{ return otherFEs; }

◆ getParentAdjacencies()

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

Get the addParentAdjacencies.

If set true add parent adjacencies

Returns
true
false
Examples
free_surface.cpp, and hanging_node_approx.cpp.

Definition at line 514 of file Simple.hpp.

514{ return addParentAdjacencies; }

◆ getProblemName() [1/2]

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

Get the Problem Name.

Returns
std::string&

Definition at line 458 of file Simple.hpp.

458{ return nameOfProblem; }

◆ getProblemName() [2/2]

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

◆ getSkeletonFEName() [1/2]

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

Get the Skeleton FE Name.

Returns
std::string&

Definition at line 431 of file Simple.hpp.

431{ 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 368 of file Simple.hpp.

368{ 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 247 of file Simple.cpp.

247 {
249 Interface &m_field = cOre;
250 if (m_field.get_comm_size() == 1)
252 else
253 CHKERR loadFile("PARALLEL=READ_PART;"
254 "PARALLEL_RESOLVE_SHARED_ENTS;"
255 "PARTITION=PARALLEL_PARTITION;",
258}
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
Definition Simple.cpp:191

◆ loadFile() [2/2]

MoFEMErrorCode MoFEM::Simple::loadFile ( const std::string options,
const std::string mesh_file_name,
LoadFileFunc loadFunc = defaultLoadFileFunc )

Load mesh file.

Parameters
optionsfile load options
mesh_file_namefile name if not set default or set by command line
loadFuncfunction for loading mesh file 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
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, dg_projection.cpp, field_evaluator.cpp, hanging_node_approx.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_method.cpp, helmholtz.cpp, higher_derivatives.cpp, loop_entities.cpp, mesh_smoothing.cpp, mixed_poisson.cpp, phase.cpp, plot_base.cpp, poisson_2d_dis_galerkin.cpp, poisson_2d_homogeneous.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 191 of file Simple.cpp.

193 {
194 Interface &m_field = cOre;
196 PetscLogEventBegin(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
197
198 if (!mesh_file_name.empty())
199 strcpy(meshFileName, mesh_file_name.c_str());
200
201 // Invoke the boost function, passing in required arguments
202 CHKERR loadFunc(m_field, meshFileName, options.c_str());
203 CHKERR m_field.rebuild_database();
204
205 // determine problem dimension
206 if (dIm == -1) {
207 int nb_ents_3d;
208 CHKERR m_field.get_moab().get_number_entities_by_dimension(
209 meshSet, 3, nb_ents_3d, true);
210 if (nb_ents_3d > 0) {
211 dIm = 3;
212 } else {
213 int nb_ents_2d;
214 CHKERR m_field.get_moab().get_number_entities_by_dimension(
215 meshSet, 2, nb_ents_2d, true);
216 if (nb_ents_2d > 0) {
217 dIm = 2;
218 } else {
219 dIm = 1;
220 }
221 }
222 }
223
224 if (!boundaryMeshset)
226 if (!skeletonMeshset)
228 if (addSkeletonFE)
230
231 if (bitLevel.any()) {
232 Range ents;
233 CHKERR m_field.get_moab().get_entities_by_dimension(meshSet, dIm, ents,
234 true);
235 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(ents, bitLevel,
236 false);
237 } else {
238 MOFEM_LOG("WORLD", Sev::warning) << "BitRefLevel is none and not set";
239 CHKERR m_field.getInterface<BitRefManager>()->addToDatabaseBitRefLevelByDim(
240 dIm, BitRefLevel().set(), BitRefLevel().set());
241 }
242
243 PetscLogEventEnd(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
245}
char mesh_file_name[255]
MoFEMErrorCode exchangeGhostCells()
Definition Simple.cpp:925
MoFEMErrorCode createBoundaryMeshset()
Definition Simple.cpp:844
MoFEMErrorCode createSkeletonMeshset()
Definition Simple.cpp:891

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

443 {
445
446 auto remove_field_from_list = [&](auto &vec) {
447 auto it = std::find(vec.begin(), vec.end(), name);
448 if (it != vec.end())
449 vec.erase(it);
450 };
451
452 remove_field_from_list(boundaryFields);
453
455}

◆ removeDomainField()

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

Remove field form domain.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 429 of file Simple.cpp.

429 {
431
432 auto remove_field_from_list = [&](auto &vec) {
433 auto it = std::find(vec.begin(), vec.end(), name);
434 if (it != vec.end())
435 vec.erase(it);
436 };
437
438 remove_field_from_list(domainFields);
439
441}

◆ removeSkeletonField()

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

Remove field form skeleton.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 457 of file Simple.cpp.

457 {
459
460 auto remove_field_from_list = [&](auto &vec) {
461 auto it = std::find(vec.begin(), vec.end(), name);
462 if (it != vec.end())
463 vec.erase(it);
464 };
465
466 remove_field_from_list(skeletonFields);
467
469}

◆ 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
Examples
child_and_parent.cpp.

Definition at line 762 of file Simple.cpp.

762 {
763 Interface &m_field = cOre;
765 PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
766
767 if (!only_dm) {
770 if (addSkeletonFE || !skeletonFields.empty()) {
772 if (addBoundaryFE || !boundaryFields.empty()) {
774 }
775 }
779 }
780
781 CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
782
783 const Problem *problem_ptr;
784 CHKERR DMMoFEMGetProblemPtr(dM, &problem_ptr);
785 const auto problem_name = problem_ptr->getName();
786 CHKERR m_field.modify_problem_ref_level_set_bit(problem_name, bitLevel);
787 CHKERR m_field.modify_problem_mask_ref_level_set_bit(problem_name,
789
790 // Set problem by the DOFs on the fields rather that by adding DOFs on the
791 // elements
792 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
794 m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
795
796 PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
798}
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition DMMoFEM.cpp:422
MoFEMErrorCode setParentAdjacency()
Definition Simple.cpp:67
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition Simple.cpp:471
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition Simple.cpp:660
MoFEMErrorCode buildFields()
Build fields.
Definition Simple.cpp:584
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
Examples
plot_base.cpp.

Definition at line 340 of file Simple.hpp.

340{ dIm = dim; };

◆ setFieldOrder()

◆ 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
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}
MoFEMErrorCode modify_finite_element_adjacency_table(const std::string &fe_name, const EntityType type, ElementAdjacencyFunct function)
modify finite element table, only for advanced user
Definition FECore.cpp:117
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 2 > > parentAdjFunctionDim2
Definition Simple.hpp:624
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
Definition Simple.hpp:622

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

◆ 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.

◆ 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>>(
29 fe_name, MBEDGE, *parentAdjSkeletonFunctionDim1);
30
32}
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 1 > > parentAdjSkeletonFunctionDim1
Definition Simple.hpp:631

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

◆ 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)

Setup problem.

Returns
error code
Examples
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, dg_projection.cpp, field_evaluator.cpp, hanging_node_approx.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_method.cpp, helmholtz.cpp, higher_derivatives.cpp, loop_entities.cpp, mixed_poisson.cpp, phase.cpp, plot_base.cpp, poisson_2d_dis_galerkin.cpp, poisson_2d_homogeneous.cpp, reaction_diffusion.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 736 of file Simple.cpp.

736 {
738
739 PetscLogEventBegin(MOFEM_EVENT_SimpleSetUP, 0, 0, 0, 0);
740
742
743 if (addSkeletonFE || !skeletonFields.empty()) {
745 if (addBoundaryFE || !boundaryFields.empty()) {
747 }
748 }
749
752
753 CHKERR defineProblem(is_partitioned);
757
758 PetscLogEventEnd(MOFEM_EVENT_SimpleSetUP, 0, 0, 0, 0);
760}
MoFEMErrorCode buildProblem()
Build problem.
Definition Simple.cpp:722
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition Simple.cpp:551

Member Data Documentation

◆ addBoundaryFE

bool MoFEM::Simple::addBoundaryFE
private

Add boundary FE.

Definition at line 576 of file Simple.hpp.

◆ addParentAdjacencies

bool MoFEM::Simple::addParentAdjacencies
private

If set to true parent adjacencies are build.

Definition at line 577 of file Simple.hpp.

◆ addSkeletonFE

bool MoFEM::Simple::addSkeletonFE
private

Add skeleton FE.

Definition at line 575 of file Simple.hpp.

◆ bitAdjEnt

BitRefLevel MoFEM::Simple::bitAdjEnt
private

bit ref level for parent

Definition at line 581 of file Simple.hpp.

◆ bitAdjEntMask

BitRefLevel MoFEM::Simple::bitAdjEntMask
private

bit ref level for parent parent

Definition at line 582 of file Simple.hpp.

◆ bitAdjParent

BitRefLevel MoFEM::Simple::bitAdjParent
private

bit ref level for parent

Definition at line 579 of file Simple.hpp.

◆ bitAdjParentMask

BitRefLevel MoFEM::Simple::bitAdjParentMask
private

bit ref level for parent parent

Definition at line 580 of file Simple.hpp.

◆ bitLevel

BitRefLevel MoFEM::Simple::bitLevel
private

BitRefLevel of the problem.

Definition at line 561 of file Simple.hpp.

◆ bitLevelMask

BitRefLevel MoFEM::Simple::bitLevelMask
private

BitRefLevel of the problem.

Definition at line 562 of file Simple.hpp.

◆ boundaryFE

std::string MoFEM::Simple::boundaryFE = "bFE"
private

boundary finite element

Definition at line 598 of file Simple.hpp.

◆ boundaryFields

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

boundary fields

Definition at line 585 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 572 of file Simple.hpp.

◆ cOre

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

Definition at line 559 of file Simple.hpp.

◆ dataFields

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

Data fields.

Definition at line 587 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 607 of file Simple.hpp.

◆ dM

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

Discrete manager (interface to PETSc/MoFEM functions)

Definition at line 610 of file Simple.hpp.

◆ domainFE

std::string MoFEM::Simple::domainFE = "dFE"
private

domain finite element

Definition at line 597 of file Simple.hpp.

◆ domainFields

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

domain fields

Definition at line 584 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 593 of file Simple.hpp.

◆ meshFileName

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

mesh file name

Definition at line 606 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 571 of file Simple.hpp.

◆ meshsetFE

std::string MoFEM::Simple::meshsetFE = "mFE"
private

meshset finite element

Definition at line 600 of file Simple.hpp.

◆ meshsetFields

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

meshset fields

Definition at line 588 of file Simple.hpp.

◆ meshsetFiniteElementEntities

std::vector<Range> MoFEM::Simple::meshsetFiniteElementEntities
private

Meshset element entities.

Definition at line 604 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 566 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 567 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 568 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 569 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 565 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleSetUP

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleSetUP
private

Definition at line 564 of file Simple.hpp.

◆ nameOfProblem

std::string MoFEM::Simple::nameOfProblem = "SimpleProblem"
private

problem name

Definition at line 596 of file Simple.hpp.

◆ noFieldDataFields

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

NOFIELD field name.

Definition at line 590 of file Simple.hpp.

◆ noFieldFields

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

NOFIELD field name.

Definition at line 589 of file Simple.hpp.

◆ otherFEs

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

Other finite elements.

Definition at line 602 of file Simple.hpp.

◆ parentAdjFunctionDim1

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

Definition at line 626 of file Simple.hpp.

◆ parentAdjFunctionDim2

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

Definition at line 624 of file Simple.hpp.

◆ parentAdjFunctionDim3

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

Definition at line 622 of file Simple.hpp.

◆ parentAdjSkeletonFunctionDim1

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

Definition at line 631 of file Simple.hpp.

◆ parentAdjSkeletonFunctionDim2

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

Definition at line 629 of file Simple.hpp.

◆ skeletonFE

std::string MoFEM::Simple::skeletonFE = "sFE"
private

skeleton finite element

Definition at line 599 of file Simple.hpp.

◆ skeletonFields

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

fields on the skeleton

Definition at line 586 of file Simple.hpp.

◆ skeletonMeshset

EntityHandle MoFEM::Simple::skeletonMeshset
private

skeleton meshset with boundary

Definition at line 573 of file Simple.hpp.


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