v0.9.0
Public Member Functions | Public Attributes | 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)
 
 ~Simple ()
 Destructor. More...
 
MoFEMErrorCode getOptions ()
 get options More...
 
MoFEMErrorCode loadFile (const std::string options="PARALLEL=READ_PART;" "PARALLEL_RESOLVE_SHARED_ENTS;" "PARTITION=PARALLEL_PARTITION;")
 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 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 buildProblem ()
 Build problem. More...
 
MoFEMErrorCode setUp (const PetscBool is_partitioned=PETSC_TRUE)
 Setup problem. More...
 
MoFEMErrorCode getDM (DM *dm)
 Get DM. More...
 
SmartPetscObj< DM > getDM ()
 Return smart DM object. More...
 
int getDim () const
 
const std::string getDomainFEName () const
 
const std::string getBoundaryFEName () const
 
const std::string getSkeletonFEName () const
 
const std::string getProblemName () const
 
std::string & getDomainFEName ()
 
std::string & getBoundaryFEName ()
 
std::string & getSkeletonFEName ()
 
std::string & getProblemName ()
 
std::vector< std::string > & getOtherFiniteElements ()
 
- 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
 

Public Attributes

MoFEM::CorecOre
 

Private Attributes

const BitRefLevel bitLevel
 
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::map< 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
 

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, hello_world.cpp, loop_entities.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 36 of file Simple.hpp.

Constructor & Destructor Documentation

◆ Simple()

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

Definition at line 34 of file Simple.cpp.

35  : cOre(const_cast<Core &>(core)), bitLevel(BitRefLevel().set(0)),
36  meshSet(0), boundaryMeshset(0), nameOfProblem("SimpleProblem"),
37  domainFE("dFE"), boundaryFE("bFE"), skeletonFE("sFE"), dIm(-1) {
38  PetscLogEventRegister("LoadMesh", 0, &MOFEM_EVENT_SimpleLoadMesh);
39  PetscLogEventRegister("buildFields", 0, &MOFEM_EVENT_SimpleBuildFields);
40  PetscLogEventRegister("buildFiniteElements", 0,
42  PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleBuildProblem);
43  PetscLogEventRegister("SimpleKSPSolve", 0, &MOFEM_EVENT_SimpleKSPSolve);
44  strcpy(meshFileName, "mesh.h5m");
45 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:253
std::string domainFE
domain finite element
Definition: Simple.hpp:252
MoFEM::Core & cOre
Definition: Simple.hpp:41
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:242
char meshFileName[255]
mesh file name
Definition: Simple.hpp:258
int dIm
dimension of problem
Definition: Simple.hpp:259
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:254
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:236
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:237
std::string nameOfProblem
problem name
Definition: Simple.hpp:251
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:243
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:238
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:239
const BitRefLevel bitLevel
Definition: Simple.hpp:234
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:240

◆ ~Simple()

MoFEM::Simple::~Simple ( )

Destructor.

Definition at line 46 of file Simple.cpp.

46 {}

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
hello_world.cpp, mesh_smoothing.cpp, and simple_interface.cpp.

Definition at line 116 of file Simple.cpp.

120  {
121  Interface &m_field = cOre;
123  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
124  verb);
125  boundaryFields.push_back(name);
127 }
MoFEM::Core & cOre
Definition: Simple.hpp:41
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:245
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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

Definition at line 145 of file Simple.cpp.

149  {
150 
151  Interface &m_field = cOre;
153  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
154  verb);
155  dataFields.push_back(name);
157 }
MoFEM::Core & cOre
Definition: Simple.hpp:41
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:247
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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, hello_world.cpp, loop_entities.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 101 of file Simple.cpp.

105  {
106 
107  Interface &m_field = cOre;
109  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
110  verb);
111  domainFields.push_back(name);
113 }
MoFEM::Core & cOre
Definition: Simple.hpp:41
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:244
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

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

Definition at line 130 of file Simple.cpp.

134  {
135 
136  Interface &m_field = cOre;
138  CHKERR m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
139  verb);
140  skeletonFields.push_back(name);
142 }
MoFEM::Core & cOre
Definition: Simple.hpp:41
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:246
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

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

264  {
265  Interface &m_field = cOre;
267  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
268  // take skin
269  {
270  Range domain_ents;
271  CHKERR m_field.get_moab().get_entities_by_dimension(meshSet, dIm,
272  domain_ents, true);
273  Skinner skin(&m_field.get_moab());
274  Range domain_skin;
275  CHKERR skin.find_skin(0, domain_ents, false, domain_skin);
276  // filter not owned entities, those are not on boundary
277  ParallelComm *pcomm =
278  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
279  Range proc_domain_skin;
280  CHKERR pcomm->filter_pstatus(domain_skin,
281  PSTATUS_SHARED | PSTATUS_MULTISHARED,
282  PSTATUS_NOT, -1, &proc_domain_skin);
283  // create boundary meshset
284  if (boundaryMeshset != 0) {
285  CHKERR m_field.get_moab().delete_entities(&boundaryMeshset, 1);
286  }
287  CHKERR m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
288  CHKERR m_field.get_moab().add_entities(boundaryMeshset, proc_domain_skin);
289  for (int dd = 0; dd != dIm - 1; dd++) {
290  Range adj;
291  CHKERR m_field.get_moab().get_adjacencies(proc_domain_skin, dd, false,
292  adj, moab::Interface::UNION);
293  CHKERR m_field.get_moab().add_entities(boundaryMeshset, adj);
294  }
295  }
296  auto comm_interface_ptr = m_field.getInterface<CommInterface>();
297  // Add entities to the fields
298  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
299  CHKERR m_field.add_ents_to_field_by_dim(meshSet, dIm, domainFields[ff]);
300  CHKERR comm_interface_ptr->synchroniseFieldEntities(domainFields[ff], 0);
301  }
302  for (unsigned int ff = 0; ff != dataFields.size(); ff++) {
303  CHKERR m_field.add_ents_to_field_by_dim(meshSet, dIm, dataFields[ff]);
304  CHKERR comm_interface_ptr->synchroniseFieldEntities(dataFields[ff], 0);
305  }
306  for (unsigned int ff = 0; ff != boundaryFields.size(); ff++) {
307  CHKERR m_field.add_ents_to_field_by_dim(boundaryMeshset, dIm - 1,
308  boundaryFields[ff]);
309  CHKERR comm_interface_ptr->synchroniseFieldEntities(boundaryFields[ff], 0);
310  }
311  for (unsigned int ff = 0; ff != skeletonFields.size(); ff++) {
312  CHKERR m_field.add_ents_to_field_by_dim(meshSet, dIm - 1,
313  skeletonFields[ff]);
314  CHKERR comm_interface_ptr->synchroniseFieldEntities(skeletonFields[ff], 0);
315  }
316  // Set order
317  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
318  if (fieldsOrder.find(domainFields[ff]) == fieldsOrder.end()) {
319  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
320  "Order for field not set %s", domainFields[ff].c_str());
321  }
322  int dds = 0;
323  const Field *field = m_field.get_field_structure(domainFields[ff]);
324  switch (field->getSpace()) {
325  case L2:
326  dds = dIm;
327  break;
328  case HDIV:
329  dds = 2;
330  break;
331  case HCURL:
332  dds = 1;
333  break;
334  case H1:
335  dds = 1;
336  break;
337  default:
338  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
339  "Huston we have a problem");
340  }
341  if (field->getSpace() == H1) {
342  CHKERR m_field.set_field_order(meshSet, MBVERTEX, domainFields[ff], 1);
343  }
344  for (int dd = dds; dd <= dIm; dd++) {
345  Range ents;
346  CHKERR m_field.get_field_entities_by_dimension(domainFields[ff], dd,
347  ents);
348  if (!fieldsOrder.at(domainFields[ff]).second.empty()) {
349  ents = intersect(ents, fieldsOrder.at(domainFields[ff]).second);
350  }
351  CHKERR m_field.set_field_order(ents, domainFields[ff],
352  fieldsOrder.at(domainFields[ff]).first);
353  }
354  }
355  // Set order to data fields
356  for (unsigned int ff = 0; ff != dataFields.size(); ff++) {
357  if (fieldsOrder.find(dataFields[ff]) == fieldsOrder.end()) {
358  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
359  "Order for field not set %s", dataFields[ff].c_str());
360  }
361  int dds = 0;
362  const Field *field = m_field.get_field_structure(dataFields[ff]);
363  switch (field->getSpace()) {
364  case L2:
365  dds = dIm;
366  break;
367  case HDIV:
368  dds = 2;
369  break;
370  case HCURL:
371  dds = 1;
372  break;
373  case H1:
374  dds = 1;
375  break;
376  default:
377  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
378  "Huston we have a problem");
379  }
380  if (field->getSpace() == H1) {
381  CHKERR m_field.set_field_order(meshSet, MBVERTEX, dataFields[ff], 1);
382  }
383  for (int dd = dds; dd <= dIm; dd++) {
384  Range ents;
385  CHKERR m_field.get_field_entities_by_dimension(dataFields[ff], dd, ents);
386  CHKERR m_field.set_field_order(ents, dataFields[ff],
387  fieldsOrder.at(dataFields[ff]).first);
388  }
389  }
390  // Set order to boundary
391  for (unsigned int ff = 0; ff != boundaryFields.size(); ff++) {
392  if (fieldsOrder.find(boundaryFields[ff]) == fieldsOrder.end()) {
393  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
394  "Order for field not set %s", boundaryFields[ff].c_str());
395  }
396  int dds = 0;
397  const Field *field = m_field.get_field_structure(boundaryFields[ff]);
398  switch (field->getSpace()) {
399  case L2:
400  dds = dIm - 1;
401  break;
402  case HDIV:
403  dds = 2;
404  break;
405  case HCURL:
406  dds = 1;
407  break;
408  case H1:
409  dds = 1;
410  break;
411  default:
412  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
413  "Huston we have a problem");
414  }
415  if (field->getSpace() == H1) {
416  CHKERR m_field.set_field_order(meshSet, MBVERTEX, boundaryFields[ff], 1);
417  }
418  for (int dd = dds; dd <= dIm - 1; dd++) {
419  Range ents;
420  CHKERR m_field.get_field_entities_by_dimension(boundaryFields[ff], dd,
421  ents);
422  CHKERR m_field.set_field_order(ents, boundaryFields[ff],
423  fieldsOrder.at(boundaryFields[ff]).first);
424  }
425  }
426  // Set order to skeleton
427  for (unsigned int ff = 0; ff != skeletonFields.size(); ff++) {
428  if (fieldsOrder.find(skeletonFields[ff]) == fieldsOrder.end()) {
429  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
430  "Order for field not set %s", skeletonFields[ff].c_str());
431  }
432  int dds = 0;
433  const Field *field = m_field.get_field_structure(skeletonFields[ff]);
434  switch (field->getSpace()) {
435  case L2:
436  dds = dIm - 1;
437  break;
438  case HDIV:
439  dds = 2;
440  break;
441  case HCURL:
442  dds = 1;
443  break;
444  case H1:
445  dds = 1;
446  break;
447  default:
448  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
449  "Huston we have a problem");
450  }
451  if (field->getSpace() == H1) {
452  CHKERR m_field.set_field_order(meshSet, MBVERTEX, skeletonFields[ff], 1);
453  }
454  for (int dd = dds; dd <= dIm - 1; dd++) {
455  Range ents;
456  CHKERR m_field.get_field_entities_by_dimension(skeletonFields[ff], dd,
457  ents);
458  CHKERR m_field.set_field_order(ents, skeletonFields[ff],
459  fieldsOrder.at(skeletonFields[ff]).first);
460  }
461  }
462  // Build fields
463  CHKERR m_field.build_fields();
464  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
466 }
MoFEM::Core & cOre
Definition: Simple.hpp:41
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:242
field with continuous normal traction
Definition: definitions.h:173
int dIm
dimension of problem
Definition: Simple.hpp:259
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:245
std::map< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:249
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:244
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:247
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:237
field with continuous tangents
Definition: definitions.h:172
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:285
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:243
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:246
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
continuous field
Definition: definitions.h:171
field with C-1 continuity
Definition: definitions.h:174

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

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

Definition at line 468 of file Simple.cpp.

468  {
469  Interface &m_field = cOre;
471  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
472  // Add finite elements
473  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
474  true);
475  CHKERR m_field.build_finite_elements(domainFE);
476  if (!boundaryFields.empty()) {
477  CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
478  boundaryFE, true);
479  CHKERR m_field.build_finite_elements(boundaryFE);
480  }
481  if (!skeletonFields.empty()) {
482  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm - 1,
483  skeletonFE, true);
484  CHKERR m_field.build_finite_elements(skeletonFE);
485  }
486  for (std::vector<std::string>::iterator fit = otherFEs.begin();
487  fit != otherFEs.end(); ++fit) {
488  CHKERR m_field.build_finite_elements(*fit);
489  }
490  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
492 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:253
std::string domainFE
domain finite element
Definition: Simple.hpp:252
MoFEM::Core & cOre
Definition: Simple.hpp:41
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:242
int dIm
dimension of problem
Definition: Simple.hpp:259
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:245
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:256
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:254
#define CHKERR
Inline error check.
Definition: definitions.h:596
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:243
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:238
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:246
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

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

Definition at line 494 of file Simple.cpp.

494  {
495  Interface &m_field = cOre;
497  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
498  CHKERR m_field.build_adjacencies(bitLevel);
499  // Set problem by the DOFs on the fields rather that by adding DOFs on the
500  // elements
501  m_field.getInterface<ProblemsManager>()->buildProblemFromFields = PETSC_TRUE;
502  CHKERR DMSetUp(dM);
503  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
505 }
MoFEM::Core & cOre
Definition: Simple.hpp:41
SmartPetscObj< DM > dM
Definition: Simple.hpp:261
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:239
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
const BitRefLevel bitLevel
Definition: Simple.hpp:234

◆ defineFiniteElements()

MoFEMErrorCode MoFEM::Simple::defineFiniteElements ( )

Define finite elements.

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

Definition at line 159 of file Simple.cpp.

159  {
160  Interface &m_field = cOre;
162  // Define finite elements
163  CHKERR m_field.add_finite_element(domainFE);
164  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
165  CHKERR m_field.modify_finite_element_add_field_row(domainFE,
166  domainFields[ff]);
167  CHKERR m_field.modify_finite_element_add_field_col(domainFE,
168  domainFields[ff]);
169  CHKERR m_field.modify_finite_element_add_field_data(domainFE,
170  domainFields[ff]);
171  }
172  for (unsigned int ff = 0; ff != dataFields.size(); ff++) {
173  CHKERR m_field.modify_finite_element_add_field_data(domainFE,
174  dataFields[ff]);
175  }
176  if (!boundaryFields.empty()) {
177  CHKERR m_field.add_finite_element(boundaryFE);
178  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
179  CHKERR m_field.modify_finite_element_add_field_row(boundaryFE,
180  domainFields[ff]);
181  CHKERR m_field.modify_finite_element_add_field_col(boundaryFE,
182  domainFields[ff]);
183  CHKERR m_field.modify_finite_element_add_field_data(boundaryFE,
184  domainFields[ff]);
185  }
186  for (unsigned int ff = 0; ff != boundaryFields.size(); ff++) {
187  CHKERR m_field.modify_finite_element_add_field_row(boundaryFE,
188  boundaryFields[ff]);
189  CHKERR m_field.modify_finite_element_add_field_col(boundaryFE,
190  boundaryFields[ff]);
191  CHKERR m_field.modify_finite_element_add_field_data(boundaryFE,
192  boundaryFields[ff]);
193  }
194  for (unsigned int ff = 0; ff != skeletonFields.size(); ff++) {
195  CHKERR m_field.modify_finite_element_add_field_row(boundaryFE,
196  skeletonFields[ff]);
197  CHKERR m_field.modify_finite_element_add_field_col(boundaryFE,
198  skeletonFields[ff]);
199  CHKERR m_field.modify_finite_element_add_field_data(boundaryFE,
200  skeletonFields[ff]);
201  }
202  }
203  if (!skeletonFields.empty()) {
204  CHKERR m_field.add_finite_element(skeletonFE);
205  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
206  CHKERR m_field.modify_finite_element_add_field_row(skeletonFE,
207  domainFields[ff]);
208  CHKERR m_field.modify_finite_element_add_field_col(skeletonFE,
209  domainFields[ff]);
210  CHKERR m_field.modify_finite_element_add_field_data(skeletonFE,
211  domainFields[ff]);
212  }
213  for (unsigned int ff = 0; ff != boundaryFields.size(); ff++) {
214  CHKERR m_field.modify_finite_element_add_field_row(skeletonFE,
215  boundaryFields[ff]);
216  CHKERR m_field.modify_finite_element_add_field_col(skeletonFE,
217  boundaryFields[ff]);
218  CHKERR m_field.modify_finite_element_add_field_data(skeletonFE,
219  boundaryFields[ff]);
220  }
221  for (unsigned int ff = 0; ff != skeletonFields.size(); ff++) {
222  CHKERR m_field.modify_finite_element_add_field_row(skeletonFE,
223  skeletonFields[ff]);
224  CHKERR m_field.modify_finite_element_add_field_col(skeletonFE,
225  skeletonFields[ff]);
226  CHKERR m_field.modify_finite_element_add_field_data(skeletonFE,
227  skeletonFields[ff]);
228  }
229  }
231 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:253
std::string domainFE
domain finite element
Definition: Simple.hpp:252
MoFEM::Core & cOre
Definition: Simple.hpp:41
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:245
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:254
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:244
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:247
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:246

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

233  {
234  Interface &m_field = cOre;
236  // Create dm instance
237  dM = createSmartDM(m_field.get_comm(), "DMMOFEM");
238  // set dm data structure which created mofem data structures
239  CHKERR DMMoFEMCreateMoFEM(dM, &m_field, nameOfProblem.c_str(), bitLevel);
240  CHKERR DMSetFromOptions(dM);
242  if (!boundaryFields.empty()) {
244  }
245  if (!skeletonFields.empty()) {
247  }
248  for (std::vector<std::string>::iterator fit = otherFEs.begin();
249  fit != otherFEs.end(); ++fit) {
250  CHKERR DMMoFEMAddElement(dM, fit->c_str());
251  }
252  CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
254 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:253
std::string domainFE
domain finite element
Definition: Simple.hpp:252
MoFEM::Core & cOre
Definition: Simple.hpp:41
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:414
SmartPetscObj< DM > dM
Definition: Simple.hpp:261
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:245
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:256
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:254
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:109
std::string nameOfProblem
problem name
Definition: Simple.hpp:251
#define CHKERR
Inline error check.
Definition: definitions.h:596
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:246
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
auto createSmartDM
Creates smart DM object.
Definition: AuxPETSc.hpp:173
const BitRefLevel bitLevel
Definition: Simple.hpp:234
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:971

◆ getBoundaryFEName() [1/2]

const std::string MoFEM::Simple::getBoundaryFEName ( ) const
Examples
hello_world.cpp, mesh_smoothing.cpp, and simple_interface.cpp.

Definition at line 222 of file Simple.hpp.

222 { return boundaryFE; }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:253

◆ getBoundaryFEName() [2/2]

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

Definition at line 227 of file Simple.hpp.

227 { return boundaryFE; }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:253

◆ getDim()

int MoFEM::Simple::getDim ( ) const

Definition at line 220 of file Simple.hpp.

220 { return dIm; }
int dIm
dimension of problem
Definition: Simple.hpp:259

◆ 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, hello_world.cpp, loop_entities.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 517 of file Simple.cpp.

517  {
519  CHKERR PetscObjectReference(getPetscObject(dM.get()));
520  *dm = dM.get();
522 }
SmartPetscObj< DM > dM
Definition: Simple.hpp:261
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscObject getPetscObject(T obj)
Definition: AuxPETSc.hpp:23
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

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

218 { return dM; }
SmartPetscObj< DM > dM
Definition: Simple.hpp:261

◆ getDomainFEName() [1/2]

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

◆ getDomainFEName() [2/2]

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

Definition at line 226 of file Simple.hpp.

226 { return domainFE; }
std::string domainFE
domain finite element
Definition: Simple.hpp:252

◆ 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, hello_world.cpp, loop_entities.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 48 of file Simple.cpp.

48  {
49  PetscBool flg = PETSC_TRUE;
51  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "Simple interface options",
52  "none");
53  CHKERRG(ierr);
54  ierr = PetscOptionsString("-file_name", "file name", "", "mesh.h5m",
55  meshFileName, 255, &flg);
56  CHKERRG(ierr);
57  ierr = PetscOptionsEnd();
58  CHKERRG(ierr);
60 }
char meshFileName[255]
mesh file name
Definition: Simple.hpp:258
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ getOtherFiniteElements()

std::vector<std::string>& MoFEM::Simple::getOtherFiniteElements ( )
Examples
mesh_smoothing.cpp.

Definition at line 231 of file Simple.hpp.

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

◆ getProblemName() [1/2]

const std::string MoFEM::Simple::getProblemName ( ) const
Examples
reaction_diffusion_equation.cpp.

Definition at line 224 of file Simple.hpp.

224 { return nameOfProblem; }
std::string nameOfProblem
problem name
Definition: Simple.hpp:251

◆ getProblemName() [2/2]

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

Definition at line 229 of file Simple.hpp.

229 { return nameOfProblem; }
std::string nameOfProblem
problem name
Definition: Simple.hpp:251

◆ getSkeletonFEName() [1/2]

const std::string MoFEM::Simple::getSkeletonFEName ( ) const
Examples
continuity_check_on_skeleton_with_simple_2d.cpp, and hello_world.cpp.

Definition at line 223 of file Simple.hpp.

223 { return skeletonFE; }
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:254

◆ getSkeletonFEName() [2/2]

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

Definition at line 228 of file Simple.hpp.

228 { return skeletonFE; }
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:254

◆ loadFile()

MoFEMErrorCode MoFEM::Simple::loadFile ( const std::string  options = "PARALLEL=READ_PART;" "PARALLEL_RESOLVE_SHARED_ENTS;" "PARTITION=PARALLEL_PARTITION;")

Load mesh file.

Parameters
field_namefile name
Returns
error code
Examples
continuity_check_on_skeleton_with_simple_2d.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hello_world.cpp, loop_entities.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 62 of file Simple.cpp.

62  {
63  Interface &m_field = cOre;
65  PetscLogEventBegin(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
66  // This is a case of distributed mesh and algebra. In that case each processor
67  // keep only part of the problem.
68  CHKERR m_field.get_moab().load_file(meshFileName, 0, options.c_str());
69  CHKERR m_field.rebuild_database();
70  // determine problem dimension
71  if (dIm == -1) {
72  int nb_ents_3d;
73  CHKERR m_field.get_moab().get_number_entities_by_dimension(
74  meshSet, 3, nb_ents_3d, true);
75  if (nb_ents_3d > 0) {
76  dIm = 3;
77  } else {
78  int nb_ents_2d;
79  CHKERR m_field.get_moab().get_number_entities_by_dimension(
80  meshSet, 2, nb_ents_2d, true);
81  if (nb_ents_2d > 0) {
82  dIm = 2;
83  } else {
84  dIm = 1;
85  }
86  }
87  }
88  Range ents;
89  CHKERR m_field.get_moab().get_entities_by_dimension(meshSet, dIm, ents, true);
90  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevel(ents, bitLevel,
91  false);
92  ParallelComm *pcomm =
93  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
94  if (pcomm == NULL)
95  pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
96  PetscLogEventEnd(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
98 }
MoFEM::Core & cOre
Definition: Simple.hpp:41
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:242
char meshFileName[255]
mesh file name
Definition: Simple.hpp:258
int dIm
dimension of problem
Definition: Simple.hpp:259
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:236
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:285
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
const BitRefLevel bitLevel
Definition: Simple.hpp:234

◆ 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:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
static const MOFEMuuid IDD_MOFEMSimple
Definition: Simple.hpp:29

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

257  {
259  fieldsOrder[field_name] =
260  std::pair<int, Range>(order, ents == NULL ? Range() : Range(*ents));
262 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
std::map< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:249

◆ 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, hello_world.cpp, loop_entities.cpp, reaction_diffusion_equation.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 507 of file Simple.cpp.

507  {
510  CHKERR defineProblem(is_partitioned);
515 }
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:468
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:159
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:494
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:264
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:233

Member Data Documentation

◆ bitLevel

const BitRefLevel MoFEM::Simple::bitLevel
private

Definition at line 234 of file Simple.hpp.

◆ boundaryFE

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

boundary finite element

Definition at line 253 of file Simple.hpp.

◆ boundaryFields

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

boundary fields

Definition at line 245 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 243 of file Simple.hpp.

◆ cOre

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

Definition at line 41 of file Simple.hpp.

◆ dataFields

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

Data fields.

Definition at line 247 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 259 of file Simple.hpp.

◆ dM

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

Definition at line 261 of file Simple.hpp.

◆ domainFE

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

domain finite element

Definition at line 252 of file Simple.hpp.

◆ domainFields

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

domain fields

Definition at line 244 of file Simple.hpp.

◆ fieldsOrder

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

fields order

Definition at line 249 of file Simple.hpp.

◆ meshFileName

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

mesh file name

Definition at line 258 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 242 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 237 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 238 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 239 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 240 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 236 of file Simple.hpp.

◆ nameOfProblem

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

problem name

Definition at line 251 of file Simple.hpp.

◆ otherFEs

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

Other finite elements.

Definition at line 256 of file Simple.hpp.

◆ skeletonFE

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

skeleton finite element

Definition at line 254 of file Simple.hpp.

◆ skeletonFields

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

fields on the skeleton

Definition at line 246 of file Simple.hpp.


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