v0.13.2
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, std::string fe_name)
 Resolve shared entities. More...
 
PetscErrorCode MoFEM::DMMoFEMGetProblemFiniteElementLayout (DM dm, std::string fe_name, PetscLayout *layout)
 Get finite elements layout in the problem. More...
 
PetscErrorCode MoFEM::DMMoFEMAddElement (DM dm, std::string fe_name)
 add element to dm More...
 
PetscErrorCode MoFEM::DMMoFEMAddElement (DM dm, std::vector< std::string > fe_name)
 add element to dm More...
 
PetscErrorCode MoFEM::DMMoFEMUnSetElement (DM dm, std::string 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[])
 
PetscErrorCode MoFEM::DMMoFEMAddSubFieldRow (DM dm, const char field_name[], boost::shared_ptr< Range > r_ptr)
 Add field to sub dm problem on rows. More...
 
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 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 DMMoFEM.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: DMMoFEM.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 1435 of file DMMoFEM.cpp.

1436 {
1437 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1439
1440 if (numFields) {
1441 *numFields = 0;
1442 }
1443 if (fieldNames) {
1444 *fieldNames = NULL;
1445 }
1446 if (fields) {
1447 *fields = NULL;
1448 }
1449
1450 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1451 auto fields_ptr = dm_field->mField_ptr->get_fields();
1452 Field_multiIndex::iterator fit, hi_fit;
1453 fit = fields_ptr->begin();
1454 hi_fit = fields_ptr->end();
1455 *numFields = std::distance(fit, hi_fit);
1456
1457 if (fieldNames) {
1458 CHKERR PetscMalloc1(*numFields, fieldNames);
1459 }
1460 if (fields) {
1461 CHKERR PetscMalloc1(*numFields, fields);
1462 }
1463
1464 for (int f = 0; fit != hi_fit; fit++, f++) {
1465 if (fieldNames) {
1466 CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1467 (char **)&(*fieldNames)[f]);
1468 }
1469 if (fields) {
1470 CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1471 ->isCreateProblemFieldAndRank(
1472 dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1473 fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1474 }
1475 }
1476
1478}
@ 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 1165 of file DMMoFEM.cpp.

1165 {
1166 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1168 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1169 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1170 dm_field->problemName, COL, g_ptr);
1171 CHKERR VecSetDM(g_ptr, dm);
1173}
@ COL
Definition: definitions.h:123

◆ 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 1155 of file DMMoFEM.cpp.

1155 {
1156 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1158 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1159 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1160 dm_field->problemName, COL, g);
1161 CHKERR VecSetDM(*g, dm);
1163}
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 1175 of file DMMoFEM.cpp.

1175 {
1176 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1178 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1179 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1180 dm_field->problemName, COL, l);
1181 CHKERR VecSetDM(*l, dm);
1183}
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: DMMoFEM.cpp:1185

◆ 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 1185 of file DMMoFEM.cpp.

1185 {
1186 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1188 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1189 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1190 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1191 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1192 M);
1193 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1194 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1195 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1196 M);
1197 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1198 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1199 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1200 dm_field->problemName, M);
1201 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1202 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1203 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1204 dm_field->problemName, M);
1205 } else {
1206 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1207 "Matrix type not implemented");
1208 }
1209 CHKERR MatSetDM(*M, dm);
1211}
@ 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 1213 of file DMMoFEM.cpp.

1213 {
1214 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1216 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1217 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1218 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1219 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1220 M);
1221 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1222 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1223 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1224 M);
1225 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1226 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1227 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1228 dm_field->problemName, M);
1229 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1230 CHKERR dm_field->mField_ptr->getInterface<MatrixManager>()
1231 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1232 dm_field->problemName, M);
1233 } else {
1234 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1235 "Matrix type not implemented");
1236 }
1237 CHKERR MatSetDM(M, dm);
1239}

◆ DMDestroy_MoFEM()

PetscErrorCode MoFEM::DMDestroy_MoFEM ( DM  dm)

Destroys dm with MoFEM data structure.

destroy the MoFEM structure

Definition at line 83 of file DMMoFEM.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 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
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308

◆ 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 1363 of file DMMoFEM.cpp.

1364 {
1365 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1367 CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1369}

◆ 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 1371 of file DMMoFEM.cpp.

1371 {
1373 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1375
1376 CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1377
1378 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1379 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1380 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1381
1382 double *array_loc, *array_glob;
1383 CHKERR VecGetArray(l, &array_loc);
1384 CHKERR VecGetArray(g, &array_glob);
1385 switch (mode) {
1386 case INSERT_VALUES:
1387 cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1388 break;
1389 case ADD_VALUES:
1390 cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1391 break;
1392 default:
1393 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1394 }
1395 CHKERR VecRestoreArray(l, &array_loc);
1396 CHKERR VecRestoreArray(g, &array_glob);
1398}

◆ 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 1400 of file DMMoFEM.cpp.

1401 {
1402
1403 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1405
1406 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1407 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1408 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1409
1410 double *array_loc, *array_glob;
1411 CHKERR VecGetArray(l, &array_loc);
1412 CHKERR VecGetArray(g, &array_glob);
1413 switch (mode) {
1414 case INSERT_VALUES:
1415 cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1416 break;
1417 case ADD_VALUES:
1418 cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1419 break;
1420 default:
1421 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1422 }
1423 CHKERR VecRestoreArray(l, &array_loc);
1424 CHKERR VecRestoreArray(g, &array_glob);
1425
1427}

◆ 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 1429 of file DMMoFEM.cpp.

1429 {
1430 //
1433}

◆ 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 373 of file DMMoFEM.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 if (!dm_field->isCompDM) {
382 dm_field->isCompDM = PETSC_TRUE;
383 }
384 if (dm_field->isSquareMatrix) {
385 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
386 "No need to add problem on column when problem block structurally "
387 "symmetric");
388 }
389 dm_field->colCompPrb.push_back(prb_name);
391}
@ MOFEM_INVALID_DATA
Definition: definitions.h:36

◆ DMMoFEMAddElement() [1/2]

PetscErrorCode MoFEM::DMMoFEMAddElement ( DM  dm,
std::string  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 485 of file DMMoFEM.cpp.

485 {
486 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
488 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
489 CHKERR dm_field->mField_ptr->modify_problem_add_finite_element(
490 dm_field->problemName, fe_name);
492}

◆ DMMoFEMAddElement() [2/2]

PetscErrorCode MoFEM::DMMoFEMAddElement ( DM  dm,
std::vector< std::string >  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 494 of file DMMoFEM.cpp.

494 {
496 for (auto fe : fe_name) {
498 }
500}
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:485

◆ 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 355 of file DMMoFEM.cpp.

355 {
356 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
358 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
359 if (!dm->data) {
360 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
361 "data structure for MoFEM not yet created");
362 }
363 if (!dm_field->isCompDM) {
364 dm_field->isCompDM = PETSC_TRUE;
365 }
366 dm_field->rowCompPrb.push_back(prb_name);
367 if (dm_field->isSquareMatrix) {
368 dm_field->colCompPrb.push_back(prb_name);
369 }
371}

◆ DMMoFEMAddSubFieldCol()

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

Add field to sub dm problem on columns

Definition at line 277 of file DMMoFEM.cpp.

277 {
278 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
280 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
281 if (!dm->data) {
282 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
283 "data structure for MoFEM not yet created");
284 }
285 if (!dm_field->isSubDM) {
286 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
287 }
288 dm_field->colFields.push_back(field_name);
289 dm_field->mapTypeCol.erase(field_name);
291}
constexpr auto field_name

◆ DMMoFEMAddSubFieldRow() [1/2]

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

Add field to sub dm problem on rows

Definition at line 244 of file DMMoFEM.cpp.

244 {
245 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
247 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
248 if (!dm->data) {
249 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
250 "data structure for MoFEM not yet created");
251 }
252 if (!dm_field->isSubDM) {
253 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
254 }
255 dm_field->rowFields.push_back(field_name);
256 dm_field->mapTypeRow.erase(field_name);
258}

◆ DMMoFEMAddSubFieldRow() [2/2]

PetscErrorCode MoFEM::DMMoFEMAddSubFieldRow ( DM  dm,
const char  field_name[],
boost::shared_ptr< Range r_ptr 
)

Add field to sub dm problem on rows.

Parameters
dm
field_name
m
Returns
PetscErrorCode

Definition at line 260 of file DMMoFEM.cpp.

261 {
262 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
264 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
265 if (!dm->data) {
266 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
267 "data structure for MoFEM not yet created");
268 }
269 if (!dm_field->isSubDM) {
270 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
271 }
272 dm_field->rowFields.push_back(field_name);
273 dm_field->mapTypeRow[field_name] = r_ptr;
275}

◆ 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 DMMoFEM.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 221 of file DMMoFEM.cpp.

221 {
223
224 if (!dm->data) {
225 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
226 "data structure for MoFEM not yet created");
227 }
228 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
229
230 CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
231 dm_field->problemPtr->getBitRefLevel(),
232 dm_field->problemPtr->getBitRefLevelMask());
233
234 DMCtx *subdm_field = (DMCtx *)subdm->data;
235 subdm_field->isSubDM = PETSC_TRUE;
236 subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
237 subdm_field->isPartitioned = dm_field->isPartitioned;
238 subdm_field->isSquareMatrix = PETSC_FALSE;
239 subdm->ops->setup = DMSubDMSetUp_MoFEM;
240
242}
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: DMMoFEM.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);
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMoFEM.cpp:1480
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76

Definition at line 1480 of file DMMoFEM.cpp.

1481 {
1482 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1484 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1485 CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1486 ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1487 field_name, 0, 1000, is);
1489}

◆ 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 393 of file DMMoFEM.cpp.

393 {
395 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
397 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
398 *is_comp_dm = dm_field->isCompDM;
400}

◆ DMMoFEMGetIsPartitioned()

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

get if read mesh is partitioned

Definition at line 1122 of file DMMoFEM.cpp.

1122 {
1123 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1125 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1126 *is_partitioned = dm_field->isPartitioned;
1128}

◆ 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 310 of file DMMoFEM.cpp.

310 {
312 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
314 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
315 *is_sub_dm = dm_field->isSubDM;
317}

◆ DMMoFEMGetKspCtx() [1/2]

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

get MoFEM::KspCtx data structure

Definition at line 1065 of file DMMoFEM.cpp.

1065 {
1066 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1068 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1069 const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
1071}

◆ DMMoFEMGetKspCtx() [2/2]

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

get MoFEM::KspCtx data structure

Definition at line 1056 of file DMMoFEM.cpp.

1056 {
1057 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1059 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1060 *ksp_ctx = dm_field->kspCtx.get();
1062}

◆ DMMoFEMGetProblemFiniteElementLayout()

PetscErrorCode MoFEM::DMMoFEMGetProblemFiniteElementLayout ( DM  dm,
std::string  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 461 of file DMMoFEM.cpp.

462 {
464 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
466 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
467
468 MPI_Comm comm;
469 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
470 CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
471 layout);
473}

◆ DMMoFEMGetProblemPtr()

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

Get pointer to problem data structure.

Definition at line 414 of file DMMoFEM.cpp.

414 {
415 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
417 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
418 if (!dm->data) {
419 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
420 "data structure for MoFEM not yet created");
421 }
422 *problem_ptr = dm_field->problemPtr;
424}

◆ DMMoFEMGetSnesCtx() [1/2]

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

get MoFEM::SnesCtx data structure

Definition at line 1091 of file DMMoFEM.cpp.

1091 {
1092 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1094 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1095 const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
1097}

◆ DMMoFEMGetSnesCtx() [2/2]

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

get MoFEM::SnesCtx data structure

Definition at line 1082 of file DMMoFEM.cpp.

1082 {
1083 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1085 DMCtx *dm_field = (DMCtx *)dm->data;
1086 *snes_ctx = dm_field->snesCtx.get();
1088}

◆ 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 475 of file DMMoFEM.cpp.

475 {
478 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
480 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
481 *square_problem = dm_field->isSquareMatrix;
483}

◆ DMMoFEMGetTsCtx() [1/2]

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

get MoFEM::TsCtx data structure

Definition at line 1138 of file DMMoFEM.cpp.

1139 {
1140 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1142 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1143 const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1145}
MoFEM::TsCtx * ts_ctx
Definition: level_set.cpp:1876

◆ DMMoFEMGetTsCtx() [2/2]

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

get MoFEM::TsCtx data structure

Definition at line 1130 of file DMMoFEM.cpp.

1130 {
1131 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1133 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1134 *ts_ctx = dm_field->tsCtx.get();
1136}

◆ 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 666 of file DMMoFEM.cpp.

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

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

◆ 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 625 of file DMMoFEM.cpp.

628 {
629 return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
631 dm, fe_name, method, pre_only, post_only);
632}
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: DMMoFEM.cpp:625

◆ 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 635 of file DMMoFEM.cpp.

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

◆ DMMoFEMResolveSharedFiniteElements()

PetscErrorCode MoFEM::DMMoFEMResolveSharedFiniteElements ( DM  dm,
std::string  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, std::string fe_name)
Resolve shared entities.
Definition: DMMoFEM.cpp:452

Definition at line 452 of file DMMoFEM.cpp.

452 {
453 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
455 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
456 CHKERR dm_field->mField_ptr->getInterface<CommInterface>()
457 ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
459}

◆ DMMoFEMSetIsPartitioned()

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

sets if read mesh is partitioned

get if read mesh is partitioned

Definition at line 1111 of file DMMoFEM.cpp.

1111 {
1112 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1114 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1115 dm_field->isPartitioned = is_partitioned;
1117}

◆ DMMoFEMSetKspCtx()

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

set MoFEM::KspCtx data structure

Definition at line 1073 of file DMMoFEM.cpp.

1074 {
1075 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1077 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1078 dm_field->kspCtx = ksp_ctx;
1080}

◆ DMMoFEMSetSnesCtx()

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

Set MoFEM::SnesCtx data structure.

Definition at line 1099 of file DMMoFEM.cpp.

1100 {
1101 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1103 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1104 dm_field->snesCtx = snes_ctx;
1106}

◆ 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 444 of file DMMoFEM.cpp.

444 {
445 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
447 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
448 dm_field->isSquareMatrix = square_problem;
450}

◆ 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 1147 of file DMMoFEM.cpp.

1147 {
1148 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1150 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1151 dm_field->tsCtx = ts_ctx;
1153}

◆ 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 706 of file DMMoFEM.cpp.

709 {
710 return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
712 dm, fe_name, method, pre_only, post_only);
713}
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: DMMoFEM.cpp:706

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

719 {
720 return DMMoFEMSNESSetFunction<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}

◆ 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 747 of file DMMoFEM.cpp.

750 {
751 return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
753 dm, fe_name, method, pre_only, post_only);
754}
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: DMMoFEM.cpp:747

◆ 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 757 of file DMMoFEM.cpp.

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

◆ 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 952 of file DMMoFEM.cpp.

955 {
956 return DMMoFEMTSSetI2Function<const char *, MoFEM::FEMethod *,
958 dm, fe_name, method, pre_only, post_only);
960}
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: DMMoFEM.cpp:963

◆ 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 963 of file DMMoFEM.cpp.

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

◆ 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 995 of file DMMoFEM.cpp.

998 {
999 return DMMoFEMTSSetI2Jacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
1000 MoFEM::BasicMethod *>(dm, fe_name, method,
1001 pre_only, post_only);
1002}
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: DMMoFEM.cpp:1005

◆ 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 1005 of file DMMoFEM.cpp.

1008 {
1009 return DMMoFEMTSSetI2Jacobian<const std::string,
1010 boost::shared_ptr<MoFEM::FEMethod>,
1011 boost::shared_ptr<MoFEM::BasicMethod>,
1012 boost::shared_ptr<MoFEM::BasicMethod>>(
1013 dm, fe_name, method, pre_only, post_only);
1014}

◆ 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 788 of file DMMoFEM.cpp.

791 {
792 return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
794 dm, fe_name, method, pre_only, post_only);
796}
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: DMMoFEM.cpp:788

◆ 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 799 of file DMMoFEM.cpp.

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

◆ 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 831 of file DMMoFEM.cpp.

834 {
835 return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
836 MoFEM::BasicMethod *>(dm, fe_name, method,
837 pre_only, post_only);
838}
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: DMMoFEM.cpp:841

◆ 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 841 of file DMMoFEM.cpp.

844 {
845 return DMMoFEMTSSetIJacobian<const std::string,
846 boost::shared_ptr<MoFEM::FEMethod>,
847 boost::shared_ptr<MoFEM::BasicMethod>,
848 boost::shared_ptr<MoFEM::BasicMethod>>(
849 dm, fe_name, method, pre_only, post_only);
850}

◆ DMMoFEMUnSetElement()

PetscErrorCode MoFEM::DMMoFEMUnSetElement ( DM  dm,
std::string  fe_name 
)

unset element from dm

Definition at line 502 of file DMMoFEM.cpp.

502 {
503 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
505 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
506 CHKERR dm_field->mField_ptr->modify_problem_unset_finite_element(
507 dm_field->problemName, fe_name);
509}

◆ 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 402 of file DMMoFEM.cpp.

402 {
403 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
405 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
406 if (!dm->data) {
407 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
408 "data structure for MoFEM not yet created");
409 }
410 *m_field_ptr = dm_field->mField_ptr;
412}

◆ DMoFEMLoopDofs()

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

execute method for dofs on field in problem

Definition at line 593 of file DMMoFEM.cpp.

594 {
595 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
597 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
598 ierr =
599 dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
600 *method, dm_field->rAnk, dm_field->rAnk);
601 CHKERRG(ierr);
603}
static PetscErrorCode ierr

◆ 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 574 of file DMMoFEM.cpp.

576 {
577 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
579 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
581 dm, fe_name, method, dm_field->rAnk, dm_field->rAnk, cache_ptr);
582 CHKERRG(ierr);
584}
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: DMMoFEM.cpp:555

◆ 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 587 of file DMMoFEM.cpp.

589 {
590 return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get(), cache_ptr);
591}
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:574

◆ 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 555 of file DMMoFEM.cpp.

557 {
559 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
560 ierr = dm_field->mField_ptr->loop_finite_elements(
561 dm_field->problemPtr, fe_name, *method, low_rank, up_rank, nullptr,
562 MF_EXIST, cache_ptr);
563 CHKERRG(ierr);
565}
@ 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 567 of file DMMoFEM.cpp.

569 {
570 return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
571 low_rank, up_rank, cache_ptr);
572}

◆ 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 523 of file DMMoFEM.cpp.

524 {
525 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
527 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
528 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
529 dm_field->problemPtr, COL, g, mode, scatter_mode);
530 CHKERRG(ierr);
532}

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

512 {
514 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
516 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
517 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
518 dm_field->problemPtr, COL, l, mode, scatter_mode);
519 CHKERRG(ierr);
521}

◆ DMoFEMPostProcessFiniteElements()

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

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

Definition at line 544 of file DMMoFEM.cpp.

544 {
545 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
547 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
548 ierr = dm_field->mField_ptr->problem_basic_method_postProcess(
549 dm_field->problemPtr, *method);
550 CHKERRG(ierr);
552}

◆ DMoFEMPreProcessFiniteElements()

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

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

Definition at line 534 of file DMMoFEM.cpp.

534 {
535 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
537 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
538 ierr = dm_field->mField_ptr->problem_basic_method_preProcess(
539 dm_field->problemPtr, *method);
540 CHKERRG(ierr);
542}

◆ 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 DMMoFEM.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 1247 of file DMMoFEM.cpp.

1247 {
1248#endif
1249
1250 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1252 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1253#if PETSC_VERSION_GE(3, 5, 3)
1254 ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1255 CHKERRG(ierr);
1256#else
1257 ierr = PetscOptionsHead("DMMoFEM Options");
1258 CHKERRG(ierr);
1259#endif
1260 ierr = PetscOptionsBool("-dm_is_partitioned",
1261 "set if mesh is partitioned (works which native MOAB "
1262 "file format, i.e. h5m",
1263 "DMSetUp", dm_field->isPartitioned,
1264 &dm_field->isPartitioned, NULL);
1265 CHKERRG(ierr);
1267}

◆ DMSetOperators_MoFEM()

PetscErrorCode MoFEM::DMSetOperators_MoFEM ( DM  dm)

Set operators for MoFEM dm.

Parameters
dm
Returns
error code

Definition at line 53 of file DMMoFEM.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: DMMoFEM.cpp:1363
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1371
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMoFEM.cpp:1175
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMoFEM.cpp:83
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1400
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMoFEM.cpp:1247
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMoFEM.cpp:1155
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMoFEM.cpp:1435
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMoFEM.cpp:1269
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1429

◆ DMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSetUp_MoFEM ( DM  dm)

sets up the MoFEM structures inside a DM object

Definition at line 1269 of file DMMoFEM.cpp.

1269 {
1270 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1271 ProblemsManager *prb_mng_ptr;
1273 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1274 CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1275
1276 if (dm_field->isCompDM) {
1277 // It is composite probelm
1278 CHKERR prb_mng_ptr->buildComposedProblem(
1279 dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1280 dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1281 } else {
1282 if (dm_field->isPartitioned) {
1283 CHKERR prb_mng_ptr->buildProblemOnDistributedMesh(
1284 dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1285 dm_field->verbosity);
1286 } else {
1287 CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1288 dm_field->isSquareMatrix == PETSC_TRUE,
1289 dm_field->verbosity);
1290 CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1291 dm_field->verbosity);
1292 }
1293 }
1294
1295 // Partition finite elements
1296 if (dm_field->isPartitioned) {
1297 CHKERR prb_mng_ptr->partitionFiniteElements(
1298 dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1299 CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1300 dm_field->problemName, dm_field->verbosity);
1301 } else {
1302 // partition finite elemnets
1303 CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1304 -1, -1, dm_field->verbosity);
1305 // Get ghost DOFs
1306 CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1307 dm_field->verbosity);
1308 }
1309
1310 // Set flag that problem is build and partitioned
1311 dm_field->isProblemBuild = PETSC_TRUE;
1312
1314}

◆ DMSubDMSetUp_MoFEM()

PetscErrorCode MoFEM::DMSubDMSetUp_MoFEM ( DM  subdm)

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

Definition at line 1316 of file DMMoFEM.cpp.

1316 {
1317 PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1318 ProblemsManager *prb_mng_ptr;
1320
1321 DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1322
1323 // build sub dm problem
1324 CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1325
1326 map<std::string, boost::shared_ptr<Range>> *entity_map_row = nullptr;
1327 map<std::string, boost::shared_ptr<Range>> *entity_map_col = nullptr;
1328
1329 if (subdm_field->mapTypeRow.size())
1330 entity_map_row = &subdm_field->mapTypeRow;
1331 if (subdm_field->mapTypeCol.size())
1332 entity_map_col = &subdm_field->mapTypeCol;
1333
1334 CHKERR prb_mng_ptr->buildSubProblem(
1335 subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1336 subdm_field->problemMainOfSubPtr->getName(),
1337 subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1338 subdm_field->verbosity);
1339
1340 // partition problem
1341 subdm_field->isPartitioned = subdm_field->isPartitioned;
1342 if (subdm_field->isPartitioned) {
1343 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1344 0, subdm_field->sIze,
1345 subdm_field->verbosity);
1346 // set ghost nodes
1347 CHKERR prb_mng_ptr->partitionGhostDofsOnDistributedMesh(
1348 subdm_field->problemName, subdm_field->verbosity);
1349 } else {
1350 // partition finite elements
1351 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1352 -1, -1, subdm_field->verbosity);
1353 // set ghost nodes
1354 CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1355 subdm_field->verbosity);
1356 }
1357
1358 subdm_field->isProblemBuild = PETSC_TRUE;
1359
1361}
if(!static_cast< bool >(ifstream(param_file)))

◆ smartCreateDMMatrix()

auto MoFEM::smartCreateDMMatrix ( DM  dm)
inline

Get smart matrix from DM.

Definition at line 983 of file DMMoFEM.hpp.

983 {
984 SmartPetscObj<Mat> a;
986 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
987 return a;
988};
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 995 of file DMMoFEM.hpp.

995 {
996 SmartPetscObj<Vec> v;
998 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
999 return v;
1000};
const double v
phase velocity of light in medium (cm/ns)