v0.8.19
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::DMMoFEMResolveSharedEntities (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::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::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[])
 
PetscErrorCode MoFEM::DMMoFEMAddSubFieldCol (DM dm, const char field_name[])
 
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:475
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ DMCreate_MoFEM()

PetscErrorCode MoFEM::DMCreate_MoFEM ( DM  dm)

Create dm data structure with MoFEM data structure.

Definition at line 120 of file DMMMoFEM.cpp.

120  {
121  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
123  dm->data = new DMCtx();
126 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:98
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

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

1112  {
1113  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1115 
1116  if (numFields) {
1117  *numFields = 0;
1118  }
1119  if (fieldNames) {
1120  *fieldNames = NULL;
1121  }
1122  if (fields) {
1123  *fields = NULL;
1124  }
1125 
1126  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1127  const Field_multiIndex *fields_ptr;
1128  CHKERR dm_field->mField_ptr->get_fields(&fields_ptr);
1129  Field_multiIndex::iterator fit, hi_fit;
1130  fit = fields_ptr->begin();
1131  hi_fit = fields_ptr->end();
1132  *numFields = std::distance(fit, hi_fit);
1133 
1134  if (fieldNames) {
1135  CHKERR PetscMalloc1(*numFields, fieldNames);
1136  }
1137  if (fields) {
1138  CHKERR PetscMalloc1(*numFields, fields);
1139  }
1140 
1141  for (int f = 0; fit != hi_fit; fit++, f++) {
1142  if (fieldNames) {
1143  CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1144  (char **)&(*fieldNames)[f]);
1145  }
1146  if (fields) {
1147  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1148  ->isCreateProblemFieldAndRank(
1149  dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1150  fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1151  }
1152  }
1153 
1155 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
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:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ DMCreateGlobalVector_MoFEM()

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.

Definition at line 890 of file DMMMoFEM.cpp.

890  {
891  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
893  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
894  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
895  dm_field->problemName, COL, g);
897 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
#define CHKERR
Inline error check.
Definition: definitions.h:594

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

899  {
900  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
902  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
903  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
904  dm_field->problemName, COL, l);
906 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ 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:475
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:908
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ DMCreateMatrix_MoFEM()

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.

Definition at line 908 of file DMMMoFEM.cpp.

908  {
909  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
911  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
912  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
913  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
914  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
915  M);
916  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
917  CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
918  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
919  M);
920  } else {
921  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
922  "Matrix type not implemented");
923  }
925 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ DMDestroy_MoFEM()

PetscErrorCode MoFEM::DMDestroy_MoFEM ( DM  dm)

Destroys dm with MoFEM data structure.

destroy the MoFEM structure

Definition at line 128 of file DMMMoFEM.cpp.

128  {
129  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
130  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
132  if (((DMCtx *)dm->data)->referenceNumber == 0) {
133  if (dm_field->destroyProblem) {
134  if (dm_field->mField_ptr->check_problem(dm_field->problemName)) {
135  dm_field->mField_ptr->delete_problem(dm_field->problemName);
136  } // else problem has to be deleted by the user
137  }
138  // cerr << "Destroy " << dm_field->problemName << endl;
139  delete ((DMCtx *)dm->data);
140  } else {
141  // cerr << "Dereference " << dm_field->problemName << " " <<
142  // ((DMCtx*)dm->data)->referenceNumber << endl;
143  (((DMCtx *)dm->data)->referenceNumber)--;
144  }
146 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

1040  {
1041  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1043  CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1045 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

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

1047  {
1049  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1051 
1052  CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1053 
1054  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1055  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1056  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1057 
1058  double *array_loc, *array_glob;
1059  CHKERR VecGetArray(l, &array_loc);
1060  CHKERR VecGetArray(g, &array_glob);
1061  switch (mode) {
1062  case INSERT_VALUES:
1063  cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1064  break;
1065  case ADD_VALUES:
1066  cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1067  break;
1068  default:
1069  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1070  }
1071  CHKERR VecRestoreArray(l, &array_loc);
1072  CHKERR VecRestoreArray(g, &array_glob);
1074 }
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:499
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:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

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

1077  {
1078 
1079  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1081 
1082  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1083  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1084  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1085 
1086  double *array_loc, *array_glob;
1087  CHKERR VecGetArray(l, &array_loc);
1088  CHKERR VecGetArray(g, &array_glob);
1089  switch (mode) {
1090  case INSERT_VALUES:
1091  cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1092  break;
1093  case ADD_VALUES:
1094  cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1095  break;
1096  default:
1097  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1098  }
1099  CHKERR VecRestoreArray(l, &array_loc);
1100  CHKERR VecRestoreArray(g, &array_glob);
1101 
1103 }
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:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

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

1105  {
1106  //
1109 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ 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:475
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ 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:475
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ 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:475
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:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
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 319 of file DMMMoFEM.cpp.

319  {
320  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
322  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
323  if (!dm->data) {
324  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
325  "data structure for MoFEM not yet created");
326  }
327  if (!dm_field->isCompDM) {
328  dm_field->isCompDM = PETSC_TRUE;
329  }
330  if (dm_field->isSquareMatrix) {
331  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
332  "No need to add problem on column when problem block structurally "
333  "symmetric");
334  }
335  dm_field->colCompPrb.push_back(prb_name);
337 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ 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, elasticity.cpp, elasticity_mixed_formulation.cpp, MagneticElement.hpp, minimal_surface_area.cpp, Remodeling.cpp, simple_elasticity.cpp, and UnsaturatedFlow.hpp.

Definition at line 433 of file DMMMoFEM.cpp.

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

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

301  {
302  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
304  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
305  if (!dm->data) {
306  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
307  "data structure for MoFEM not yet created");
308  }
309  if (!dm_field->isCompDM) {
310  dm_field->isCompDM = PETSC_TRUE;
311  }
312  dm_field->rowCompPrb.push_back(prb_name);
313  if (dm_field->isSquareMatrix) {
314  dm_field->colCompPrb.push_back(prb_name);
315  }
317 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ DMMoFEMAddSubFieldCol()

PetscErrorCode MoFEM::DMMoFEMAddSubFieldCol ( DM  dm,
const char  field_name[] 
)

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

241  {
242  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
244  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
245  if (!dm->data) {
246  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
247  "data structure for MoFEM not yet created");
248  }
249  if (!dm_field->isSubDM) {
250  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
251  }
252  dm_field->colFields.push_back(field_name);
254 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMAddSubFieldRow()

PetscErrorCode MoFEM::DMMoFEMAddSubFieldRow ( DM  dm,
const char  field_name[] 
)

Add field to sub dm problem on rows

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

Definition at line 226 of file DMMMoFEM.cpp.

226  {
227  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
229  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
230  if (!dm->data) {
231  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
232  "data structure for MoFEM not yet created");
233  }
234  if (!dm_field->isSubDM) {
235  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
236  }
237  dm_field->rowFields.push_back(field_name);
239 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ 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, elasticity.cpp, elasticity_mixed_formulation.cpp, MagneticElement.hpp, minimal_surface_area.cpp, Remodeling.cpp, and UnsaturatedFlow.hpp.

Definition at line 148 of file DMMMoFEM.cpp.

151  {
153 
154  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
155  if (!dm->data) {
156  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
157  "data structure for MoFEM not yet created");
158  }
159  if (!m_field_ptr) {
160  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
161  "DM function not implemented into MoFEM");
162  }
163  dm_field->mField_ptr = m_field_ptr;
164  dm_field->problemName = problem_name;
165  if (!m_field_ptr->check_problem(dm_field->problemName)) {
166  // problem is not defined, declare problem internally set bool to
167  // destroyProblem problem with DM
168  dm_field->destroyProblem = PETSC_TRUE;
169  CHKERR dm_field->mField_ptr->add_problem(dm_field->problemName, MF_EXCL,
170  dm_field->verbosity);
171  } else {
172  dm_field->destroyProblem = PETSC_FALSE;
173  }
174  CHKERR dm_field->mField_ptr->modify_problem_ref_level_add_bit(
175  dm_field->problemName, bit_level);
176  CHKERR dm_field->mField_ptr->modify_problem_mask_ref_level_add_bit(
177  dm_field->problemName, bit_mask);
178  dm_field->kspCtx =
179  boost::shared_ptr<KspCtx>(new KspCtx(*m_field_ptr, problem_name));
180  dm_field->snesCtx =
181  boost::shared_ptr<SnesCtx>(new SnesCtx(*m_field_ptr, problem_name));
182  dm_field->tsCtx =
183  boost::shared_ptr<TsCtx>(new TsCtx(*m_field_ptr, problem_name));
184 
185  MPI_Comm comm;
186  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
187  int result = 0;
188  MPI_Comm_compare(comm, m_field_ptr->get_comm(), &result);
189  // std::cerr << result << " " << MPI_IDENT << " " << MPI_CONGRUENT << " " <<
190  // MPI_SIMILAR << " " << MPI_UNEQUAL << std::endl;
191  if (result > MPI_CONGRUENT) {
192  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
193  "MoFEM and DM using different communicators");
194  }
195  MPI_Comm_size(comm, &dm_field->sIze);
196  MPI_Comm_rank(comm, &dm_field->rAnk);
197 
198  // problem structure
199  CHKERR dm_field->mField_ptr->get_problem(dm_field->problemName,
200  &dm_field->problemPtr);
201 
203 }
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:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
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, and dm_create_subdm.cpp.

Definition at line 205 of file DMMMoFEM.cpp.

205  {
207 
208  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
209  if (!dm->data) {
210  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
211  "data structure for MoFEM not yet created");
212  }
213  CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
214  dm_field->problemPtr->getBitRefLevel());
215 
216  DMCtx *subdm_field = (DMCtx *)subdm->data;
217  subdm_field->isSubDM = PETSC_TRUE;
218  subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
219  subdm_field->isPartitioned = dm_field->isPartitioned;
220  subdm_field->isSquareMatrix = PETSC_FALSE;
221  subdm->ops->setup = DMSubDMSetUp_MoFEM;
222 
224 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1002
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:148
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

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

1158  {
1159  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1161  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1162  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1163  ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1164  field_name, 0, 1000, is);
1166 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

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

339  {
341  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
343  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
344  *is_comp_dm = dm_field->isCompDM;
346 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMGetIsPartitioned()

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

get if read mesh is partitioned

Definition at line 857 of file DMMMoFEM.cpp.

857  {
858  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
860  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
861  *is_partitioned = dm_field->isPartitioned;
863 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

256  {
258  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
260  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
261  *is_sub_dm = dm_field->isSubDM;
263 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMGetKspCtx() [1/2]

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

get MoFEM::KspCtx data structure

Definition at line 791 of file DMMMoFEM.cpp.

791  {
792  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
794  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
795  *ksp_ctx = dm_field->kspCtx.get();
797 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMGetKspCtx() [2/2]

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

get MoFEM::KspCtx data structure

Definition at line 800 of file DMMMoFEM.cpp.

800  {
801  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
803  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
804  const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
806 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

410  {
412  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
414  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
415 
416  MPI_Comm comm;
417  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
418  CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
419  layout);
421 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ DMMoFEMGetProblemPtr()

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

Get pointer to problem data structure.

Examples:
MagneticElement.hpp, and Remodeling.cpp.

Definition at line 360 of file DMMMoFEM.cpp.

360  {
361  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
363  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
364  if (!dm->data) {
365  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
366  "data structure for MoFEM not yet created");
367  }
368  *problem_ptr = dm_field->problemPtr;
370 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMGetSnesCtx() [1/2]

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

get MoFEM::SnesCtx data structure

Examples:
minimal_surface_area.cpp, and testing_jacobian_of_hook_element.cpp.

Definition at line 817 of file DMMMoFEM.cpp.

817  {
818  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
820  DMCtx *dm_field = (DMCtx *)dm->data;
821  *snes_ctx = dm_field->snesCtx.get();
823 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMGetSnesCtx() [2/2]

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

get MoFEM::SnesCtx data structure

Definition at line 826 of file DMMMoFEM.cpp.

826  {
827  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
829  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
830  const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
832 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

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

◆ DMMoFEMGetTsCtx() [1/2]

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

get MoFEM::TsCtx data structure

Examples:
Remodeling.cpp, and UnsaturatedFlow.hpp.

Definition at line 865 of file DMMMoFEM.cpp.

865  {
866  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
868  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
869  *ts_ctx = dm_field->tsCtx.get();
871 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMGetTsCtx() [2/2]

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

get MoFEM::TsCtx data structure

Definition at line 873 of file DMMMoFEM.cpp.

874  {
875  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
877  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
878  const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
880 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

608  {
609  return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
612  dm, fe_name, method, pre_only, post_only);
613 }
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:586
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 616 of file DMMMoFEM.cpp.

619  {
620  return DMMoFEMKSPSetComputeOperators<const std::string &,
621  boost::shared_ptr<MoFEM::FEMethod>>(
622  dm, fe_name, method, pre_only, post_only);
623 }
static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:586

◆ 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 code set KSP right hand side evaluation function
Examples:
analytical_poisson.cpp, and simple_elasticity.cpp.

Definition at line 564 of file DMMMoFEM.cpp.

567  {
568  return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
570  dm, fe_name, method, pre_only, post_only);
571 }
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:545
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 574 of file DMMMoFEM.cpp.

577  {
578  return DMMoFEMKSPSetComputeRHS<const std::string &,
579  boost::shared_ptr<MoFEM::FEMethod>,
580  boost::shared_ptr<MoFEM::BasicMethod>,
581  boost::shared_ptr<MoFEM::BasicMethod>>(
582  dm, fe_name, method, pre_only, post_only);
583 }
static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:545

◆ DMMoFEMResolveSharedEntities()

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

Resolve shared entities.

Parameters
dmdm
fe_namefinite element for which shared entities are resolved
Returns
error code

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

ierr = DMMoFEMGetSquareProblem(dm,"SHELL_ELEMENT"); CHKERRG(ierr);
Tag th;
rval = mField.get_moab().tag_get_handle("ADAPT_ORDER",th); CHKERRQ_MOAB(rval);
ParallelComm* pcomm =
ParallelComm::get_pcomm(&mField.get_moab(),MYPCOMM_INDEX);
// rval = pcomm->reduce_tags(th,MPI_SUM,prisms);
rval = pcomm->exchange_tags(th,prisms);

Definition at line 399 of file DMMMoFEM.cpp.

399  {
401  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
403  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
404  CHKERR dm_field->mField_ptr->resolve_shared_ents(dm_field->problemPtr,
405  fe_name);
407 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

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

Definition at line 846 of file DMMMoFEM.cpp.

846  {
847  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
849  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
850  dm_field->isPartitioned = is_partitioned;
852 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMSetKspCtx()

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

set MoFEM::KspCtx data structure

Definition at line 808 of file DMMMoFEM.cpp.

809  {
810  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
812  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
813  dm_field->kspCtx = ksp_ctx;
815 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMMoFEMSetSnesCtx()

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

Set MoFEM::SnesCtx data structure.

Definition at line 834 of file DMMMoFEM.cpp.

835  {
836  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
838  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
839  dm_field->snesCtx = snes_ctx;
841 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

Definition at line 390 of file DMMMoFEM.cpp.

390  {
392  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
394  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
395  dm_field->isSquareMatrix = square_problem;
397 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

882  {
883  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
885  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
886  dm_field->tsCtx = ts_ctx;
888 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

Definition at line 645 of file DMMMoFEM.cpp.

648  {
649  return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
651  dm, fe_name, method, pre_only, post_only);
652 }
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:626
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 655 of file DMMMoFEM.cpp.

658  {
659  return DMMoFEMSNESSetFunction<const std::string &,
660  boost::shared_ptr<MoFEM::FEMethod>,
661  boost::shared_ptr<MoFEM::BasicMethod>,
662  boost::shared_ptr<MoFEM::BasicMethod>>(
663  dm, fe_name, method, pre_only, post_only);
664 }
static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:626

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

Definition at line 686 of file DMMMoFEM.cpp.

689  {
690  return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
692  dm, fe_name, method, pre_only, post_only);
693 }
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:667
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 696 of file DMMMoFEM.cpp.

699  {
700  return DMMoFEMSNESSetJacobian<const std::string &,
701  boost::shared_ptr<MoFEM::FEMethod>,
702  boost::shared_ptr<MoFEM::BasicMethod>,
703  boost::shared_ptr<MoFEM::BasicMethod>>(
704  dm, fe_name, method, pre_only, post_only);
705 }
static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:667

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

Definition at line 727 of file DMMMoFEM.cpp.

730  {
731  return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
733  dm, fe_name, method, pre_only, post_only);
735 }
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:506
static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:708
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 738 of file DMMMoFEM.cpp.

741  {
742  return DMMoFEMTSSetIFunction<const std::string,
743  boost::shared_ptr<MoFEM::FEMethod>,
744  boost::shared_ptr<MoFEM::BasicMethod>,
745  boost::shared_ptr<MoFEM::BasicMethod>>(
746  dm, fe_name, method, pre_only, post_only);
748 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:708

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

Definition at line 780 of file DMMMoFEM.cpp.

783  {
784  return DMMoFEMTSSetIJacobian<const std::string &,
785  boost::shared_ptr<MoFEM::FEMethod>,
786  boost::shared_ptr<MoFEM::BasicMethod>,
787  boost::shared_ptr<MoFEM::BasicMethod>>(
788  dm, fe_name, method, pre_only, post_only);
789 }
static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:751

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

773  {
774  return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
775  MoFEM::BasicMethod *>(dm, fe_name, method,
776  pre_only, post_only);
777 }
static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method, T1 pre_only, T2 post_only)
Definition: DMMMoFEM.cpp:751
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 443 of file DMMMoFEM.cpp.

443  {
444  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
446  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
447  ierr = dm_field->mField_ptr->modify_problem_unset_finite_element(
448  dm_field->problemName, fe_name);
449  CHKERRG(ierr);
451 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

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

◆ DMoFEMLoopDofs()

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

execute method for dofs on field in problem

Definition at line 532 of file DMMMoFEM.cpp.

533  {
534  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
536  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
537  ierr =
538  dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
539  *method, dm_field->rAnk, dm_field->rAnk);
540  CHKERRG(ierr);
542 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ 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, cell_forces.cpp, elasticity.cpp, elasticity_mixed_formulation.cpp, hello_world.cpp, HookeElement.cpp, MagneticElement.hpp, mesh_smoothing.cpp, minimal_surface_area.cpp, Remodeling.cpp, simple_elasticity.cpp, simple_interface.cpp, and UnsaturatedFlow.hpp.

Definition at line 515 of file DMMMoFEM.cpp.

516  {
517  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
519  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
520  ierr = DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name, method,
521  dm_field->rAnk, dm_field->rAnk);
522  CHKERRG(ierr);
524 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
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:496

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

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

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

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

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

510  {
511  return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
512  low_rank, up_rank);
513 }
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:496

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

466  {
467  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
469  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
470  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
471  dm_field->problemPtr, COL, g, mode, scatter_mode);
472  CHKERRG(ierr);
474 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

Definition at line 453 of file DMMMoFEM.cpp.

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

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

486  {
487  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
489  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
490  ierr = dm_field->mField_ptr->problem_basic_method_postProcess(
491  dm_field->problemPtr, *method);
492  CHKERRG(ierr);
494 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

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

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

◆ 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:475
MoFEMErrorCode DMCreate_MGViaApproxOrders(DM dm)
Create DM data structure for Multi-Grid via approximation orders.
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ 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, cell_forces.cpp, dm_build_partitioned_mesh.cpp, dm_create_subdm.cpp, elasticity_mixed_formulation.cpp, hello_world.cpp, MagneticElement.hpp, mesh_smoothing.cpp, minimal_surface_area.cpp, Remodeling.cpp, simple_elasticity.cpp, simple_interface.cpp, testing_jacobian_of_hook_element.cpp, and unsaturated_transport.cpp.

Definition at line 92 of file DMMMoFEM.cpp.

92  {
94  CHKERR DMRegister(sname, DMCreate_MoFEM);
96 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:120
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ DMSetFromOptions_MoFEM()

PetscErrorCode MoFEM::DMSetFromOptions_MoFEM ( DM  dm)

Set options for MoFEM DM

Definition at line 933 of file DMMMoFEM.cpp.

933  {
934 #endif
935 
936  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
938  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
939 #if PETSC_VERSION_GE(3, 5, 3)
940  ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
941  CHKERRG(ierr);
942 #else
943  ierr = PetscOptionsHead("DMMoFEM Options");
944  CHKERRG(ierr);
945 #endif
946  ierr = PetscOptionsBool("-dm_is_partitioned",
947  "set if mesh is partitioned (works which native MOAB "
948  "file format, i.e. h5m",
949  "DMSetUp", dm_field->isPartitioned,
950  &dm_field->isPartitioned, NULL);
951  CHKERRG(ierr);
953 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506

◆ DMSetOperators_MoFEM()

PetscErrorCode MoFEM::DMSetOperators_MoFEM ( DM  dm)

Set operators for MoFEM dm.

Parameters
dm
Returns
error code

Definition at line 98 of file DMMMoFEM.cpp.

98  {
99  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
101 
102  dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
103  dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
104  dm->ops->creatematrix = DMCreateMatrix_MoFEM;
105  dm->ops->setup = DMSetUp_MoFEM;
106  dm->ops->destroy = DMDestroy_MoFEM;
107  dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
108  dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
109  dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
110  dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
111  dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
112  dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
113 
114  // Default matrix type
115  CHKERR DMSetMatType(dm, MATMPIAIJ);
116 
118 }
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1039
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:955
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1076
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:908
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:128
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM...
Definition: DMMMoFEM.cpp:890
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1047
#define CHKERR
Inline error check.
Definition: definitions.h:594
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1111
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1105
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:933
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM...
Definition: DMMMoFEM.cpp:899

◆ DMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSetUp_MoFEM ( DM  dm)

sets up the MoFEM structures inside a DM object

Definition at line 955 of file DMMMoFEM.cpp.

955  {
956  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
957  ProblemsManager *prb_mng_ptr;
959  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
960  CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
961 
962  if (dm_field->isCompDM) {
963  // It is composite probelm
964  CHKERR prb_mng_ptr->buildCompsedProblem(
965  dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
966  dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
967  } else {
968  if (dm_field->isPartitioned) {
969  CHKERR prb_mng_ptr->buildProblemOnDistributedMesh(
970  dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
971  dm_field->verbosity);
972  } else {
973  CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
974  dm_field->isSquareMatrix == PETSC_TRUE,
975  dm_field->verbosity);
976  CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
977  dm_field->verbosity);
978  }
979  }
980 
981  // Partition finite elements
982  if (dm_field->isPartitioned) {
983  CHKERR prb_mng_ptr->partitionFiniteElements(
984  dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
985  CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
986  dm_field->problemName, dm_field->verbosity);
987  } else {
988  // partition finite elemnets
989  CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
990  -1, -1, dm_field->verbosity);
991  // Get ghost DOFs
992  CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
993  dm_field->verbosity);
994  }
995 
996  // Set flag that problem is build and partitioned
997  dm_field->isProblemBuild = PETSC_TRUE;
998 
1000 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ DMSubDMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSubDMSetUp_MoFEM ( DM  subdm)

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

Definition at line 1002 of file DMMMoFEM.cpp.

1002  {
1003  PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1004  ProblemsManager *prb_mng_ptr;
1006 
1007  DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1008 
1009  // build sub dm problem
1010  CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1011  CHKERR prb_mng_ptr->buildSubProblem(
1012  subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1013  subdm_field->problemMainOfSubPtr->getName(),
1014  subdm_field->isSquareMatrix == PETSC_TRUE, subdm_field->verbosity);
1015 
1016  // partition problem
1017  subdm_field->isPartitioned = subdm_field->isPartitioned;
1018  if (subdm_field->isPartitioned) {
1019  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1020  0, subdm_field->sIze,
1021  subdm_field->verbosity);
1022  // set ghost nodes
1023  CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1024  subdm_field->problemName, subdm_field->verbosity);
1025  } else {
1026  // partition finite elements
1027  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1028  -1, -1, subdm_field->verbosity);
1029  // set ghost nodes
1030  CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1031  subdm_field->verbosity);
1032  }
1033 
1034  subdm_field->isProblemBuild = PETSC_TRUE;
1035 
1037 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405