v0.9.1
Classes | Functions | Variables
Distributed mesh manager

Implementation of PETSc DM, managing interactions between mesh data structures and vectors and matrices. More...

Collaboration diagram for Distributed mesh manager:

Classes

struct  MoFEM::DMCtx
 PETSc Discrete Manager data structure. More...
 
struct  DMMGViaApproxOrdersCtx
 Structure for DM for multi-grid via approximation orders. More...
 

Functions

PetscErrorCode MoFEM::DMRegister_MoFEM (const char sname[])
 Register MoFEM problem. More...
 
PetscErrorCode MoFEM::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. More...
 
PetscErrorCode MoFEM::DMMoFEMCreateSubDM (DM subdm, DM dm, const char problem_name[])
 Must be called by user to set Sub DM MoFEM data structures. More...
 
PetscErrorCode MoFEM::DMoFEMGetInterfacePtr (DM dm, MoFEM::Interface **m_field_ptr)
 Get pointer to MoFEM::Interface. More...
 
PetscErrorCode MoFEM::DMMoFEMGetProblemPtr (DM dm, const MoFEM::Problem **problem_ptr)
 Get pointer to problem data structure. More...
 
PetscErrorCode MoFEM::DMMoFEMSetSquareProblem (DM dm, PetscBool square_problem)
 set squared problemIt if true is assumed that matrix has the same indexing on rows and columns. This reduces interprocessor communication. More...
 
PetscErrorCode MoFEM::DMMoFEMGetSquareProblem (DM dm, PetscBool *square_problem)
 get squared problemIt if true is assumed that matrix has the same indexing on rows and columns. This reduces interprocessor communication. More...
 
PetscErrorCode MoFEM::DMMoFEMResolveSharedFiniteElements (DM dm, const char fe_name[])
 Resolve shared entities. More...
 
PetscErrorCode MoFEM::DMMoFEMGetProblemFiniteElementLayout (DM dm, const char fe_name[], PetscLayout *layout)
 Get finite elements layout in the problem. More...
 
PetscErrorCode MoFEM::DMMoFEMAddElement (DM dm, const char fe_name[])
 add element to dm More...
 
PetscErrorCode MoFEM::DMMoFEMUnSetElement (DM dm, const char fe_name[])
 unset element from dm More...
 
PetscErrorCode MoFEM::DMoFEMMeshToLocalVector (DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
 set local (or ghosted) vector values on mesh for partition only More...
 
PetscErrorCode MoFEM::DMoFEMMeshToGlobalVector (DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
 set ghosted vector values on all existing mesh entities More...
 
PetscErrorCode MoFEM::DMoFEMPreProcessFiniteElements (DM dm, MoFEM::FEMethod *method)
 execute finite element method for each element in dm (problem) More...
 
PetscErrorCode MoFEM::DMoFEMPostProcessFiniteElements (DM dm, MoFEM::FEMethod *method)
 execute finite element method for each element in dm (problem) More...
 
PetscErrorCode MoFEM::DMoFEMLoopFiniteElementsUpAndLowRank (DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank)
 Executes FEMethod for finite elements in DM. More...
 
PetscErrorCode MoFEM::DMoFEMLoopFiniteElementsUpAndLowRank (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, int low_rank, int up_rank)
 Executes FEMethod for finite elements in DM. More...
 
PetscErrorCode MoFEM::DMoFEMLoopFiniteElements (DM dm, const char fe_name[], MoFEM::FEMethod *method)
 Executes FEMethod for finite elements in DM. More...
 
PetscErrorCode MoFEM::DMoFEMLoopFiniteElements (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method)
 Executes FEMethod for finite elements in DM. More...
 
PetscErrorCode MoFEM::DMoFEMLoopDofs (DM dm, const char field_name[], MoFEM::DofMethod *method)
 execute method for dofs on field in problem More...
 
PetscErrorCode MoFEM::DMMoFEMKSPSetComputeRHS (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 Set compute operator for KSP solver via sub-matrix and IS. More...
 
PetscErrorCode MoFEM::DMMoFEMKSPSetComputeRHS (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set KSP right hand side evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMKSPSetComputeOperators (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 Set KSP operators and push mofem finite element methods. More...
 
PetscErrorCode MoFEM::DMMoFEMKSPSetComputeOperators (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 Set KSP operators and push mofem finite element methods. More...
 
PetscErrorCode MoFEM::DMMoFEMSNESSetFunction (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 set SNES residual evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMSNESSetFunction (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set SNES residual evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMSNESSetJacobian (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 set SNES Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMSNESSetJacobian (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set SNES Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetIFunction (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 set TS implicit function evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetIFunction (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set TS implicit function evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetIJacobian (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set TS Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetIJacobian (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 set TS Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetI2Function (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set TS implicit function evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetI2Jacobian (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set TS Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMGetKspCtx (DM dm, MoFEM::KspCtx **ksp_ctx)
 get MoFEM::KspCtx data structure More...
 
PetscErrorCode MoFEM::DMMoFEMGetKspCtx (DM dm, const boost::shared_ptr< MoFEM::KspCtx > &ksp_ctx)
 get MoFEM::KspCtx data structure More...
 
PetscErrorCode MoFEM::DMMoFEMSetKspCtx (DM dm, boost::shared_ptr< MoFEM::KspCtx > &ksp_ctx)
 set MoFEM::KspCtx data structure More...
 
PetscErrorCode MoFEM::DMMoFEMGetSnesCtx (DM dm, MoFEM::SnesCtx **snes_ctx)
 get MoFEM::SnesCtx data structure More...
 
PetscErrorCode MoFEM::DMMoFEMGetSnesCtx (DM dm, const boost::shared_ptr< MoFEM::SnesCtx > &snes_ctx)
 get MoFEM::SnesCtx data structure More...
 
PetscErrorCode MoFEM::DMMoFEMSetSnesCtx (DM dm, boost::shared_ptr< MoFEM::SnesCtx > &snes_ctx)
 Set MoFEM::SnesCtx data structure. More...
 
PetscErrorCode MoFEM::DMMoFEMGetTsCtx (DM dm, MoFEM::TsCtx **ts_ctx)
 get MoFEM::TsCtx data structure More...
 
PetscErrorCode MoFEM::DMMoFEMGetTsCtx (DM dm, const boost::shared_ptr< MoFEM::TsCtx > &ts_ctx)
 get MoFEM::TsCtx data structure More...
 
PetscErrorCode MoFEM::DMMoFEMSetTsCtx (DM dm, boost::shared_ptr< MoFEM::TsCtx > &ts_ctx)
 Set MoFEM::TsCtx data structureIt take over pointer, do not delete it, DM will destroy pointer when is destroyed. More...
 
PetscErrorCode MoFEM::DMMoFEMSetIsPartitioned (DM dm, PetscBool is_partitioned)
 
PetscErrorCode MoFEM::DMMoFEMGetIsPartitioned (DM dm, PetscBool *is_partitioned)
 
PetscErrorCode MoFEM::DMSetOperators_MoFEM (DM dm)
 Set operators for MoFEM dm. More...
 
PetscErrorCode MoFEM::DMCreate_MoFEM (DM dm)
 Create dm data structure with MoFEM data structure. More...
 
PetscErrorCode MoFEM::DMDestroy_MoFEM (DM dm)
 Destroys dm with MoFEM data structure. More...
 
PetscErrorCode MoFEM::DMCreateGlobalVector_MoFEM (DM dm, Vec *g)
 DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM. More...
 
PetscErrorCode MoFEM::DMCreateGlobalVector_MoFEM (DM dm, SmartPetscObj< Vec > &g_ptr)
 DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM. More...
 
PetscErrorCode MoFEM::DMCreateLocalVector_MoFEM (DM dm, Vec *l)
 DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM. More...
 
PetscErrorCode MoFEM::DMCreateMatrix_MoFEM (DM dm, Mat *M)
 
PetscErrorCode MoFEM::DMCreateMatrix_MoFEM (DM dm, SmartPetscObj< Mat > &M)
 
PetscErrorCode MoFEM::DMSetFromOptions_MoFEM (DM dm)
 
PetscErrorCode MoFEM::DMSetUp_MoFEM (DM dm)
 
PetscErrorCode MoFEM::DMSubDMSetUp_MoFEM (DM subdm)
 
PetscErrorCode MoFEM::DMMoFEMAddSubFieldRow (DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
 
PetscErrorCode MoFEM::DMMoFEMAddSubFieldCol (DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
 
PetscErrorCode MoFEM::DMMoFEMGetIsSubDM (DM dm, PetscBool *is_sub_dm)
 
PetscErrorCode MoFEM::DMMoFEMAddRowCompositeProblem (DM dm, const char prb_name[])
 Add problem to composite DM on rowThis create block on row with DOFs from problem of given name. More...
 
PetscErrorCode MoFEM::DMMoFEMAddColCompositeProblem (DM dm, const char prb_name[])
 Add problem to composite DM on colThis create block on col with DOFs from problem of given name. More...
 
PetscErrorCode MoFEM::DMMoFEMGetIsCompDM (DM dm, PetscBool *is_comp_dm)
 Get if this DM is composite DM. More...
 
PetscErrorCode MoFEM::DMGlobalToLocalBegin_MoFEM (DM dm, Vec, InsertMode, Vec)
 
PetscErrorCode MoFEM::DMGlobalToLocalEnd_MoFEM (DM dm, Vec, InsertMode, Vec)
 
PetscErrorCode MoFEM::DMLocalToGlobalBegin_MoFEM (DM, Vec, InsertMode, Vec)
 
PetscErrorCode MoFEM::DMLocalToGlobalEnd_MoFEM (DM, Vec, InsertMode, Vec)
 
PetscErrorCode MoFEM::DMCreateFieldIS_MoFEM (DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
 
PetscErrorCode MoFEM::DMMoFEMGetFieldIS (DM dm, RowColData rc, const char field_name[], IS *is)
 get field is in the problem More...
 
MoFEMErrorCode DMMGViaApproxOrdersPushBackCoarseningIS (DM, IS is, Mat A, Mat *subA, bool create_sub_matrix, bool shell_sub_a)
 Push back coarsening level to MG via approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS (DM)
 Pop is form MG via approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS (DM)
 Clear approximation orders. More...
 
MoFEMErrorCode DMRegister_MGViaApproxOrders (const char sname[])
 Register DM for Multi-Grid via approximation orders. More...
 
MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders (DM dm, Mat *M)
 Create matrix for Multi-Grid via approximation orders. More...
 
MoFEMErrorCode DMCoarsen_MGViaApproxOrders (DM dm, MPI_Comm comm, DM *dmc)
 Coarsen DM. More...
 

Variables

auto MoFEM::smartCreateDMMatrix
 Get smart matrix from DM. More...
 
auto MoFEM::smartCreateDMVector
 Get smart vector from DM. More...
 

Detailed Description

Implementation of PETSc DM, managing interactions between mesh data structures and vectors and matrices.

DM objects are used to manage communication between the algebraic structures in PETSc (Vec and Mat) and mesh data structures in PDE-based (or other) simulations.

DM is abstract interface, here is it particular implementation for MoFEM code.

Function Documentation

◆ DMCoarsen_MGViaApproxOrders()

MoFEMErrorCode DMCoarsen_MGViaApproxOrders ( DM  dm,
MPI_Comm  comm,
DM *  dmc 
)

Coarsen DM.

Not used directly by user

Parameters
dmDistributed mesh data structure
commCommunicator
dmcCoarse distributed mesh data structure
Returns
Error code

Definition at line 485 of file PCMGSetUpViaApproxOrders.cpp.

485  {
486  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
488  GET_DM_FIELD(dm);
489  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
490  CHKERR DMCreate(comm, dmc);
491  (*dmc)->data = dm->data;
492  DMType type;
493  CHKERR DMGetType(dm, &type);
494  CHKERR DMSetType(*dmc, type);
495  CHKERR PetscObjectReference((PetscObject)(*dmc));
496  PetscInfo1(dm, "Coarsen DMMGViaApproxOrders leveldown = %d\n", dm->leveldown);
498 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMCreate_MoFEM()

PetscErrorCode MoFEM::DMCreate_MoFEM ( DM  dm)

Create dm data structure with MoFEM data structure.

Definition at line 76 of file DMMMoFEM.cpp.

76  {
77  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
79  dm->data = new DMCtx();
82 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:54
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMCreateFieldIS_MoFEM()

PetscErrorCode MoFEM::DMCreateFieldIS_MoFEM ( DM  dm,
PetscInt *  numFields,
char ***  fieldNames,
IS **  fields 
)

Creates a set of IS objects with the global indices of dofs for each field

Parameters
dmThe number of fields (or NULL if not requested)

Output:

Parameters
numFieldsThe number of fields (or NULL if not requested)
fieldNamesThe name for each field (or NULL if not requested)
fieldsThe global indices for each field (or NULL if not requested)
Returns
error code
Note
The user is responsible for freeing all requested arrays. In particular, every entry of names should be freed with PetscFree(), every entry of fields should be destroyed with ISDestroy(), and both arrays should be freed with PetscFree().

Definition at line 1285 of file DMMMoFEM.cpp.

1286  {
1287  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1289 
1290  if (numFields) {
1291  *numFields = 0;
1292  }
1293  if (fieldNames) {
1294  *fieldNames = NULL;
1295  }
1296  if (fields) {
1297  *fields = NULL;
1298  }
1299 
1300  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1301  auto fields_ptr = dm_field->mField_ptr->get_fields();
1302  Field_multiIndex::iterator fit, hi_fit;
1303  fit = fields_ptr->begin();
1304  hi_fit = fields_ptr->end();
1305  *numFields = std::distance(fit, hi_fit);
1306 
1307  if (fieldNames) {
1308  CHKERR PetscMalloc1(*numFields, fieldNames);
1309  }
1310  if (fields) {
1311  CHKERR PetscMalloc1(*numFields, fields);
1312  }
1313 
1314  for (int f = 0; fit != hi_fit; fit++, f++) {
1315  if (fieldNames) {
1316  CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1317  (char **)&(*fieldNames)[f]);
1318  }
1319  if (fields) {
1320  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1321  ->isCreateProblemFieldAndRank(
1322  dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1323  fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1324  }
1325  }
1326 
1328 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMCreateGlobalVector_MoFEM() [1/2]

PetscErrorCode MoFEM::DMCreateGlobalVector_MoFEM ( DM  dm,
Vec *  g 
)

DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM.

Examples
elasticity_mixed_formulation.cpp, ep.cpp, and reaction_diffusion_equation.cpp.

Definition at line 1026 of file DMMMoFEM.cpp.

1026  {
1027  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1029  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1030  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1031  dm_field->problemName, COL, g);
1033 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
#define CHKERR
Inline error check.
Definition: definitions.h:602

◆ DMCreateGlobalVector_MoFEM() [2/2]

PetscErrorCode MoFEM::DMCreateGlobalVector_MoFEM ( DM  dm,
SmartPetscObj< Vec > &  g_ptr 
)

DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM.

Definition at line 1035 of file DMMMoFEM.cpp.

1035  {
1036  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1038  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1039  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1040  dm_field->problemName, COL, g_ptr);
1042 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
#define CHKERR
Inline error check.
Definition: definitions.h:602

◆ DMCreateLocalVector_MoFEM()

PetscErrorCode MoFEM::DMCreateLocalVector_MoFEM ( DM  dm,
Vec *  l 
)

DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM.

Definition at line 1044 of file DMMMoFEM.cpp.

1044  {
1045  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1047  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1048  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1049  dm_field->problemName, COL, l);
1051 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
FTensor::Index< 'l', 2 > l
Definition: ContactOps.hpp:29

◆ DMCreateMatrix_MGViaApproxOrders()

MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders ( DM  dm,
Mat *  M 
)

Create matrix for Multi-Grid via approximation orders.

Not used directly by user

Parameters
dmDistributed mesh data structure
MMatrix
Returns
Error code

Definition at line 454 of file PCMGSetUpViaApproxOrders.cpp.

454  {
455 
456  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
458  GET_DM_FIELD(dm);
459 
460  int leveldown = dm->leveldown;
461 
462  if (dm_field->kspOperators.empty()) {
464  } else {
465  MPI_Comm comm;
466  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
467  if (dm_field->kspOperators.empty()) {
468  SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
469  "data inconsistency, operator can not be set");
470  }
471  if (static_cast<int>(dm_field->kspOperators.size()) < leveldown) {
472  SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
473  "data inconsistency, no IS for that level");
474  }
475  *M = dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown];
476  CHKERR PetscObjectReference((PetscObject)*M);
477  }
478 
479  PetscInfo1(dm, "Create Matrix DMMGViaApproxOrders leveldown = %d\n",
480  leveldown);
481 
483 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1053
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:602
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:72
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMCreateMatrix_MoFEM() [1/2]

PetscErrorCode MoFEM::DMCreateMatrix_MoFEM ( DM  dm,
Mat *  M 
)

DMShellSetCreateMatrix

sets the routine to create a matrix associated with the shell DM

Examples
dm_build_partitioned_mesh.cpp, elasticity_mixed_formulation.cpp, ep.cpp, and reaction_diffusion_equation.cpp.

Definition at line 1053 of file DMMMoFEM.cpp.

1053  {
1054  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1056  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1057  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1058  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1059  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1060  M);
1061  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1062  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1063  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1064  M);
1065  } else {
1066  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1067  "Matrix type not implemented");
1068  }
1070 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:72
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMCreateMatrix_MoFEM() [2/2]

PetscErrorCode MoFEM::DMCreateMatrix_MoFEM ( DM  dm,
SmartPetscObj< Mat > &  M 
)

DMShellSetCreateMatrix

sets the routine to create a matrix associated with the shell DM

Definition at line 1072 of file DMMMoFEM.cpp.

1072  {
1073  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1075  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1076  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1077  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1078  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1079  M);
1080  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1081  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1082  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1083  M);
1084  } else {
1085  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1086  "Matrix type not implemented");
1087  }
1089 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:72
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMDestroy_MoFEM()

PetscErrorCode MoFEM::DMDestroy_MoFEM ( DM  dm)

Destroys dm with MoFEM data structure.

destroy the MoFEM structure

Definition at line 84 of file DMMMoFEM.cpp.

84  {
85  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
86  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
88  if (static_cast<DMCtx *>(dm->data)->referenceNumber == 0) {
89  if (dm_field->destroyProblem) {
90 
91  if (dm_field->mField_ptr->check_problem(dm_field->problemName)) {
92  dm_field->mField_ptr->delete_problem(dm_field->problemName);
93  } // else problem has to be deleted by the user
94 
95  }
96 
97  delete static_cast<DMCtx *>(dm->data);
98 
99  } else
100  --static_cast<DMCtx *>(dm->data)->referenceNumber;
101 
103 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMGlobalToLocalBegin_MoFEM()

PetscErrorCode MoFEM::DMGlobalToLocalBegin_MoFEM ( DM  dm,
Vec  g,
InsertMode  mode,
Vec  l 
)

DMShellSetGlobalToLocal

the routine that begins the global to local scatter

Definition at line 1213 of file DMMMoFEM.cpp.

1214  {
1215  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1217  CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1219 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMGlobalToLocalEnd_MoFEM()

PetscErrorCode MoFEM::DMGlobalToLocalEnd_MoFEM ( DM  dm,
Vec  g,
InsertMode  mode,
Vec  l 
)

DMShellSetGlobalToLocal

the routine that begins the global to local scatter

Definition at line 1221 of file DMMMoFEM.cpp.

1221  {
1223  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1225 
1226  CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1227 
1228  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1229  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1230  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1231 
1232  double *array_loc, *array_glob;
1233  CHKERR VecGetArray(l, &array_loc);
1234  CHKERR VecGetArray(g, &array_glob);
1235  switch (mode) {
1236  case INSERT_VALUES:
1237  cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1238  break;
1239  case ADD_VALUES:
1240  cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1241  break;
1242  default:
1243  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1244  }
1245  CHKERR VecRestoreArray(l, &array_loc);
1246  CHKERR VecRestoreArray(g, &array_glob);
1248 }
void cblas_daxpy(const int N, const double alpha, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_daxpy.c:11
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
void cblas_dcopy(const int N, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_dcopy.c:11
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
FTensor::Index< 'l', 2 > l
Definition: ContactOps.hpp:29

◆ DMLocalToGlobalBegin_MoFEM()

PetscErrorCode MoFEM::DMLocalToGlobalBegin_MoFEM ( DM  dm,
Vec  l,
InsertMode  mode,
Vec  g 
)

DMShellSetLocalToGlobal

the routine that begins the local to global scatter

DMShellSetLocalToGlobal

the routine that ends the local to global scatter

Definition at line 1250 of file DMMMoFEM.cpp.

1251  {
1252 
1253  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1255 
1256  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1257  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1258  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1259 
1260  double *array_loc, *array_glob;
1261  CHKERR VecGetArray(l, &array_loc);
1262  CHKERR VecGetArray(g, &array_glob);
1263  switch (mode) {
1264  case INSERT_VALUES:
1265  cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1266  break;
1267  case ADD_VALUES:
1268  cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1269  break;
1270  default:
1271  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1272  }
1273  CHKERR VecRestoreArray(l, &array_loc);
1274  CHKERR VecRestoreArray(g, &array_glob);
1275 
1277 }
void cblas_daxpy(const int N, const double alpha, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_daxpy.c:11
void cblas_dcopy(const int N, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_dcopy.c:11
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
FTensor::Index< 'l', 2 > l
Definition: ContactOps.hpp:29

◆ DMLocalToGlobalEnd_MoFEM()

PetscErrorCode MoFEM::DMLocalToGlobalEnd_MoFEM ( DM  ,
Vec  l,
InsertMode  mode,
Vec  g 
)

DMShellSetLocalToGlobal

the routine that ends the local to global scatter

Definition at line 1279 of file DMMMoFEM.cpp.

1279  {
1280  //
1283 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMGViaApproxOrdersClearCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS ( DM  )

Clear approximation orders.

Parameters
DMdm
Returns
Error code

Definition at line 300 of file PCMGSetUpViaApproxOrders.cpp.

300  {
301  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
303  GET_DM_FIELD(dm);
304  CHKERR dm_field->destroyCoarseningIS();
305  PetscInfo(dm, "Clear DMs data structures\n");
307 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMGViaApproxOrdersPopBackCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS ( DM  )

Pop is form MG via approximation orders.

Parameters
DMdm
ispop back IS
Returns
error code

Definition at line 284 of file PCMGSetUpViaApproxOrders.cpp.

284  {
285  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
287  GET_DM_FIELD(dm);
288  if (dm_field->coarseningIS.back()) {
289  CHKERR ISDestroy(&dm_field->coarseningIS.back());
290  dm_field->coarseningIS.pop_back();
291  }
292  if (dm_field->kspOperators.back()) {
293  CHKERR MatDestroy(&dm_field->kspOperators.back());
294  }
295  dm_field->kspOperators.pop_back();
296  PetscInfo(dm, "Pop back IS to DMMGViaApproxOrders\n");
298 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMGViaApproxOrdersPushBackCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersPushBackCoarseningIS ( DM  ,
IS  is,
Mat  A,
Mat *  subA,
bool  create_sub_matrix,
bool  shell_sub_a 
)

Push back coarsening level to MG via approximation orders.

Parameters
DMdiscrete manager
isPush back IS used for coarsening
AGet sub-matrix of A using is (that sets operators for coarsening levels)
subAReturning pointer to created sub matrix
subAIf true create sub matrix, otherwise in subA has to be valid pointer to subA
Returns
Error code

Definition at line 230 of file PCMGSetUpViaApproxOrders.cpp.

233  {
234  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
236  GET_DM_FIELD(dm);
237  dm_field->coarseningIS.push_back(is);
238  dm_field->shellMatrixCtxPtr.push_back(new PCMGSubMatrixCtx_private(A, is));
239  if (is) {
240  CHKERR PetscObjectReference((PetscObject)is);
241  }
242  if (is) {
243  IS is2 = is;
244  if (dm_field->aO) {
245  CHKERR ISDuplicate(is, &is2);
246  CHKERR ISCopy(is, is2);
247  CHKERR AOApplicationToPetscIS(dm_field->aO, is2);
248  }
249  if (create_sub_matrix) {
250  if (shell_sub_a) {
251  int n, N;
252  CHKERR ISGetSize(is, &N);
253  CHKERR ISGetLocalSize(is, &n);
254  MPI_Comm comm;
255  CHKERR PetscObjectGetComm((PetscObject)A, &comm);
256  CHKERR MatCreateShell(comm, n, n, N, N,
257  &(dm_field->shellMatrixCtxPtr.back()), subA);
258  CHKERR MatShellSetOperation(*subA, MATOP_MULT,
259  (void (*)(void))sub_mat_mult);
260  CHKERR MatShellSetOperation(*subA, MATOP_MULT_ADD,
261  (void (*)(void))sub_mat_mult_add);
262  CHKERR MatShellSetOperation(*subA, MATOP_SOR,
263  (void (*)(void))sub_mat_sor);
264  } else {
265  #if PETSC_VERSION_GE(3, 8, 0)
266  CHKERR MatCreateSubMatrix(A, is2, is2, MAT_INITIAL_MATRIX, subA);
267  #else
268  CHKERR MatGetSubMatrix(A, is2, is2, MAT_INITIAL_MATRIX, subA);
269  #endif
270  }
271  }
272  if (dm_field->aO) {
273  CHKERR ISDestroy(&is2);
274  }
275  dm_field->kspOperators.push_back(*subA);
276  CHKERR PetscObjectReference((PetscObject)(*subA));
277  } else {
278  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
279  }
280  PetscInfo(dm, "Push back IS to DMMGViaApproxOrders\n");
282 }
MoFEMErrorCode sub_mat_mult_add(Mat a, Vec x, Vec f)
FTensor::Index< 'n', 2 > n
Definition: PlasticOps.hpp:68
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
MoFEMErrorCode sub_mat_sor(Mat mat, Vec b, PetscReal omega, MatSORType flag, PetscReal shift, PetscInt its, PetscInt lits, Vec x)
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
MoFEMErrorCode sub_mat_mult(Mat a, Vec x, Vec f)
const int N
Definition: speed_test.cpp:3

◆ DMMoFEMAddColCompositeProblem()

PetscErrorCode MoFEM::DMMoFEMAddColCompositeProblem ( DM  dm,
const char  prb_name[] 
)

Add problem to composite DM on colThis create block on col with DOFs from problem of given name.

Parameters
dmthe DM object
prb_nameadd problem name
Returns
error code

Definition at line 307 of file DMMMoFEM.cpp.

307  {
308  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
310  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
311  if (!dm->data) {
312  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
313  "data structure for MoFEM not yet created");
314  }
315  if (!dm_field->isCompDM) {
316  dm_field->isCompDM = PETSC_TRUE;
317  }
318  if (dm_field->isSquareMatrix) {
319  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
320  "No need to add problem on column when problem block structurally "
321  "symmetric");
322  }
323  dm_field->colCompPrb.push_back(prb_name);
325 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMoFEMAddElement()

PetscErrorCode MoFEM::DMMoFEMAddElement ( DM  dm,
const char  fe_name[] 
)

add element to dm

Note
add_file is a collective, should be executed on all processors. Otherwise could lead to deadlock.
Examples
analytical_poisson_field_split.cpp, cell_forces.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, dm_partitioned_no_field.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, EshelbianPlasticity.cpp, lorentz_force.cpp, MagneticElement.hpp, minimal_surface_area.cpp, Remodeling.cpp, simple_contact.cpp, simple_elasticity.cpp, test_jacobian_of_simple_contact_element.cpp, and UnsaturatedFlow.hpp.

Definition at line 425 of file DMMMoFEM.cpp.

425  {
426  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
428  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
429  ierr = dm_field->mField_ptr->modify_problem_add_finite_element(
430  dm_field->problemName, fe_name);
431  CHKERRG(ierr);
433 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMAddRowCompositeProblem()

PetscErrorCode MoFEM::DMMoFEMAddRowCompositeProblem ( DM  dm,
const char  prb_name[] 
)

Add problem to composite DM on rowThis create block on row with DOFs from problem of given name.

Parameters
dmthe DM object
prb_nameadd problem name
Returns
error code

Definition at line 289 of file DMMMoFEM.cpp.

289  {
290  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
292  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
293  if (!dm->data) {
294  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
295  "data structure for MoFEM not yet created");
296  }
297  if (!dm_field->isCompDM) {
298  dm_field->isCompDM = PETSC_TRUE;
299  }
300  dm_field->rowCompPrb.push_back(prb_name);
301  if (dm_field->isSquareMatrix) {
302  dm_field->colCompPrb.push_back(prb_name);
303  }
305 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMoFEMAddSubFieldCol()

PetscErrorCode MoFEM::DMMoFEMAddSubFieldCol ( DM  dm,
const char  field_name[],
EntityType  lo_type = MBVERTEX,
EntityType  hi_type = MBMAXTYPE 
)

Add field to sub dm problem on columns

Examples
analytical_poisson_field_split.cpp, cell_forces.cpp, and dm_create_subdm.cpp.

Definition at line 221 of file DMMMoFEM.cpp.

222  {
223  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
225  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
226  if (!dm->data) {
227  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
228  "data structure for MoFEM not yet created");
229  }
230  if (!dm_field->isSubDM) {
231  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
232  }
233  dm_field->colFields.push_back(field_name);
234  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
235  if (!dm_field->mapTypeCol)
236  dm_field->mapTypeCol = boost::make_shared<
237  std::map<std::string, std::pair<EntityType, EntityType>>>();
238  (*dm_field->mapTypeCol)[field_name] =
239  std::pair<EntityType, EntityType>(lo_type, hi_type);
240  }
242 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMAddSubFieldRow()

PetscErrorCode MoFEM::DMMoFEMAddSubFieldRow ( DM  dm,
const char  field_name[],
EntityType  lo_type = MBVERTEX,
EntityType  hi_type = MBMAXTYPE 
)

Add field to sub dm problem on rows

Examples
analytical_poisson_field_split.cpp, cell_forces.cpp, dm_create_subdm.cpp, and EshelbianPlasticity.cpp.

Definition at line 198 of file DMMMoFEM.cpp.

199  {
200  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
202  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
203  if (!dm->data) {
204  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
205  "data structure for MoFEM not yet created");
206  }
207  if (!dm_field->isSubDM) {
208  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
209  }
210  dm_field->rowFields.push_back(field_name);
211  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
212  if (!dm_field->mapTypeRow)
213  dm_field->mapTypeRow = boost::make_shared<
214  std::map<std::string, std::pair<EntityType, EntityType>>>();
215  (*dm_field->mapTypeRow)[field_name] =
216  std::pair<EntityType, EntityType>(lo_type, hi_type);
217  }
219 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMCreateMoFEM()

PetscErrorCode MoFEM::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.

Examples
cell_forces.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, dm_partitioned_no_field.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, EshelbianPlasticity.cpp, lorentz_force.cpp, MagneticElement.hpp, minimal_surface_area.cpp, Remodeling.cpp, simple_contact.cpp, test_jacobian_of_simple_contact_element.cpp, and UnsaturatedFlow.hpp.

Definition at line 105 of file DMMMoFEM.cpp.

108  {
110 
111  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
112  if (!dm->data) {
113  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
114  "data structure for MoFEM not yet created");
115  }
116  if (!m_field_ptr) {
117  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
118  "DM function not implemented into MoFEM");
119  }
120  dm_field->mField_ptr = m_field_ptr;
121  dm_field->problemName = problem_name;
122  if (!m_field_ptr->check_problem(dm_field->problemName)) {
123  // problem is not defined, declare problem internally set bool to
124  // destroyProblem problem with DM
125  dm_field->destroyProblem = PETSC_TRUE;
126  CHKERR dm_field->mField_ptr->add_problem(dm_field->problemName, MF_EXCL,
127  dm_field->verbosity);
128  } else {
129  dm_field->destroyProblem = PETSC_FALSE;
130  }
131  CHKERR dm_field->mField_ptr->modify_problem_ref_level_add_bit(
132  dm_field->problemName, bit_level);
133  CHKERR dm_field->mField_ptr->modify_problem_mask_ref_level_add_bit(
134  dm_field->problemName, bit_mask);
135  dm_field->kspCtx =
136  boost::shared_ptr<KspCtx>(new KspCtx(*m_field_ptr, problem_name));
137  dm_field->snesCtx =
138  boost::shared_ptr<SnesCtx>(new SnesCtx(*m_field_ptr, problem_name));
139  dm_field->tsCtx =
140  boost::shared_ptr<TsCtx>(new TsCtx(*m_field_ptr, problem_name));
141 
142  MPI_Comm comm;
143  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
144  int result = 0;
145  MPI_Comm_compare(comm, m_field_ptr->get_comm(), &result);
146  if (result > MPI_CONGRUENT) {
147  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
148  "MoFEM and DM using different communicators");
149  }
150  MPI_Comm_size(comm, &dm_field->sIze);
151  MPI_Comm_rank(comm, &dm_field->rAnk);
152 
153  // problem structure
154  CHKERR dm_field->mField_ptr->get_problem(dm_field->problemName,
155  &dm_field->problemPtr);
156 
158 }
virtual bool check_problem(const std::string name)=0
check if problem exist
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
virtual MPI_Comm & get_comm() const =0

◆ DMMoFEMCreateSubDM()

PetscErrorCode MoFEM::DMMoFEMCreateSubDM ( DM  subdm,
DM  dm,
const char  problem_name[] 
)

Must be called by user to set Sub DM MoFEM data structures.

Examples
analytical_poisson_field_split.cpp, cell_forces.cpp, dm_create_subdm.cpp, and EshelbianPlasticity.cpp.

Definition at line 177 of file DMMMoFEM.cpp.

177  {
179 
180  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
181  if (!dm->data) {
182  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
183  "data structure for MoFEM not yet created");
184  }
185  CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
186  dm_field->problemPtr->getBitRefLevel());
187 
188  DMCtx *subdm_field = (DMCtx *)subdm->data;
189  subdm_field->isSubDM = PETSC_TRUE;
190  subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
191  subdm_field->isPartitioned = dm_field->isPartitioned;
192  subdm_field->isSquareMatrix = PETSC_FALSE;
193  subdm->ops->setup = DMSubDMSetUp_MoFEM;
194 
196 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1166
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:105
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMoFEMGetFieldIS()

PetscErrorCode MoFEM::DMMoFEMGetFieldIS ( DM  dm,
RowColData  rc,
const char  field_name[],
IS *  is 
)

get field is in the problem

Parameters
dmthe DM objects
rcROW or COL (no difference is problem is squared)
field_namename of the field
isreturned the IS object
Returns
error code
IS is;
ierr = DMMoFEMGetFieldIS(dm,ROW,"DISP",&is_disp); CHKERRG(ierr);

Definition at line 1330 of file DMMMoFEM.cpp.

1331  {
1332  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1334  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1335  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1336  ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1337  field_name, 0, 1000, is);
1339 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMoFEMGetIsCompDM()

PetscErrorCode MoFEM::DMMoFEMGetIsCompDM ( DM  dm,
PetscBool *  is_comp_dm 
)

Get if this DM is composite DM.

Parameters
dmthe DM object
is_comp_dmreturn true if composite problem here
Returns
error code

Definition at line 327 of file DMMMoFEM.cpp.

327  {
329  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
331  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
332  *is_comp_dm = dm_field->isCompDM;
334 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetIsPartitioned()

PetscErrorCode MoFEM::DMMoFEMGetIsPartitioned ( DM  dm,
PetscBool *  is_partitioned 
)

get if read mesh is partitioned

Definition at line 993 of file DMMMoFEM.cpp.

993  {
994  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
996  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
997  *is_partitioned = dm_field->isPartitioned;
999 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetIsSubDM()

PetscErrorCode MoFEM::DMMoFEMGetIsSubDM ( DM  dm,
PetscBool *  is_sub_dm 
)

Return true if this DM is sub problem

Parameters
dmthe DM object
is_subproblemtrue if subproblem
Returns
error code

Definition at line 244 of file DMMMoFEM.cpp.

244  {
246  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
248  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
249  *is_sub_dm = dm_field->isSubDM;
251 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetKspCtx() [1/2]

PetscErrorCode MoFEM::DMMoFEMGetKspCtx ( DM  dm,
MoFEM::KspCtx **  ksp_ctx 
)

get MoFEM::KspCtx data structure

Definition at line 927 of file DMMMoFEM.cpp.

927  {
928  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
930  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
931  *ksp_ctx = dm_field->kspCtx.get();
933 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetKspCtx() [2/2]

PetscErrorCode MoFEM::DMMoFEMGetKspCtx ( DM  dm,
const boost::shared_ptr< MoFEM::KspCtx > &  ksp_ctx 
)

get MoFEM::KspCtx data structure

Definition at line 936 of file DMMMoFEM.cpp.

936  {
937  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
939  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
940  const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
942 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetProblemFiniteElementLayout()

PetscErrorCode MoFEM::DMMoFEMGetProblemFiniteElementLayout ( DM  dm,
const char  fe_name[],
PetscLayout *  layout 
)

Get finite elements layout in the problem.

In layout is stored information how many elements is on each processor, for more information look int petsc documentation http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/IS/PetscLayoutCreate.html#PetscLayoutCreate

Parameters
dmdiscrete manager for this problem
fe_namefinite element name
layoutpointer to layout, for created layout user takes responsibility for destroying it.
Returns
error code

Definition at line 401 of file DMMMoFEM.cpp.

402  {
404  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
406  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
407 
408  MPI_Comm comm;
409  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
410  CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
411  layout);
413 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMoFEMGetProblemPtr()

PetscErrorCode MoFEM::DMMoFEMGetProblemPtr ( DM  dm,
const MoFEM::Problem **  problem_ptr 
)

Get pointer to problem data structure.

Examples
EshelbianPlasticity.cpp, field_evaluator.cpp, loop_entities.cpp, lorentz_force.cpp, MagneticElement.hpp, and Remodeling.cpp.

Definition at line 348 of file DMMMoFEM.cpp.

348  {
349  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
351  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
352  if (!dm->data) {
353  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
354  "data structure for MoFEM not yet created");
355  }
356  *problem_ptr = dm_field->problemPtr;
358 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetSnesCtx() [1/2]

PetscErrorCode MoFEM::DMMoFEMGetSnesCtx ( DM  dm,
MoFEM::SnesCtx **  snes_ctx 
)

get MoFEM::SnesCtx data structure

Examples
minimal_surface_area.cpp, simple_contact.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 953 of file DMMMoFEM.cpp.

953  {
954  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
956  DMCtx *dm_field = (DMCtx *)dm->data;
957  *snes_ctx = dm_field->snesCtx.get();
959 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetSnesCtx() [2/2]

PetscErrorCode MoFEM::DMMoFEMGetSnesCtx ( DM  dm,
const boost::shared_ptr< MoFEM::SnesCtx > &  snes_ctx 
)

get MoFEM::SnesCtx data structure

Definition at line 962 of file DMMMoFEM.cpp.

962  {
963  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
965  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
966  const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
968 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetSquareProblem()

PetscErrorCode MoFEM::DMMoFEMGetSquareProblem ( DM  dm,
PetscBool *  square_problem 
)

get squared problemIt if true is assumed that matrix has the same indexing on rows and columns. This reduces interprocessor communication.

Definition at line 415 of file DMMMoFEM.cpp.

415  {
418  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
420  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
421  *square_problem = dm_field->isSquareMatrix;
423 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetTsCtx() [1/2]

PetscErrorCode MoFEM::DMMoFEMGetTsCtx ( DM  dm,
MoFEM::TsCtx **  ts_ctx 
)

get MoFEM::TsCtx data structure

Examples
EshelbianPlasticity.cpp, Remodeling.cpp, and UnsaturatedFlow.hpp.

Definition at line 1001 of file DMMMoFEM.cpp.

1001  {
1002  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1004  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1005  *ts_ctx = dm_field->tsCtx.get();
1007 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMGetTsCtx() [2/2]

PetscErrorCode MoFEM::DMMoFEMGetTsCtx ( DM  dm,
const boost::shared_ptr< MoFEM::TsCtx > &  ts_ctx 
)

get MoFEM::TsCtx data structure

Definition at line 1009 of file DMMMoFEM.cpp.

1010  {
1011  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1013  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1014  const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1016 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMKSPSetComputeOperators() [1/2]

PetscErrorCode MoFEM::DMMoFEMKSPSetComputeOperators ( DM  dm,
const char  fe_name[],
MoFEM::FEMethod method,
MoFEM::BasicMethod pre_only,
MoFEM::BasicMethod post_only 
)

Set KSP operators and push mofem finite element methods.

Parameters
dmDM
fe_namefinite element name
methodmethod on the element (executed for each element in the problem which given name)
pre_onlymethod for pre-process before element method
post_onlymethod for post-process after element method
Returns
error code
Examples
analytical_poisson.cpp, elasticity_mixed_formulation.cpp, and simple_elasticity.cpp.

Definition at line 597 of file DMMMoFEM.cpp.

600  {
601  return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
604  dm, fe_name, method, pre_only, post_only);
605 }
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:578
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...

◆ DMMoFEMKSPSetComputeOperators() [2/2]

PetscErrorCode MoFEM::DMMoFEMKSPSetComputeOperators ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
boost::shared_ptr< MoFEM::BasicMethod pre_only,
boost::shared_ptr< MoFEM::BasicMethod post_only 
)

Set KSP operators and push mofem finite element methods.

Parameters
dmDM
fe_namefinite element name
methodmethod on the element (executed for each element in the problem which given name)
pre_onlymethod for pre-process before element method
post_onlymethod for post-process after element method
Returns
error code

Definition at line 608 of file DMMMoFEM.cpp.

611  {
612  return DMMoFEMKSPSetComputeOperators<const std::string,
613  boost::shared_ptr<MoFEM::FEMethod>>(
614  dm, fe_name, method, pre_only, post_only);
615 }
static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:578

◆ DMMoFEMKSPSetComputeRHS() [1/2]

PetscErrorCode MoFEM::DMMoFEMKSPSetComputeRHS ( DM  dm,
const char  fe_name[],
MoFEM::FEMethod method,
MoFEM::BasicMethod pre_only,
MoFEM::BasicMethod post_only 
)

Set compute operator for KSP solver via sub-matrix and IS.

Parameters
dmDM
Returns
error codeset KSP right hand side evaluation function
Examples
analytical_poisson.cpp, and simple_elasticity.cpp.

Definition at line 556 of file DMMMoFEM.cpp.

559  {
560  return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
562  dm, fe_name, method, pre_only, post_only);
563 }
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:537
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...

◆ DMMoFEMKSPSetComputeRHS() [2/2]

PetscErrorCode MoFEM::DMMoFEMKSPSetComputeRHS ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
boost::shared_ptr< MoFEM::BasicMethod pre_only,
boost::shared_ptr< MoFEM::BasicMethod post_only 
)

set KSP right hand side evaluation function

Definition at line 566 of file DMMMoFEM.cpp.

569  {
570  return DMMoFEMKSPSetComputeRHS<const std::string,
571  boost::shared_ptr<MoFEM::FEMethod>,
572  boost::shared_ptr<MoFEM::BasicMethod>,
573  boost::shared_ptr<MoFEM::BasicMethod>>(
574  dm, fe_name, method, pre_only, post_only);
575 }
static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:537

◆ DMMoFEMResolveSharedFiniteElements()

PetscErrorCode MoFEM::DMMoFEMResolveSharedFiniteElements ( DM  dm,
const char  fe_name[] 
)

Resolve shared entities.

Parameters
dmdm
fe_namefinite element for which shared entities are resolved
Returns
error code
Note
This function is valid for parallel algebra and serial mesh. It should be run collectively, i.e. on all processors.

This allows for tag reduction or tag exchange, f.e.

Tag th;
CHKERR mField.get_moab().tag_get_handle("ADAPT_ORDER",th);
ParallelComm* pcomm =
ParallelComm::get_pcomm(&mField.get_moab(),MYPCOMM_INDEX);
// CHKERR pcomm->reduce_tags(th,MPI_SUM,prisms);
CHKERR pcomm->exchange_tags(th,prisms);

Definition at line 387 of file DMMMoFEM.cpp.

387  {
389  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
391  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
392  CHKERR dm_field->mField_ptr->getInterface<CommInterface>()
393  ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
395 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMoFEMSetIsPartitioned()

PetscErrorCode MoFEM::DMMoFEMSetIsPartitioned ( DM  dm,
PetscBool  is_partitioned 
)

sets if read mesh is partitioned

get if read mesh is partitioned

Examples
cell_forces.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, EshelbianPlasticity.cpp, Remodeling.cpp, simple_contact.cpp, simple_elasticity.cpp, test_jacobian_of_simple_contact_element.cpp, and UnsaturatedFlow.hpp.

Definition at line 982 of file DMMMoFEM.cpp.

982  {
983  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
985  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
986  dm_field->isPartitioned = is_partitioned;
988 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMSetKspCtx()

PetscErrorCode MoFEM::DMMoFEMSetKspCtx ( DM  dm,
boost::shared_ptr< MoFEM::KspCtx > &  ksp_ctx 
)

set MoFEM::KspCtx data structure

Definition at line 944 of file DMMMoFEM.cpp.

945  {
946  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
948  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
949  dm_field->kspCtx = ksp_ctx;
951 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMSetSnesCtx()

PetscErrorCode MoFEM::DMMoFEMSetSnesCtx ( DM  dm,
boost::shared_ptr< MoFEM::SnesCtx > &  snes_ctx 
)

Set MoFEM::SnesCtx data structure.

Definition at line 970 of file DMMMoFEM.cpp.

971  {
972  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
974  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
975  dm_field->snesCtx = snes_ctx;
977 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMSetSquareProblem()

PetscErrorCode MoFEM::DMMoFEMSetSquareProblem ( DM  dm,
PetscBool  square_problem 
)

set squared problemIt if true is assumed that matrix has the same indexing on rows and columns. This reduces interprocessor communication.

Examples
analytical_poisson_field_split.cpp, cell_forces.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, dm_partitioned_no_field.cpp, and EshelbianPlasticity.cpp.

Definition at line 378 of file DMMMoFEM.cpp.

378  {
380  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
382  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
383  dm_field->isSquareMatrix = square_problem;
385 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMSetTsCtx()

PetscErrorCode MoFEM::DMMoFEMSetTsCtx ( DM  dm,
boost::shared_ptr< MoFEM::TsCtx > &  ts_ctx 
)

Set MoFEM::TsCtx data structureIt take over pointer, do not delete it, DM will destroy pointer when is destroyed.

Definition at line 1018 of file DMMMoFEM.cpp.

1018  {
1019  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1021  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1022  dm_field->tsCtx = ts_ctx;
1024 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMMoFEMSNESSetFunction() [1/2]

PetscErrorCode MoFEM::DMMoFEMSNESSetFunction ( DM  dm,
const char  fe_name[],
MoFEM::FEMethod method,
MoFEM::BasicMethod pre_only,
MoFEM::BasicMethod post_only 
)

set SNES residual evaluation function

Examples
analytical_nonlinear_poisson.cpp, mesh_smoothing.cpp, minimal_surface_area.cpp, simple_contact.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 637 of file DMMMoFEM.cpp.

640  {
641  return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
643  dm, fe_name, method, pre_only, post_only);
644 }
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:618
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...

◆ DMMoFEMSNESSetFunction() [2/2]

PetscErrorCode MoFEM::DMMoFEMSNESSetFunction ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
boost::shared_ptr< MoFEM::BasicMethod pre_only,
boost::shared_ptr< MoFEM::BasicMethod post_only 
)

set SNES residual evaluation function

Definition at line 647 of file DMMMoFEM.cpp.

650  {
651  return DMMoFEMSNESSetFunction<const std::string,
652  boost::shared_ptr<MoFEM::FEMethod>,
653  boost::shared_ptr<MoFEM::BasicMethod>,
654  boost::shared_ptr<MoFEM::BasicMethod>>(
655  dm, fe_name, method, pre_only, post_only);
656 }
static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:618

◆ DMMoFEMSNESSetJacobian() [1/2]

PetscErrorCode MoFEM::DMMoFEMSNESSetJacobian ( DM  dm,
const char  fe_name[],
MoFEM::FEMethod method,
MoFEM::BasicMethod pre_only,
MoFEM::BasicMethod post_only 
)

set SNES Jacobian evaluation function

Examples
analytical_nonlinear_poisson.cpp, mesh_smoothing.cpp, minimal_surface_area.cpp, simple_contact.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 678 of file DMMMoFEM.cpp.

681  {
682  return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
684  dm, fe_name, method, pre_only, post_only);
685 }
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:659
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...

◆ DMMoFEMSNESSetJacobian() [2/2]

PetscErrorCode MoFEM::DMMoFEMSNESSetJacobian ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
boost::shared_ptr< MoFEM::BasicMethod pre_only,
boost::shared_ptr< MoFEM::BasicMethod post_only 
)

set SNES Jacobian evaluation function

Definition at line 688 of file DMMMoFEM.cpp.

691  {
692  return DMMoFEMSNESSetJacobian<const std::string,
693  boost::shared_ptr<MoFEM::FEMethod>,
694  boost::shared_ptr<MoFEM::BasicMethod>,
695  boost::shared_ptr<MoFEM::BasicMethod>>(
696  dm, fe_name, method, pre_only, post_only);
697 }
static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:659

◆ DMMoFEMTSSetI2Function()

PetscErrorCode MoFEM::DMMoFEMTSSetI2Function ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
boost::shared_ptr< MoFEM::BasicMethod pre_only,
boost::shared_ptr< MoFEM::BasicMethod post_only 
)

set TS implicit function evaluation function

Definition at line 844 of file DMMMoFEM.cpp.

847  {
848  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
850  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
851  if (pre_only) {
852  dm_field->tsCtx->get_preProcess_to_do_IFunction().push_back(pre_only);
853  }
854  if (method) {
855  dm_field->tsCtx->get_loops_to_do_IFunction().push_back(
856  PairNameFEMethodPtr(fe_name, method));
857  }
858  if (post_only) {
859  dm_field->tsCtx->get_postProcess_to_do_IFunction().push_back(post_only);
860  }
861  CHKERR DMTSSetI2Function(dm, TsSetI2Function, dm_field->tsCtx.get());
863 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, Vec F, void *ctx)
Calculation the right hand side for second order PDE in time.
Definition: TsCtx.cpp:525
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMoFEMTSSetI2Jacobian()

PetscErrorCode MoFEM::DMMoFEMTSSetI2Jacobian ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
boost::shared_ptr< MoFEM::BasicMethod pre_only,
boost::shared_ptr< MoFEM::BasicMethod post_only 
)

set TS Jacobian evaluation function

Definition at line 866 of file DMMMoFEM.cpp.

869  {
870  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
872  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
873  if (pre_only) {
874  dm_field->tsCtx->get_preProcess_to_do_IJacobian().push_back(pre_only);
875  }
876  if (method) {
877  dm_field->tsCtx->get_loops_to_do_IJacobian().push_back(
878  PairNameFEMethodPtr(fe_name, method));
879  }
880  if (post_only) {
881  dm_field->tsCtx->get_postProcess_to_do_IJacobian().push_back(post_only);
882  }
883  CHKERR DMTSSetI2Jacobian(dm, TsSetI2Jacobian, dm_field->tsCtx.get());
885 }
PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, PetscReal v, PetscReal a, Mat A, Mat B, void *ctx)
Calculation Jaconian for second order PDE in time.
Definition: TsCtx.cpp:434
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMMoFEMTSSetIFunction() [1/2]

PetscErrorCode MoFEM::DMMoFEMTSSetIFunction ( DM  dm,
const char  fe_name[],
MoFEM::FEMethod method,
MoFEM::BasicMethod pre_only,
MoFEM::BasicMethod post_only 
)

set TS implicit function evaluation function

Examples
bone_adaptation.cpp, EshelbianPlasticity.cpp, reaction_diffusion_equation.cpp, and UnsaturatedFlow.hpp.

Definition at line 719 of file DMMMoFEM.cpp.

722  {
723  return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
725  dm, fe_name, method, pre_only, post_only);
727 }
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:700
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...

◆ DMMoFEMTSSetIFunction() [2/2]

PetscErrorCode MoFEM::DMMoFEMTSSetIFunction ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
boost::shared_ptr< MoFEM::BasicMethod pre_only,
boost::shared_ptr< MoFEM::BasicMethod post_only 
)

set TS implicit function evaluation function

Definition at line 730 of file DMMMoFEM.cpp.

733  {
734  return DMMoFEMTSSetIFunction<const std::string,
735  boost::shared_ptr<MoFEM::FEMethod>,
736  boost::shared_ptr<MoFEM::BasicMethod>,
737  boost::shared_ptr<MoFEM::BasicMethod>>(
738  dm, fe_name, method, pre_only, post_only);
740 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:700

◆ DMMoFEMTSSetIJacobian() [1/2]

PetscErrorCode MoFEM::DMMoFEMTSSetIJacobian ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
boost::shared_ptr< MoFEM::BasicMethod pre_only,
boost::shared_ptr< MoFEM::BasicMethod post_only 
)

set TS Jacobian evaluation function

Examples
bone_adaptation.cpp, EshelbianPlasticity.cpp, reaction_diffusion_equation.cpp, and UnsaturatedFlow.hpp.

Definition at line 772 of file DMMMoFEM.cpp.

775  {
776  return DMMoFEMTSSetIJacobian<const std::string,
777  boost::shared_ptr<MoFEM::FEMethod>,
778  boost::shared_ptr<MoFEM::BasicMethod>,
779  boost::shared_ptr<MoFEM::BasicMethod>>(
780  dm, fe_name, method, pre_only, post_only);
781 }
static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:743

◆ DMMoFEMTSSetIJacobian() [2/2]

PetscErrorCode MoFEM::DMMoFEMTSSetIJacobian ( DM  dm,
const char  fe_name[],
MoFEM::FEMethod method,
MoFEM::BasicMethod pre_only,
MoFEM::BasicMethod post_only 
)

set TS Jacobian evaluation function

Definition at line 762 of file DMMMoFEM.cpp.

765  {
766  return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
767  MoFEM::BasicMethod *>(dm, fe_name, method,
768  pre_only, post_only);
769 }
static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:743
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...

◆ DMMoFEMUnSetElement()

PetscErrorCode MoFEM::DMMoFEMUnSetElement ( DM  dm,
const char  fe_name[] 
)

unset element from dm

Definition at line 435 of file DMMMoFEM.cpp.

435  {
436  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
438  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
439  ierr = dm_field->mField_ptr->modify_problem_unset_finite_element(
440  dm_field->problemName, fe_name);
441  CHKERRG(ierr);
443 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMoFEMGetInterfacePtr()

PetscErrorCode MoFEM::DMoFEMGetInterfacePtr ( DM  dm,
MoFEM::Interface **  m_field_ptr 
)

Get pointer to MoFEM::Interface.

Parameters
dmDistributed mesh manager
m_field_ptrPointer to pointer of field interface
Returns
Error code
Examples
HookeElement.cpp, and mesh_smoothing.cpp.

Definition at line 336 of file DMMMoFEM.cpp.

336  {
337  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
339  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
340  if (!dm->data) {
341  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
342  "data structure for MoFEM not yet created");
343  }
344  *m_field_ptr = dm_field->mField_ptr;
346 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMoFEMLoopDofs()

PetscErrorCode MoFEM::DMoFEMLoopDofs ( DM  dm,
const char  field_name[],
MoFEM::DofMethod method 
)

execute method for dofs on field in problem

Definition at line 524 of file DMMMoFEM.cpp.

525  {
526  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
528  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
529  ierr =
530  dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
531  *method, dm_field->rAnk, dm_field->rAnk);
532  CHKERRG(ierr);
534 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMoFEMLoopFiniteElements() [1/2]

PetscErrorCode MoFEM::DMoFEMLoopFiniteElements ( DM  dm,
const char  fe_name[],
MoFEM::FEMethod method 
)

Executes FEMethod for finite elements in DM.

Parameters
dmMoFEM discrete manager
fe_namename of element
methodpointer to MOFEM::FEMethod
Returns
Error code
Examples
analytical_nonlinear_poisson.cpp, analytical_poisson.cpp, analytical_poisson_field_split.cpp, boundary_marker.cpp, cell_forces.cpp, continuity_check_on_skeleton_with_simple_2d.cpp, dm_build_partitioned_mesh.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, EshelbianPlasticity.cpp, hello_world.cpp, HookeElement.cpp, MagneticElement.hpp, mesh_smoothing.cpp, minimal_surface_area.cpp, reaction_diffusion_equation.cpp, Remodeling.cpp, simple_contact.cpp, simple_elasticity.cpp, simple_interface.cpp, and UnsaturatedFlow.hpp.

Definition at line 507 of file DMMMoFEM.cpp.

508  {
509  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
511  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
512  ierr = DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name, method,
513  dm_field->rAnk, dm_field->rAnk);
514  CHKERRG(ierr);
516 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:488

◆ DMoFEMLoopFiniteElements() [2/2]

PetscErrorCode MoFEM::DMoFEMLoopFiniteElements ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method 
)

Executes FEMethod for finite elements in DM.

Parameters
dmMoFEM discrete manager
fe_namename of element
methodpointer to MOFEM::FEMethod
Returns
Error code

Definition at line 519 of file DMMMoFEM.cpp.

520  {
521  return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get());
522 }
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:507

◆ DMoFEMLoopFiniteElementsUpAndLowRank() [1/2]

PetscErrorCode MoFEM::DMoFEMLoopFiniteElementsUpAndLowRank ( DM  dm,
const char  fe_name[],
MoFEM::FEMethod method,
int  low_rank,
int  up_rank 
)

Executes FEMethod for finite elements in DM.

Parameters
dmMoFEM discrete manager
fe_namename of finite element
methodpointer to MoFEM::FEMethod
low_ranklowest rank of processor
up_rankupper run of processor
Returns
Error code

Definition at line 488 of file DMMMoFEM.cpp.

490  {
492  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
493  ierr = dm_field->mField_ptr->loop_finite_elements(
494  dm_field->problemPtr, fe_name, *method, low_rank, up_rank);
495  CHKERRG(ierr);
497 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMoFEMLoopFiniteElementsUpAndLowRank() [2/2]

PetscErrorCode MoFEM::DMoFEMLoopFiniteElementsUpAndLowRank ( DM  dm,
const std::string  fe_name,
boost::shared_ptr< MoFEM::FEMethod method,
int  low_rank,
int  up_rank 
)

Executes FEMethod for finite elements in DM.

Parameters
dmMoFEM discrete manager
fe_namename of finite element
methodpointer to MoFEM::FEMethod
low_ranklowest rank of processor
up_rankupper run of processor
Returns
Error code

Definition at line 500 of file DMMMoFEM.cpp.

502  {
503  return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
504  low_rank, up_rank);
505 }
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:488

◆ DMoFEMMeshToGlobalVector()

PetscErrorCode MoFEM::DMoFEMMeshToGlobalVector ( DM  dm,
Vec  g,
InsertMode  mode,
ScatterMode  scatter_mode 
)

set ghosted vector values on all existing mesh entities

Parameters
gvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities from V vector.

SCATTER_FORWARD set vector V from data field entities

Examples
analytical_nonlinear_poisson.cpp, analytical_poisson.cpp, analytical_poisson_field_split.cpp, elasticity_mixed_formulation.cpp, mesh_smoothing.cpp, minimal_surface_area.cpp, simple_contact.cpp, and simple_elasticity.cpp.

Definition at line 457 of file DMMMoFEM.cpp.

458  {
459  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
461  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
462  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
463  dm_field->problemPtr, COL, g, mode, scatter_mode);
464  CHKERRG(ierr);
466 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMoFEMMeshToLocalVector()

PetscErrorCode MoFEM::DMoFEMMeshToLocalVector ( DM  dm,
Vec  l,
InsertMode  mode,
ScatterMode  scatter_mode 
)

set local (or ghosted) vector values on mesh for partition only

Parameters
lvector
modesee petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
scatter_modesee petsc manual for ScatterMode (The available modes are: SCATTER_FORWARD or SCATTER_REVERSE)

SCATTER_REVERSE set data to field entities from V vector.

SCATTER_FORWARD set vector V from data field entities

Examples
bone_adaptation.cpp, cell_forces.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, ep.cpp, EshelbianPlasticity.cpp, hcurl_check_approx_in_2d.cpp, lesson2_approximaton.cpp, lesson3_poisson.cpp, lesson4_elastic.cpp, lesson5_helmholtz.cpp, lesson6_radiation.cpp, lesson7_plastic.cpp, lesson8_contact.cpp, MagneticElement.hpp, mesh_smoothing.cpp, minimal_surface_area.cpp, reaction_diffusion_equation.cpp, scalar_check_approximation_2d.cpp, simple_contact.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and UnsaturatedFlow.hpp.

Definition at line 445 of file DMMMoFEM.cpp.

446  {
448  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
450  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
451  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
452  dm_field->problemPtr, COL, l, mode, scatter_mode);
453  CHKERRG(ierr);
455 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
FTensor::Index< 'l', 2 > l
Definition: ContactOps.hpp:29

◆ DMoFEMPostProcessFiniteElements()

PetscErrorCode MoFEM::DMoFEMPostProcessFiniteElements ( DM  dm,
MoFEM::FEMethod method 
)

execute finite element method for each element in dm (problem)

Examples
cell_forces.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, and simple_elasticity.cpp.

Definition at line 478 of file DMMMoFEM.cpp.

478  {
479  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
481  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
482  ierr = dm_field->mField_ptr->problem_basic_method_postProcess(
483  dm_field->problemPtr, *method);
484  CHKERRG(ierr);
486 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMoFEMPreProcessFiniteElements()

PetscErrorCode MoFEM::DMoFEMPreProcessFiniteElements ( DM  dm,
MoFEM::FEMethod method 
)

execute finite element method for each element in dm (problem)

Examples
cell_forces.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, and simple_contact.cpp.

Definition at line 468 of file DMMMoFEM.cpp.

468  {
469  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
471  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
472  ierr = dm_field->mField_ptr->problem_basic_method_preProcess(
473  dm_field->problemPtr, *method);
474  CHKERRG(ierr);
476 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMRegister_MGViaApproxOrders()

MoFEMErrorCode DMRegister_MGViaApproxOrders ( const char  sname[])

Register DM for Multi-Grid via approximation orders.

Parameters
snameproblem/dm registered name
Returns
error code

Definition at line 400 of file PCMGSetUpViaApproxOrders.cpp.

400  {
402  CHKERR DMRegister(sname, DMCreate_MGViaApproxOrders);
404 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
MoFEMErrorCode DMCreate_MGViaApproxOrders(DM dm)
Create DM data structure for Multi-Grid via approximation orders.
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMRegister_MoFEM()

PetscErrorCode MoFEM::DMRegister_MoFEM ( const char  sname[])

Register MoFEM problem.

Examples
analytical_nonlinear_poisson.cpp, analytical_poisson.cpp, analytical_poisson_field_split.cpp, boundary_marker.cpp, cell_forces.cpp, continuity_check_on_skeleton_with_simple_2d.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, dm_partitioned_no_field.cpp, elasticity_mixed_formulation.cpp, ep.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, hello_world.cpp, lesson1_moment_of_inertia.cpp, lesson2_approximaton.cpp, lesson3_poisson.cpp, lesson4_elastic.cpp, lesson5_helmholtz.cpp, lesson6_radiation.cpp, lesson7_plastic.cpp, lesson8_contact.cpp, loop_entities.cpp, lorentz_force.cpp, MagneticElement.hpp, mesh_smoothing.cpp, minimal_surface_area.cpp, reaction_diffusion_equation.cpp, Remodeling.cpp, scalar_check_approximation_2d.cpp, simple_contact.cpp, simple_elasticity.cpp, simple_interface.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and unsaturated_transport.cpp.

Definition at line 48 of file DMMMoFEM.cpp.

48  {
50  CHKERR DMRegister(sname, DMCreate_MoFEM);
52 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:76
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMSetFromOptions_MoFEM()

PetscErrorCode MoFEM::DMSetFromOptions_MoFEM ( DM  dm)

Set options for MoFEM DM

Definition at line 1097 of file DMMMoFEM.cpp.

1097  {
1098 #endif
1099 
1100  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1102  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1103 #if PETSC_VERSION_GE(3, 5, 3)
1104  ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1105  CHKERRG(ierr);
1106 #else
1107  ierr = PetscOptionsHead("DMMoFEM Options");
1108  CHKERRG(ierr);
1109 #endif
1110  ierr = PetscOptionsBool("-dm_is_partitioned",
1111  "set if mesh is partitioned (works which native MOAB "
1112  "file format, i.e. h5m",
1113  "DMSetUp", dm_field->isPartitioned,
1114  &dm_field->isPartitioned, NULL);
1115  CHKERRG(ierr);
1117 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514

◆ DMSetOperators_MoFEM()

PetscErrorCode MoFEM::DMSetOperators_MoFEM ( DM  dm)

Set operators for MoFEM dm.

Parameters
dm
Returns
error code

Definition at line 54 of file DMMMoFEM.cpp.

54  {
55  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
57 
58  dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
59  dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
60  dm->ops->creatematrix = DMCreateMatrix_MoFEM;
61  dm->ops->setup = DMSetUp_MoFEM;
62  dm->ops->destroy = DMDestroy_MoFEM;
63  dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
64  dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
65  dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
66  dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
67  dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
68  dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
69 
70  // Default matrix type
71  CHKERR DMSetMatType(dm, MATMPIAIJ);
72 
74 }
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1213
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1119
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1250
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1053
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:84
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM.
Definition: DMMMoFEM.cpp:1026
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1221
#define CHKERR
Inline error check.
Definition: definitions.h:602
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1285
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1279
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1097
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM.
Definition: DMMMoFEM.cpp:1044

◆ DMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSetUp_MoFEM ( DM  dm)

sets up the MoFEM structures inside a DM object

Definition at line 1119 of file DMMMoFEM.cpp.

1119  {
1120  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1121  ProblemsManager *prb_mng_ptr;
1123  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1124  CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1125 
1126  if (dm_field->isCompDM) {
1127  // It is composite probelm
1128  CHKERR prb_mng_ptr->buildCompsedProblem(
1129  dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1130  dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1131  } else {
1132  if (dm_field->isPartitioned) {
1133  CHKERR prb_mng_ptr->buildProblemOnDistributedMesh(
1134  dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1135  dm_field->verbosity);
1136  } else {
1137  CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1138  dm_field->isSquareMatrix == PETSC_TRUE,
1139  dm_field->verbosity);
1140  CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1141  dm_field->verbosity);
1142  }
1143  }
1144 
1145  // Partition finite elements
1146  if (dm_field->isPartitioned) {
1147  CHKERR prb_mng_ptr->partitionFiniteElements(
1148  dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1149  CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1150  dm_field->problemName, dm_field->verbosity);
1151  } else {
1152  // partition finite elemnets
1153  CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1154  -1, -1, dm_field->verbosity);
1155  // Get ghost DOFs
1156  CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1157  dm_field->verbosity);
1158  }
1159 
1160  // Set flag that problem is build and partitioned
1161  dm_field->isProblemBuild = PETSC_TRUE;
1162 
1164 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ DMSubDMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSubDMSetUp_MoFEM ( DM  subdm)

Sets up the MoFEM structures inside a DM object for sub dm

Definition at line 1166 of file DMMMoFEM.cpp.

1166  {
1167  PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1168  ProblemsManager *prb_mng_ptr;
1170 
1171  DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1172 
1173  // build sub dm problem
1174  CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1175 
1176  map<std::string, std::pair<EntityType, EntityType>> *entity_map_row = nullptr;
1177  map<std::string, std::pair<EntityType, EntityType>> *entity_map_col = nullptr;
1178 
1179  if (subdm_field->mapTypeRow)
1180  entity_map_row = subdm_field->mapTypeRow.get();
1181  if (subdm_field->mapTypeCol)
1182  entity_map_row = subdm_field->mapTypeCol.get();
1183 
1184  CHKERR prb_mng_ptr->buildSubProblem(
1185  subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1186  subdm_field->problemMainOfSubPtr->getName(),
1187  subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1188  subdm_field->verbosity);
1189 
1190  // partition problem
1191  subdm_field->isPartitioned = subdm_field->isPartitioned;
1192  if (subdm_field->isPartitioned) {
1193  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1194  0, subdm_field->sIze,
1195  subdm_field->verbosity);
1196  // set ghost nodes
1197  CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1198  subdm_field->problemName, subdm_field->verbosity);
1199  } else {
1200  // partition finite elements
1201  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1202  -1, -1, subdm_field->verbosity);
1203  // set ghost nodes
1204  CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1205  subdm_field->verbosity);
1206  }
1207 
1208  subdm_field->isProblemBuild = PETSC_TRUE;
1209 
1211 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

Variable Documentation

◆ smartCreateDMMatrix

auto MoFEM::smartCreateDMMatrix
Initial value:
= [](DM dm) {
SmartPetscObj<Mat> a;
CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
return a;
}
MPI_Comm getCommFromPetscObject(PetscObject obj)
Get the Comm From Petsc Object object.
Definition: AuxPETSc.hpp:189
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, SmartPetscObj< Mat > &M)
Definition: DMMMoFEM.cpp:1072

Get smart matrix from DM.

Examples
simple_contact.cpp, and test_jacobian_of_simple_contact_element.cpp.

Definition at line 927 of file DMMoFEM.hpp.

◆ smartCreateDMVector

auto MoFEM::smartCreateDMVector
Initial value:
= [](DM dm) {
SmartPetscObj<Vec> v;
CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
return v;
}
MPI_Comm getCommFromPetscObject(PetscObject obj)
Get the Comm From Petsc Object object.
Definition: AuxPETSc.hpp:189
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj< Vec > &g_ptr)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM.
Definition: DMMMoFEM.cpp:1035

Get smart vector from DM.

Examples
hcurl_check_approx_in_2d.cpp, lesson2_approximaton.cpp, lesson3_poisson.cpp, lesson4_elastic.cpp, lesson5_helmholtz.cpp, lesson6_radiation.cpp, lesson7_plastic.cpp, lesson8_contact.cpp, scalar_check_approximation_2d.cpp, simple_contact.cpp, and test_jacobian_of_simple_contact_element.cpp.

Definition at line 939 of file DMMoFEM.hpp.