v0.8.16
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 #include <Includes.hpp>
16 #include <version.h>
17 #include <definitions.h>
18 #include <Common.hpp>
19 
20 #include <h1_hdiv_hcurl_l2.h>
21 
22 #include <UnknownInterface.hpp>
23 
24 #include <MaterialBlocks.hpp>
25 #include <BCData.hpp>
26 #include <TagMultiIndices.hpp>
27 #include <CoordSysMultiIndices.hpp>
28 #include <FieldMultiIndices.hpp>
29 #include <EntsMultiIndices.hpp>
30 #include <DofsMultiIndices.hpp>
31 #include <FEMultiIndices.hpp>
32 #include <ProblemsMultiIndices.hpp>
34 #include <BCMultiIndices.hpp>
35 #include <CoreDataStructures.hpp>
36 #include <SeriesMultiIndices.hpp>
37 
38 #include <LoopMethods.hpp>
39 #include <Interface.hpp>
40 #include <MeshRefinement.hpp>
41 #include <PrismInterface.hpp>
42 #include <SeriesRecorder.hpp>
43 #include <ProblemsManager.hpp>
44 #include <ISManager.hpp>
45 #include <VecManager.hpp>
46 #include <Core.hpp>
47 
48 #include <AuxPETSc.hpp>
49 #include <KspCtx.hpp>
50 #include <SnesCtx.hpp>
51 #include <TsCtx.hpp>
52 
53 // #undef PETSC_VERSION_RELEASE
54 // #define PETSC_VERSION_RELEASE 1
55 
56 #if PETSC_VERSION_GE(3,6,0)
57  #include <petsc/private/dmimpl.h> /*I "petscdm.h" I*/
58  // #include <petsc/private/vecimpl.h> /*I "petscdm.h" I*/
59 #else
60  #include <petsc-private/dmimpl.h> /*I "petscdm.h" I*/
61  #include <petsc-private/vecimpl.h> /*I "petscdm.h" I*/
62 #endif
63 
64 #include <DMMoFEM.hpp>
65 
66 namespace MoFEM {
67 
69  : mField_ptr(PETSC_NULL), isProblemBuild(PETSC_FALSE),
70  isPartitioned(PETSC_FALSE), isSquareMatrix(PETSC_TRUE),
71  isSubDM(PETSC_FALSE), isCompDM(PETSC_FALSE), destroyProblem(PETSC_FALSE),
72  verbosity(VERBOSE), referenceNumber(0) {}
73 
75  // cerr << "Snes " << snesCtx.use_count() << endl;
76  // cerr << "Destroy DMCtx" << endl;
77 }
78 
80  UnknownInterface **iface) const {
82  *iface = NULL;
83  if (uuid == IDD_DMCTX) {
84  *iface = const_cast<DMCtx *>(this);
86  }
87  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
89 }
90 
91 PetscErrorCode DMRegister_MoFEM(const char sname[]) {
93  CHKERR DMRegister(sname, DMCreate_MoFEM);
95 }
96 
97 PetscErrorCode DMSetOperators_MoFEM(DM dm) {
98  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
100 
101  dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
102  dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
103  dm->ops->creatematrix = DMCreateMatrix_MoFEM;
104  dm->ops->setup = DMSetUp_MoFEM;
105  dm->ops->destroy = DMDestroy_MoFEM;
106  dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
107  dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
108  dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
109  dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
110  dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
111  dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
112 
113  // Default matrix type
114  CHKERR DMSetMatType(dm, MATMPIAIJ);
115 
117 }
118 
119 PetscErrorCode DMCreate_MoFEM(DM dm) {
120  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
122  dm->data = new DMCtx();
125 }
126 
127 PetscErrorCode DMDestroy_MoFEM(DM dm) {
128  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
129  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
131  if (((DMCtx *)dm->data)->referenceNumber == 0) {
132  if (dm_field->destroyProblem) {
133  if (dm_field->mField_ptr->check_problem(dm_field->problemName)) {
134  dm_field->mField_ptr->delete_problem(dm_field->problemName);
135  } // else problem has to be deleted by the user
136  }
137  // cerr << "Destroy " << dm_field->problemName << endl;
138  delete ((DMCtx *)dm->data);
139  } else {
140  // cerr << "Dereference " << dm_field->problemName << " " <<
141  // ((DMCtx*)dm->data)->referenceNumber << endl;
142  (((DMCtx *)dm->data)->referenceNumber)--;
143  }
145 }
146 
147 PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr,
148  const char problem_name[],
149  const MoFEM::BitRefLevel bit_level,
150  const MoFEM::BitRefLevel bit_mask) {
152 
153  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
154  if (!dm->data) {
155  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
156  "data structure for MoFEM not yet created");
157  }
158  if (!m_field_ptr) {
159  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
160  "DM function not implemented into MoFEM");
161  }
162  dm_field->mField_ptr = m_field_ptr;
163  dm_field->problemName = problem_name;
164  if (!m_field_ptr->check_problem(dm_field->problemName)) {
165  // problem is not defined, declare problem internally set bool to
166  // destroyProblem problem with DM
167  dm_field->destroyProblem = PETSC_TRUE;
168  CHKERR dm_field->mField_ptr->add_problem(dm_field->problemName, MF_EXCL,
169  dm_field->verbosity);
170  } else {
171  dm_field->destroyProblem = PETSC_FALSE;
172  }
174  dm_field->problemName, bit_level);
176  dm_field->problemName, bit_mask);
177  dm_field->kspCtx =
178  boost::shared_ptr<KspCtx>(new KspCtx(*m_field_ptr, problem_name));
179  dm_field->snesCtx =
180  boost::shared_ptr<SnesCtx>(new SnesCtx(*m_field_ptr, problem_name));
181  dm_field->tsCtx =
182  boost::shared_ptr<TsCtx>(new TsCtx(*m_field_ptr, problem_name));
183 
184  MPI_Comm comm;
185  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
186  int result = 0;
187  MPI_Comm_compare(comm, m_field_ptr->get_comm(), &result);
188  // std::cerr << result << " " << MPI_IDENT << " " << MPI_CONGRUENT << " " <<
189  // MPI_SIMILAR << " " << MPI_UNEQUAL << std::endl;
190  if (result > MPI_CONGRUENT) {
191  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
192  "MoFEM and DM using different communicators");
193  }
194  MPI_Comm_size(comm, &dm_field->sIze);
195  MPI_Comm_rank(comm, &dm_field->rAnk);
196 
197  // problem structure
198  CHKERR dm_field->mField_ptr->get_problem(dm_field->problemName,
199  &dm_field->problemPtr);
200 
202 }
203 
204 PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]) {
206 
207  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
208  if (!dm->data) {
209  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
210  "data structure for MoFEM not yet created");
211  }
212  CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
213  dm_field->problemPtr->getBitRefLevel());
214 
215  DMCtx *subdm_field = (DMCtx *)subdm->data;
216  subdm_field->isSubDM = PETSC_TRUE;
217  subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
218  subdm_field->isPartitioned = dm_field->isPartitioned;
219  subdm_field->isSquareMatrix = PETSC_FALSE;
220  subdm->ops->setup = DMSubDMSetUp_MoFEM;
221 
223 }
224 
225 PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[]) {
226  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
228  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
229  if (!dm->data) {
230  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
231  "data structure for MoFEM not yet created");
232  }
233  if (!dm_field->isSubDM) {
234  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
235  }
236  dm_field->rowFields.push_back(field_name);
238 }
239 
240 PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[]) {
241  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
243  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
244  if (!dm->data) {
245  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
246  "data structure for MoFEM not yet created");
247  }
248  if (!dm_field->isSubDM) {
249  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
250  }
251  dm_field->colFields.push_back(field_name);
253 }
254 
255 PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm) {
257  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
259  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
260  *is_sub_dm = dm_field->isSubDM;
262 }
263 
264 PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is) {
266  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
268  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
269  if (dm_field->isSubDM != PETSC_TRUE) {
270  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
271  "This DM is not created as a SubDM");
272  }
273  if (dm_field->isProblemBuild != PETSC_TRUE) {
274  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
275  }
276  boost::shared_ptr<Problem::SubProblemData> sub_data =
277  dm_field->problemPtr->getSubData();
278  CHKERR sub_data->getRowIs(is);
280 }
281 
282 PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is) {
284  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
286  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
287  if (dm_field->isSubDM != PETSC_TRUE) {
288  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
289  "This DM is not created as a SubDM");
290  }
291  if (dm_field->isProblemBuild != PETSC_TRUE) {
292  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
293  }
294  boost::shared_ptr<Problem::SubProblemData> sub_data =
295  dm_field->problemPtr->getSubData();
296  CHKERR sub_data->getColIs(is);
298 }
299 
300 PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]) {
301  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
303  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
304  if (!dm->data) {
305  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
306  "data structure for MoFEM not yet created");
307  }
308  if (!dm_field->isCompDM) {
309  dm_field->isCompDM = PETSC_TRUE;
310  }
311  dm_field->rowCompPrb.push_back(prb_name);
312  if (dm_field->isSquareMatrix) {
313  dm_field->colCompPrb.push_back(prb_name);
314  }
316 }
317 
318 PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]) {
319  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
321  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
322  if (!dm->data) {
323  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
324  "data structure for MoFEM not yet created");
325  }
326  if (!dm_field->isCompDM) {
327  dm_field->isCompDM = PETSC_TRUE;
328  }
329  if (dm_field->isSquareMatrix) {
330  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
331  "No need to add problem on column when problem block structurally "
332  "symmetric");
333  }
334  dm_field->colCompPrb.push_back(prb_name);
336 }
337 
338 PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm) {
340  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
342  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
343  *is_comp_dm = dm_field->isCompDM;
345 }
346 
347 PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr) {
348  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
350  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
351  if (!dm->data) {
352  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
353  "data structure for MoFEM not yet created");
354  }
355  *m_field_ptr = dm_field->mField_ptr;
357 }
358 
359 PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr) {
360  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
362  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
363  if (!dm->data) {
364  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
365  "data structure for MoFEM not yet created");
366  }
367  *problem_ptr = dm_field->problemPtr;
369 }
370 
371 PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem) {
373  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
375  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
376  dm_field->destroyProblem = destroy_problem;
378 }
379 
380 PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem) {
382  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
384  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
385  *destroy_problem = dm_field->destroyProblem;
387 }
388 
389 PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem) {
391  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
393  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
394  dm_field->isSquareMatrix = square_problem;
396 }
397 
398 PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[]) {
400  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
402  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
403  CHKERR dm_field->mField_ptr->resolve_shared_ents(dm_field->problemPtr,
404  fe_name);
406 }
407 
408 PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[],
409  PetscLayout *layout) {
411  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
413  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
414 
415  MPI_Comm comm;
416  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
417  CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
418  layout);
420 }
421 
422 PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem) {
425  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
427  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
428  *square_problem = dm_field->isSquareMatrix;
430 }
431 
432 PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[]) {
433  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
435  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
437  dm_field->problemName, fe_name);
438  CHKERRG(ierr);
440 }
441 
442 PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[]) {
443  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
445  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
447  dm_field->problemName, fe_name);
448  CHKERRG(ierr);
450 }
451 
452 PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
453  ScatterMode scatter_mode) {
455  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
457  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
458  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
459  dm_field->problemPtr, COL, l, mode, scatter_mode);
460  CHKERRG(ierr);
462 }
463 
464 PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
465  ScatterMode scatter_mode) {
466  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
468  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
469  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
470  dm_field->problemPtr, COL, g, mode, scatter_mode);
471  CHKERRG(ierr);
473 }
474 
475 PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
476  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
478  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
480  dm_field->problemPtr, *method);
481  CHKERRG(ierr);
483 }
484 
485 PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
486  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
488  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
490  dm_field->problemPtr, *method);
491  CHKERRG(ierr);
493 }
494 
495 PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
496  MoFEM::FEMethod *method,
497  int low_rank, int up_rank) {
499  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
500  ierr = dm_field->mField_ptr->loop_finite_elements(
501  dm_field->problemPtr, fe_name, *method, low_rank, up_rank);
502  CHKERRG(ierr);
504 }
505 
506 PetscErrorCode
507 DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const std::string &fe_name,
508  boost::shared_ptr<MoFEM::FEMethod> method,
509  int low_rank, int up_rank) {
510  return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
511  low_rank, up_rank);
512 }
513 
514 PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
515  MoFEM::FEMethod *method) {
516  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
518  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
519  ierr = DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name, method,
520  dm_field->rAnk, dm_field->rAnk);
521  CHKERRG(ierr);
523 }
524 
525 PetscErrorCode
526 DMoFEMLoopFiniteElements(DM dm, const std::string &fe_name,
527  boost::shared_ptr<MoFEM::FEMethod> method) {
528  return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get());
529 }
530 
531 PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
532  MoFEM::DofMethod *method) {
533  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
535  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
536  ierr =
537  dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
538  *method, dm_field->rAnk, dm_field->rAnk);
539  CHKERRG(ierr);
541 }
542 
543 template <class S, class T0, class T1, class T2 >
544 static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method,
545  T1 pre_only, T2 post_only) {
546  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
548  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
549  if (pre_only) {
550  dm_field->kspCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
551  }
552  if (method) {
553  dm_field->kspCtx->get_loops_to_do_Rhs().push_back(
554  PairNameFEMethodPtr(fe_name, method));
555  }
556  if (post_only) {
557  dm_field->kspCtx->get_postProcess_to_do_Rhs().push_back(post_only);
558  }
559  CHKERR DMKSPSetComputeRHS(dm, KspRhs, dm_field->kspCtx.get());
561 }
562 
563 PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
564  MoFEM::FEMethod *method,
565  MoFEM::BasicMethod *pre_only,
566  MoFEM::BasicMethod *post_only) {
567  return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
569  dm, fe_name, method, pre_only, post_only);
570 }
571 
572 PetscErrorCode
573 DMMoFEMKSPSetComputeRHS(DM dm, const std::string &fe_name,
574  boost::shared_ptr<MoFEM::FEMethod> method,
575  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
576  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
577  return DMMoFEMKSPSetComputeRHS<const std::string &,
578  boost::shared_ptr<MoFEM::FEMethod>,
579  boost::shared_ptr<MoFEM::BasicMethod>,
580  boost::shared_ptr<MoFEM::BasicMethod>>(
581  dm, fe_name, method, pre_only, post_only);
582 }
583 
584 template <class S, class T0, class T1, class T2>
585 static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method,
586  T1 pre_only, T2 post_only) {
587  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
589  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
590  if (pre_only) {
591  dm_field->kspCtx->get_preProcess_to_do_Mat().push_back(pre_only);
592  }
593  if (method) {
594  dm_field->kspCtx->get_loops_to_do_Mat().push_back(
595  PairNameFEMethodPtr(fe_name, method));
596  }
597  if (post_only) {
598  dm_field->kspCtx->get_postProcess_to_do_Mat().push_back(post_only);
599  }
600  CHKERR DMKSPSetComputeOperators(dm, KspMat, dm_field->kspCtx.get());
602 }
603 
604 PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
605  MoFEM::FEMethod *method,
606  MoFEM::BasicMethod *pre_only,
607  MoFEM::BasicMethod *post_only) {
608  return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
611  dm, fe_name, method, pre_only, post_only);
612 }
613 
614 PetscErrorCode
615 DMMoFEMKSPSetComputeOperators(DM dm, const std::string &fe_name,
616  boost::shared_ptr<MoFEM::FEMethod> method,
617  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
618  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
619  return DMMoFEMKSPSetComputeOperators<const std::string &,
620  boost::shared_ptr<MoFEM::FEMethod> >(
621  dm, fe_name, method, pre_only, post_only);
622 }
623 
624 template <class S, class T0, class T1, class T2>
625 static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method,
626  T1 pre_only, T2 post_only) {
627  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
629  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
630  if (pre_only) {
631  dm_field->snesCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
632  }
633  if (method) {
634  dm_field->snesCtx->get_loops_to_do_Rhs().push_back(
635  PairNameFEMethodPtr(fe_name, method));
636  }
637  if (post_only) {
638  dm_field->snesCtx->get_postProcess_to_do_Rhs().push_back(post_only);
639  }
640  CHKERR DMSNESSetFunction(dm, SnesRhs, dm_field->snesCtx.get());
642 }
643 
644 PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
645  MoFEM::FEMethod *method,
646  MoFEM::BasicMethod *pre_only,
647  MoFEM::BasicMethod *post_only) {
648  return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
650  dm, fe_name, method, pre_only, post_only);
651 }
652 
653 PetscErrorCode
654 DMMoFEMSNESSetFunction(DM dm, const std::string &fe_name,
655  boost::shared_ptr<MoFEM::FEMethod> method,
656  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
657  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
658  return DMMoFEMSNESSetFunction<const std::string &,
659  boost::shared_ptr<MoFEM::FEMethod>,
660  boost::shared_ptr<MoFEM::BasicMethod>,
661  boost::shared_ptr<MoFEM::BasicMethod>>(
662  dm, fe_name, method, pre_only, post_only);
663 }
664 
665 template <class S, class T0, class T1, class T2>
666 static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method,
667  T1 pre_only, T2 post_only) {
668  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
670  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
671  if (pre_only) {
672  dm_field->snesCtx->get_preProcess_to_do_Mat().push_back(pre_only);
673  }
674  if (method) {
675  dm_field->snesCtx->get_loops_to_do_Mat().push_back(
676  PairNameFEMethodPtr(fe_name, method));
677  }
678  if (post_only) {
679  dm_field->snesCtx->get_postProcess_to_do_Mat().push_back(post_only);
680  }
681  CHKERR DMSNESSetJacobian(dm, SnesMat, dm_field->snesCtx.get());
683 }
684 
685 PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
686  MoFEM::FEMethod *method,
687  MoFEM::BasicMethod *pre_only,
688  MoFEM::BasicMethod *post_only) {
689  return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
691  dm, fe_name, method, pre_only, post_only);
692 }
693 
694 PetscErrorCode
695 DMMoFEMSNESSetJacobian(DM dm, const std::string &fe_name,
696  boost::shared_ptr<MoFEM::FEMethod> method,
697  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
698  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
699  return DMMoFEMSNESSetJacobian<const std::string &,
700  boost::shared_ptr<MoFEM::FEMethod>,
701  boost::shared_ptr<MoFEM::BasicMethod>,
702  boost::shared_ptr<MoFEM::BasicMethod>>(
703  dm, fe_name, method, pre_only, post_only);
704 }
705 
706 template <class S, class T0, class T1, class T2>
707 static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method,
708  T1 pre_only, T2 post_only) {
709  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
711  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
712  if (pre_only) {
713  dm_field->tsCtx->get_preProcess_to_do_IFunction().push_back(pre_only);
714  }
715  if (method) {
716  dm_field->tsCtx->get_loops_to_do_IFunction().push_back(
717  PairNameFEMethodPtr(fe_name, method));
718  }
719  if (post_only) {
720  dm_field->tsCtx->get_postProcess_to_do_IFunction().push_back(post_only);
721  }
722  CHKERR DMTSSetIFunction(dm, f_TSSetIFunction, dm_field->tsCtx.get());
724 }
725 
726 PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
727  MoFEM::FEMethod *method,
728  MoFEM::BasicMethod *pre_only,
729  MoFEM::BasicMethod *post_only) {
730  return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
732  dm, fe_name, method, pre_only, post_only);
734 }
735 
736 PetscErrorCode
737 DMMoFEMTSSetIFunction(DM dm, const std::string &fe_name,
738  boost::shared_ptr<MoFEM::FEMethod> method,
739  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
740  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
741  return DMMoFEMTSSetIFunction<const std::string,
742  boost::shared_ptr<MoFEM::FEMethod>,
743  boost::shared_ptr<MoFEM::BasicMethod>,
744  boost::shared_ptr<MoFEM::BasicMethod>>(
745  dm, fe_name, method, pre_only, post_only);
747 }
748 
749 template <class S, class T0, class T1, class T2>
750 static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method,
751  T1 pre_only, T2 post_only) {
752  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
754  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
755  if (pre_only) {
756  dm_field->tsCtx->get_preProcess_to_do_IJacobian().push_back(pre_only);
757  }
758  if (method) {
759  dm_field->tsCtx->get_loops_to_do_IJacobian().push_back(
760  PairNameFEMethodPtr(fe_name, method));
761  }
762  if (post_only) {
763  dm_field->tsCtx->get_postProcess_to_do_IJacobian().push_back(post_only);
764  }
765  CHKERR DMTSSetIJacobian(dm, f_TSSetIJacobian, dm_field->tsCtx.get());
767 }
768 
769 PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
770  MoFEM::FEMethod *method,
771  MoFEM::BasicMethod *pre_only,
772  MoFEM::BasicMethod *post_only) {
773  return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
774  MoFEM::BasicMethod *>(dm, fe_name, method,
775  pre_only, post_only);
776 }
777 
778 PetscErrorCode
779 DMMoFEMTSSetIJacobian(DM dm, const std::string &fe_name,
780  boost::shared_ptr<MoFEM::FEMethod> method,
781  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
782  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
783  return DMMoFEMTSSetIJacobian<const std::string &,
784  boost::shared_ptr<MoFEM::FEMethod>,
785  boost::shared_ptr<MoFEM::BasicMethod>,
786  boost::shared_ptr<MoFEM::BasicMethod>>(
787  dm, fe_name, method, pre_only, post_only);
788 }
789 
790 PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx) {
791  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
793  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
794  *ksp_ctx = dm_field->kspCtx.get();
796 }
797 
798 PetscErrorCode
799 DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
800  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
802  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
803  const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
805 }
806 
807 PetscErrorCode DMMoFEMSetKspCtx(DM dm,
808  boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
809  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
811  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
812  dm_field->kspCtx = ksp_ctx;
814 }
815 
816 PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx) {
817  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
819  DMCtx *dm_field = (DMCtx *)dm->data;
820  *snes_ctx = dm_field->snesCtx.get();
822 }
823 
824 PetscErrorCode
825 DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
826  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
828  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
829  const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
831 }
832 
833 PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
834  boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
835  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
837  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
838  dm_field->snesCtx = snes_ctx;
840 }
841 
842 /** get if read mesh is partitioned
843  * \ingroup dm
844  */
845 PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned) {
846  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
848  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
849  dm_field->isPartitioned = is_partitioned;
851 }
852 
853 /** get if read mesh is partitioned
854  * \ingroup dm
855  */
856 PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned) {
857  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
859  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
860  *is_partitioned = dm_field->isPartitioned;
862 }
863 
864 PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx) {
865  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
867  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
868  *ts_ctx = dm_field->tsCtx.get();
870 }
871 
872 PetscErrorCode DMMoFEMGetTsCtx(DM dm,
873  const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
874  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
876  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
877  const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
879 }
880 
881 PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
882  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
884  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
885  dm_field->tsCtx = ts_ctx;
887 }
888 
889 PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g) {
890  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
892  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
893  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
894  dm_field->problemName, COL, g);
896 }
897 
898 PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l) {
899  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
901  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
902  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
903  dm_field->problemName, COL, l);
905 }
906 
907 PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M) {
908  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
910  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
911  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
913  dm_field->problemName, M);
914  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
915  PetscInt *i;
916  PetscInt *j;
917  PetscScalar *v;
918 #if PETSC_VERSION_GE(3, 7, 0)
920  dm_field->problemName, M, &i, &j, &v);
921  CHKERR MatConvert(*M, MATAIJ, MAT_INPLACE_MATRIX, M);
922 #else
923  Mat N;
925  dm_field->problemName, &N, &i, &j, &v);
926  CHKERR MatConvert(N, MATAIJ, MAT_INITIAL_MATRIX, M);
927  CHKERR MatDestroy(&N);
928 #endif
929  CHKERR PetscFree(i);
930  CHKERR PetscFree(j);
931  CHKERR PetscFree(v);
932  } else {
933  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
934  "Matrix type not implemented");
935  }
937 }
938 
939 #if PETSC_VERSION_GE(3, 7, 0)
940 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
941  DM dm) {
942 #elif PETSC_VERSION_GE(3, 5, 3)
943 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm) {
944 #else
945 PetscErrorCode DMSetFromOptions_MoFEM(DM dm) {
946 #endif
947 
948  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
950  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
951 #if PETSC_VERSION_GE(3, 5, 3)
952  ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
953  CHKERRG(ierr);
954 #else
955  ierr = PetscOptionsHead("DMMoFEM Options");
956  CHKERRG(ierr);
957 #endif
958  ierr = PetscOptionsBool("-dm_is_partitioned",
959  "set if mesh is partitioned (works which native MOAB "
960  "file format, i.e. h5m",
961  "DMSetUp", dm_field->isPartitioned,
962  &dm_field->isPartitioned, NULL);
963  CHKERRG(ierr);
965 }
966 
967 PetscErrorCode DMSetUp_MoFEM(DM dm) {
968  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
969  ProblemsManager *prb_mng_ptr;
971  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
972  CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
973 
974  if (dm_field->isCompDM) {
975  // It is composite probelm
976  CHKERR prb_mng_ptr->buildCompsedProblem(
977  dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
978  dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
979  } else {
980  if (dm_field->isPartitioned) {
982  dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
983  dm_field->verbosity);
984  } else {
985  CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
986  dm_field->isSquareMatrix == PETSC_TRUE,
987  dm_field->verbosity);
988  CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
989  dm_field->verbosity);
990  }
991  }
992 
993  // Partition finite elements
994  if (dm_field->isPartitioned) {
995  CHKERR prb_mng_ptr->partitionFiniteElements(
996  dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
998  dm_field->problemName, dm_field->verbosity);
999  } else {
1000  // partition finite elemnets
1001  CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
1002  -1, -1, dm_field->verbosity);
1003  // Get ghost DOFs
1004  CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
1005  dm_field->verbosity);
1006  }
1007 
1008  // Set flag that problem is build and partitioned
1009  dm_field->isProblemBuild = PETSC_TRUE;
1010 
1012 }
1013 
1014 PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm) {
1015  PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1016  ProblemsManager *prb_mng_ptr;
1018 
1019  DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1020 
1021  // build sub dm problem
1022  CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1023  CHKERR prb_mng_ptr->buildSubProblem(
1024  subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1025  subdm_field->problemMainOfSubPtr->getName(),
1026  subdm_field->isSquareMatrix == PETSC_TRUE, subdm_field->verbosity);
1027 
1028  // partition problem
1029  subdm_field->isPartitioned = subdm_field->isPartitioned;
1030  if (subdm_field->isPartitioned) {
1031  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1032  0, subdm_field->sIze,
1033  subdm_field->verbosity);
1034  // set ghost nodes
1036  subdm_field->problemName, subdm_field->verbosity);
1037  } else {
1038  // partition finite elements
1039  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1040  -1, -1, subdm_field->verbosity);
1041  // set ghost nodes
1042  CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1043  subdm_field->verbosity);
1044  }
1045 
1046  subdm_field->isProblemBuild = PETSC_TRUE;
1047 
1049 }
1050 
1051 PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec g, InsertMode mode,
1052  Vec l) {
1053  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1055  CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1057 }
1058 
1059 PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec g, InsertMode mode, Vec l) {
1061  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1063 
1064  CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1065 
1066  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1067  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1068  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1069 
1070  double *array_loc, *array_glob;
1071  CHKERR VecGetArray(l, &array_loc);
1072  CHKERR VecGetArray(g, &array_glob);
1073  switch (mode) {
1074  case INSERT_VALUES:
1075  cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1076  break;
1077  case ADD_VALUES:
1078  cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1079  break;
1080  default:
1081  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1082  }
1083  CHKERR VecRestoreArray(l, &array_loc);
1084  CHKERR VecRestoreArray(g, &array_glob);
1086 }
1087 
1088 PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM dm, Vec l, InsertMode mode,
1089  Vec g) {
1090 
1091  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1093 
1094  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1095  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1096  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1097 
1098  double *array_loc, *array_glob;
1099  CHKERR VecGetArray(l, &array_loc);
1100  CHKERR VecGetArray(g, &array_glob);
1101  switch (mode) {
1102  case INSERT_VALUES:
1103  cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1104  break;
1105  case ADD_VALUES:
1106  cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1107  break;
1108  default:
1109  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1110  }
1111  CHKERR VecRestoreArray(l, &array_loc);
1112  CHKERR VecRestoreArray(g, &array_glob);
1113 
1115 }
1116 
1117 PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec l, InsertMode mode, Vec g) {
1118  //
1121 }
1122 
1123 PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
1124  char ***fieldNames, IS **fields) {
1125  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1127 
1128  if (numFields) {
1129  *numFields = 0;
1130  }
1131  if (fieldNames) {
1132  *fieldNames = NULL;
1133  }
1134  if (fields) {
1135  *fields = NULL;
1136  }
1137 
1138  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1139  const Field_multiIndex *fields_ptr;
1140  CHKERR dm_field->mField_ptr->get_fields(&fields_ptr);
1141  Field_multiIndex::iterator fit, hi_fit;
1142  fit = fields_ptr->begin();
1143  hi_fit = fields_ptr->end();
1144  *numFields = std::distance(fit, hi_fit);
1145 
1146  if (fieldNames) {
1147  CHKERR PetscMalloc1(*numFields, fieldNames);
1148  }
1149  if (fields) {
1150  CHKERR PetscMalloc1(*numFields, fields);
1151  }
1152 
1153  for (int f = 0; fit != hi_fit; fit++, f++) {
1154  if (fieldNames) {
1155  CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1156  (char **)&(*fieldNames)[f]);
1157  }
1158  if (fields) {
1159  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1160  ->isCreateProblemFieldAndRank(
1161  dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1162  fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1163  }
1164  }
1165 
1167 }
1168 
1169 PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
1170  IS *is) {
1171  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1173  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1174  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1175  ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1176  field_name, 0, 1000, is);
1178 }
1179 
1180 PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb) {
1181  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1183  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1184  dm_field->verbosity = verb;
1186 }
1187 
1188 } // namespace MoFEM
Interface for mesh refinement.
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:744
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:264
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices...
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:644
MoFEM interface unique ID.
Problem manager is used to build and partition problems .
PetscBool isSquareMatrix
true if rows equals to cols
Definition: DMMoFEM.hpp:748
virtual bool check_problem(const std::string name)=0
check if problem exist
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition: DMMMoFEM.cpp:1180
int verbosity
verbosity
Definition: DMMoFEM.hpp:770
MoFEMErrorCode buildProblem(const std::string &name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
MoFEMErrorCode partitionGhostDofs(const std::string &name, int verb=VERBOSE)
determine ghost nodes
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Common.hpp:60
PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set compute operator for KSP solver via sub-matrix and IS.
Definition: DMMMoFEM.cpp:563
Multi-index containers, data structures for problems and other low-level functions.
virtual MoFEMErrorCode delete_problem(const std::string name)=0
Delete problem.
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
Definition: DMMMoFEM.cpp:240
Myltindex containers, data structures and other low-level functions.
Multi-index contains, for mofem entities data structures and other low-level functions.
Data structures for Meshset/Blocsk with material data.
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1051
boost::shared_ptr< SubProblemData > & getSubData() const
Get main problem of sub-problem is.
MoFEMErrorCode partitionFiniteElements(const std::string &name, bool part_from_moab=false, int low_proc=-1, int hi_proc=-1, int verb=VERBOSE)
partition finite elementsFunction which partition finite elements based on dofs partitioning. In addition it sets information about local row and cols dofs at given element on partition.
Context for PETSc KSP, i.e. nonlinear solver.
MoFEMErrorCode partitionProblem(const std::string &name, int verb=VERBOSE)
partition problem dofs (collective)
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
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:282
PetscBool isSubDM
Definition: DMMoFEM.hpp:756
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:33
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
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:452
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problemIt if true is assumed that matrix has the same indexing on rows and columns...
Definition: DMMMoFEM.cpp:389
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMMoFEM.cpp:255
PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:398
PetscErrorCode f_TSSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Definition: TsCtx.cpp:119
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:432
void cblas_dcopy(const int N, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_dcopy.c:11
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:128
PetscErrorCode f_TSSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Definition: TsCtx.cpp:56
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 ...
Interface * mField_ptr
MoFEM interface.
Definition: DMMoFEM.hpp:742
base class for all interface classes
PetscBool isPartitioned
true if read mesh is on parts
Definition: DMMoFEM.hpp:747
virtual MoFEMErrorCode MatCreateSeqAIJWithArrays(const std::string &name, Mat *Aij, PetscInt **i, PetscInt **j, PetscScalar **v, int verb=DEFAULT_VERBOSITY)=0
create Mat (AIJ) for problem
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:967
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:526
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on colThis create block on col with DOFs from problem of given name...
Definition: DMMMoFEM.cpp:318
Basic structures and data.
multi_index_container< boost::shared_ptr< Field >, indexed_by< hashed_unique< tag< BitFieldId_mi_tag >, const_mem_fun< Field, const BitFieldId &, &Field::getId >, HashBit< BitFieldId >, EqBit< BitFieldId > >, ordered_unique< tag< Meshset_mi_tag >, member< Field, EntityHandle, &Field::meshSet > >, ordered_unique< tag< FieldName_mi_tag >, const_mem_fun< Field, boost::string_ref, &Field::getNameRef > >, ordered_non_unique< tag< BitFieldId_space_mi_tag >, const_mem_fun< Field, FieldSpace, &Field::getSpace > > > > Field_multiIndex
Field_multiIndex for Field.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: DMMMoFEM.cpp:79
DofIdx getNbLocalDofsRow() const
Tags for Multi-index containers.
const Problem * problemPtr
pinter to problem data structure
Definition: DMMoFEM.hpp:753
std::vector< std::string > colCompPrb
Definition: DMMoFEM.hpp:763
MoFEM interface.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
Data structure to exchange data between mofem and User Loop Methods on entities.It allows to exchange...
virtual ~DMCtx()
Definition: DMMMoFEM.cpp:74
keeps basic data about problemThis is low level structure with information about problem, what elements compose problem and what DOFs are on rows and columns.
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1088
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:907
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:604
DofIdx getNbGhostDofsRow() const
MoFEMErrorCode buildProblemOnDistributedMesh(const std::string &name, const bool square_matrix, int verb=VERBOSE)
build problem data structures, assuming that mesh is distributed (collective)Mesh is distributed...
RowColData
RowColData.
Definition: definitions.h:183
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:514
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
Discrete manager interface for MoFEM.
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition: DMMMoFEM.cpp:338
PetscBool destroyProblem
If true destroy problem with DM.
Definition: DMMoFEM.hpp:765
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:127
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM...
Definition: DMMMoFEM.cpp:889
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[], PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMMoFEM.cpp:408
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition: DMMoFEM.hpp:775
MoFEM interface.
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:204
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on rowThis create block on row with DOFs from problem of given name...
Definition: DMMMoFEM.cpp:300
PetscErrorCode KspRhs(KSP ksp, Vec f, void *ctx)
Run over elements in the lists.
Definition: KspCtx.cpp:54
std::vector< std::string > rowFields
Definition: DMMoFEM.hpp:757
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Common.hpp:147
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > &ts_ctx)
Set MoFEM::TsCtx data structureIt take over pointer, do not delete it, DM will destroy pointer when i...
Definition: DMMMoFEM.cpp:881
virtual MoFEMErrorCode loop_finite_elements(const std::string &problem_name, const std::string &fe_name, FEMethod &method, MoFEMTypes bh=MF_EXIST, int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
Myltindex containes, data structures for mofem adjacencies and other low-level functions.
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1014
Coordinate systems attached to DOFs.
MoFEMErrorCode partitionGhostDofsOnDistributedMesh(const std::string &name, int verb=VERBOSE)
determine ghost nodes on distributed meshes
virtual MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const =0
Get problem database (data structure)
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problemIt if true is assumed that matrix has the same indexing on rows and columns...
Definition: DMMMoFEM.cpp:422
useful compiler directives and definitions
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMMoFEM.cpp:816
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:464
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
Interface managing problemsManaging problems, build and partitioning.
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition: DMMoFEM.hpp:773
PetscBool isProblemBuild
True if problem is build.
Definition: DMMoFEM.hpp:743
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:27
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
Definition: DMMoFEM.hpp:759
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:91
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMMoFEM.cpp:380
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition: DMMoFEM.hpp:774
Auxuliary MoFEM-PETSc structures.
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:790
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...
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:685
virtual MoFEMErrorCode resolve_shared_ents(const Problem *problem_ptr, const std::string &fe_name, int verb=DEFAULT_VERBOSITY)=0
resolve shared entities for finite elements in the problem
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.
PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[])
unset element from dm
Definition: DMMMoFEM.cpp:442
PetscErrorCode KspMat(KSP ksp, Mat A, Mat B, void *ctx)
Run over elenents in the list.
Definition: KspCtx.cpp:90
virtual MoFEMErrorCode modify_problem_ref_level_add_bit(const std::string &name_problem, const BitRefLevel &bit)=0
add ref level to problem
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:475
Multi-index contains, data structures for mofem finite elements and other low-level functions...
Vector manager is used to create vectors .
Definition: VecManager.hpp:35
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1059
Multi-index containers, data boundary data structures and other low-level functions.
Interface managing ISManaging problems, build and partitioning.
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:485
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:147
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:119
virtual MoFEMErrorCode problem_basic_method_preProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethod.
MoFEM interface.
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:97
#define CHKERR
Inline error check.
Definition: definitions.h:578
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1123
Field data structure storing information about space, approximation base, coordinate systems...
virtual MoFEMErrorCode MatCreateMPIAIJWithArrays(const std::string &name, Mat *Aij, int verb=DEFAULT_VERBOSITY)=0
create Mat (MPIAIJ) for problem (collective)
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:758
Functions to approximate hierarchical spaces.
BitRefLevel getBitRefLevel() const
PetscBool isCompDM
Definition: DMMoFEM.hpp:761
MoFEM interface.
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:762
Myltindex containers, for mofem fields data structures and other low-level functions.
Interface for linear (KSP) solver.
Definition: KspCtx.hpp:27
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition: DMMMoFEM.cpp:864
virtual MoFEMErrorCode add_problem(const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add problem.
Core interface class for user interface.
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMMoFEM.cpp:1169
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:347
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:737
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:371
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:726
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
Definition: DMMMoFEM.cpp:225
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:359
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMMoFEM.cpp:856
static const MOFEMuuid IDD_DMCTX
Definition: DMMoFEM.hpp:722
virtual MPI_Comm & get_comm() const =0
std::string getName() const
virtual MoFEMErrorCode get_fields(const Field_multiIndex **fields_ptr) const =0
Get fields multi-index from database.
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMMoFEM.cpp:531
MoFEMErrorCode getNumberOfElementsByNameAndPart(MPI_Comm comm, const std::string name, PetscLayout *layout) const
Get number of finite elements by name on processors.
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:495
const int N
Definition: speed_test.cpp:3
MoFEM interface.
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:55
virtual MoFEMErrorCode problem_basic_method_postProcess(const Problem *problem_ptr, BasicMethod &method, int verb=DEFAULT_VERBOSITY)=0
Set data for BasicMethodThis function set data about problem, adjacencies and other multi-indices in ...
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > &ksp_ctx)
set MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:807
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:845
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:779
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...
Context for PETSc Time Stepping.
Interface managing vectorsManaging problems, build and partitioning.
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1117
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
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:945
Context for PETSc SNES, i.e. nonlinear solver.
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, int verb=VERBOSE)
build sub problem
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM...
Definition: DMMMoFEM.cpp:898
Data structure with Cubit native blocks/meshet with boundary conditions.
Multi-Index contains, data structures for mofem dofs and other low-level functions.
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > &snes_ctx)
Set MoFEM::SnesCtx data structure.
Definition: DMMMoFEM.cpp:833