v0.14.0
Public Types | Public Member Functions | Static 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 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 More...
 
MoFEMErrorCode loadFile (const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
 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 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. 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 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. 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 & getMeshsetFEName ()
 Get the Skeleton FE Name. More...
 
auto & getMeshsetFiniteElementEntities ()
 Get the Domain Fields. 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
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface reference 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
 

Static Public Member Functions

static MoFEMErrorCode defaultLoadFileFunc (Interface &m_field, const char *file_name, const char *options)
 Set function for loading mesh file. More...
 
- 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...
 

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 > meshsetFields
 meshset 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 = "SimpleProblem"
 problem name More...
 
std::string domainFE = "dFE"
 domain finite element More...
 
std::string boundaryFE = "bFE"
 boundary finite element More...
 
std::string skeletonFE = "sFE"
 skeleton finite element More...
 
std::string meshsetFE = "mFE"
 meshset finite element More...
 
std::vector< std::string > otherFEs
 Other finite elements. More...
 
std::vector< RangemeshsetFiniteElementEntities
 Meshset element entities. 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
 

Detailed Description

Simple interface for fast problem set-up.

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, free_surface.cpp, hanging_node_approx.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_check_approx_in_3d.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, 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, reaction_diffusion.cpp, scalar_check_approximation.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, thermo_elastic.cpp, and wave_equation.cpp.

Definition at line 27 of file Simple.hpp.

Member Typedef Documentation

◆ LoadFileFunc

typedef 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),
169  bitAdjEnt(BitRefLevel().set()), bitAdjEntMask(BitRefLevel().set()) {
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 }

◆ ~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, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, mesh_smoothing.cpp, and simple_interface.cpp.

Definition at line 358 of file Simple.cpp.

362  {
363  Interface &m_field = cOre;
365  CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
366  verb);
367  boundaryFields.push_back(name);
368  if (space == NOFIELD)
369  SETERRQ(
370  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
371  "NOFIELD space for boundary filed not implemented in Simple interface");
373 }

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

400  {
401 
402  Interface &m_field = cOre;
404  CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
405  verb);
406  if (space == NOFIELD)
407  noFieldDataFields.push_back(name);
408  else
409  dataFields.push_back(name);
411 }

◆ 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

Definition at line 285 of file Simple.cpp.

289  {
290 
291  Interface &m_field = cOre;
293 
294  auto get_side_map_hcurl = [&]() {
295  return std::vector<
296 
297  std::pair<EntityType,
298  std::function<MoFEMErrorCode(BaseFunction::DofsSideMap &)>
299 
300  >>{
301 
302  {MBTRI,
303  [&](BaseFunction::DofsSideMap &dofs_side_map) -> MoFEMErrorCode {
305  dofs_side_map);
306  }}
307 
308  };
309  };
310 
311  auto get_side_map_hdiv = [&]() {
312  return std::vector<
313 
314  std::pair<EntityType,
315  std::function<MoFEMErrorCode(BaseFunction::DofsSideMap &)>
316 
317  >>{
318 
319  {MBTET,
320  [&](BaseFunction::DofsSideMap &dofs_side_map) -> MoFEMErrorCode {
322  dofs_side_map);
323  }}
324 
325  };
326  };
327 
328  auto get_side_map = [&](auto space) {
329  switch (space) {
330  case HCURL:
331  return get_side_map_hcurl();
332  case HDIV:
333  return get_side_map_hdiv();
334  default:
335  MOFEM_LOG_CHANNEL("WORLD");
336  MOFEM_LOG("WORLD", Sev::warning)
337  << "Side dof map for broken space <" << space << "> not implemented";
338  }
339  return std::vector<
340 
341  std::pair<EntityType,
342  std::function<MoFEMErrorCode(BaseFunction::DofsSideMap &)>
343 
344  >>{};
345  };
346 
347  CHKERR m_field.add_broken_field(name, space, base, nb_of_coefficients,
348  get_side_map(space), tag_type, bh, verb);
349  domainFields.push_back(name);
350  if (space == NOFIELD)
351  SETERRQ(
352  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
353  "NOFIELD space for boundary filed not implemented in Simple interface");
355 }

◆ 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
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, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, mesh_smoothing.cpp, mixed_poisson.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 264 of file Simple.cpp.

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

◆ 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

Definition at line 869 of file Simple.cpp.

870  {
871  Interface &m_field = cOre;
873  CHKERR m_field.getInterface<ProblemsManager>()->addFieldToEmptyFieldBlocks(
874  getProblemName(), row_field, col_field);
876 }

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

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

418  {
419 
420  Interface &m_field = cOre;
422  CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
423  verb);
424  if (space != NOFIELD) {
425  meshsetFields.push_back(name);
426  } else {
427  noFieldFields.push_back(name);
428  }
430 }

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

380  {
381 
382  Interface &m_field = cOre;
384  CHKERR m_field.add_field(name, space, base, nb_of_coefficients, tag_type, bh,
385  verb);
386  skeletonFields.push_back(name);
387  if (space == NOFIELD)
388  SETERRQ(
389  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
390  "NOFIELD space for boundary filed not implemented in Simple interface");
391 
393 }

◆ buildFields()

MoFEMErrorCode MoFEM::Simple::buildFields ( )

Build fields.

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

Definition at line 613 of file Simple.cpp.

613  {
614  Interface &m_field = cOre;
616  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
617 
618  auto comm_interface_ptr = m_field.getInterface<CommInterface>();
619  auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
620 
621  // Add entities to the fields
622  auto add_ents_to_field = [&](auto meshset, auto dim, auto &fields) {
624  for (auto &field : fields) {
625  CHKERR m_field.add_ents_to_field_by_dim(meshset, dim, field);
626  CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
627  }
629  };
630 
631  auto make_no_field_ents = [&](auto &fields) {
633  for (auto &field : fields) {
634  std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
635  CHKERR m_field.create_vertices_and_add_to_field(field, coords.data(), 1);
636  CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
637  CHKERR bit_ref_ptr->setFieldEntitiesBitRefLevel(field, bitLevel);
638  }
640  };
641 
642  CHKERR add_ents_to_field(meshSet, dIm, domainFields);
643  CHKERR add_ents_to_field(meshSet, dIm, dataFields);
644  CHKERR add_ents_to_field(boundaryMeshset, dIm - 1, boundaryFields);
645  CHKERR add_ents_to_field(skeletonMeshset, dIm - 1, skeletonFields);
646 
647  std::set<std::string> nofield_fields;
648  for (auto &f : noFieldFields)
649  nofield_fields.insert(f);
650  for (auto &f : noFieldDataFields)
651  nofield_fields.insert(f);
652 
653  CHKERR make_no_field_ents(nofield_fields);
654 
655  // Set order
656 
657  auto get_field_ptr = [&](auto &f) { return m_field.get_field_structure(f); };
658  for (auto &t : fieldsOrder) {
659  const auto f = std::get<0>(t);
660  const auto order = std::get<1>(t);
661 
662  MOFEM_TAG_AND_LOG("WORLD", Sev::inform, "Simple")
663  << "Set order to field " << f << " order " << order;
664  MOFEM_LOG_CHANNEL("SYNC");
665  if (std::get<3>(t)) {
666  MOFEM_TAG_AND_LOG("SYNC", Sev::verbose, "Simple")
667  << "To ents: " << std::endl
668  << std::get<2>(t) << std::endl;
669  }
670  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::verbose);
671 
672  if (std::get<3>(t)) {
673 
674  CHKERR m_field.set_field_order(std::get<2>(t), f, order);
675 
676  } else {
677  auto f_ptr = get_field_ptr(f);
678 
679  if (f_ptr->getSpace() == H1) {
680  if (f_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
681  CHKERR m_field.set_field_order(meshSet, MBVERTEX, f, order);
682  } else {
683  CHKERR m_field.set_field_order(meshSet, MBVERTEX, f, 1);
684  }
685  }
686 
687  for (auto d = 1; d <= dIm; ++d) {
688  for (EntityType t = CN::TypeDimensionMap[d].first;
689  t <= CN::TypeDimensionMap[d].second; ++t) {
690  CHKERR m_field.set_field_order(meshSet, t, f, order);
691  }
692  }
693  }
694  }
695  MOFEM_LOG_CHANNEL("WORLD");
696  // Build fields
697  CHKERR m_field.build_fields();
698  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
700 }

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

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

Definition at line 702 of file Simple.cpp.

702  {
703  Interface &m_field = cOre;
705  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
706  // Add finite elements
707  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
708  true);
709  CHKERR m_field.build_finite_elements(domainFE);
710  if (addBoundaryFE || boundaryFields.size()) {
711  CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
712  boundaryFE, true);
713  CHKERR m_field.build_finite_elements(boundaryFE);
714  }
715  if (addSkeletonFE || skeletonFields.size()) {
716  CHKERR m_field.add_ents_to_finite_element_by_dim(skeletonMeshset, dIm - 1,
717  skeletonFE, true);
718  CHKERR m_field.build_finite_elements(skeletonFE);
719  }
720 
721  if (noFieldFields.size() || meshsetFields.size()) {
722  auto add_fields_ents = [&](auto list) {
724  auto fe_meshset = m_field.get_finite_element_meshset(meshsetFE);
725  for (auto &fe_ents : meshsetFiniteElementEntities) {
726  EntityHandle meshset_entity_fe;
727  CHKERR m_field.get_moab().create_meshset(MESHSET_SET,
728  meshset_entity_fe);
729  CHKERR m_field.getInterface<BitRefManager>()->setBitLevelToMeshset(
730  meshset_entity_fe, BitRefLevel().set());
731  for (auto f : list) {
732  auto field_meshset = m_field.get_field_meshset(f);
733  Range ents;
734  CHKERR m_field.get_moab().get_entities_by_handle(field_meshset, ents,
735  true);
736  CHKERR m_field.get_moab().add_entities(meshset_entity_fe, ents);
737  }
738  CHKERR m_field.get_moab().add_entities(meshset_entity_fe, fe_ents);
739  CHKERR m_field.get_moab().add_entities(fe_meshset, &meshset_entity_fe,
740  1);
741  }
743  };
744  CHKERR add_fields_ents(noFieldFields);
745  CHKERR m_field.build_finite_elements(meshsetFE);
746  }
747 
748  for (std::vector<std::string>::iterator fit = otherFEs.begin();
749  fit != otherFEs.end(); ++fit) {
750  CHKERR m_field.build_finite_elements(*fit);
751  }
752  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
754 }

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

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

Definition at line 756 of file Simple.cpp.

756  {
757  Interface &m_field = cOre;
759  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
760  CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
761  // Set problem by the DOFs on the fields rather that by adding DOFs on the
762  // elements
763  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
765  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
766  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
768 }

◆ createBoundaryMeshset()

MoFEMErrorCode MoFEM::Simple::createBoundaryMeshset ( )
private

Definition at line 878 of file Simple.cpp.

878  {
879  Interface &m_field = cOre;
881  ParallelComm *pcomm =
882  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
883 
884  auto get_skin = [&](auto meshset) {
885  // filter not owned entities, those are not on boundary
886 
887  Range domain_ents;
888  CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
889  domain_ents, true);
890  CHKERR pcomm->filter_pstatus(domain_ents,
891  PSTATUS_SHARED | PSTATUS_MULTISHARED,
892  PSTATUS_NOT, -1, nullptr);
893 
894  Skinner skin(&m_field.get_moab());
895  Range domain_skin;
896  CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
897  CHKERR pcomm->filter_pstatus(domain_skin,
898  PSTATUS_SHARED | PSTATUS_MULTISHARED,
899  PSTATUS_NOT, -1, nullptr);
900  return domain_skin;
901  };
902 
903  auto create_boundary_meshset = [&](auto &&domain_skin) {
905  // create boundary meshset
906  if (boundaryMeshset != 0) {
908  }
909  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
910  CHKERR m_field.get_moab().add_entities(boundaryMeshset, domain_skin);
911  for (int dd = 0; dd != dIm - 1; dd++) {
912  Range adj;
913  CHKERR m_field.get_moab().get_adjacencies(domain_skin, dd, false, adj,
914  moab::Interface::UNION);
915  CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
916  }
918  };
919 
920  CHKERR create_boundary_meshset(get_skin(meshSet));
921 
923 }

◆ createSkeletonMeshset()

MoFEMErrorCode MoFEM::Simple::createSkeletonMeshset ( )
private

Definition at line 925 of file Simple.cpp.

925  {
926  Interface &m_field = cOre;
928 
929  ParallelComm *pcomm =
930  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
931 
932  auto create_skeleton_meshset = [&](auto meshset) {
934  // create boundary meshset
935  if (skeletonMeshset != 0) {
937  }
938  Range boundary_ents, skeleton_ents;
939  CHKERR m_field.get_moab().get_entities_by_dimension(boundaryMeshset,
940  dIm - 1, boundary_ents);
941  Range domain_ents;
942  CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
943  domain_ents, true);
944  CHKERR m_field.get_moab().get_adjacencies(
945  domain_ents, dIm - 1, true, skeleton_ents, moab::Interface::UNION);
946  skeleton_ents = subtract(skeleton_ents, boundary_ents);
947  CHKERR pcomm->filter_pstatus(skeleton_ents, PSTATUS_NOT_OWNED, PSTATUS_NOT,
948  -1, nullptr);
949  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, skeletonMeshset);
950  CHKERR m_field.get_moab().add_entities(skeletonMeshset, skeleton_ents);
952  };
953 
954  CHKERR create_skeleton_meshset(meshSet);
955 
957 }

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

474  {
475  Interface &m_field = cOre;
477 
478  auto clear_rows_and_cols = [&](auto &fe_name) {
480  auto fe_ptr = m_field.get_finite_elements();
481  auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
482  ->get<FiniteElement_name_mi_tag>();
483  auto it_fe = fe_by_name.find(fe_name);
484  if (it_fe != fe_by_name.end()) {
485 
486  if (!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
487  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
488  "modification unsuccessful");
489 
490  if (!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
491  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
492  "modification unsuccessful");
493  }
495  };
496  CHKERR clear_rows_and_cols(domainFE);
497  CHKERR clear_rows_and_cols(boundaryFE);
498  CHKERR clear_rows_and_cols(skeletonFE);
499 
500  // Define finite elements
501  CHKERR m_field.add_finite_element(domainFE, MF_ZERO);
502 
503  auto add_fields = [&](auto &fe_name, auto &fields) {
505  for (auto &field : fields) {
506  CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
507  CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
508  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
509  }
511  };
512 
513  auto add_data_fields = [&](auto &fe_name, auto &fields) {
515  for (auto &field : fields)
516  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
518  };
519 
520  CHKERR add_fields(domainFE, domainFields);
521  CHKERR add_data_fields(domainFE, dataFields);
522  CHKERR add_fields(domainFE, noFieldFields);
523  CHKERR add_data_fields(domainFE, noFieldDataFields);
524 
525  if (addBoundaryFE || !boundaryFields.empty()) {
526  CHKERR m_field.add_finite_element(boundaryFE, MF_ZERO);
527  CHKERR add_fields(boundaryFE, domainFields);
528  if (!boundaryFields.empty())
529  CHKERR add_fields(boundaryFE, boundaryFields);
530  CHKERR add_data_fields(boundaryFE, dataFields);
531  CHKERR add_data_fields(boundaryFE, noFieldDataFields);
532  CHKERR add_fields(boundaryFE, noFieldFields);
533  }
534  if (addSkeletonFE || !skeletonFields.empty()) {
535  CHKERR m_field.add_finite_element(skeletonFE, MF_ZERO);
536  CHKERR add_fields(skeletonFE, domainFields);
537  if (!skeletonFields.empty())
538  CHKERR add_fields(skeletonFE, skeletonFields);
539  CHKERR add_data_fields(skeletonFE, dataFields);
540  CHKERR add_data_fields(skeletonFE, noFieldDataFields);
541  CHKERR add_fields(skeletonFE, noFieldFields);
542  }
543 
544  if (noFieldFields.size() || meshsetFields.size()) {
545  CHKERR m_field.add_finite_element(meshsetFE, MF_ZERO);
546 
547  auto create_meshset = [&]() {
548  EntityHandle fe_meshset;
549  auto create = [&]() {
551  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, fe_meshset);
552  {
553  EntityHandle meshset;
554  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, meshset);
555  CHKERR m_field.getInterface<BitRefManager>()->setBitLevelToMeshset(
556  meshset, BitRefLevel().set());
557  ParallelComm *pcomm =
558  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
559  Tag part_tag = pcomm->part_tag();
560  int rank = pcomm->rank();
561  CHKERR m_field.get_moab().tag_set_data(part_tag, &meshset, 1, &rank);
562  }
564  };
565  CHK_THROW_MESSAGE(create(), "Failed to create meshset");
566  return fe_meshset;
567  };
568 
569  CHKERR m_field.add_ents_to_finite_element_by_MESHSET(create_meshset(),
570  meshsetFE, false);
571 
572  CHKERR add_fields(meshsetFE, meshsetFields);
573  CHKERR add_fields(meshsetFE, noFieldFields);
574  CHKERR add_data_fields(meshsetFE, noFieldDataFields);
575  }
576 
578 }

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

580  {
581  Interface &m_field = cOre;
583  // Create dm instance
584  dM = createDM(m_field.get_comm(), "DMMOFEM");
585  // set dm data structure which created mofem data structures
586  CHKERR DMMoFEMCreateMoFEM(dM, &m_field, nameOfProblem.c_str(), bitLevel,
587  bitLevelMask);
588  CHKERR DMSetFromOptions(dM);
590  if (addBoundaryFE || !boundaryFields.empty()) {
592  }
593  if (addSkeletonFE || !skeletonFields.empty()) {
595  }
596  if (noFieldFields.size() || meshsetFields.size()) {
598  }
600  CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
602 }

◆ deleteDM()

MoFEMErrorCode MoFEM::Simple::deleteDM ( )

Delete dm.

Returns
MoFEMErrorCode

Definition at line 846 of file Simple.cpp.

846  {
848  dM.reset();
850 }

◆ deleteFiniteElements()

MoFEMErrorCode MoFEM::Simple::deleteFiniteElements ( )

Delete finite elements.

Returns
MoFEMErrorCode

Definition at line 857 of file Simple.cpp.

857  {
858  Interface &m_field = cOre;
860  for (auto fe : {domainFE, boundaryFE, skeletonFE}) {
861  if (m_field.check_finite_element(fe)) {
862  CHKERR m_field.delete_finite_element(fe);
863  }
864  }
866 }

◆ exchangeGhostCells()

MoFEMErrorCode MoFEM::Simple::exchangeGhostCells ( )
private

get all adjacent ghosted entities

Definition at line 959 of file Simple.cpp.

959  {
960  Interface &m_field = cOre;
962  MOFEM_LOG("WORLD", Sev::verbose) << "Exchange ghost cells";
963 
964  ParallelComm *pcomm =
965  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
966  if (pcomm == NULL)
967  pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
968 
969  CHKERR pcomm->exchange_ghost_cells(getDim(), getDim() - 1, 1,
970  3 /**get all adjacent ghosted entities */,
971  true, true, meshSet ? &meshSet : nullptr);
972 
973  Range shared;
974  CHKERR m_field.get_moab().get_entities_by_dimension(0, dIm, shared);
975  for (auto d = dIm - 1; d >= 0; --d) {
976  CHKERR m_field.get_moab().get_adjacencies(shared, d, false, shared,
977  moab::Interface::UNION);
978  }
979  CHKERR pcomm->filter_pstatus(shared, PSTATUS_SHARED | PSTATUS_MULTISHARED,
980  PSTATUS_OR, -1, &shared);
981  Tag part_tag = pcomm->part_tag();
982  CHKERR pcomm->exchange_tags(part_tag, shared);
983  CHKERR m_field.getInterface<CommInterface>()->resolveParentEntities(shared,
984  VERBOSE);
985 
987 }

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

498 { 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 487 of file Simple.hpp.

487 { return addSkeletonFE; }

◆ getBitAdjEnt()

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

bit ref level for parent

Returns
auto&

Definition at line 529 of file Simple.hpp.

529 { return bitAdjEnt; }

◆ getBitAdjEntMask()

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

bit ref level for parent parent

Returns
auto&

Definition at line 536 of file Simple.hpp.

536 { return bitAdjEntMask; }

◆ getBitAdjParent()

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

bit ref level for parent

Returns
auto&

Definition at line 515 of file Simple.hpp.

515 { return bitAdjParent; }

◆ getBitAdjParentMask()

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

bit ref level for parent parent

Returns
auto&

Definition at line 522 of file Simple.hpp.

522 { return bitAdjParentMask; }

◆ getBitRefLevel()

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

Get the BitRefLevel.

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

Definition at line 373 of file Simple.hpp.

373 { return bitLevel; }

◆ getBitRefLevelMask()

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

Get the BitRefLevel.

Returns
BitRefLevel
Examples
child_and_parent.cpp.

Definition at line 380 of file Simple.hpp.

380 { return bitLevelMask; }

◆ getBoundaryFEName() [1/2]

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

Get the Boundary FE Name.

Returns
std::string&

Definition at line 422 of file Simple.hpp.

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

394 { return boundaryFE; }

◆ getBoundaryMeshSet()

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

Get the BoundaryMeshSet object.

Returns
EntityHandle&

Definition at line 359 of file Simple.hpp.

359 { 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, and field_evaluator.cpp.

Definition at line 331 of file Simple.hpp.

331 { 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 321 of file Simple.hpp.

321 { 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 415 of file Simple.hpp.

415 { 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 345 of file Simple.hpp.

345 { return meshSet; }

◆ getMeshset()

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

Get the MeshSet object.

Returns
EntityHandle&

Definition at line 352 of file Simple.hpp.

352 { return meshSet; }

◆ getMeshsetFEName()

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

Get the Skeleton FE Name.

Returns
std::string&
Examples
simple_interface.cpp.

Definition at line 436 of file Simple.hpp.

436 { return meshsetFE; }

◆ getMeshsetFiniteElementEntities()

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

Get the Domain Fields.

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

Definition at line 443 of file Simple.hpp.

443  {
445  }

◆ getOptions()

MoFEMErrorCode MoFEM::Simple::getOptions ( )

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

462 { return otherFEs; }

◆ getParentAdjacencies()

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

Get the addParentAdjacencies.

If set true add parent adjacencies

Returns
true
false

Definition at line 508 of file Simple.hpp.

508 { return addParentAdjacencies; }

◆ getProblemName() [1/2]

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

Get the Problem Name.

Returns
std::string&

Definition at line 452 of file Simple.hpp.

452 { return nameOfProblem; }

◆ getProblemName() [2/2]

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

Get the Problem Name.

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

Definition at line 408 of file Simple.hpp.

408 { return nameOfProblem; }

◆ getSkeletonFEName() [1/2]

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

Get the Skeleton FE Name.

Returns
std::string&

Definition at line 429 of file Simple.hpp.

429 { 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 366 of file Simple.hpp.

366 { 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;",
259  meshFileName);
261 }

◆ 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, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, mesh_smoothing.cpp, mixed_poisson.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 194 of file Simple.cpp.

196  {
197  Interface &m_field = cOre;
199  PetscLogEventBegin(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
200 
201  if (!mesh_file_name.empty())
202  strcpy(meshFileName, mesh_file_name.c_str());
203 
204  // Invoke the boost function, passing in required arguments
205  CHKERR loadFunc(m_field, meshFileName, 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 }

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

◆ removeBoundaryField()

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

Remove field form boundary.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 446 of file Simple.cpp.

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

◆ removeDomainField()

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

Remove field form domain.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 432 of file Simple.cpp.

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

◆ removeSkeletonField()

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

Remove field form skeleton.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 460 of file Simple.cpp.

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

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

796  {
797  Interface &m_field = cOre;
799  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
800 
801  if (!only_dm) {
804  if (addSkeletonFE || !skeletonFields.empty()) {
806  if (addBoundaryFE || !boundaryFields.empty()) {
808  }
809  }
813  }
814 
815  CHKERR m_field.build_adjacencies(bitLevel, bitLevelMask);
816 
817  const Problem *problem_ptr;
818  CHKERR DMMoFEMGetProblemPtr(dM, &problem_ptr);
819  const auto problem_name = problem_ptr->getName();
820  CHKERR m_field.modify_problem_ref_level_set_bit(problem_name, bitLevel);
821  CHKERR m_field.modify_problem_mask_ref_level_set_bit(problem_name,
822  bitLevelMask);
823 
824  // Set problem by the DOFs on the fields rather that by adding DOFs on the
825  // elements
826  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
828  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
829 
830  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
832 }

◆ setDim()

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

Set the problem dimension.

Returns
int

Definition at line 338 of file Simple.hpp.

338 { 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 }

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

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

◆ 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");
69  return MOFEM_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 }

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

◆ 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");
17  return MOFEM_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 570 of file Simple.hpp.

◆ addParentAdjacencies

bool MoFEM::Simple::addParentAdjacencies
private

If set to true parent adjacencies are build.

Definition at line 571 of file Simple.hpp.

◆ addSkeletonFE

bool MoFEM::Simple::addSkeletonFE
private

Add skeleton FE.

Definition at line 569 of file Simple.hpp.

◆ bitAdjEnt

BitRefLevel MoFEM::Simple::bitAdjEnt
private

bit ref level for parent

Definition at line 575 of file Simple.hpp.

◆ bitAdjEntMask

BitRefLevel MoFEM::Simple::bitAdjEntMask
private

bit ref level for parent parent

Definition at line 576 of file Simple.hpp.

◆ bitAdjParent

BitRefLevel MoFEM::Simple::bitAdjParent
private

bit ref level for parent

Definition at line 573 of file Simple.hpp.

◆ bitAdjParentMask

BitRefLevel MoFEM::Simple::bitAdjParentMask
private

bit ref level for parent parent

Definition at line 574 of file Simple.hpp.

◆ bitLevel

BitRefLevel MoFEM::Simple::bitLevel
private

BitRefLevel of the problem.

Definition at line 555 of file Simple.hpp.

◆ bitLevelMask

BitRefLevel MoFEM::Simple::bitLevelMask
private

BitRefLevel of the problem.

Definition at line 556 of file Simple.hpp.

◆ boundaryFE

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

boundary finite element

Definition at line 592 of file Simple.hpp.

◆ boundaryFields

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

boundary fields

Definition at line 579 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 566 of file Simple.hpp.

◆ cOre

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

Definition at line 553 of file Simple.hpp.

◆ dataFields

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

Data fields.

Definition at line 581 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 601 of file Simple.hpp.

◆ dM

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

Discrete manager (interface to PETSc/MoFEM functions)

Definition at line 604 of file Simple.hpp.

◆ domainFE

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

domain finite element

Definition at line 591 of file Simple.hpp.

◆ domainFields

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

domain fields

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

◆ meshFileName

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

mesh file name

Definition at line 600 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 565 of file Simple.hpp.

◆ meshsetFE

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

meshset finite element

Definition at line 594 of file Simple.hpp.

◆ meshsetFields

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

meshset fields

Definition at line 582 of file Simple.hpp.

◆ meshsetFiniteElementEntities

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

Meshset element entities.

Definition at line 598 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 560 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 561 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 562 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 563 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 559 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleSetUP

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleSetUP
private

Definition at line 558 of file Simple.hpp.

◆ nameOfProblem

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

problem name

Definition at line 590 of file Simple.hpp.

◆ noFieldDataFields

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

NOFIELD field name.

Definition at line 584 of file Simple.hpp.

◆ noFieldFields

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

NOFIELD field name.

Definition at line 583 of file Simple.hpp.

◆ otherFEs

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

Other finite elements.

Definition at line 596 of file Simple.hpp.

◆ parentAdjFunctionDim1

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

Definition at line 620 of file Simple.hpp.

◆ parentAdjFunctionDim2

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

Definition at line 618 of file Simple.hpp.

◆ parentAdjFunctionDim3

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

Definition at line 616 of file Simple.hpp.

◆ parentAdjSkeletonFunctionDim1

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

Definition at line 625 of file Simple.hpp.

◆ parentAdjSkeletonFunctionDim2

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

Definition at line 623 of file Simple.hpp.

◆ skeletonFE

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

skeleton finite element

Definition at line 593 of file Simple.hpp.

◆ skeletonFields

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

fields on the skeleton

Definition at line 580 of file Simple.hpp.

◆ skeletonMeshset

EntityHandle MoFEM::Simple::skeletonMeshset
private

skeleton meshset with boundary

Definition at line 567 of file Simple.hpp.


The documentation for this struct was generated from the following files:
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::Simple::boundaryFields
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:579
MoFEM::Simple::parentAdjFunctionDim1
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 1 > > parentAdjFunctionDim1
Definition: Simple.hpp:620
MoFEM::Simple::setSkeletonAdjacency
MoFEMErrorCode setSkeletonAdjacency(int dim=-1, std::string fe_name="")
Set the skeleton adjacency object.
Definition: Simple.cpp:143
MYPCOMM_INDEX
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:228
MoFEM::Simple::otherFEs
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:596
H1
@ H1
continuous field
Definition: definitions.h:85
MoFEM::Simple::addParentAdjacencies
bool addParentAdjacencies
If set to true parent adjacencies are build.
Definition: Simple.hpp:571
MoFEM::Simple::addBoundaryFE
bool addBoundaryFE
Add boundary FE.
Definition: Simple.hpp:570
MOFEM_LOG_SEVERITY_SYNC
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
Definition: LogManager.hpp:352
MoFEM::Simple::exchangeGhostCells
MoFEMErrorCode exchangeGhostCells()
Definition: Simple.cpp:959
MoFEM::Simple::fieldsOrder
std::list< std::tuple< std::string, int, Range, bool > > fieldsOrder
Definition: Simple.hpp:587
MOFEM_LOG_CHANNEL
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:284
EntityHandle
MoFEM::Simple::Simple
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:162
MoFEM::Simple::parentAdjFunctionDim2
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 2 > > parentAdjFunctionDim2
Definition: Simple.hpp:618
MoFEM::Simple::noFieldFields
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:583
MoFEM::BaseFunction::DofsSideMap
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.
Definition: BaseFunction.hpp:73
MoFEM::Simple::meshSet
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:565
CHK_THROW_MESSAGE
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:609
MoFEM::Simple::createBoundaryMeshset
MoFEMErrorCode createBoundaryMeshset()
Definition: Simple.cpp:878
MoFEM::Simple::meshFileName
char meshFileName[255]
mesh file name
Definition: Simple.hpp:600
MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:563
MoFEM::Simple::buildProblem
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:756
MoFEM::Simple::meshsetFields
std::vector< std::string > meshsetFields
meshset fields
Definition: Simple.hpp:582
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:560
MoFEM::getPetscObject
PetscObject getPetscObject(T obj)
Definition: PetscSmartObj.hpp:9
MoFEM::Simple::loadFile
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
Definition: Simple.cpp:194
MoFEM::DMSetUp_MoFEM
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMoFEM.cpp:1303
MoFEM::Simple::dIm
int dIm
dimension of problem
Definition: Simple.hpp:601
MoFEM::Simple::addSkeletonFE
bool addSkeletonFE
Add skeleton FE.
Definition: Simple.hpp:569
MoFEM::Simple::nameOfProblem
std::string nameOfProblem
problem name
Definition: Simple.hpp:590
MoFEM::Simple::parentAdjSkeletonFunctionDim2
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 2 > > parentAdjSkeletonFunctionDim2
Definition: Simple.hpp:623
THROW_MESSAGE
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:574
MoFEM::DMMoFEMAddElement
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:497
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:559
MoFEM::Simple::bitLevelMask
BitRefLevel bitLevelMask
BitRefLevel of the problem.
Definition: Simple.hpp:556
MoFEM::Simple::dataFields
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:581
MoFEM::Simple::meshsetFiniteElementEntities
std::vector< Range > meshsetFiniteElementEntities
Meshset element entities.
Definition: Simple.hpp:598
VERBOSE
@ VERBOSE
Definition: definitions.h:222
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::Simple::bitAdjEntMask
BitRefLevel bitAdjEntMask
bit ref level for parent parent
Definition: Simple.hpp:576
MoFEM::Simple::meshsetFE
std::string meshsetFE
meshset finite element
Definition: Simple.hpp:594
MoFEM::createDM
auto createDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.
Definition: PetscSmartObj.hpp:141
MoFEM::Simple::domainFE
std::string domainFE
domain finite element
Definition: Simple.hpp:591
MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:561
MoFEM::Simple::createSkeletonMeshset
MoFEMErrorCode createSkeletonMeshset()
Definition: Simple.cpp:925
MoFEM::TriPolynomialBase::setDofsSideMap
static MoFEMErrorCode setDofsSideMap(const FieldSpace space, const FieldContinuity continuity, const FieldApproximationBase base, DofsSideMap &)
Set map of dof to side number.
Definition: TriPolynomialBase.cpp:1229
MoFEM::Simple::boundaryFE
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:592
MoFEM::Simple::parentAdjFunctionDim3
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
Definition: Simple.hpp:616
MoFEM::Simple::bitAdjParentMask
BitRefLevel bitAdjParentMask
bit ref level for parent parent
Definition: Simple.hpp:574
MoFEM::Simple::skeletonFE
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:593
get_skin
static auto get_skin(MoFEM::Interface &m_field, Range body_ents)
Definition: EshelbianPlasticity.cpp:163
MOFEM_OPERATION_UNSUCCESSFUL
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition: definitions.h:34
MoFEM::Simple::bitAdjEnt
BitRefLevel bitAdjEnt
bit ref level for parent
Definition: Simple.hpp:575
t
constexpr double t
plate stiffness
Definition: plate.cpp:58
MoFEM::DMMoFEMCreateMoFEM
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
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
AINSWORTH_BERNSTEIN_BEZIER_BASE
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:64
mesh_file_name
char mesh_file_name[255]
Definition: mesh_smoothing.cpp:22
MoFEM::Simple::buildFiniteElements
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:702
Range
MoFEM::TetPolynomialBase::setDofsSideMap
static MoFEMErrorCode setDofsSideMap(const FieldSpace space, const FieldContinuity continuity, const FieldApproximationBase base, DofsSideMap &)
Set map of dof to side number.
Definition: TetPolynomialBase.cpp:2129
FTensor::dd
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
MF_ZERO
@ MF_ZERO
Definition: definitions.h:111
MOFEM_TAG_AND_LOG
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:362
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::Simple::domainFields
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:578
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
HenckyOps::f
auto f
Definition: HenckyOps.hpp:15
DISCONTINUOUS
@ DISCONTINUOUS
Broken continuity (No effect on L2 space)
Definition: definitions.h:101
MoFEM::Simple::skeletonFields
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:580
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
HCURL
@ HCURL
field with continuous tangents
Definition: definitions.h:86
MoFEM::DMMoFEMGetProblemPtr
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMoFEM.cpp:426
MoFEM::Simple::bitLevel
BitRefLevel bitLevel
BitRefLevel of the problem.
Definition: Simple.hpp:555
MoFEM::Simple::defineFiniteElements
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:474
MoFEM::Simple::getDim
int getDim() const
Get the problem dimension.
Definition: Simple.hpp:331
MoFEM::Simple::buildFields
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:613
MoFEM::Simple::bitAdjParent
BitRefLevel bitAdjParent
bit ref level for parent
Definition: Simple.hpp:573
MoFEM::Simple::noFieldDataFields
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:584
sdf_hertz_2d_axisymm_plane.d
float d
Definition: sdf_hertz_2d_axisymm_plane.py:4
FiniteElement_multiIndex
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.
Definition: FEMultiIndices.hpp:849
MoFEM::Simple::boundaryMeshset
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:566
MoFEM::Simple::parentAdjSkeletonFunctionDim1
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 1 > > parentAdjSkeletonFunctionDim1
Definition: Simple.hpp:625
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
MoFEM::Simple::defineProblem
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:580
MoFEM::Simple::dM
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:604
MoFEM::Simple::setParentAdjacency
MoFEMErrorCode setParentAdjacency()
Definition: Simple.cpp:67
MoFEM::Simple::getProblemName
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:408
MoFEM::Simple::addFieldToEmptyFieldBlocks
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
add empty block to problem
Definition: Simple.cpp:869
MoFEM::Simple::cOre
MoFEM::Core & cOre
Definition: Simple.hpp:553
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEM::Simple::skeletonMeshset
EntityHandle skeletonMeshset
skeleton meshset with boundary
Definition: Simple.hpp:567
MoFEM::Simple::MOFEM_EVENT_SimpleSetUP
PetscLogEvent MOFEM_EVENT_SimpleSetUP
Definition: Simple.hpp:558
HDIV
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
CHKERRG
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
NOFIELD
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:84
MoFEM::DMMoFEMSetIsPartitioned
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMoFEM.cpp:1123
MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:562