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="PARALLEL=READ_PART;" "PARALLEL_RESOLVE_SHARED_ENTS;" "PARTITION=PARALLEL_PARTITION;", const std::string mesh_file_name="")
 Load mesh file. 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
basic_approx.cpp, basic_contact.cpp, basic_elastic.cpp, basic_helmholtz.cpp, basic_moment_of_inertia.cpp, basic_plastic.cpp, basic_poisson.cpp, basic_radiation.cpp, 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, 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:363
std::string domainFE
domain finite element
Definition: Simple.hpp:362
MoFEM::Core & cOre
Definition: Simple.hpp:340
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:350
char meshFileName[255]
mesh file name
Definition: Simple.hpp:368
int dIm
dimension of problem
Definition: Simple.hpp:369
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:344
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:345
std::string nameOfProblem
problem name
Definition: Simple.hpp:361
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:351
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:346
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:347
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:342
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:348

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

278  {
279  Interface &m_field = cOre;
281  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
282  verb);
283  boundaryFields.push_back(name);
284  if (space == NOFIELD)
285  SETERRQ(
286  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
287  "NOFIELD space for boundary filed not implemented in Simple interface");
289 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:353
#define CHKERR
Inline error check.
Definition: definitions.h:601
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

315  {
316 
317  Interface &m_field = cOre;
319  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
320  verb);
321  if (space == NOFIELD)
322  noFieldDataFields.push_back(name);
323  else
324  dataFields.push_back(name);
326 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:357
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:355
#define CHKERR
Inline error check.
Definition: definitions.h:601
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

260  {
261 
262  Interface &m_field = cOre;
264  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
265  verb);
266  if (space == NOFIELD)
267  noFieldFields.push_back(name);
268  else
269  domainFields.push_back(name);
271 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:356
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:352
#define CHKERR
Inline error check.
Definition: definitions.h:601
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

296  {
297 
298  Interface &m_field = cOre;
300  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
301  verb);
302  skeletonFields.push_back(name);
303  if (space == NOFIELD)
304  SETERRQ(
305  PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
306  "NOFIELD space for boundary filed not implemented in Simple interface");
308 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define CHKERR
Inline error check.
Definition: definitions.h:601
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

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

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

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

Definition at line 620 of file Simple.cpp.

620  {
621  Interface &m_field = cOre;
623  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
624  // Add finite elements
625  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
626  true);
627  CHKERR m_field.build_finite_elements(domainFE);
628  if (!boundaryFields.empty()) {
629  CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
630  boundaryFE, true);
631  CHKERR m_field.build_finite_elements(boundaryFE);
632  }
633  if (!skeletonFields.empty()) {
634  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm - 1,
635  skeletonFE, true);
636  CHKERR m_field.build_finite_elements(skeletonFE);
637  }
638  for (std::vector<std::string>::iterator fit = otherFEs.begin();
639  fit != otherFEs.end(); ++fit) {
640  CHKERR m_field.build_finite_elements(*fit);
641  }
642  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
644 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:363
std::string domainFE
domain finite element
Definition: Simple.hpp:362
MoFEM::Core & cOre
Definition: Simple.hpp:340
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:350
int dIm
dimension of problem
Definition: Simple.hpp:369
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:353
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:366
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364
#define CHKERR
Inline error check.
Definition: definitions.h:601
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:351
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:346
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

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

Definition at line 646 of file Simple.cpp.

646  {
647  Interface &m_field = cOre;
649  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
650  CHKERR m_field.build_adjacencies(bitLevel);
651  // Set problem by the DOFs on the fields rather that by adding DOFs on the
652  // elements
653  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
654  CHKERR DMSetUp(dM);
655  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
656  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
658 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:372
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define CHKERR
Inline error check.
Definition: definitions.h:601
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:347
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:342

◆ defineFiniteElements()

MoFEMErrorCode MoFEM::Simple::defineFiniteElements ( )

Define finite elements.

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

Definition at line 374 of file Simple.cpp.

374  {
375  Interface &m_field = cOre;
377 
378  auto clear_rows_and_cols = [&](auto &fe_name) {
380  const FiniteElement_multiIndex *fe_ptr;
381  CHKERR m_field.get_finite_elements(&fe_ptr);
382  auto &fe_by_name = const_cast<FiniteElement_multiIndex *>(fe_ptr)
383  ->get<FiniteElement_name_mi_tag>();
384  auto it_fe = fe_by_name.find(fe_name);
385  if (it_fe != fe_by_name.end()) {
386 
387  if(!fe_by_name.modify(it_fe, FiniteElement_row_change_bit_reset()))
388  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
389  "modification unsuccessful");
390 
391  if(!fe_by_name.modify(it_fe, FiniteElement_col_change_bit_reset()))
392  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
393  "modification unsuccessful");
394  }
396  };
397  CHKERR clear_rows_and_cols(domainFE);
398  CHKERR clear_rows_and_cols(boundaryFE);
399  CHKERR clear_rows_and_cols(skeletonFE);
400 
401  // Define finite elements
402  CHKERR m_field.add_finite_element(domainFE, MF_ZERO);
403 
404  auto add_fields = [&](auto &fe_name, auto &fields) {
406  for (auto &field : fields) {
407  CHKERR m_field.modify_finite_element_add_field_row(fe_name, field);
408  CHKERR m_field.modify_finite_element_add_field_col(fe_name, field);
409  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
410  }
412  };
413 
414  auto add_data_fields = [&](auto &fe_name, auto &fields) {
416  for (auto &field : fields)
417  CHKERR m_field.modify_finite_element_add_field_data(fe_name, field);
419  };
420 
421  CHKERR add_fields(domainFE, domainFields);
422  CHKERR add_data_fields(domainFE, dataFields);
423  CHKERR add_fields(domainFE, noFieldFields);
424  CHKERR add_data_fields(domainFE, noFieldDataFields);
425 
426  if (!boundaryFields.empty()) {
427  CHKERR m_field.add_finite_element(boundaryFE, MF_ZERO);
428  CHKERR add_fields(boundaryFE, domainFields);
429  CHKERR add_fields(boundaryFE, boundaryFields);
430  CHKERR add_fields(boundaryFE, skeletonFields);
431  CHKERR add_data_fields(boundaryFE, dataFields);
432  CHKERR add_data_fields(boundaryFE, noFieldDataFields);
433  CHKERR add_fields(boundaryFE, noFieldFields);
434  }
435  if (!skeletonFields.empty()) {
436  CHKERR m_field.add_finite_element(skeletonFE, MF_ZERO);
437  CHKERR add_fields(skeletonFE, domainFields);
438  CHKERR add_fields(skeletonFE, boundaryFields);
439  CHKERR add_fields(skeletonFE, skeletonFields);
440  CHKERR add_data_fields(skeletonFE, dataFields);
441  CHKERR add_data_fields(skeletonFE, noFieldDataFields);
442  CHKERR add_fields(skeletonFE, noFieldFields);
443  }
445 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:363
std::string domainFE
domain finite element
Definition: Simple.hpp:362
MoFEM::Core & cOre
Definition: Simple.hpp:340
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.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:353
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:356
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:352
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:357
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:355
#define CHKERR
Inline error check.
Definition: definitions.h:601
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791

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

447  {
448  Interface &m_field = cOre;
450  // Create dm instance
451  dM = createSmartDM(m_field.get_comm(), "DMMOFEM");
452  // set dm data structure which created mofem data structures
453  CHKERR DMMoFEMCreateMoFEM(dM, &m_field, nameOfProblem.c_str(), bitLevel);
454  CHKERR DMSetFromOptions(dM);
456  if (!boundaryFields.empty()) {
458  }
459  if (!skeletonFields.empty()) {
461  }
462  for (std::vector<std::string>::iterator fit = otherFEs.begin();
463  fit != otherFEs.end(); ++fit) {
464  CHKERR DMMoFEMAddElement(dM, fit->c_str());
465  }
466  CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
468 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:363
std::string domainFE
domain finite element
Definition: Simple.hpp:362
MoFEM::Core & cOre
Definition: Simple.hpp:340
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:372
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:353
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:366
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364
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:361
#define CHKERR
Inline error check.
Definition: definitions.h:601
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
auto createSmartDM
Creates smart DM object.
Definition: AuxPETSc.hpp:174
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:342
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 285 of file Simple.hpp.

285 { return boundaryFE; }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:363

◆ getBoundaryFEName() [2/2]

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

Get the Boundary FE Name.

Returns
std::string&

Definition at line 313 of file Simple.hpp.

313 { return boundaryFE; }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:363

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

271 { return dIm; }
int dIm
dimension of problem
Definition: Simple.hpp:369

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

694  {
696  CHKERR PetscObjectReference(getPetscObject(dM.get()));
697  *dm = dM.get();
699 }
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:372
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
PetscObject getPetscObject(T obj)
Definition: AuxPETSc.hpp:23
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

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

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

278 { return domainFE; }
std::string domainFE
domain finite element
Definition: Simple.hpp:362

◆ getDomainFEName() [2/2]

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

Get the Domain FE Name.

Returns
std::string&

Definition at line 306 of file Simple.hpp.

306 { return domainFE; }
std::string domainFE
domain finite element
Definition: Simple.hpp:362

◆ 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:368
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:549
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513

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

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

◆ getProblemName() [1/2]

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

Get the Problem Name.

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

Definition at line 299 of file Simple.hpp.

299 { return nameOfProblem; }
std::string nameOfProblem
problem name
Definition: Simple.hpp:361

◆ getProblemName() [2/2]

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

Get the Problem Name.

Returns
std::string&

Definition at line 327 of file Simple.hpp.

327 { return nameOfProblem; }
std::string nameOfProblem
problem name
Definition: Simple.hpp:361

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

292 { return skeletonFE; }
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364

◆ getSkeletonFEName() [2/2]

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

Get the Skeleton FE Name.

Returns
std::string&

Definition at line 320 of file Simple.hpp.

320 { return skeletonFE; }
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364

◆ loadFile()

MoFEMErrorCode MoFEM::Simple::loadFile ( const std::string  options = "PARALLEL=READ_PART;" "PARALLEL_RESOLVE_SHARED_ENTS;" "PARTITION=PARALLEL_PARTITION;",
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:340
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:350
char meshFileName[255]
mesh file name
Definition: Simple.hpp:368
int dIm
dimension of problem
Definition: Simple.hpp:369
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
char mesh_file_name[255]
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:344
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:290
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:342

◆ 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:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
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 344 of file Simple.cpp.

344  {
345  Interface &m_field = cOre;
347 
348  auto remove_field_from_list = [&](auto &vec) {
349  auto it = std::find(vec.begin(), vec.end(), name);
350  if (it != vec.end())
351  vec.erase(it);
352  };
353 
354  remove_field_from_list(boundaryFields);
355 
357 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:353
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ removeDomainField()

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

Remove field form domain.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 328 of file Simple.cpp.

328  {
329  Interface &m_field = cOre;
331 
332  auto remove_field_from_list = [&](auto &vec) {
333  auto it = std::find(vec.begin(), vec.end(), name);
334  if (it != vec.end())
335  vec.erase(it);
336  };
337 
338  remove_field_from_list(noFieldFields);
339  remove_field_from_list(domainFields);
340 
342 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:356
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:352
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ removeSkeletonField()

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

Remove field form skeleton.

Parameters
name
Returns
MoFEMErrorCode

Definition at line 359 of file Simple.cpp.

359  {
360  Interface &m_field = cOre;
362 
363  auto remove_field_from_list = [&](auto &vec) {
364  auto it = std::find(vec.begin(), vec.end(), name);
365  if (it != vec.end())
366  vec.erase(it);
367  };
368 
369  remove_field_from_list(skeletonFields);
370 
372 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

672  {
673  Interface &m_field = cOre;
675 
677  if (!skeletonFields.empty())
681 
682  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
683  CHKERR m_field.build_adjacencies(bitLevel);
684  // Set problem by the DOFs on the fields rather that by adding DOFs on the
685  // elements
686  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
688  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_FALSE;
689  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
690 
692 }
MoFEM::Core & cOre
Definition: Simple.hpp:340
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:620
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:374
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:372
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1119
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:478
#define CHKERR
Inline error check.
Definition: definitions.h:601
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:347
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:342
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 470 of file Simple.cpp.

471  {
473  fieldsOrder.insert(
474  {field_name, {order, ents == NULL ? Range() : Range(*ents)}});
476 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:359
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:340
field with continuous normal traction
Definition: definitions.h:178
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:352
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
field with continuous tangents
Definition: definitions.h:177
#define CHKERR
Inline error check.
Definition: definitions.h:601
static MoFEMErrorCode defaultEdge(Interface &moab, const Field &field, const EntFiniteElement &fe, Range &adjacency)
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
continuous field
Definition: definitions.h:176
field with C-1 continuity
Definition: definitions.h:179

◆ 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:340
field with continuous normal traction
Definition: definitions.h:178
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
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:482
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:352
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
field with continuous tangents
Definition: definitions.h:177
#define CHKERR
Inline error check.
Definition: definitions.h:601
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1791
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
continuous field
Definition: definitions.h:176
field with C-1 continuity
Definition: definitions.h:179

◆ 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:271
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:625
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ 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:271
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
const int dim
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

660  {
663  if (!skeletonFields.empty())
665  CHKERR defineProblem(is_partitioned);
670 }
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:620
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:374
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:646
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:478
#define CHKERR
Inline error check.
Definition: definitions.h:601
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:447
MoFEMErrorCode setSkeletonAdjacency()
Definition: Simple.cpp:35

Member Data Documentation

◆ bitLevel

const BitRefLevel MoFEM::Simple::bitLevel
private

BitRefLevel of the probelm.

Definition at line 342 of file Simple.hpp.

◆ boundaryFE

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

boundary finite element

Definition at line 363 of file Simple.hpp.

◆ boundaryFields

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

boundary fields

Definition at line 353 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 351 of file Simple.hpp.

◆ cOre

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

Definition at line 340 of file Simple.hpp.

◆ dataFields

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

Data fields.

Definition at line 355 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 369 of file Simple.hpp.

◆ dM

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

Discrete manager (interface to PETSc/MoFEM functions)

Definition at line 372 of file Simple.hpp.

◆ domainFE

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

domain finite element

Definition at line 362 of file Simple.hpp.

◆ domainFields

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

domain fields

Definition at line 352 of file Simple.hpp.

◆ fieldsOrder

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

fields order

Definition at line 359 of file Simple.hpp.

◆ meshFileName

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

mesh file name

Definition at line 368 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 350 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 345 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 346 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 347 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 348 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 344 of file Simple.hpp.

◆ nameOfProblem

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

problem name

Definition at line 361 of file Simple.hpp.

◆ noFieldDataFields

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

NOFIELD field name.

Definition at line 357 of file Simple.hpp.

◆ noFieldFields

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

NOFIELD field name.

Definition at line 356 of file Simple.hpp.

◆ otherFEs

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

Other finite elements.

Definition at line 366 of file Simple.hpp.

◆ skeletonFE

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

skeleton finite element

Definition at line 364 of file Simple.hpp.

◆ skeletonFields

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

fields on the skeleton

Definition at line 354 of file Simple.hpp.


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