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

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

#include <src/interfaces/Simple.hpp>

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

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 Simple (const MoFEM::Core &core)
 
virtual ~Simple ()=default
 
MoFEMErrorCode getOptions ()
 get options More...
 
MoFEMErrorCode loadFile (const std::string options, const std::string mesh_file_name)
 Load mesh file. More...
 
MoFEMErrorCode loadFile (const std::string mesh_file_name="")
 Load mesh file with parallel options if number of cores > 1. More...
 
MoFEMErrorCode addDomainField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
 Add field on domain. More...
 
MoFEMErrorCode addBoundaryField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
 Add field on boundary. More...
 
MoFEMErrorCode addSkeletonField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
 Add field on skeleton. More...
 
MoFEMErrorCode addDataField (const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
 Add field on domain. More...
 
MoFEMErrorCode removeDomainField (const std::string &name)
 Remove field form domain. More...
 
MoFEMErrorCode removeBoundaryField (const std::string &name)
 Remove field form boundary. More...
 
MoFEMErrorCode removeSkeletonField (const std::string &name)
 Remove field form skeleton. More...
 
MoFEMErrorCode defineFiniteElements ()
 Define finite elements. More...
 
MoFEMErrorCode defineProblem (const PetscBool is_partitioned=PETSC_TRUE)
 define problem More...
 
MoFEMErrorCode setFieldOrder (const std::string field_name, const int order, const Range *ents=NULL)
 Set field order. More...
 
MoFEMErrorCode buildFields ()
 Build fields. More...
 
MoFEMErrorCode buildFiniteElements ()
 Build finite elements. More...
 
MoFEMErrorCode setSkeletonAdjacency (int dim=-1)
 Set the skeleton adjacency object. More...
 
MoFEMErrorCode buildProblem ()
 Build problem. More...
 
MoFEMErrorCode setUp (const PetscBool is_partitioned=PETSC_TRUE)
 Setup problem. More...
 
MoFEMErrorCode reSetUp ()
 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 dimensio. More...
 
const std::string getDomainFEName () const
 Get the Domain FE Name. More...
 
const std::string getBoundaryFEName () const
 Get the Boundary FE Name. More...
 
const std::string getSkeletonFEName () const
 Get the Skeleton FE Name. More...
 
const std::string getProblemName () const
 Get the Problem Name. More...
 
std::string & getDomainFEName ()
 Get the Domain FE Name. More...
 
std::string & getBoundaryFEName ()
 Get the Boundary FE Name. More...
 
std::string & getSkeletonFEName ()
 Get the Skeleton FE Name. More...
 
std::string & getProblemName ()
 Get the Problem Name. More...
 
std::vector< std::string > & getOtherFiniteElements ()
 Get the Other Finite Elements. More...
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

Private Member Functions

template<int DIM = -1>
MoFEMErrorCode setSkeletonAdjacency ()
 
template<>
MoFEMErrorCode setSkeletonAdjacency ()
 
template<>
MoFEMErrorCode setSkeletonAdjacency ()
 
template<>
MoFEMErrorCode setSkeletonAdjacency ()
 

Private Attributes

MoFEM::CorecOre
 
const BitRefLevel bitLevel
 BitRefLevel of the probelm. More...
 
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...
 
std::vector< std::string > domainFields
 domain fields More...
 
std::vector< std::string > boundaryFields
 boundary fields More...
 
std::vector< std::string > skeletonFields
 fields on the skeleton More...
 
std::vector< std::string > dataFields
 Data fields. More...
 
std::vector< std::string > noFieldFields
 NOFIELD field name. More...
 
std::vector< std::string > noFieldDataFields
 NOFIELD field name. More...
 
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
 fields order More...
 
std::string nameOfProblem
 problem name More...
 
std::string domainFE
 domain finite element More...
 
std::string boundaryFE
 boundary finite element More...
 
std::string skeletonFE
 skeleton finite element More...
 
std::vector< std::string > otherFEs
 Other finite elements. More...
 
char meshFileName [255]
 mesh file name More...
 
int dIm
 dimension of problem More...
 
SmartPetscObj< DM > dM
 Discrete manager (interface to PETSc/MoFEM functions) More...
 

Additional Inherited Members

- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Simple interface for fast problem set-up.

Examples
boundary_marker.cpp, continuity_check_on_skeleton_with_simple_2d.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, hello_world.cpp, lesson1_moment_of_inertia.cpp, lesson2_approximaton.cpp, lesson3_poisson.cpp, lesson4_elastic.cpp, lesson5_helmholtz.cpp, lesson6_radiation.cpp, lesson7_plastic.cpp, lesson8_contact.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion_equation.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 36 of file Simple.hpp.

Constructor & Destructor Documentation

◆ Simple()

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

Definition at line 185 of file Simple.cpp.

186  : cOre(const_cast<Core &>(core)), bitLevel(BitRefLevel().set(0)),
187  meshSet(0), boundaryMeshset(0), nameOfProblem("SimpleProblem"),
188  domainFE("dFE"), boundaryFE("bFE"), skeletonFE("sFE"), dIm(-1) {
189  PetscLogEventRegister("LoadMesh", 0, &MOFEM_EVENT_SimpleLoadMesh);
190  PetscLogEventRegister("buildFields", 0, &MOFEM_EVENT_SimpleBuildFields);
191  PetscLogEventRegister("buildFiniteElements", 0,
193  PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleBuildProblem);
194  PetscLogEventRegister("SimpleKSPSolve", 0, &MOFEM_EVENT_SimpleKSPSolve);
195  strcpy(meshFileName, "mesh.h5m");
196 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:366
std::string domainFE
domain finite element
Definition: Simple.hpp:365
MoFEM::Core & cOre
Definition: Simple.hpp:343
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:353
char meshFileName[255]
mesh file name
Definition: Simple.hpp:371
int dIm
dimension of problem
Definition: Simple.hpp:372
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:347
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:348
std::string nameOfProblem
problem name
Definition: Simple.hpp:364
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:50
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:354
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:349
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:350
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:345
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:351

◆ ~Simple()

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

Member Function Documentation

◆ addBoundaryField()

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

Add field on boundary.

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

Definition at line 287 of file Simple.cpp.

291  {
292  Interface &m_field = cOre;
294  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
295  verb);
296  boundaryFields.push_back(name);
297  if (space == NOFIELD)
298  SETERRQ(
299  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
300  "NOFIELD space for boundary filed not implemented in Simple interface");
302 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:356
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ addDataField()

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

Add field on domain.

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

Definition at line 324 of file Simple.cpp.

328  {
329 
330  Interface &m_field = cOre;
332  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
333  verb);
334  if (space == NOFIELD)
335  noFieldDataFields.push_back(name);
336  else
337  dataFields.push_back(name);
339 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:360
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:358
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ addDomainField()

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

Add field on domain.

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

Definition at line 269 of file Simple.cpp.

273  {
274 
275  Interface &m_field = cOre;
277  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
278  verb);
279  if (space == NOFIELD)
280  noFieldFields.push_back(name);
281  else
282  domainFields.push_back(name);
284 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:359
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:355
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ addSkeletonField()

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

Add field on skeleton.

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

Definition at line 305 of file Simple.cpp.

309  {
310 
311  Interface &m_field = cOre;
313  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
314  verb);
315  skeletonFields.push_back(name);
316  if (space == NOFIELD)
317  SETERRQ(
318  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
319  "NOFIELD space for boundary filed not implemented in Simple interface");
321 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ buildFields()

MoFEMErrorCode MoFEM::Simple::buildFields ( )

Build fields.

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

Definition at line 490 of file Simple.cpp.

490  {
491  Interface &m_field = cOre;
493  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
494 
495  auto get_skin = [&](auto meshset) {
496  Range domain_ents;
497  CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dIm,
498  domain_ents, true);
499  Skinner skin(&m_field.get_moab());
500  Range domain_skin;
501  CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
502  // filter not owned entities, those are not on boundary
503  ParallelComm *pcomm =
504  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
505  Range proc_domain_skin;
506  CHKERR pcomm->filter_pstatus(domain_skin,
507  PSTATUS_SHARED | PSTATUS_MULTISHARED,
508  PSTATUS_NOT, -1, &proc_domain_skin);
509  return proc_domain_skin;
510  };
511 
512  auto create_boundary_meshset = [&](auto &&proc_domain_skin) {
514  // create boundary meshset
515  if (boundaryMeshset != 0) {
516  CHKERR m_field.get_moab().delete_entities(&boundaryMeshset, 1);
517  }
518  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
519  CHKERR m_field.get_moab().add_entities(boundaryMeshset, proc_domain_skin);
520  for (int dd = 0; dd != dIm - 1; dd++) {
521  Range adj;
522  CHKERR m_field.get_moab().get_adjacencies(proc_domain_skin, dd, false,
523  adj, moab::Interface::UNION);
524  CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
525  }
527  };
528 
529  CHKERR create_boundary_meshset(get_skin(meshSet));
530 
531  auto comm_interface_ptr = m_field.getInterface<CommInterface>();
532  auto bit_ref_ptr = m_field.getInterface<BitRefManager>();
533 
534  // Add entities to the fields
535  auto add_ents_to_field = [&](auto meshset, auto dim, auto &fields) {
537  for (auto &field : fields) {
538  CHKERR m_field.add_ents_to_field_by_dim(meshset, dim, field);
539  CHKERR comm_interface_ptr->synchroniseFieldEntities(field, 0);
540  }
542  };
543 
544  auto make_no_field_ents = [&](auto &fields) {
546  for (auto &field : fields) {
547  std::array<double, 6> coords = {0, 0, 0, 0, 0, 0};
548  CHKERR m_field.create_vertices_and_add_to_field(field, coords.data(), 2);
549  CHKERR comm_interface_ptr->makeFieldEntitiesMultishared(field, 0);
550  CHKERR bit_ref_ptr->setFieldEntitiesBitRefLevel(field, bitLevel);
551  }
553  };
554 
555  CHKERR add_ents_to_field(meshSet, dIm, domainFields);
556  CHKERR add_ents_to_field(meshSet, dIm, dataFields);
557  CHKERR add_ents_to_field(boundaryMeshset, dIm - 1, boundaryFields);
558  CHKERR add_ents_to_field(meshSet, dIm - 1, skeletonFields);
559  CHKERR make_no_field_ents(noFieldFields);
560  CHKERR make_no_field_ents(noFieldDataFields);
561 
562  // Set order
563  auto set_order = [&](auto meshset, auto dim, auto &fields) {
565  for (auto &field : fields) {
566  if (fieldsOrder.find(field) == fieldsOrder.end()) {
567  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
568  "Order for field not set %s", field.c_str());
569  }
570  int dds = 0;
571  const Field *field_ptr = m_field.get_field_structure(field);
572  switch (field_ptr->getSpace()) {
573  case L2:
574  dds = dim;
575  break;
576  case HDIV:
577  dds = 2;
578  break;
579  case HCURL:
580  dds = 1;
581  break;
582  case H1:
583  dds = 1;
584  break;
585  default:
586  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
587  "Glasgow we have a problem");
588  }
589 
590  auto set_order = [&](auto field, auto &ents) {
592 
593  auto range = fieldsOrder.equal_range(field);
594  for (auto o = range.first; o != range.second; ++o) {
595  if (!o->second.second.empty())
596  ents = intersect(ents, o->second.second);
597  CHKERR m_field.set_field_order(ents, field, o->second.first);
598  }
599 
601  };
602 
603  if (field_ptr->getSpace() == H1) {
604  if (field_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
605  Range ents;
606  CHKERR m_field.get_field_entities_by_dimension(field, 0, ents);
607  CHKERR set_order(field, ents);
608  } else {
609  CHKERR m_field.set_field_order(meshSet, MBVERTEX, field, 1);
610  }
611  }
612  for (int dd = dds; dd <= dim; dd++) {
613  Range ents;
614  CHKERR m_field.get_field_entities_by_dimension(field, dd, ents);
615  CHKERR set_order(field, ents);
616  }
617  }
619  };
620 
621  CHKERR set_order(meshSet, dIm, domainFields);
622  CHKERR set_order(meshSet, dIm, dataFields);
623  CHKERR set_order(meshSet, dIm - 1, boundaryFields);
624  CHKERR set_order(meshSet, dIm - 1, skeletonFields);
625 
626  // Build fields
627  CHKERR m_field.build_fields();
628  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
630 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:353
field with continuous normal traction
Definition: definitions.h:179
int dIm
dimension of problem
Definition: Simple.hpp:372
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:356
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:359
const int dim
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:355
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:360
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:362
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:358
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:348
field with continuous tangents
Definition: definitions.h:178
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:291
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:354
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
continuous field
Definition: definitions.h:177
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:345
field with C-1 continuity
Definition: definitions.h:180

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

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

Definition at line 632 of file Simple.cpp.

632  {
633  Interface &m_field = cOre;
635  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
636  // Add finite elements
637  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
638  true);
639  CHKERR m_field.build_finite_elements(domainFE);
640  if (!boundaryFields.empty()) {
641  CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
642  boundaryFE, true);
643  CHKERR m_field.build_finite_elements(boundaryFE);
644  }
645  if (!skeletonFields.empty()) {
646  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm - 1,
647  skeletonFE, true);
648  CHKERR m_field.build_finite_elements(skeletonFE);
649  }
650  for (std::vector<std::string>::iterator fit = otherFEs.begin();
651  fit != otherFEs.end(); ++fit) {
652  CHKERR m_field.build_finite_elements(*fit);
653  }
654  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
656 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:366
std::string domainFE
domain finite element
Definition: Simple.hpp:365
MoFEM::Core & cOre
Definition: Simple.hpp:343
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:353
int dIm
dimension of problem
Definition: Simple.hpp:372
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:356
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:369
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367
#define CHKERR
Inline error check.
Definition: definitions.h:602
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:354
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:349
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

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

Definition at line 658 of file Simple.cpp.

658  {
659  Interface &m_field = cOre;
661  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
662  CHKERR m_field.build_adjacencies(bitLevel);
663  // Set problem by the DOFs on the fields rather that by adding DOFs on the
664  // elements
665  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
666  CHKERR DMSetUp(dM);
667  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
668  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
670 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:375
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:350
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:345

◆ defineFiniteElements()

MoFEMErrorCode MoFEM::Simple::defineFiniteElements ( )

Define finite elements.

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

Definition at line 387 of file Simple.cpp.

387  {
388  Interface &m_field = cOre;
390 
391  auto clear_rows_and_cols = [&](auto &fe_name) {
393  auto fe_ptr = m_field.get_finite_elements();
394  auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
395  ->get<FiniteElement_name_mi_tag>();
396  auto it_fe = fe_by_name.find(fe_name);
397  if (it_fe != fe_by_name.end()) {
398 
399  if(!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
400  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
401  "modification unsuccessful");
402 
403  if(!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
404  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
405  "modification unsuccessful");
406  }
408  };
409  CHKERR clear_rows_and_cols(domainFE);
410  CHKERR clear_rows_and_cols(boundaryFE);
411  CHKERR clear_rows_and_cols(skeletonFE);
412 
413  // Define finite elements
414  CHKERR m_field.add_finite_element(domainFE, MF_ZERO);
415 
416  auto add_fields = [&](auto &fe_name, auto &fields) {
418  for (auto &field : fields) {
419  CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
420  CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
421  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
422  }
424  };
425 
426  auto add_data_fields = [&](auto &fe_name, auto &fields) {
428  for (auto &field : fields)
429  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
431  };
432 
433  CHKERR add_fields(domainFE, domainFields);
434  CHKERR add_data_fields(domainFE, dataFields);
435  CHKERR add_fields(domainFE, noFieldFields);
436  CHKERR add_data_fields(domainFE, noFieldDataFields);
437 
438  if (!boundaryFields.empty()) {
439  CHKERR m_field.add_finite_element(boundaryFE, MF_ZERO);
440  CHKERR add_fields(boundaryFE, domainFields);
441  CHKERR add_fields(boundaryFE, boundaryFields);
442  CHKERR add_fields(boundaryFE, skeletonFields);
443  CHKERR add_data_fields(boundaryFE, dataFields);
444  CHKERR add_data_fields(boundaryFE, noFieldDataFields);
445  CHKERR add_fields(boundaryFE, noFieldFields);
446  }
447  if (!skeletonFields.empty()) {
448  CHKERR m_field.add_finite_element(skeletonFE, MF_ZERO);
449  CHKERR add_fields(skeletonFE, domainFields);
450  CHKERR add_fields(skeletonFE, boundaryFields);
451  CHKERR add_fields(skeletonFE, skeletonFields);
452  CHKERR add_data_fields(skeletonFE, dataFields);
453  CHKERR add_data_fields(skeletonFE, noFieldDataFields);
454  CHKERR add_fields(skeletonFE, noFieldFields);
455  }
457 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:366
std::string domainFE
domain finite element
Definition: Simple.hpp:365
MoFEM::Core & cOre
Definition: Simple.hpp:343
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:356
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:359
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:355
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:360
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:358
#define CHKERR
Inline error check.
Definition: definitions.h:602
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879

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

459  {
460  Interface &m_field = cOre;
462  // Create dm instance
463  dM = createSmartDM(m_field.get_comm(), "DMMOFEM");
464  // set dm data structure which created mofem data structures
465  CHKERR DMMoFEMCreateMoFEM(dM, &m_field, nameOfProblem.c_str(), bitLevel);
466  CHKERR DMSetFromOptions(dM);
468  if (!boundaryFields.empty()) {
470  }
471  if (!skeletonFields.empty()) {
473  }
474  for (std::vector<std::string>::iterator fit = otherFEs.begin();
475  fit != otherFEs.end(); ++fit) {
476  CHKERR DMMoFEMAddElement(dM, fit->c_str());
477  }
478  CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
480 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:366
std::string domainFE
domain finite element
Definition: Simple.hpp:365
MoFEM::Core & cOre
Definition: Simple.hpp:343
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:425
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:375
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:356
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:369
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:105
std::string nameOfProblem
problem name
Definition: Simple.hpp:364
#define CHKERR
Inline error check.
Definition: definitions.h:602
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
auto createSmartDM
Creates smart DM object.
Definition: AuxPETSc.hpp:174
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:345
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:982

◆ getBoundaryFEName() [1/2]

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

Get the Boundary FE Name.

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

Definition at line 288 of file Simple.hpp.

288 { return boundaryFE; }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:366

◆ getBoundaryFEName() [2/2]

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

Get the Boundary FE Name.

Returns
std::string&

Definition at line 316 of file Simple.hpp.

316 { return boundaryFE; }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:366

◆ getDim()

int MoFEM::Simple::getDim ( ) const

Get the problem dimensio.

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

Returns
int

Definition at line 274 of file Simple.hpp.

274 { return dIm; }
int dIm
dimension of problem
Definition: Simple.hpp:372

◆ getDM() [1/2]

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

Get DM.

Parameters
dmdiscrete manager
Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion_equation.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 706 of file Simple.cpp.

706  {
708  CHKERR PetscObjectReference(getPetscObject(dM.get()));
709  *dm = dM.get();
711 }
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:375
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscObject getPetscObject(T obj)
Definition: AuxPETSc.hpp:23
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ getDM() [2/2]

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

Return smart DM object.

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

Definition at line 264 of file Simple.hpp.

264 { return dM; }
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:375

◆ getDomainFEName() [1/2]

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

Get the Domain FE Name.

Returns
const std::string
Examples
field_evaluator.cpp, mesh_smoothing.cpp, reaction_diffusion_equation.cpp, simple_elasticity.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 281 of file Simple.hpp.

281 { return domainFE; }
std::string domainFE
domain finite element
Definition: Simple.hpp:365

◆ getDomainFEName() [2/2]

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

Get the Domain FE Name.

Returns
std::string&

Definition at line 309 of file Simple.hpp.

309 { return domainFE; }
std::string domainFE
domain finite element
Definition: Simple.hpp:365

◆ getOptions()

MoFEMErrorCode MoFEM::Simple::getOptions ( )

get options

Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion_equation.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 198 of file Simple.cpp.

198  {
199  PetscBool flg = PETSC_TRUE;
201  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "Simple interface options",
202  "none");
203  CHKERRG(ierr);
204  ierr = PetscOptionsString("-file_name", "file name", "", "mesh.h5m",
205  meshFileName, 255, &flg);
206  CHKERRG(ierr);
207  ierr = PetscOptionsEnd();
208  CHKERRG(ierr);
210 }
char meshFileName[255]
mesh file name
Definition: Simple.hpp:371
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ getOtherFiniteElements()

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

Get the Other Finite Elements.

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

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

Definition at line 340 of file Simple.hpp.

340 { return otherFEs; }
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:369

◆ getProblemName() [1/2]

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

Get the Problem Name.

Returns
const std::string
Examples
lesson8_contact.cpp, and reaction_diffusion_equation.cpp.

Definition at line 302 of file Simple.hpp.

302 { return nameOfProblem; }
std::string nameOfProblem
problem name
Definition: Simple.hpp:364

◆ getProblemName() [2/2]

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

Get the Problem Name.

Returns
std::string&

Definition at line 330 of file Simple.hpp.

330 { return nameOfProblem; }
std::string nameOfProblem
problem name
Definition: Simple.hpp:364

◆ getSkeletonFEName() [1/2]

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

Get the Skeleton FE Name.

Returns
const std::string
Examples
continuity_check_on_skeleton_with_simple_2d.cpp.

Definition at line 295 of file Simple.hpp.

295 { return skeletonFE; }
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367

◆ getSkeletonFEName() [2/2]

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

Get the Skeleton FE Name.

Returns
std::string&

Definition at line 323 of file Simple.hpp.

323 { return skeletonFE; }
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367

◆ loadFile() [1/2]

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

Load mesh file.

Parameters
optionsfile load options
mesh_file_namefile name if not set default or set by command line is used.
Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion_equation.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 212 of file Simple.cpp.

213  {
214  Interface &m_field = cOre;
216  PetscLogEventBegin(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
217 
218  if (!mesh_file_name.empty())
219  strcpy(meshFileName, mesh_file_name.c_str());
220 
221  // This is a case of distributed mesh and algebra. In that case each processor
222  // keep only part of the problem.
223  CHKERR m_field.get_moab().load_file(meshFileName, 0, options.c_str());
224  CHKERR m_field.rebuild_database();
225  // determine problem dimension
226  if (dIm == -1) {
227  int nb_ents_3d;
228  CHKERR m_field.get_moab().get_number_entities_by_dimension(
229  meshSet, 3, nb_ents_3d, true);
230  if (nb_ents_3d > 0) {
231  dIm = 3;
232  } else {
233  int nb_ents_2d;
234  CHKERR m_field.get_moab().get_number_entities_by_dimension(
235  meshSet, 2, nb_ents_2d, true);
236  if (nb_ents_2d > 0) {
237  dIm = 2;
238  } else {
239  dIm = 1;
240  }
241  }
242  }
243  Range ents;
244  CHKERR m_field.get_moab().get_entities_by_dimension(meshSet, dIm, ents, true);
245  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(ents, bitLevel,
246  false);
247  ParallelComm *pcomm =
248  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
249  if (pcomm == NULL)
250  pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
251  PetscLogEventEnd(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
253 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:353
char meshFileName[255]
mesh file name
Definition: Simple.hpp:371
int dIm
dimension of problem
Definition: Simple.hpp:372
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
char mesh_file_name[255]
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:347
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:291
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:345

◆ loadFile() [2/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 255 of file Simple.cpp.

255  {
257  Interface &m_field = cOre;
258  if (m_field.get_comm_size() == 1)
260  else
261  CHKERR loadFile("PARALLEL=READ_PART;"
262  "PARALLEL_RESOLVE_SHARED_ENTS;"
263  "PARTITION=PARALLEL_PARTITION;",
264  meshFileName);
266 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
char meshFileName[255]
mesh file name
Definition: Simple.hpp:371
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:212
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ query_interface()

MoFEMErrorCode MoFEM::Simple::query_interface ( const MOFEMuuid uuid,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 22 of file Simple.cpp.

23  {
25  *iface = NULL;
26  if (uuid == IDD_MOFEMSimple) {
27  *iface = const_cast<Simple *>(this);
29  }
30  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
32 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
static const MOFEMuuid IDD_MOFEMSimple
Definition: Simple.hpp:29

◆ removeBoundaryField()

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

Remove field form boundary.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 357 of file Simple.cpp.

357  {
358  Interface &m_field = cOre;
360 
361  auto remove_field_from_list = [&](auto &vec) {
362  auto it = std::find(vec.begin(), vec.end(), name);
363  if (it != vec.end())
364  vec.erase(it);
365  };
366 
367  remove_field_from_list(boundaryFields);
368 
370 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:356
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ removeDomainField()

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

Remove field form domain.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 341 of file Simple.cpp.

341  {
342  Interface &m_field = cOre;
344 
345  auto remove_field_from_list = [&](auto &vec) {
346  auto it = std::find(vec.begin(), vec.end(), name);
347  if (it != vec.end())
348  vec.erase(it);
349  };
350 
351  remove_field_from_list(noFieldFields);
352  remove_field_from_list(domainFields);
353 
355 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:359
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:355
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ removeSkeletonField()

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

Remove field form skeleton.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 372 of file Simple.cpp.

372  {
373  Interface &m_field = cOre;
375 
376  auto remove_field_from_list = [&](auto &vec) {
377  auto it = std::find(vec.begin(), vec.end(), name);
378  if (it != vec.end())
379  vec.erase(it);
380  };
381 
382  remove_field_from_list(skeletonFields);
383 
385 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ reSetUp()

MoFEMErrorCode MoFEM::Simple::reSetUp ( )

Rebuild internal MoFEM data structures.

Call this function after you add field or remove it.

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

Definition at line 684 of file Simple.cpp.

684  {
685  Interface &m_field = cOre;
687 
689  if (!skeletonFields.empty())
693 
694  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
695  CHKERR m_field.build_adjacencies(bitLevel);
696  // Set problem by the DOFs on the fields rather that by adding DOFs on the
697  // elements
698  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
700  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
701  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
702 
704 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:632
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:387
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:375
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1119
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:490
#define CHKERR
Inline error check.
Definition: definitions.h:602
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:350
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:345
MoFEMErrorCode setSkeletonAdjacency()
Definition: Simple.cpp:35

◆ setFieldOrder()

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

Set field order.

Parameters
std::field_namefield name
orderorder
rangeof entities to which order is set (If null it sat to all entities)
Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, mesh_smoothing.cpp, reaction_diffusion_equation.cpp, scalar_check_approximation_2d.cpp, simple_elasticity.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 482 of file Simple.cpp.

483  {
485  fieldsOrder.insert(
486  {field_name, {order, ents == NULL ? Range() : Range(*ents)}});
488 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:362
constexpr int order

◆ setSkeletonAdjacency() [1/5]

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

Definition at line 41 of file Simple.cpp.

41  {
42  Interface &m_field = cOre;
44 
45  auto defaultSkeletonEdge = [&](moab::Interface &moab, const Field &field,
46  const EntFiniteElement &fe,
47  Range &adjacency) -> MoFEMErrorCode {
49 
50  CHKERR DefaultElementAdjacency::defaultEdge(moab, field, fe, adjacency);
51 
52  if (std::find(domainFields.begin(), domainFields.end(),
53  field.getName()) != domainFields.end()) {
54 
55  const EntityHandle fe_ent = fe.getEnt();
56  Range bride_adjacency_edge, bride_adjacency;
57  CHKERR moab.get_adjacencies(&fe_ent, 1, 2, false, bride_adjacency_edge);
58 
59  switch (field.getSpace()) {
60  case H1:
61  CHKERR moab.get_connectivity(bride_adjacency_edge, bride_adjacency,
62  true);
63  case L2:
64  case HCURL:
65  case HDIV:
66  CHKERR moab.get_adjacencies(bride_adjacency_edge, 1, false,
67  bride_adjacency, moab::Interface::UNION);
68  bride_adjacency.merge(bride_adjacency_edge);
69  break;
70  case NOFIELD:
71  break;
72  default:
73  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
74  "this field is not implemented for TRI finite element");
75  }
76 
77  bride_adjacency = subtract(bride_adjacency, adjacency);
78 
79  for (auto e : bride_adjacency)
80  const_cast<SideNumber_multiIndex &>(fe.getSideNumberTable())
81  .insert(boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
82 
83  adjacency.merge(bride_adjacency);
84  }
85 
87  };
88 
89  CHKERR m_field.modify_finite_element_adjacency_table(skeletonFE, MBEDGE,
90  defaultSkeletonEdge);
91 
93 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
field with continuous normal traction
Definition: definitions.h:179
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:355
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
field with continuous tangents
Definition: definitions.h:178
#define CHKERR
Inline error check.
Definition: definitions.h:602
static MoFEMErrorCode defaultEdge(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
continuous field
Definition: definitions.h:177
field with C-1 continuity
Definition: definitions.h:180

◆ setSkeletonAdjacency() [2/5]

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

Definition at line 96 of file Simple.cpp.

96  {
97  Interface &m_field = cOre;
99 
100  auto defaultSkeletonEdge = [&](moab::Interface &moab, const Field &field,
101  const EntFiniteElement &fe,
102  Range &adjacency) -> MoFEMErrorCode {
104 
105  CHKERR DefaultElementAdjacency::defaultFace(moab, field, fe, adjacency);
106 
107  if (std::find(domainFields.begin(), domainFields.end(),
108  field.getName()) != domainFields.end()) {
109 
110  const EntityHandle fe_ent = fe.getEnt();
111  Range bride_adjacency_face, bride_adjacency;
112  CHKERR moab.get_adjacencies(&fe_ent, 1, 3, false, bride_adjacency_face);
113 
114  switch (field.getSpace()) {
115  case H1:
116  CHKERR moab.get_connectivity(bride_adjacency_face, bride_adjacency,
117  true);
118  case HCURL:
119  CHKERR moab.get_adjacencies(bride_adjacency_face, 1, false,
120  bride_adjacency, moab::Interface::UNION);
121  case L2:
122  case HDIV:
123  CHKERR moab.get_adjacencies(bride_adjacency_face, 2, false,
124  bride_adjacency, moab::Interface::UNION);
125  bride_adjacency.merge(bride_adjacency_face);
126  break;
127  case NOFIELD:
128  break;
129  default:
130  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
131  "this field is not implemented for TRI finite element");
132  }
133 
134  bride_adjacency = subtract(bride_adjacency, adjacency);
135 
136  for (auto e : bride_adjacency)
137  const_cast<SideNumber_multiIndex &>(fe.getSideNumberTable())
138  .insert(boost::shared_ptr<SideNumber>(new SideNumber(e, -1, 0, 0)));
139 
140  adjacency.merge(bride_adjacency);
141  }
142 
144  };
145 
146  CHKERR m_field.modify_finite_element_adjacency_table(skeletonFE, MBTRI,
147  defaultSkeletonEdge);
148  CHKERR m_field.modify_finite_element_adjacency_table(skeletonFE, MBQUAD,
149  defaultSkeletonEdge);
150 
152 }
MoFEM::Core & cOre
Definition: Simple.hpp:343
field with continuous normal traction
Definition: definitions.h:179
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
static MoFEMErrorCode defaultFace(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:355
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
field with continuous tangents
Definition: definitions.h:178
#define CHKERR
Inline error check.
Definition: definitions.h:602
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1879
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
continuous field
Definition: definitions.h:177
field with C-1 continuity
Definition: definitions.h:180

◆ setSkeletonAdjacency() [3/5]

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

Definition at line 155 of file Simple.cpp.

155  {
157  switch (getDim()) {
158  case 1:
159  THROW_MESSAGE("Not implemented");
160  case 2:
161  return setSkeletonAdjacency<2>();
162  case 3:
163  return setSkeletonAdjacency<3>();
164  default:
165  THROW_MESSAGE("Not implemented");
166  }
168 }
int getDim() const
Get the problem dimensio.
Definition: Simple.hpp:274
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:626
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ setSkeletonAdjacency() [4/5]

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

Set the skeleton adjacency object.

Parameters
dim
Returns
MoFEMErrorCode

Definition at line 170 of file Simple.cpp.

170  {
172  if(dim == -1)
173  dim = getDim();
174  switch(dim) {
175  case 2:
176  return setSkeletonAdjacency<2>();
177  case 3:
178  return setSkeletonAdjacency<3>();
179  default:
180  SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_IMPLEMENTED, "Not implemented");
181  }
183 }
int getDim() const
Get the problem dimensio.
Definition: Simple.hpp:274
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
const int dim
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ setSkeletonAdjacency() [5/5]

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

Definition at line 35 of file Simple.cpp.

35  {
36  static_assert(DIM == 2 || DIM == 3, "not implemented");
37  return MOFEM_NOT_IMPLEMENTED;
38 }

◆ setUp()

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

Setup problem.

Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, loop_entities.cpp, reaction_diffusion_equation.cpp, scalar_check_approximation_2d.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 672 of file Simple.cpp.

672  {
675  if (!skeletonFields.empty())
677  CHKERR defineProblem(is_partitioned);
682 }
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:632
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:387
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:658
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:490
#define CHKERR
Inline error check.
Definition: definitions.h:602
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:459
MoFEMErrorCode setSkeletonAdjacency()
Definition: Simple.cpp:35

Member Data Documentation

◆ bitLevel

const BitRefLevel MoFEM::Simple::bitLevel
private

BitRefLevel of the probelm.

Definition at line 345 of file Simple.hpp.

◆ boundaryFE

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

boundary finite element

Definition at line 366 of file Simple.hpp.

◆ boundaryFields

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

boundary fields

Definition at line 356 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 354 of file Simple.hpp.

◆ cOre

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

Definition at line 343 of file Simple.hpp.

◆ dataFields

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

Data fields.

Definition at line 358 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 372 of file Simple.hpp.

◆ dM

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

Discrete manager (interface to PETSc/MoFEM functions)

Definition at line 375 of file Simple.hpp.

◆ domainFE

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

domain finite element

Definition at line 365 of file Simple.hpp.

◆ domainFields

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

domain fields

Definition at line 355 of file Simple.hpp.

◆ fieldsOrder

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

fields order

Definition at line 362 of file Simple.hpp.

◆ meshFileName

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

mesh file name

Definition at line 371 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 353 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 348 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 349 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 350 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 351 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 347 of file Simple.hpp.

◆ nameOfProblem

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

problem name

Definition at line 364 of file Simple.hpp.

◆ noFieldDataFields

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

NOFIELD field name.

Definition at line 360 of file Simple.hpp.

◆ noFieldFields

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

NOFIELD field name.

Definition at line 359 of file Simple.hpp.

◆ otherFEs

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

Other finite elements.

Definition at line 369 of file Simple.hpp.

◆ skeletonFE

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

skeleton finite element

Definition at line 367 of file Simple.hpp.

◆ skeletonFields

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

fields on the skeleton

Definition at line 357 of file Simple.hpp.


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