v0.8.23
Classes | Functions
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::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...
 

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 477 of file PCMGSetUpViaApproxOrders.cpp.

477  {
478  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
480  GET_DM_FIELD(dm);
481  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
482  CHKERR DMCreate(comm, dmc);
483  (*dmc)->data = dm->data;
484  DMType type;
485  CHKERR DMGetType(dm, &type);
486  CHKERR DMSetType(*dmc, type);
487  CHKERR PetscObjectReference((PetscObject)(*dmc));
488  PetscInfo1(dm, "Coarsen DMMGViaApproxOrders leveldown = %d\n", dm->leveldown);
490 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ DMCreate_MoFEM()

PetscErrorCode MoFEM::DMCreate_MoFEM ( DM  dm)

Create dm data structure with MoFEM data structure.

Definition at line 81 of file DMMMoFEM.cpp.

81  {
82  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
84  dm->data = new DMCtx();
87 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:59
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 1250 of file DMMMoFEM.cpp.

1251  {
1252  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1254 
1255  if (numFields) {
1256  *numFields = 0;
1257  }
1258  if (fieldNames) {
1259  *fieldNames = NULL;
1260  }
1261  if (fields) {
1262  *fields = NULL;
1263  }
1264 
1265  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1266  const Field_multiIndex *fields_ptr;
1267  CHKERR dm_field->mField_ptr->get_fields(&fields_ptr);
1268  Field_multiIndex::iterator fit, hi_fit;
1269  fit = fields_ptr->begin();
1270  hi_fit = fields_ptr->end();
1271  *numFields = std::distance(fit, hi_fit);
1272 
1273  if (fieldNames) {
1274  CHKERR PetscMalloc1(*numFields, fieldNames);
1275  }
1276  if (fields) {
1277  CHKERR PetscMalloc1(*numFields, fields);
1278  }
1279 
1280  for (int f = 0; fit != hi_fit; fit++, f++) {
1281  if (fieldNames) {
1282  CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1283  (char **)&(*fieldNames)[f]);
1284  }
1285  if (fields) {
1286  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1287  ->isCreateProblemFieldAndRank(
1288  dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1289  fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1290  }
1291  }
1292 
1294 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
multi_index_container< boost::shared_ptr< Field >, indexed_by< hashed_unique< tag< BitFieldId_mi_tag >, const_mem_fun< Field, const BitFieldId &, &Field::getId >, HashBit< BitFieldId >, EqBit< BitFieldId > >, ordered_unique< tag< Meshset_mi_tag >, member< Field, EntityHandle, &Field::meshSet > >, ordered_unique< tag< FieldName_mi_tag >, const_mem_fun< Field, boost::string_ref, &Field::getNameRef > >, ordered_non_unique< tag< BitFieldId_space_mi_tag >, const_mem_fun< Field, FieldSpace, &Field::getSpace > > > > Field_multiIndex
Field_multiIndex for Field.
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

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

Definition at line 991 of file DMMMoFEM.cpp.

991  {
992  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
994  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
995  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
996  dm_field->problemName, COL, g);
998 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
#define CHKERR
Inline error check.
Definition: definitions.h:596

◆ 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 1000 of file DMMMoFEM.cpp.

1000  {
1001  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1003  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1004  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1005  dm_field->problemName, COL, g_ptr);
1007 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
#define CHKERR
Inline error check.
Definition: definitions.h:596

◆ 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 1009 of file DMMMoFEM.cpp.

1009  {
1010  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1012  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1013  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1014  dm_field->problemName, COL, l);
1016 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 446 of file PCMGSetUpViaApproxOrders.cpp.

446  {
447 
448  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
450  GET_DM_FIELD(dm);
451 
452  int leveldown = dm->leveldown;
453 
454  if (dm_field->kspOperators.empty()) {
456  } else {
457  MPI_Comm comm;
458  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
459  if (dm_field->kspOperators.empty()) {
460  SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
461  "data inconsistency, operator can not be set");
462  }
463  if (static_cast<int>(dm_field->kspOperators.size()) < leveldown) {
464  SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
465  "data inconsistency, no IS for that level");
466  }
467  *M = dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown];
468  CHKERR PetscObjectReference((PetscObject)*M);
469  }
470 
471  PetscInfo1(dm, "Create Matrix DMMGViaApproxOrders leveldown = %d\n",
472  leveldown);
473 
475 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1018
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

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

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  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1023  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1024  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1025  M);
1026  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1027  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1028  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1029  M);
1030  } else {
1031  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1032  "Matrix type not implemented");
1033  }
1035 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 1037 of file DMMMoFEM.cpp.

1037  {
1038  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1040  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1041  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1042  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1043  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1044  M);
1045  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1046  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1047  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1048  M);
1049  } else {
1050  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1051  "Matrix type not implemented");
1052  }
1054 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ DMDestroy_MoFEM()

PetscErrorCode MoFEM::DMDestroy_MoFEM ( DM  dm)

Destroys dm with MoFEM data structure.

destroy the MoFEM structure

Definition at line 89 of file DMMMoFEM.cpp.

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

◆ 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 1178 of file DMMMoFEM.cpp.

1179  {
1180  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1182  CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1184 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 1186 of file DMMMoFEM.cpp.

1186  {
1188  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1190 
1191  CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1192 
1193  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1194  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1195  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1196 
1197  double *array_loc, *array_glob;
1198  CHKERR VecGetArray(l, &array_loc);
1199  CHKERR VecGetArray(g, &array_glob);
1200  switch (mode) {
1201  case INSERT_VALUES:
1202  cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1203  break;
1204  case ADD_VALUES:
1205  cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1206  break;
1207  default:
1208  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1209  }
1210  CHKERR VecRestoreArray(l, &array_loc);
1211  CHKERR VecRestoreArray(g, &array_glob);
1213 }
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:501
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:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 1215 of file DMMMoFEM.cpp.

1216  {
1217 
1218  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1220 
1221  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1222  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1223  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1224 
1225  double *array_loc, *array_glob;
1226  CHKERR VecGetArray(l, &array_loc);
1227  CHKERR VecGetArray(g, &array_glob);
1228  switch (mode) {
1229  case INSERT_VALUES:
1230  cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1231  break;
1232  case ADD_VALUES:
1233  cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1234  break;
1235  default:
1236  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1237  }
1238  CHKERR VecRestoreArray(l, &array_loc);
1239  CHKERR VecRestoreArray(g, &array_glob);
1240 
1242 }
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:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 1244 of file DMMMoFEM.cpp.

1244  {
1245  //
1248 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMMGViaApproxOrdersClearCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS ( DM  )

Clear approximation orders.

Parameters
DMdm
Returns
Error code

Definition at line 296 of file PCMGSetUpViaApproxOrders.cpp.

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

◆ DMMGViaApproxOrdersPopBackCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS ( DM  )

Pop is form MG via approximation orders.

Parameters
DMdm
ispop back IS
Returns
error code

Definition at line 280 of file PCMGSetUpViaApproxOrders.cpp.

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

◆ 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  CHKERR MatGetSubMatrix(A, is2, is2, MAT_INITIAL_MATRIX, subA);
266  }
267  }
268  if (dm_field->aO) {
269  CHKERR ISDestroy(&is2);
270  }
271  dm_field->kspOperators.push_back(*subA);
272  CHKERR PetscObjectReference((PetscObject)(*subA));
273  } else {
274  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
275  }
276  PetscInfo(dm, "Push back IS to DMMGViaApproxOrders\n");
278 }
MoFEMErrorCode sub_mat_mult_add(Mat a, Vec x, Vec f)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
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:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
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 296 of file DMMMoFEM.cpp.

296  {
297  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
299  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
300  if (!dm->data) {
301  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
302  "data structure for MoFEM not yet created");
303  }
304  if (!dm_field->isCompDM) {
305  dm_field->isCompDM = PETSC_TRUE;
306  }
307  if (dm_field->isSquareMatrix) {
308  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
309  "No need to add problem on column when problem block structurally "
310  "symmetric");
311  }
312  dm_field->colCompPrb.push_back(prb_name);
314 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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_elasticity.cpp, and UnsaturatedFlow.hpp.

Definition at line 414 of file DMMMoFEM.cpp.

414  {
415  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
417  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
418  ierr = dm_field->mField_ptr->modify_problem_add_finite_element(
419  dm_field->problemName, fe_name);
420  CHKERRG(ierr);
422 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 278 of file DMMMoFEM.cpp.

278  {
279  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
281  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
282  if (!dm->data) {
283  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
284  "data structure for MoFEM not yet created");
285  }
286  if (!dm_field->isCompDM) {
287  dm_field->isCompDM = PETSC_TRUE;
288  }
289  dm_field->rowCompPrb.push_back(prb_name);
290  if (dm_field->isSquareMatrix) {
291  dm_field->colCompPrb.push_back(prb_name);
292  }
294 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 210 of file DMMMoFEM.cpp.

211  {
212  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
214  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
215  if (!dm->data) {
216  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
217  "data structure for MoFEM not yet created");
218  }
219  if (!dm_field->isSubDM) {
220  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
221  }
222  dm_field->colFields.push_back(field_name);
223  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
224  if (!dm_field->mapTypeCol)
225  dm_field->mapTypeCol = boost::make_shared<
226  std::map<std::string, std::pair<EntityType, EntityType>>>();
227  (*dm_field->mapTypeCol)[field_name] =
228  std::pair<EntityType, EntityType>(lo_type, hi_type);
229  }
231 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 187 of file DMMMoFEM.cpp.

188  {
189  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
191  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
192  if (!dm->data) {
193  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
194  "data structure for MoFEM not yet created");
195  }
196  if (!dm_field->isSubDM) {
197  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
198  }
199  dm_field->rowFields.push_back(field_name);
200  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
201  if (!dm_field->mapTypeRow)
202  dm_field->mapTypeRow = boost::make_shared<
203  std::map<std::string, std::pair<EntityType, EntityType>>>();
204  (*dm_field->mapTypeRow)[field_name] =
205  std::pair<EntityType, EntityType>(lo_type, hi_type);
206  }
208 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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, and UnsaturatedFlow.hpp.

Definition at line 109 of file DMMMoFEM.cpp.

112  {
114 
115  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
116  if (!dm->data) {
117  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
118  "data structure for MoFEM not yet created");
119  }
120  if (!m_field_ptr) {
121  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
122  "DM function not implemented into MoFEM");
123  }
124  dm_field->mField_ptr = m_field_ptr;
125  dm_field->problemName = problem_name;
126  if (!m_field_ptr->check_problem(dm_field->problemName)) {
127  // problem is not defined, declare problem internally set bool to
128  // destroyProblem problem with DM
129  dm_field->destroyProblem = PETSC_TRUE;
130  CHKERR dm_field->mField_ptr->add_problem(dm_field->problemName, MF_EXCL,
131  dm_field->verbosity);
132  } else {
133  dm_field->destroyProblem = PETSC_FALSE;
134  }
135  CHKERR dm_field->mField_ptr->modify_problem_ref_level_add_bit(
136  dm_field->problemName, bit_level);
137  CHKERR dm_field->mField_ptr->modify_problem_mask_ref_level_add_bit(
138  dm_field->problemName, bit_mask);
139  dm_field->kspCtx =
140  boost::shared_ptr<KspCtx>(new KspCtx(*m_field_ptr, problem_name));
141  dm_field->snesCtx =
142  boost::shared_ptr<SnesCtx>(new SnesCtx(*m_field_ptr, problem_name));
143  dm_field->tsCtx =
144  boost::shared_ptr<TsCtx>(new TsCtx(*m_field_ptr, problem_name));
145 
146  MPI_Comm comm;
147  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
148  int result = 0;
149  MPI_Comm_compare(comm, m_field_ptr->get_comm(), &result);
150  // std::cerr << result << " " << MPI_IDENT << " " << MPI_CONGRUENT << " " <<
151  // MPI_SIMILAR << " " << MPI_UNEQUAL << std::endl;
152  if (result > MPI_CONGRUENT) {
153  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
154  "MoFEM and DM using different communicators");
155  }
156  MPI_Comm_size(comm, &dm_field->sIze);
157  MPI_Comm_rank(comm, &dm_field->rAnk);
158 
159  // problem structure
160  CHKERR dm_field->mField_ptr->get_problem(dm_field->problemName,
161  &dm_field->problemPtr);
162 
164 }
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:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
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 166 of file DMMMoFEM.cpp.

166  {
168 
169  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
170  if (!dm->data) {
171  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
172  "data structure for MoFEM not yet created");
173  }
174  CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
175  dm_field->problemPtr->getBitRefLevel());
176 
177  DMCtx *subdm_field = (DMCtx *)subdm->data;
178  subdm_field->isSubDM = PETSC_TRUE;
179  subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
180  subdm_field->isPartitioned = dm_field->isPartitioned;
181  subdm_field->isSquareMatrix = PETSC_FALSE;
182  subdm->ops->setup = DMSubDMSetUp_MoFEM;
183 
185 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1131
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:109
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 1296 of file DMMMoFEM.cpp.

1297  {
1298  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1300  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1301  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1302  ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1303  field_name, 0, 1000, is);
1305 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 316 of file DMMMoFEM.cpp.

316  {
318  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
320  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
321  *is_comp_dm = dm_field->isCompDM;
323 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMMoFEMGetIsPartitioned()

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

get if read mesh is partitioned

Definition at line 958 of file DMMMoFEM.cpp.

958  {
959  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
961  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
962  *is_partitioned = dm_field->isPartitioned;
964 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

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

233  {
235  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
237  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
238  *is_sub_dm = dm_field->isSubDM;
240 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMMoFEMGetKspCtx() [1/2]

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

get MoFEM::KspCtx data structure

Definition at line 892 of file DMMMoFEM.cpp.

892  {
893  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
895  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
896  *ksp_ctx = dm_field->kspCtx.get();
898 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMMoFEMGetKspCtx() [2/2]

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

get MoFEM::KspCtx data structure

Definition at line 901 of file DMMMoFEM.cpp.

901  {
902  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
904  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
905  const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
907 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 390 of file DMMMoFEM.cpp.

391  {
393  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
395  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
396 
397  MPI_Comm comm;
398  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
399  CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
400  layout);
402 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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 337 of file DMMMoFEM.cpp.

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

◆ DMMoFEMGetSnesCtx() [1/2]

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

get MoFEM::SnesCtx data structure

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

Definition at line 918 of file DMMMoFEM.cpp.

918  {
919  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
921  DMCtx *dm_field = (DMCtx *)dm->data;
922  *snes_ctx = dm_field->snesCtx.get();
924 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMMoFEMGetSnesCtx() [2/2]

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

get MoFEM::SnesCtx 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  const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
933 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 404 of file DMMMoFEM.cpp.

404  {
407  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
409  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
410  *square_problem = dm_field->isSquareMatrix;
412 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 966 of file DMMMoFEM.cpp.

966  {
967  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
969  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
970  *ts_ctx = dm_field->tsCtx.get();
972 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMMoFEMGetTsCtx() [2/2]

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

get MoFEM::TsCtx data structure

Definition at line 974 of file DMMMoFEM.cpp.

975  {
976  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
978  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
979  const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
981 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 586 of file DMMMoFEM.cpp.

589  {
590  return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
593  dm, fe_name, method, pre_only, post_only);
594 }
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:567
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 597 of file DMMMoFEM.cpp.

600  {
601  return DMMoFEMKSPSetComputeOperators<const std::string,
602  boost::shared_ptr<MoFEM::FEMethod>>(
603  dm, fe_name, method, pre_only, post_only);
604 }
static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:567

◆ 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 545 of file DMMMoFEM.cpp.

548  {
549  return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
551  dm, fe_name, method, pre_only, post_only);
552 }
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:526
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 555 of file DMMMoFEM.cpp.

558  {
559  return DMMoFEMKSPSetComputeRHS<const std::string,
560  boost::shared_ptr<MoFEM::FEMethod>,
561  boost::shared_ptr<MoFEM::BasicMethod>,
562  boost::shared_ptr<MoFEM::BasicMethod>>(
563  dm, fe_name, method, pre_only, post_only);
564 }
static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:526

◆ 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 376 of file DMMMoFEM.cpp.

376  {
378  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
380  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
381  CHKERR dm_field->mField_ptr->getInterface<CommInterface>()
382  ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
384 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ 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_elasticity.cpp, and UnsaturatedFlow.hpp.

Definition at line 947 of file DMMMoFEM.cpp.

947  {
948  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
950  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
951  dm_field->isPartitioned = is_partitioned;
953 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMMoFEMSetKspCtx()

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

set MoFEM::KspCtx data structure

Definition at line 909 of file DMMMoFEM.cpp.

910  {
911  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
913  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
914  dm_field->kspCtx = ksp_ctx;
916 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMMoFEMSetSnesCtx()

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

Set MoFEM::SnesCtx data structure.

Definition at line 935 of file DMMMoFEM.cpp.

936  {
937  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
939  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
940  dm_field->snesCtx = snes_ctx;
942 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 367 of file DMMMoFEM.cpp.

367  {
369  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
371  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
372  dm_field->isSquareMatrix = square_problem;
374 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 983 of file DMMMoFEM.cpp.

983  {
984  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
986  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
987  dm_field->tsCtx = ts_ctx;
989 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

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

Definition at line 626 of file DMMMoFEM.cpp.

629  {
630  return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
632  dm, fe_name, method, pre_only, post_only);
633 }
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:607
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 636 of file DMMMoFEM.cpp.

639  {
640  return DMMoFEMSNESSetFunction<const std::string,
641  boost::shared_ptr<MoFEM::FEMethod>,
642  boost::shared_ptr<MoFEM::BasicMethod>,
643  boost::shared_ptr<MoFEM::BasicMethod>>(
644  dm, fe_name, method, pre_only, post_only);
645 }
static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:607

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

Definition at line 667 of file DMMMoFEM.cpp.

670  {
671  return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
673  dm, fe_name, method, pre_only, post_only);
674 }
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:648
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 677 of file DMMMoFEM.cpp.

680  {
681  return DMMoFEMSNESSetJacobian<const std::string,
682  boost::shared_ptr<MoFEM::FEMethod>,
683  boost::shared_ptr<MoFEM::BasicMethod>,
684  boost::shared_ptr<MoFEM::BasicMethod>>(
685  dm, fe_name, method, pre_only, post_only);
686 }
static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:648

◆ 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 708 of file DMMMoFEM.cpp.

711  {
712  return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
714  dm, fe_name, method, pre_only, post_only);
716 }
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:508
static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:689
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 719 of file DMMMoFEM.cpp.

722  {
723  return DMMoFEMTSSetIFunction<const std::string,
724  boost::shared_ptr<MoFEM::FEMethod>,
725  boost::shared_ptr<MoFEM::BasicMethod>,
726  boost::shared_ptr<MoFEM::BasicMethod>>(
727  dm, fe_name, method, pre_only, post_only);
729 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:689

◆ 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 761 of file DMMMoFEM.cpp.

764  {
765  return DMMoFEMTSSetIJacobian<const std::string,
766  boost::shared_ptr<MoFEM::FEMethod>,
767  boost::shared_ptr<MoFEM::BasicMethod>,
768  boost::shared_ptr<MoFEM::BasicMethod>>(
769  dm, fe_name, method, pre_only, post_only);
770 }
static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:732

◆ 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 751 of file DMMMoFEM.cpp.

754  {
755  return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
756  MoFEM::BasicMethod *>(dm, fe_name, method,
757  pre_only, post_only);
758 }
static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:732
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 424 of file DMMMoFEM.cpp.

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

◆ 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 325 of file DMMMoFEM.cpp.

325  {
326  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
328  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
329  if (!dm->data) {
330  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
331  "data structure for MoFEM not yet created");
332  }
333  *m_field_ptr = dm_field->mField_ptr;
335 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMoFEMLoopDofs()

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

execute method for dofs on field in problem

Definition at line 513 of file DMMMoFEM.cpp.

514  {
515  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
517  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
518  ierr =
519  dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
520  *method, dm_field->rAnk, dm_field->rAnk);
521  CHKERRG(ierr);
523 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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, 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_elasticity.cpp, simple_interface.cpp, and UnsaturatedFlow.hpp.

Definition at line 496 of file DMMMoFEM.cpp.

497  {
498  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
500  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
501  ierr = DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name, method,
502  dm_field->rAnk, dm_field->rAnk);
503  CHKERRG(ierr);
505 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
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:477

◆ 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 508 of file DMMMoFEM.cpp.

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

◆ 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 477 of file DMMMoFEM.cpp.

479  {
481  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
482  ierr = dm_field->mField_ptr->loop_finite_elements(
483  dm_field->problemPtr, fe_name, *method, low_rank, up_rank);
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:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 489 of file DMMMoFEM.cpp.

491  {
492  return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
493  low_rank, up_rank);
494 }
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:477

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

Definition at line 446 of file DMMMoFEM.cpp.

447  {
448  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
450  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
451  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
452  dm_field->problemPtr, COL, g, 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:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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, EshelbianPlasticity.cpp, MagneticElement.hpp, mesh_smoothing.cpp, minimal_surface_area.cpp, reaction_diffusion_equation.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and UnsaturatedFlow.hpp.

Definition at line 434 of file DMMMoFEM.cpp.

435  {
437  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
439  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
440  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
441  dm_field->problemPtr, COL, l, mode, scatter_mode);
442  CHKERRG(ierr);
444 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 467 of file DMMMoFEM.cpp.

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

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

Definition at line 457 of file DMMMoFEM.cpp.

457  {
458  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
460  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
461  ierr = dm_field->mField_ptr->problem_basic_method_preProcess(
462  dm_field->problemPtr, *method);
463  CHKERRG(ierr);
465 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ 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 392 of file PCMGSetUpViaApproxOrders.cpp.

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

◆ 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, hello_world.cpp, loop_entities.cpp, lorentz_force.cpp, MagneticElement.hpp, mesh_smoothing.cpp, minimal_surface_area.cpp, reaction_diffusion_equation.cpp, Remodeling.cpp, simple_elasticity.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, and unsaturated_transport.cpp.

Definition at line 53 of file DMMMoFEM.cpp.

53  {
55  CHKERR DMRegister(sname, DMCreate_MoFEM);
57 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:81
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ DMSetFromOptions_MoFEM()

PetscErrorCode MoFEM::DMSetFromOptions_MoFEM ( DM  dm)

Set options for MoFEM DM

Definition at line 1062 of file DMMMoFEM.cpp.

1062  {
1063 #endif
1064 
1065  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1067  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1068 #if PETSC_VERSION_GE(3, 5, 3)
1069  ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1070  CHKERRG(ierr);
1071 #else
1072  ierr = PetscOptionsHead("DMMoFEM Options");
1073  CHKERRG(ierr);
1074 #endif
1075  ierr = PetscOptionsBool("-dm_is_partitioned",
1076  "set if mesh is partitioned (works which native MOAB "
1077  "file format, i.e. h5m",
1078  "DMSetUp", dm_field->isPartitioned,
1079  &dm_field->isPartitioned, NULL);
1080  CHKERRG(ierr);
1082 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508

◆ DMSetOperators_MoFEM()

PetscErrorCode MoFEM::DMSetOperators_MoFEM ( DM  dm)

Set operators for MoFEM dm.

Parameters
dm
Returns
error code

Definition at line 59 of file DMMMoFEM.cpp.

59  {
60  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
62 
63  dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
64  dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
65  dm->ops->creatematrix = DMCreateMatrix_MoFEM;
66  dm->ops->setup = DMSetUp_MoFEM;
67  dm->ops->destroy = DMDestroy_MoFEM;
68  dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
69  dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
70  dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
71  dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
72  dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
73  dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
74 
75  // Default matrix type
76  CHKERR DMSetMatType(dm, MATMPIAIJ);
77 
79 }
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1178
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1084
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1215
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1018
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:89
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM.
Definition: DMMMoFEM.cpp:991
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1186
#define CHKERR
Inline error check.
Definition: definitions.h:596
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1250
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1244
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1062
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM.
Definition: DMMMoFEM.cpp:1009

◆ DMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSetUp_MoFEM ( DM  dm)

sets up the MoFEM structures inside a DM object

Definition at line 1084 of file DMMMoFEM.cpp.

1084  {
1085  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1086  ProblemsManager *prb_mng_ptr;
1088  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1089  CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1090 
1091  if (dm_field->isCompDM) {
1092  // It is composite probelm
1093  CHKERR prb_mng_ptr->buildCompsedProblem(
1094  dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1095  dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1096  } else {
1097  if (dm_field->isPartitioned) {
1098  CHKERR prb_mng_ptr->buildProblemOnDistributedMesh(
1099  dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1100  dm_field->verbosity);
1101  } else {
1102  CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1103  dm_field->isSquareMatrix == PETSC_TRUE,
1104  dm_field->verbosity);
1105  CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1106  dm_field->verbosity);
1107  }
1108  }
1109 
1110  // Partition finite elements
1111  if (dm_field->isPartitioned) {
1112  CHKERR prb_mng_ptr->partitionFiniteElements(
1113  dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1114  CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1115  dm_field->problemName, dm_field->verbosity);
1116  } else {
1117  // partition finite elemnets
1118  CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1119  -1, -1, dm_field->verbosity);
1120  // Get ghost DOFs
1121  CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1122  dm_field->verbosity);
1123  }
1124 
1125  // Set flag that problem is build and partitioned
1126  dm_field->isProblemBuild = PETSC_TRUE;
1127 
1129 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

◆ DMSubDMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSubDMSetUp_MoFEM ( DM  subdm)

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

Definition at line 1131 of file DMMMoFEM.cpp.

1131  {
1132  PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1133  ProblemsManager *prb_mng_ptr;
1135 
1136  DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1137 
1138  // build sub dm problem
1139  CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1140 
1141  map<std::string, std::pair<EntityType, EntityType>> *entity_map_row = nullptr;
1142  map<std::string, std::pair<EntityType, EntityType>> *entity_map_col = nullptr;
1143 
1144  if (subdm_field->mapTypeRow)
1145  entity_map_row = subdm_field->mapTypeRow.get();
1146  if (subdm_field->mapTypeCol)
1147  entity_map_row = subdm_field->mapTypeCol.get();
1148 
1149  CHKERR prb_mng_ptr->buildSubProblem(
1150  subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1151  subdm_field->problemMainOfSubPtr->getName(),
1152  subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1153  subdm_field->verbosity);
1154 
1155  // partition problem
1156  subdm_field->isPartitioned = subdm_field->isPartitioned;
1157  if (subdm_field->isPartitioned) {
1158  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1159  0, subdm_field->sIze,
1160  subdm_field->verbosity);
1161  // set ghost nodes
1162  CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1163  subdm_field->problemName, subdm_field->verbosity);
1164  } else {
1165  // partition finite elements
1166  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1167  -1, -1, subdm_field->verbosity);
1168  // set ghost nodes
1169  CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1170  subdm_field->verbosity);
1171  }
1172 
1173  subdm_field->isProblemBuild = PETSC_TRUE;
1174 
1176 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407