v0.8.13
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 ()
 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 ()
 Setup problem. More...
 
MoFEMErrorCode getDM (DM *dm)
 Get DM. More...
 
int getDim () const
 
const std::string & getDomainFEName () const
 
const std::string & getBoundaryFEName () const
 
const std::string & getSkeletonFEName () const
 
std::string & getDomainFEName ()
 
std::string & getBoundaryFEName ()
 
std::string & getSkeletonFEName ()
 
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 ()
 
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::string nameOfProblem
 problem name 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 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...
 
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:
hello_world.cpp, mesh_smoothing.cpp, simple_elasticity.cpp, and simple_interface.cpp.

Definition at line 33 of file Simple.hpp.

Constructor & Destructor Documentation

◆ Simple()

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

Definition at line 40 of file Simple.cpp.

41  : cOre(const_cast<Core &>(core)), bitLevel(BitRefLevel().set(0)),
42  meshSet(0), boundaryMeshset(0), domainFE("dFE"), boundaryFE("bFE"),
43  skeletonFE("sFE"), dIm(-1), dM(PETSC_NULL) {
44  PetscLogEventRegister("LoadMesh", 0, &MOFEM_EVENT_SimpleLoadMesh);
45  PetscLogEventRegister("buildFields", 0, &MOFEM_EVENT_SimpleBuildFields);
46  PetscLogEventRegister("buildFiniteElements", 0,
48  PetscLogEventRegister("SimpleSetUp", 0, &MOFEM_EVENT_SimpleBuildProblem);
49  PetscLogEventRegister("SimpleKSPSolve", 0, &MOFEM_EVENT_SimpleKSPSolve);
50  strcpy(meshFileName, "mesh.h5m");
51 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:230
std::string domainFE
domain finite element
Definition: Simple.hpp:229
MoFEM::Core & cOre
Definition: Simple.hpp:37
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:219
char meshFileName[255]
mesh file name
Definition: Simple.hpp:235
int dIm
dimension of problem
Definition: Simple.hpp:236
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:231
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Common.hpp:147
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:213
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:214
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:220
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:215
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:216
const BitRefLevel bitLevel
Definition: Simple.hpp:211
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:217

◆ ~Simple()

MoFEM::Simple::~Simple ( )

Destructor.

Definition at line 52 of file Simple.cpp.

52  {
53 
54  if (dM != PETSC_NULL) {
55  ierr = DMDestroy(&dM);
56  CHKERRABORT(PETSC_COMM_WORLD, ierr);
57  }
58 }
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80

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

143  {
144 
145  Interface &m_field = cOre;
147  ierr = m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
148  verb);
149  CHKERRG(ierr);
150  boundaryFields.push_back(name);
152 }
MoFEM::Core & cOre
Definition: Simple.hpp:37
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:223
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80

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

175  {
176 
177  Interface &m_field = cOre;
179  ierr = m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
180  verb);
181  CHKERRG(ierr);
182  dataFields.push_back(name);
184 }
MoFEM::Core & cOre
Definition: Simple.hpp:37
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:225

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

Definition at line 123 of file Simple.cpp.

127  {
128 
129  Interface &m_field = cOre;
131  ierr = m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
132  verb);
133  CHKERRG(ierr);
134  domainFields.push_back(name);
136 }
MoFEM::Core & cOre
Definition: Simple.hpp:37
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:222

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

Definition at line 155 of file Simple.cpp.

159  {
160 
161  Interface &m_field = cOre;
163  ierr = m_field.add_field(name, space, base, nb_of_cooficients, tag_type, bh,
164  verb);
165  CHKERRG(ierr);
166  skeletonFields.push_back(name);
168 }
MoFEM::Core & cOre
Definition: Simple.hpp:37
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:224

◆ buildFields()

MoFEMErrorCode MoFEM::Simple::buildFields ( )

Build fields.

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

Definition at line 320 of file Simple.cpp.

320  {
321  Interface &m_field = cOre;
323  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
324  // take skin
325  {
326  Range domain_ents;
327  rval = m_field.get_moab().get_entities_by_dimension(meshSet, dIm,
328  domain_ents, true);
330  Skinner skin(&m_field.get_moab());
331  Range domain_skin;
332  rval = skin.find_skin(0, domain_ents, false, domain_skin);
334  // filter not owned entities, those are not on boundary
335  ParallelComm *pcomm =
336  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
337  Range proc_domain_skin;
338  rval =
339  pcomm->filter_pstatus(domain_skin, PSTATUS_SHARED | PSTATUS_MULTISHARED,
340  PSTATUS_NOT, -1, &proc_domain_skin);
342  // cerr << proc_domain_skin << endl;
343  // create boundary meshset
344  if (boundaryMeshset != 0) {
345  rval = m_field.get_moab().delete_entities(&boundaryMeshset, 1);
347  }
348  rval = m_field.get_moab().create_meshset(MESHSET_SET, boundaryMeshset);
350  rval = m_field.get_moab().add_entities(boundaryMeshset, proc_domain_skin);
352  for (int dd = 0; dd != dIm - 1; dd++) {
353  Range adj;
354  rval = m_field.get_moab().get_adjacencies(proc_domain_skin, dd, false,
355  adj, moab::Interface::UNION);
357  rval = m_field.get_moab().add_entities(boundaryMeshset, adj);
359  }
360  }
361  // Add entities to the fields
362  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
363  ierr = m_field.add_ents_to_field_by_dim(meshSet, dIm, domainFields[ff]);
364  CHKERRG(ierr);
365  ierr = m_field.synchronise_field_entities(domainFields[ff], 0);
366  CHKERRG(ierr);
367  }
368  for (unsigned int ff = 0; ff != dataFields.size(); ff++) {
369  ierr = m_field.add_ents_to_field_by_dim(meshSet, dIm, dataFields[ff]);
370  CHKERRG(ierr);
371  ierr = m_field.synchronise_field_entities(dataFields[ff], 0);
372  CHKERRG(ierr);
373  }
374  for (unsigned int ff = 0; ff != boundaryFields.size(); ff++) {
375  ierr = m_field.add_ents_to_field_by_dim(boundaryMeshset, dIm - 1,
376  boundaryFields[ff]);
377  CHKERRG(ierr);
378  ierr = m_field.synchronise_field_entities(boundaryFields[ff], 0);
379  CHKERRG(ierr);
380  }
381  for (unsigned int ff = 0; ff != skeletonFields.size(); ff++) {
382  ierr =
383  m_field.add_ents_to_field_by_dim(meshSet, dIm - 1, skeletonFields[ff]);
384  CHKERRG(ierr);
385  ierr = m_field.synchronise_field_entities(skeletonFields[ff], 0);
386  CHKERRG(ierr);
387  }
388  // Set order
389  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
390  if (fieldsOrder.find(domainFields[ff]) == fieldsOrder.end()) {
391  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
392  "Order for field not set %s", domainFields[ff].c_str());
393  }
394  int dds = 0;
395  const Field *field = m_field.get_field_structure(domainFields[ff]);
396  switch (field->getSpace()) {
397  case L2:
398  dds = dIm;
399  break;
400  case HDIV:
401  dds = 2;
402  break;
403  case HCURL:
404  dds = 1;
405  break;
406  case H1:
407  dds = 1;
408  break;
409  default:
410  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
411  "Huston we have a problem");
412  }
413  if (field->getSpace() == H1) {
414  ierr = m_field.set_field_order(meshSet, MBVERTEX, domainFields[ff], 1);
415  CHKERRG(ierr);
416  }
417  for (int dd = dds; dd <= dIm; dd++) {
418  Range ents;
419  ierr =
420  m_field.get_field_entities_by_dimension(domainFields[ff], dd, ents);
421  CHKERRG(ierr);
422  if (!fieldsOrder.at(domainFields[ff]).second.empty()) {
423  ents = intersect(ents, fieldsOrder.at(domainFields[ff]).second);
424  }
425  ierr = m_field.set_field_order(ents, domainFields[ff],
426  fieldsOrder.at(domainFields[ff]).first);
427  CHKERRG(ierr);
428  }
429  }
430  // Set order to data fields
431  for (unsigned int ff = 0; ff != dataFields.size(); ff++) {
432  if (fieldsOrder.find(dataFields[ff]) == fieldsOrder.end()) {
433  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
434  "Order for field not set %s", dataFields[ff].c_str());
435  }
436  int dds = 0;
437  const Field *field = m_field.get_field_structure(dataFields[ff]);
438  switch (field->getSpace()) {
439  case L2:
440  dds = dIm;
441  break;
442  case HDIV:
443  dds = 2;
444  break;
445  case HCURL:
446  dds = 1;
447  break;
448  case H1:
449  dds = 1;
450  break;
451  default:
452  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
453  "Huston we have a problem");
454  }
455  if (field->getSpace() == H1) {
456  ierr = m_field.set_field_order(meshSet, MBVERTEX, dataFields[ff], 1);
457  CHKERRG(ierr);
458  }
459  for (int dd = dds; dd <= dIm; dd++) {
460  Range ents;
461  ierr = m_field.get_field_entities_by_dimension(dataFields[ff], dd, ents);
462  CHKERRG(ierr);
463  ierr = m_field.set_field_order(ents, dataFields[ff],
464  fieldsOrder.at(dataFields[ff]).first);
465  CHKERRG(ierr);
466  }
467  }
468  // Set order to boundary
469  for (unsigned int ff = 0; ff != boundaryFields.size(); ff++) {
470  if (fieldsOrder.find(boundaryFields[ff]) == fieldsOrder.end()) {
471  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
472  "Order for field not set %s", boundaryFields[ff].c_str());
473  }
474  int dds = 0;
475  const Field *field = m_field.get_field_structure(boundaryFields[ff]);
476  switch (field->getSpace()) {
477  case L2:
478  dds = dIm - 1;
479  break;
480  case HDIV:
481  dds = 2;
482  break;
483  case HCURL:
484  dds = 1;
485  break;
486  case H1:
487  dds = 1;
488  break;
489  default:
490  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
491  "Huston we have a problem");
492  }
493  if (field->getSpace() == H1) {
494  ierr = m_field.set_field_order(meshSet, MBVERTEX, boundaryFields[ff], 1);
495  CHKERRG(ierr);
496  }
497  for (int dd = dds; dd <= dIm - 1; dd++) {
498  Range ents;
499  ierr =
500  m_field.get_field_entities_by_dimension(boundaryFields[ff], dd, ents);
501  CHKERRG(ierr);
502  ierr = m_field.set_field_order(ents, boundaryFields[ff],
503  fieldsOrder.at(boundaryFields[ff]).first);
504  CHKERRG(ierr);
505  }
506  }
507  // Set order to skeleton
508  for (unsigned int ff = 0; ff != skeletonFields.size(); ff++) {
509  if (fieldsOrder.find(skeletonFields[ff]) == fieldsOrder.end()) {
510  SETERRQ1(PETSC_COMM_WORLD, MOFEM_INVALID_DATA,
511  "Order for field not set %s", skeletonFields[ff].c_str());
512  }
513  int dds = 0;
514  const Field *field = m_field.get_field_structure(skeletonFields[ff]);
515  switch (field->getSpace()) {
516  case L2:
517  dds = dIm - 1;
518  break;
519  case HDIV:
520  dds = 2;
521  break;
522  case HCURL:
523  dds = 1;
524  break;
525  case H1:
526  dds = 1;
527  break;
528  default:
529  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
530  "Huston we have a problem");
531  }
532  if (field->getSpace() == H1) {
533  ierr = m_field.set_field_order(meshSet, MBVERTEX, skeletonFields[ff], 1);
534  CHKERRG(ierr);
535  }
536  for (int dd = dds; dd <= dIm - 1; dd++) {
537  Range ents;
538  ierr =
539  m_field.get_field_entities_by_dimension(skeletonFields[ff], dd, ents);
540  CHKERRG(ierr);
541  ierr = m_field.set_field_order(ents, skeletonFields[ff],
542  fieldsOrder.at(skeletonFields[ff]).first);
543  CHKERRG(ierr);
544  }
545  }
546  // Build fields
547  ierr = m_field.build_fields();
548  CHKERRG(ierr);
549  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFields, 0, 0, 0, 0);
551 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:533
MoFEM::Core & cOre
Definition: Simple.hpp:37
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:219
field with continuous normal traction
Definition: definitions.h:181
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
int dIm
dimension of problem
Definition: Simple.hpp:236
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
std::map< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:227
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:223
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:222
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:225
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:214
field with continuous tangents
Definition: definitions.h:180
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:318
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:220
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:224
continuous field
Definition: definitions.h:179
field with C-1 continuity
Definition: definitions.h:182

◆ buildFiniteElements()

MoFEMErrorCode MoFEM::Simple::buildFiniteElements ( )

Build finite elements.

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

Definition at line 553 of file Simple.cpp.

553  {
554  Interface &m_field = cOre;
556  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
557  // Add finite elements
558  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm, domainFE,
559  true);
560  CHKERR m_field.build_finite_elements(domainFE);
561  if (!boundaryFields.empty()) {
562  CHKERR m_field.add_ents_to_finite_element_by_dim(boundaryMeshset, dIm - 1,
563  boundaryFE, true);
564  CHKERR m_field.build_finite_elements(boundaryFE);
565  }
566  if (!skeletonFields.empty()) {
567  CHKERR m_field.add_ents_to_finite_element_by_dim(meshSet, dIm - 1,
568  skeletonFE, true);
569  CHKERR m_field.build_finite_elements(skeletonFE);
570  }
571  for (std::vector<std::string>::iterator fit = otherFEs.begin();
572  fit != otherFEs.end(); ++fit) {
573  CHKERR m_field.build_finite_elements(*fit);
574  }
575  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildFiniteElements, 0, 0, 0, 0);
577 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:230
std::string domainFE
domain finite element
Definition: Simple.hpp:229
MoFEM::Core & cOre
Definition: Simple.hpp:37
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:219
int dIm
dimension of problem
Definition: Simple.hpp:236
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:495
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:223
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:233
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:231
#define CHKERR
Inline error check.
Definition: definitions.h:614
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:220
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:215
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:224
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:439

◆ buildProblem()

MoFEMErrorCode MoFEM::Simple::buildProblem ( )

Build problem.

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

Definition at line 579 of file Simple.cpp.

579  {
580  Interface &m_field = cOre;
582  PetscLogEventBegin(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
583  ierr = m_field.build_adjacencies(bitLevel);
584  CHKERRG(ierr);
585  // Set problem by the DOFs on the fields rather that by adding DOFs on the
586  // elements
587  ierr = m_field.getInterface<ProblemsManager>()->buildProblemFromFields =
588  PETSC_TRUE;
589  ierr = DMSetUp(dM);
590  CHKERRG(ierr);
591  PetscLogEventEnd(MOFEM_EVENT_SimpleBuildProblem, 0, 0, 0, 0);
593 }
MoFEM::Core & cOre
Definition: Simple.hpp:37
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:216
const BitRefLevel bitLevel
Definition: Simple.hpp:211

◆ defineFiniteElements()

MoFEMErrorCode MoFEM::Simple::defineFiniteElements ( )

Define finite elements.

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

Definition at line 186 of file Simple.cpp.

186  {
187  Interface &m_field = cOre;
189  // Define finite elements
190  ierr = m_field.add_finite_element(domainFE);
191  CHKERRG(ierr);
192  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
193  ierr =
194  m_field.modify_finite_element_add_field_row(domainFE, domainFields[ff]);
195  CHKERRG(ierr);
196  ierr =
197  m_field.modify_finite_element_add_field_col(domainFE, domainFields[ff]);
198  CHKERRG(ierr);
199  ierr = m_field.modify_finite_element_add_field_data(domainFE,
200  domainFields[ff]);
201  CHKERRG(ierr);
202  }
203  for (unsigned int ff = 0; ff != dataFields.size(); ff++) {
204  ierr =
205  m_field.modify_finite_element_add_field_data(domainFE, dataFields[ff]);
206  CHKERRG(ierr);
207  }
208  if (!boundaryFields.empty()) {
209  ierr = m_field.add_finite_element(boundaryFE);
210  CHKERRG(ierr);
211  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
212  ierr = m_field.modify_finite_element_add_field_row(boundaryFE,
213  domainFields[ff]);
214  CHKERRG(ierr);
215  ierr = m_field.modify_finite_element_add_field_col(boundaryFE,
216  domainFields[ff]);
217  CHKERRG(ierr);
218  ierr = m_field.modify_finite_element_add_field_data(boundaryFE,
219  domainFields[ff]);
220  CHKERRG(ierr);
221  }
222  for (unsigned int ff = 0; ff != boundaryFields.size(); ff++) {
223  ierr = m_field.modify_finite_element_add_field_row(boundaryFE,
224  boundaryFields[ff]);
225  CHKERRG(ierr);
226  ierr = m_field.modify_finite_element_add_field_col(boundaryFE,
227  boundaryFields[ff]);
228  CHKERRG(ierr);
229  ierr = m_field.modify_finite_element_add_field_data(boundaryFE,
230  boundaryFields[ff]);
231  CHKERRG(ierr);
232  }
233  for (unsigned int ff = 0; ff != skeletonFields.size(); ff++) {
234  ierr = m_field.modify_finite_element_add_field_row(boundaryFE,
235  skeletonFields[ff]);
236  CHKERRG(ierr);
237  ierr = m_field.modify_finite_element_add_field_col(boundaryFE,
238  skeletonFields[ff]);
239  CHKERRG(ierr);
240  ierr = m_field.modify_finite_element_add_field_data(boundaryFE,
241  skeletonFields[ff]);
242  CHKERRG(ierr);
243  }
244  }
245  if (!skeletonFields.empty()) {
246  ierr = m_field.add_finite_element(skeletonFE);
247  CHKERRG(ierr);
248  for (unsigned int ff = 0; ff != domainFields.size(); ff++) {
249  ierr = m_field.modify_finite_element_add_field_row(skeletonFE,
250  domainFields[ff]);
251  CHKERRG(ierr);
252  ierr = m_field.modify_finite_element_add_field_col(skeletonFE,
253  domainFields[ff]);
254  CHKERRG(ierr);
255  ierr = m_field.modify_finite_element_add_field_data(skeletonFE,
256  domainFields[ff]);
257  CHKERRG(ierr);
258  }
259  for (unsigned int ff = 0; ff != boundaryFields.size(); ff++) {
260  ierr = m_field.modify_finite_element_add_field_row(skeletonFE,
261  boundaryFields[ff]);
262  CHKERRG(ierr);
263  ierr = m_field.modify_finite_element_add_field_col(skeletonFE,
264  boundaryFields[ff]);
265  CHKERRG(ierr);
266  ierr = m_field.modify_finite_element_add_field_data(skeletonFE,
267  boundaryFields[ff]);
268  CHKERRG(ierr);
269  }
270  for (unsigned int ff = 0; ff != skeletonFields.size(); ff++) {
271  ierr = m_field.modify_finite_element_add_field_row(skeletonFE,
272  skeletonFields[ff]);
273  CHKERRG(ierr);
274  ierr = m_field.modify_finite_element_add_field_col(skeletonFE,
275  skeletonFields[ff]);
276  CHKERRG(ierr);
277  ierr = m_field.modify_finite_element_add_field_data(skeletonFE,
278  skeletonFields[ff]);
279  CHKERRG(ierr);
280  }
281  }
283 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:230
std::string domainFE
domain finite element
Definition: Simple.hpp:229
MoFEM::Core & cOre
Definition: Simple.hpp:37
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:223
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:231
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:222
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:225
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:224

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

285  {
286  Interface &m_field = cOre;
288  if (dM != PETSC_NULL) {
289  CHKERR DMDestroy(&dM);
290  }
291  // Create dm instance
292  CHKERR DMCreate(m_field.get_comm(), &dM);
293  CHKERR DMSetType(dM, "DMMOFEM");
294  // set dm data structure which created mofem data structures
295  CHKERR DMMoFEMCreateMoFEM(dM, &m_field, "SimpleProblem", bitLevel);
296  CHKERR DMSetFromOptions(dM);
298  if (!boundaryFields.empty()) {
300  }
301  if (!skeletonFields.empty()) {
303  }
304  for (std::vector<std::string>::iterator fit = otherFEs.begin();
305  fit != otherFEs.end(); ++fit) {
306  CHKERR DMMoFEMAddElement(dM, fit->c_str());
307  }
308  CHKERR DMMoFEMSetIsPartitioned(dM, is_partitioned);
310 }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:230
std::string domainFE
domain finite element
Definition: Simple.hpp:229
MoFEM::Core & cOre
Definition: Simple.hpp:37
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:432
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:495
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:223
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:233
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:231
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:147
#define CHKERR
Inline error check.
Definition: definitions.h:614
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:224
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:439
const BitRefLevel bitLevel
Definition: Simple.hpp:211
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:845

◆ getBoundaryFEName() [1/2]

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

Definition at line 198 of file Simple.hpp.

198 { return boundaryFE; }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:230

◆ getBoundaryFEName() [2/2]

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

Definition at line 202 of file Simple.hpp.

202 { return boundaryFE; }
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:230

◆ getDim()

int MoFEM::Simple::getDim ( ) const
Examples:
simple_elasticity.cpp.

Definition at line 196 of file Simple.hpp.

196 { return dIm; }
int dIm
dimension of problem
Definition: Simple.hpp:236

◆ getDM()

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

Get DM.

Parameters
dmdiscrete manager
Returns
error code
Examples:
hello_world.cpp, mesh_smoothing.cpp, simple_elasticity.cpp, and simple_interface.cpp.

Definition at line 605 of file Simple.cpp.

605  {
607  PetscObjectReference((PetscObject)dM);
608  *dm = dM;
610 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526

◆ getDomainFEName() [1/2]

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

Definition at line 197 of file Simple.hpp.

197 { return domainFE; }
std::string domainFE
domain finite element
Definition: Simple.hpp:229

◆ getDomainFEName() [2/2]

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

Definition at line 201 of file Simple.hpp.

201 { return domainFE; }
std::string domainFE
domain finite element
Definition: Simple.hpp:229

◆ getOptions()

MoFEMErrorCode MoFEM::Simple::getOptions ( )

get options

Returns
error code
Examples:
hello_world.cpp, mesh_smoothing.cpp, simple_elasticity.cpp, and simple_interface.cpp.

Definition at line 60 of file Simple.cpp.

60  {
61 
62  PetscBool flg = PETSC_TRUE;
64  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "Simple interface options",
65  "none");
66  CHKERRG(ierr);
67  ierr = PetscOptionsString("-file_name", "file name", "", "mesh.h5m",
68  meshFileName, 255, &flg);
69  CHKERRG(ierr);
70  ierr = PetscOptionsEnd();
71  CHKERRG(ierr);
73 }
char meshFileName[255]
mesh file name
Definition: Simple.hpp:235
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80

◆ getOtherFiniteElements()

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

Definition at line 205 of file Simple.hpp.

205  {
206  return otherFEs;
207  }
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:233

◆ getSkeletonFEName() [1/2]

const std::string& MoFEM::Simple::getSkeletonFEName ( ) const
Examples:
hello_world.cpp.

Definition at line 199 of file Simple.hpp.

199 { return skeletonFE; }
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:231

◆ getSkeletonFEName() [2/2]

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

Definition at line 203 of file Simple.hpp.

203 { return skeletonFE; }
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:231

◆ loadFile()

MoFEMErrorCode MoFEM::Simple::loadFile ( )

Load mesh file.

Parameters
field_namefile name
Returns
error code
Examples:
hello_world.cpp, mesh_smoothing.cpp, simple_elasticity.cpp, and simple_interface.cpp.

Definition at line 75 of file Simple.cpp.

75  {
76  Interface &m_field = cOre;
78  PetscLogEventBegin(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
79  // This is a case of distributed mesh and algebra. In that case each processor
80  // keep only part of the problem.
81  const char *option;
82  option = "PARALLEL=READ_PART;"
83  "PARALLEL_RESOLVE_SHARED_ENTS;"
84  "PARTITION=PARALLEL_PARTITION;";
85  rval = m_field.get_moab().load_file(meshFileName, 0, option);
87  ierr = m_field.rebuild_database();
88  CHKERRG(ierr);
89  // determine problem dimension
90  if (dIm == -1) {
91  int nb_ents_3d;
92  rval = m_field.get_moab().get_number_entities_by_dimension(
93  meshSet, 3, nb_ents_3d, true);
95  if (nb_ents_3d > 0) {
96  dIm = 3;
97  } else {
98  int nb_ents_2d;
99  rval = m_field.get_moab().get_number_entities_by_dimension(
100  meshSet, 2, nb_ents_2d, true);
102  if (nb_ents_2d > 0) {
103  dIm = 2;
104  } else {
105  dIm = 1;
106  }
107  }
108  }
109  Range ents;
110  ierr = m_field.get_moab().get_entities_by_dimension(meshSet, dIm, ents, true);
111  ierr = m_field.getInterface<BitRefManager>()->setBitRefLevel(ents, bitLevel,
112  false);
113  CHKERRG(ierr);
114  ParallelComm *pcomm =
115  ParallelComm::get_pcomm(&m_field.get_moab(), MYPCOMM_INDEX);
116  if (pcomm == NULL)
117  pcomm = new ParallelComm(&m_field.get_moab(), m_field.get_comm());
118  PetscLogEventEnd(MOFEM_EVENT_SimpleLoadMesh, 0, 0, 0, 0);
120 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:533
MoFEM::Core & cOre
Definition: Simple.hpp:37
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:219
char meshFileName[255]
mesh file name
Definition: Simple.hpp:235
int dIm
dimension of problem
Definition: Simple.hpp:236
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:213
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:318
const BitRefLevel bitLevel
Definition: Simple.hpp:211

◆ query_interface()

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

Implements MoFEM::UnknownInterface.

Definition at line 28 of file Simple.cpp.

29  {
31  *iface = NULL;
32  if (uuid == IDD_MOFEMSimple) {
33  *iface = const_cast<Simple *>(this);
35  }
36  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
38 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
static const MOFEMuuid IDD_MOFEMSimple
Definition: Simple.hpp:27
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:40

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

Definition at line 312 of file Simple.cpp.

313  {
315  fieldsOrder[field_name] =
316  std::pair<int, Range>(order, ents == NULL ? Range() : Range(*ents));
318 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:519
std::map< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:227
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526

◆ setUp()

MoFEMErrorCode MoFEM::Simple::setUp ( )

Setup problem.

Returns
error code
Examples:
hello_world.cpp, and simple_interface.cpp.

Definition at line 595 of file Simple.cpp.

595  {
603 }
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:553
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:186
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:495
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:579
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:320
#define CHKERR
Inline error check.
Definition: definitions.h:614
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:439
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:285

Member Data Documentation

◆ bitLevel

const BitRefLevel MoFEM::Simple::bitLevel
private

Definition at line 211 of file Simple.hpp.

◆ boundaryFE

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

boundary finite element

Definition at line 230 of file Simple.hpp.

◆ boundaryFields

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

boundary fields

Definition at line 223 of file Simple.hpp.

◆ boundaryMeshset

EntityHandle MoFEM::Simple::boundaryMeshset
private

meshset with boundary

Definition at line 220 of file Simple.hpp.

◆ cOre

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

Definition at line 37 of file Simple.hpp.

◆ dataFields

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

Data fields.

Definition at line 225 of file Simple.hpp.

◆ dIm

int MoFEM::Simple::dIm
private

dimension of problem

Definition at line 236 of file Simple.hpp.

◆ dM

DM MoFEM::Simple::dM
private

Definition at line 238 of file Simple.hpp.

◆ domainFE

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

domain finite element

Definition at line 229 of file Simple.hpp.

◆ domainFields

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

domain fields

Definition at line 222 of file Simple.hpp.

◆ fieldsOrder

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

fields order

Definition at line 227 of file Simple.hpp.

◆ meshFileName

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

mesh file name

Definition at line 235 of file Simple.hpp.

◆ meshSet

EntityHandle MoFEM::Simple::meshSet
private

domain meshset

Definition at line 219 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFields

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFields
private

Definition at line 214 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildFiniteElements

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildFiniteElements
private

Definition at line 215 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleBuildProblem

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleBuildProblem
private

Definition at line 216 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleKSPSolve

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleKSPSolve
private

Definition at line 217 of file Simple.hpp.

◆ MOFEM_EVENT_SimpleLoadMesh

PetscLogEvent MoFEM::Simple::MOFEM_EVENT_SimpleLoadMesh
private

Definition at line 213 of file Simple.hpp.

◆ nameOfProblem

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

problem name

Definition at line 221 of file Simple.hpp.

◆ otherFEs

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

Other finite elements.

Definition at line 233 of file Simple.hpp.

◆ skeletonFE

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

skeleton finite element

Definition at line 231 of file Simple.hpp.

◆ skeletonFields

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

fields on the skeleton

Definition at line 224 of file Simple.hpp.


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