v0.13.1
Loading...
Searching...
No Matches
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 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...
 
auto MoFEM::smartCreateDMMatrix (DM dm)
 Get smart matrix from DM. More...
 
auto MoFEM::smartCreateDMVector (DM dm)
 Get smart vector from DM. 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 462 of file PCMGSetUpViaApproxOrders.cpp.

462 {
463 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
465 GET_DM_FIELD(dm);
466 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
467 CHKERR DMCreate(comm, dmc);
468 (*dmc)->data = dm->data;
469 DMType type;
470 CHKERR DMGetType(dm, &type);
471 CHKERR DMSetType(*dmc, type);
472 CHKERR PetscObjectReference((PetscObject)(*dmc));
473 PetscInfo1(dm, "Coarsen DMMGViaApproxOrders leveldown = %d\n", dm->leveldown);
475}
#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:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535

◆ DMCreate_MoFEM()

PetscErrorCode MoFEM::DMCreate_MoFEM ( DM  dm)

Create dm data structure with MoFEM data structure.

Definition at line 75 of file DMMMoFEM.cpp.

75 {
76 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
78 dm->data = new DMCtx();
81}
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:53

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

1395 {
1396 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1398
1399 if (numFields) {
1400 *numFields = 0;
1401 }
1402 if (fieldNames) {
1403 *fieldNames = NULL;
1404 }
1405 if (fields) {
1406 *fields = NULL;
1407 }
1408
1409 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1410 auto fields_ptr = dm_field->mField_ptr->get_fields();
1411 Field_multiIndex::iterator fit, hi_fit;
1412 fit = fields_ptr->begin();
1413 hi_fit = fields_ptr->end();
1414 *numFields = std::distance(fit, hi_fit);
1415
1416 if (fieldNames) {
1417 CHKERR PetscMalloc1(*numFields, fieldNames);
1418 }
1419 if (fields) {
1420 CHKERR PetscMalloc1(*numFields, fields);
1421 }
1422
1423 for (int f = 0; fit != hi_fit; fit++, f++) {
1424 if (fieldNames) {
1425 CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1426 (char **)&(*fieldNames)[f]);
1427 }
1428 if (fields) {
1429 CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1430 ->isCreateProblemFieldAndRank(
1431 dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1432 fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1433 }
1434 }
1435
1437}
@ ROW
Definition: definitions.h:123
auto f
Definition: HenckyOps.hpp:5

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

1124 {
1125 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1127 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1128 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1129 dm_field->problemName, COL, g_ptr);
1130 CHKERR VecSetDM(g_ptr, dm);
1132}
@ COL
Definition: definitions.h:123
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440

◆ 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

Definition at line 1114 of file DMMMoFEM.cpp.

1114 {
1115 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1117 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1118 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1119 dm_field->problemName, COL, g);
1120 CHKERR VecSetDM(*g, dm);
1122}
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 1134 of file DMMMoFEM.cpp.

1134 {
1135 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1137 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1138 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1139 dm_field->problemName, COL, l);
1140 CHKERR VecSetDM(*l, dm);
1142}
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 430 of file PCMGSetUpViaApproxOrders.cpp.

430 {
431
432 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
434 GET_DM_FIELD(dm);
435
436 int leveldown = dm->leveldown;
437
438 if (dm_field->kspOperators.empty()) {
440 } else {
441 MPI_Comm comm;
442 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
443 if (dm_field->kspOperators.empty()) {
444 SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
445 "data inconsistency, operator can not be set");
446 }
447 if (static_cast<int>(dm_field->kspOperators.size()) < leveldown) {
448 SETERRQ(comm, MOFEM_DATA_INCONSISTENCY,
449 "data inconsistency, no IS for that level");
450 }
451 *M = dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown];
452 CHKERR PetscObjectReference((PetscObject)*M);
453 CHKERR MatSetDM(*M, dm);
454 }
455
456 PetscInfo1(dm, "Create Matrix DMMGViaApproxOrders leveldown = %d\n",
457 leveldown);
458
460}
static Index< 'M', 3 > M
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1144

◆ 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

Definition at line 1144 of file DMMMoFEM.cpp.

1144 {
1145 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1147 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1148 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1149 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1150 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1151 M);
1152 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1153 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1154 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1155 M);
1156 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1157 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1158 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1159 dm_field->problemName, M);
1160 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1161 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1162 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1163 dm_field->problemName, M);
1164 } else {
1165 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1166 "Matrix type not implemented");
1167 }
1168 CHKERR MatSetDM(*M, dm);
1170}
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32

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

1172 {
1173 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1175 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1176 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1177 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1178 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1179 M);
1180 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1181 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1182 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1183 M);
1184 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1185 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1186 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1187 dm_field->problemName, M);
1188 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1189 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1190 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1191 dm_field->problemName, M);
1192 } else {
1193 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1194 "Matrix type not implemented");
1195 }
1196 CHKERR MatSetDM(M, dm);
1198}

◆ DMDestroy_MoFEM()

PetscErrorCode MoFEM::DMDestroy_MoFEM ( DM  dm)

Destroys dm with MoFEM data structure.

destroy the MoFEM structure

Definition at line 83 of file DMMMoFEM.cpp.

83 {
84 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
85 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
87
88 MPI_Comm comm;
89 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
90
91 int result;
92 MPI_Comm_compare(comm, PETSC_COMM_SELF, &result);
93 if (result == MPI_IDENT)
94 MOFEM_LOG("DMSELF", Sev::noisy)
95 << "MoFEM DM destroy for problem " << dm_field->problemName
96 << " referenceNumber " << dm_field->referenceNumber;
97 else
98 MOFEM_LOG("DMWORLD", Sev::noisy)
99 << "MoFEM DM destroy for problem " << dm_field->problemName
100 << " referenceNumber " << dm_field->referenceNumber;
101
102 if (dm_field->referenceNumber == 0) {
103 if (dm_field->destroyProblem) {
104
105 if (dm_field->mField_ptr->check_problem(dm_field->problemName)) {
106 dm_field->mField_ptr->delete_problem(dm_field->problemName);
107 } // else problem has to be deleted by the user
108 }
109
110 delete static_cast<DMCtx *>(dm->data);
111
112 } else
113 --dm_field->referenceNumber;
114
116}
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301

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

1323 {
1324 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1326 CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1328}

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

1330 {
1332 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1334
1335 CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
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_glob, 1, array_loc, 1);
1347 break;
1348 case ADD_VALUES:
1349 cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 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);
1357}

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

1360 {
1361
1362 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1364
1365 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1366 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1367 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1368
1369 double *array_loc, *array_glob;
1370 CHKERR VecGetArray(l, &array_loc);
1371 CHKERR VecGetArray(g, &array_glob);
1372 switch (mode) {
1373 case INSERT_VALUES:
1374 cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1375 break;
1376 case ADD_VALUES:
1377 cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1378 break;
1379 default:
1380 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1381 }
1382 CHKERR VecRestoreArray(l, &array_loc);
1383 CHKERR VecRestoreArray(g, &array_glob);
1384
1386}

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

1388 {
1389 //
1392}

◆ DMMGViaApproxOrdersClearCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS ( DM  dm)

Clear approximation orders.

Parameters
DMdm
Returns
Error code

Definition at line 276 of file PCMGSetUpViaApproxOrders.cpp.

276 {
277 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
279 GET_DM_FIELD(dm);
280 CHKERR dm_field->destroyCoarseningIS();
281 PetscInfo(dm, "Clear DMs data structures\n");
283}

◆ DMMGViaApproxOrdersPopBackCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS ( DM  dm)

Pop is form MG via approximation orders.

Parameters
DMdm
ispop back IS
Returns
error code

Definition at line 260 of file PCMGSetUpViaApproxOrders.cpp.

260 {
261 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
263 GET_DM_FIELD(dm);
264 if (dm_field->coarseningIS.back()) {
265 CHKERR ISDestroy(&dm_field->coarseningIS.back());
266 dm_field->coarseningIS.pop_back();
267 }
268 if (dm_field->kspOperators.back()) {
269 CHKERR MatDestroy(&dm_field->kspOperators.back());
270 }
271 dm_field->kspOperators.pop_back();
272 PetscInfo(dm, "Pop back IS to DMMGViaApproxOrders\n");
274}

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

209 {
210 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
212 GET_DM_FIELD(dm);
213 dm_field->coarseningIS.push_back(is);
214 dm_field->shellMatrixCtxPtr.push_back(new PCMGSubMatrixCtx_private(A, is));
215 if (is) {
216 CHKERR PetscObjectReference((PetscObject)is);
217 }
218 if (is) {
219 IS is2 = is;
220 if (dm_field->aO) {
221 CHKERR ISDuplicate(is, &is2);
222 CHKERR ISCopy(is, is2);
223 CHKERR AOApplicationToPetscIS(dm_field->aO, is2);
224 }
225 if (create_sub_matrix) {
226 if (shell_sub_a) {
227 int n, N;
228 CHKERR ISGetSize(is, &N);
229 CHKERR ISGetLocalSize(is, &n);
230 MPI_Comm comm;
231 CHKERR PetscObjectGetComm((PetscObject)A, &comm);
232 CHKERR MatCreateShell(comm, n, n, N, N,
233 &(dm_field->shellMatrixCtxPtr.back()), subA);
234 CHKERR MatShellSetOperation(*subA, MATOP_MULT,
235 (void (*)(void))sub_mat_mult);
236 CHKERR MatShellSetOperation(*subA, MATOP_MULT_ADD,
237 (void (*)(void))sub_mat_mult_add);
238 CHKERR MatShellSetOperation(*subA, MATOP_SOR,
239 (void (*)(void))sub_mat_sor);
240 } else {
241#if PETSC_VERSION_GE(3, 8, 0)
242 CHKERR MatCreateSubMatrix(A, is2, is2, MAT_INITIAL_MATRIX, subA);
243#else
244 CHKERR MatGetSubMatrix(A, is2, is2, MAT_INITIAL_MATRIX, subA);
245#endif
246 }
247 }
248 if (dm_field->aO) {
249 CHKERR ISDestroy(&is2);
250 }
251 dm_field->kspOperators.push_back(*subA);
252 CHKERR PetscObjectReference((PetscObject)(*subA));
253 } else {
254 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
255 }
256 PetscInfo(dm, "Push back IS to DMMGViaApproxOrders\n");
258}
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)
FTensor::Index< 'n', SPACE_DIM > n
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 332 of file DMMMoFEM.cpp.

332 {
333 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
335 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
336 if (!dm->data) {
337 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
338 "data structure for MoFEM not yet created");
339 }
340 if (!dm_field->isCompDM) {
341 dm_field->isCompDM = PETSC_TRUE;
342 }
343 if (dm_field->isSquareMatrix) {
344 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
345 "No need to add problem on column when problem block structurally "
346 "symmetric");
347 }
348 dm_field->colCompPrb.push_back(prb_name);
350}
@ MOFEM_INVALID_DATA
Definition: definitions.h:36

◆ 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.

Definition at line 450 of file DMMMoFEM.cpp.

450 {
451 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
453 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
454 ierr = dm_field->mField_ptr->modify_problem_add_finite_element(
455 dm_field->problemName, fe_name);
456 CHKERRG(ierr);
458}
static PetscErrorCode ierr
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483

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

314 {
315 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
317 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
318 if (!dm->data) {
319 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
320 "data structure for MoFEM not yet created");
321 }
322 if (!dm_field->isCompDM) {
323 dm_field->isCompDM = PETSC_TRUE;
324 }
325 dm_field->rowCompPrb.push_back(prb_name);
326 if (dm_field->isSquareMatrix) {
327 dm_field->colCompPrb.push_back(prb_name);
328 }
330}

◆ 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

Definition at line 246 of file DMMMoFEM.cpp.

247 {
248 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
250 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
251 if (!dm->data) {
252 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
253 "data structure for MoFEM not yet created");
254 }
255 if (!dm_field->isSubDM) {
256 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
257 }
258 dm_field->colFields.push_back(field_name);
259 if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
260 if (!dm_field->mapTypeCol)
261 dm_field->mapTypeCol = boost::make_shared<
262 std::map<std::string, std::pair<EntityType, EntityType>>>();
263 (*dm_field->mapTypeCol)[field_name] =
264 std::pair<EntityType, EntityType>(lo_type, hi_type);
265 }
267}
constexpr auto field_name

◆ 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

Definition at line 223 of file DMMMoFEM.cpp.

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

◆ 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.

Definition at line 118 of file DMMMoFEM.cpp.

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

Definition at line 201 of file DMMMoFEM.cpp.

201 {
203
204 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
205 if (!dm->data) {
206 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
207 "data structure for MoFEM not yet created");
208 }
209 CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
210 dm_field->problemPtr->getBitRefLevel(),
211 dm_field->problemPtr->getBitRefLevelMask());
212
213 DMCtx *subdm_field = (DMCtx *)subdm->data;
214 subdm_field->isSubDM = PETSC_TRUE;
215 subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
216 subdm_field->isPartitioned = dm_field->isPartitioned;
217 subdm_field->isSquareMatrix = PETSC_FALSE;
218 subdm->ops->setup = DMSubDMSetUp_MoFEM;
219
221}
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:118

◆ 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:1439
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76

Definition at line 1439 of file DMMMoFEM.cpp.

1440 {
1441 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1443 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1444 CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1445 ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1446 field_name, 0, 1000, is);
1448}

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

352 {
354 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
356 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
357 *is_comp_dm = dm_field->isCompDM;
359}

◆ DMMoFEMGetIsPartitioned()

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

get if read mesh is partitioned

Definition at line 1081 of file DMMMoFEM.cpp.

1081 {
1082 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1084 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1085 *is_partitioned = dm_field->isPartitioned;
1087}

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

269 {
271 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
273 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
274 *is_sub_dm = dm_field->isSubDM;
276}

◆ DMMoFEMGetKspCtx() [1/2]

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

get MoFEM::KspCtx data structure

Definition at line 1024 of file DMMMoFEM.cpp.

1024 {
1025 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1027 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1028 const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
1030}

◆ DMMoFEMGetKspCtx() [2/2]

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

get MoFEM::KspCtx data structure

Definition at line 1015 of file DMMMoFEM.cpp.

1015 {
1016 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1018 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1019 *ksp_ctx = dm_field->kspCtx.get();
1021}

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

427 {
429 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
431 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
432
433 MPI_Comm comm;
434 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
435 CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
436 layout);
438}

◆ DMMoFEMGetProblemPtr()

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

Get pointer to problem data structure.

Definition at line 373 of file DMMMoFEM.cpp.

373 {
374 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
376 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
377 if (!dm->data) {
378 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
379 "data structure for MoFEM not yet created");
380 }
381 *problem_ptr = dm_field->problemPtr;
383}

◆ DMMoFEMGetSnesCtx() [1/2]

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

get MoFEM::SnesCtx data structure

Definition at line 1050 of file DMMMoFEM.cpp.

1050 {
1051 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1053 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1054 const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
1056}

◆ DMMoFEMGetSnesCtx() [2/2]

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

get MoFEM::SnesCtx data structure

Definition at line 1041 of file DMMMoFEM.cpp.

1041 {
1042 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1044 DMCtx *dm_field = (DMCtx *)dm->data;
1045 *snes_ctx = dm_field->snesCtx.get();
1047}

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

440 {
443 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
445 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
446 *square_problem = dm_field->isSquareMatrix;
448}

◆ DMMoFEMGetTsCtx() [1/2]

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

get MoFEM::TsCtx data structure

Definition at line 1097 of file DMMMoFEM.cpp.

1098 {
1099 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1101 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1102 const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1104}

◆ DMMoFEMGetTsCtx() [2/2]

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

get MoFEM::TsCtx data structure

Definition at line 1089 of file DMMMoFEM.cpp.

1089 {
1090 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1092 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1093 *ts_ctx = dm_field->tsCtx.get();
1095}

◆ 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

Definition at line 625 of file DMMMoFEM.cpp.

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

639 {
640 return DMMoFEMKSPSetComputeOperators<const std::string,
641 boost::shared_ptr<MoFEM::FEMethod>>(
642 dm, fe_name, method, pre_only, post_only);
643}

◆ 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

Definition at line 584 of file DMMMoFEM.cpp.

587 {
588 return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
590 dm, fe_name, method, pre_only, post_only);
591}
PetscErrorCode 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
Definition: DMMMoFEM.cpp:584

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

597 {
598 return DMMoFEMKSPSetComputeRHS<const std::string,
599 boost::shared_ptr<MoFEM::FEMethod>,
600 boost::shared_ptr<MoFEM::BasicMethod>,
601 boost::shared_ptr<MoFEM::BasicMethod>>(
602 dm, fe_name, method, pre_only, post_only);
603}

◆ 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:215
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:412

Definition at line 412 of file DMMMoFEM.cpp.

412 {
414 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
416 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
417 CHKERR dm_field->mField_ptr->getInterface<CommInterface>()
418 ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
420}

◆ DMMoFEMSetIsPartitioned()

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

sets if read mesh is partitioned

get if read mesh is partitioned

Definition at line 1070 of file DMMMoFEM.cpp.

1070 {
1071 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1073 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1074 dm_field->isPartitioned = is_partitioned;
1076}

◆ DMMoFEMSetKspCtx()

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

set MoFEM::KspCtx data structure

Definition at line 1032 of file DMMMoFEM.cpp.

1033 {
1034 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1036 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1037 dm_field->kspCtx = ksp_ctx;
1039}

◆ DMMoFEMSetSnesCtx()

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

Set MoFEM::SnesCtx data structure.

Definition at line 1058 of file DMMMoFEM.cpp.

1059 {
1060 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1062 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1063 dm_field->snesCtx = snes_ctx;
1065}

◆ 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.

Definition at line 403 of file DMMMoFEM.cpp.

403 {
405 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
407 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
408 dm_field->isSquareMatrix = square_problem;
410}

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

1106 {
1107 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1109 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1110 dm_field->tsCtx = ts_ctx;
1112}

◆ 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

Definition at line 665 of file DMMMoFEM.cpp.

668 {
669 return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
671 dm, fe_name, method, pre_only, post_only);
672}
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMMoFEM.cpp:665

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

678 {
679 return DMMoFEMSNESSetFunction<const std::string,
680 boost::shared_ptr<MoFEM::FEMethod>,
681 boost::shared_ptr<MoFEM::BasicMethod>,
682 boost::shared_ptr<MoFEM::BasicMethod>>(
683 dm, fe_name, method, pre_only, post_only);
684}

◆ 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

Definition at line 706 of file DMMMoFEM.cpp.

709 {
710 return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
712 dm, fe_name, method, pre_only, post_only);
713}
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMMoFEM.cpp:706

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

719 {
720 return DMMoFEMSNESSetJacobian<const std::string,
721 boost::shared_ptr<MoFEM::FEMethod>,
722 boost::shared_ptr<MoFEM::BasicMethod>,
723 boost::shared_ptr<MoFEM::BasicMethod>>(
724 dm, fe_name, method, pre_only, post_only);
725}

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

914 {
915 return DMMoFEMTSSetI2Function<const char *, MoFEM::FEMethod *,
917 dm, fe_name, method, pre_only, post_only);
919}
PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:922

◆ 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

Definition at line 922 of file DMMMoFEM.cpp.

925 {
926 return DMMoFEMTSSetI2Function<const std::string,
927 boost::shared_ptr<MoFEM::FEMethod>,
928 boost::shared_ptr<MoFEM::BasicMethod>,
929 boost::shared_ptr<MoFEM::BasicMethod>>(
930 dm, fe_name, method, pre_only, post_only);
932}

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

957 {
958 return DMMoFEMTSSetI2Jacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
959 MoFEM::BasicMethod *>(dm, fe_name, method,
960 pre_only, post_only);
961}
PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:964

◆ 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

Definition at line 964 of file DMMMoFEM.cpp.

967 {
968 return DMMoFEMTSSetI2Jacobian<const std::string,
969 boost::shared_ptr<MoFEM::FEMethod>,
970 boost::shared_ptr<MoFEM::BasicMethod>,
971 boost::shared_ptr<MoFEM::BasicMethod>>(
972 dm, fe_name, method, pre_only, post_only);
973}

◆ 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

Definition at line 747 of file DMMMoFEM.cpp.

750 {
751 return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
753 dm, fe_name, method, pre_only, post_only);
755}
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:747

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

761 {
762 return DMMoFEMTSSetIFunction<const std::string,
763 boost::shared_ptr<MoFEM::FEMethod>,
764 boost::shared_ptr<MoFEM::BasicMethod>,
765 boost::shared_ptr<MoFEM::BasicMethod>>(
766 dm, fe_name, method, pre_only, post_only);
768}

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

793 {
794 return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
795 MoFEM::BasicMethod *>(dm, fe_name, method,
796 pre_only, post_only);
797}
PetscErrorCode 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
Definition: DMMMoFEM.cpp:800

◆ 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

Definition at line 800 of file DMMMoFEM.cpp.

803 {
804 return DMMoFEMTSSetIJacobian<const std::string,
805 boost::shared_ptr<MoFEM::FEMethod>,
806 boost::shared_ptr<MoFEM::BasicMethod>,
807 boost::shared_ptr<MoFEM::BasicMethod>>(
808 dm, fe_name, method, pre_only, post_only);
809}

◆ DMMoFEMUnSetElement()

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

unset element from dm

Definition at line 460 of file DMMMoFEM.cpp.

460 {
461 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
463 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
464 ierr = dm_field->mField_ptr->modify_problem_unset_finite_element(
465 dm_field->problemName, fe_name);
466 CHKERRG(ierr);
468}

◆ 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

Definition at line 361 of file DMMMoFEM.cpp.

361 {
362 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
364 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
365 if (!dm->data) {
366 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
367 "data structure for MoFEM not yet created");
368 }
369 *m_field_ptr = dm_field->mField_ptr;
371}

◆ DMoFEMLoopDofs()

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

execute method for dofs on field in problem

Definition at line 552 of file DMMMoFEM.cpp.

553 {
554 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
556 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
557 ierr =
558 dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
559 *method, dm_field->rAnk, dm_field->rAnk);
560 CHKERRG(ierr);
562}

◆ 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

Definition at line 533 of file DMMMoFEM.cpp.

535 {
536 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
538 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
540 dm, fe_name, method, dm_field->rAnk, dm_field->rAnk, cache_ptr);
541 CHKERRG(ierr);
543}
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:514

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

548 {
549 return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get(), cache_ptr);
550}
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:533

◆ 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

Definition at line 514 of file DMMMoFEM.cpp.

516 {
518 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
519 ierr = dm_field->mField_ptr->loop_finite_elements(
520 dm_field->problemPtr, fe_name, *method, low_rank, up_rank, nullptr,
521 MF_EXIST, cache_ptr);
522 CHKERRG(ierr);
524}
@ MF_EXIST
Definition: definitions.h:100

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

528 {
529 return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
530 low_rank, up_rank, cache_ptr);
531}

◆ 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

Definition at line 482 of file DMMMoFEM.cpp.

483 {
484 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
486 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
487 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
488 dm_field->problemPtr, COL, g, mode, scatter_mode);
489 CHKERRG(ierr);
491}

◆ 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

Definition at line 470 of file DMMMoFEM.cpp.

471 {
473 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
475 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
476 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
477 dm_field->problemPtr, COL, l, mode, scatter_mode);
478 CHKERRG(ierr);
480}

◆ DMoFEMPostProcessFiniteElements()

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

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

Definition at line 503 of file DMMMoFEM.cpp.

503 {
504 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
506 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
507 ierr = dm_field->mField_ptr->problem_basic_method_postProcess(
508 dm_field->problemPtr, *method);
509 CHKERRG(ierr);
511}

◆ DMoFEMPreProcessFiniteElements()

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

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

Definition at line 493 of file DMMMoFEM.cpp.

493 {
494 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
496 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
497 ierr = dm_field->mField_ptr->problem_basic_method_preProcess(
498 dm_field->problemPtr, *method);
499 CHKERRG(ierr);
501}

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

376 {
378 CHKERR DMRegister(sname, DMCreate_MGViaApproxOrders);
380}
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.

Definition at line 47 of file DMMMoFEM.cpp.

47 {
49 CHKERR DMRegister(sname, DMCreate_MoFEM);
51}

◆ DMSetFromOptions_MoFEM()

PetscErrorCode MoFEM::DMSetFromOptions_MoFEM ( DM  dm)

Set options for MoFEM DM

Definition at line 1206 of file DMMMoFEM.cpp.

1206 {
1207#endif
1208
1209 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1211 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1212#if PETSC_VERSION_GE(3, 5, 3)
1213 ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1214 CHKERRG(ierr);
1215#else
1216 ierr = PetscOptionsHead("DMMoFEM Options");
1217 CHKERRG(ierr);
1218#endif
1219 ierr = PetscOptionsBool("-dm_is_partitioned",
1220 "set if mesh is partitioned (works which native MOAB "
1221 "file format, i.e. h5m",
1222 "DMSetUp", dm_field->isPartitioned,
1223 &dm_field->isPartitioned, NULL);
1224 CHKERRG(ierr);
1226}

◆ DMSetOperators_MoFEM()

PetscErrorCode MoFEM::DMSetOperators_MoFEM ( DM  dm)

Set operators for MoFEM dm.

Parameters
dm
Returns
error code

Definition at line 53 of file DMMMoFEM.cpp.

53 {
54 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
56
57 dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
58 dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
59 dm->ops->creatematrix = DMCreateMatrix_MoFEM;
60 dm->ops->setup = DMSetUp_MoFEM;
61 dm->ops->destroy = DMDestroy_MoFEM;
62 dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
63 dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
64 dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
65 dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
66 dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
67 dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
68
69 // Default matrix type
70 CHKERR DMSetMatType(dm, MATMPIAIJ);
71
73}
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1322
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1330
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMMoFEM.cpp:1134
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:83
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1359
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1206
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1114
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1394
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1228
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1388

◆ DMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSetUp_MoFEM ( DM  dm)

sets up the MoFEM structures inside a DM object

Definition at line 1228 of file DMMMoFEM.cpp.

1228 {
1229 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1230 ProblemsManager *prb_mng_ptr;
1232 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1233 CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1234
1235 if (dm_field->isCompDM) {
1236 // It is composite probelm
1237 CHKERR prb_mng_ptr->buildComposedProblem(
1238 dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1239 dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1240 } else {
1241 if (dm_field->isPartitioned) {
1242 CHKERR prb_mng_ptr->buildProblemOnDistributedMesh(
1243 dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1244 dm_field->verbosity);
1245 } else {
1246 CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1247 dm_field->isSquareMatrix == PETSC_TRUE,
1248 dm_field->verbosity);
1249 CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1250 dm_field->verbosity);
1251 }
1252 }
1253
1254 // Partition finite elements
1255 if (dm_field->isPartitioned) {
1256 CHKERR prb_mng_ptr->partitionFiniteElements(
1257 dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1258 CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1259 dm_field->problemName, dm_field->verbosity);
1260 } else {
1261 // partition finite elemnets
1262 CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1263 -1, -1, dm_field->verbosity);
1264 // Get ghost DOFs
1265 CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1266 dm_field->verbosity);
1267 }
1268
1269 // Set flag that problem is build and partitioned
1270 dm_field->isProblemBuild = PETSC_TRUE;
1271
1273}

◆ DMSubDMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSubDMSetUp_MoFEM ( DM  subdm)

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

Definition at line 1275 of file DMMMoFEM.cpp.

1275 {
1276 PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1277 ProblemsManager *prb_mng_ptr;
1279
1280 DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1281
1282 // build sub dm problem
1283 CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1284
1285 map<std::string, std::pair<EntityType, EntityType>> *entity_map_row = nullptr;
1286 map<std::string, std::pair<EntityType, EntityType>> *entity_map_col = nullptr;
1287
1288 if (subdm_field->mapTypeRow)
1289 entity_map_row = subdm_field->mapTypeRow.get();
1290 if (subdm_field->mapTypeCol)
1291 entity_map_row = subdm_field->mapTypeCol.get();
1292
1293 CHKERR prb_mng_ptr->buildSubProblem(
1294 subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1295 subdm_field->problemMainOfSubPtr->getName(),
1296 subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1297 subdm_field->verbosity);
1298
1299 // partition problem
1300 subdm_field->isPartitioned = subdm_field->isPartitioned;
1301 if (subdm_field->isPartitioned) {
1302 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1303 0, subdm_field->sIze,
1304 subdm_field->verbosity);
1305 // set ghost nodes
1306 CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1307 subdm_field->problemName, subdm_field->verbosity);
1308 } else {
1309 // partition finite elements
1310 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1311 -1, -1, subdm_field->verbosity);
1312 // set ghost nodes
1313 CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1314 subdm_field->verbosity);
1315 }
1316
1317 subdm_field->isProblemBuild = PETSC_TRUE;
1318
1320}

◆ smartCreateDMMatrix()

auto MoFEM::smartCreateDMMatrix ( DM  dm)
inline

Get smart matrix from DM.

Definition at line 953 of file DMMoFEM.hpp.

953 {
954 SmartPetscObj<Mat> a;
956 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
957 return a;
958};
constexpr double a
MPI_Comm getCommFromPetscObject(PetscObject obj)
Get the Comm From Petsc Object object.

◆ smartCreateDMVector()

auto MoFEM::smartCreateDMVector ( DM  dm)
inline

Get smart vector from DM.

Definition at line 965 of file DMMoFEM.hpp.

965 {
966 SmartPetscObj<Vec> v;
968 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
969 return v;
970};
const double v
phase velocity of light in medium (cm/ns)