v0.13.2
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 auto *dm_field = static_cast<DMCtx *>(dm->data);
192
193 if (static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber == 0)
194 delete static_cast<DMCtx *>(dm_duplicate->data);
195
196 dm_duplicate->data = dm->data;
197 ++(static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber);
198
200}
201
202PetscErrorCode DMMoFEMSwapDMCtx(DM dm, DM dm_swap) {
204 if (!dm->data)
205 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
206 "data structure for MoFEM not yet created on dm");
207 if (!dm_swap->data)
208 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
209 "data structure for MoFEM not yet created on swap dm");
210
211 auto *dm_field = static_cast<DMCtx *>(dm->data);
212 auto *dm_field_swap = static_cast<DMCtx *>(dm_swap->data);
213
214 auto tmp_field = dm_field;
215 dm_field = dm_field_swap;
216 dm_field_swap = tmp_field;
217
219}
220
221PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]) {
223
224 if (!dm->data) {
225 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
226 "data structure for MoFEM not yet created");
227 }
228 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
229
230 CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
231 dm_field->problemPtr->getBitRefLevel(),
232 dm_field->problemPtr->getBitRefLevelMask());
233
234 DMCtx *subdm_field = (DMCtx *)subdm->data;
235 subdm_field->isSubDM = PETSC_TRUE;
236 subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
237 subdm_field->isPartitioned = dm_field->isPartitioned;
238 subdm_field->isSquareMatrix = PETSC_FALSE;
239 subdm->ops->setup = DMSubDMSetUp_MoFEM;
240
242}
243
244PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[]) {
245 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
247 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
248 if (!dm->data) {
249 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
250 "data structure for MoFEM not yet created");
251 }
252 if (!dm_field->isSubDM) {
253 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
254 }
255 dm_field->rowFields.push_back(field_name);
256 dm_field->mapTypeRow.erase(field_name);
258}
259
260PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
261 boost::shared_ptr<Range> r_ptr) {
262 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
264 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
265 if (!dm->data) {
266 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
267 "data structure for MoFEM not yet created");
268 }
269 if (!dm_field->isSubDM) {
270 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
271 }
272 dm_field->rowFields.push_back(field_name);
273 dm_field->mapTypeRow[field_name] = r_ptr;
275}
276
277PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[]) {
278 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
280 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
281 if (!dm->data) {
282 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
283 "data structure for MoFEM not yet created");
284 }
285 if (!dm_field->isSubDM) {
286 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
287 }
288 dm_field->colFields.push_back(field_name);
289 dm_field->mapTypeCol.erase(field_name);
291}
292
293PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
294 boost::shared_ptr<Range> r_ptr) {
295 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
297 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
298 if (!dm->data) {
299 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
300 "data structure for MoFEM not yet created");
301 }
302 if (!dm_field->isSubDM) {
303 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
304 }
305 dm_field->colFields.push_back(field_name);
306 dm_field->mapTypeCol[field_name] = r_ptr;
308}
309
310PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm) {
312 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
314 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
315 *is_sub_dm = dm_field->isSubDM;
317}
318
319PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is) {
321 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
323 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
324 if (dm_field->isSubDM != PETSC_TRUE) {
325 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
326 "This DM is not created as a SubDM");
327 }
328 if (dm_field->isProblemBuild != PETSC_TRUE) {
329 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
330 }
331 boost::shared_ptr<Problem::SubProblemData> sub_data =
332 dm_field->problemPtr->getSubData();
333 CHKERR sub_data->getRowIs(is);
335}
336
337PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is) {
339 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
341 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
342 if (dm_field->isSubDM != PETSC_TRUE) {
343 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
344 "This DM is not created as a SubDM");
345 }
346 if (dm_field->isProblemBuild != PETSC_TRUE) {
347 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
348 }
349 boost::shared_ptr<Problem::SubProblemData> sub_data =
350 dm_field->problemPtr->getSubData();
351 CHKERR sub_data->getColIs(is);
353}
354
355PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]) {
356 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
358 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
359 if (!dm->data) {
360 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
361 "data structure for MoFEM not yet created");
362 }
363 if (!dm_field->isCompDM) {
364 dm_field->isCompDM = PETSC_TRUE;
365 }
366 dm_field->rowCompPrb.push_back(prb_name);
367 if (dm_field->isSquareMatrix) {
368 dm_field->colCompPrb.push_back(prb_name);
369 }
371}
372
373PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]) {
374 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
376 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
377 if (!dm->data) {
378 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
379 "data structure for MoFEM not yet created");
380 }
381 if (!dm_field->isCompDM) {
382 dm_field->isCompDM = PETSC_TRUE;
383 }
384 if (dm_field->isSquareMatrix) {
385 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
386 "No need to add problem on column when problem block structurally "
387 "symmetric");
388 }
389 dm_field->colCompPrb.push_back(prb_name);
391}
392
393PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm) {
395 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
397 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
398 *is_comp_dm = dm_field->isCompDM;
400}
401
402PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr) {
403 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
405 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
406 if (!dm->data) {
407 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
408 "data structure for MoFEM not yet created");
409 }
410 *m_field_ptr = dm_field->mField_ptr;
412}
413
414PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr) {
415 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
417 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
418 if (!dm->data) {
419 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
420 "data structure for MoFEM not yet created");
421 }
422 *problem_ptr = dm_field->problemPtr;
424}
425
426PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem) {
428 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
430 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
431 dm_field->destroyProblem = destroy_problem;
433}
434
435PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem) {
437 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
439 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
440 *destroy_problem = dm_field->destroyProblem;
442}
443
444PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem) {
445 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
447 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
448 dm_field->isSquareMatrix = square_problem;
450}
451
452PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, std::string fe_name) {
453 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
455 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
457 ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
459}
460
461PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, std::string fe_name,
462 PetscLayout *layout) {
464 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
466 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
467
468 MPI_Comm comm;
469 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
470 CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
471 layout);
473}
474
475PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem) {
478 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
480 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
481 *square_problem = dm_field->isSquareMatrix;
483}
484
485PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name) {
486 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
488 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
490 dm_field->problemName, fe_name);
492}
493
494PetscErrorCode DMMoFEMAddElement(DM dm, std::vector<std::string> fe_name) {
496 for (auto fe : fe_name) {
498 }
500}
501
502PetscErrorCode DMMoFEMUnSetElement(DM dm, std::string fe_name) {
503 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
505 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
507 dm_field->problemName, fe_name);
509}
510
511PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
512 ScatterMode scatter_mode) {
514 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
516 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
517 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
518 dm_field->problemPtr, COL, l, mode, scatter_mode);
519 CHKERRG(ierr);
521}
522
523PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
524 ScatterMode scatter_mode) {
525 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
527 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
528 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
529 dm_field->problemPtr, COL, g, mode, scatter_mode);
530 CHKERRG(ierr);
532}
533
534PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
535 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
537 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
539 dm_field->problemPtr, *method);
540 CHKERRG(ierr);
542}
543
544PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
545 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
547 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
549 dm_field->problemPtr, *method);
550 CHKERRG(ierr);
552}
553
554PetscErrorCode
555DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
556 MoFEM::FEMethod *method, int low_rank,
557 int up_rank, CacheTupleWeakPtr cache_ptr) {
559 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
561 dm_field->problemPtr, fe_name, *method, low_rank, up_rank, nullptr,
562 MF_EXIST, cache_ptr);
563 CHKERRG(ierr);
565}
566
568 DM dm, const std::string fe_name, boost::shared_ptr<MoFEM::FEMethod> method,
569 int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr) {
570 return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
571 low_rank, up_rank, cache_ptr);
572}
573
574PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
575 MoFEM::FEMethod *method,
576 CacheTupleWeakPtr cache_ptr) {
577 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
579 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
581 dm, fe_name, method, dm_field->rAnk, dm_field->rAnk, cache_ptr);
582 CHKERRG(ierr);
584}
585
586PetscErrorCode
587DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
588 boost::shared_ptr<MoFEM::FEMethod> method,
589 CacheTupleWeakPtr cache_ptr) {
590 return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get(), cache_ptr);
591}
592
593PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
594 MoFEM::DofMethod *method) {
595 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
597 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
598 ierr =
599 dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
600 *method, dm_field->rAnk, dm_field->rAnk);
601 CHKERRG(ierr);
603}
604
605template <class S, class T0, class T1, class T2>
606static PetscErrorCode DMMoFEMKSPSetComputeRHS(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_Rhs().push_back(pre_only);
613 }
614 if (method) {
615 dm_field->kspCtx->get_loops_to_do_Rhs().push_back(
616 PairNameFEMethodPtr(fe_name, method));
617 }
618 if (post_only) {
619 dm_field->kspCtx->get_postProcess_to_do_Rhs().push_back(post_only);
620 }
621 CHKERR DMKSPSetComputeRHS(dm, KspRhs, dm_field->kspCtx.get());
623}
624
625PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
626 MoFEM::FEMethod *method,
627 MoFEM::BasicMethod *pre_only,
628 MoFEM::BasicMethod *post_only) {
629 return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
631 dm, fe_name, method, pre_only, post_only);
632}
633
634PetscErrorCode
635DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
636 boost::shared_ptr<MoFEM::FEMethod> method,
637 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
638 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
639 return DMMoFEMKSPSetComputeRHS<const std::string,
640 boost::shared_ptr<MoFEM::FEMethod>,
641 boost::shared_ptr<MoFEM::BasicMethod>,
642 boost::shared_ptr<MoFEM::BasicMethod>>(
643 dm, fe_name, method, pre_only, post_only);
644}
645
646template <class S, class T0, class T1, class T2>
647static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method,
648 T1 pre_only, T2 post_only) {
649 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
651 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
652 if (pre_only) {
653 dm_field->kspCtx->get_preProcess_to_do_Mat().push_back(pre_only);
654 }
655 if (method) {
656 dm_field->kspCtx->get_loops_to_do_Mat().push_back(
657 PairNameFEMethodPtr(fe_name, method));
658 }
659 if (post_only) {
660 dm_field->kspCtx->get_postProcess_to_do_Mat().push_back(post_only);
661 }
662 CHKERR DMKSPSetComputeOperators(dm, KspMat, dm_field->kspCtx.get());
664}
665
666PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
667 MoFEM::FEMethod *method,
668 MoFEM::BasicMethod *pre_only,
669 MoFEM::BasicMethod *post_only) {
670 return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
673 dm, fe_name, method, pre_only, post_only);
674}
675
676PetscErrorCode
677DMMoFEMKSPSetComputeOperators(DM dm, const std::string fe_name,
678 boost::shared_ptr<MoFEM::FEMethod> method,
679 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
680 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
681 return DMMoFEMKSPSetComputeOperators<const std::string,
682 boost::shared_ptr<MoFEM::FEMethod>>(
683 dm, fe_name, method, pre_only, post_only);
684}
685
686template <class S, class T0, class T1, class T2>
687static PetscErrorCode DMMoFEMSNESSetFunction(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_Rhs().push_back(pre_only);
694 }
695 if (method) {
696 dm_field->snesCtx->get_loops_to_do_Rhs().push_back(
697 PairNameFEMethodPtr(fe_name, method));
698 }
699 if (post_only) {
700 dm_field->snesCtx->get_postProcess_to_do_Rhs().push_back(post_only);
701 }
702 CHKERR DMSNESSetFunction(dm, SnesRhs, dm_field->snesCtx.get());
704}
705
706PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
707 MoFEM::FEMethod *method,
708 MoFEM::BasicMethod *pre_only,
709 MoFEM::BasicMethod *post_only) {
710 return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
712 dm, fe_name, method, pre_only, post_only);
713}
714
715PetscErrorCode
716DMMoFEMSNESSetFunction(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 DMMoFEMSNESSetFunction<const std::string,
721 boost::shared_ptr<MoFEM::FEMethod>,
722 boost::shared_ptr<MoFEM::BasicMethod>,
723 boost::shared_ptr<MoFEM::BasicMethod>>(
724 dm, fe_name, method, pre_only, post_only);
725}
726
727template <class S, class T0, class T1, class T2>
728static PetscErrorCode DMMoFEMSNESSetJacobian(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->snesCtx->get_preProcess_to_do_Mat().push_back(pre_only);
735 }
736 if (method) {
737 dm_field->snesCtx->get_loops_to_do_Mat().push_back(
738 PairNameFEMethodPtr(fe_name, method));
739 }
740 if (post_only) {
741 dm_field->snesCtx->get_postProcess_to_do_Mat().push_back(post_only);
742 }
743 CHKERR DMSNESSetJacobian(dm, SnesMat, dm_field->snesCtx.get());
745}
746
747PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
748 MoFEM::FEMethod *method,
749 MoFEM::BasicMethod *pre_only,
750 MoFEM::BasicMethod *post_only) {
751 return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
753 dm, fe_name, method, pre_only, post_only);
754}
755
756PetscErrorCode
757DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
758 boost::shared_ptr<MoFEM::FEMethod> method,
759 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
760 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
761 return DMMoFEMSNESSetJacobian<const std::string,
762 boost::shared_ptr<MoFEM::FEMethod>,
763 boost::shared_ptr<MoFEM::BasicMethod>,
764 boost::shared_ptr<MoFEM::BasicMethod>>(
765 dm, fe_name, method, pre_only, post_only);
766}
767
768template <class S, class T0, class T1, class T2>
769static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method,
770 T1 pre_only, T2 post_only) {
771 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
773 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
774 if (pre_only) {
775 dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
776 }
777 if (method) {
778 dm_field->tsCtx->getLoopsIFunction().push_back(
779 PairNameFEMethodPtr(fe_name, method));
780 }
781 if (post_only) {
782 dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
783 }
784 CHKERR DMTSSetIFunction(dm, TsSetIFunction, dm_field->tsCtx.get());
786}
787
788PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
789 MoFEM::FEMethod *method,
790 MoFEM::BasicMethod *pre_only,
791 MoFEM::BasicMethod *post_only) {
792 return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
794 dm, fe_name, method, pre_only, post_only);
796}
797
798PetscErrorCode
799DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
800 boost::shared_ptr<MoFEM::FEMethod> method,
801 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
802 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
803 return DMMoFEMTSSetIFunction<const std::string,
804 boost::shared_ptr<MoFEM::FEMethod>,
805 boost::shared_ptr<MoFEM::BasicMethod>,
806 boost::shared_ptr<MoFEM::BasicMethod>>(
807 dm, fe_name, method, pre_only, post_only);
809}
810
811template <class S, class T0, class T1, class T2>
812static PetscErrorCode DMMoFEMTSSetIJacobian(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->getPreProcessIJacobian().push_back(pre_only);
819 }
820 if (method) {
821 dm_field->tsCtx->getLoopsIJacobian().push_back(
822 PairNameFEMethodPtr(fe_name, method));
823 }
824 if (post_only) {
825 dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
826 }
827 CHKERR DMTSSetIJacobian(dm, TsSetIJacobian, dm_field->tsCtx.get());
829}
830
831PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
832 MoFEM::FEMethod *method,
833 MoFEM::BasicMethod *pre_only,
834 MoFEM::BasicMethod *post_only) {
835 return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
836 MoFEM::BasicMethod *>(dm, fe_name, method,
837 pre_only, post_only);
838}
839
840PetscErrorCode
841DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
842 boost::shared_ptr<MoFEM::FEMethod> method,
843 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
844 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
845 return DMMoFEMTSSetIJacobian<const std::string,
846 boost::shared_ptr<MoFEM::FEMethod>,
847 boost::shared_ptr<MoFEM::BasicMethod>,
848 boost::shared_ptr<MoFEM::BasicMethod>>(
849 dm, fe_name, method, pre_only, post_only);
850}
851
852template <class S, class T0, class T1, class T2>
853static PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, S fe_name, T0 method,
854 T1 pre_only, T2 post_only) {
855 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
857 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
858 if (pre_only)
859 dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
860 if (method)
861 dm_field->tsCtx->getLoopsRHSFunction().push_back(
862 PairNameFEMethodPtr(fe_name, method));
863 if (post_only)
864 dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
865 CHKERR DMTSSetRHSFunction(dm, TsSetRHSFunction, dm_field->tsCtx.get());
867}
868
869PetscErrorCode
870DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name,
871 boost::shared_ptr<MoFEM::FEMethod> method,
872 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
873 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
874 return DMMoFEMTSSetRHSFunction<const std::string,
875 boost::shared_ptr<MoFEM::FEMethod>,
876 boost::shared_ptr<MoFEM::BasicMethod>,
877 boost::shared_ptr<MoFEM::BasicMethod>>(
878 dm, fe_name, method, pre_only, post_only);
880}
881
882PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const char fe_name[],
883 MoFEM::FEMethod *method,
884 MoFEM::BasicMethod *pre_only,
885 MoFEM::BasicMethod *post_only) {
886 return DMMoFEMTSSetRHSFunction<const char *, MoFEM::FEMethod *,
888 dm, fe_name, method, pre_only, post_only);
890}
891
892template <class S, class T0, class T1, class T2>
893static PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, S fe_name, T0 method,
894 T1 pre_only, T2 post_only) {
895 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
897 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
898 if (pre_only)
899 dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
900 if (method)
901 dm_field->tsCtx->getLoopsRHSFunction().push_back(
902 PairNameFEMethodPtr(fe_name, method));
903 if (post_only)
904 dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
905 CHKERR DMTSSetRHSJacobian(dm, TsSetRHSJacobian, dm_field->tsCtx.get());
907}
908
909PetscErrorCode
910DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
911 boost::shared_ptr<MoFEM::FEMethod> method,
912 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
913 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
914 return DMMoFEMTSSetRHSJacobian<const std::string,
915 boost::shared_ptr<MoFEM::FEMethod>,
916 boost::shared_ptr<MoFEM::BasicMethod>,
917 boost::shared_ptr<MoFEM::BasicMethod>>(
918 dm, fe_name, method, pre_only, post_only);
920}
921
922PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const char fe_name[],
923 MoFEM::FEMethod *method,
924 MoFEM::BasicMethod *pre_only,
925 MoFEM::BasicMethod *post_only) {
926 return DMMoFEMTSSetRHSJacobian<const char *, MoFEM::FEMethod *,
928 dm, fe_name, method, pre_only, post_only);
930}
931
932template <class S, class T0, class T1, class T2>
933static PetscErrorCode DMMoFEMTSSetI2Function(DM dm, S fe_name, T0 method,
934 T1 pre_only, T2 post_only) {
935 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
937 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
938 if (pre_only) {
939 dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
940 }
941 if (method) {
942 dm_field->tsCtx->getLoopsIFunction().push_back(
943 PairNameFEMethodPtr(fe_name, method));
944 }
945 if (post_only) {
946 dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
947 }
948 CHKERR DMTSSetI2Function(dm, TsSetI2Function, dm_field->tsCtx.get());
950}
951
952PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const char fe_name[],
953 MoFEM::FEMethod *method,
954 MoFEM::BasicMethod *pre_only,
955 MoFEM::BasicMethod *post_only) {
956 return DMMoFEMTSSetI2Function<const char *, MoFEM::FEMethod *,
958 dm, fe_name, method, pre_only, post_only);
960}
961
962PetscErrorCode
963DMMoFEMTSSetI2Function(DM dm, const std::string fe_name,
964 boost::shared_ptr<MoFEM::FEMethod> method,
965 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
966 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
967 return DMMoFEMTSSetI2Function<const std::string,
968 boost::shared_ptr<MoFEM::FEMethod>,
969 boost::shared_ptr<MoFEM::BasicMethod>,
970 boost::shared_ptr<MoFEM::BasicMethod>>(
971 dm, fe_name, method, pre_only, post_only);
973}
974
975template <class S, class T0, class T1, class T2>
976static PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, 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->getPreProcessIJacobian().push_back(pre_only);
983 }
984 if (method) {
985 dm_field->tsCtx->getLoopsIJacobian().push_back(
986 PairNameFEMethodPtr(fe_name, method));
987 }
988 if (post_only) {
989 dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
990 }
991 CHKERR DMTSSetI2Jacobian(dm, TsSetI2Jacobian, dm_field->tsCtx.get());
993}
994
995PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const char fe_name[],
996 MoFEM::FEMethod *method,
997 MoFEM::BasicMethod *pre_only,
998 MoFEM::BasicMethod *post_only) {
999 return DMMoFEMTSSetI2Jacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
1000 MoFEM::BasicMethod *>(dm, fe_name, method,
1001 pre_only, post_only);
1002}
1003
1004PetscErrorCode
1005DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name,
1006 boost::shared_ptr<MoFEM::FEMethod> method,
1007 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
1008 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
1009 return DMMoFEMTSSetI2Jacobian<const std::string,
1010 boost::shared_ptr<MoFEM::FEMethod>,
1011 boost::shared_ptr<MoFEM::BasicMethod>,
1012 boost::shared_ptr<MoFEM::BasicMethod>>(
1013 dm, fe_name, method, pre_only, post_only);
1014}
1015
1016template <class S, class T0, class T1, class T2>
1017static PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, S fe_name, T0 method,
1018 T1 pre_only, T2 post_only) {
1019 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1021 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1022 if (pre_only)
1023 dm_field->tsCtx->getPreProcessMonitor().push_back(pre_only);
1024 if (method)
1025 dm_field->tsCtx->getLoopsMonitor().push_back(
1026 PairNameFEMethodPtr(fe_name, method));
1027 if (post_only)
1028 dm_field->tsCtx->getPostProcessMonitor().push_back(post_only);
1029 CHKERR TSMonitorSet(ts, TsMonitorSet, dm_field->tsCtx.get(), PETSC_NULL);
1031}
1032
1033PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const char fe_name[],
1034 MoFEM::FEMethod *method,
1035 MoFEM::BasicMethod *pre_only,
1036 MoFEM::BasicMethod *post_only) {
1037 return DMMoFEMTSSetMonitor<const char *, MoFEM::FEMethod *,
1039 dm, ts, fe_name, method, pre_only, post_only);
1041}
1042
1043PetscErrorCode
1044DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name,
1045 boost::shared_ptr<MoFEM::FEMethod> method,
1046 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
1047 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
1048 return DMMoFEMTSSetMonitor<const std::string,
1049 boost::shared_ptr<MoFEM::FEMethod>,
1050 boost::shared_ptr<MoFEM::BasicMethod>,
1051 boost::shared_ptr<MoFEM::BasicMethod>>(
1052 dm, ts, fe_name, method, pre_only, post_only);
1054}
1055
1056PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx) {
1057 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1059 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1060 *ksp_ctx = dm_field->kspCtx.get();
1062}
1063
1064PetscErrorCode
1065DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
1066 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1068 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1069 const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
1071}
1072
1073PetscErrorCode DMMoFEMSetKspCtx(DM dm,
1074 boost::shared_ptr<MoFEM::KspCtx> ksp_ctx) {
1075 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1077 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1078 dm_field->kspCtx = ksp_ctx;
1080}
1081
1082PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx) {
1083 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1085 DMCtx *dm_field = (DMCtx *)dm->data;
1086 *snes_ctx = dm_field->snesCtx.get();
1088}
1089
1090PetscErrorCode
1091DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
1092 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1094 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1095 const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
1097}
1098
1099PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
1100 boost::shared_ptr<MoFEM::SnesCtx> snes_ctx) {
1101 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1103 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1104 dm_field->snesCtx = snes_ctx;
1106}
1107
1108/** get if read mesh is partitioned
1109 * \ingroup dm
1110 */
1111PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned) {
1112 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1114 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1115 dm_field->isPartitioned = is_partitioned;
1117}
1118
1119/** get if read mesh is partitioned
1120 * \ingroup dm
1121 */
1122PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned) {
1123 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1125 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1126 *is_partitioned = dm_field->isPartitioned;
1128}
1129
1130PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx) {
1131 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1133 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1134 *ts_ctx = dm_field->tsCtx.get();
1136}
1137
1138PetscErrorCode DMMoFEMGetTsCtx(DM dm,
1139 const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
1140 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1142 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1143 const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1145}
1146
1147PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> ts_ctx) {
1148 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1150 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1151 dm_field->tsCtx = ts_ctx;
1153}
1154
1155PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g) {
1156 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1158 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1159 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1160 dm_field->problemName, COL, g);
1161 CHKERR VecSetDM(*g, dm);
1163}
1164
1165PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj<Vec> &g_ptr) {
1166 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1168 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1169 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1170 dm_field->problemName, COL, g_ptr);
1171 CHKERR VecSetDM(g_ptr, dm);
1173}
1174
1175PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l) {
1176 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1178 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1179 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1180 dm_field->problemName, COL, l);
1181 CHKERR VecSetDM(*l, dm);
1183}
1184
1185PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M) {
1186 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1188 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1189 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1191 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1192 M);
1193 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1195 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1196 M);
1197 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1199 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1200 dm_field->problemName, M);
1201 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1203 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1204 dm_field->problemName, M);
1205 } else {
1206 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1207 "Matrix type not implemented");
1208 }
1209 CHKERR MatSetDM(*M, dm);
1211}
1212
1214 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1216 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1217 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1219 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1220 M);
1221 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1223 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1224 M);
1225 } else if (strcmp(dm->mattype, MATAIJCUSPARSE) == 0) {
1227 ->createMPIAIJCUSPARSEWithArrays<PetscGlobalIdx_mi_tag>(
1228 dm_field->problemName, M);
1229 } else if (strcmp(dm->mattype, MATSEQAIJCUSPARSE) == 0) {
1231 ->createSeqAIJCUSPARSEWithArrays<PetscLocalIdx_mi_tag>(
1232 dm_field->problemName, M);
1233 } else {
1234 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1235 "Matrix type not implemented");
1236 }
1237 CHKERR MatSetDM(M, dm);
1239}
1240
1241#if PETSC_VERSION_GE(3, 7, 0)
1242PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
1243 DM dm) {
1244#elif PETSC_VERSION_GE(3, 5, 3)
1245PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm) {
1246#else
1247PetscErrorCode DMSetFromOptions_MoFEM(DM dm) {
1248#endif
1249
1250 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1252 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1253#if PETSC_VERSION_GE(3, 5, 3)
1254 ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1255 CHKERRG(ierr);
1256#else
1257 ierr = PetscOptionsHead("DMMoFEM Options");
1258 CHKERRG(ierr);
1259#endif
1260 ierr = PetscOptionsBool("-dm_is_partitioned",
1261 "set if mesh is partitioned (works which native MOAB "
1262 "file format, i.e. h5m",
1263 "DMSetUp", dm_field->isPartitioned,
1264 &dm_field->isPartitioned, NULL);
1265 CHKERRG(ierr);
1267}
1268
1269PetscErrorCode DMSetUp_MoFEM(DM dm) {
1270 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1271 ProblemsManager *prb_mng_ptr;
1273 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1274 CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1275
1276 if (dm_field->isCompDM) {
1277 // It is composite probelm
1278 CHKERR prb_mng_ptr->buildComposedProblem(
1279 dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1280 dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1281 } else {
1282 if (dm_field->isPartitioned) {
1284 dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1285 dm_field->verbosity);
1286 } else {
1287 CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1288 dm_field->isSquareMatrix == PETSC_TRUE,
1289 dm_field->verbosity);
1290 CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1291 dm_field->verbosity);
1292 }
1293 }
1294
1295 // Partition finite elements
1296 if (dm_field->isPartitioned) {
1297 CHKERR prb_mng_ptr->partitionFiniteElements(
1298 dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1300 dm_field->problemName, dm_field->verbosity);
1301 } else {
1302 // partition finite elemnets
1303 CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1304 -1, -1, dm_field->verbosity);
1305 // Get ghost DOFs
1306 CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1307 dm_field->verbosity);
1308 }
1309
1310 // Set flag that problem is build and partitioned
1311 dm_field->isProblemBuild = PETSC_TRUE;
1312
1314}
1315
1316PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm) {
1317 PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1318 ProblemsManager *prb_mng_ptr;
1320
1321 DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1322
1323 // build sub dm problem
1324 CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1325
1326 map<std::string, boost::shared_ptr<Range>> *entity_map_row = nullptr;
1327 map<std::string, boost::shared_ptr<Range>> *entity_map_col = nullptr;
1328
1329 if (subdm_field->mapTypeRow.size())
1330 entity_map_row = &subdm_field->mapTypeRow;
1331 if (subdm_field->mapTypeCol.size())
1332 entity_map_col = &subdm_field->mapTypeCol;
1333
1334 CHKERR prb_mng_ptr->buildSubProblem(
1335 subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1336 subdm_field->problemMainOfSubPtr->getName(),
1337 subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1338 subdm_field->verbosity);
1339
1340 // partition problem
1341 subdm_field->isPartitioned = subdm_field->isPartitioned;
1342 if (subdm_field->isPartitioned) {
1343 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1344 0, subdm_field->sIze,
1345 subdm_field->verbosity);
1346 // set ghost nodes
1348 subdm_field->problemName, subdm_field->verbosity);
1349 } else {
1350 // partition finite elements
1351 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1352 -1, -1, subdm_field->verbosity);
1353 // set ghost nodes
1354 CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1355 subdm_field->verbosity);
1356 }
1357
1358 subdm_field->isProblemBuild = PETSC_TRUE;
1359
1361}
1362
1363PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec g, InsertMode mode,
1364 Vec l) {
1365 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1367 CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1369}
1370
1371PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec g, InsertMode mode, Vec l) {
1373 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1375
1376 CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1377
1378 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1379 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1380 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1381
1382 double *array_loc, *array_glob;
1383 CHKERR VecGetArray(l, &array_loc);
1384 CHKERR VecGetArray(g, &array_glob);
1385 switch (mode) {
1386 case INSERT_VALUES:
1387 cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1388 break;
1389 case ADD_VALUES:
1390 cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1391 break;
1392 default:
1393 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1394 }
1395 CHKERR VecRestoreArray(l, &array_loc);
1396 CHKERR VecRestoreArray(g, &array_glob);
1398}
1399
1400PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM dm, Vec l, InsertMode mode,
1401 Vec g) {
1402
1403 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1405
1406 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1407 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1408 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1409
1410 double *array_loc, *array_glob;
1411 CHKERR VecGetArray(l, &array_loc);
1412 CHKERR VecGetArray(g, &array_glob);
1413 switch (mode) {
1414 case INSERT_VALUES:
1415 cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1416 break;
1417 case ADD_VALUES:
1418 cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1419 break;
1420 default:
1421 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1422 }
1423 CHKERR VecRestoreArray(l, &array_loc);
1424 CHKERR VecRestoreArray(g, &array_glob);
1425
1427}
1428
1429PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec l, InsertMode mode, Vec g) {
1430 //
1433}
1434
1435PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
1436 char ***fieldNames, IS **fields) {
1437 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1439
1440 if (numFields) {
1441 *numFields = 0;
1442 }
1443 if (fieldNames) {
1444 *fieldNames = NULL;
1445 }
1446 if (fields) {
1447 *fields = NULL;
1448 }
1449
1450 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1451 auto fields_ptr = dm_field->mField_ptr->get_fields();
1452 Field_multiIndex::iterator fit, hi_fit;
1453 fit = fields_ptr->begin();
1454 hi_fit = fields_ptr->end();
1455 *numFields = std::distance(fit, hi_fit);
1456
1457 if (fieldNames) {
1458 CHKERR PetscMalloc1(*numFields, fieldNames);
1459 }
1460 if (fields) {
1461 CHKERR PetscMalloc1(*numFields, fields);
1462 }
1463
1464 for (int f = 0; fit != hi_fit; fit++, f++) {
1465 if (fieldNames) {
1466 CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1467 (char **)&(*fieldNames)[f]);
1468 }
1469 if (fields) {
1471 ->isCreateProblemFieldAndRank(
1472 dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1473 fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1474 }
1475 }
1476
1478}
1479
1480PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
1481 IS *is) {
1482 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1484 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1486 ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1487 field_name, 0, 1000, is);
1489}
1490
1491PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb) {
1492 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1494 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1495 dm_field->verbosity = verb;
1497}
1498
1499} // 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:1363
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1371
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMoFEM.cpp:1111
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:221
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMoFEM.cpp:1122
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition: DMMoFEM.cpp:485
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMoFEM.cpp:706
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, std::string fe_name, PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMoFEM.cpp:461
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMoFEM.cpp:1175
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition: DMMoFEM.cpp:444
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMoFEM.cpp:788
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:544
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:511
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on col.
Definition: DMMoFEM.cpp:373
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1400
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMoFEM.cpp:414
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:1099
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > ksp_ctx)
set MoFEM::KspCtx data structure
Definition: DMMoFEM.cpp:1073
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
Definition: DMMoFEM.cpp:244
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMoFEM.cpp:310
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on row.
Definition: DMMoFEM.cpp:355
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMoFEM.cpp:747
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMoFEM.cpp:1185
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMoFEM.cpp:1247
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMoFEM.cpp:1130
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:1056
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problem
Definition: DMMoFEM.cpp:475
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMoFEM.cpp:593
PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set KSP right hand side evaluation function
Definition: DMMoFEM.cpp:625
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, std::string fe_name)
Resolve shared entities.
Definition: DMMoFEM.cpp:452
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > ts_ctx)
Set MoFEM::TsCtx data structure.
Definition: DMMoFEM.cpp:1147
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:574
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMoFEM.cpp:1480
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMoFEM.cpp:841
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMoFEM.cpp:402
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
Definition: DMMoFEM.cpp:277
PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMoFEM.cpp:1005
PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS implicit function evaluation function
Definition: DMMoFEM.cpp:963
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMoFEM.cpp:1155
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMoFEM.cpp:1082
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMoFEM.cpp:523
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition: DMMoFEM.cpp:393
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMoFEM.cpp:1316
PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set KSP operators and push mofem finite element methods.
Definition: DMMoFEM.cpp:666
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMoFEM.cpp:1435
PetscErrorCode DMMoFEMUnSetElement(DM dm, std::string fe_name)
unset element from dm
Definition: DMMoFEM.cpp:502
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMoFEM.cpp:1269
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMoFEM.cpp:1429
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMoFEM.cpp:555
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMoFEM.cpp:534
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:391
#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:1876
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:1044
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: DMMoFEM.cpp:870
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:426
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMoFEM.cpp:435
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: DMMoFEM.cpp:319
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:337
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: DMMoFEM.cpp:910
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition: DMMoFEM.cpp:1491
PetscErrorCode DMMoFEMSwapDMCtx(DM dm, DM dm_swap)
Swap internal data struture.
Definition: DMMoFEM.cpp:202
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:924
PetscBool isCompDM
Definition: DMMoFEM.hpp:948
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:931
const Problem * problemPtr
pinter to problem data structure
Definition: DMMoFEM.hpp:940
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:949
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition: DMMoFEM.hpp:963
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition: DMMoFEM.hpp:964
PetscBool isSquareMatrix
true if rows equals to cols
Definition: DMMoFEM.hpp:935
PetscBool isPartitioned
true if read mesh is on parts
Definition: DMMoFEM.hpp:934
PetscBool isSubDM
Definition: DMMoFEM.hpp:943
int referenceNumber
Definition: DMMoFEM.hpp:960
int verbosity
verbosity
Definition: DMMoFEM.hpp:959
PetscBool isProblemBuild
True if problem is build.
Definition: DMMoFEM.hpp:930
Interface * mField_ptr
MoFEM interface.
Definition: DMMoFEM.hpp:929
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition: DMMoFEM.hpp:962
std::vector< std::string > colCompPrb
Definition: DMMoFEM.hpp:950
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: DMMoFEM.cpp:41
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:945
PetscBool destroyProblem
If true destroy problem with DM.
Definition: DMMoFEM.hpp:954
std::vector< std::string > rowFields
Definition: DMMoFEM.hpp:944
std::map< std::string, boost::shared_ptr< Range > > mapTypeRow
Definition: DMMoFEM.hpp:951
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
Definition: DMMoFEM.hpp:946
std::map< std::string, boost::shared_ptr< Range > > mapTypeCol
Definition: DMMoFEM.hpp:952
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:300
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:346
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:350
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:401
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:342
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