v0.13.1
DMMMoFEM.cpp
Go to the documentation of this file.
1/* This file is part of MoFEM.
2 * MoFEM is free software: you can redistribute it and/or modify it under
3 * the terms of the GNU Lesser General Public License as published by the
4 * Free Software Foundation, either version 3 of the License, or (at your
5 * option) any later version.
6 *
7 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
8 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
9 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
10 * License for more details.
11 *
12 * You should have received a copy of the GNU Lesser General Public
13 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
14
15// #undef PETSC_VERSION_RELEASE
16// #define PETSC_VERSION_RELEASE 1
17
18#if PETSC_VERSION_GE(3, 6, 0)
19#include <petsc/private/dmimpl.h> /*I "petscdm.h" I*/
20// #include <petsc/private/vecimpl.h> /*I "petscdm.h" I*/
21#else
22#include <petsc-private/dmimpl.h> /*I "petscdm.h" I*/
23#include <petsc-private/vecimpl.h> /*I "petscdm.h" I*/
24#endif
25
26#include <DMMoFEM.hpp>
27
28namespace MoFEM {
29
31 : mField_ptr(PETSC_NULL), isProblemBuild(PETSC_FALSE),
32 isPartitioned(PETSC_FALSE), isSquareMatrix(PETSC_TRUE),
33 isSubDM(PETSC_FALSE), isCompDM(PETSC_FALSE), destroyProblem(PETSC_FALSE),
34 verbosity(VERBOSE), referenceNumber(0) {
35
36 if (!LogManager::checkIfChannelExist("DMWORLD")) {
37 auto core_log = logging::core::get();
38 core_log->add_sink(
40 core_log->add_sink(
42 core_log->add_sink(
44 LogManager::setLog("DMWORLD");
45 LogManager::setLog("DMSYNC");
46 LogManager::setLog("DMSELF");
47 MOFEM_LOG_TAG("DMWORLD", "DM");
48 MOFEM_LOG_TAG("DMSYNC", "DM");
49 MOFEM_LOG_TAG("DMSELF", "DM");
50 }
51}
52
53MoFEMErrorCode DMCtx::query_interface(boost::typeindex::type_index type_index,
54 UnknownInterface **iface) const {
55 *iface = const_cast<DMCtx *>(this);
56 return 0;
57}
58
59PetscErrorCode DMRegister_MoFEM(const char sname[]) {
61 CHKERR DMRegister(sname, DMCreate_MoFEM);
63}
64
65PetscErrorCode DMSetOperators_MoFEM(DM dm) {
66 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
68
69 dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
70 dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
71 dm->ops->creatematrix = DMCreateMatrix_MoFEM;
72 dm->ops->setup = DMSetUp_MoFEM;
73 dm->ops->destroy = DMDestroy_MoFEM;
74 dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
75 dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
76 dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
77 dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
78 dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
79 dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
80
81 // Default matrix type
82 CHKERR DMSetMatType(dm, MATMPIAIJ);
83
85}
86
87PetscErrorCode DMCreate_MoFEM(DM dm) {
88 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
90 dm->data = new DMCtx();
93}
94
95PetscErrorCode DMDestroy_MoFEM(DM dm) {
96 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
97 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
99
100 MPI_Comm comm;
101 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
102
103 int result;
104 MPI_Comm_compare(comm, PETSC_COMM_SELF, &result);
105 if (result == MPI_IDENT)
106 MOFEM_LOG("DMSELF", Sev::noisy)
107 << "MoFEM DM destroy for problem " << dm_field->problemName
108 << " referenceNumber " << dm_field->referenceNumber;
109 else
110 MOFEM_LOG("DMWORLD", Sev::noisy)
111 << "MoFEM DM destroy for problem " << dm_field->problemName
112 << " referenceNumber " << dm_field->referenceNumber;
113
114 if (dm_field->referenceNumber == 0) {
115 if (dm_field->destroyProblem) {
116
117 if (dm_field->mField_ptr->check_problem(dm_field->problemName)) {
118 dm_field->mField_ptr->delete_problem(dm_field->problemName);
119 } // else problem has to be deleted by the user
120 }
121
122 delete static_cast<DMCtx *>(dm->data);
123
124 } else
125 --dm_field->referenceNumber;
126
128}
129
130PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr,
131 const char problem_name[],
132 const MoFEM::BitRefLevel bit_level,
133 const MoFEM::BitRefLevel bit_mask) {
135
136 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
137 if (!dm->data) {
138 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
139 "data structure for MoFEM not yet created");
140 }
141 if (!m_field_ptr) {
142 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
143 "DM function not implemented into MoFEM");
144 }
145 dm_field->mField_ptr = m_field_ptr;
146 dm_field->problemName = problem_name;
147 if (!m_field_ptr->check_problem(dm_field->problemName)) {
148 // problem is not defined, declare problem internally set bool to
149 // destroyProblem problem with DM
150 dm_field->destroyProblem = PETSC_TRUE;
151 CHKERR dm_field->mField_ptr->add_problem(dm_field->problemName, MF_EXCL,
152 dm_field->verbosity);
153 } else {
154 dm_field->destroyProblem = PETSC_FALSE;
155 }
157 dm_field->problemName, bit_level);
159 dm_field->problemName, bit_mask);
160 dm_field->kspCtx =
161 boost::shared_ptr<KspCtx>(new KspCtx(*m_field_ptr, problem_name));
162 dm_field->snesCtx =
163 boost::shared_ptr<SnesCtx>(new SnesCtx(*m_field_ptr, problem_name));
164 dm_field->tsCtx =
165 boost::shared_ptr<TsCtx>(new TsCtx(*m_field_ptr, problem_name));
166
167 MPI_Comm comm;
168 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
169 int result = 0;
170 MPI_Comm_compare(comm, m_field_ptr->get_comm(), &result);
171 if (result > MPI_CONGRUENT) {
172 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
173 "MoFEM and DM using different communicators");
174 }
175 MPI_Comm_size(comm, &dm_field->sIze);
176 MPI_Comm_rank(comm, &dm_field->rAnk);
177
178 // problem structure
179 CHKERR dm_field->mField_ptr->get_problem(dm_field->problemName,
180 &dm_field->problemPtr);
181
182 MPI_Comm_compare(comm, PETSC_COMM_SELF, &result);
183 if (result == MPI_IDENT) {
184 MOFEM_LOG("DMSELF", Sev::noisy)
185 << "MoFEM DM created for problem " << dm_field->problemName;
186 MOFEM_LOG("DMSELF", Sev::noisy) << *dm_field->problemPtr;
187 } else {
188 MOFEM_LOG("DMWORLD", Sev::noisy)
189 << "MoFEM DM created for problem " << dm_field->problemName;
190 MOFEM_LOG("DMWORLD", Sev::noisy) << *dm_field->problemPtr;
191 }
192
194}
195
196PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate) {
198
199 auto *dm_field = static_cast<DMCtx *>(dm->data);
200 if (!dm->data)
201 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
202 "data structure for MoFEM not yet created");
203
204 if (static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber == 0)
205 delete static_cast<DMCtx *>(dm_duplicate->data);
206
207 dm_duplicate->data = dm->data;
208 ++(static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber);
209
211}
212
213PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]) {
215
216 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
217 if (!dm->data) {
218 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
219 "data structure for MoFEM not yet created");
220 }
221 CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
222 dm_field->problemPtr->getBitRefLevel(),
223 dm_field->problemPtr->getBitRefLevelMask());
224
225 DMCtx *subdm_field = (DMCtx *)subdm->data;
226 subdm_field->isSubDM = PETSC_TRUE;
227 subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
228 subdm_field->isPartitioned = dm_field->isPartitioned;
229 subdm_field->isSquareMatrix = PETSC_FALSE;
230 subdm->ops->setup = DMSubDMSetUp_MoFEM;
231
233}
234
235PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
236 EntityType lo_type, EntityType hi_type) {
237 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
239 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
240 if (!dm->data) {
241 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
242 "data structure for MoFEM not yet created");
243 }
244 if (!dm_field->isSubDM) {
245 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
246 }
247 dm_field->rowFields.push_back(field_name);
248 if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
249 if (!dm_field->mapTypeRow)
250 dm_field->mapTypeRow = boost::make_shared<
251 std::map<std::string, std::pair<EntityType, EntityType>>>();
252 (*dm_field->mapTypeRow)[field_name] =
253 std::pair<EntityType, EntityType>(lo_type, hi_type);
254 }
256}
257
258PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
259 EntityType lo_type, EntityType hi_type) {
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->colFields.push_back(field_name);
271 if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
272 if (!dm_field->mapTypeCol)
273 dm_field->mapTypeCol = boost::make_shared<
274 std::map<std::string, std::pair<EntityType, EntityType>>>();
275 (*dm_field->mapTypeCol)[field_name] =
276 std::pair<EntityType, EntityType>(lo_type, hi_type);
277 }
279}
280
281PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm) {
283 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
285 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
286 *is_sub_dm = dm_field->isSubDM;
288}
289
290PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is) {
292 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
294 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
295 if (dm_field->isSubDM != PETSC_TRUE) {
296 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
297 "This DM is not created as a SubDM");
298 }
299 if (dm_field->isProblemBuild != PETSC_TRUE) {
300 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
301 }
302 boost::shared_ptr<Problem::SubProblemData> sub_data =
303 dm_field->problemPtr->getSubData();
304 CHKERR sub_data->getRowIs(is);
306}
307
308PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is) {
310 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
312 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
313 if (dm_field->isSubDM != PETSC_TRUE) {
314 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
315 "This DM is not created as a SubDM");
316 }
317 if (dm_field->isProblemBuild != PETSC_TRUE) {
318 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
319 }
320 boost::shared_ptr<Problem::SubProblemData> sub_data =
321 dm_field->problemPtr->getSubData();
322 CHKERR sub_data->getColIs(is);
324}
325
326PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]) {
327 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
329 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
330 if (!dm->data) {
331 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
332 "data structure for MoFEM not yet created");
333 }
334 if (!dm_field->isCompDM) {
335 dm_field->isCompDM = PETSC_TRUE;
336 }
337 dm_field->rowCompPrb.push_back(prb_name);
338 if (dm_field->isSquareMatrix) {
339 dm_field->colCompPrb.push_back(prb_name);
340 }
342}
343
344PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]) {
345 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
347 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
348 if (!dm->data) {
349 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
350 "data structure for MoFEM not yet created");
351 }
352 if (!dm_field->isCompDM) {
353 dm_field->isCompDM = PETSC_TRUE;
354 }
355 if (dm_field->isSquareMatrix) {
356 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
357 "No need to add problem on column when problem block structurally "
358 "symmetric");
359 }
360 dm_field->colCompPrb.push_back(prb_name);
362}
363
364PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm) {
366 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
368 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
369 *is_comp_dm = dm_field->isCompDM;
371}
372
373PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr) {
374 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
376 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
377 if (!dm->data) {
378 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
379 "data structure for MoFEM not yet created");
380 }
381 *m_field_ptr = dm_field->mField_ptr;
383}
384
385PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr) {
386 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
388 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
389 if (!dm->data) {
390 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
391 "data structure for MoFEM not yet created");
392 }
393 *problem_ptr = dm_field->problemPtr;
395}
396
397PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem) {
399 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
401 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
402 dm_field->destroyProblem = destroy_problem;
404}
405
406PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem) {
408 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
410 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
411 *destroy_problem = dm_field->destroyProblem;
413}
414
415PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem) {
417 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
419 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
420 dm_field->isSquareMatrix = square_problem;
422}
423
424PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[]) {
426 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
428 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
430 ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
432}
433
434PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[]) {
435 return DMMoFEMResolveSharedFiniteElements(dm, fe_name);
436}
437
438PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[],
439 PetscLayout *layout) {
441 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
443 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
444
445 MPI_Comm comm;
446 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
447 CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
448 layout);
450}
451
452PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem) {
455 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
457 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
458 *square_problem = dm_field->isSquareMatrix;
460}
461
462PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[]) {
463 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
465 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
467 dm_field->problemName, fe_name);
468 CHKERRG(ierr);
470}
471
472PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[]) {
473 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
475 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
477 dm_field->problemName, fe_name);
478 CHKERRG(ierr);
480}
481
482PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
483 ScatterMode scatter_mode) {
485 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
487 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
488 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
489 dm_field->problemPtr, COL, l, mode, scatter_mode);
490 CHKERRG(ierr);
492}
493
494PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
495 ScatterMode scatter_mode) {
496 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
498 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
499 ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
500 dm_field->problemPtr, COL, g, mode, scatter_mode);
501 CHKERRG(ierr);
503}
504
505PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
506 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
508 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
510 dm_field->problemPtr, *method);
511 CHKERRG(ierr);
513}
514
515PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
516 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
518 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
520 dm_field->problemPtr, *method);
521 CHKERRG(ierr);
523}
524
525PetscErrorCode
526DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
527 MoFEM::FEMethod *method, int low_rank,
528 int up_rank, CacheTupleWeakPtr cache_ptr) {
530 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
532 dm_field->problemPtr, fe_name, *method, low_rank, up_rank, nullptr,
533 MF_EXIST, cache_ptr);
534 CHKERRG(ierr);
536}
537
539 DM dm, const std::string fe_name, boost::shared_ptr<MoFEM::FEMethod> method,
540 int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr) {
541 return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
542 low_rank, up_rank, cache_ptr);
543}
544
545PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
546 MoFEM::FEMethod *method,
547 CacheTupleWeakPtr cache_ptr) {
548 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
550 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
552 dm, fe_name, method, dm_field->rAnk, dm_field->rAnk, cache_ptr);
553 CHKERRG(ierr);
555}
556
557PetscErrorCode
558DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
559 boost::shared_ptr<MoFEM::FEMethod> method,
560 CacheTupleWeakPtr cache_ptr) {
561 return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get(), cache_ptr);
562}
563
564PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
565 MoFEM::DofMethod *method) {
566 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
568 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
569 ierr =
570 dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
571 *method, dm_field->rAnk, dm_field->rAnk);
572 CHKERRG(ierr);
574}
575
576template <class S, class T0, class T1, class T2>
577static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method,
578 T1 pre_only, T2 post_only) {
579 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
581 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
582 if (pre_only) {
583 dm_field->kspCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
584 }
585 if (method) {
586 dm_field->kspCtx->get_loops_to_do_Rhs().push_back(
587 PairNameFEMethodPtr(fe_name, method));
588 }
589 if (post_only) {
590 dm_field->kspCtx->get_postProcess_to_do_Rhs().push_back(post_only);
591 }
592 CHKERR DMKSPSetComputeRHS(dm, KspRhs, dm_field->kspCtx.get());
594}
595
596PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
597 MoFEM::FEMethod *method,
598 MoFEM::BasicMethod *pre_only,
599 MoFEM::BasicMethod *post_only) {
600 return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
602 dm, fe_name, method, pre_only, post_only);
603}
604
605PetscErrorCode
606DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
607 boost::shared_ptr<MoFEM::FEMethod> method,
608 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
609 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
610 return DMMoFEMKSPSetComputeRHS<const std::string,
611 boost::shared_ptr<MoFEM::FEMethod>,
612 boost::shared_ptr<MoFEM::BasicMethod>,
613 boost::shared_ptr<MoFEM::BasicMethod>>(
614 dm, fe_name, method, pre_only, post_only);
615}
616
617template <class S, class T0, class T1, class T2>
618static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method,
619 T1 pre_only, T2 post_only) {
620 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
622 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
623 if (pre_only) {
624 dm_field->kspCtx->get_preProcess_to_do_Mat().push_back(pre_only);
625 }
626 if (method) {
627 dm_field->kspCtx->get_loops_to_do_Mat().push_back(
628 PairNameFEMethodPtr(fe_name, method));
629 }
630 if (post_only) {
631 dm_field->kspCtx->get_postProcess_to_do_Mat().push_back(post_only);
632 }
633 CHKERR DMKSPSetComputeOperators(dm, KspMat, dm_field->kspCtx.get());
635}
636
637PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
638 MoFEM::FEMethod *method,
639 MoFEM::BasicMethod *pre_only,
640 MoFEM::BasicMethod *post_only) {
641 return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
644 dm, fe_name, method, pre_only, post_only);
645}
646
647PetscErrorCode
648DMMoFEMKSPSetComputeOperators(DM dm, const std::string fe_name,
649 boost::shared_ptr<MoFEM::FEMethod> method,
650 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
651 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
652 return DMMoFEMKSPSetComputeOperators<const std::string,
653 boost::shared_ptr<MoFEM::FEMethod>>(
654 dm, fe_name, method, pre_only, post_only);
655}
656
657template <class S, class T0, class T1, class T2>
658static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method,
659 T1 pre_only, T2 post_only) {
660 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
662 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
663 if (pre_only) {
664 dm_field->snesCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
665 }
666 if (method) {
667 dm_field->snesCtx->get_loops_to_do_Rhs().push_back(
668 PairNameFEMethodPtr(fe_name, method));
669 }
670 if (post_only) {
671 dm_field->snesCtx->get_postProcess_to_do_Rhs().push_back(post_only);
672 }
673 CHKERR DMSNESSetFunction(dm, SnesRhs, dm_field->snesCtx.get());
675}
676
677PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
678 MoFEM::FEMethod *method,
679 MoFEM::BasicMethod *pre_only,
680 MoFEM::BasicMethod *post_only) {
681 return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
683 dm, fe_name, method, pre_only, post_only);
684}
685
686PetscErrorCode
687DMMoFEMSNESSetFunction(DM dm, const std::string fe_name,
688 boost::shared_ptr<MoFEM::FEMethod> method,
689 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
690 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
691 return DMMoFEMSNESSetFunction<const std::string,
692 boost::shared_ptr<MoFEM::FEMethod>,
693 boost::shared_ptr<MoFEM::BasicMethod>,
694 boost::shared_ptr<MoFEM::BasicMethod>>(
695 dm, fe_name, method, pre_only, post_only);
696}
697
698template <class S, class T0, class T1, class T2>
699static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method,
700 T1 pre_only, T2 post_only) {
701 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
703 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
704 if (pre_only) {
705 dm_field->snesCtx->get_preProcess_to_do_Mat().push_back(pre_only);
706 }
707 if (method) {
708 dm_field->snesCtx->get_loops_to_do_Mat().push_back(
709 PairNameFEMethodPtr(fe_name, method));
710 }
711 if (post_only) {
712 dm_field->snesCtx->get_postProcess_to_do_Mat().push_back(post_only);
713 }
714 CHKERR DMSNESSetJacobian(dm, SnesMat, dm_field->snesCtx.get());
716}
717
718PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
719 MoFEM::FEMethod *method,
720 MoFEM::BasicMethod *pre_only,
721 MoFEM::BasicMethod *post_only) {
722 return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
724 dm, fe_name, method, pre_only, post_only);
725}
726
727PetscErrorCode
728DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
729 boost::shared_ptr<MoFEM::FEMethod> method,
730 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
731 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
732 return DMMoFEMSNESSetJacobian<const std::string,
733 boost::shared_ptr<MoFEM::FEMethod>,
734 boost::shared_ptr<MoFEM::BasicMethod>,
735 boost::shared_ptr<MoFEM::BasicMethod>>(
736 dm, fe_name, method, pre_only, post_only);
737}
738
739template <class S, class T0, class T1, class T2>
740static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method,
741 T1 pre_only, T2 post_only) {
742 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
744 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
745 if (pre_only) {
746 dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
747 }
748 if (method) {
749 dm_field->tsCtx->getLoopsIFunction().push_back(
750 PairNameFEMethodPtr(fe_name, method));
751 }
752 if (post_only) {
753 dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
754 }
755 CHKERR DMTSSetIFunction(dm, TsSetIFunction, dm_field->tsCtx.get());
757}
758
759PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
760 MoFEM::FEMethod *method,
761 MoFEM::BasicMethod *pre_only,
762 MoFEM::BasicMethod *post_only) {
763 return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
765 dm, fe_name, method, pre_only, post_only);
767}
768
769PetscErrorCode
770DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
771 boost::shared_ptr<MoFEM::FEMethod> method,
772 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
773 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
774 return DMMoFEMTSSetIFunction<const std::string,
775 boost::shared_ptr<MoFEM::FEMethod>,
776 boost::shared_ptr<MoFEM::BasicMethod>,
777 boost::shared_ptr<MoFEM::BasicMethod>>(
778 dm, fe_name, method, pre_only, post_only);
780}
781
782template <class S, class T0, class T1, class T2>
783static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method,
784 T1 pre_only, T2 post_only) {
785 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
787 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
788 if (pre_only) {
789 dm_field->tsCtx->getPreProcessIJacobian().push_back(pre_only);
790 }
791 if (method) {
792 dm_field->tsCtx->getLoopsIJacobian().push_back(
793 PairNameFEMethodPtr(fe_name, method));
794 }
795 if (post_only) {
796 dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
797 }
798 CHKERR DMTSSetIJacobian(dm, TsSetIJacobian, dm_field->tsCtx.get());
800}
801
802PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
803 MoFEM::FEMethod *method,
804 MoFEM::BasicMethod *pre_only,
805 MoFEM::BasicMethod *post_only) {
806 return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
807 MoFEM::BasicMethod *>(dm, fe_name, method,
808 pre_only, post_only);
809}
810
811PetscErrorCode
812DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
813 boost::shared_ptr<MoFEM::FEMethod> method,
814 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
815 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
816 return DMMoFEMTSSetIJacobian<const std::string,
817 boost::shared_ptr<MoFEM::FEMethod>,
818 boost::shared_ptr<MoFEM::BasicMethod>,
819 boost::shared_ptr<MoFEM::BasicMethod>>(
820 dm, fe_name, method, pre_only, post_only);
821}
822
823template <class S, class T0, class T1, class T2>
824static PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, S fe_name, T0 method,
825 T1 pre_only, T2 post_only) {
826 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
828 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
829 if (pre_only)
830 dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
831 if (method)
832 dm_field->tsCtx->getLoopsRHSFunction().push_back(
833 PairNameFEMethodPtr(fe_name, method));
834 if (post_only)
835 dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
836 CHKERR DMTSSetRHSFunction(dm, TsSetRHSFunction, dm_field->tsCtx.get());
838}
839
840PetscErrorCode
841DMMoFEMTSSetRHSFunction(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 DMMoFEMTSSetRHSFunction<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);
851}
852
853PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const char fe_name[],
854 MoFEM::FEMethod *method,
855 MoFEM::BasicMethod *pre_only,
856 MoFEM::BasicMethod *post_only) {
857 return DMMoFEMTSSetRHSFunction<const char *, MoFEM::FEMethod *,
859 dm, fe_name, method, pre_only, post_only);
861}
862
863template <class S, class T0, class T1, class T2>
864static PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, S fe_name, T0 method,
865 T1 pre_only, T2 post_only) {
866 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
868 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
869 if (pre_only)
870 dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
871 if (method)
872 dm_field->tsCtx->getLoopsRHSFunction().push_back(
873 PairNameFEMethodPtr(fe_name, method));
874 if (post_only)
875 dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
876 CHKERR DMTSSetRHSJacobian(dm, TsSetRHSJacobian, dm_field->tsCtx.get());
878}
879
880PetscErrorCode
881DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
882 boost::shared_ptr<MoFEM::FEMethod> method,
883 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
884 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
885 return DMMoFEMTSSetRHSJacobian<const std::string,
886 boost::shared_ptr<MoFEM::FEMethod>,
887 boost::shared_ptr<MoFEM::BasicMethod>,
888 boost::shared_ptr<MoFEM::BasicMethod>>(
889 dm, fe_name, method, pre_only, post_only);
891}
892
893PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const char fe_name[],
894 MoFEM::FEMethod *method,
895 MoFEM::BasicMethod *pre_only,
896 MoFEM::BasicMethod *post_only) {
897 return DMMoFEMTSSetRHSJacobian<const char *, MoFEM::FEMethod *,
899 dm, fe_name, method, pre_only, post_only);
901}
902
903template <class S, class T0, class T1, class T2>
904static PetscErrorCode DMMoFEMTSSetI2Function(DM dm, S fe_name, T0 method,
905 T1 pre_only, T2 post_only) {
906 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
908 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
909 if (pre_only) {
910 dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
911 }
912 if (method) {
913 dm_field->tsCtx->getLoopsIFunction().push_back(
914 PairNameFEMethodPtr(fe_name, method));
915 }
916 if (post_only) {
917 dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
918 }
919 CHKERR DMTSSetI2Function(dm, TsSetI2Function, dm_field->tsCtx.get());
921}
922
923PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const char fe_name[],
924 MoFEM::FEMethod *method,
925 MoFEM::BasicMethod *pre_only,
926 MoFEM::BasicMethod *post_only) {
927 return DMMoFEMTSSetI2Function<const char *, MoFEM::FEMethod *,
929 dm, fe_name, method, pre_only, post_only);
931}
932
933PetscErrorCode
934DMMoFEMTSSetI2Function(DM dm, const std::string fe_name,
935 boost::shared_ptr<MoFEM::FEMethod> method,
936 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
937 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
938 return DMMoFEMTSSetI2Function<const std::string,
939 boost::shared_ptr<MoFEM::FEMethod>,
940 boost::shared_ptr<MoFEM::BasicMethod>,
941 boost::shared_ptr<MoFEM::BasicMethod>>(
942 dm, fe_name, method, pre_only, post_only);
944}
945
946template <class S, class T0, class T1, class T2>
947static PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, S fe_name, T0 method,
948 T1 pre_only, T2 post_only) {
949 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
951 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
952 if (pre_only) {
953 dm_field->tsCtx->getPreProcessIJacobian().push_back(pre_only);
954 }
955 if (method) {
956 dm_field->tsCtx->getLoopsIJacobian().push_back(
957 PairNameFEMethodPtr(fe_name, method));
958 }
959 if (post_only) {
960 dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
961 }
962 CHKERR DMTSSetI2Jacobian(dm, TsSetI2Jacobian, dm_field->tsCtx.get());
964}
965
966PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const char fe_name[],
967 MoFEM::FEMethod *method,
968 MoFEM::BasicMethod *pre_only,
969 MoFEM::BasicMethod *post_only) {
970 return DMMoFEMTSSetI2Jacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
971 MoFEM::BasicMethod *>(dm, fe_name, method,
972 pre_only, post_only);
973}
974
975PetscErrorCode
976DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name,
977 boost::shared_ptr<MoFEM::FEMethod> method,
978 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
979 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
980 return DMMoFEMTSSetI2Jacobian<const std::string,
981 boost::shared_ptr<MoFEM::FEMethod>,
982 boost::shared_ptr<MoFEM::BasicMethod>,
983 boost::shared_ptr<MoFEM::BasicMethod>>(
984 dm, fe_name, method, pre_only, post_only);
985}
986
987template <class S, class T0, class T1, class T2>
988static PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, S fe_name, T0 method,
989 T1 pre_only, T2 post_only) {
990 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
992 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
993 if (pre_only)
994 dm_field->tsCtx->getPreProcessMonitor().push_back(pre_only);
995 if (method)
996 dm_field->tsCtx->getLoopsMonitor().push_back(
997 PairNameFEMethodPtr(fe_name, method));
998 if (post_only)
999 dm_field->tsCtx->getPostProcessMonitor().push_back(post_only);
1000 CHKERR TSMonitorSet(ts, TsMonitorSet, dm_field->tsCtx.get(), PETSC_NULL);
1002}
1003
1004PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const char fe_name[],
1005 MoFEM::FEMethod *method,
1006 MoFEM::BasicMethod *pre_only,
1007 MoFEM::BasicMethod *post_only) {
1008 return DMMoFEMTSSetMonitor<const char *, MoFEM::FEMethod *,
1010 dm, ts, fe_name, method, pre_only, post_only);
1012}
1013
1014PetscErrorCode
1015DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name,
1016 boost::shared_ptr<MoFEM::FEMethod> method,
1017 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
1018 boost::shared_ptr<MoFEM::BasicMethod> post_only) {
1019 return DMMoFEMTSSetMonitor<const std::string,
1020 boost::shared_ptr<MoFEM::FEMethod>,
1021 boost::shared_ptr<MoFEM::BasicMethod>,
1022 boost::shared_ptr<MoFEM::BasicMethod>>(
1023 dm, ts, fe_name, method, pre_only, post_only);
1025}
1026
1027PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx) {
1028 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1030 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1031 *ksp_ctx = dm_field->kspCtx.get();
1033}
1034
1035PetscErrorCode
1036DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
1037 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1039 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1040 const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
1042}
1043
1044PetscErrorCode DMMoFEMSetKspCtx(DM dm,
1045 boost::shared_ptr<MoFEM::KspCtx> ksp_ctx) {
1046 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1048 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1049 dm_field->kspCtx = ksp_ctx;
1051}
1052
1053PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx) {
1054 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1056 DMCtx *dm_field = (DMCtx *)dm->data;
1057 *snes_ctx = dm_field->snesCtx.get();
1059}
1060
1061PetscErrorCode
1062DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
1063 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1065 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1066 const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
1068}
1069
1070PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
1071 boost::shared_ptr<MoFEM::SnesCtx> snes_ctx) {
1072 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1074 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1075 dm_field->snesCtx = snes_ctx;
1077}
1078
1079/** get if read mesh is partitioned
1080 * \ingroup dm
1081 */
1082PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned) {
1083 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1085 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1086 dm_field->isPartitioned = is_partitioned;
1088}
1089
1090/** get if read mesh is partitioned
1091 * \ingroup dm
1092 */
1093PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned) {
1094 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1096 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1097 *is_partitioned = dm_field->isPartitioned;
1099}
1100
1101PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx) {
1102 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1104 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1105 *ts_ctx = dm_field->tsCtx.get();
1107}
1108
1109PetscErrorCode DMMoFEMGetTsCtx(DM dm,
1110 const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
1111 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1113 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1114 const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1116}
1117
1118PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> ts_ctx) {
1119 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1121 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1122 dm_field->tsCtx = ts_ctx;
1124}
1125
1126PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g) {
1127 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1129 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1130 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1131 dm_field->problemName, COL, g);
1132 CHKERR VecSetDM(*g, dm);
1134}
1135
1136PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj<Vec> &g_ptr) {
1137 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1139 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1140 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1141 dm_field->problemName, COL, g_ptr);
1142 CHKERR VecSetDM(g_ptr, dm);
1144}
1145
1146PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l) {
1147 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1149 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1150 CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1151 dm_field->problemName, COL, l);
1152 CHKERR VecSetDM(*l, dm);
1154}
1155
1156PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M) {
1157 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1159 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1160 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1162 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1163 M);
1164 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1166 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1167 M);
1168 } else {
1169 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1170 "Matrix type not implemented");
1171 }
1172 CHKERR MatSetDM(*M, dm);
1174}
1175
1177 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1179 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1180 if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1182 ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1183 M);
1184 } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1186 ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1187 M);
1188 } else {
1189 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1190 "Matrix type not implemented");
1191 }
1192 CHKERR MatSetDM(M, dm);
1194}
1195
1196#if PETSC_VERSION_GE(3, 7, 0)
1197PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
1198 DM dm) {
1199#elif PETSC_VERSION_GE(3, 5, 3)
1200PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm) {
1201#else
1202PetscErrorCode DMSetFromOptions_MoFEM(DM dm) {
1203#endif
1204
1205 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1207 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1208#if PETSC_VERSION_GE(3, 5, 3)
1209 ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1210 CHKERRG(ierr);
1211#else
1212 ierr = PetscOptionsHead("DMMoFEM Options");
1213 CHKERRG(ierr);
1214#endif
1215 ierr = PetscOptionsBool("-dm_is_partitioned",
1216 "set if mesh is partitioned (works which native MOAB "
1217 "file format, i.e. h5m",
1218 "DMSetUp", dm_field->isPartitioned,
1219 &dm_field->isPartitioned, NULL);
1220 CHKERRG(ierr);
1222}
1223
1224PetscErrorCode DMSetUp_MoFEM(DM dm) {
1225 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1226 ProblemsManager *prb_mng_ptr;
1228 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1229 CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1230
1231 if (dm_field->isCompDM) {
1232 // It is composite probelm
1233 CHKERR prb_mng_ptr->buildComposedProblem(
1234 dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1235 dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1236 } else {
1237 if (dm_field->isPartitioned) {
1239 dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1240 dm_field->verbosity);
1241 } else {
1242 CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1243 dm_field->isSquareMatrix == PETSC_TRUE,
1244 dm_field->verbosity);
1245 CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1246 dm_field->verbosity);
1247 }
1248 }
1249
1250 // Partition finite elements
1251 if (dm_field->isPartitioned) {
1252 CHKERR prb_mng_ptr->partitionFiniteElements(
1253 dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1255 dm_field->problemName, dm_field->verbosity);
1256 } else {
1257 // partition finite elemnets
1258 CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1259 -1, -1, dm_field->verbosity);
1260 // Get ghost DOFs
1261 CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1262 dm_field->verbosity);
1263 }
1264
1265 // Set flag that problem is build and partitioned
1266 dm_field->isProblemBuild = PETSC_TRUE;
1267
1269}
1270
1271PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm) {
1272 PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1273 ProblemsManager *prb_mng_ptr;
1275
1276 DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1277
1278 // build sub dm problem
1279 CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1280
1281 map<std::string, std::pair<EntityType, EntityType>> *entity_map_row = nullptr;
1282 map<std::string, std::pair<EntityType, EntityType>> *entity_map_col = nullptr;
1283
1284 if (subdm_field->mapTypeRow)
1285 entity_map_row = subdm_field->mapTypeRow.get();
1286 if (subdm_field->mapTypeCol)
1287 entity_map_row = subdm_field->mapTypeCol.get();
1288
1289 CHKERR prb_mng_ptr->buildSubProblem(
1290 subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1291 subdm_field->problemMainOfSubPtr->getName(),
1292 subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1293 subdm_field->verbosity);
1294
1295 // partition problem
1296 subdm_field->isPartitioned = subdm_field->isPartitioned;
1297 if (subdm_field->isPartitioned) {
1298 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1299 0, subdm_field->sIze,
1300 subdm_field->verbosity);
1301 // set ghost nodes
1303 subdm_field->problemName, subdm_field->verbosity);
1304 } else {
1305 // partition finite elements
1306 CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1307 -1, -1, subdm_field->verbosity);
1308 // set ghost nodes
1309 CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1310 subdm_field->verbosity);
1311 }
1312
1313 subdm_field->isProblemBuild = PETSC_TRUE;
1314
1316}
1317
1318PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec g, InsertMode mode,
1319 Vec l) {
1320 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1322 CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1324}
1325
1326PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec g, InsertMode mode, Vec l) {
1328 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1330
1331 CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1332
1333 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1334 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1335 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1336
1337 double *array_loc, *array_glob;
1338 CHKERR VecGetArray(l, &array_loc);
1339 CHKERR VecGetArray(g, &array_glob);
1340 switch (mode) {
1341 case INSERT_VALUES:
1342 cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1343 break;
1344 case ADD_VALUES:
1345 cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1346 break;
1347 default:
1348 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1349 }
1350 CHKERR VecRestoreArray(l, &array_loc);
1351 CHKERR VecRestoreArray(g, &array_glob);
1353}
1354
1355PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM dm, Vec l, InsertMode mode,
1356 Vec g) {
1357
1358 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1360
1361 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1362 int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1363 int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1364
1365 double *array_loc, *array_glob;
1366 CHKERR VecGetArray(l, &array_loc);
1367 CHKERR VecGetArray(g, &array_glob);
1368 switch (mode) {
1369 case INSERT_VALUES:
1370 cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1371 break;
1372 case ADD_VALUES:
1373 cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1374 break;
1375 default:
1376 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1377 }
1378 CHKERR VecRestoreArray(l, &array_loc);
1379 CHKERR VecRestoreArray(g, &array_glob);
1380
1382}
1383
1384PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec l, InsertMode mode, Vec g) {
1385 //
1388}
1389
1390PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
1391 char ***fieldNames, IS **fields) {
1392 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1394
1395 if (numFields) {
1396 *numFields = 0;
1397 }
1398 if (fieldNames) {
1399 *fieldNames = NULL;
1400 }
1401 if (fields) {
1402 *fields = NULL;
1403 }
1404
1405 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1406 auto fields_ptr = dm_field->mField_ptr->get_fields();
1407 Field_multiIndex::iterator fit, hi_fit;
1408 fit = fields_ptr->begin();
1409 hi_fit = fields_ptr->end();
1410 *numFields = std::distance(fit, hi_fit);
1411
1412 if (fieldNames) {
1413 CHKERR PetscMalloc1(*numFields, fieldNames);
1414 }
1415 if (fields) {
1416 CHKERR PetscMalloc1(*numFields, fields);
1417 }
1418
1419 for (int f = 0; fit != hi_fit; fit++, f++) {
1420 if (fieldNames) {
1421 CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1422 (char **)&(*fieldNames)[f]);
1423 }
1424 if (fields) {
1426 ->isCreateProblemFieldAndRank(
1427 dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1428 fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1429 }
1430 }
1431
1433}
1434
1435PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
1436 IS *is) {
1437 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1439 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1441 ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1442 field_name, 0, 1000, is);
1444}
1445
1446PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb) {
1447 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1449 DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1450 dm_field->verbosity = verb;
1452}
1453
1454} // namespace MoFEM
static Index< 'M', 3 > M
Discrete manager interface for MoFEM.
@ VERBOSE
Definition: definitions.h:222
RowColData
RowColData.
Definition: definitions.h:136
@ COL
Definition: definitions.h:136
@ ROW
Definition: definitions.h:136
@ MF_EXIST
Definition: definitions.h:113
@ MF_EXCL
Definition: definitions.h:112
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:496
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
@ MOFEM_INVALID_DATA
Definition: definitions.h:49
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1318
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1326
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:424
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:1082
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
Definition: DMMMoFEM.cpp:213
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMMoFEM.cpp:1093
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[], PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMMoFEM.cpp:438
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMMoFEM.cpp:677
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMMoFEM.cpp:1146
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition: DMMMoFEM.cpp:415
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:759
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:235
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:95
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:130
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:462
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:515
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:482
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on col.
Definition: DMMMoFEM.cpp:344
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1355
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:385
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:87
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > snes_ctx)
Set MoFEM::SnesCtx data structure.
Definition: DMMMoFEM.cpp:1070
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > ksp_ctx)
set MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:1044
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMMoFEM.cpp:281
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on row.
Definition: DMMMoFEM.cpp:326
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMMoFEM.cpp:718
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1156
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1202
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMMoFEM.cpp:1101
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:258
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[])
unset element from dm
Definition: DMMMoFEM.cpp:472
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:1027
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problem
Definition: DMMMoFEM.cpp:452
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMMoFEM.cpp:564
PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set KSP right hand side evaluation function
Definition: DMMMoFEM.cpp:596
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > ts_ctx)
Set MoFEM::TsCtx data structure.
Definition: DMMMoFEM.cpp:1118
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:545
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMMoFEM.cpp:1435
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:812
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:373
PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:976
PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:934
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1126
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMMoFEM.cpp:1053
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:494
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition: DMMMoFEM.cpp:364
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1271
PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set KSP operators and push mofem finite element methods.
Definition: DMMMoFEM.cpp:637
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1390
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1224
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1384
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:526
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:505
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:65
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:364
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:342
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
virtual MoFEMErrorCode problem_basic_method_postProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
virtual MoFEMErrorCode loop_finite_elements(const std::string &problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
MoFEMErrorCode buildProblemOnDistributedMesh(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures, assuming that mesh is distributed (collective)
MoFEMErrorCode partitionGhostDofs(const std::string name, int verb=VERBOSE)
determine ghost nodes
MoFEMErrorCode buildComposedProblem(const std::string out_name, const std::vector< std::string > add_row_problems, const std::vector< std::string > add_col_problems, const bool square_matrix=true, int verb=1)
build composite problem
MoFEMErrorCode buildSubProblem(const std::string out_name, const std::vector< std::string > &fields_row, const std::vector< std::string > &fields_col, const std::string main_problem, const bool square_matrix=true, const map< std::string, std::pair< EntityType, EntityType > > *entityMapRow=nullptr, const map< std::string, std::pair< EntityType, EntityType > > *entityMapCol=nullptr, int verb=VERBOSE)
build sub problem
MoFEMErrorCode buildProblem(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
MoFEMErrorCode partitionProblem(const std::string name, int verb=VERBOSE)
partition problem dofs (collective)
MoFEMErrorCode partitionGhostDofsOnDistributedMesh(const std::string name, int verb=VERBOSE)
determine ghost nodes on distributed meshes
MoFEMErrorCode partitionFiniteElements(const std::string name, bool part_from_moab=false, int low_proc=-1, int hi_proc=-1, int verb=VERBOSE)
partition finite elements
virtual MoFEMErrorCode modify_problem_unset_finite_element(const std::string &name_problem, const std::string &fe_name)=0
unset finite element from problem, this remove entities assigned to finite element to a particular pr...
virtual MoFEMErrorCode modify_problem_add_finite_element(const std::string &name_problem, const std::string &fe_name)=0
add finite element to problem, this add entities assigned to finite element to a particular problem
virtual MoFEMErrorCode add_problem(const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add problem.
virtual MoFEMErrorCode delete_problem(const std::string name)=0
Delete problem.
virtual MoFEMErrorCode modify_problem_ref_level_add_bit(const std::string &name_problem, const BitRefLevel &bit)=0
add ref level to problem
virtual bool check_problem(const std::string name)=0
check if problem exist
virtual MoFEMErrorCode modify_problem_mask_ref_level_add_bit(const std::string &name_problem, const BitRefLevel &bit)=0
set dof mask ref level for problem
FTensor::Index< 'l', 3 > l
const FTensor::Tensor2< T, Dim, Dim > Vec
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
Definition: DMMMoFEM.cpp:1015
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:140
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:229
PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS the right hand side function
Definition: DMMMoFEM.cpp:841
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:148
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:397
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMMoFEM.cpp:406
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:569
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:287
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:39
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:290
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:37
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:308
PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx)
Run over elements in the lists.
Definition: KspCtx.cpp:33
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:473
PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx)
Run over elenents in the list.
Definition: KspCtx.cpp:97
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:385
PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS the right hand side jacobian
Definition: DMMMoFEM.cpp:881
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition: DMMMoFEM.cpp:1446
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
Definition: DMMMoFEM.cpp:196
DEPRECATED PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[])
Definition: DMMMoFEM.cpp:434
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:904
PetscBool isCompDM
Definition: DMMoFEM.hpp:928
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:911
const Problem * problemPtr
pinter to problem data structure
Definition: DMMoFEM.hpp:920
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:929
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition: DMMoFEM.hpp:945
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition: DMMoFEM.hpp:946
PetscBool isSquareMatrix
true if rows equals to cols
Definition: DMMoFEM.hpp:915
PetscBool isPartitioned
true if read mesh is on parts
Definition: DMMoFEM.hpp:914
PetscBool isSubDM
Definition: DMMoFEM.hpp:923
int referenceNumber
Definition: DMMoFEM.hpp:942
int verbosity
verbosity
Definition: DMMoFEM.hpp:941
PetscBool isProblemBuild
True if problem is build.
Definition: DMMoFEM.hpp:910
Interface * mField_ptr
MoFEM interface.
Definition: DMMoFEM.hpp:909
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition: DMMoFEM.hpp:944
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeRow
Definition: DMMoFEM.hpp:932
std::vector< std::string > colCompPrb
Definition: DMMoFEM.hpp:930
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: DMMMoFEM.cpp:53
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:925
PetscBool destroyProblem
If true destroy problem with DM.
Definition: DMMoFEM.hpp:936
std::vector< std::string > rowFields
Definition: DMMoFEM.hpp:924
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
Definition: DMMoFEM.hpp:926
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeCol
Definition: DMMoFEM.hpp:934
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:33
Interface for linear (KSP) solver.
Definition: KspCtx.hpp:27
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:327
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:374
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:319
Matrix manager is used to build and partition problems.
keeps basic data about problem
DofIdx getNbLocalDofsRow() const
MoFEMErrorCode getNumberOfElementsByNameAndPart(MPI_Comm comm, const std::string name, PetscLayout *layout) const
Get number of finite elements by name on processors.
BitRefLevel getBitRefLevel() const
BitRefLevel getBitRefLevelMask() const
boost::shared_ptr< SubProblemData > & getSubData() const
Get main problem of sub-problem is.
DofIdx getNbGhostDofsRow() const
Problem manager is used to build and partition problems.
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:27
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
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:33