v0.12.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 
28 namespace 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 
53 MoFEMErrorCode DMCtx::query_interface(boost::typeindex::type_index type_index,
54  UnknownInterface **iface) const {
55  *iface = const_cast<DMCtx *>(this);
56  return 0;
57 }
58 
59 PetscErrorCode DMRegister_MoFEM(const char sname[]) {
61  CHKERR DMRegister(sname, DMCreate_MoFEM);
63 }
64 
65 PetscErrorCode 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 
87 PetscErrorCode DMCreate_MoFEM(DM dm) {
88  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
90  dm->data = new DMCtx();
93 }
94 
95 PetscErrorCode 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 
130 PetscErrorCode 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  else
187  MOFEM_LOG("DMWORLD", Sev::noisy)
188  << "MoFEM DM created for problem " << dm_field->problemName;
189 
191 }
192 
193 PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate) {
195 
196  auto *dm_field = static_cast<DMCtx *>(dm->data);
197  if (!dm->data)
198  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
199  "data structure for MoFEM not yet created");
200 
201  if (static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber == 0)
202  delete static_cast<DMCtx *>(dm_duplicate->data);
203 
204  dm_duplicate->data = dm->data;
205  ++(static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber);
206 
208 }
209 
210 PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]) {
212 
213  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
214  if (!dm->data) {
215  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
216  "data structure for MoFEM not yet created");
217  }
218  CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
219  dm_field->problemPtr->getBitRefLevel());
220 
221  DMCtx *subdm_field = (DMCtx *)subdm->data;
222  subdm_field->isSubDM = PETSC_TRUE;
223  subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
224  subdm_field->isPartitioned = dm_field->isPartitioned;
225  subdm_field->isSquareMatrix = PETSC_FALSE;
226  subdm->ops->setup = DMSubDMSetUp_MoFEM;
227 
229 }
230 
231 PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
232  EntityType lo_type, EntityType hi_type) {
233  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
235  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
236  if (!dm->data) {
237  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
238  "data structure for MoFEM not yet created");
239  }
240  if (!dm_field->isSubDM) {
241  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
242  }
243  dm_field->rowFields.push_back(field_name);
244  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
245  if (!dm_field->mapTypeRow)
246  dm_field->mapTypeRow = boost::make_shared<
247  std::map<std::string, std::pair<EntityType, EntityType>>>();
248  (*dm_field->mapTypeRow)[field_name] =
249  std::pair<EntityType, EntityType>(lo_type, hi_type);
250  }
252 }
253 
254 PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
255  EntityType lo_type, EntityType hi_type) {
256  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
258  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
259  if (!dm->data) {
260  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
261  "data structure for MoFEM not yet created");
262  }
263  if (!dm_field->isSubDM) {
264  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
265  }
266  dm_field->colFields.push_back(field_name);
267  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
268  if (!dm_field->mapTypeCol)
269  dm_field->mapTypeCol = boost::make_shared<
270  std::map<std::string, std::pair<EntityType, EntityType>>>();
271  (*dm_field->mapTypeCol)[field_name] =
272  std::pair<EntityType, EntityType>(lo_type, hi_type);
273  }
275 }
276 
277 PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm) {
279  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
281  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
282  *is_sub_dm = dm_field->isSubDM;
284 }
285 
286 PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is) {
288  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
290  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
291  if (dm_field->isSubDM != PETSC_TRUE) {
292  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
293  "This DM is not created as a SubDM");
294  }
295  if (dm_field->isProblemBuild != PETSC_TRUE) {
296  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
297  }
298  boost::shared_ptr<Problem::SubProblemData> sub_data =
299  dm_field->problemPtr->getSubData();
300  CHKERR sub_data->getRowIs(is);
302 }
303 
304 PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is) {
306  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
308  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
309  if (dm_field->isSubDM != PETSC_TRUE) {
310  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
311  "This DM is not created as a SubDM");
312  }
313  if (dm_field->isProblemBuild != PETSC_TRUE) {
314  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
315  }
316  boost::shared_ptr<Problem::SubProblemData> sub_data =
317  dm_field->problemPtr->getSubData();
318  CHKERR sub_data->getColIs(is);
320 }
321 
322 PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]) {
323  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
325  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
326  if (!dm->data) {
327  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
328  "data structure for MoFEM not yet created");
329  }
330  if (!dm_field->isCompDM) {
331  dm_field->isCompDM = PETSC_TRUE;
332  }
333  dm_field->rowCompPrb.push_back(prb_name);
334  if (dm_field->isSquareMatrix) {
335  dm_field->colCompPrb.push_back(prb_name);
336  }
338 }
339 
340 PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]) {
341  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
343  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
344  if (!dm->data) {
345  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
346  "data structure for MoFEM not yet created");
347  }
348  if (!dm_field->isCompDM) {
349  dm_field->isCompDM = PETSC_TRUE;
350  }
351  if (dm_field->isSquareMatrix) {
352  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
353  "No need to add problem on column when problem block structurally "
354  "symmetric");
355  }
356  dm_field->colCompPrb.push_back(prb_name);
358 }
359 
360 PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm) {
362  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
364  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
365  *is_comp_dm = dm_field->isCompDM;
367 }
368 
369 PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr) {
370  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
372  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
373  if (!dm->data) {
374  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
375  "data structure for MoFEM not yet created");
376  }
377  *m_field_ptr = dm_field->mField_ptr;
379 }
380 
381 PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr) {
382  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
384  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
385  if (!dm->data) {
386  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
387  "data structure for MoFEM not yet created");
388  }
389  *problem_ptr = dm_field->problemPtr;
391 }
392 
393 PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem) {
395  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
397  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
398  dm_field->destroyProblem = destroy_problem;
400 }
401 
402 PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem) {
404  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
406  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
407  *destroy_problem = dm_field->destroyProblem;
409 }
410 
411 PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem) {
413  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
415  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
416  dm_field->isSquareMatrix = square_problem;
418 }
419 
420 PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[]) {
422  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
424  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
426  ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
428 }
429 
430 PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[]) {
431  return DMMoFEMResolveSharedFiniteElements(dm, fe_name);
432 }
433 
434 PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[],
435  PetscLayout *layout) {
437  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
439  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
440 
441  MPI_Comm comm;
442  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
443  CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
444  layout);
446 }
447 
448 PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem) {
451  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
453  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
454  *square_problem = dm_field->isSquareMatrix;
456 }
457 
458 PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[]) {
459  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
461  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
463  dm_field->problemName, fe_name);
464  CHKERRG(ierr);
466 }
467 
468 PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[]) {
469  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
471  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
473  dm_field->problemName, fe_name);
474  CHKERRG(ierr);
476 }
477 
478 PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
479  ScatterMode scatter_mode) {
481  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
483  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
484  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
485  dm_field->problemPtr, COL, l, mode, scatter_mode);
486  CHKERRG(ierr);
488 }
489 
490 PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
491  ScatterMode scatter_mode) {
492  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
494  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
495  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
496  dm_field->problemPtr, COL, g, mode, scatter_mode);
497  CHKERRG(ierr);
499 }
500 
501 PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
502  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
504  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
506  dm_field->problemPtr, *method);
507  CHKERRG(ierr);
509 }
510 
511 PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
512  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
514  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
516  dm_field->problemPtr, *method);
517  CHKERRG(ierr);
519 }
520 
521 PetscErrorCode
522 DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
523  MoFEM::FEMethod *method, int low_rank,
524  int up_rank, CacheTupleWeakPtr cache_ptr) {
526  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
527  ierr = dm_field->mField_ptr->loop_finite_elements(
528  dm_field->problemPtr, fe_name, *method, low_rank, up_rank, nullptr,
529  MF_EXIST, cache_ptr);
530  CHKERRG(ierr);
532 }
533 
535  DM dm, const std::string fe_name, boost::shared_ptr<MoFEM::FEMethod> method,
536  int low_rank, int up_rank, CacheTupleWeakPtr cache_ptr) {
537  return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
538  low_rank, up_rank, cache_ptr);
539 }
540 
541 PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
542  MoFEM::FEMethod *method,
543  CacheTupleWeakPtr cache_ptr) {
544  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
546  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
548  dm, fe_name, method, dm_field->rAnk, dm_field->rAnk, cache_ptr);
549  CHKERRG(ierr);
551 }
552 
553 PetscErrorCode
554 DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
555  boost::shared_ptr<MoFEM::FEMethod> method,
556  CacheTupleWeakPtr cache_ptr) {
557  return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get(), cache_ptr);
558 }
559 
560 PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
561  MoFEM::DofMethod *method) {
562  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
564  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
565  ierr =
566  dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
567  *method, dm_field->rAnk, dm_field->rAnk);
568  CHKERRG(ierr);
570 }
571 
572 template <class S, class T0, class T1, class T2>
573 static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method,
574  T1 pre_only, T2 post_only) {
575  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
577  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
578  if (pre_only) {
579  dm_field->kspCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
580  }
581  if (method) {
582  dm_field->kspCtx->get_loops_to_do_Rhs().push_back(
583  PairNameFEMethodPtr(fe_name, method));
584  }
585  if (post_only) {
586  dm_field->kspCtx->get_postProcess_to_do_Rhs().push_back(post_only);
587  }
588  CHKERR DMKSPSetComputeRHS(dm, KspRhs, dm_field->kspCtx.get());
590 }
591 
592 PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
593  MoFEM::FEMethod *method,
594  MoFEM::BasicMethod *pre_only,
595  MoFEM::BasicMethod *post_only) {
596  return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
598  dm, fe_name, method, pre_only, post_only);
599 }
600 
601 PetscErrorCode
602 DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
603  boost::shared_ptr<MoFEM::FEMethod> method,
604  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
605  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
606  return DMMoFEMKSPSetComputeRHS<const std::string,
607  boost::shared_ptr<MoFEM::FEMethod>,
608  boost::shared_ptr<MoFEM::BasicMethod>,
609  boost::shared_ptr<MoFEM::BasicMethod>>(
610  dm, fe_name, method, pre_only, post_only);
611 }
612 
613 template <class S, class T0, class T1, class T2>
614 static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method,
615  T1 pre_only, T2 post_only) {
616  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
618  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
619  if (pre_only) {
620  dm_field->kspCtx->get_preProcess_to_do_Mat().push_back(pre_only);
621  }
622  if (method) {
623  dm_field->kspCtx->get_loops_to_do_Mat().push_back(
624  PairNameFEMethodPtr(fe_name, method));
625  }
626  if (post_only) {
627  dm_field->kspCtx->get_postProcess_to_do_Mat().push_back(post_only);
628  }
629  CHKERR DMKSPSetComputeOperators(dm, KspMat, dm_field->kspCtx.get());
631 }
632 
633 PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
634  MoFEM::FEMethod *method,
635  MoFEM::BasicMethod *pre_only,
636  MoFEM::BasicMethod *post_only) {
637  return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
640  dm, fe_name, method, pre_only, post_only);
641 }
642 
643 PetscErrorCode
644 DMMoFEMKSPSetComputeOperators(DM dm, const std::string fe_name,
645  boost::shared_ptr<MoFEM::FEMethod> method,
646  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
647  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
648  return DMMoFEMKSPSetComputeOperators<const std::string,
649  boost::shared_ptr<MoFEM::FEMethod>>(
650  dm, fe_name, method, pre_only, post_only);
651 }
652 
653 template <class S, class T0, class T1, class T2>
654 static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method,
655  T1 pre_only, T2 post_only) {
656  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
658  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
659  if (pre_only) {
660  dm_field->snesCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
661  }
662  if (method) {
663  dm_field->snesCtx->get_loops_to_do_Rhs().push_back(
664  PairNameFEMethodPtr(fe_name, method));
665  }
666  if (post_only) {
667  dm_field->snesCtx->get_postProcess_to_do_Rhs().push_back(post_only);
668  }
669  CHKERR DMSNESSetFunction(dm, SnesRhs, dm_field->snesCtx.get());
671 }
672 
673 PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
674  MoFEM::FEMethod *method,
675  MoFEM::BasicMethod *pre_only,
676  MoFEM::BasicMethod *post_only) {
677  return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
679  dm, fe_name, method, pre_only, post_only);
680 }
681 
682 PetscErrorCode
683 DMMoFEMSNESSetFunction(DM dm, const std::string fe_name,
684  boost::shared_ptr<MoFEM::FEMethod> method,
685  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
686  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
687  return DMMoFEMSNESSetFunction<const std::string,
688  boost::shared_ptr<MoFEM::FEMethod>,
689  boost::shared_ptr<MoFEM::BasicMethod>,
690  boost::shared_ptr<MoFEM::BasicMethod>>(
691  dm, fe_name, method, pre_only, post_only);
692 }
693 
694 template <class S, class T0, class T1, class T2>
695 static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method,
696  T1 pre_only, T2 post_only) {
697  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
699  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
700  if (pre_only) {
701  dm_field->snesCtx->get_preProcess_to_do_Mat().push_back(pre_only);
702  }
703  if (method) {
704  dm_field->snesCtx->get_loops_to_do_Mat().push_back(
705  PairNameFEMethodPtr(fe_name, method));
706  }
707  if (post_only) {
708  dm_field->snesCtx->get_postProcess_to_do_Mat().push_back(post_only);
709  }
710  CHKERR DMSNESSetJacobian(dm, SnesMat, dm_field->snesCtx.get());
712 }
713 
714 PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
715  MoFEM::FEMethod *method,
716  MoFEM::BasicMethod *pre_only,
717  MoFEM::BasicMethod *post_only) {
718  return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
720  dm, fe_name, method, pre_only, post_only);
721 }
722 
723 PetscErrorCode
724 DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
725  boost::shared_ptr<MoFEM::FEMethod> method,
726  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
727  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
728  return DMMoFEMSNESSetJacobian<const std::string,
729  boost::shared_ptr<MoFEM::FEMethod>,
730  boost::shared_ptr<MoFEM::BasicMethod>,
731  boost::shared_ptr<MoFEM::BasicMethod>>(
732  dm, fe_name, method, pre_only, post_only);
733 }
734 
735 template <class S, class T0, class T1, class T2>
736 static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method,
737  T1 pre_only, T2 post_only) {
738  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
740  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
741  if (pre_only) {
742  dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
743  }
744  if (method) {
745  dm_field->tsCtx->getLoopsIFunction().push_back(
746  PairNameFEMethodPtr(fe_name, method));
747  }
748  if (post_only) {
749  dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
750  }
751  CHKERR DMTSSetIFunction(dm, TsSetIFunction, dm_field->tsCtx.get());
753 }
754 
755 PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
756  MoFEM::FEMethod *method,
757  MoFEM::BasicMethod *pre_only,
758  MoFEM::BasicMethod *post_only) {
759  return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
761  dm, fe_name, method, pre_only, post_only);
763 }
764 
765 PetscErrorCode
766 DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
767  boost::shared_ptr<MoFEM::FEMethod> method,
768  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
769  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
770  return DMMoFEMTSSetIFunction<const std::string,
771  boost::shared_ptr<MoFEM::FEMethod>,
772  boost::shared_ptr<MoFEM::BasicMethod>,
773  boost::shared_ptr<MoFEM::BasicMethod>>(
774  dm, fe_name, method, pre_only, post_only);
776 }
777 
778 template <class S, class T0, class T1, class T2>
779 static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method,
780  T1 pre_only, T2 post_only) {
781  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
783  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
784  if (pre_only) {
785  dm_field->tsCtx->getPreProcessIJacobian().push_back(pre_only);
786  }
787  if (method) {
788  dm_field->tsCtx->getLoopsIJacobian().push_back(
789  PairNameFEMethodPtr(fe_name, method));
790  }
791  if (post_only) {
792  dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
793  }
794  CHKERR DMTSSetIJacobian(dm, TsSetIJacobian, dm_field->tsCtx.get());
796 }
797 
798 PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
799  MoFEM::FEMethod *method,
800  MoFEM::BasicMethod *pre_only,
801  MoFEM::BasicMethod *post_only) {
802  return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
803  MoFEM::BasicMethod *>(dm, fe_name, method,
804  pre_only, post_only);
805 }
806 
807 PetscErrorCode
808 DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
809  boost::shared_ptr<MoFEM::FEMethod> method,
810  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
811  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
812  return DMMoFEMTSSetIJacobian<const std::string,
813  boost::shared_ptr<MoFEM::FEMethod>,
814  boost::shared_ptr<MoFEM::BasicMethod>,
815  boost::shared_ptr<MoFEM::BasicMethod>>(
816  dm, fe_name, method, pre_only, post_only);
817 }
818 
819 template <class S, class T0, class T1, class T2>
820 static PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, S fe_name, T0 method,
821  T1 pre_only, T2 post_only) {
822  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
824  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
825  if (pre_only)
826  dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
827  if (method)
828  dm_field->tsCtx->getLoopsRHSFunction().push_back(
829  PairNameFEMethodPtr(fe_name, method));
830  if (post_only)
831  dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
832  CHKERR DMTSSetRHSFunction(dm, TsSetRHSFunction, dm_field->tsCtx.get());
834 }
835 
836 PetscErrorCode
837 DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name,
838  boost::shared_ptr<MoFEM::FEMethod> method,
839  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
840  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
841  return DMMoFEMTSSetRHSFunction<const std::string,
842  boost::shared_ptr<MoFEM::FEMethod>,
843  boost::shared_ptr<MoFEM::BasicMethod>,
844  boost::shared_ptr<MoFEM::BasicMethod>>(
845  dm, fe_name, method, pre_only, post_only);
847 }
848 
849 template <class S, class T0, class T1, class T2>
850 static PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, S fe_name, T0 method,
851  T1 pre_only, T2 post_only) {
852  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
854  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
855  if (pre_only)
856  dm_field->tsCtx->getPreProcessRHSFunction().push_back(pre_only);
857  if (method)
858  dm_field->tsCtx->getLoopsRHSFunction().push_back(
859  PairNameFEMethodPtr(fe_name, method));
860  if (post_only)
861  dm_field->tsCtx->getPostProcessRHSFunction().push_back(post_only);
862  CHKERR DMTSSetRHSJacobian(dm, TsSetRHSJacobian, dm_field->tsCtx.get());
864 }
865 
866 PetscErrorCode
867 DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
868  boost::shared_ptr<MoFEM::FEMethod> method,
869  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
870  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
871  return DMMoFEMTSSetRHSJacobian<const std::string,
872  boost::shared_ptr<MoFEM::FEMethod>,
873  boost::shared_ptr<MoFEM::BasicMethod>,
874  boost::shared_ptr<MoFEM::BasicMethod>>(
875  dm, fe_name, method, pre_only, post_only);
877 }
878 
879 template <class S, class T0, class T1, class T2>
880 static PetscErrorCode DMMoFEMTSSetI2Function(DM dm, S fe_name, T0 method,
881  T1 pre_only, T2 post_only) {
882  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
884  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
885  if (pre_only) {
886  dm_field->tsCtx->getPreProcessIFunction().push_back(pre_only);
887  }
888  if (method) {
889  dm_field->tsCtx->getLoopsIFunction().push_back(
890  PairNameFEMethodPtr(fe_name, method));
891  }
892  if (post_only) {
893  dm_field->tsCtx->getPostProcessIFunction().push_back(post_only);
894  }
895  CHKERR DMTSSetI2Function(dm, TsSetI2Function, dm_field->tsCtx.get());
897 }
898 
899 PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const char fe_name[],
900  MoFEM::FEMethod *method,
901  MoFEM::BasicMethod *pre_only,
902  MoFEM::BasicMethod *post_only) {
903  return DMMoFEMTSSetI2Function<const char *, MoFEM::FEMethod *,
905  dm, fe_name, method, pre_only, post_only);
907 }
908 
909 PetscErrorCode
910 DMMoFEMTSSetI2Function(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 DMMoFEMTSSetI2Function<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 
922 template <class S, class T0, class T1, class T2>
923 static PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, S fe_name, T0 method,
924  T1 pre_only, T2 post_only) {
925  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
927  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
928  if (pre_only) {
929  dm_field->tsCtx->getPreProcessIJacobian().push_back(pre_only);
930  }
931  if (method) {
932  dm_field->tsCtx->getLoopsIJacobian().push_back(
933  PairNameFEMethodPtr(fe_name, method));
934  }
935  if (post_only) {
936  dm_field->tsCtx->getPostProcessIJacobian().push_back(post_only);
937  }
938  CHKERR DMTSSetI2Jacobian(dm, TsSetI2Jacobian, dm_field->tsCtx.get());
940 }
941 
942 PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const char fe_name[],
943  MoFEM::FEMethod *method,
944  MoFEM::BasicMethod *pre_only,
945  MoFEM::BasicMethod *post_only) {
946  return DMMoFEMTSSetI2Jacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
947  MoFEM::BasicMethod *>(dm, fe_name, method,
948  pre_only, post_only);
949 }
950 
951 PetscErrorCode
952 DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name,
953  boost::shared_ptr<MoFEM::FEMethod> method,
954  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
955  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
956  return DMMoFEMTSSetI2Jacobian<const std::string,
957  boost::shared_ptr<MoFEM::FEMethod>,
958  boost::shared_ptr<MoFEM::BasicMethod>,
959  boost::shared_ptr<MoFEM::BasicMethod>>(
960  dm, fe_name, method, pre_only, post_only);
961 }
962 
963 template <class S, class T0, class T1, class T2>
964 static PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, S fe_name, T0 method,
965  T1 pre_only, T2 post_only) {
966  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
968  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
969  if (pre_only)
970  dm_field->tsCtx->getPreProcessMonitor().push_back(pre_only);
971  if (method)
972  dm_field->tsCtx->getLoopsMonitor().push_back(
973  PairNameFEMethodPtr(fe_name, method));
974  if (post_only)
975  dm_field->tsCtx->getPostProcessMonitor().push_back(post_only);
976  CHKERR TSMonitorSet(ts, TsMonitorSet, dm_field->tsCtx.get(), PETSC_NULL);
978 }
979 
980 PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const char fe_name[],
981  MoFEM::FEMethod *method,
982  MoFEM::BasicMethod *pre_only,
983  MoFEM::BasicMethod *post_only) {
984  return DMMoFEMTSSetMonitor<const char *, MoFEM::FEMethod *,
986  dm, ts, fe_name, method, pre_only, post_only);
988 }
989 
990 PetscErrorCode
991 DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name,
992  boost::shared_ptr<MoFEM::FEMethod> method,
993  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
994  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
995  return DMMoFEMTSSetMonitor<const std::string,
996  boost::shared_ptr<MoFEM::FEMethod>,
997  boost::shared_ptr<MoFEM::BasicMethod>,
998  boost::shared_ptr<MoFEM::BasicMethod>>(
999  dm, ts, fe_name, method, pre_only, post_only);
1001 }
1002 
1003 PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx) {
1004  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1006  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1007  *ksp_ctx = dm_field->kspCtx.get();
1009 }
1010 
1011 PetscErrorCode
1012 DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
1013  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1015  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1016  const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
1018 }
1019 
1020 PetscErrorCode DMMoFEMSetKspCtx(DM dm,
1021  boost::shared_ptr<MoFEM::KspCtx> ksp_ctx) {
1022  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1024  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1025  dm_field->kspCtx = ksp_ctx;
1027 }
1028 
1029 PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx) {
1030  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1032  DMCtx *dm_field = (DMCtx *)dm->data;
1033  *snes_ctx = dm_field->snesCtx.get();
1035 }
1036 
1037 PetscErrorCode
1038 DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
1039  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1041  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1042  const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
1044 }
1045 
1046 PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
1047  boost::shared_ptr<MoFEM::SnesCtx> snes_ctx) {
1048  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1050  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1051  dm_field->snesCtx = snes_ctx;
1053 }
1054 
1055 /** get if read mesh is partitioned
1056  * \ingroup dm
1057  */
1058 PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned) {
1059  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1061  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1062  dm_field->isPartitioned = is_partitioned;
1064 }
1065 
1066 /** get if read mesh is partitioned
1067  * \ingroup dm
1068  */
1069 PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned) {
1070  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1072  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1073  *is_partitioned = dm_field->isPartitioned;
1075 }
1076 
1077 PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx) {
1078  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1080  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1081  *ts_ctx = dm_field->tsCtx.get();
1083 }
1084 
1085 PetscErrorCode DMMoFEMGetTsCtx(DM dm,
1086  const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
1087  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1089  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1090  const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1092 }
1093 
1094 PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> ts_ctx) {
1095  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1097  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1098  dm_field->tsCtx = ts_ctx;
1100 }
1101 
1102 PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g) {
1103  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1105  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1106  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1107  dm_field->problemName, COL, g);
1108  CHKERR VecSetDM(*g, dm);
1110 }
1111 
1112 PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj<Vec> &g_ptr) {
1113  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1115  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1116  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1117  dm_field->problemName, COL, g_ptr);
1118  CHKERR VecSetDM(g_ptr, dm);
1120 }
1121 
1122 PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l) {
1123  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1125  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1126  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1127  dm_field->problemName, COL, l);
1128  CHKERR VecSetDM(*l, dm);
1130 }
1131 
1132 PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M) {
1133  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1135  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1136  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1138  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1139  M);
1140  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1142  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1143  M);
1144  } else {
1145  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1146  "Matrix type not implemented");
1147  }
1148  CHKERR MatSetDM(*M, dm);
1150 }
1151 
1152 PetscErrorCode DMCreateMatrix_MoFEM(DM dm, SmartPetscObj<Mat> &M) {
1153  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1155  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1156  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1158  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1159  M);
1160  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1162  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1163  M);
1164  } else {
1165  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1166  "Matrix type not implemented");
1167  }
1168  CHKERR MatSetDM(M, dm);
1170 }
1171 
1172 #if PETSC_VERSION_GE(3, 7, 0)
1173 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
1174  DM dm) {
1175 #elif PETSC_VERSION_GE(3, 5, 3)
1176 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm) {
1177 #else
1178 PetscErrorCode DMSetFromOptions_MoFEM(DM dm) {
1179 #endif
1180 
1181  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1183  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1184 #if PETSC_VERSION_GE(3, 5, 3)
1185  ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1186  CHKERRG(ierr);
1187 #else
1188  ierr = PetscOptionsHead("DMMoFEM Options");
1189  CHKERRG(ierr);
1190 #endif
1191  ierr = PetscOptionsBool("-dm_is_partitioned",
1192  "set if mesh is partitioned (works which native MOAB "
1193  "file format, i.e. h5m",
1194  "DMSetUp", dm_field->isPartitioned,
1195  &dm_field->isPartitioned, NULL);
1196  CHKERRG(ierr);
1198 }
1199 
1200 PetscErrorCode DMSetUp_MoFEM(DM dm) {
1201  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1202  ProblemsManager *prb_mng_ptr;
1204  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1205  CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1206 
1207  if (dm_field->isCompDM) {
1208  // It is composite probelm
1209  CHKERR prb_mng_ptr->buildCompsedProblem(
1210  dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1211  dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1212  } else {
1213  if (dm_field->isPartitioned) {
1215  dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1216  dm_field->verbosity);
1217  } else {
1218  CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1219  dm_field->isSquareMatrix == PETSC_TRUE,
1220  dm_field->verbosity);
1221  CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1222  dm_field->verbosity);
1223  }
1224  }
1225 
1226  // Partition finite elements
1227  if (dm_field->isPartitioned) {
1228  CHKERR prb_mng_ptr->partitionFiniteElements(
1229  dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1231  dm_field->problemName, dm_field->verbosity);
1232  } else {
1233  // partition finite elemnets
1234  CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1235  -1, -1, dm_field->verbosity);
1236  // Get ghost DOFs
1237  CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1238  dm_field->verbosity);
1239  }
1240 
1241  // Set flag that problem is build and partitioned
1242  dm_field->isProblemBuild = PETSC_TRUE;
1243 
1245 }
1246 
1247 PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm) {
1248  PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1249  ProblemsManager *prb_mng_ptr;
1251 
1252  DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1253 
1254  // build sub dm problem
1255  CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1256 
1257  map<std::string, std::pair<EntityType, EntityType>> *entity_map_row = nullptr;
1258  map<std::string, std::pair<EntityType, EntityType>> *entity_map_col = nullptr;
1259 
1260  if (subdm_field->mapTypeRow)
1261  entity_map_row = subdm_field->mapTypeRow.get();
1262  if (subdm_field->mapTypeCol)
1263  entity_map_row = subdm_field->mapTypeCol.get();
1264 
1265  CHKERR prb_mng_ptr->buildSubProblem(
1266  subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1267  subdm_field->problemMainOfSubPtr->getName(),
1268  subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1269  subdm_field->verbosity);
1270 
1271  // partition problem
1272  subdm_field->isPartitioned = subdm_field->isPartitioned;
1273  if (subdm_field->isPartitioned) {
1274  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1275  0, subdm_field->sIze,
1276  subdm_field->verbosity);
1277  // set ghost nodes
1279  subdm_field->problemName, subdm_field->verbosity);
1280  } else {
1281  // partition finite elements
1282  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1283  -1, -1, subdm_field->verbosity);
1284  // set ghost nodes
1285  CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1286  subdm_field->verbosity);
1287  }
1288 
1289  subdm_field->isProblemBuild = PETSC_TRUE;
1290 
1292 }
1293 
1294 PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec g, InsertMode mode,
1295  Vec l) {
1296  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1298  CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1300 }
1301 
1302 PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec g, InsertMode mode, Vec l) {
1304  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1306 
1307  CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1308 
1309  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1310  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1311  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1312 
1313  double *array_loc, *array_glob;
1314  CHKERR VecGetArray(l, &array_loc);
1315  CHKERR VecGetArray(g, &array_glob);
1316  switch (mode) {
1317  case INSERT_VALUES:
1318  cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1319  break;
1320  case ADD_VALUES:
1321  cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1322  break;
1323  default:
1324  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1325  }
1326  CHKERR VecRestoreArray(l, &array_loc);
1327  CHKERR VecRestoreArray(g, &array_glob);
1329 }
1330 
1331 PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM dm, Vec l, InsertMode mode,
1332  Vec g) {
1333 
1334  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1336 
1337  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1338  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1339  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1340 
1341  double *array_loc, *array_glob;
1342  CHKERR VecGetArray(l, &array_loc);
1343  CHKERR VecGetArray(g, &array_glob);
1344  switch (mode) {
1345  case INSERT_VALUES:
1346  cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1347  break;
1348  case ADD_VALUES:
1349  cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1350  break;
1351  default:
1352  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1353  }
1354  CHKERR VecRestoreArray(l, &array_loc);
1355  CHKERR VecRestoreArray(g, &array_glob);
1356 
1358 }
1359 
1360 PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec l, InsertMode mode, Vec g) {
1361  //
1364 }
1365 
1366 PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
1367  char ***fieldNames, IS **fields) {
1368  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1370 
1371  if (numFields) {
1372  *numFields = 0;
1373  }
1374  if (fieldNames) {
1375  *fieldNames = NULL;
1376  }
1377  if (fields) {
1378  *fields = NULL;
1379  }
1380 
1381  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1382  auto fields_ptr = dm_field->mField_ptr->get_fields();
1383  Field_multiIndex::iterator fit, hi_fit;
1384  fit = fields_ptr->begin();
1385  hi_fit = fields_ptr->end();
1386  *numFields = std::distance(fit, hi_fit);
1387 
1388  if (fieldNames) {
1389  CHKERR PetscMalloc1(*numFields, fieldNames);
1390  }
1391  if (fields) {
1392  CHKERR PetscMalloc1(*numFields, fields);
1393  }
1394 
1395  for (int f = 0; fit != hi_fit; fit++, f++) {
1396  if (fieldNames) {
1397  CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1398  (char **)&(*fieldNames)[f]);
1399  }
1400  if (fields) {
1401  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1402  ->isCreateProblemFieldAndRank(
1403  dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1404  fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1405  }
1406  }
1407 
1409 }
1410 
1411 PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
1412  IS *is) {
1413  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1415  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1416  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1417  ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1418  field_name, 0, 1000, is);
1420 }
1421 
1422 PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb) {
1423  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1425  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1426  dm_field->verbosity = verb;
1428 }
1429 
1430 } // 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:1294
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1302
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:420
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:1058
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:210
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMMoFEM.cpp:1069
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[], PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMMoFEM.cpp:434
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:673
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMMoFEM.cpp:1122
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition: DMMMoFEM.cpp:411
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:755
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:231
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:458
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:511
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:478
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on col.
Definition: DMMMoFEM.cpp:340
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1331
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:381
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:1046
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > ksp_ctx)
set MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:1020
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMMoFEM.cpp:277
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on row.
Definition: DMMMoFEM.cpp:322
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:714
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1132
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1178
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMMoFEM.cpp:1077
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:254
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:468
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:1003
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problem
Definition: DMMMoFEM.cpp:448
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMMoFEM.cpp:560
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:592
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > ts_ctx)
Set MoFEM::TsCtx data structure.
Definition: DMMMoFEM.cpp:1094
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:541
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMMoFEM.cpp:1411
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:808
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:369
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:952
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:910
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1102
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMMoFEM.cpp:1029
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:490
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition: DMMMoFEM.cpp:360
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1247
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:633
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1366
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1200
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1360
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:522
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:501
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:65
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem 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:309
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:340
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.
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 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 buildCompsedProblem(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 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:991
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:120
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:209
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:837
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:126
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:393
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMMoFEM.cpp:402
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:549
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:267
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:17
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:286
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:17
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:304
PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx)
Run over elements in the lists.
Definition: KspCtx.cpp:23
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:453
PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx)
Run over elenents in the list.
Definition: KspCtx.cpp:87
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:365
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:867
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition: DMMMoFEM.cpp:1422
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
Definition: DMMMoFEM.cpp:193
DEPRECATED PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[])
Definition: DMMMoFEM.cpp:430
boost::weak_ptr< CacheTuple > CacheTupleWeakPtr
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:894
PetscBool isCompDM
Definition: DMMoFEM.hpp:918
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:901
const Problem * problemPtr
pinter to problem data structure
Definition: DMMoFEM.hpp:910
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:919
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition: DMMoFEM.hpp:935
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition: DMMoFEM.hpp:936
PetscBool isSquareMatrix
true if rows equals to cols
Definition: DMMoFEM.hpp:905
PetscBool isPartitioned
true if read mesh is on parts
Definition: DMMoFEM.hpp:904
PetscBool isSubDM
Definition: DMMoFEM.hpp:913
int referenceNumber
Definition: DMMoFEM.hpp:932
int verbosity
verbosity
Definition: DMMoFEM.hpp:931
PetscBool isProblemBuild
True if problem is build.
Definition: DMMoFEM.hpp:900
Interface * mField_ptr
MoFEM interface.
Definition: DMMoFEM.hpp:899
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition: DMMoFEM.hpp:934
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeRow
Definition: DMMoFEM.hpp:922
std::vector< std::string > colCompPrb
Definition: DMMoFEM.hpp:920
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: DMMMoFEM.cpp:53
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:915
PetscBool destroyProblem
If true destroy problem with DM.
Definition: DMMoFEM.hpp:926
std::vector< std::string > rowFields
Definition: DMMoFEM.hpp:914
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
Definition: DMMoFEM.hpp:916
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeCol
Definition: DMMoFEM.hpp:924
Deprecated interface functions.
DEPRECATED MoFEMErrorCode loop_finite_elements(const Problem *problem_ptr, const std::string &fe_name, FEMethod &method, int lower_rank, int upper_rank, MoFEMTypes bh, int verb=DEFAULT_VERBOSITY)
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
std::string getName() const
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
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