v0.14.0
Loading...
Searching...
No Matches
DMMoFEM.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 if (!dm->data)
188 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
189 "data structure for MoFEM not yet created");
190
191 if (static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber == 0)
192 delete static_cast<DMCtx *>(dm_duplicate->data);
193
194 dm_duplicate->data = dm->data;
195 ++(static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber);
196
198}
199
200PetscErrorCode DMMoFEMSwapDMCtx(DM dm, DM dm_swap) {
202 if (!dm->data)
203 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
204 "data structure for MoFEM not yet created on dm");
205 if (!dm_swap->data)
206 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
207 "data structure for MoFEM not yet created on swap dm");
208
209 auto *dm_field = static_cast<DMCtx *>(dm->data);
210 auto *dm_field_swap = static_cast<DMCtx *>(dm_swap->data);
211
212 auto tmp_field = dm_field;
213 dm_field = dm_field_swap;
214 dm_field_swap = tmp_field;
215
217}
218
219PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]) {
221
222 if (!dm->data) {
223 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
224 "data structure for MoFEM not yet created");
225 }
226 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
227
228 CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
229 dm_field->problemPtr->getBitRefLevel(),
230 dm_field->problemPtr->getBitRefLevelMask());
231
232 DMCtx *subdm_field = (DMCtx *)subdm->data;
233 subdm_field->isSubDM = PETSC_TRUE;
234 subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
235 subdm_field->isPartitioned = dm_field->isPartitioned;
236 subdm_field->isSquareMatrix = PETSC_FALSE;
237 subdm->ops->setup = DMSubDMSetUp_MoFEM;
238
240}
241
242PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[]) {
243 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
245 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
246 if (!dm->data) {
247 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
248 "data structure for MoFEM not yet created");
249 }
250 if (!dm_field->isSubDM) {
251 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
252 }
253 dm_field->rowFields.push_back(field_name);
254 dm_field->mapTypeRow.erase(field_name);
256}
257
258PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
259 boost::shared_ptr<Range> r_ptr) {
260 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
262 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
263 if (!dm->data) {
264 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
265 "data structure for MoFEM not yet created");
266 }
267 if (!dm_field->isSubDM) {
268 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
269 }
270 dm_field->rowFields.push_back(field_name);
271 dm_field->mapTypeRow[field_name] = r_ptr;
273}
274
275PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[]) {
276 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
278 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
279 if (!dm->data) {
280 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
281 "data structure for MoFEM not yet created");
282 }
283 if (!dm_field->isSubDM) {
284 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
285 }
286 dm_field->colFields.push_back(field_name);
287 dm_field->mapTypeCol.erase(field_name);
289}
290
291PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
292 boost::shared_ptr<Range> r_ptr) {
293 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
295 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
296 if (!dm->data) {
297 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
298 "data structure for MoFEM not yet created");
299 }
300 if (!dm_field->isSubDM) {
301 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
302 }
303 dm_field->colFields.push_back(field_name);
304 dm_field->mapTypeCol[field_name] = r_ptr;
306}
307
308PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm) {
310 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
312 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
313 *is_sub_dm = dm_field->isSubDM;
315}
316
317PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is) {
319 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
321 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
322 if (dm_field->isSubDM != PETSC_TRUE) {
323 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
324 "This DM is not created as a SubDM");
325 }
326 if (dm_field->isProblemBuild != PETSC_TRUE) {
327 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
328 }
329 boost::shared_ptr<Problem::SubProblemData> sub_data =
330 dm_field->problemPtr->getSubData();
331 CHKERR sub_data->getRowIs(is);
333}
334
335PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is) {
337 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
339 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
340 if (dm_field->isSubDM != PETSC_TRUE) {
341 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
342 "This DM is not created as a SubDM");
343 }
344 if (dm_field->isProblemBuild != PETSC_TRUE) {
345 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
346 }
347 boost::shared_ptr<Problem::SubProblemData> sub_data =
348 dm_field->problemPtr->getSubData();
349 CHKERR sub_data->getColIs(is);
351}
352
353PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]) {
354 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
356 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
357 if (!dm->data) {
358 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
359 "data structure for MoFEM not yet created");
360 }
361 if (!dm_field->isCompDM) {
362 dm_field->isCompDM = PETSC_TRUE;
363 }
364 dm_field->rowCompPrb.push_back(prb_name);
365 if (dm_field->isSquareMatrix) {
366 dm_field->colCompPrb.push_back(prb_name);
367 }
369}
370
371PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]) {
372 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
374 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
375 if (!dm->data) {
376 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
377 "data structure for MoFEM not yet created");
378 }
379 if (!dm_field->isCompDM) {
380 dm_field->isCompDM = PETSC_TRUE;
381 }
382 if (dm_field->isSquareMatrix) {
383 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
384 "No need to add problem on column when problem block structurally "
385 "symmetric");
386 }
387 dm_field->colCompPrb.push_back(prb_name);
389}
390
391PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm) {
393 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
395 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
396 *is_comp_dm = dm_field->isCompDM;
398}
399
400PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr) {
401 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
403 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
404 if (!dm->data) {
405 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
406 "data structure for MoFEM not yet created");
407 }
408 *m_field_ptr = dm_field->mField_ptr;
410}
411
412PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr) {
413 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
415 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
416 if (!dm->data) {
417 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
418 "data structure for MoFEM not yet created");
419 }
420 *problem_ptr = dm_field->problemPtr;
422}
423
424PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem) {
426 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
428 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
429 dm_field->destroyProblem = destroy_problem;
431}
432
433PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem) {
435 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
437 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
438 *destroy_problem = dm_field->destroyProblem;
440}
441
442PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem) {
443 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
445 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
446 dm_field->isSquareMatrix = square_problem;
448}
449
450PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, std::string fe_name) {
451 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
453 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
455 ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
457}
458
459PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, std::string fe_name,
460 PetscLayout *layout) {
462 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
464 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
465
466 MPI_Comm comm;
467 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
468 CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
469 layout);
471}
472
473PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem) {
476 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
478 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
479 *square_problem = dm_field->isSquareMatrix;
481}
482
483PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name) {
484 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
486 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
488 dm_field->problemName, fe_name);
490}
491
492PetscErrorCode DMMoFEMAddElement(DM dm, std::vector<std::string> fe_name) {
494 for (auto fe : fe_name) {
496 }
498}
499
500PetscErrorCode DMMoFEMUnSetElement(DM dm, std::string fe_name) {
501 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
503 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
505 dm_field->problemName, fe_name);
507}
508
509PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
510 ScatterMode scatter_mode) {
512 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
514 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
515 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
516 dm_field->problemPtr, COL, l, mode, scatter_mode);
517 CHKERRG(ierr);
519}
520
521PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
522 ScatterMode scatter_mode) {
523 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
525 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
526 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
527 dm_field->problemPtr, COL, g, mode, scatter_mode);
528 CHKERRG(ierr);
530}
531
532PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
533 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
535 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
537 dm_field->problemPtr, *method);
538 CHKERRG(ierr);
540}
541
542PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
543 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
545 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
547 dm_field->problemPtr, *method);
548 CHKERRG(ierr);
550}
551
552PetscErrorCode
553DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
554 MoFEM::FEMethod *method, int low_rank,
555 int up_rank, CacheTupleWeakPtr cache_ptr) {
557 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
559 dm_field->problemPtr, fe_name, *method, low_rank, up_rank, nullptr,
560 MF_EXIST, cache_ptr);
561 CHKERRG(ierr);
563}
564
566 DM dm, const std::string fe_name, boost::shared_ptr<MoFEM::FEMethod> method,
567 int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr) {
568 return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
569 low_rank, up_rank, cache_ptr);
570}
571
572PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
573 MoFEM::FEMethod *method,
574 CacheTupleWeakPtr cache_ptr) {
575 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
577 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
579 dm, fe_name, method, dm_field->rAnk, dm_field->rAnk, cache_ptr);
580 CHKERRG(ierr);
582}
583
584PetscErrorCode
585DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
586 boost::shared_ptr<MoFEM::FEMethod> method,
587 CacheTupleWeakPtr cache_ptr) {
588 return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get(), cache_ptr);
589}
590
591PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
592 MoFEM::DofMethod *method) {
593 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
595 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
596 ierr =
597 dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
598 *method, dm_field->rAnk, dm_field->rAnk);
599 CHKERRG(ierr);
601}
602
603template <class S, class T0, class T1, class T2>
604static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method,
605 T1 pre_only, T2 post_only) {
606 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
608 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
609 if (pre_only) {
610 dm_field->kspCtx->getPreProcComputeRhs().push_back(pre_only);
611 }
612 if (method) {
613 dm_field->kspCtx->getComputeRhs().push_back(
614 PairNameFEMethodPtr(fe_name, method));
615 }
616 if (post_only) {
617 dm_field->kspCtx->getPostProcComputeRhs().push_back(post_only);
618 }
619 CHKERR DMKSPSetComputeRHS(dm, KspRhs, dm_field->kspCtx.get());
621}
622
623PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
624 MoFEM::FEMethod *method,
625 MoFEM::BasicMethod *pre_only,
626 MoFEM::BasicMethod *post_only) {
627 return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
629 dm, fe_name, method, pre_only, post_only);
630}
631
632PetscErrorCode
633DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
634 boost::shared_ptr<MoFEM::FEMethod> method,
635 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
636 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
637 return DMMoFEMKSPSetComputeRHS<const std::string,
638 boost::shared_ptr<MoFEM::FEMethod>,
639 boost::shared_ptr<MoFEM::BasicMethod>,
640 boost::shared_ptr<MoFEM::BasicMethod>>(
641 dm, fe_name, method, pre_only, post_only);
642}
643
644template <class S, class T0, class T1, class T2>
645static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method,
646 T1 pre_only, T2 post_only) {
647 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
649 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
650 if (pre_only) {
651 dm_field->kspCtx->getPreProcSetOperators().push_back(pre_only);
652 }
653 if (method) {
654 dm_field->kspCtx->getSetOperators().push_back(
655 PairNameFEMethodPtr(fe_name, method));
656 }
657 if (post_only) {
658 dm_field->kspCtx->getPostProcSetOperators().push_back(post_only);
659 }
660 CHKERR DMKSPSetComputeOperators(dm, KspMat, dm_field->kspCtx.get());
662}
663
664PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
665 MoFEM::FEMethod *method,
666 MoFEM::BasicMethod *pre_only,
667 MoFEM::BasicMethod *post_only) {
668 return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
671 dm, fe_name, method, pre_only, post_only);
672}
673
674PetscErrorCode
675DMMoFEMKSPSetComputeOperators(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 DMMoFEMKSPSetComputeOperators<const std::string,
680 boost::shared_ptr<MoFEM::FEMethod>>(
681 dm, fe_name, method, pre_only, post_only);
682}
683
684template <class S, class T0, class T1, class T2>
685static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method,
686 T1 pre_only, T2 post_only) {
687 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
689 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
690 if (pre_only) {
691 dm_field->snesCtx->getPreProcComputeRhs().push_back(pre_only);
692 }
693 if (method) {
694 dm_field->snesCtx->getComputeRhs().push_back(
695 PairNameFEMethodPtr(fe_name, method));
696 }
697 if (post_only) {
698 dm_field->snesCtx->getPostProcComputeRhs().push_back(post_only);
699 }
700 CHKERR DMSNESSetFunction(dm, SnesRhs, dm_field->snesCtx.get());
702}
703
704PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
705 MoFEM::FEMethod *method,
706 MoFEM::BasicMethod *pre_only,
707 MoFEM::BasicMethod *post_only) {
708 return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
710 dm, fe_name, method, pre_only, post_only);
711}
712
713PetscErrorCode
714DMMoFEMSNESSetFunction(DM dm, const std::string fe_name,
715 boost::shared_ptr<MoFEM::FEMethod> method,
716 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
717 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
718 return DMMoFEMSNESSetFunction<const std::string,
719 boost::shared_ptr<MoFEM::FEMethod>,
720 boost::shared_ptr<MoFEM::BasicMethod>,
721 boost::shared_ptr<MoFEM::BasicMethod>>(
722 dm, fe_name, method, pre_only, post_only);
723}
724
725template <class S, class T0, class T1, class T2>
726static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method,
727 T1 pre_only, T2 post_only) {
728 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
730 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
731 if (pre_only) {
732 dm_field->snesCtx->getPreProcSetOperators().push_back(pre_only);
733 }
734 if (method) {
735 dm_field->snesCtx->getSetOperators().push_back(
736 PairNameFEMethodPtr(fe_name, method));
737 }
738 if (post_only) {
739 dm_field->snesCtx->getPostProcSetOperators().push_back(post_only);
740 }
741 CHKERR DMSNESSetJacobian(dm, SnesMat, dm_field->snesCtx.get());
743}
744
745PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
746 MoFEM::FEMethod *method,
747 MoFEM::BasicMethod *pre_only,
748 MoFEM::BasicMethod *post_only) {
749 return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
751 dm, fe_name, method, pre_only, post_only);
752}
753
754PetscErrorCode
755DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
756 boost::shared_ptr<MoFEM::FEMethod> method,
757 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
758 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
759 return DMMoFEMSNESSetJacobian<const std::string,
760 boost::shared_ptr<MoFEM::FEMethod>,
761 boost::shared_ptr<MoFEM::BasicMethod>,
762 boost::shared_ptr<MoFEM::BasicMethod>>(
763 dm, fe_name, method, pre_only, post_only);
764}
765
766template <class S, class T0, class T1, class T2>
767static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method,
768 T1 pre_only, T2 post_only) {
769 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
771 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
772 if (pre_only) {
773 dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
774 }
775 if (method) {
776 dm_field->tsCtx->getLoopsIFunction().push_back(
777 PairNameFEMethodPtr(fe_name, method));
778 }
779 if (post_only) {
780 dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
781 }
782 CHKERR DMTSSetIFunction(dm, TsSetIFunction, dm_field->tsCtx.get());
784}
785
786PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
787 MoFEM::FEMethod *method,
788 MoFEM::BasicMethod *pre_only,
789 MoFEM::BasicMethod *post_only) {
790 return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
792 dm, fe_name, method, pre_only, post_only);
794}
795
796PetscErrorCode
797DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
798 boost::shared_ptr<MoFEM::FEMethod> method,
799 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
800 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
801 return DMMoFEMTSSetIFunction<const std::string,
802 boost::shared_ptr<MoFEM::FEMethod>,
803 boost::shared_ptr<MoFEM::BasicMethod>,
804 boost::shared_ptr<MoFEM::BasicMethod>>(
805 dm, fe_name, method, pre_only, post_only);
807}
808
809template <class S, class T0, class T1, class T2>
810static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method,
811 T1 pre_only, T2 post_only) {
812 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
814 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
815 if (pre_only) {
816 dm_field->tsCtx->getPreProcessIJacobian().push_back(pre_only);
817 }
818 if (method) {
819 dm_field->tsCtx->getLoopsIJacobian().push_back(
820 PairNameFEMethodPtr(fe_name, method));
821 }
822 if (post_only) {
823 dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
824 }
825 CHKERR DMTSSetIJacobian(dm, TsSetIJacobian, dm_field->tsCtx.get());
827}
828
829PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
830 MoFEM::FEMethod *method,
831 MoFEM::BasicMethod *pre_only,
832 MoFEM::BasicMethod *post_only) {
833 return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
834 MoFEM::BasicMethod *>(dm, fe_name, method,
835 pre_only, post_only);
836}
837
838PetscErrorCode
839DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
840 boost::shared_ptr<MoFEM::FEMethod> method,
841 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
842 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
843 return DMMoFEMTSSetIJacobian<const std::string,
844 boost::shared_ptr<MoFEM::FEMethod>,
845 boost::shared_ptr<MoFEM::BasicMethod>,
846 boost::shared_ptr<MoFEM::BasicMethod>>(
847 dm, fe_name, method, pre_only, post_only);
848}
849
850template <class S, class T0, class T1, class T2>
851static PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, S fe_name, T0 method,
852 T1 pre_only, T2 post_only) {
853 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
855 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
856 if (pre_only)
857 dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
858 if (method)
859 dm_field->tsCtx->getLoopsRHSFunction().push_back(
860 PairNameFEMethodPtr(fe_name, method));
861 if (post_only)
862 dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
863 CHKERR DMTSSetRHSFunction(dm, TsSetRHSFunction, dm_field->tsCtx.get());
865}
866
867PetscErrorCode
868DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name,
869 boost::shared_ptr<MoFEM::FEMethod> method,
870 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
871 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
872 return DMMoFEMTSSetRHSFunction<const std::string,
873 boost::shared_ptr<MoFEM::FEMethod>,
874 boost::shared_ptr<MoFEM::BasicMethod>,
875 boost::shared_ptr<MoFEM::BasicMethod>>(
876 dm, fe_name, method, pre_only, post_only);
878}
879
880PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const char fe_name[],
881 MoFEM::FEMethod *method,
882 MoFEM::BasicMethod *pre_only,
883 MoFEM::BasicMethod *post_only) {
884 return DMMoFEMTSSetRHSFunction<const char *, MoFEM::FEMethod *,
886 dm, fe_name, method, pre_only, post_only);
888}
889
890template <class S, class T0, class T1, class T2>
891static PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, S fe_name, T0 method,
892 T1 pre_only, T2 post_only) {
893 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
895 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
896 if (pre_only)
897 dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
898 if (method)
899 dm_field->tsCtx->getLoopsRHSFunction().push_back(
900 PairNameFEMethodPtr(fe_name, method));
901 if (post_only)
902 dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
903 CHKERR DMTSSetRHSJacobian(dm, TsSetRHSJacobian, dm_field->tsCtx.get());
905}
906
907PetscErrorCode
908DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
909 boost::shared_ptr<MoFEM::FEMethod> method,
910 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
911 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
912 return DMMoFEMTSSetRHSJacobian<const std::string,
913 boost::shared_ptr<MoFEM::FEMethod>,
914 boost::shared_ptr<MoFEM::BasicMethod>,
915 boost::shared_ptr<MoFEM::BasicMethod>>(
916 dm, fe_name, method, pre_only, post_only);
918}
919
920PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const char fe_name[],
921 MoFEM::FEMethod *method,
922 MoFEM::BasicMethod *pre_only,
923 MoFEM::BasicMethod *post_only) {
924 return DMMoFEMTSSetRHSJacobian<const char *, MoFEM::FEMethod *,
926 dm, fe_name, method, pre_only, post_only);
928}
929
930template <class S, class T0, class T1, class T2>
931static PetscErrorCode DMMoFEMTSSetI2Function(DM dm, S fe_name, T0 method,
932 T1 pre_only, T2 post_only) {
933 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
935 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
936 if (pre_only) {
937 dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
938 }
939 if (method) {
940 dm_field->tsCtx->getLoopsIFunction().push_back(
941 PairNameFEMethodPtr(fe_name, method));
942 }
943 if (post_only) {
944 dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
945 }
946 CHKERR DMTSSetI2Function(dm, TsSetI2Function, dm_field->tsCtx.get());
948}
949
950PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const char fe_name[],
951 MoFEM::FEMethod *method,
952 MoFEM::BasicMethod *pre_only,
953 MoFEM::BasicMethod *post_only) {
954 return DMMoFEMTSSetI2Function<const char *, MoFEM::FEMethod *,
956 dm, fe_name, method, pre_only, post_only);
958}
959
960PetscErrorCode
961DMMoFEMTSSetI2Function(DM dm, const std::string fe_name,
962 boost::shared_ptr<MoFEM::FEMethod> method,
963 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
964 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
965 return DMMoFEMTSSetI2Function<const std::string,
966 boost::shared_ptr<MoFEM::FEMethod>,
967 boost::shared_ptr<MoFEM::BasicMethod>,
968 boost::shared_ptr<MoFEM::BasicMethod>>(
969 dm, fe_name, method, pre_only, post_only);
971}
972
973template <class S, class T0, class T1, class T2>
974static PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, S fe_name, T0 method,
975 T1 pre_only, T2 post_only) {
976 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
978 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
979 if (pre_only) {
980 dm_field->tsCtx->getPreProcessIJacobian().push_back(pre_only);
981 }
982 if (method) {
983 dm_field->tsCtx->getLoopsIJacobian().push_back(
984 PairNameFEMethodPtr(fe_name, method));
985 }
986 if (post_only) {
987 dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
988 }
989 CHKERR DMTSSetI2Jacobian(dm, TsSetI2Jacobian, dm_field->tsCtx.get());
991}
992
993PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const char fe_name[],
994 MoFEM::FEMethod *method,
995 MoFEM::BasicMethod *pre_only,
996 MoFEM::BasicMethod *post_only) {
997 return DMMoFEMTSSetI2Jacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
998 MoFEM::BasicMethod *>(dm, fe_name, method,
999 pre_only, post_only);
1000}
1001
1002PetscErrorCode
1003DMMoFEMTSSetI2Jacobian(DM dm, 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 DMMoFEMTSSetI2Jacobian<const std::string,
1008 boost::shared_ptr<MoFEM::FEMethod>,
1009 boost::shared_ptr<MoFEM::BasicMethod>,
1010 boost::shared_ptr<MoFEM::BasicMethod>>(
1011 dm, fe_name, method, pre_only, post_only);
1012}
1013
1014template <class S, class T0, class T1, class T2>
1015static PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, S fe_name, T0 method,
1016 T1 pre_only, T2 post_only) {
1017 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1019 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1020 if (pre_only)
1021 dm_field->tsCtx->getPreProcessMonitor().push_back(pre_only);
1022 if (method)
1023 dm_field->tsCtx->getLoopsMonitor().push_back(
1024 PairNameFEMethodPtr(fe_name, method));
1025 if (post_only)
1026 dm_field->tsCtx->getPostProcessMonitor().push_back(post_only);
1027 CHKERR TSMonitorSet(ts, TsMonitorSet, dm_field->tsCtx.get(), PETSC_NULL);
1029}
1030
1031PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const char fe_name[],
1032 MoFEM::FEMethod *method,
1033 MoFEM::BasicMethod *pre_only,
1034 MoFEM::BasicMethod *post_only) {
1035 return DMMoFEMTSSetMonitor<const char *, MoFEM::FEMethod *,
1037 dm, ts, fe_name, method, pre_only, post_only);
1039}
1040
1041PetscErrorCode
1042DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name,
1043 boost::shared_ptr<MoFEM::FEMethod> method,
1044 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
1045 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
1046 return DMMoFEMTSSetMonitor<const std::string,
1047 boost::shared_ptr<MoFEM::FEMethod>,
1048 boost::shared_ptr<MoFEM::BasicMethod>,
1049 boost::shared_ptr<MoFEM::BasicMethod>>(
1050 dm, ts, fe_name, method, pre_only, post_only);
1052}
1053
1054PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx) {
1055 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1057 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1058 *ksp_ctx = dm_field->kspCtx.get();
1060}
1061
1062PetscErrorCode
1063DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
1064 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1066 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1067 const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
1069}
1070
1071PetscErrorCode DMMoFEMSetKspCtx(DM dm,
1072 boost::shared_ptr<MoFEM::KspCtx> ksp_ctx) {
1073 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1075 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1076 dm_field->kspCtx = ksp_ctx;
1078}
1079
1080PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx) {
1081 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1083 DMCtx *dm_field = (DMCtx *)dm->data;
1084 *snes_ctx = dm_field->snesCtx.get();
1086}
1087
1088PetscErrorCode
1089DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
1090 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1092 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1093 const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
1095}
1096
1097PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
1098 boost::shared_ptr<MoFEM::SnesCtx> snes_ctx) {
1099 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1101 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1102 dm_field->snesCtx = snes_ctx;
1104}
1105
1106/** get if read mesh is partitioned
1107 * \ingroup dm
1108 */
1109PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned) {
1110 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1112 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1113 dm_field->isPartitioned = is_partitioned;
1115}
1116
1117/** get if read mesh is partitioned
1118 * \ingroup dm
1119 */
1120PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned) {
1121 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1123 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1124 *is_partitioned = dm_field->isPartitioned;
1126}
1127
1128PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx) {
1129 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1131 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1132 *ts_ctx = dm_field->tsCtx.get();
1134}
1135
1136PetscErrorCode DMMoFEMGetTsCtx(DM dm,
1137 const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
1138 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1140 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1141 const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1143}
1144
1145PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> ts_ctx) {
1146 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1148 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1149 dm_field->tsCtx = ts_ctx;
1151}
1152
1153PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g) {
1154 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1156 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1157 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1158 dm_field->problemName, COL, g);
1159 CHKERR VecSetDM(*g, dm);
1161}
1162
1163PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj<Vec> &g_ptr) {
1164 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1166 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1167 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1168 dm_field->problemName, COL, g_ptr);
1169 CHKERR VecSetDM(g_ptr, dm);
1171}
1172
1173PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l) {
1174 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1176 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1177 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1178 dm_field->problemName, COL, l);
1179 CHKERR VecSetDM(*l, dm);
1181}
1182
1183PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M) {
1184 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1186 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1187 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1189 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1190 M);
1191 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1193 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1194 M);
1195 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1197 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1198 dm_field->problemName, M);
1199 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1201 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1202 dm_field->problemName, M);
1203 } else {
1204 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1205 "Matrix type not implemented");
1206 }
1207 CHKERR MatSetDM(*M, dm);
1209}
1210
1212 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1214 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1215 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1217 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1218 M);
1219 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1221 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1222 M);
1223 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1225 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1226 dm_field->problemName, M);
1227 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1229 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1230 dm_field->problemName, M);
1231 } else {
1232 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1233 "Matrix type not implemented");
1234 }
1235 CHKERR MatSetDM(M, dm);
1237}
1238
1239#if PETSC_VERSION_GE(3, 7, 0)
1240PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
1241 DM dm) {
1242#elif PETSC_VERSION_GE(3, 5, 3)
1243PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm) {
1244#else
1245PetscErrorCode DMSetFromOptions_MoFEM(DM dm) {
1246#endif
1247
1248 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1250 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1251#if PETSC_VERSION_GE(3, 5, 3)
1252 ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1253 CHKERRG(ierr);
1254#else
1255 ierr = PetscOptionsHead("DMMoFEM Options");
1256 CHKERRG(ierr);
1257#endif
1258 ierr = PetscOptionsBool("-dm_is_partitioned",
1259 "set if mesh is partitioned (works which native MOAB "
1260 "file format, i.e. h5m",
1261 "DMSetUp", dm_field->isPartitioned,
1262 &dm_field->isPartitioned, NULL);
1263 CHKERRG(ierr);
1265}
1266
1267PetscErrorCode DMSetUp_MoFEM(DM dm) {
1268 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1269 ProblemsManager *prb_mng_ptr;
1271 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1272 CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1273
1274 if (dm_field->isCompDM) {
1275 // It is composite probelm
1276 CHKERR prb_mng_ptr->buildComposedProblem(
1277 dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1278 dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1279 } else {
1280 if (dm_field->isPartitioned) {
1282 dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1283 dm_field->verbosity);
1284 } else {
1285 CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1286 dm_field->isSquareMatrix == PETSC_TRUE,
1287 dm_field->verbosity);
1288 CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1289 dm_field->verbosity);
1290 }
1291 }
1292
1293 // Partition finite elements
1294 if (dm_field->isPartitioned) {
1295 CHKERR prb_mng_ptr->partitionFiniteElements(
1296 dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1298 dm_field->problemName, dm_field->verbosity);
1299 } else {
1300 // partition finite elemnets
1301 CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1302 -1, -1, dm_field->verbosity);
1303 // Get ghost DOFs
1304 CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1305 dm_field->verbosity);
1306 }
1307
1308 // Set flag that problem is build and partitioned
1309 dm_field->isProblemBuild = PETSC_TRUE;
1310
1312}
1313
1314PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm) {
1315 PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1316 ProblemsManager *prb_mng_ptr;
1318
1319 DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1320
1321 // build sub dm problem
1322 CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1323
1324 map<std::string, boost::shared_ptr<Range>> *entity_map_row = nullptr;
1325 map<std::string, boost::shared_ptr<Range>> *entity_map_col = nullptr;
1326
1327 if (subdm_field->mapTypeRow.size())
1328 entity_map_row = &subdm_field->mapTypeRow;
1329 if (subdm_field->mapTypeCol.size())
1330 entity_map_col = &subdm_field->mapTypeCol;
1331
1332 CHKERR prb_mng_ptr->buildSubProblem(
1333 subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1334 subdm_field->problemMainOfSubPtr->getName(),
1335 subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1336 subdm_field->verbosity);
1337
1338 // partition problem
1339 subdm_field->isPartitioned = subdm_field->isPartitioned;
1340 if (subdm_field->isPartitioned) {
1341 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1342 0, subdm_field->sIze,
1343 subdm_field->verbosity);
1344 // set ghost nodes
1346 subdm_field->problemName, subdm_field->verbosity);
1347 } else {
1348 // partition finite elements
1349 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1350 -1, -1, subdm_field->verbosity);
1351 // set ghost nodes
1352 CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1353 subdm_field->verbosity);
1354 }
1355
1356 subdm_field->isProblemBuild = PETSC_TRUE;
1357
1359}
1360
1361PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec g, InsertMode mode,
1362 Vec l) {
1363 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1365 CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1367}
1368
1369PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec g, InsertMode mode, Vec l) {
1371 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1373
1374 CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1375
1376 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1377 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1378 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1379
1380 double *array_loc, *array_glob;
1381 CHKERR VecGetArray(l, &array_loc);
1382 CHKERR VecGetArray(g, &array_glob);
1383 switch (mode) {
1384 case INSERT_VALUES:
1385 cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1386 break;
1387 case ADD_VALUES:
1388 cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1389 break;
1390 default:
1391 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1392 }
1393 CHKERR VecRestoreArray(l, &array_loc);
1394 CHKERR VecRestoreArray(g, &array_glob);
1396}
1397
1398PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM dm, Vec l, InsertMode mode,
1399 Vec g) {
1400
1401 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1403
1404 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1405 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1406 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1407
1408 double *array_loc, *array_glob;
1409 CHKERR VecGetArray(l, &array_loc);
1410 CHKERR VecGetArray(g, &array_glob);
1411 switch (mode) {
1412 case INSERT_VALUES:
1413 cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1414 break;
1415 case ADD_VALUES:
1416 cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1417 break;
1418 default:
1419 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1420 }
1421 CHKERR VecRestoreArray(l, &array_loc);
1422 CHKERR VecRestoreArray(g, &array_glob);
1423
1425}
1426
1427PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec l, InsertMode mode, Vec g) {
1428 //
1431}
1432
1433PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
1434 char ***fieldNames, IS **fields) {
1435 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1437
1438 if (numFields) {
1439 *numFields = 0;
1440 }
1441 if (fieldNames) {
1442 *fieldNames = NULL;
1443 }
1444 if (fields) {
1445 *fields = NULL;
1446 }
1447
1448 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1449 auto fields_ptr = dm_field->mField_ptr->get_fields();
1450 Field_multiIndex::iterator fit, hi_fit;
1451 fit = fields_ptr->begin();
1452 hi_fit = fields_ptr->end();
1453 *numFields = std::distance(fit, hi_fit);
1454
1455 if (fieldNames) {
1456 CHKERR PetscMalloc1(*numFields, fieldNames);
1457 }
1458 if (fields) {
1459 CHKERR PetscMalloc1(*numFields, fields);
1460 }
1461
1462 for (int f = 0; fit != hi_fit; fit++, f++) {
1463 if (fieldNames) {
1464 CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1465 (char **)&(*fieldNames)[f]);
1466 }
1467 if (fields) {
1469 ->isCreateProblemFieldAndRank(
1470 dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1471 fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1472 }
1473 }
1474
1476}
1477
1478PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
1479 IS *is) {
1480 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1482 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1484 ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1485 field_name, 0, 1000, is);
1487}
1488
1489PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb) {
1490 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1492 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1493 dm_field->verbosity = verb;
1495}
1496
1497} // namespace MoFEM
static Index< 'M', 3 > M
Discrete manager interface for MoFEM.
if(!static_cast< bool >(ifstream(param_file)))
@ 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: DMMoFEM.cpp:1361
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1369
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMoFEM.cpp:1109
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
Definition: DMMoFEM.cpp:219
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMoFEM.cpp:1120
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:483
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMoFEM.cpp:704
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, std::string fe_name, PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMoFEM.cpp:459
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMoFEM.cpp:1173
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition: DMMoFEM.cpp:442
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMoFEM.cpp:786
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMoFEM.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: DMMoFEM.cpp:118
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMoFEM.cpp:542
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMoFEM.cpp:509
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on col.
Definition: DMMoFEM.cpp:371
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1398
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMoFEM.cpp:412
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMoFEM.cpp:75
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > snes_ctx)
Set MoFEM::SnesCtx data structure.
Definition: DMMoFEM.cpp:1097
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > ksp_ctx)
set MoFEM::KspCtx data structure
Definition: DMMoFEM.cpp:1071
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
Definition: DMMoFEM.cpp:242
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMoFEM.cpp:308
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on row.
Definition: DMMoFEM.cpp:353
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMoFEM.cpp:745
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMoFEM.cpp:1183
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMoFEM.cpp:1245
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMoFEM.cpp:1128
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:47
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
Definition: DMMoFEM.cpp:1054
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problem
Definition: DMMoFEM.cpp:473
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMoFEM.cpp:591
PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set KSP right hand side evaluation function
Definition: DMMoFEM.cpp:623
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, std::string fe_name)
Resolve shared entities.
Definition: DMMoFEM.cpp:450
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > ts_ctx)
Set MoFEM::TsCtx data structure.
Definition: DMMoFEM.cpp:1145
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:572
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMoFEM.cpp:1478
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMoFEM.cpp:839
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMoFEM.cpp:400
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
Definition: DMMoFEM.cpp:275
PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMoFEM.cpp:1003
PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS implicit function evaluation function
Definition: DMMoFEM.cpp:961
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMoFEM.cpp:1153
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMoFEM.cpp:1080
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMoFEM.cpp:521
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition: DMMoFEM.cpp:391
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMoFEM.cpp:1314
PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set KSP operators and push mofem finite element methods.
Definition: DMMoFEM.cpp:664
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMoFEM.cpp:1433
PetscErrorCode DMMoFEMUnSetElement(DM dm, std::string fe_name)
unset element from dm
Definition: DMMoFEM.cpp:500
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMoFEM.cpp:1267
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1427
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:553
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMoFEM.cpp:532
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMoFEM.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:389
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:339
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
MoFEM::TsCtx * ts_ctx
Definition: level_set.cpp:1884
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: DMMoFEM.cpp:1042
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:131
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:221
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: DMMoFEM.cpp:868
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: DMMoFEM.cpp:424
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMoFEM.cpp:433
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:566
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:281
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: DMMoFEM.cpp:317
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: DMMoFEM.cpp:335
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:469
PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx)
Run over elements in the list.
Definition: KspCtx.cpp:86
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:380
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: DMMoFEM.cpp:908
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition: DMMoFEM.cpp:1489
PetscErrorCode DMMoFEMSwapDMCtx(DM dm, DM dm_swap)
Swap internal data struture.
Definition: DMMoFEM.cpp:200
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
Definition: DMMoFEM.cpp:184
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:929
PetscBool isCompDM
Definition: DMMoFEM.hpp:953
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:936
const Problem * problemPtr
pinter to problem data structure
Definition: DMMoFEM.hpp:945
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:954
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition: DMMoFEM.hpp:968
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition: DMMoFEM.hpp:969
PetscBool isSquareMatrix
true if rows equals to cols
Definition: DMMoFEM.hpp:940
PetscBool isPartitioned
true if read mesh is on parts
Definition: DMMoFEM.hpp:939
PetscBool isSubDM
Definition: DMMoFEM.hpp:948
int referenceNumber
Definition: DMMoFEM.hpp:965
int verbosity
verbosity
Definition: DMMoFEM.hpp:964
PetscBool isProblemBuild
True if problem is build.
Definition: DMMoFEM.hpp:935
Interface * mField_ptr
MoFEM interface.
Definition: DMMoFEM.hpp:934
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition: DMMoFEM.hpp:967
std::vector< std::string > colCompPrb
Definition: DMMoFEM.hpp:955
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: DMMoFEM.cpp:41
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:950
PetscBool destroyProblem
If true destroy problem with DM.
Definition: DMMoFEM.hpp:959
std::vector< std::string > rowFields
Definition: DMMoFEM.hpp:949
std::map< std::string, boost::shared_ptr< Range > > mapTypeRow
Definition: DMMoFEM.hpp:956
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
Definition: DMMoFEM.hpp:951
std::map< std::string, boost::shared_ptr< Range > > mapTypeCol
Definition: DMMoFEM.hpp:957
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:13
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:298
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:344
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:348
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:399
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:340
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:13
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