v0.13.1
Loading...
Searching...
No Matches
DMMMoFEM.cpp
Go to the documentation of this file.
1
2
3// #undef PETSC_VERSION_RELEASE
4// #define PETSC_VERSION_RELEASE 1
5
6#if PETSC_VERSION_GE(3, 6, 0)
7#include <petsc/private/dmimpl.h> /*I "petscdm.h" I*/
8// #include <petsc/private/vecimpl.h> /*I "petscdm.h" I*/
9#else
10#include <petsc-private/dmimpl.h> /*I "petscdm.h" I*/
11#include <petsc-private/vecimpl.h> /*I "petscdm.h" I*/
12#endif
13
14#include <DMMoFEM.hpp>
15
16namespace MoFEM {
17
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) {
23
24 if (!LogManager::checkIfChannelExist("DMWORLD")) {
25 auto core_log = logging::core::get();
26 core_log->add_sink(
28 core_log->add_sink(
30 core_log->add_sink(
32 LogManager::setLog("DMWORLD");
33 LogManager::setLog("DMSYNC");
34 LogManager::setLog("DMSELF");
35 MOFEM_LOG_TAG("DMWORLD", "DM");
36 MOFEM_LOG_TAG("DMSYNC", "DM");
37 MOFEM_LOG_TAG("DMSELF", "DM");
38 }
39}
40
41MoFEMErrorCode DMCtx::query_interface(boost::typeindex::type_index type_index,
42 UnknownInterface **iface) const {
43 *iface = const_cast<DMCtx *>(this);
44 return 0;
45}
46
47PetscErrorCode DMRegister_MoFEM(const char sname[]) {
49 CHKERR DMRegister(sname, DMCreate_MoFEM);
51}
52
53PetscErrorCode DMSetOperators_MoFEM(DM dm) {
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}
74
75PetscErrorCode DMCreate_MoFEM(DM dm) {
76 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
78 dm->data = new DMCtx();
81}
82
83PetscErrorCode DMDestroy_MoFEM(DM dm) {
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}
117
118PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr,
119 const char problem_name[],
120 const MoFEM::BitRefLevel bit_level,
121 const MoFEM::BitRefLevel bit_mask) {
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 }
145 dm_field->problemName, bit_level);
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}
183
184PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate) {
186
187 auto *dm_field = static_cast<DMCtx *>(dm->data);
188 if (!dm->data)
189 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
190 "data structure for MoFEM not yet created");
191
192 if (static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber == 0)
193 delete static_cast<DMCtx *>(dm_duplicate->data);
194
195 dm_duplicate->data = dm->data;
196 ++(static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber);
197
199}
200
201PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]) {
203
204 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
205 if (!dm->data) {
206 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
207 "data structure for MoFEM not yet created");
208 }
209 CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
210 dm_field->problemPtr->getBitRefLevel(),
211 dm_field->problemPtr->getBitRefLevelMask());
212
213 DMCtx *subdm_field = (DMCtx *)subdm->data;
214 subdm_field->isSubDM = PETSC_TRUE;
215 subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
216 subdm_field->isPartitioned = dm_field->isPartitioned;
217 subdm_field->isSquareMatrix = PETSC_FALSE;
218 subdm->ops->setup = DMSubDMSetUp_MoFEM;
219
221}
222
223PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
224 EntityType lo_type, EntityType hi_type) {
225 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
227 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
228 if (!dm->data) {
229 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
230 "data structure for MoFEM not yet created");
231 }
232 if (!dm_field->isSubDM) {
233 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
234 }
235 dm_field->rowFields.push_back(field_name);
236 if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
237 if (!dm_field->mapTypeRow)
238 dm_field->mapTypeRow = boost::make_shared<
239 std::map<std::string, std::pair<EntityType, EntityType>>>();
240 (*dm_field->mapTypeRow)[field_name] =
241 std::pair<EntityType, EntityType>(lo_type, hi_type);
242 }
244}
245
246PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
247 EntityType lo_type, EntityType hi_type) {
248 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
250 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
251 if (!dm->data) {
252 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
253 "data structure for MoFEM not yet created");
254 }
255 if (!dm_field->isSubDM) {
256 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
257 }
258 dm_field->colFields.push_back(field_name);
259 if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
260 if (!dm_field->mapTypeCol)
261 dm_field->mapTypeCol = boost::make_shared<
262 std::map<std::string, std::pair<EntityType, EntityType>>>();
263 (*dm_field->mapTypeCol)[field_name] =
264 std::pair<EntityType, EntityType>(lo_type, hi_type);
265 }
267}
268
269PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm) {
271 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
273 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
274 *is_sub_dm = dm_field->isSubDM;
276}
277
278PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is) {
280 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
282 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
283 if (dm_field->isSubDM != PETSC_TRUE) {
284 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
285 "This DM is not created as a SubDM");
286 }
287 if (dm_field->isProblemBuild != PETSC_TRUE) {
288 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
289 }
290 boost::shared_ptr<Problem::SubProblemData> sub_data =
291 dm_field->problemPtr->getSubData();
292 CHKERR sub_data->getRowIs(is);
294}
295
296PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is) {
298 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
300 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
301 if (dm_field->isSubDM != PETSC_TRUE) {
302 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
303 "This DM is not created as a SubDM");
304 }
305 if (dm_field->isProblemBuild != PETSC_TRUE) {
306 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
307 }
308 boost::shared_ptr<Problem::SubProblemData> sub_data =
309 dm_field->problemPtr->getSubData();
310 CHKERR sub_data->getColIs(is);
312}
313
314PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]) {
315 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
317 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
318 if (!dm->data) {
319 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
320 "data structure for MoFEM not yet created");
321 }
322 if (!dm_field->isCompDM) {
323 dm_field->isCompDM = PETSC_TRUE;
324 }
325 dm_field->rowCompPrb.push_back(prb_name);
326 if (dm_field->isSquareMatrix) {
327 dm_field->colCompPrb.push_back(prb_name);
328 }
330}
331
332PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]) {
333 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
335 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
336 if (!dm->data) {
337 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
338 "data structure for MoFEM not yet created");
339 }
340 if (!dm_field->isCompDM) {
341 dm_field->isCompDM = PETSC_TRUE;
342 }
343 if (dm_field->isSquareMatrix) {
344 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
345 "No need to add problem on column when problem block structurally "
346 "symmetric");
347 }
348 dm_field->colCompPrb.push_back(prb_name);
350}
351
352PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm) {
354 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
356 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
357 *is_comp_dm = dm_field->isCompDM;
359}
360
361PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr) {
362 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
364 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
365 if (!dm->data) {
366 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
367 "data structure for MoFEM not yet created");
368 }
369 *m_field_ptr = dm_field->mField_ptr;
371}
372
373PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr) {
374 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
376 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
377 if (!dm->data) {
378 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
379 "data structure for MoFEM not yet created");
380 }
381 *problem_ptr = dm_field->problemPtr;
383}
384
385PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem) {
387 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
389 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
390 dm_field->destroyProblem = destroy_problem;
392}
393
394PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem) {
396 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
398 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
399 *destroy_problem = dm_field->destroyProblem;
401}
402
403PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem) {
405 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
407 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
408 dm_field->isSquareMatrix = square_problem;
410}
411
412PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[]) {
414 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
416 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
418 ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
420}
421
422PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[]) {
423 return DMMoFEMResolveSharedFiniteElements(dm, fe_name);
424}
425
426PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[],
427 PetscLayout *layout) {
429 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
431 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
432
433 MPI_Comm comm;
434 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
435 CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
436 layout);
438}
439
440PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem) {
443 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
445 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
446 *square_problem = dm_field->isSquareMatrix;
448}
449
450PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[]) {
451 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
453 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
455 dm_field->problemName, fe_name);
456 CHKERRG(ierr);
458}
459
460PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[]) {
461 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
463 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
465 dm_field->problemName, fe_name);
466 CHKERRG(ierr);
468}
469
470PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
471 ScatterMode scatter_mode) {
473 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
475 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
476 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
477 dm_field->problemPtr, COL, l, mode, scatter_mode);
478 CHKERRG(ierr);
480}
481
482PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
483 ScatterMode scatter_mode) {
484 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
486 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
487 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
488 dm_field->problemPtr, COL, g, mode, scatter_mode);
489 CHKERRG(ierr);
491}
492
493PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
494 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
496 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
498 dm_field->problemPtr, *method);
499 CHKERRG(ierr);
501}
502
503PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
504 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
506 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
508 dm_field->problemPtr, *method);
509 CHKERRG(ierr);
511}
512
513PetscErrorCode
514DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
515 MoFEM::FEMethod *method, int low_rank,
516 int up_rank, CacheTupleWeakPtr cache_ptr) {
518 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
520 dm_field->problemPtr, fe_name, *method, low_rank, up_rank, nullptr,
521 MF_EXIST, cache_ptr);
522 CHKERRG(ierr);
524}
525
527 DM dm, const std::string fe_name, boost::shared_ptr<MoFEM::FEMethod> method,
528 int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr) {
529 return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
530 low_rank, up_rank, cache_ptr);
531}
532
533PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
534 MoFEM::FEMethod *method,
535 CacheTupleWeakPtr cache_ptr) {
536 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
538 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
540 dm, fe_name, method, dm_field->rAnk, dm_field->rAnk, cache_ptr);
541 CHKERRG(ierr);
543}
544
545PetscErrorCode
546DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
547 boost::shared_ptr<MoFEM::FEMethod> method,
548 CacheTupleWeakPtr cache_ptr) {
549 return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get(), cache_ptr);
550}
551
552PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
553 MoFEM::DofMethod *method) {
554 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
556 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
557 ierr =
558 dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
559 *method, dm_field->rAnk, dm_field->rAnk);
560 CHKERRG(ierr);
562}
563
564template <class S, class T0, class T1, class T2>
565static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method,
566 T1 pre_only, T2 post_only) {
567 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
569 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
570 if (pre_only) {
571 dm_field->kspCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
572 }
573 if (method) {
574 dm_field->kspCtx->get_loops_to_do_Rhs().push_back(
575 PairNameFEMethodPtr(fe_name, method));
576 }
577 if (post_only) {
578 dm_field->kspCtx->get_postProcess_to_do_Rhs().push_back(post_only);
579 }
580 CHKERR DMKSPSetComputeRHS(dm, KspRhs, dm_field->kspCtx.get());
582}
583
584PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
585 MoFEM::FEMethod *method,
586 MoFEM::BasicMethod *pre_only,
587 MoFEM::BasicMethod *post_only) {
588 return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
590 dm, fe_name, method, pre_only, post_only);
591}
592
593PetscErrorCode
594DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
595 boost::shared_ptr<MoFEM::FEMethod> method,
596 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
597 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
598 return DMMoFEMKSPSetComputeRHS<const std::string,
599 boost::shared_ptr<MoFEM::FEMethod>,
600 boost::shared_ptr<MoFEM::BasicMethod>,
601 boost::shared_ptr<MoFEM::BasicMethod>>(
602 dm, fe_name, method, pre_only, post_only);
603}
604
605template <class S, class T0, class T1, class T2>
606static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method,
607 T1 pre_only, T2 post_only) {
608 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
610 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
611 if (pre_only) {
612 dm_field->kspCtx->get_preProcess_to_do_Mat().push_back(pre_only);
613 }
614 if (method) {
615 dm_field->kspCtx->get_loops_to_do_Mat().push_back(
616 PairNameFEMethodPtr(fe_name, method));
617 }
618 if (post_only) {
619 dm_field->kspCtx->get_postProcess_to_do_Mat().push_back(post_only);
620 }
621 CHKERR DMKSPSetComputeOperators(dm, KspMat, dm_field->kspCtx.get());
623}
624
625PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
626 MoFEM::FEMethod *method,
627 MoFEM::BasicMethod *pre_only,
628 MoFEM::BasicMethod *post_only) {
629 return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
632 dm, fe_name, method, pre_only, post_only);
633}
634
635PetscErrorCode
636DMMoFEMKSPSetComputeOperators(DM dm, const std::string fe_name,
637 boost::shared_ptr<MoFEM::FEMethod> method,
638 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
639 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
640 return DMMoFEMKSPSetComputeOperators<const std::string,
641 boost::shared_ptr<MoFEM::FEMethod>>(
642 dm, fe_name, method, pre_only, post_only);
643}
644
645template <class S, class T0, class T1, class T2>
646static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method,
647 T1 pre_only, T2 post_only) {
648 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
650 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
651 if (pre_only) {
652 dm_field->snesCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
653 }
654 if (method) {
655 dm_field->snesCtx->get_loops_to_do_Rhs().push_back(
656 PairNameFEMethodPtr(fe_name, method));
657 }
658 if (post_only) {
659 dm_field->snesCtx->get_postProcess_to_do_Rhs().push_back(post_only);
660 }
661 CHKERR DMSNESSetFunction(dm, SnesRhs, dm_field->snesCtx.get());
663}
664
665PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
666 MoFEM::FEMethod *method,
667 MoFEM::BasicMethod *pre_only,
668 MoFEM::BasicMethod *post_only) {
669 return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
671 dm, fe_name, method, pre_only, post_only);
672}
673
674PetscErrorCode
675DMMoFEMSNESSetFunction(DM dm, const std::string fe_name,
676 boost::shared_ptr<MoFEM::FEMethod> method,
677 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
678 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
679 return DMMoFEMSNESSetFunction<const std::string,
680 boost::shared_ptr<MoFEM::FEMethod>,
681 boost::shared_ptr<MoFEM::BasicMethod>,
682 boost::shared_ptr<MoFEM::BasicMethod>>(
683 dm, fe_name, method, pre_only, post_only);
684}
685
686template <class S, class T0, class T1, class T2>
687static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method,
688 T1 pre_only, T2 post_only) {
689 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
691 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
692 if (pre_only) {
693 dm_field->snesCtx->get_preProcess_to_do_Mat().push_back(pre_only);
694 }
695 if (method) {
696 dm_field->snesCtx->get_loops_to_do_Mat().push_back(
697 PairNameFEMethodPtr(fe_name, method));
698 }
699 if (post_only) {
700 dm_field->snesCtx->get_postProcess_to_do_Mat().push_back(post_only);
701 }
702 CHKERR DMSNESSetJacobian(dm, SnesMat, dm_field->snesCtx.get());
704}
705
706PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
707 MoFEM::FEMethod *method,
708 MoFEM::BasicMethod *pre_only,
709 MoFEM::BasicMethod *post_only) {
710 return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
712 dm, fe_name, method, pre_only, post_only);
713}
714
715PetscErrorCode
716DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
717 boost::shared_ptr<MoFEM::FEMethod> method,
718 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
719 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
720 return DMMoFEMSNESSetJacobian<const std::string,
721 boost::shared_ptr<MoFEM::FEMethod>,
722 boost::shared_ptr<MoFEM::BasicMethod>,
723 boost::shared_ptr<MoFEM::BasicMethod>>(
724 dm, fe_name, method, pre_only, post_only);
725}
726
727template <class S, class T0, class T1, class T2>
728static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method,
729 T1 pre_only, T2 post_only) {
730 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
732 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
733 if (pre_only) {
734 dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
735 }
736 if (method) {
737 dm_field->tsCtx->getLoopsIFunction().push_back(
738 PairNameFEMethodPtr(fe_name, method));
739 }
740 if (post_only) {
741 dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
742 }
743 CHKERR DMTSSetIFunction(dm, TsSetIFunction, dm_field->tsCtx.get());
745}
746
747PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
748 MoFEM::FEMethod *method,
749 MoFEM::BasicMethod *pre_only,
750 MoFEM::BasicMethod *post_only) {
751 return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
753 dm, fe_name, method, pre_only, post_only);
755}
756
757PetscErrorCode
758DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
759 boost::shared_ptr<MoFEM::FEMethod> method,
760 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
761 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
762 return DMMoFEMTSSetIFunction<const std::string,
763 boost::shared_ptr<MoFEM::FEMethod>,
764 boost::shared_ptr<MoFEM::BasicMethod>,
765 boost::shared_ptr<MoFEM::BasicMethod>>(
766 dm, fe_name, method, pre_only, post_only);
768}
769
770template <class S, class T0, class T1, class T2>
771static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method,
772 T1 pre_only, T2 post_only) {
773 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
775 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
776 if (pre_only) {
777 dm_field->tsCtx->getPreProcessIJacobian().push_back(pre_only);
778 }
779 if (method) {
780 dm_field->tsCtx->getLoopsIJacobian().push_back(
781 PairNameFEMethodPtr(fe_name, method));
782 }
783 if (post_only) {
784 dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
785 }
786 CHKERR DMTSSetIJacobian(dm, TsSetIJacobian, dm_field->tsCtx.get());
788}
789
790PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
791 MoFEM::FEMethod *method,
792 MoFEM::BasicMethod *pre_only,
793 MoFEM::BasicMethod *post_only) {
794 return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
795 MoFEM::BasicMethod *>(dm, fe_name, method,
796 pre_only, post_only);
797}
798
799PetscErrorCode
800DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
801 boost::shared_ptr<MoFEM::FEMethod> method,
802 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
803 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
804 return DMMoFEMTSSetIJacobian<const std::string,
805 boost::shared_ptr<MoFEM::FEMethod>,
806 boost::shared_ptr<MoFEM::BasicMethod>,
807 boost::shared_ptr<MoFEM::BasicMethod>>(
808 dm, fe_name, method, pre_only, post_only);
809}
810
811template <class S, class T0, class T1, class T2>
812static PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, S fe_name, T0 method,
813 T1 pre_only, T2 post_only) {
814 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
816 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
817 if (pre_only)
818 dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
819 if (method)
820 dm_field->tsCtx->getLoopsRHSFunction().push_back(
821 PairNameFEMethodPtr(fe_name, method));
822 if (post_only)
823 dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
824 CHKERR DMTSSetRHSFunction(dm, TsSetRHSFunction, dm_field->tsCtx.get());
826}
827
828PetscErrorCode
829DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name,
830 boost::shared_ptr<MoFEM::FEMethod> method,
831 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
832 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
833 return DMMoFEMTSSetRHSFunction<const std::string,
834 boost::shared_ptr<MoFEM::FEMethod>,
835 boost::shared_ptr<MoFEM::BasicMethod>,
836 boost::shared_ptr<MoFEM::BasicMethod>>(
837 dm, fe_name, method, pre_only, post_only);
839}
840
841PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const char fe_name[],
842 MoFEM::FEMethod *method,
843 MoFEM::BasicMethod *pre_only,
844 MoFEM::BasicMethod *post_only) {
845 return DMMoFEMTSSetRHSFunction<const char *, MoFEM::FEMethod *,
847 dm, fe_name, method, pre_only, post_only);
849}
850
851template <class S, class T0, class T1, class T2>
852static PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, S fe_name, T0 method,
853 T1 pre_only, T2 post_only) {
854 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
856 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
857 if (pre_only)
858 dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
859 if (method)
860 dm_field->tsCtx->getLoopsRHSFunction().push_back(
861 PairNameFEMethodPtr(fe_name, method));
862 if (post_only)
863 dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
864 CHKERR DMTSSetRHSJacobian(dm, TsSetRHSJacobian, dm_field->tsCtx.get());
866}
867
868PetscErrorCode
869DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
870 boost::shared_ptr<MoFEM::FEMethod> method,
871 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
872 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
873 return DMMoFEMTSSetRHSJacobian<const std::string,
874 boost::shared_ptr<MoFEM::FEMethod>,
875 boost::shared_ptr<MoFEM::BasicMethod>,
876 boost::shared_ptr<MoFEM::BasicMethod>>(
877 dm, fe_name, method, pre_only, post_only);
879}
880
881PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const char fe_name[],
882 MoFEM::FEMethod *method,
883 MoFEM::BasicMethod *pre_only,
884 MoFEM::BasicMethod *post_only) {
885 return DMMoFEMTSSetRHSJacobian<const char *, MoFEM::FEMethod *,
887 dm, fe_name, method, pre_only, post_only);
889}
890
891template <class S, class T0, class T1, class T2>
892static PetscErrorCode DMMoFEMTSSetI2Function(DM dm, S fe_name, T0 method,
893 T1 pre_only, T2 post_only) {
894 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
896 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
897 if (pre_only) {
898 dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
899 }
900 if (method) {
901 dm_field->tsCtx->getLoopsIFunction().push_back(
902 PairNameFEMethodPtr(fe_name, method));
903 }
904 if (post_only) {
905 dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
906 }
907 CHKERR DMTSSetI2Function(dm, TsSetI2Function, dm_field->tsCtx.get());
909}
910
911PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const char fe_name[],
912 MoFEM::FEMethod *method,
913 MoFEM::BasicMethod *pre_only,
914 MoFEM::BasicMethod *post_only) {
915 return DMMoFEMTSSetI2Function<const char *, MoFEM::FEMethod *,
917 dm, fe_name, method, pre_only, post_only);
919}
920
921PetscErrorCode
922DMMoFEMTSSetI2Function(DM dm, const std::string fe_name,
923 boost::shared_ptr<MoFEM::FEMethod> method,
924 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
925 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
926 return DMMoFEMTSSetI2Function<const std::string,
927 boost::shared_ptr<MoFEM::FEMethod>,
928 boost::shared_ptr<MoFEM::BasicMethod>,
929 boost::shared_ptr<MoFEM::BasicMethod>>(
930 dm, fe_name, method, pre_only, post_only);
932}
933
934template <class S, class T0, class T1, class T2>
935static PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, S fe_name, T0 method,
936 T1 pre_only, T2 post_only) {
937 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
939 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
940 if (pre_only) {
941 dm_field->tsCtx->getPreProcessIJacobian().push_back(pre_only);
942 }
943 if (method) {
944 dm_field->tsCtx->getLoopsIJacobian().push_back(
945 PairNameFEMethodPtr(fe_name, method));
946 }
947 if (post_only) {
948 dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
949 }
950 CHKERR DMTSSetI2Jacobian(dm, TsSetI2Jacobian, dm_field->tsCtx.get());
952}
953
954PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const char fe_name[],
955 MoFEM::FEMethod *method,
956 MoFEM::BasicMethod *pre_only,
957 MoFEM::BasicMethod *post_only) {
958 return DMMoFEMTSSetI2Jacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
959 MoFEM::BasicMethod *>(dm, fe_name, method,
960 pre_only, post_only);
961}
962
963PetscErrorCode
964DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name,
965 boost::shared_ptr<MoFEM::FEMethod> method,
966 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
967 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
968 return DMMoFEMTSSetI2Jacobian<const std::string,
969 boost::shared_ptr<MoFEM::FEMethod>,
970 boost::shared_ptr<MoFEM::BasicMethod>,
971 boost::shared_ptr<MoFEM::BasicMethod>>(
972 dm, fe_name, method, pre_only, post_only);
973}
974
975template <class S, class T0, class T1, class T2>
976static PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, S fe_name, T0 method,
977 T1 pre_only, T2 post_only) {
978 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
980 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
981 if (pre_only)
982 dm_field->tsCtx->getPreProcessMonitor().push_back(pre_only);
983 if (method)
984 dm_field->tsCtx->getLoopsMonitor().push_back(
985 PairNameFEMethodPtr(fe_name, method));
986 if (post_only)
987 dm_field->tsCtx->getPostProcessMonitor().push_back(post_only);
988 CHKERR TSMonitorSet(ts, TsMonitorSet, dm_field->tsCtx.get(), PETSC_NULL);
990}
991
992PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const char fe_name[],
993 MoFEM::FEMethod *method,
994 MoFEM::BasicMethod *pre_only,
995 MoFEM::BasicMethod *post_only) {
996 return DMMoFEMTSSetMonitor<const char *, MoFEM::FEMethod *,
998 dm, ts, fe_name, method, pre_only, post_only);
1000}
1001
1002PetscErrorCode
1003DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name,
1004 boost::shared_ptr<MoFEM::FEMethod> method,
1005 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
1006 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
1007 return DMMoFEMTSSetMonitor<const std::string,
1008 boost::shared_ptr<MoFEM::FEMethod>,
1009 boost::shared_ptr<MoFEM::BasicMethod>,
1010 boost::shared_ptr<MoFEM::BasicMethod>>(
1011 dm, ts, fe_name, method, pre_only, post_only);
1013}
1014
1015PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx) {
1016 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1018 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1019 *ksp_ctx = dm_field->kspCtx.get();
1021}
1022
1023PetscErrorCode
1024DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
1025 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1027 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1028 const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
1030}
1031
1032PetscErrorCode DMMoFEMSetKspCtx(DM dm,
1033 boost::shared_ptr<MoFEM::KspCtx> ksp_ctx) {
1034 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1036 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1037 dm_field->kspCtx = ksp_ctx;
1039}
1040
1041PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx) {
1042 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1044 DMCtx *dm_field = (DMCtx *)dm->data;
1045 *snes_ctx = dm_field->snesCtx.get();
1047}
1048
1049PetscErrorCode
1050DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
1051 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1053 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1054 const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
1056}
1057
1058PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
1059 boost::shared_ptr<MoFEM::SnesCtx> snes_ctx) {
1060 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1062 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1063 dm_field->snesCtx = snes_ctx;
1065}
1066
1067/** get if read mesh is partitioned
1068 * \ingroup dm
1069 */
1070PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned) {
1071 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1073 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1074 dm_field->isPartitioned = is_partitioned;
1076}
1077
1078/** get if read mesh is partitioned
1079 * \ingroup dm
1080 */
1081PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned) {
1082 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1084 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1085 *is_partitioned = dm_field->isPartitioned;
1087}
1088
1089PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx) {
1090 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1092 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1093 *ts_ctx = dm_field->tsCtx.get();
1095}
1096
1097PetscErrorCode DMMoFEMGetTsCtx(DM dm,
1098 const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
1099 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1101 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1102 const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1104}
1105
1106PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> ts_ctx) {
1107 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1109 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1110 dm_field->tsCtx = ts_ctx;
1112}
1113
1114PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g) {
1115 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1117 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1118 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1119 dm_field->problemName, COL, g);
1120 CHKERR VecSetDM(*g, dm);
1122}
1123
1124PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj<Vec> &g_ptr) {
1125 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1127 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1128 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1129 dm_field->problemName, COL, g_ptr);
1130 CHKERR VecSetDM(g_ptr, dm);
1132}
1133
1134PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l) {
1135 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1137 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1138 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1139 dm_field->problemName, COL, l);
1140 CHKERR VecSetDM(*l, dm);
1142}
1143
1144PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M) {
1145 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1147 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1148 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1150 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1151 M);
1152 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1154 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1155 M);
1156 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1158 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1159 dm_field->problemName, M);
1160 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1162 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1163 dm_field->problemName, M);
1164 } else {
1165 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1166 "Matrix type not implemented");
1167 }
1168 CHKERR MatSetDM(*M, dm);
1170}
1171
1173 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1175 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1176 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1178 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1179 M);
1180 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1182 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1183 M);
1184 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1186 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1187 dm_field->problemName, M);
1188 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1190 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1191 dm_field->problemName, M);
1192 } else {
1193 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1194 "Matrix type not implemented");
1195 }
1196 CHKERR MatSetDM(M, dm);
1198}
1199
1200#if PETSC_VERSION_GE(3, 7, 0)
1201PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
1202 DM dm) {
1203#elif PETSC_VERSION_GE(3, 5, 3)
1204PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm) {
1205#else
1206PetscErrorCode DMSetFromOptions_MoFEM(DM dm) {
1207#endif
1208
1209 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1211 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1212#if PETSC_VERSION_GE(3, 5, 3)
1213 ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1214 CHKERRG(ierr);
1215#else
1216 ierr = PetscOptionsHead("DMMoFEM Options");
1217 CHKERRG(ierr);
1218#endif
1219 ierr = PetscOptionsBool("-dm_is_partitioned",
1220 "set if mesh is partitioned (works which native MOAB "
1221 "file format, i.e. h5m",
1222 "DMSetUp", dm_field->isPartitioned,
1223 &dm_field->isPartitioned, NULL);
1224 CHKERRG(ierr);
1226}
1227
1228PetscErrorCode DMSetUp_MoFEM(DM dm) {
1229 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1230 ProblemsManager *prb_mng_ptr;
1232 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1233 CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1234
1235 if (dm_field->isCompDM) {
1236 // It is composite probelm
1237 CHKERR prb_mng_ptr->buildComposedProblem(
1238 dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1239 dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1240 } else {
1241 if (dm_field->isPartitioned) {
1243 dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1244 dm_field->verbosity);
1245 } else {
1246 CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1247 dm_field->isSquareMatrix == PETSC_TRUE,
1248 dm_field->verbosity);
1249 CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1250 dm_field->verbosity);
1251 }
1252 }
1253
1254 // Partition finite elements
1255 if (dm_field->isPartitioned) {
1256 CHKERR prb_mng_ptr->partitionFiniteElements(
1257 dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1259 dm_field->problemName, dm_field->verbosity);
1260 } else {
1261 // partition finite elemnets
1262 CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1263 -1, -1, dm_field->verbosity);
1264 // Get ghost DOFs
1265 CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1266 dm_field->verbosity);
1267 }
1268
1269 // Set flag that problem is build and partitioned
1270 dm_field->isProblemBuild = PETSC_TRUE;
1271
1273}
1274
1275PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm) {
1276 PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1277 ProblemsManager *prb_mng_ptr;
1279
1280 DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1281
1282 // build sub dm problem
1283 CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1284
1285 map<std::string, std::pair<EntityType, EntityType>> *entity_map_row = nullptr;
1286 map<std::string, std::pair<EntityType, EntityType>> *entity_map_col = nullptr;
1287
1288 if (subdm_field->mapTypeRow)
1289 entity_map_row = subdm_field->mapTypeRow.get();
1290 if (subdm_field->mapTypeCol)
1291 entity_map_row = subdm_field->mapTypeCol.get();
1292
1293 CHKERR prb_mng_ptr->buildSubProblem(
1294 subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1295 subdm_field->problemMainOfSubPtr->getName(),
1296 subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1297 subdm_field->verbosity);
1298
1299 // partition problem
1300 subdm_field->isPartitioned = subdm_field->isPartitioned;
1301 if (subdm_field->isPartitioned) {
1302 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1303 0, subdm_field->sIze,
1304 subdm_field->verbosity);
1305 // set ghost nodes
1307 subdm_field->problemName, subdm_field->verbosity);
1308 } else {
1309 // partition finite elements
1310 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1311 -1, -1, subdm_field->verbosity);
1312 // set ghost nodes
1313 CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1314 subdm_field->verbosity);
1315 }
1316
1317 subdm_field->isProblemBuild = PETSC_TRUE;
1318
1320}
1321
1322PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec g, InsertMode mode,
1323 Vec l) {
1324 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1326 CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1328}
1329
1330PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec g, InsertMode mode, Vec l) {
1332 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1334
1335 CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1336
1337 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1338 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1339 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1340
1341 double *array_loc, *array_glob;
1342 CHKERR VecGetArray(l, &array_loc);
1343 CHKERR VecGetArray(g, &array_glob);
1344 switch (mode) {
1345 case INSERT_VALUES:
1346 cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1347 break;
1348 case ADD_VALUES:
1349 cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1350 break;
1351 default:
1352 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1353 }
1354 CHKERR VecRestoreArray(l, &array_loc);
1355 CHKERR VecRestoreArray(g, &array_glob);
1357}
1358
1359PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM dm, Vec l, InsertMode mode,
1360 Vec g) {
1361
1362 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1364
1365 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1366 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1367 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1368
1369 double *array_loc, *array_glob;
1370 CHKERR VecGetArray(l, &array_loc);
1371 CHKERR VecGetArray(g, &array_glob);
1372 switch (mode) {
1373 case INSERT_VALUES:
1374 cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1375 break;
1376 case ADD_VALUES:
1377 cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1378 break;
1379 default:
1380 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1381 }
1382 CHKERR VecRestoreArray(l, &array_loc);
1383 CHKERR VecRestoreArray(g, &array_glob);
1384
1386}
1387
1388PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec l, InsertMode mode, Vec g) {
1389 //
1392}
1393
1394PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
1395 char ***fieldNames, IS **fields) {
1396 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1398
1399 if (numFields) {
1400 *numFields = 0;
1401 }
1402 if (fieldNames) {
1403 *fieldNames = NULL;
1404 }
1405 if (fields) {
1406 *fields = NULL;
1407 }
1408
1409 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1410 auto fields_ptr = dm_field->mField_ptr->get_fields();
1411 Field_multiIndex::iterator fit, hi_fit;
1412 fit = fields_ptr->begin();
1413 hi_fit = fields_ptr->end();
1414 *numFields = std::distance(fit, hi_fit);
1415
1416 if (fieldNames) {
1417 CHKERR PetscMalloc1(*numFields, fieldNames);
1418 }
1419 if (fields) {
1420 CHKERR PetscMalloc1(*numFields, fields);
1421 }
1422
1423 for (int f = 0; fit != hi_fit; fit++, f++) {
1424 if (fieldNames) {
1425 CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1426 (char **)&(*fieldNames)[f]);
1427 }
1428 if (fields) {
1430 ->isCreateProblemFieldAndRank(
1431 dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1432 fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1433 }
1434 }
1435
1437}
1438
1439PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
1440 IS *is) {
1441 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1443 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1445 ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1446 field_name, 0, 1000, is);
1448}
1449
1450PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb) {
1451 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1453 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1454 dm_field->verbosity = verb;
1456}
1457
1458} // namespace MoFEM
static Index< 'M', 3 > M
Discrete manager interface for MoFEM.
@ VERBOSE
Definition: definitions.h:209
RowColData
RowColData.
Definition: definitions.h:123
@ COL
Definition: definitions.h:123
@ ROW
Definition: definitions.h:123
@ MF_EXIST
Definition: definitions.h:100
@ MF_EXCL
Definition: definitions.h:99
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
@ MOFEM_INVALID_DATA
Definition: definitions.h:36
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#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
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1322
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1330
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:412
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:1070
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
Definition: DMMMoFEM.cpp:201
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMMoFEM.cpp:1081
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[], PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMMoFEM.cpp:426
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMMoFEM.cpp:665
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMMoFEM.cpp:1134
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition: DMMMoFEM.cpp:403
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:747
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:223
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:83
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:118
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:450
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:503
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:470
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on col.
Definition: DMMMoFEM.cpp:332
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1359
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:373
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:75
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > snes_ctx)
Set MoFEM::SnesCtx data structure.
Definition: DMMMoFEM.cpp:1058
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > ksp_ctx)
set MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:1032
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMMoFEM.cpp:269
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on row.
Definition: DMMMoFEM.cpp:314
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMMoFEM.cpp:706
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1144
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1206
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMMoFEM.cpp:1089
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:246
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:47
PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[])
unset element from dm
Definition: DMMMoFEM.cpp:460
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:1015
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problem
Definition: DMMMoFEM.cpp:440
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMMoFEM.cpp:552
PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set KSP right hand side evaluation function
Definition: DMMMoFEM.cpp:584
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > ts_ctx)
Set MoFEM::TsCtx data structure.
Definition: DMMMoFEM.cpp:1106
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:533
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMMoFEM.cpp:1439
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:800
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:361
PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:964
PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:922
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1114
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMMoFEM.cpp:1041
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:482
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition: DMMMoFEM.cpp:352
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1275
PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set KSP operators and push mofem finite element methods.
Definition: DMMMoFEM.cpp:625
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1394
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1228
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1388
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:514
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:493
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:53
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.
Definition: LogManager.cpp:364
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:332
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 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, std::pair< EntityType, EntityType > > *entityMapRow=nullptr, const map< std::string, std::pair< EntityType, EntityType > > *entityMapCol=nullptr, int verb=VERBOSE)
build sub problem
MoFEMErrorCode buildProblem(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
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 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 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
virtual MoFEMErrorCode add_problem(const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add problem.
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
FTensor::Index< 'l', 3 > l
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
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.
Definition: DMMMoFEM.cpp:1003
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.
Definition: TsCtx.cpp:129
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:219
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
Definition: DMMMoFEM.cpp:829
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.
Definition: SnesCtx.cpp:136
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:385
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMMoFEM.cpp:394
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.
Definition: TsCtx.cpp:563
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:278
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.
Definition: SnesCtx.cpp:27
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:278
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:25
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:296
PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx)
Run over elements in the lists.
Definition: KspCtx.cpp:21
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.
Definition: TsCtx.cpp:466
PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx)
Run over elenents in the list.
Definition: KspCtx.cpp:85
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:377
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
Definition: DMMMoFEM.cpp:869
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition: DMMMoFEM.cpp:1450
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
Definition: DMMMoFEM.cpp:184
DEPRECATED PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[])
Definition: DMMMoFEM.cpp:422
boost::weak_ptr< CacheTuple > CacheTupleWeakPtr
constexpr auto field_name
constexpr double g
Data structure to exchange data between mofem and User Loop Methods.
Managing BitRefLevels.
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.
Definition: DMMoFEM.hpp:892
PetscBool isCompDM
Definition: DMMoFEM.hpp:916
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:899
const Problem * problemPtr
pinter to problem data structure
Definition: DMMoFEM.hpp:908
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:917
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition: DMMoFEM.hpp:933
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition: DMMoFEM.hpp:934
PetscBool isSquareMatrix
true if rows equals to cols
Definition: DMMoFEM.hpp:903
PetscBool isPartitioned
true if read mesh is on parts
Definition: DMMoFEM.hpp:902
PetscBool isSubDM
Definition: DMMoFEM.hpp:911
int referenceNumber
Definition: DMMoFEM.hpp:930
int verbosity
verbosity
Definition: DMMoFEM.hpp:929
PetscBool isProblemBuild
True if problem is build.
Definition: DMMoFEM.hpp:898
Interface * mField_ptr
MoFEM interface.
Definition: DMMoFEM.hpp:897
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition: DMMoFEM.hpp:932
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeRow
Definition: DMMoFEM.hpp:920
std::vector< std::string > colCompPrb
Definition: DMMoFEM.hpp:918
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: DMMMoFEM.cpp:41
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:913
PetscBool destroyProblem
If true destroy problem with DM.
Definition: DMMoFEM.hpp:924
std::vector< std::string > rowFields
Definition: DMMoFEM.hpp:912
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
Definition: DMMoFEM.hpp:914
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeCol
Definition: DMMoFEM.hpp:922
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.
Definition: ISManager.hpp:23
Interface for linear (KSP) solver.
Definition: KspCtx.hpp:15
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:327
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:374
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:319
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.
Definition: SnesCtx.hpp:15
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:15
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.
Definition: VecManager.hpp:23