6#if PETSC_VERSION_GE(3, 6, 0)
7#include <petsc/private/dmimpl.h>
10#include <petsc-private/dmimpl.h>
11#include <petsc-private/vecimpl.h>
19 : mField_ptr(PETSC_NULL), isProblemBuild(PETSC_FALSE),
20 isPartitioned(PETSC_FALSE), isSquareMatrix(PETSC_TRUE),
21 isSubDM(PETSC_FALSE), isCompDM(PETSC_FALSE), destroyProblem(PETSC_FALSE),
22 verbosity(
VERBOSE), referenceNumber(0) {
25 auto core_log = logging::core::get();
43 *iface =
const_cast<DMCtx *
>(
this);
54 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
70 CHKERR DMSetMatType(dm, MATMPIAIJ);
76 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
78 dm->data =
new DMCtx();
84 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
85 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
89 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
92 MPI_Comm_compare(comm, PETSC_COMM_SELF, &result);
93 if (result == MPI_IDENT)
95 <<
"MoFEM DM destroy for problem " << dm_field->
problemName
99 <<
"MoFEM DM destroy for problem " << dm_field->
problemName
110 delete static_cast<DMCtx *
>(dm->data);
119 const char problem_name[],
124 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
127 "data structure for MoFEM not yet created");
131 "DM function not implemented into MoFEM");
149 boost::shared_ptr<KspCtx>(
new KspCtx(*m_field_ptr, problem_name));
151 boost::shared_ptr<SnesCtx>(
new SnesCtx(*m_field_ptr, problem_name));
153 boost::shared_ptr<TsCtx>(
new TsCtx(*m_field_ptr, problem_name));
156 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
158 MPI_Comm_compare(comm, m_field_ptr->
get_comm(), &result);
159 if (result > MPI_CONGRUENT) {
161 "MoFEM and DM using different communicators");
163 MPI_Comm_size(comm, &dm_field->
sIze);
164 MPI_Comm_rank(comm, &dm_field->
rAnk);
170 MPI_Comm_compare(comm, PETSC_COMM_SELF, &result);
171 if (result == MPI_IDENT) {
173 <<
"MoFEM DM created for problem " << dm_field->
problemName;
177 <<
"MoFEM DM created for problem " << dm_field->
problemName;
189 "data structure for MoFEM not yet created");
191 auto *dm_field =
static_cast<DMCtx *
>(dm->data);
194 delete static_cast<DMCtx *
>(dm_duplicate->data);
196 dm_duplicate->data = dm->data;
197 ++(
static_cast<DMCtx *
>(dm_duplicate->data)->referenceNumber);
206 "data structure for MoFEM not yet created on dm");
209 "data structure for MoFEM not yet created on swap dm");
211 auto *dm_field =
static_cast<DMCtx *
>(dm->data);
212 auto *dm_field_swap =
static_cast<DMCtx *
>(dm_swap->data);
214 auto tmp_field = dm_field;
215 dm_field = dm_field_swap;
216 dm_field_swap = tmp_field;
226 "data structure for MoFEM not yet created");
228 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
235 subdm_field->
isSubDM = PETSC_TRUE;
245 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
247 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
250 "data structure for MoFEM not yet created");
261 boost::shared_ptr<Range> r_ptr) {
262 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
264 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
267 "data structure for MoFEM not yet created");
278 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
280 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
283 "data structure for MoFEM not yet created");
294 boost::shared_ptr<Range> r_ptr) {
295 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
297 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
300 "data structure for MoFEM not yet created");
312 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
314 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
315 *is_sub_dm = dm_field->
isSubDM;
321 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
323 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
324 if (dm_field->
isSubDM != PETSC_TRUE) {
326 "This DM is not created as a SubDM");
331 boost::shared_ptr<Problem::SubProblemData> sub_data =
333 CHKERR sub_data->getRowIs(is);
339 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
341 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
342 if (dm_field->
isSubDM != PETSC_TRUE) {
344 "This DM is not created as a SubDM");
349 boost::shared_ptr<Problem::SubProblemData> sub_data =
351 CHKERR sub_data->getColIs(is);
356 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
358 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
361 "data structure for MoFEM not yet created");
374 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
376 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
379 "data structure for MoFEM not yet created");
386 "No need to add problem on column when problem block structurally "
395 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
397 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
403 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
405 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
408 "data structure for MoFEM not yet created");
415 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
417 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
420 "data structure for MoFEM not yet created");
428 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
430 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
437 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
439 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
445 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
447 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
453 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
455 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
457 ->resolveSharedFiniteElements(dm_field->
problemPtr, fe_name);
462 PetscLayout *layout) {
464 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
466 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
469 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
478 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
480 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
486 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
488 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
496 for (
auto fe : fe_name) {
503 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
505 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
512 ScatterMode scatter_mode) {
514 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
516 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
524 ScatterMode scatter_mode) {
525 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
527 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
535 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
537 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
545 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
547 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
559 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
561 dm_field->
problemPtr, fe_name, *method, low_rank, up_rank,
nullptr,
568 DM dm,
const std::string fe_name, boost::shared_ptr<MoFEM::FEMethod> method,
571 low_rank, up_rank, cache_ptr);
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);
588 boost::shared_ptr<MoFEM::FEMethod> method,
595 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
597 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
600 *method, dm_field->
rAnk, dm_field->
rAnk);
605template <
class S,
class T0,
class T1,
class T2>
607 T1 pre_only, T2 post_only) {
608 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
610 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
612 dm_field->
kspCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
615 dm_field->
kspCtx->get_loops_to_do_Rhs().push_back(
619 dm_field->
kspCtx->get_postProcess_to_do_Rhs().push_back(post_only);
631 dm, fe_name, method, pre_only, post_only);
636 boost::shared_ptr<MoFEM::FEMethod> method,
637 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
638 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
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);
646template <
class S,
class T0,
class T1,
class T2>
648 T1 pre_only, T2 post_only) {
649 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
651 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
653 dm_field->
kspCtx->get_preProcess_to_do_Mat().push_back(pre_only);
656 dm_field->
kspCtx->get_loops_to_do_Mat().push_back(
660 dm_field->
kspCtx->get_postProcess_to_do_Mat().push_back(post_only);
673 dm, fe_name, method, pre_only, post_only);
678 boost::shared_ptr<MoFEM::FEMethod> method,
679 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
680 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
682 boost::shared_ptr<MoFEM::FEMethod>>(
683 dm, fe_name, method, pre_only, post_only);
686template <
class S,
class T0,
class T1,
class T2>
688 T1 pre_only, T2 post_only) {
689 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
691 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
693 dm_field->
snesCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
696 dm_field->
snesCtx->get_loops_to_do_Rhs().push_back(
700 dm_field->
snesCtx->get_postProcess_to_do_Rhs().push_back(post_only);
712 dm, fe_name, method, pre_only, post_only);
717 boost::shared_ptr<MoFEM::FEMethod> method,
718 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
719 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
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);
727template <
class S,
class T0,
class T1,
class T2>
729 T1 pre_only, T2 post_only) {
730 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
732 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
734 dm_field->
snesCtx->get_preProcess_to_do_Mat().push_back(pre_only);
737 dm_field->
snesCtx->get_loops_to_do_Mat().push_back(
741 dm_field->
snesCtx->get_postProcess_to_do_Mat().push_back(post_only);
753 dm, fe_name, method, pre_only, post_only);
758 boost::shared_ptr<MoFEM::FEMethod> method,
759 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
760 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
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);
768template <
class S,
class T0,
class T1,
class T2>
770 T1 pre_only, T2 post_only) {
771 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
773 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
775 dm_field->
tsCtx->getPreProcessIFunction().push_back(pre_only);
778 dm_field->
tsCtx->getLoopsIFunction().push_back(
782 dm_field->
tsCtx->getPostProcessIFunction().push_back(post_only);
794 dm, fe_name, method, pre_only, post_only);
800 boost::shared_ptr<MoFEM::FEMethod> method,
801 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
802 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
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);
811template <
class S,
class T0,
class T1,
class T2>
813 T1 pre_only, T2 post_only) {
814 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
816 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
818 dm_field->
tsCtx->getPreProcessIJacobian().push_back(pre_only);
821 dm_field->
tsCtx->getLoopsIJacobian().push_back(
825 dm_field->
tsCtx->getPostProcessIJacobian().push_back(post_only);
837 pre_only, post_only);
842 boost::shared_ptr<MoFEM::FEMethod> method,
843 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
844 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
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);
852template <
class S,
class T0,
class T1,
class T2>
854 T1 pre_only, T2 post_only) {
855 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
857 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
859 dm_field->
tsCtx->getPreProcessRHSFunction().push_back(pre_only);
861 dm_field->
tsCtx->getLoopsRHSFunction().push_back(
864 dm_field->
tsCtx->getPostProcessRHSFunction().push_back(post_only);
871 boost::shared_ptr<MoFEM::FEMethod> method,
872 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
873 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
875 boost::shared_ptr<MoFEM::FEMethod>,
876 boost::shared_ptr<MoFEM::BasicMethod>,
877 boost::shared_ptr<MoFEM::BasicMethod>>(
878 dm, fe_name, method, pre_only, post_only);
888 dm, fe_name, method, pre_only, post_only);
892template <
class S,
class T0,
class T1,
class T2>
894 T1 pre_only, T2 post_only) {
895 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
897 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
899 dm_field->
tsCtx->getPreProcessRHSFunction().push_back(pre_only);
901 dm_field->
tsCtx->getLoopsRHSFunction().push_back(
904 dm_field->
tsCtx->getPostProcessRHSFunction().push_back(post_only);
911 boost::shared_ptr<MoFEM::FEMethod> method,
912 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
913 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
915 boost::shared_ptr<MoFEM::FEMethod>,
916 boost::shared_ptr<MoFEM::BasicMethod>,
917 boost::shared_ptr<MoFEM::BasicMethod>>(
918 dm, fe_name, method, pre_only, post_only);
928 dm, fe_name, method, pre_only, post_only);
932template <
class S,
class T0,
class T1,
class T2>
934 T1 pre_only, T2 post_only) {
935 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
937 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
939 dm_field->
tsCtx->getPreProcessIFunction().push_back(pre_only);
942 dm_field->
tsCtx->getLoopsIFunction().push_back(
946 dm_field->
tsCtx->getPostProcessIFunction().push_back(post_only);
958 dm, fe_name, method, pre_only, post_only);
964 boost::shared_ptr<MoFEM::FEMethod> method,
965 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
966 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
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);
975template <
class S,
class T0,
class T1,
class T2>
977 T1 pre_only, T2 post_only) {
978 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
980 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
982 dm_field->
tsCtx->getPreProcessIJacobian().push_back(pre_only);
985 dm_field->
tsCtx->getLoopsIJacobian().push_back(
989 dm_field->
tsCtx->getPostProcessIJacobian().push_back(post_only);
1001 pre_only, post_only);
1006 boost::shared_ptr<MoFEM::FEMethod> method,
1007 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
1008 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
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);
1016template <
class S,
class T0,
class T1,
class T2>
1018 T1 pre_only, T2 post_only) {
1019 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1021 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1023 dm_field->
tsCtx->getPreProcessMonitor().push_back(pre_only);
1025 dm_field->
tsCtx->getLoopsMonitor().push_back(
1028 dm_field->
tsCtx->getPostProcessMonitor().push_back(post_only);
1039 dm, ts, fe_name, method, pre_only, post_only);
1045 boost::shared_ptr<MoFEM::FEMethod> method,
1046 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
1047 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
1049 boost::shared_ptr<MoFEM::FEMethod>,
1050 boost::shared_ptr<MoFEM::BasicMethod>,
1051 boost::shared_ptr<MoFEM::BasicMethod>>(
1052 dm, ts, fe_name, method, pre_only, post_only);
1057 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1059 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1060 *ksp_ctx = dm_field->
kspCtx.get();
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;
1074 boost::shared_ptr<MoFEM::KspCtx> ksp_ctx) {
1075 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1077 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1078 dm_field->
kspCtx = ksp_ctx;
1083 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1086 *snes_ctx = dm_field->
snesCtx.get();
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;
1100 boost::shared_ptr<MoFEM::SnesCtx> snes_ctx) {
1101 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1103 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1112 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1114 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1123 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1125 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1131 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1133 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1139 const boost::shared_ptr<MoFEM::TsCtx> &
ts_ctx) {
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;
1148 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1150 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1156 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1158 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1166 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1168 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1171 CHKERR VecSetDM(g_ptr, dm);
1176 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1178 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1186 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1188 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1189 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1191 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->
problemName,
1193 }
else if (strcmp(dm->mattype, MATAIJ) == 0) {
1195 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->
problemName,
1197 }
else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1199 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1201 }
else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1203 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1207 "Matrix type not implemented");
1214 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1216 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1217 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1219 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->
problemName,
1221 }
else if (strcmp(dm->mattype, MATAIJ) == 0) {
1223 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->
problemName,
1225 }
else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1227 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1229 }
else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1231 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1235 "Matrix type not implemented");
1241#if PETSC_VERSION_GE(3, 7, 0)
1244#elif PETSC_VERSION_GE(3, 5, 3)
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");
1257 ierr = PetscOptionsHead(
"DMMoFEM Options");
1260 ierr = PetscOptionsBool(
"-dm_is_partitioned",
1261 "set if mesh is partitioned (works which native MOAB "
1262 "file format, i.e. h5m",
1270 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1273 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1317 PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1321 DMCtx *subdm_field =
static_cast<DMCtx *
>(subdm->data);
1326 map<std::string, boost::shared_ptr<Range>> *entity_map_row =
nullptr;
1327 map<std::string, boost::shared_ptr<Range>> *entity_map_col =
nullptr;
1337 subdm_field->
isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1344 0, subdm_field->
sIze,
1365 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1367 CHKERR VecGhostUpdateBegin(
g, INSERT_VALUES, SCATTER_FORWARD);
1373 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1376 CHKERR VecGhostUpdateEnd(
g, INSERT_VALUES, SCATTER_FORWARD);
1378 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1382 double *array_loc, *array_glob;
1383 CHKERR VecGetArray(
l, &array_loc);
1384 CHKERR VecGetArray(
g, &array_glob);
1387 cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1390 cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1395 CHKERR VecRestoreArray(
l, &array_loc);
1396 CHKERR VecRestoreArray(
g, &array_glob);
1403 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1406 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1410 double *array_loc, *array_glob;
1411 CHKERR VecGetArray(
l, &array_loc);
1412 CHKERR VecGetArray(
g, &array_glob);
1415 cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1418 cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1423 CHKERR VecRestoreArray(
l, &array_loc);
1424 CHKERR VecRestoreArray(
g, &array_glob);
1436 char ***fieldNames, IS **fields) {
1437 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1450 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
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);
1458 CHKERR PetscMalloc1(*numFields, fieldNames);
1461 CHKERR PetscMalloc1(*numFields, fields);
1464 for (
int f = 0; fit != hi_fit; fit++, f++) {
1466 CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1467 (
char **)&(*fieldNames)[f]);
1471 ->isCreateProblemFieldAndRank(
1473 fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1482 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1484 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
1492 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1494 DMCtx *dm_field =
static_cast<DMCtx *
>(dm->data);
Discrete manager interface for MoFEM.
if(!static_cast< bool >(ifstream(param_file)))
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, std::string fe_name, PetscLayout *layout)
Get finite elements layout in the problem.
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
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
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
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.
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on col.
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > snes_ctx)
Set MoFEM::SnesCtx data structure.
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > ksp_ctx)
set MoFEM::KspCtx data structure
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on row.
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problem
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
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
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, std::string fe_name)
Resolve shared entities.
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > ts_ctx)
Set MoFEM::TsCtx data structure.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
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
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
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
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
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
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.
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
PetscErrorCode DMMoFEMUnSetElement(DM dm, std::string fe_name)
unset element from dm
PetscErrorCode DMSetUp_MoFEM(DM dm)
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
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.
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
virtual MoFEMErrorCode problem_basic_method_postProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
virtual MoFEMErrorCode loop_finite_elements(const std::string problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
MoFEMErrorCode buildProblemOnDistributedMesh(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures, assuming that mesh is distributed (collective)
MoFEMErrorCode partitionGhostDofs(const std::string name, int verb=VERBOSE)
determine ghost nodes
MoFEMErrorCode buildComposedProblem(const std::string out_name, const std::vector< std::string > add_row_problems, const std::vector< std::string > add_col_problems, const bool square_matrix=true, int verb=1)
build composite problem
MoFEMErrorCode buildProblem(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
MoFEMErrorCode buildSubProblem(const std::string out_name, const std::vector< std::string > &fields_row, const std::vector< std::string > &fields_col, const std::string main_problem, const bool square_matrix=true, const map< std::string, boost::shared_ptr< Range > > *entityMapRow=nullptr, const map< std::string, boost::shared_ptr< Range > > *entityMapCol=nullptr, int verb=VERBOSE)
build sub problem
MoFEMErrorCode partitionProblem(const std::string name, int verb=VERBOSE)
partition problem dofs (collective)
MoFEMErrorCode partitionGhostDofsOnDistributedMesh(const std::string name, int verb=VERBOSE)
determine ghost nodes on distributed meshes
MoFEMErrorCode partitionFiniteElements(const std::string name, bool part_from_moab=false, int low_proc=-1, int hi_proc=-1, int verb=VERBOSE)
partition finite elements
virtual MoFEMErrorCode add_problem(const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add problem.
virtual MoFEMErrorCode modify_problem_unset_finite_element(const std::string name_problem, const std::string &fe_name)=0
unset finite element from problem, this remove entities assigned to finite element to a particular pr...
virtual MoFEMErrorCode delete_problem(const std::string name)=0
Delete problem.
virtual MoFEMErrorCode modify_problem_ref_level_add_bit(const std::string &name_problem, const BitRefLevel &bit)=0
add ref level to problem
virtual bool check_problem(const std::string name)=0
check if problem exist
virtual MoFEMErrorCode modify_problem_mask_ref_level_add_bit(const std::string &name_problem, const BitRefLevel &bit)=0
set dof mask ref level for problem
virtual MoFEMErrorCode modify_problem_add_finite_element(const std::string name_problem, const std::string &fe_name)=0
add finite element to problem, this add entities assigned to finite element to a particular problem
FTensor::Index< 'l', 3 > l
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
implementation of Data Operators for Forces and Sources
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, 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 Monitor To TS solver.
PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Set function evaluating jacobina in TS solver.
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
PetscErrorCode DMMoFEMTSSetRHSFunction(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 the right hand side function
PetscErrorCode SnesMat(SNES snes, Vec x, Mat A, Mat B, void *ctx)
This is MoFEM implementation for the left hand side (tangent matrix) evaluation in SNES solver.
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
PetscErrorCode TsSetI2Function(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, Vec F, void *ctx)
Calculation the right hand side for second order PDE in time.
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
PetscErrorCode SnesRhs(SNES snes, Vec x, Vec f, void *ctx)
This is MoFEM implementation for the right hand side (residual vector) evaluation in SNES solver.
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx)
Run over elements in the lists.
PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, PetscReal a, PetscReal aa, Mat A, Mat B, void *ctx)
Calculation Jaconian for second order PDE in time.
PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx)
Run over elenents in the list.
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
PetscErrorCode DMMoFEMTSSetRHSJacobian(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 the right hand side jacobian
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
PetscErrorCode DMMoFEMSwapDMCtx(DM dm, DM dm_swap)
Swap internal data struture.
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
boost::weak_ptr< CacheTuple > CacheTupleWeakPtr
constexpr auto field_name
Data structure to exchange data between mofem and User Loop Methods.
virtual MoFEMErrorCode problem_basic_method_preProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
virtual MPI_Comm & get_comm() const =0
PETSc Discrete Manager data structure.
std::string problemName
Problem name.
const Problem * problemPtr
pinter to problem data structure
std::vector< std::string > rowCompPrb
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
PetscBool isSquareMatrix
true if rows equals to cols
PetscBool isPartitioned
true if read mesh is on parts
PetscBool isProblemBuild
True if problem is build.
Interface * mField_ptr
MoFEM interface.
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
std::vector< std::string > colCompPrb
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
std::vector< std::string > colFields
PetscBool destroyProblem
If true destroy problem with DM.
std::vector< std::string > rowFields
std::map< std::string, boost::shared_ptr< Range > > mapTypeRow
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
std::map< std::string, boost::shared_ptr< Range > > mapTypeCol
Deprecated interface functions.
Data structure to exchange data between mofem and User Loop Methods on entities.
structure for User Loop Methods on finite elements
Section manager is used to create indexes and sections.
Interface for linear (KSP) solver.
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Matrix manager is used to build and partition problems.
keeps basic data about problem
DofIdx getNbLocalDofsRow() const
MoFEMErrorCode getNumberOfElementsByNameAndPart(MPI_Comm comm, const std::string name, PetscLayout *layout) const
Get number of finite elements by name on processors.
BitRefLevel getBitRefLevel() const
BitRefLevel getBitRefLevelMask() const
boost::shared_ptr< SubProblemData > & getSubData() const
Get main problem of sub-problem is.
DofIdx getNbGhostDofsRow() const
Problem manager is used to build and partition problems.
intrusive_ptr for managing petsc objects
Interface for nonlinear (SNES) solver.
Interface for Time Stepping (TS) solver.
base class for all interface classes
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Vector manager is used to create vectors \mofem_vectors.