v0.10.0
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 
37  UnknownInterface **iface) const {
39  *iface = NULL;
40  if (uuid == IDD_DMCTX) {
41  *iface = const_cast<DMCtx *>(this);
43  }
44  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
46 }
47 
48 PetscErrorCode DMRegister_MoFEM(const char sname[]) {
50  CHKERR DMRegister(sname, DMCreate_MoFEM);
52 }
53 
54 PetscErrorCode DMSetOperators_MoFEM(DM dm) {
55  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
57 
58  dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
59  dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
60  dm->ops->creatematrix = DMCreateMatrix_MoFEM;
61  dm->ops->setup = DMSetUp_MoFEM;
62  dm->ops->destroy = DMDestroy_MoFEM;
63  dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
64  dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
65  dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
66  dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
67  dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
68  dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
69 
70  // Default matrix type
71  CHKERR DMSetMatType(dm, MATMPIAIJ);
72 
74 }
75 
76 PetscErrorCode DMCreate_MoFEM(DM dm) {
77  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
79  dm->data = new DMCtx();
82 }
83 
84 PetscErrorCode DMDestroy_MoFEM(DM dm) {
85  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
86  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
88  if (static_cast<DMCtx *>(dm->data)->referenceNumber == 0) {
89  if (dm_field->destroyProblem) {
90 
91  if (dm_field->mField_ptr->check_problem(dm_field->problemName)) {
92  dm_field->mField_ptr->delete_problem(dm_field->problemName);
93  } // else problem has to be deleted by the user
94 
95  }
96 
97  delete static_cast<DMCtx *>(dm->data);
98 
99  } else
100  --static_cast<DMCtx *>(dm->data)->referenceNumber;
101 
103 }
104 
105 PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr,
106  const char problem_name[],
107  const MoFEM::BitRefLevel bit_level,
108  const MoFEM::BitRefLevel bit_mask) {
110 
111  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
112  if (!dm->data) {
113  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
114  "data structure for MoFEM not yet created");
115  }
116  if (!m_field_ptr) {
117  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
118  "DM function not implemented into MoFEM");
119  }
120  dm_field->mField_ptr = m_field_ptr;
121  dm_field->problemName = problem_name;
122  if (!m_field_ptr->check_problem(dm_field->problemName)) {
123  // problem is not defined, declare problem internally set bool to
124  // destroyProblem problem with DM
125  dm_field->destroyProblem = PETSC_TRUE;
126  CHKERR dm_field->mField_ptr->add_problem(dm_field->problemName, MF_EXCL,
127  dm_field->verbosity);
128  } else {
129  dm_field->destroyProblem = PETSC_FALSE;
130  }
132  dm_field->problemName, bit_level);
134  dm_field->problemName, bit_mask);
135  dm_field->kspCtx =
136  boost::shared_ptr<KspCtx>(new KspCtx(*m_field_ptr, problem_name));
137  dm_field->snesCtx =
138  boost::shared_ptr<SnesCtx>(new SnesCtx(*m_field_ptr, problem_name));
139  dm_field->tsCtx =
140  boost::shared_ptr<TsCtx>(new TsCtx(*m_field_ptr, problem_name));
141 
142  MPI_Comm comm;
143  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
144  int result = 0;
145  MPI_Comm_compare(comm, m_field_ptr->get_comm(), &result);
146  if (result > MPI_CONGRUENT) {
147  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
148  "MoFEM and DM using different communicators");
149  }
150  MPI_Comm_size(comm, &dm_field->sIze);
151  MPI_Comm_rank(comm, &dm_field->rAnk);
152 
153  // problem structure
154  CHKERR dm_field->mField_ptr->get_problem(dm_field->problemName,
155  &dm_field->problemPtr);
156 
158 }
159 
160 PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate) {
162 
163  auto *dm_field = static_cast<DMCtx *>(dm->data);
164  if (!dm->data)
165  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
166  "data structure for MoFEM not yet created");
167 
168  if (static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber == 0)
169  delete static_cast<DMCtx *>(dm_duplicate->data);
170 
171  dm_duplicate->data = dm->data;
172  ++(static_cast<DMCtx *>(dm_duplicate->data)->referenceNumber);
173 
175 }
176 
177 PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]) {
179 
180  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
181  if (!dm->data) {
182  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
183  "data structure for MoFEM not yet created");
184  }
185  CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
186  dm_field->problemPtr->getBitRefLevel());
187 
188  DMCtx *subdm_field = (DMCtx *)subdm->data;
189  subdm_field->isSubDM = PETSC_TRUE;
190  subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
191  subdm_field->isPartitioned = dm_field->isPartitioned;
192  subdm_field->isSquareMatrix = PETSC_FALSE;
193  subdm->ops->setup = DMSubDMSetUp_MoFEM;
194 
196 }
197 
198 PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
199  EntityType lo_type, EntityType hi_type) {
200  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
202  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
203  if (!dm->data) {
204  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
205  "data structure for MoFEM not yet created");
206  }
207  if (!dm_field->isSubDM) {
208  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
209  }
210  dm_field->rowFields.push_back(field_name);
211  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
212  if (!dm_field->mapTypeRow)
213  dm_field->mapTypeRow = boost::make_shared<
214  std::map<std::string, std::pair<EntityType, EntityType>>>();
215  (*dm_field->mapTypeRow)[field_name] =
216  std::pair<EntityType, EntityType>(lo_type, hi_type);
217  }
219 }
220 
221 PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
222  EntityType lo_type, EntityType hi_type) {
223  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
225  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
226  if (!dm->data) {
227  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
228  "data structure for MoFEM not yet created");
229  }
230  if (!dm_field->isSubDM) {
231  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
232  }
233  dm_field->colFields.push_back(field_name);
234  if (lo_type != MBVERTEX || hi_type != MBMAXTYPE) {
235  if (!dm_field->mapTypeCol)
236  dm_field->mapTypeCol = boost::make_shared<
237  std::map<std::string, std::pair<EntityType, EntityType>>>();
238  (*dm_field->mapTypeCol)[field_name] =
239  std::pair<EntityType, EntityType>(lo_type, hi_type);
240  }
242 }
243 
244 PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm) {
246  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
248  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
249  *is_sub_dm = dm_field->isSubDM;
251 }
252 
253 PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is) {
255  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
257  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
258  if (dm_field->isSubDM != PETSC_TRUE) {
259  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
260  "This DM is not created as a SubDM");
261  }
262  if (dm_field->isProblemBuild != PETSC_TRUE) {
263  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
264  }
265  boost::shared_ptr<Problem::SubProblemData> sub_data =
266  dm_field->problemPtr->getSubData();
267  CHKERR sub_data->getRowIs(is);
269 }
270 
271 PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is) {
273  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
275  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
276  if (dm_field->isSubDM != PETSC_TRUE) {
277  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
278  "This DM is not created as a SubDM");
279  }
280  if (dm_field->isProblemBuild != PETSC_TRUE) {
281  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
282  }
283  boost::shared_ptr<Problem::SubProblemData> sub_data =
284  dm_field->problemPtr->getSubData();
285  CHKERR sub_data->getColIs(is);
287 }
288 
289 PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]) {
290  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
292  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
293  if (!dm->data) {
294  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
295  "data structure for MoFEM not yet created");
296  }
297  if (!dm_field->isCompDM) {
298  dm_field->isCompDM = PETSC_TRUE;
299  }
300  dm_field->rowCompPrb.push_back(prb_name);
301  if (dm_field->isSquareMatrix) {
302  dm_field->colCompPrb.push_back(prb_name);
303  }
305 }
306 
307 PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]) {
308  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
310  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
311  if (!dm->data) {
312  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
313  "data structure for MoFEM not yet created");
314  }
315  if (!dm_field->isCompDM) {
316  dm_field->isCompDM = PETSC_TRUE;
317  }
318  if (dm_field->isSquareMatrix) {
319  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
320  "No need to add problem on column when problem block structurally "
321  "symmetric");
322  }
323  dm_field->colCompPrb.push_back(prb_name);
325 }
326 
327 PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm) {
329  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
331  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
332  *is_comp_dm = dm_field->isCompDM;
334 }
335 
336 PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr) {
337  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
339  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
340  if (!dm->data) {
341  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
342  "data structure for MoFEM not yet created");
343  }
344  *m_field_ptr = dm_field->mField_ptr;
346 }
347 
348 PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr) {
349  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
351  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
352  if (!dm->data) {
353  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
354  "data structure for MoFEM not yet created");
355  }
356  *problem_ptr = dm_field->problemPtr;
358 }
359 
360 PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem) {
362  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
364  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
365  dm_field->destroyProblem = destroy_problem;
367 }
368 
369 PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem) {
371  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
373  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
374  *destroy_problem = dm_field->destroyProblem;
376 }
377 
378 PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem) {
380  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
382  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
383  dm_field->isSquareMatrix = square_problem;
385 }
386 
387 PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[]) {
389  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
391  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
393  ->resolveSharedFiniteElements(dm_field->problemPtr, fe_name);
395 }
396 
397 PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[]) {
398  return DMMoFEMResolveSharedFiniteElements(dm, fe_name);
399 }
400 
401 PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[],
402  PetscLayout *layout) {
404  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
406  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
407 
408  MPI_Comm comm;
409  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
410  CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
411  layout);
413 }
414 
415 PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem) {
418  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
420  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
421  *square_problem = dm_field->isSquareMatrix;
423 }
424 
425 PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[]) {
426  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
428  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
430  dm_field->problemName, fe_name);
431  CHKERRG(ierr);
433 }
434 
435 PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[]) {
436  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
438  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
440  dm_field->problemName, fe_name);
441  CHKERRG(ierr);
443 }
444 
445 PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
446  ScatterMode scatter_mode) {
448  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
450  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
451  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
452  dm_field->problemPtr, COL, l, mode, scatter_mode);
453  CHKERRG(ierr);
455 }
456 
457 PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
458  ScatterMode scatter_mode) {
459  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
461  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
462  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
463  dm_field->problemPtr, COL, g, mode, scatter_mode);
464  CHKERRG(ierr);
466 }
467 
468 PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
469  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
471  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
473  dm_field->problemPtr, *method);
474  CHKERRG(ierr);
476 }
477 
478 PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
479  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
481  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
483  dm_field->problemPtr, *method);
484  CHKERRG(ierr);
486 }
487 
488 PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
489  MoFEM::FEMethod *method,
490  int low_rank, int up_rank) {
492  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
493  ierr = dm_field->mField_ptr->loop_finite_elements(
494  dm_field->problemPtr, fe_name, *method, low_rank, up_rank);
495  CHKERRG(ierr);
497 }
498 
499 PetscErrorCode
500 DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const std::string fe_name,
501  boost::shared_ptr<MoFEM::FEMethod> method,
502  int low_rank, int up_rank) {
503  return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
504  low_rank, up_rank);
505 }
506 
507 PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
508  MoFEM::FEMethod *method) {
509  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
511  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
512  ierr = DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name, method,
513  dm_field->rAnk, dm_field->rAnk);
514  CHKERRG(ierr);
516 }
517 
518 PetscErrorCode
519 DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
520  boost::shared_ptr<MoFEM::FEMethod> method) {
521  return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get());
522 }
523 
524 PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
525  MoFEM::DofMethod *method) {
526  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
528  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
529  ierr =
530  dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
531  *method, dm_field->rAnk, dm_field->rAnk);
532  CHKERRG(ierr);
534 }
535 
536 template <class S, class T0, class T1, class T2>
537 static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method,
538  T1 pre_only, T2 post_only) {
539  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
541  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
542  if (pre_only) {
543  dm_field->kspCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
544  }
545  if (method) {
546  dm_field->kspCtx->get_loops_to_do_Rhs().push_back(
547  PairNameFEMethodPtr(fe_name, method));
548  }
549  if (post_only) {
550  dm_field->kspCtx->get_postProcess_to_do_Rhs().push_back(post_only);
551  }
552  CHKERR DMKSPSetComputeRHS(dm, KspRhs, dm_field->kspCtx.get());
554 }
555 
556 PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
557  MoFEM::FEMethod *method,
558  MoFEM::BasicMethod *pre_only,
559  MoFEM::BasicMethod *post_only) {
560  return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
562  dm, fe_name, method, pre_only, post_only);
563 }
564 
565 PetscErrorCode
566 DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
567  boost::shared_ptr<MoFEM::FEMethod> method,
568  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
569  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
570  return DMMoFEMKSPSetComputeRHS<const std::string,
571  boost::shared_ptr<MoFEM::FEMethod>,
572  boost::shared_ptr<MoFEM::BasicMethod>,
573  boost::shared_ptr<MoFEM::BasicMethod>>(
574  dm, fe_name, method, pre_only, post_only);
575 }
576 
577 template <class S, class T0, class T1, class T2>
578 static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method,
579  T1 pre_only, T2 post_only) {
580  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
582  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
583  if (pre_only) {
584  dm_field->kspCtx->get_preProcess_to_do_Mat().push_back(pre_only);
585  }
586  if (method) {
587  dm_field->kspCtx->get_loops_to_do_Mat().push_back(
588  PairNameFEMethodPtr(fe_name, method));
589  }
590  if (post_only) {
591  dm_field->kspCtx->get_postProcess_to_do_Mat().push_back(post_only);
592  }
593  CHKERR DMKSPSetComputeOperators(dm, KspMat, dm_field->kspCtx.get());
595 }
596 
597 PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
598  MoFEM::FEMethod *method,
599  MoFEM::BasicMethod *pre_only,
600  MoFEM::BasicMethod *post_only) {
601  return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
604  dm, fe_name, method, pre_only, post_only);
605 }
606 
607 PetscErrorCode
608 DMMoFEMKSPSetComputeOperators(DM dm, const std::string fe_name,
609  boost::shared_ptr<MoFEM::FEMethod> method,
610  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
611  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
612  return DMMoFEMKSPSetComputeOperators<const std::string,
613  boost::shared_ptr<MoFEM::FEMethod>>(
614  dm, fe_name, method, pre_only, post_only);
615 }
616 
617 template <class S, class T0, class T1, class T2>
618 static PetscErrorCode DMMoFEMSNESSetFunction(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->snesCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
625  }
626  if (method) {
627  dm_field->snesCtx->get_loops_to_do_Rhs().push_back(
628  PairNameFEMethodPtr(fe_name, method));
629  }
630  if (post_only) {
631  dm_field->snesCtx->get_postProcess_to_do_Rhs().push_back(post_only);
632  }
633  CHKERR DMSNESSetFunction(dm, SnesRhs, dm_field->snesCtx.get());
635 }
636 
637 PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
638  MoFEM::FEMethod *method,
639  MoFEM::BasicMethod *pre_only,
640  MoFEM::BasicMethod *post_only) {
641  return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
643  dm, fe_name, method, pre_only, post_only);
644 }
645 
646 PetscErrorCode
647 DMMoFEMSNESSetFunction(DM dm, const std::string fe_name,
648  boost::shared_ptr<MoFEM::FEMethod> method,
649  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
650  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
651  return DMMoFEMSNESSetFunction<const std::string,
652  boost::shared_ptr<MoFEM::FEMethod>,
653  boost::shared_ptr<MoFEM::BasicMethod>,
654  boost::shared_ptr<MoFEM::BasicMethod>>(
655  dm, fe_name, method, pre_only, post_only);
656 }
657 
658 template <class S, class T0, class T1, class T2>
659 static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method,
660  T1 pre_only, T2 post_only) {
661  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
663  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
664  if (pre_only) {
665  dm_field->snesCtx->get_preProcess_to_do_Mat().push_back(pre_only);
666  }
667  if (method) {
668  dm_field->snesCtx->get_loops_to_do_Mat().push_back(
669  PairNameFEMethodPtr(fe_name, method));
670  }
671  if (post_only) {
672  dm_field->snesCtx->get_postProcess_to_do_Mat().push_back(post_only);
673  }
674  CHKERR DMSNESSetJacobian(dm, SnesMat, dm_field->snesCtx.get());
676 }
677 
678 PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
679  MoFEM::FEMethod *method,
680  MoFEM::BasicMethod *pre_only,
681  MoFEM::BasicMethod *post_only) {
682  return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
684  dm, fe_name, method, pre_only, post_only);
685 }
686 
687 PetscErrorCode
688 DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
689  boost::shared_ptr<MoFEM::FEMethod> method,
690  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
691  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
692  return DMMoFEMSNESSetJacobian<const std::string,
693  boost::shared_ptr<MoFEM::FEMethod>,
694  boost::shared_ptr<MoFEM::BasicMethod>,
695  boost::shared_ptr<MoFEM::BasicMethod>>(
696  dm, fe_name, method, pre_only, post_only);
697 }
698 
699 template <class S, class T0, class T1, class T2>
700 static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method,
701  T1 pre_only, T2 post_only) {
702  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
704  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
705  if (pre_only) {
706  dm_field->tsCtx->get_preProcess_to_do_IFunction().push_back(pre_only);
707  }
708  if (method) {
709  dm_field->tsCtx->get_loops_to_do_IFunction().push_back(
710  PairNameFEMethodPtr(fe_name, method));
711  }
712  if (post_only) {
713  dm_field->tsCtx->get_postProcess_to_do_IFunction().push_back(post_only);
714  }
715  CHKERR DMTSSetIFunction(dm, TsSetIFunction, dm_field->tsCtx.get());
717 }
718 
719 PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
720  MoFEM::FEMethod *method,
721  MoFEM::BasicMethod *pre_only,
722  MoFEM::BasicMethod *post_only) {
723  return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
725  dm, fe_name, method, pre_only, post_only);
727 }
728 
729 PetscErrorCode
730 DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
731  boost::shared_ptr<MoFEM::FEMethod> method,
732  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
733  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
734  return DMMoFEMTSSetIFunction<const std::string,
735  boost::shared_ptr<MoFEM::FEMethod>,
736  boost::shared_ptr<MoFEM::BasicMethod>,
737  boost::shared_ptr<MoFEM::BasicMethod>>(
738  dm, fe_name, method, pre_only, post_only);
740 }
741 
742 template <class S, class T0, class T1, class T2>
743 static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method,
744  T1 pre_only, T2 post_only) {
745  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
747  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
748  if (pre_only) {
749  dm_field->tsCtx->get_preProcess_to_do_IJacobian().push_back(pre_only);
750  }
751  if (method) {
752  dm_field->tsCtx->get_loops_to_do_IJacobian().push_back(
753  PairNameFEMethodPtr(fe_name, method));
754  }
755  if (post_only) {
756  dm_field->tsCtx->get_postProcess_to_do_IJacobian().push_back(post_only);
757  }
758  CHKERR DMTSSetIJacobian(dm, TsSetIJacobian, dm_field->tsCtx.get());
760 }
761 
762 PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
763  MoFEM::FEMethod *method,
764  MoFEM::BasicMethod *pre_only,
765  MoFEM::BasicMethod *post_only) {
766  return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
767  MoFEM::BasicMethod *>(dm, fe_name, method,
768  pre_only, post_only);
769 }
770 
771 PetscErrorCode
772 DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
773  boost::shared_ptr<MoFEM::FEMethod> method,
774  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
775  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
776  return DMMoFEMTSSetIJacobian<const std::string,
777  boost::shared_ptr<MoFEM::FEMethod>,
778  boost::shared_ptr<MoFEM::BasicMethod>,
779  boost::shared_ptr<MoFEM::BasicMethod>>(
780  dm, fe_name, method, pre_only, post_only);
781 }
782 
783 template <class S, class T0, class T1, class T2>
784 static PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, S fe_name, T0 method,
785  T1 pre_only, T2 post_only) {
786  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
788  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
789  if (pre_only)
790  dm_field->tsCtx->get_preProcess_to_do_RHSFunction().push_back(pre_only);
791  if (method)
792  dm_field->tsCtx->get_loops_to_do_RHSFunction().push_back(
793  PairNameFEMethodPtr(fe_name, method));
794  if (post_only)
795  dm_field->tsCtx->get_postProcess_to_do_RHSFunction().push_back(post_only);
796  CHKERR DMTSSetRHSFunction(dm, TsSetRHSFunction, dm_field->tsCtx.get());
798 }
799 
800 PetscErrorCode
801 DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name,
802  boost::shared_ptr<MoFEM::FEMethod> method,
803  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
804  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
805  return DMMoFEMTSSetRHSFunction<const std::string,
806  boost::shared_ptr<MoFEM::FEMethod>,
807  boost::shared_ptr<MoFEM::BasicMethod>,
808  boost::shared_ptr<MoFEM::BasicMethod>>(
809  dm, fe_name, method, pre_only, post_only);
811 }
812 
813 template <class S, class T0, class T1, class T2>
814 static PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, S fe_name, T0 method,
815  T1 pre_only, T2 post_only) {
816  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
818  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
819  if (pre_only)
820  dm_field->tsCtx->get_preProcess_to_do_RHSFunction().push_back(pre_only);
821  if (method)
822  dm_field->tsCtx->get_loops_to_do_RHSFunction().push_back(
823  PairNameFEMethodPtr(fe_name, method));
824  if (post_only)
825  dm_field->tsCtx->get_postProcess_to_do_RHSFunction().push_back(post_only);
826  CHKERR DMTSSetRHSJacobian(dm, TsSetRHSJacobian, dm_field->tsCtx.get());
828 }
829 
830 PetscErrorCode
831 DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
832  boost::shared_ptr<MoFEM::FEMethod> method,
833  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
834  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
835  return DMMoFEMTSSetRHSJacobian<const std::string,
836  boost::shared_ptr<MoFEM::FEMethod>,
837  boost::shared_ptr<MoFEM::BasicMethod>,
838  boost::shared_ptr<MoFEM::BasicMethod>>(
839  dm, fe_name, method, pre_only, post_only);
841 }
842 
843 PetscErrorCode
844 DMMoFEMTSSetI2Function(DM dm, const std::string fe_name,
845  boost::shared_ptr<MoFEM::FEMethod> method,
846  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
847  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
848  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
850  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
851  if (pre_only) {
852  dm_field->tsCtx->get_preProcess_to_do_IFunction().push_back(pre_only);
853  }
854  if (method) {
855  dm_field->tsCtx->get_loops_to_do_IFunction().push_back(
856  PairNameFEMethodPtr(fe_name, method));
857  }
858  if (post_only) {
859  dm_field->tsCtx->get_postProcess_to_do_IFunction().push_back(post_only);
860  }
861  CHKERR DMTSSetI2Function(dm, TsSetI2Function, dm_field->tsCtx.get());
863 }
864 
865 PetscErrorCode
866 DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name,
867  boost::shared_ptr<MoFEM::FEMethod> method,
868  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
869  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
870  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
872  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
873  if (pre_only) {
874  dm_field->tsCtx->get_preProcess_to_do_IJacobian().push_back(pre_only);
875  }
876  if (method) {
877  dm_field->tsCtx->get_loops_to_do_IJacobian().push_back(
878  PairNameFEMethodPtr(fe_name, method));
879  }
880  if (post_only) {
881  dm_field->tsCtx->get_postProcess_to_do_IJacobian().push_back(post_only);
882  }
883  CHKERR DMTSSetI2Jacobian(dm, TsSetI2Jacobian, dm_field->tsCtx.get());
885 }
886 
887 template <class S, class T0, class T1, class T2>
888 static PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, S fe_name, T0 method,
889  T1 pre_only, T2 post_only) {
890  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
892  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
893  if (pre_only)
894  dm_field->tsCtx->get_preProcess_to_do_Monitor().push_back(pre_only);
895  if (method)
896  dm_field->tsCtx->get_loops_to_do_Monitor().push_back(
897  PairNameFEMethodPtr(fe_name, method));
898  if (post_only)
899  dm_field->tsCtx->get_postProcess_to_do_Monitor().push_back(post_only);
900  CHKERR TSMonitorSet(ts, TsMonitorSet, dm_field->tsCtx.get(), PETSC_NULL);
902 }
903 
904 PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const char fe_name[],
905  MoFEM::FEMethod *method,
906  MoFEM::BasicMethod *pre_only,
907  MoFEM::BasicMethod *post_only) {
908  return DMMoFEMTSSetMonitor<const char *, MoFEM::FEMethod *,
910  dm, ts, fe_name, method, pre_only, post_only);
912 }
913 
914 PetscErrorCode
915 DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name,
916  boost::shared_ptr<MoFEM::FEMethod> method,
917  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
918  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
919  return DMMoFEMTSSetMonitor<const std::string,
920  boost::shared_ptr<MoFEM::FEMethod>,
921  boost::shared_ptr<MoFEM::BasicMethod>,
922  boost::shared_ptr<MoFEM::BasicMethod>>(
923  dm, ts, fe_name, method, pre_only, post_only);
925 }
926 
927 PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx) {
928  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
930  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
931  *ksp_ctx = dm_field->kspCtx.get();
933 }
934 
935 PetscErrorCode
936 DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
937  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
939  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
940  const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
942 }
943 
944 PetscErrorCode DMMoFEMSetKspCtx(DM dm,
945  boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
946  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
948  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
949  dm_field->kspCtx = ksp_ctx;
951 }
952 
953 PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx) {
954  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
956  DMCtx *dm_field = (DMCtx *)dm->data;
957  *snes_ctx = dm_field->snesCtx.get();
959 }
960 
961 PetscErrorCode
962 DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
963  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
965  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
966  const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
968 }
969 
970 PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
971  boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
972  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
974  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
975  dm_field->snesCtx = snes_ctx;
977 }
978 
979 /** get if read mesh is partitioned
980  * \ingroup dm
981  */
982 PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned) {
983  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
985  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
986  dm_field->isPartitioned = is_partitioned;
988 }
989 
990 /** get if read mesh is partitioned
991  * \ingroup dm
992  */
993 PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned) {
994  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
996  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
997  *is_partitioned = dm_field->isPartitioned;
999 }
1000 
1001 PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx) {
1002  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1004  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1005  *ts_ctx = dm_field->tsCtx.get();
1007 }
1008 
1009 PetscErrorCode DMMoFEMGetTsCtx(DM dm,
1010  const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
1011  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1013  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1014  const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
1016 }
1017 
1018 PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
1019  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1021  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1022  dm_field->tsCtx = ts_ctx;
1024 }
1025 
1026 PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g) {
1027  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1029  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1030  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1031  dm_field->problemName, COL, g);
1033 }
1034 
1035 PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj<Vec> &g_ptr) {
1036  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1038  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1039  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
1040  dm_field->problemName, COL, g_ptr);
1042 }
1043 
1044 PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l) {
1045  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1047  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1048  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
1049  dm_field->problemName, COL, l);
1051 }
1052 
1053 PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M) {
1054  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1056  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1057  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1059  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1060  M);
1061  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1063  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1064  M);
1065  } else {
1066  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1067  "Matrix type not implemented");
1068  }
1070 }
1071 
1072 PetscErrorCode DMCreateMatrix_MoFEM(DM dm, SmartPetscObj<Mat> &M) {
1073  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1075  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1076  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
1078  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
1079  M);
1080  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
1082  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
1083  M);
1084  } else {
1085  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1086  "Matrix type not implemented");
1087  }
1089 }
1090 
1091 #if PETSC_VERSION_GE(3, 7, 0)
1092 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
1093  DM dm) {
1094 #elif PETSC_VERSION_GE(3, 5, 3)
1095 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm) {
1096 #else
1097 PetscErrorCode DMSetFromOptions_MoFEM(DM dm) {
1098 #endif
1099 
1100  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1102  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1103 #if PETSC_VERSION_GE(3, 5, 3)
1104  ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
1105  CHKERRG(ierr);
1106 #else
1107  ierr = PetscOptionsHead("DMMoFEM Options");
1108  CHKERRG(ierr);
1109 #endif
1110  ierr = PetscOptionsBool("-dm_is_partitioned",
1111  "set if mesh is partitioned (works which native MOAB "
1112  "file format, i.e. h5m",
1113  "DMSetUp", dm_field->isPartitioned,
1114  &dm_field->isPartitioned, NULL);
1115  CHKERRG(ierr);
1117 }
1118 
1119 PetscErrorCode DMSetUp_MoFEM(DM dm) {
1120  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1121  ProblemsManager *prb_mng_ptr;
1123  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1124  CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
1125 
1126  if (dm_field->isCompDM) {
1127  // It is composite probelm
1128  CHKERR prb_mng_ptr->buildCompsedProblem(
1129  dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
1130  dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
1131  } else {
1132  if (dm_field->isPartitioned) {
1134  dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
1135  dm_field->verbosity);
1136  } else {
1137  CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
1138  dm_field->isSquareMatrix == PETSC_TRUE,
1139  dm_field->verbosity);
1140  CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
1141  dm_field->verbosity);
1142  }
1143  }
1144 
1145  // Partition finite elements
1146  if (dm_field->isPartitioned) {
1147  CHKERR prb_mng_ptr->partitionFiniteElements(
1148  dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
1150  dm_field->problemName, dm_field->verbosity);
1151  } else {
1152  // partition finite elemnets
1153  CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1154  -1, -1, dm_field->verbosity);
1155  // Get ghost DOFs
1156  CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1157  dm_field->verbosity);
1158  }
1159 
1160  // Set flag that problem is build and partitioned
1161  dm_field->isProblemBuild = PETSC_TRUE;
1162 
1164 }
1165 
1166 PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm) {
1167  PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1168  ProblemsManager *prb_mng_ptr;
1170 
1171  DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1172 
1173  // build sub dm problem
1174  CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1175 
1176  map<std::string, std::pair<EntityType, EntityType>> *entity_map_row = nullptr;
1177  map<std::string, std::pair<EntityType, EntityType>> *entity_map_col = nullptr;
1178 
1179  if (subdm_field->mapTypeRow)
1180  entity_map_row = subdm_field->mapTypeRow.get();
1181  if (subdm_field->mapTypeCol)
1182  entity_map_row = subdm_field->mapTypeCol.get();
1183 
1184  CHKERR prb_mng_ptr->buildSubProblem(
1185  subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1186  subdm_field->problemMainOfSubPtr->getName(),
1187  subdm_field->isSquareMatrix == PETSC_TRUE, entity_map_row, entity_map_col,
1188  subdm_field->verbosity);
1189 
1190  // partition problem
1191  subdm_field->isPartitioned = subdm_field->isPartitioned;
1192  if (subdm_field->isPartitioned) {
1193  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1194  0, subdm_field->sIze,
1195  subdm_field->verbosity);
1196  // set ghost nodes
1198  subdm_field->problemName, subdm_field->verbosity);
1199  } else {
1200  // partition finite elements
1201  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1202  -1, -1, subdm_field->verbosity);
1203  // set ghost nodes
1204  CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1205  subdm_field->verbosity);
1206  }
1207 
1208  subdm_field->isProblemBuild = PETSC_TRUE;
1209 
1211 }
1212 
1213 PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec g, InsertMode mode,
1214  Vec l) {
1215  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1217  CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1219 }
1220 
1221 PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec g, InsertMode mode, Vec l) {
1223  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1225 
1226  CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1227 
1228  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1229  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1230  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1231 
1232  double *array_loc, *array_glob;
1233  CHKERR VecGetArray(l, &array_loc);
1234  CHKERR VecGetArray(g, &array_glob);
1235  switch (mode) {
1236  case INSERT_VALUES:
1237  cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1238  break;
1239  case ADD_VALUES:
1240  cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1241  break;
1242  default:
1243  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1244  }
1245  CHKERR VecRestoreArray(l, &array_loc);
1246  CHKERR VecRestoreArray(g, &array_glob);
1248 }
1249 
1250 PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM dm, Vec l, InsertMode mode,
1251  Vec g) {
1252 
1253  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1255 
1256  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1257  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1258  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1259 
1260  double *array_loc, *array_glob;
1261  CHKERR VecGetArray(l, &array_loc);
1262  CHKERR VecGetArray(g, &array_glob);
1263  switch (mode) {
1264  case INSERT_VALUES:
1265  cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1266  break;
1267  case ADD_VALUES:
1268  cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1269  break;
1270  default:
1271  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1272  }
1273  CHKERR VecRestoreArray(l, &array_loc);
1274  CHKERR VecRestoreArray(g, &array_glob);
1275 
1277 }
1278 
1279 PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec l, InsertMode mode, Vec g) {
1280  //
1283 }
1284 
1285 PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
1286  char ***fieldNames, IS **fields) {
1287  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1289 
1290  if (numFields) {
1291  *numFields = 0;
1292  }
1293  if (fieldNames) {
1294  *fieldNames = NULL;
1295  }
1296  if (fields) {
1297  *fields = NULL;
1298  }
1299 
1300  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1301  auto fields_ptr = dm_field->mField_ptr->get_fields();
1302  Field_multiIndex::iterator fit, hi_fit;
1303  fit = fields_ptr->begin();
1304  hi_fit = fields_ptr->end();
1305  *numFields = std::distance(fit, hi_fit);
1306 
1307  if (fieldNames) {
1308  CHKERR PetscMalloc1(*numFields, fieldNames);
1309  }
1310  if (fields) {
1311  CHKERR PetscMalloc1(*numFields, fields);
1312  }
1313 
1314  for (int f = 0; fit != hi_fit; fit++, f++) {
1315  if (fieldNames) {
1316  CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1317  (char **)&(*fieldNames)[f]);
1318  }
1319  if (fields) {
1320  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1321  ->isCreateProblemFieldAndRank(
1322  dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1323  fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1324  }
1325  }
1326 
1328 }
1329 
1330 PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
1331  IS *is) {
1332  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1334  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1335  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1336  ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1337  field_name, 0, 1000, is);
1339 }
1340 
1341 PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb) {
1342  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1344  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1345  dm_field->verbosity = verb;
1347 }
1348 
1349 } // namespace MoFEM
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
MoFEM::DMCtx::query_interface
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: DMMMoFEM.cpp:36
MoFEM::TsCtx
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
cblas_daxpy
void cblas_daxpy(const int N, const double alpha, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_daxpy.c:11
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:26
MoFEM::DMSetFromOptions_MoFEM
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1097
MoFEM::CoreInterface::problem_basic_method_postProcess
virtual MoFEMErrorCode problem_basic_method_postProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
MoFEM::KspCtx
Interface for linear (KSP) solver.
Definition: KspCtx.hpp:27
MoFEM::CoreInterface::loop_dofs
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.
MoFEM::ProblemsManager::buildProblem
MoFEMErrorCode buildProblem(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
Definition: ProblemsManager.cpp:428
MoFEM::DMMoFEMSetSnesCtx
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > &snes_ctx)
Set MoFEM::SnesCtx data structure.
Definition: DMMMoFEM.cpp:970
MoFEM::DMMoFEMKSPSetComputeRHS
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:556
MoFEM::BasicMethod
Data structure to exchange data between mofem and User Loop Methods.
Definition: LoopMethods.hpp:225
MoFEM::VecManager
Vector manager is used to create vectors \mofem_vectors.
Definition: VecManager.hpp:36
MoFEM::Problem::getBitRefLevel
BitRefLevel getBitRefLevel() const
Definition: ProblemsMultiIndices.hpp:538
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
MoFEM::IDD_DMCTX
static const MOFEMuuid IDD_DMCTX
Definition: DMMoFEM.hpp:862
MoFEM::DMCtx::isCompDM
PetscBool isCompDM
Definition: DMMoFEM.hpp:901
MoFEM::DMMoFEMResolveSharedFiniteElements
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:387
MoFEM::CommInterface
Managing BitRefLevels.
Definition: CommInterface.hpp:33
MoFEM::SnesMat
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
MoFEM::DMCtx::tsCtx
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition: DMMoFEM.hpp:919
MoFEM::DMMoFEMSetSquareProblem
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition: DMMMoFEM.cpp:378
MoFEM::DMMoFEMTSSetRHSJacobian
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:831
MoFEM::DMMoFEMAddSubFieldCol
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:221
MoFEM::DMCtx::rowFields
std::vector< std::string > rowFields
Definition: DMMoFEM.hpp:897
MoFEM::DMMoFEMAddElement
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:425
MoFEM::FEMethod
structure for User Loop Methods on finite elements
Definition: LoopMethods.hpp:369
MoFEM::TsSetIJacobian
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
MoFEM::DMGlobalToLocalBegin_MoFEM
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1213
MoFEM::DMMoFEMAddColCompositeProblem
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on col.
Definition: DMMMoFEM.cpp:307
MoFEM::DMCtx::isProblemBuild
PetscBool isProblemBuild
True if problem is build.
Definition: DMMoFEM.hpp:883
MoFEM::Problem::getName
std::string getName() const
Definition: ProblemsMultiIndices.hpp:527
MoFEM::CoreInterface::get_problem
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
MoFEM::DMCtx::kspCtx
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition: DMMoFEM.hpp:917
MoFEM::DMSetUp_MoFEM
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1119
MoFEM::DMoFEMMeshToLocalVector
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:445
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition: MatrixFunction.hpp:70
MoFEM::DMCtx::mapTypeCol
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeCol
Definition: DMMoFEM.hpp:907
MoFEM::DMMoFEMGetIsSubDM
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMMoFEM.cpp:244
MoFEM::CoreInterface::modify_problem_add_finite_element
virtual MoFEMErrorCode modify_problem_add_finite_element(const std::string &name_problem, const std::string &fe_name)=0
add finite element to problem, this add entities assigned to finite element to a particular problem
MoFEM::DMMoFEMKSPSetComputeOperators
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:597
MoFEM::DMLocalToGlobalBegin_MoFEM
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1250
MoFEM::TsSetIFunction
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:17
MoFEM::CoreInterface::check_problem
virtual bool check_problem(const std::string name)=0
check if problem exist
ROW
@ ROW
Definition: definitions.h:192
MoFEM::DMMoFEMGetIsCompDM
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition: DMMMoFEM.cpp:327
MoFEM::DMCreateGlobalVector_MoFEM
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1026
MoFEM::DMCtx::snesCtx
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition: DMMoFEM.hpp:918
MoFEM::DMCtx::colFields
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:898
VERBOSE
@ VERBOSE
Definition: definitions.h:278
MoFEM::DMCreateMatrix_MoFEM
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1053
MoFEM::DMMoFEMAddRowCompositeProblem
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on row.
Definition: DMMMoFEM.cpp:289
MoFEM::DMMoFEMGetSubColIS
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:271
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:604
MoFEM::TsSetI2Jacobian
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
MoFEM::DMMoFEMAddSubFieldRow
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:198
MoFEM::DMMoFEMResolveSharedEntities
DEPRECATED PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[])
Definition: DMMMoFEM.cpp:397
MoFEM::DMMoFEMSetTsCtx
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > &ts_ctx)
Set MoFEM::TsCtx data structure.
Definition: DMMMoFEM.cpp:1018
MoFEM::DMMoFEMTSSetIJacobian
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:772
MoFEM::Problem::getSubData
boost::shared_ptr< SubProblemData > & getSubData() const
Get main problem of sub-problem is.
Definition: ProblemsMultiIndices.hpp:234
MoFEM::DMoFEMLoopFiniteElements
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:507
HenkyOps::f
auto f
Definition: HenkyOps.hpp:19
MoFEM::DMCtx::isSubDM
PetscBool isSubDM
Definition: DMMoFEM.hpp:896
MoFEM::DMCtx::isSquareMatrix
PetscBool isSquareMatrix
true if rows equals to cols
Definition: DMMoFEM.hpp:888
MoFEM::DMMoFEMTSSetI2Jacobian
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:866
MoFEM::DMDestroy_MoFEM
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:84
MoFEM::DMMoFEMGetProblemFiniteElementLayout
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[], PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMMoFEM.cpp:401
MoFEM::ProblemsManager::partitionGhostDofsOnDistributedMesh
MoFEMErrorCode partitionGhostDofsOnDistributedMesh(const std::string name, int verb=VERBOSE)
determine ghost nodes on distributed meshes
Definition: ProblemsManager.cpp:2733
MoFEM::DMSubDMSetUp_MoFEM
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1166
MoFEM::TsSetRHSFunction
PetscErrorCode TsSetRHSFunction(TS ts, PetscReal t, Vec u, Vec F, void *ctx)
TS solver function.
Definition: TsCtx.cpp:267
MoFEM::DMMoFEMTSSetI2Function
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:844
MoFEM::DMMoFEMSetVerbosity
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition: DMMMoFEM.cpp:1341
MoFEM::DMMoFEMCreateSubDM
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:177
RowColData
RowColData
RowColData.
Definition: definitions.h:192
M
static Index< 'M', 3 > M
Definition: BasicFeTools.hpp:86
MoFEM::DMMoFEMSetDestroyProblem
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:360
MoFEM::SnesRhs
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
MoFEM::DMCtx::rowCompPrb
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:902
MoFEM::DMCtx::isPartitioned
PetscBool isPartitioned
true if read mesh is on parts
Definition: DMMoFEM.hpp:887
MoFEM::DMoFEMMeshToGlobalVector
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:457
COL
@ COL
Definition: definitions.h:192
MoFEM::DMMoFEMGetSnesCtx
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMMoFEM.cpp:953
MoFEM::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:48
MoFEM::CoreInterface::problem_basic_method_preProcess
virtual MoFEMErrorCode problem_basic_method_preProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
MoFEM::Problem::getNumberOfElementsByNameAndPart
MoFEMErrorCode getNumberOfElementsByNameAndPart(MPI_Comm comm, const std::string name, PetscLayout *layout) const
Get number of finite elements by name on processors.
Definition: ProblemsMultiIndices.cpp:119
MoFEM::MOFEMuuid
MoFEM interface unique ID.
Definition: UnknownInterface.hpp:26
MoFEM::SnesCtx
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:27
MoFEM::DofMethod
Data structure to exchange data between mofem and User Loop Methods on entities.
Definition: LoopMethods.hpp:554
MoFEM::SmartPetscObj< Vec >
MoFEM::ProblemsManager
Problem manager is used to build and partition problems.
Definition: ProblemsManager.hpp:34
MoFEM::DMCtx::DMCtx
DMCtx()
Definition: DMMMoFEM.cpp:30
MoFEM::ProblemsManager::buildCompsedProblem
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
Definition: ProblemsManager.cpp:1521
MoFEM::DMMoFEMSNESSetJacobian
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:678
MoFEM::DMCtx::problemPtr
const Problem * problemPtr
pinter to problem data structure
Definition: DMMoFEM.hpp:893
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MoFEM::DMoFEMPreProcessFiniteElements
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:468
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
MoFEM::DMMoFEMUnSetElement
PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[])
unset element from dm
Definition: DMMMoFEM.cpp:435
MoFEM::DMMoFEMGetSquareProblem
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problem
Definition: DMMMoFEM.cpp:415
MoFEM::DMGlobalToLocalEnd_MoFEM
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1221
MoFEM::DMoFEMPostProcessFiniteElements
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:478
MoFEM::DMMoFEMCreateMoFEM
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:105
MoFEM::ISManager
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:36
DMMoFEM.hpp
Discrete manager interface for MoFEM.
MoFEM::ProblemsManager::partitionFiniteElements
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
Definition: ProblemsManager.cpp:2416
MoFEM::ProblemsManager::buildSubProblem
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
Definition: ProblemsManager.cpp:1265
MoFEM::DMCtx::destroyProblem
PetscBool destroyProblem
If true destroy problem with DM.
Definition: DMMoFEM.hpp:909
MoFEM::DMSetOperators_MoFEM
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:54
MoFEM::CoreInterface::modify_problem_ref_level_add_bit
virtual MoFEMErrorCode modify_problem_ref_level_add_bit(const std::string &name_problem, const BitRefLevel &bit)=0
add ref level to problem
MoFEM::DMCreateFieldIS_MoFEM
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1285
MoFEM::Problem::getNbLocalDofsRow
DofIdx getNbLocalDofsRow() const
Definition: ProblemsMultiIndices.hpp:533
MoFEM::DMCtx::colCompPrb
std::vector< std::string > colCompPrb
Definition: DMMoFEM.hpp:903
MoFEM::DMCtx::rAnk
int rAnk
Definition: DMMoFEM.hpp:890
MoFEM::DMCtx
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:877
MoFEM::DMMoFEMGetKspCtx
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:927
MoFEM::CoreInterface::get_fields
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
MoFEM::Problem::getNbGhostDofsRow
DofIdx getNbGhostDofsRow() const
Definition: ProblemsMultiIndices.hpp:535
MoFEM::DeprecatedCoreInterface::loop_finite_elements
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)
Definition: DeprecatedCoreInterface.cpp:579
MoFEM::DMCreate_MoFEM
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:76
MoFEM::TsSetI2Function
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
MoFEM::DMMoFEMGetSubRowIS
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:253
MoFEM::DMMoFEMGetFieldIS
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMMoFEM.cpp:1330
MoFEM::DMoFEMGetInterfacePtr
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:336
MoFEM::DMCtx::problemName
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:884
MoFEM::DMCtx::sIze
int sIze
Definition: DMMoFEM.hpp:890
MoFEM::DMMoFEMTSSetIFunction
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:719
MoFEM::DMCtx::mapTypeRow
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeRow
Definition: DMMoFEM.hpp:905
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:123
MoFEM::KspMat
PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx)
Run over elenents in the list.
Definition: KspCtx.cpp:87
MoFEM::ProblemsManager::buildProblemOnDistributedMesh
MoFEMErrorCode buildProblemOnDistributedMesh(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures, assuming that mesh is distributed (collective)
Definition: ProblemsManager.cpp:631
MoFEM::DMMoFEMGetIsPartitioned
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMMoFEM.cpp:993
MoFEM::DMMoFEMGetProblemPtr
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:348
MoFEM::DMMoFEMGetDestroyProblem
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMMoFEM.cpp:369
MoFEM::KspRhs
PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx)
Run over elements in the lists.
Definition: KspCtx.cpp:23
MoFEM::DMCtx::verbosity
int verbosity
verbosity
Definition: DMMoFEM.hpp:914
MoFEM::DMMoFEMGetTsCtx
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMMoFEM.cpp:1001
MoFEM::DMoFEMLoopFiniteElementsUpAndLowRank
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:488
MoFEM::DMCtx::problemMainOfSubPtr
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
Definition: DMMoFEM.hpp:899
MoFEM::CoreInterface::modify_problem_unset_finite_element
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...
MoFEM::TsMonitorSet
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:209
MoFEM::PairNameFEMethodPtr
Definition: AuxPETSc.hpp:24
MoFEM::DMMoFEMTSSetMonitor
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:915
MoFEM::TsSetRHSJacobian
PetscErrorCode TsSetRHSJacobian(TS ts, PetscReal t, Vec u, Mat A, Mat B, void *ctx)
TS solver function.
Definition: TsCtx.cpp:365
MoFEM::CoreInterface::modify_problem_mask_ref_level_add_bit
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
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
MoFEM::Problem
keeps basic data about problem
Definition: ProblemsMultiIndices.hpp:79
MoFEM::MatrixManager
Matrix manager is used to build and partition problems.
Definition: MatrixManager.hpp:34
MoFEM::ProblemsManager::partitionGhostDofs
MoFEMErrorCode partitionGhostDofs(const std::string name, int verb=VERBOSE)
determine ghost nodes
Definition: ProblemsManager.cpp:2599
MF_EXCL
@ MF_EXCL
Definition: definitions.h:189
MoFEM::DMLocalToGlobalEnd_MoFEM
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1279
MoFEM::CoreInterface::add_problem
virtual MoFEMErrorCode add_problem(const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add problem.
MoFEM::DMMoFEMDuplicateDMCtx
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
Definition: DMMMoFEM.cpp:160
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
MoFEM::DMCtx::referenceNumber
int referenceNumber
Definition: DMMoFEM.hpp:915
MoFEM::DMoFEMLoopDofs
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMMoFEM.cpp:524
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:124
MoFEM::CoreInterface::delete_problem
virtual MoFEMErrorCode delete_problem(const std::string name)=0
Delete problem.
MoFEM::DMCreateLocalVector_MoFEM
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition: DMMMoFEM.cpp:1044
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
CHKERRG
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:552
MoFEM::DMCtx::mField_ptr
Interface * mField_ptr
MoFEM interface.
Definition: DMMoFEM.hpp:882
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
Definition: UnknownInterface.hpp:130
cblas_dcopy
void cblas_dcopy(const int N, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_dcopy.c:11
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MOFEM_INVALID_DATA
@ MOFEM_INVALID_DATA
Definition: definitions.h:128
MoFEM::DMMoFEMSetKspCtx
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > &ksp_ctx)
set MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:944
MoFEM::UnknownInterface
base class for all interface classes
Definition: UnknownInterface.hpp:84
MoFEM::ProblemsManager::partitionProblem
MoFEMErrorCode partitionProblem(const std::string name, int verb=VERBOSE)
partition problem dofs (collective)
Definition: ProblemsManager.cpp:1954
MoFEM::DMMoFEMTSSetRHSFunction
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:801
MoFEM::DMMoFEMSetIsPartitioned
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:982
MoFEM::DMMoFEMSNESSetFunction
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:637