v0.12.1
Functions | Variables
Distributed mesh manager

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

Collaboration diagram for Distributed mesh manager:

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 problem More...
 
PetscErrorCode MoFEM::DMMoFEMGetSquareProblem (DM dm, PetscBool *square_problem)
 get squared problem 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, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
 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, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
 Executes FEMethod for finite elements in DM. More...
 
PetscErrorCode MoFEM::DMoFEMLoopFiniteElements (DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
 Executes FEMethod for finite elements in DM. More...
 
PetscErrorCode MoFEM::DMoFEMLoopFiniteElements (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
 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 KSP right hand side evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMKSPSetComputeRHS (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set KSP right hand side evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMKSPSetComputeOperators (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 Set KSP operators and push mofem finite element methods. More...
 
PetscErrorCode MoFEM::DMMoFEMKSPSetComputeOperators (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 Set KSP operators and push mofem finite element methods. More...
 
PetscErrorCode MoFEM::DMMoFEMSNESSetFunction (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 set SNES residual evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMSNESSetFunction (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set SNES residual evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMSNESSetJacobian (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 set SNES Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMSNESSetJacobian (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set SNES Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetIFunction (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 set TS implicit function evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetIFunction (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set TS implicit function evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetIJacobian (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set TS Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetIJacobian (DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
 set TS Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetI2Function (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set TS implicit function evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetI2Function (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::DMMoFEMTSSetI2Jacobian (DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
 set TS Jacobian evaluation function More...
 
PetscErrorCode MoFEM::DMMoFEMTSSetI2Jacobian (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 structure. 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)
 DMShellSetCreateGlobalVector. More...
 
PetscErrorCode MoFEM::DMCreateGlobalVector_MoFEM (DM dm, SmartPetscObj< Vec > &g_ptr)
 DMShellSetCreateGlobalVector. More...
 
PetscErrorCode MoFEM::DMCreateLocalVector_MoFEM (DM dm, Vec *l)
 DMShellSetCreateLocalVector. 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 row. More...
 
PetscErrorCode MoFEM::DMMoFEMAddColCompositeProblem (DM dm, const char prb_name[])
 Add problem to composite DM on col. More...
 
PetscErrorCode MoFEM::DMMoFEMGetIsCompDM (DM dm, PetscBool *is_comp_dm)
 Get if this DM is composite DM. More...
 
PetscErrorCode MoFEM::DMGlobalToLocalBegin_MoFEM (DM dm, Vec, InsertMode, Vec)
 
PetscErrorCode MoFEM::DMGlobalToLocalEnd_MoFEM (DM dm, Vec, InsertMode, Vec)
 
PetscErrorCode MoFEM::DMLocalToGlobalBegin_MoFEM (DM, Vec, InsertMode, Vec)
 
PetscErrorCode MoFEM::DMLocalToGlobalEnd_MoFEM (DM, Vec, InsertMode, Vec)
 
PetscErrorCode MoFEM::DMCreateFieldIS_MoFEM (DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
 
PetscErrorCode MoFEM::DMMoFEMGetFieldIS (DM dm, RowColData rc, const char field_name[], IS *is)
 get field is in the problem More...
 
MoFEMErrorCode DMMGViaApproxOrdersPushBackCoarseningIS (DM, IS is, Mat A, Mat *subA, bool create_sub_matrix, bool shell_sub_a)
 Push back coarsening level to MG via approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS (DM)
 Pop is form MG via approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS (DM)
 Clear approximation orders. More...
 
MoFEMErrorCode DMRegister_MGViaApproxOrders (const char sname[])
 Register DM for Multi-Grid via approximation orders. More...
 
MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders (DM dm, Mat *M)
 Create matrix for Multi-Grid via approximation orders. More...
 
MoFEMErrorCode DMCoarsen_MGViaApproxOrders (DM dm, MPI_Comm comm, DM *dmc)
 Coarsen DM. More...
 

Variables

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

Detailed Description

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

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

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

Function Documentation

◆ DMCoarsen_MGViaApproxOrders()

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

Coarsen DM.

Not used directly by user

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

Definition at line 473 of file PCMGSetUpViaApproxOrders.cpp.

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

◆ DMCreate_MoFEM()

PetscErrorCode MoFEM::DMCreate_MoFEM ( DM  dm)

Create dm data structure with MoFEM data structure.

Definition at line 87 of file DMMMoFEM.cpp.

87  {
88  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
90  dm->data = new DMCtx();
93 }
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:65

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

1367  {
1368  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1370 
1371  if (numFields) {
1372  *numFields = 0;
1373  }
1374  if (fieldNames) {
1375  *fieldNames = NULL;
1376  }
1377  if (fields) {
1378  *fields = NULL;
1379  }
1380 
1381  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1382  auto fields_ptr = dm_field->mField_ptr->get_fields();
1383  Field_multiIndex::iterator fit, hi_fit;
1384  fit = fields_ptr->begin();
1385  hi_fit = fields_ptr->end();
1386  *numFields = std::distance(fit, hi_fit);
1387 
1388  if (fieldNames) {
1389  CHKERR PetscMalloc1(*numFields, fieldNames);
1390  }
1391  if (fields) {
1392  CHKERR PetscMalloc1(*numFields, fields);
1393  }
1394 
1395  for (int f = 0; fit != hi_fit; fit++, f++) {
1396  if (fieldNames) {
1397  CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1398  (char **)&(*fieldNames)[f]);
1399  }
1400  if (fields) {
1401  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1402  ->isCreateProblemFieldAndRank(
1403  dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1404  fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1405  }
1406  }
1407 
1409 }
@ ROW
Definition: definitions.h:136

◆ DMCreateGlobalVector_MoFEM() [1/2]

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

DMShellSetCreateGlobalVector.

sets the routine to create a global vector associated with the shell DM

Definition at line 1112 of file DMMMoFEM.cpp.

1112  {
1113  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1115  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1116  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1117  dm_field->problemName, COL, g_ptr);
1118  CHKERR VecSetDM(g_ptr, dm);
1120 }
@ COL
Definition: definitions.h:136
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453

◆ DMCreateGlobalVector_MoFEM() [2/2]

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

DMShellSetCreateGlobalVector.

sets the routine to create a global vector associated with the shell DM

Examples
elasticity.cpp, elasticity_mixed_formulation.cpp, ep.cpp, reaction_diffusion.cpp, and test_cache_on_entities.cpp.

Definition at line 1102 of file DMMMoFEM.cpp.

1102  {
1103  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1105  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1106  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1107  dm_field->problemName, COL, g);
1108  CHKERR VecSetDM(*g, dm);
1110 }
constexpr double g

◆ DMCreateLocalVector_MoFEM()

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

DMShellSetCreateLocalVector.

sets the routine to create a local vector associated with the shell DM

Definition at line 1122 of file DMMMoFEM.cpp.

1122  {
1123  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1125  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1126  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1127  dm_field->problemName, COL, l);
1128  CHKERR VecSetDM(*l, dm);
1130 }
FTensor::Index< 'l', 3 > l

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

441  {
442 
443  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
445  GET_DM_FIELD(dm);
446 
447  int leveldown = dm->leveldown;
448 
449  if (dm_field->kspOperators.empty()) {
451  } else {
452  MPI_Comm comm;
453  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
454  if (dm_field->kspOperators.empty()) {
455  SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
456  "data inconsistency, operator can not be set");
457  }
458  if (static_cast<int>(dm_field->kspOperators.size()) < leveldown) {
459  SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
460  "data inconsistency, no IS for that level");
461  }
462  *M = dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown];
463  CHKERR PetscObjectReference((PetscObject)*M);
464  CHKERR MatSetDM(*M, dm);
465  }
466 
467  PetscInfo1(dm, "Create Matrix DMMGViaApproxOrders leveldown = %d\n",
468  leveldown);
469 
471 }
static Index< 'M', 3 > M
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1132

◆ 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
EshelbianPlasticity.cpp, dm_build_partitioned_mesh.cpp, eigen_elastic.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, ep.cpp, and reaction_diffusion.cpp.

Definition at line 1132 of file DMMMoFEM.cpp.

1132  {
1133  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1135  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1136  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1137  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1138  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1139  M);
1140  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1141  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1142  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1143  M);
1144  } else {
1145  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1146  "Matrix type not implemented");
1147  }
1148  CHKERR MatSetDM(*M, dm);
1150 }
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45

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

1152  {
1153  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1155  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1156  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1157  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1158  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1159  M);
1160  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1161  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1162  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1163  M);
1164  } else {
1165  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1166  "Matrix type not implemented");
1167  }
1168  CHKERR MatSetDM(M, dm);
1170 }

◆ DMDestroy_MoFEM()

PetscErrorCode MoFEM::DMDestroy_MoFEM ( DM  dm)

Destroys dm with MoFEM data structure.

destroy the MoFEM structure

Definition at line 95 of file DMMMoFEM.cpp.

95  {
96  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
97  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
99 
100  MPI_Comm comm;
101  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
102 
103  int result;
104  MPI_Comm_compare(comm, PETSC_COMM_SELF, &result);
105  if (result == MPI_IDENT)
106  MOFEM_LOG("DMSELF", Sev::noisy)
107  << "MoFEM DM destroy for problem " << dm_field->problemName
108  << " referenceNumber " << dm_field->referenceNumber;
109  else
110  MOFEM_LOG("DMWORLD", Sev::noisy)
111  << "MoFEM DM destroy for problem " << dm_field->problemName
112  << " referenceNumber " << dm_field->referenceNumber;
113 
114  if (dm_field->referenceNumber == 0) {
115  if (dm_field->destroyProblem) {
116 
117  if (dm_field->mField_ptr->check_problem(dm_field->problemName)) {
118  dm_field->mField_ptr->delete_problem(dm_field->problemName);
119  } // else problem has to be deleted by the user
120  }
121 
122  delete static_cast<DMCtx *>(dm->data);
123 
124  } else
125  --dm_field->referenceNumber;
126 
128 }
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:309

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

1295  {
1296  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1298  CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1300 }

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

1302  {
1304  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1306 
1307  CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1308 
1309  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1310  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1311  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1312 
1313  double *array_loc, *array_glob;
1314  CHKERR VecGetArray(l, &array_loc);
1315  CHKERR VecGetArray(g, &array_glob);
1316  switch (mode) {
1317  case INSERT_VALUES:
1318  cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1319  break;
1320  case ADD_VALUES:
1321  cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1322  break;
1323  default:
1324  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1325  }
1326  CHKERR VecRestoreArray(l, &array_loc);
1327  CHKERR VecRestoreArray(g, &array_glob);
1329 }

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

1332  {
1333 
1334  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1336 
1337  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1338  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1339  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1340 
1341  double *array_loc, *array_glob;
1342  CHKERR VecGetArray(l, &array_loc);
1343  CHKERR VecGetArray(g, &array_glob);
1344  switch (mode) {
1345  case INSERT_VALUES:
1346  cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1347  break;
1348  case ADD_VALUES:
1349  cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1350  break;
1351  default:
1352  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1353  }
1354  CHKERR VecRestoreArray(l, &array_loc);
1355  CHKERR VecRestoreArray(g, &array_glob);
1356 
1358 }

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

1360  {
1361  //
1364 }

◆ DMMGViaApproxOrdersClearCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS ( DM  dm)

Clear approximation orders.

Parameters
DMdm
Returns
Error code

Definition at line 287 of file PCMGSetUpViaApproxOrders.cpp.

287  {
288  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
290  GET_DM_FIELD(dm);
291  CHKERR dm_field->destroyCoarseningIS();
292  PetscInfo(dm, "Clear DMs data structures\n");
294 }

◆ DMMGViaApproxOrdersPopBackCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS ( DM  dm)

Pop is form MG via approximation orders.

Parameters
DMdm
ispop back IS
Returns
error code

Definition at line 271 of file PCMGSetUpViaApproxOrders.cpp.

271  {
272  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
274  GET_DM_FIELD(dm);
275  if (dm_field->coarseningIS.back()) {
276  CHKERR ISDestroy(&dm_field->coarseningIS.back());
277  dm_field->coarseningIS.pop_back();
278  }
279  if (dm_field->kspOperators.back()) {
280  CHKERR MatDestroy(&dm_field->kspOperators.back());
281  }
282  dm_field->kspOperators.pop_back();
283  PetscInfo(dm, "Pop back IS to DMMGViaApproxOrders\n");
285 }

◆ DMMGViaApproxOrdersPushBackCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersPushBackCoarseningIS ( DM  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 217 of file PCMGSetUpViaApproxOrders.cpp.

220  {
221  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
223  GET_DM_FIELD(dm);
224  dm_field->coarseningIS.push_back(is);
225  dm_field->shellMatrixCtxPtr.push_back(new PCMGSubMatrixCtx_private(A, is));
226  if (is) {
227  CHKERR PetscObjectReference((PetscObject)is);
228  }
229  if (is) {
230  IS is2 = is;
231  if (dm_field->aO) {
232  CHKERR ISDuplicate(is, &is2);
233  CHKERR ISCopy(is, is2);
234  CHKERR AOApplicationToPetscIS(dm_field->aO, is2);
235  }
236  if (create_sub_matrix) {
237  if (shell_sub_a) {
238  int n, N;
239  CHKERR ISGetSize(is, &N);
240  CHKERR ISGetLocalSize(is, &n);
241  MPI_Comm comm;
242  CHKERR PetscObjectGetComm((PetscObject)A, &comm);
243  CHKERR MatCreateShell(comm, n, n, N, N,
244  &(dm_field->shellMatrixCtxPtr.back()), subA);
245  CHKERR MatShellSetOperation(*subA, MATOP_MULT,
246  (void (*)(void))sub_mat_mult);
247  CHKERR MatShellSetOperation(*subA, MATOP_MULT_ADD,
248  (void (*)(void))sub_mat_mult_add);
249  CHKERR MatShellSetOperation(*subA, MATOP_SOR,
250  (void (*)(void))sub_mat_sor);
251  } else {
252 #if PETSC_VERSION_GE(3, 8, 0)
253  CHKERR MatCreateSubMatrix(A, is2, is2, MAT_INITIAL_MATRIX, subA);
254 #else
255  CHKERR MatGetSubMatrix(A, is2, is2, MAT_INITIAL_MATRIX, subA);
256 #endif
257  }
258  }
259  if (dm_field->aO) {
260  CHKERR ISDestroy(&is2);
261  }
262  dm_field->kspOperators.push_back(*subA);
263  CHKERR PetscObjectReference((PetscObject)(*subA));
264  } else {
265  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
266  }
267  PetscInfo(dm, "Push back IS to DMMGViaApproxOrders\n");
269 }
static Index< 'n', 3 > n
MoFEMErrorCode sub_mat_mult(Mat a, Vec x, Vec f)
MoFEMErrorCode sub_mat_sor(Mat mat, Vec b, PetscReal omega, MatSORType flag, PetscReal shift, PetscInt its, PetscInt lits, Vec x)
MoFEMErrorCode sub_mat_mult_add(Mat a, Vec x, Vec f)
double A
const int N
Definition: speed_test.cpp:3

◆ DMMoFEMAddColCompositeProblem()

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

Add problem to composite DM on col.

This 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 340 of file DMMMoFEM.cpp.

340  {
341  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
343  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
344  if (!dm->data) {
345  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
346  "data structure for MoFEM not yet created");
347  }
348  if (!dm_field->isCompDM) {
349  dm_field->isCompDM = PETSC_TRUE;
350  }
351  if (dm_field->isSquareMatrix) {
352  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
353  "No need to add problem on column when problem block structurally "
354  "symmetric");
355  }
356  dm_field->colCompPrb.push_back(prb_name);
358 }
@ MOFEM_INVALID_DATA
Definition: definitions.h:49

◆ 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
EshelbianPlasticity.cpp, MagneticElement.hpp, Remodeling.cpp, UnsaturatedFlow.hpp, 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, lorentz_force.cpp, minimal_surface_area.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, photon_diffusion.cpp, simple_contact.cpp, simple_contact_thermal.cpp, simple_elasticity.cpp, and test_jacobian_of_simple_contact_element.cpp.

Definition at line 458 of file DMMMoFEM.cpp.

458  {
459  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
461  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
462  ierr = dm_field->mField_ptr->modify_problem_add_finite_element(
463  dm_field->problemName, fe_name);
464  CHKERRG(ierr);
466 }
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87

◆ DMMoFEMAddRowCompositeProblem()

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

Add problem to composite DM on row.

This 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 322 of file DMMMoFEM.cpp.

322  {
323  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
325  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
326  if (!dm->data) {
327  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
328  "data structure for MoFEM not yet created");
329  }
330  if (!dm_field->isCompDM) {
331  dm_field->isCompDM = PETSC_TRUE;
332  }
333  dm_field->rowCompPrb.push_back(prb_name);
334  if (dm_field->isSquareMatrix) {
335  dm_field->colCompPrb.push_back(prb_name);
336  }
338 }

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

255  {
256  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
258  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
259  if (!dm->data) {
260  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
261  "data structure for MoFEM not yet created");
262  }
263  if (!dm_field->isSubDM) {
264  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
265  }
266  dm_field->colFields.push_back(field_name);
267  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
268  if (!dm_field->mapTypeCol)
269  dm_field->mapTypeCol = boost::make_shared<
270  std::map<std::string, std::pair<EntityType, EntityType>>>();
271  (*dm_field->mapTypeCol)[field_name] =
272  std::pair<EntityType, EntityType>(lo_type, hi_type);
273  }
275 }

◆ 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
EshelbianPlasticity.cpp, analytical_poisson_field_split.cpp, cell_forces.cpp, and dm_create_subdm.cpp.

Definition at line 231 of file DMMMoFEM.cpp.

232  {
233  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
235  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
236  if (!dm->data) {
237  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
238  "data structure for MoFEM not yet created");
239  }
240  if (!dm_field->isSubDM) {
241  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
242  }
243  dm_field->rowFields.push_back(field_name);
244  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
245  if (!dm_field->mapTypeRow)
246  dm_field->mapTypeRow = boost::make_shared<
247  std::map<std::string, std::pair<EntityType, EntityType>>>();
248  (*dm_field->mapTypeRow)[field_name] =
249  std::pair<EntityType, EntityType>(lo_type, hi_type);
250  }
252 }

◆ 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
EshelbianPlasticity.cpp, MagneticElement.hpp, Remodeling.cpp, UnsaturatedFlow.hpp, cell_forces.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, dm_partitioned_no_field.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, lorentz_force.cpp, minimal_surface_area.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, simple_contact.cpp, simple_contact_thermal.cpp, and test_jacobian_of_simple_contact_element.cpp.

Definition at line 130 of file DMMMoFEM.cpp.

133  {
135 
136  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
137  if (!dm->data) {
138  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
139  "data structure for MoFEM not yet created");
140  }
141  if (!m_field_ptr) {
142  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
143  "DM function not implemented into MoFEM");
144  }
145  dm_field->mField_ptr = m_field_ptr;
146  dm_field->problemName = problem_name;
147  if (!m_field_ptr->check_problem(dm_field->problemName)) {
148  // problem is not defined, declare problem internally set bool to
149  // destroyProblem problem with DM
150  dm_field->destroyProblem = PETSC_TRUE;
151  CHKERR dm_field->mField_ptr->add_problem(dm_field->problemName, MF_EXCL,
152  dm_field->verbosity);
153  } else {
154  dm_field->destroyProblem = PETSC_FALSE;
155  }
156  CHKERR dm_field->mField_ptr->modify_problem_ref_level_add_bit(
157  dm_field->problemName, bit_level);
158  CHKERR dm_field->mField_ptr->modify_problem_mask_ref_level_add_bit(
159  dm_field->problemName, bit_mask);
160  dm_field->kspCtx =
161  boost::shared_ptr<KspCtx>(new KspCtx(*m_field_ptr, problem_name));
162  dm_field->snesCtx =
163  boost::shared_ptr<SnesCtx>(new SnesCtx(*m_field_ptr, problem_name));
164  dm_field->tsCtx =
165  boost::shared_ptr<TsCtx>(new TsCtx(*m_field_ptr, problem_name));
166 
167  MPI_Comm comm;
168  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
169  int result = 0;
170  MPI_Comm_compare(comm, m_field_ptr->get_comm(), &result);
171  if (result > MPI_CONGRUENT) {
172  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
173  "MoFEM and DM using different communicators");
174  }
175  MPI_Comm_size(comm, &dm_field->sIze);
176  MPI_Comm_rank(comm, &dm_field->rAnk);
177 
178  // problem structure
179  CHKERR dm_field->mField_ptr->get_problem(dm_field->problemName,
180  &dm_field->problemPtr);
181 
182  MPI_Comm_compare(comm, PETSC_COMM_SELF, &result);
183  if (result == MPI_IDENT)
184  MOFEM_LOG("DMSELF", Sev::noisy)
185  << "MoFEM DM created for problem " << dm_field->problemName;
186  else
187  MOFEM_LOG("DMWORLD", Sev::noisy)
188  << "MoFEM DM created for problem " << dm_field->problemName;
189 
191 }
@ MF_EXCL
Definition: definitions.h:112
virtual bool check_problem(const std::string name)=0
check if problem exist
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
EshelbianPlasticity.cpp, analytical_poisson_field_split.cpp, cell_forces.cpp, and dm_create_subdm.cpp.

Definition at line 210 of file DMMMoFEM.cpp.

210  {
212 
213  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
214  if (!dm->data) {
215  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
216  "data structure for MoFEM not yet created");
217  }
218  CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
219  dm_field->problemPtr->getBitRefLevel());
220 
221  DMCtx *subdm_field = (DMCtx *)subdm->data;
222  subdm_field->isSubDM = PETSC_TRUE;
223  subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
224  subdm_field->isPartitioned = dm_field->isPartitioned;
225  subdm_field->isSquareMatrix = PETSC_FALSE;
226  subdm->ops->setup = DMSubDMSetUp_MoFEM;
227 
229 }
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:130
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1247

◆ 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);
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMMoFEM.cpp:1411

Definition at line 1411 of file DMMMoFEM.cpp.

1412  {
1413  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1415  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1416  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1417  ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1418  field_name, 0, 1000, is);
1420 }

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

360  {
362  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
364  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
365  *is_comp_dm = dm_field->isCompDM;
367 }

◆ DMMoFEMGetIsPartitioned()

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

get if read mesh is partitioned

Definition at line 1069 of file DMMMoFEM.cpp.

1069  {
1070  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1072  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1073  *is_partitioned = dm_field->isPartitioned;
1075 }

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

277  {
279  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
281  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
282  *is_sub_dm = dm_field->isSubDM;
284 }

◆ DMMoFEMGetKspCtx() [1/2]

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

get MoFEM::KspCtx data structure

Definition at line 1012 of file DMMMoFEM.cpp.

1012  {
1013  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1015  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1016  const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
1018 }

◆ DMMoFEMGetKspCtx() [2/2]

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

get MoFEM::KspCtx data structure

Definition at line 1003 of file DMMMoFEM.cpp.

1003  {
1004  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1006  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1007  *ksp_ctx = dm_field->kspCtx.get();
1009 }

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

435  {
437  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
439  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
440 
441  MPI_Comm comm;
442  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
443  CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
444  layout);
446 }

◆ DMMoFEMGetProblemPtr()

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

Get pointer to problem data structure.

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

Definition at line 381 of file DMMMoFEM.cpp.

381  {
382  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
384  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
385  if (!dm->data) {
386  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
387  "data structure for MoFEM not yet created");
388  }
389  *problem_ptr = dm_field->problemPtr;
391 }

◆ DMMoFEMGetSnesCtx() [1/2]

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

get MoFEM::SnesCtx data structure

Definition at line 1038 of file DMMMoFEM.cpp.

1038  {
1039  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1041  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1042  const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
1044 }

◆ DMMoFEMGetSnesCtx() [2/2]

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

get MoFEM::SnesCtx data structure

Examples
minimal_surface_area.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, simple_contact.cpp, simple_contact_thermal.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 1029 of file DMMMoFEM.cpp.

1029  {
1030  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1032  DMCtx *dm_field = (DMCtx *)dm->data;
1033  *snes_ctx = dm_field->snesCtx.get();
1035 }

◆ DMMoFEMGetSquareProblem()

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

get squared problem

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

Definition at line 448 of file DMMMoFEM.cpp.

448  {
451  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
453  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
454  *square_problem = dm_field->isSquareMatrix;
456 }

◆ DMMoFEMGetTsCtx() [1/2]

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

get MoFEM::TsCtx data structure

Definition at line 1085 of file DMMMoFEM.cpp.

1086  {
1087  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1089  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1090  const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1092 }

◆ DMMoFEMGetTsCtx() [2/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 1077 of file DMMMoFEM.cpp.

1077  {
1078  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1080  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1081  *ts_ctx = dm_field->tsCtx.get();
1083 }

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

636  {
637  return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
640  dm, fe_name, method, pre_only, post_only);
641 }
static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:614
Data structure to exchange data between mofem and User Loop Methods.
structure for User Loop Methods on finite elements

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

647  {
648  return DMMoFEMKSPSetComputeOperators<const std::string,
649  boost::shared_ptr<MoFEM::FEMethod>>(
650  dm, fe_name, method, pre_only, post_only);
651 }

◆ DMMoFEMKSPSetComputeRHS() [1/2]

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

set KSP right hand side evaluation function

Examples
analytical_poisson.cpp, and simple_elasticity.cpp.

Definition at line 592 of file DMMMoFEM.cpp.

595  {
596  return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
598  dm, fe_name, method, pre_only, post_only);
599 }
static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:573

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

605  {
606  return DMMoFEMKSPSetComputeRHS<const std::string,
607  boost::shared_ptr<MoFEM::FEMethod>,
608  boost::shared_ptr<MoFEM::BasicMethod>,
609  boost::shared_ptr<MoFEM::BasicMethod>>(
610  dm, fe_name, method, pre_only, post_only);
611 }

◆ 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);
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:228
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:420

Definition at line 420 of file DMMMoFEM.cpp.

420  {
422  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
424  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
425  CHKERR dm_field->mField_ptr->getInterface<CommInterface>()
426  ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
428 }

◆ DMMoFEMSetIsPartitioned()

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

sets if read mesh is partitioned

get if read mesh is partitioned

Examples
EshelbianPlasticity.cpp, MagneticElement.hpp, Remodeling.cpp, UnsaturatedFlow.hpp, cell_forces.cpp, dm_build_partitioned_mesh.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, simple_contact.cpp, simple_contact_thermal.cpp, simple_elasticity.cpp, and test_jacobian_of_simple_contact_element.cpp.

Definition at line 1058 of file DMMMoFEM.cpp.

1058  {
1059  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1061  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1062  dm_field->isPartitioned = is_partitioned;
1064 }

◆ DMMoFEMSetKspCtx()

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

set MoFEM::KspCtx data structure

Definition at line 1020 of file DMMMoFEM.cpp.

1021  {
1022  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1024  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1025  dm_field->kspCtx = ksp_ctx;
1027 }

◆ DMMoFEMSetSnesCtx()

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

Set MoFEM::SnesCtx data structure.

Definition at line 1046 of file DMMMoFEM.cpp.

1047  {
1048  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1050  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1051  dm_field->snesCtx = snes_ctx;
1053 }

◆ DMMoFEMSetSquareProblem()

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

set squared problem

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

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

Definition at line 411 of file DMMMoFEM.cpp.

411  {
413  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
415  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
416  dm_field->isSquareMatrix = square_problem;
418 }

◆ DMMoFEMSetTsCtx()

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

Set MoFEM::TsCtx data structure.

It take over pointer, do not delete it, DM will destroy pointer when is destroyed.

Definition at line 1094 of file DMMMoFEM.cpp.

1094  {
1095  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1097  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1098  dm_field->tsCtx = ts_ctx;
1100 }

◆ 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, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, simple_contact.cpp, simple_contact_thermal.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 673 of file DMMMoFEM.cpp.

676  {
677  return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
679  dm, fe_name, method, pre_only, post_only);
680 }
static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:654

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

686  {
687  return DMMoFEMSNESSetFunction<const std::string,
688  boost::shared_ptr<MoFEM::FEMethod>,
689  boost::shared_ptr<MoFEM::BasicMethod>,
690  boost::shared_ptr<MoFEM::BasicMethod>>(
691  dm, fe_name, method, pre_only, post_only);
692 }

◆ 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, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, simple_contact.cpp, simple_contact_thermal.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, and testing_jacobian_of_hook_scaled_with_density_element.cpp.

Definition at line 714 of file DMMMoFEM.cpp.

717  {
718  return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
720  dm, fe_name, method, pre_only, post_only);
721 }
static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:695

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

727  {
728  return DMMoFEMSNESSetJacobian<const std::string,
729  boost::shared_ptr<MoFEM::FEMethod>,
730  boost::shared_ptr<MoFEM::BasicMethod>,
731  boost::shared_ptr<MoFEM::BasicMethod>>(
732  dm, fe_name, method, pre_only, post_only);
733 }

◆ DMMoFEMTSSetI2Function() [1/2]

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

set TS implicit function evaluation function

Definition at line 899 of file DMMMoFEM.cpp.

902  {
903  return DMMoFEMTSSetI2Function<const char *, MoFEM::FEMethod *,
905  dm, fe_name, method, pre_only, post_only);
907 }
static PetscErrorCode DMMoFEMTSSetI2Function(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:880

◆ DMMoFEMTSSetI2Function() [2/2]

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

set TS implicit function evaluation function

Examples
EshelbianPlasticity.cpp.

Definition at line 910 of file DMMMoFEM.cpp.

913  {
914  return DMMoFEMTSSetI2Function<const std::string,
915  boost::shared_ptr<MoFEM::FEMethod>,
916  boost::shared_ptr<MoFEM::BasicMethod>,
917  boost::shared_ptr<MoFEM::BasicMethod>>(
918  dm, fe_name, method, pre_only, post_only);
920 }

◆ DMMoFEMTSSetI2Jacobian() [1/2]

PetscErrorCode MoFEM::DMMoFEMTSSetI2Jacobian ( 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 942 of file DMMMoFEM.cpp.

945  {
946  return DMMoFEMTSSetI2Jacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
947  MoFEM::BasicMethod *>(dm, fe_name, method,
948  pre_only, post_only);
949 }
static PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:923

◆ DMMoFEMTSSetI2Jacobian() [2/2]

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

set TS Jacobian evaluation function

Examples
EshelbianPlasticity.cpp.

Definition at line 952 of file DMMMoFEM.cpp.

955  {
956  return DMMoFEMTSSetI2Jacobian<const std::string,
957  boost::shared_ptr<MoFEM::FEMethod>,
958  boost::shared_ptr<MoFEM::BasicMethod>,
959  boost::shared_ptr<MoFEM::BasicMethod>>(
960  dm, fe_name, method, pre_only, post_only);
961 }

◆ 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
EshelbianPlasticity.cpp, UnsaturatedFlow.hpp, bone_adaptation.cpp, plastic.cpp, reaction_diffusion.cpp, and thermo_plastic.cpp.

Definition at line 755 of file DMMMoFEM.cpp.

758  {
759  return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
761  dm, fe_name, method, pre_only, post_only);
763 }
static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:736

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

769  {
770  return DMMoFEMTSSetIFunction<const std::string,
771  boost::shared_ptr<MoFEM::FEMethod>,
772  boost::shared_ptr<MoFEM::BasicMethod>,
773  boost::shared_ptr<MoFEM::BasicMethod>>(
774  dm, fe_name, method, pre_only, post_only);
776 }

◆ DMMoFEMTSSetIJacobian() [1/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 798 of file DMMMoFEM.cpp.

801  {
802  return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
803  MoFEM::BasicMethod *>(dm, fe_name, method,
804  pre_only, post_only);
805 }
static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:779

◆ DMMoFEMTSSetIJacobian() [2/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
EshelbianPlasticity.cpp, UnsaturatedFlow.hpp, bone_adaptation.cpp, plastic.cpp, reaction_diffusion.cpp, and thermo_plastic.cpp.

Definition at line 808 of file DMMMoFEM.cpp.

811  {
812  return DMMoFEMTSSetIJacobian<const std::string,
813  boost::shared_ptr<MoFEM::FEMethod>,
814  boost::shared_ptr<MoFEM::BasicMethod>,
815  boost::shared_ptr<MoFEM::BasicMethod>>(
816  dm, fe_name, method, pre_only, post_only);
817 }

◆ DMMoFEMUnSetElement()

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

unset element from dm

Definition at line 468 of file DMMMoFEM.cpp.

468  {
469  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
471  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
472  ierr = dm_field->mField_ptr->modify_problem_unset_finite_element(
473  dm_field->problemName, fe_name);
474  CHKERRG(ierr);
476 }

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

369  {
370  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
372  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
373  if (!dm->data) {
374  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
375  "data structure for MoFEM not yet created");
376  }
377  *m_field_ptr = dm_field->mField_ptr;
379 }

◆ DMoFEMLoopDofs()

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

execute method for dofs on field in problem

Definition at line 560 of file DMMMoFEM.cpp.

561  {
562  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
564  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
565  ierr =
566  dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
567  *method, dm_field->rAnk, dm_field->rAnk);
568  CHKERRG(ierr);
570 }

◆ DMoFEMLoopFiniteElements() [1/2]

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

Executes FEMethod for finite elements in DM.

Parameters
dmMoFEM discrete manager
fe_namename of element
methodpointer to MOFEM::FEMethod
Returns
Error code
Examples
EshelbianPlasticity.cpp, HookeElement.cpp, MagneticElement.hpp, PlasticOpsMonitor.hpp, Remodeling.cpp, UnsaturatedFlow.hpp, analytical_nonlinear_poisson.cpp, analytical_poisson.cpp, analytical_poisson_field_split.cpp, approx_sphere.cpp, boundary_marker.cpp, cell_forces.cpp, continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, dm_build_partitioned_mesh.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, heat_equation.cpp, mesh_smoothing.cpp, minimal_surface_area.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, nonlinear_elastic.cpp, photon_diffusion.cpp, reaction_diffusion.cpp, shallow_wave.cpp, simple_contact.cpp, simple_contact_thermal.cpp, simple_elasticity.cpp, simple_interface.cpp, test_cache_on_entities.cpp, and wave_equation.cpp.

Definition at line 541 of file DMMMoFEM.cpp.

543  {
544  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
546  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
548  dm, fe_name, method, dm_field->rAnk, dm_field->rAnk, cache_ptr);
549  CHKERRG(ierr);
551 }
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:522

◆ DMoFEMLoopFiniteElements() [2/2]

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

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

556  {
557  return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get(), cache_ptr);
558 }
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:541

◆ DMoFEMLoopFiniteElementsUpAndLowRank() [1/2]

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

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
Examples
mortar_contact_thermal.cpp, and simple_contact_thermal.cpp.

Definition at line 522 of file DMMMoFEM.cpp.

524  {
526  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
527  ierr = dm_field->mField_ptr->loop_finite_elements(
528  dm_field->problemPtr, fe_name, *method, low_rank, up_rank, nullptr,
529  MF_EXIST, cache_ptr);
530  CHKERRG(ierr);
532 }
@ MF_EXIST
Definition: definitions.h:113

◆ 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,
CacheTupleWeakPtr  cache_ptr = CacheTupleSharedPtr() 
)

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

536  {
537  return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
538  low_rank, up_rank, cache_ptr);
539 }

◆ 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, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, simple_contact.cpp, simple_contact_thermal.cpp, and simple_elasticity.cpp.

Definition at line 490 of file DMMMoFEM.cpp.

491  {
492  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
494  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
495  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
496  dm_field->problemPtr, COL, g, mode, scatter_mode);
497  CHKERRG(ierr);
499 }

◆ 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
EshelbianPlasticity.cpp, MagneticElement.hpp, UnsaturatedFlow.hpp, approx_sphere.cpp, bone_adaptation.cpp, cell_forces.cpp, contact.cpp, eigen_elastic.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, ep.cpp, hcurl_check_approx_in_2d.cpp, heat_equation.cpp, heat_method.cpp, helmholtz.cpp, mesh_smoothing.cpp, minimal_surface_area.cpp, mixed_poisson.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, nonlinear_elastic.cpp, photon_diffusion.cpp, plastic.cpp, reaction_diffusion.cpp, scalar_check_approximation.cpp, shallow_wave.cpp, simple_contact.cpp, simple_contact_thermal.cpp, test_cache_on_entities.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, thermo_plastic.cpp, and wave_equation.cpp.

Definition at line 478 of file DMMMoFEM.cpp.

479  {
481  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
483  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
484  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
485  dm_field->problemPtr, COL, l, mode, scatter_mode);
486  CHKERRG(ierr);
488 }

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

511  {
512  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
514  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
515  ierr = dm_field->mField_ptr->problem_basic_method_postProcess(
516  dm_field->problemPtr, *method);
517  CHKERRG(ierr);
519 }

◆ DMoFEMPreProcessFiniteElements()

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

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

Examples
cell_forces.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, simple_contact.cpp, and simple_contact_thermal.cpp.

Definition at line 501 of file DMMMoFEM.cpp.

501  {
502  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
504  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
505  ierr = dm_field->mField_ptr->problem_basic_method_preProcess(
506  dm_field->problemPtr, *method);
507  CHKERRG(ierr);
509 }

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

387  {
389  CHKERR DMRegister(sname, DMCreate_MGViaApproxOrders);
391 }
MoFEMErrorCode DMCreate_MGViaApproxOrders(DM dm)
Create DM data structure for Multi-Grid via approximation orders.

◆ DMRegister_MoFEM()

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

Register MoFEM problem.

Examples
MagneticElement.hpp, Remodeling.cpp, analytical_nonlinear_poisson.cpp, analytical_poisson.cpp, analytical_poisson_field_split.cpp, approx_sphere.cpp, boundary_marker.cpp, cell_forces.cpp, contact.cpp, continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, dm_partitioned_no_field.cpp, eigen_elastic.cpp, elasticity_mixed_formulation.cpp, ep.cpp, field_blas_set_vertex_dofs.cpp, field_evaluator.cpp, hcurl_check_approx_in_2d.cpp, hcurl_curl_operator.cpp, hdiv_divergence_operator.cpp, heat_equation.cpp, heat_method.cpp, hello_world.cpp, helmholtz.cpp, loop_entities.cpp, lorentz_force.cpp, mesh_smoothing.cpp, minimal_surface_area.cpp, mixed_poisson.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, nonlinear_elastic.cpp, partition_mesh.cpp, photon_diffusion.cpp, plastic.cpp, plot_base.cpp, reaction_diffusion.cpp, scalar_check_approximation.cpp, shallow_wave.cpp, simple_contact.cpp, simple_contact_thermal.cpp, simple_elasticity.cpp, simple_interface.cpp, test_cache_on_entities.cpp, test_jacobian_of_simple_contact_element.cpp, testing_jacobian_of_hook_element.cpp, testing_jacobian_of_hook_scaled_with_density_element.cpp, thermo_plastic.cpp, unsaturated_transport.cpp, and wave_equation.cpp.

Definition at line 59 of file DMMMoFEM.cpp.

59  {
61  CHKERR DMRegister(sname, DMCreate_MoFEM);
63 }
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:87

◆ DMSetFromOptions_MoFEM()

PetscErrorCode MoFEM::DMSetFromOptions_MoFEM ( DM  dm)

Set options for MoFEM DM

Definition at line 1178 of file DMMMoFEM.cpp.

1178  {
1179 #endif
1180 
1181  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1183  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1184 #if PETSC_VERSION_GE(3, 5, 3)
1185  ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1186  CHKERRG(ierr);
1187 #else
1188  ierr = PetscOptionsHead("DMMoFEM Options");
1189  CHKERRG(ierr);
1190 #endif
1191  ierr = PetscOptionsBool("-dm_is_partitioned",
1192  "set if mesh is partitioned (works which native MOAB "
1193  "file format, i.e. h5m",
1194  "DMSetUp", dm_field->isPartitioned,
1195  &dm_field->isPartitioned, NULL);
1196  CHKERRG(ierr);
1198 }

◆ DMSetOperators_MoFEM()

PetscErrorCode MoFEM::DMSetOperators_MoFEM ( DM  dm)

Set operators for MoFEM dm.

Parameters
dm
Returns
error code

Definition at line 65 of file DMMMoFEM.cpp.

65  {
66  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
68 
69  dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
70  dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
71  dm->ops->creatematrix = DMCreateMatrix_MoFEM;
72  dm->ops->setup = DMSetUp_MoFEM;
73  dm->ops->destroy = DMDestroy_MoFEM;
74  dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
75  dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
76  dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
77  dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
78  dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
79  dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
80 
81  // Default matrix type
82  CHKERR DMSetMatType(dm, MATMPIAIJ);
83 
85 }
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1294
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1302
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMMoFEM.cpp:1122
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:95
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1331
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1178
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1102
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1366
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1200
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1360

◆ DMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSetUp_MoFEM ( DM  dm)

sets up the MoFEM structures inside a DM object

Examples
mixed_poisson.cpp.

Definition at line 1200 of file DMMMoFEM.cpp.

1200  {
1201  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1202  ProblemsManager *prb_mng_ptr;
1204  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1205  CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1206 
1207  if (dm_field->isCompDM) {
1208  // It is composite probelm
1209  CHKERR prb_mng_ptr->buildCompsedProblem(
1210  dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1211  dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1212  } else {
1213  if (dm_field->isPartitioned) {
1214  CHKERR prb_mng_ptr->buildProblemOnDistributedMesh(
1215  dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1216  dm_field->verbosity);
1217  } else {
1218  CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1219  dm_field->isSquareMatrix == PETSC_TRUE,
1220  dm_field->verbosity);
1221  CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1222  dm_field->verbosity);
1223  }
1224  }
1225 
1226  // Partition finite elements
1227  if (dm_field->isPartitioned) {
1228  CHKERR prb_mng_ptr->partitionFiniteElements(
1229  dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1230  CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1231  dm_field->problemName, dm_field->verbosity);
1232  } else {
1233  // partition finite elemnets
1234  CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1235  -1, -1, dm_field->verbosity);
1236  // Get ghost DOFs
1237  CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1238  dm_field->verbosity);
1239  }
1240 
1241  // Set flag that problem is build and partitioned
1242  dm_field->isProblemBuild = PETSC_TRUE;
1243 
1245 }

◆ DMSubDMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSubDMSetUp_MoFEM ( DM  subdm)

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

Definition at line 1247 of file DMMMoFEM.cpp.

1247  {
1248  PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1249  ProblemsManager *prb_mng_ptr;
1251 
1252  DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1253 
1254  // build sub dm problem
1255  CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1256 
1257  map<std::string, std::pair<EntityType, EntityType>> *entity_map_row = nullptr;
1258  map<std::string, std::pair<EntityType, EntityType>> *entity_map_col = nullptr;
1259 
1260  if (subdm_field->mapTypeRow)
1261  entity_map_row = subdm_field->mapTypeRow.get();
1262  if (subdm_field->mapTypeCol)
1263  entity_map_row = subdm_field->mapTypeCol.get();
1264 
1265  CHKERR prb_mng_ptr->buildSubProblem(
1266  subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1267  subdm_field->problemMainOfSubPtr->getName(),
1268  subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1269  subdm_field->verbosity);
1270 
1271  // partition problem
1272  subdm_field->isPartitioned = subdm_field->isPartitioned;
1273  if (subdm_field->isPartitioned) {
1274  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1275  0, subdm_field->sIze,
1276  subdm_field->verbosity);
1277  // set ghost nodes
1278  CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1279  subdm_field->problemName, subdm_field->verbosity);
1280  } else {
1281  // partition finite elements
1282  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1283  -1, -1, subdm_field->verbosity);
1284  // set ghost nodes
1285  CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1286  subdm_field->verbosity);
1287  }
1288 
1289  subdm_field->isProblemBuild = PETSC_TRUE;
1290 
1292 }

Variable Documentation

◆ smartCreateDMMatrix

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

Get smart matrix from DM.

Examples
mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, simple_contact.cpp, simple_contact_thermal.cpp, and test_jacobian_of_simple_contact_element.cpp.

Definition at line 944 of file DMMoFEM.hpp.

◆ smartCreateDMVector

auto MoFEM::smartCreateDMVector
Initial value:
= [](DM dm) {
SmartPetscObj<Vec> v;
CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
return v;
}
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj< Vec > &g_ptr)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1112
double v
phase velocity of light in medium (cm/ns)

Get smart vector from DM.

Examples
PlasticOpsMonitor.hpp, approx_sphere.cpp, contact.cpp, eigen_elastic.cpp, hcurl_check_approx_in_2d.cpp, heat_equation.cpp, heat_method.cpp, helmholtz.cpp, mixed_poisson.cpp, mortar_contact.cpp, mortar_contact_thermal.cpp, navier_stokes.cpp, nonlinear_elastic.cpp, photon_diffusion.cpp, plastic.cpp, scalar_check_approximation.cpp, shallow_wave.cpp, simple_contact.cpp, simple_contact_thermal.cpp, test_jacobian_of_simple_contact_element.cpp, thermo_plastic.cpp, and wave_equation.cpp.

Definition at line 956 of file DMMoFEM.hpp.