v0.8.19
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 <MatrixManager.hpp>
45 #include <ISManager.hpp>
46 #include <VecManager.hpp>
47 #include <Core.hpp>
48 
49 #include <AuxPETSc.hpp>
50 #include <KspCtx.hpp>
51 #include <SnesCtx.hpp>
52 #include <TsCtx.hpp>
53 
54 // #undef PETSC_VERSION_RELEASE
55 // #define PETSC_VERSION_RELEASE 1
56 
57 #if PETSC_VERSION_GE(3, 6, 0)
58 #include <petsc/private/dmimpl.h> /*I "petscdm.h" I*/
59 // #include <petsc/private/vecimpl.h> /*I "petscdm.h" I*/
60 #else
61 #include <petsc-private/dmimpl.h> /*I "petscdm.h" I*/
62 #include <petsc-private/vecimpl.h> /*I "petscdm.h" I*/
63 #endif
64 
65 #include <DMMoFEM.hpp>
66 
67 namespace MoFEM {
68 
70  : mField_ptr(PETSC_NULL), isProblemBuild(PETSC_FALSE),
71  isPartitioned(PETSC_FALSE), isSquareMatrix(PETSC_TRUE),
72  isSubDM(PETSC_FALSE), isCompDM(PETSC_FALSE), destroyProblem(PETSC_FALSE),
73  verbosity(VERBOSE), referenceNumber(0) {}
74 
76  // cerr << "Snes " << snesCtx.use_count() << endl;
77  // cerr << "Destroy DMCtx" << endl;
78 }
79 
81  UnknownInterface **iface) const {
83  *iface = NULL;
84  if (uuid == IDD_DMCTX) {
85  *iface = const_cast<DMCtx *>(this);
87  }
88  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
90 }
91 
92 PetscErrorCode DMRegister_MoFEM(const char sname[]) {
94  CHKERR DMRegister(sname, DMCreate_MoFEM);
96 }
97 
98 PetscErrorCode DMSetOperators_MoFEM(DM dm) {
99  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
101 
102  dm->ops->createglobalvector = DMCreateGlobalVector_MoFEM;
103  dm->ops->createlocalvector = DMCreateLocalVector_MoFEM;
104  dm->ops->creatematrix = DMCreateMatrix_MoFEM;
105  dm->ops->setup = DMSetUp_MoFEM;
106  dm->ops->destroy = DMDestroy_MoFEM;
107  dm->ops->setfromoptions = DMSetFromOptions_MoFEM;
108  dm->ops->globaltolocalbegin = DMGlobalToLocalBegin_MoFEM;
109  dm->ops->globaltolocalend = DMGlobalToLocalEnd_MoFEM;
110  dm->ops->localtoglobalbegin = DMLocalToGlobalBegin_MoFEM;
111  dm->ops->localtoglobalend = DMLocalToGlobalEnd_MoFEM;
112  dm->ops->createfieldis = DMCreateFieldIS_MoFEM;
113 
114  // Default matrix type
115  CHKERR DMSetMatType(dm, MATMPIAIJ);
116 
118 }
119 
120 PetscErrorCode DMCreate_MoFEM(DM dm) {
121  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
123  dm->data = new DMCtx();
126 }
127 
128 PetscErrorCode DMDestroy_MoFEM(DM dm) {
129  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
130  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
132  if (((DMCtx *)dm->data)->referenceNumber == 0) {
133  if (dm_field->destroyProblem) {
134  if (dm_field->mField_ptr->check_problem(dm_field->problemName)) {
135  dm_field->mField_ptr->delete_problem(dm_field->problemName);
136  } // else problem has to be deleted by the user
137  }
138  // cerr << "Destroy " << dm_field->problemName << endl;
139  delete ((DMCtx *)dm->data);
140  } else {
141  // cerr << "Dereference " << dm_field->problemName << " " <<
142  // ((DMCtx*)dm->data)->referenceNumber << endl;
143  (((DMCtx *)dm->data)->referenceNumber)--;
144  }
146 }
147 
148 PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr,
149  const char problem_name[],
150  const MoFEM::BitRefLevel bit_level,
151  const MoFEM::BitRefLevel bit_mask) {
153 
154  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
155  if (!dm->data) {
156  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
157  "data structure for MoFEM not yet created");
158  }
159  if (!m_field_ptr) {
160  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
161  "DM function not implemented into MoFEM");
162  }
163  dm_field->mField_ptr = m_field_ptr;
164  dm_field->problemName = problem_name;
165  if (!m_field_ptr->check_problem(dm_field->problemName)) {
166  // problem is not defined, declare problem internally set bool to
167  // destroyProblem problem with DM
168  dm_field->destroyProblem = PETSC_TRUE;
169  CHKERR dm_field->mField_ptr->add_problem(dm_field->problemName, MF_EXCL,
170  dm_field->verbosity);
171  } else {
172  dm_field->destroyProblem = PETSC_FALSE;
173  }
175  dm_field->problemName, bit_level);
177  dm_field->problemName, bit_mask);
178  dm_field->kspCtx =
179  boost::shared_ptr<KspCtx>(new KspCtx(*m_field_ptr, problem_name));
180  dm_field->snesCtx =
181  boost::shared_ptr<SnesCtx>(new SnesCtx(*m_field_ptr, problem_name));
182  dm_field->tsCtx =
183  boost::shared_ptr<TsCtx>(new TsCtx(*m_field_ptr, problem_name));
184 
185  MPI_Comm comm;
186  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
187  int result = 0;
188  MPI_Comm_compare(comm, m_field_ptr->get_comm(), &result);
189  // std::cerr << result << " " << MPI_IDENT << " " << MPI_CONGRUENT << " " <<
190  // MPI_SIMILAR << " " << MPI_UNEQUAL << std::endl;
191  if (result > MPI_CONGRUENT) {
192  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
193  "MoFEM and DM using different communicators");
194  }
195  MPI_Comm_size(comm, &dm_field->sIze);
196  MPI_Comm_rank(comm, &dm_field->rAnk);
197 
198  // problem structure
199  CHKERR dm_field->mField_ptr->get_problem(dm_field->problemName,
200  &dm_field->problemPtr);
201 
203 }
204 
205 PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]) {
207 
208  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
209  if (!dm->data) {
210  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
211  "data structure for MoFEM not yet created");
212  }
213  CHKERR DMMoFEMCreateMoFEM(subdm, dm_field->mField_ptr, problem_name,
214  dm_field->problemPtr->getBitRefLevel());
215 
216  DMCtx *subdm_field = (DMCtx *)subdm->data;
217  subdm_field->isSubDM = PETSC_TRUE;
218  subdm_field->problemMainOfSubPtr = dm_field->problemPtr;
219  subdm_field->isPartitioned = dm_field->isPartitioned;
220  subdm_field->isSquareMatrix = PETSC_FALSE;
221  subdm->ops->setup = DMSubDMSetUp_MoFEM;
222 
224 }
225 
226 PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[]) {
227  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
229  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
230  if (!dm->data) {
231  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
232  "data structure for MoFEM not yet created");
233  }
234  if (!dm_field->isSubDM) {
235  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
236  }
237  dm_field->rowFields.push_back(field_name);
239 }
240 
241 PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[]) {
242  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
244  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
245  if (!dm->data) {
246  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
247  "data structure for MoFEM not yet created");
248  }
249  if (!dm_field->isSubDM) {
250  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "this is not sub-dm");
251  }
252  dm_field->colFields.push_back(field_name);
254 }
255 
256 PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm) {
258  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
260  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
261  *is_sub_dm = dm_field->isSubDM;
263 }
264 
265 PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is) {
267  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
269  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
270  if (dm_field->isSubDM != PETSC_TRUE) {
271  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
272  "This DM is not created as a SubDM");
273  }
274  if (dm_field->isProblemBuild != PETSC_TRUE) {
275  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
276  }
277  boost::shared_ptr<Problem::SubProblemData> sub_data =
278  dm_field->problemPtr->getSubData();
279  CHKERR sub_data->getRowIs(is);
281 }
282 
283 PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is) {
285  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
287  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
288  if (dm_field->isSubDM != PETSC_TRUE) {
289  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
290  "This DM is not created as a SubDM");
291  }
292  if (dm_field->isProblemBuild != PETSC_TRUE) {
293  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Problem is not build");
294  }
295  boost::shared_ptr<Problem::SubProblemData> sub_data =
296  dm_field->problemPtr->getSubData();
297  CHKERR sub_data->getColIs(is);
299 }
300 
301 PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]) {
302  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
304  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
305  if (!dm->data) {
306  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
307  "data structure for MoFEM not yet created");
308  }
309  if (!dm_field->isCompDM) {
310  dm_field->isCompDM = PETSC_TRUE;
311  }
312  dm_field->rowCompPrb.push_back(prb_name);
313  if (dm_field->isSquareMatrix) {
314  dm_field->colCompPrb.push_back(prb_name);
315  }
317 }
318 
319 PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]) {
320  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
322  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
323  if (!dm->data) {
324  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
325  "data structure for MoFEM not yet created");
326  }
327  if (!dm_field->isCompDM) {
328  dm_field->isCompDM = PETSC_TRUE;
329  }
330  if (dm_field->isSquareMatrix) {
331  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
332  "No need to add problem on column when problem block structurally "
333  "symmetric");
334  }
335  dm_field->colCompPrb.push_back(prb_name);
337 }
338 
339 PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm) {
341  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
343  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
344  *is_comp_dm = dm_field->isCompDM;
346 }
347 
348 PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_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  *m_field_ptr = dm_field->mField_ptr;
358 }
359 
360 PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr) {
361  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
363  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
364  if (!dm->data) {
365  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
366  "data structure for MoFEM not yet created");
367  }
368  *problem_ptr = dm_field->problemPtr;
370 }
371 
372 PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem) {
374  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
376  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
377  dm_field->destroyProblem = destroy_problem;
379 }
380 
381 PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem) {
383  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
385  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
386  *destroy_problem = dm_field->destroyProblem;
388 }
389 
390 PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem) {
392  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
394  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
395  dm_field->isSquareMatrix = square_problem;
397 }
398 
399 PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[]) {
401  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
403  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
404  CHKERR dm_field->mField_ptr->resolve_shared_ents(dm_field->problemPtr,
405  fe_name);
407 }
408 
409 PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[],
410  PetscLayout *layout) {
412  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
414  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
415 
416  MPI_Comm comm;
417  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
418  CHKERR dm_field->problemPtr->getNumberOfElementsByNameAndPart(comm, fe_name,
419  layout);
421 }
422 
423 PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem) {
426  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
428  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
429  *square_problem = dm_field->isSquareMatrix;
431 }
432 
433 PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[]) {
434  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
436  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
438  dm_field->problemName, fe_name);
439  CHKERRG(ierr);
441 }
442 
443 PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[]) {
444  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
446  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
448  dm_field->problemName, fe_name);
449  CHKERRG(ierr);
451 }
452 
453 PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
454  ScatterMode scatter_mode) {
456  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
458  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
459  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setLocalGhostVector(
460  dm_field->problemPtr, COL, l, mode, scatter_mode);
461  CHKERRG(ierr);
463 }
464 
465 PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
466  ScatterMode scatter_mode) {
467  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
469  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
470  ierr = dm_field->mField_ptr->getInterface<VecManager>()->setGlobalGhostVector(
471  dm_field->problemPtr, COL, g, mode, scatter_mode);
472  CHKERRG(ierr);
474 }
475 
476 PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
477  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
479  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
481  dm_field->problemPtr, *method);
482  CHKERRG(ierr);
484 }
485 
486 PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method) {
487  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
489  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
491  dm_field->problemPtr, *method);
492  CHKERRG(ierr);
494 }
495 
496 PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
497  MoFEM::FEMethod *method,
498  int low_rank, int up_rank) {
500  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
501  ierr = dm_field->mField_ptr->loop_finite_elements(
502  dm_field->problemPtr, fe_name, *method, low_rank, up_rank);
503  CHKERRG(ierr);
505 }
506 
507 PetscErrorCode
508 DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const std::string &fe_name,
509  boost::shared_ptr<MoFEM::FEMethod> method,
510  int low_rank, int up_rank) {
511  return DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name.c_str(), method.get(),
512  low_rank, up_rank);
513 }
514 
515 PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
516  MoFEM::FEMethod *method) {
517  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
519  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
520  ierr = DMoFEMLoopFiniteElementsUpAndLowRank(dm, fe_name, method,
521  dm_field->rAnk, dm_field->rAnk);
522  CHKERRG(ierr);
524 }
525 
526 PetscErrorCode
527 DMoFEMLoopFiniteElements(DM dm, const std::string &fe_name,
528  boost::shared_ptr<MoFEM::FEMethod> method) {
529  return DMoFEMLoopFiniteElements(dm, fe_name.c_str(), method.get());
530 }
531 
532 PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
533  MoFEM::DofMethod *method) {
534  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
536  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
537  ierr =
538  dm_field->mField_ptr->loop_dofs(dm_field->problemPtr, field_name, COL,
539  *method, dm_field->rAnk, dm_field->rAnk);
540  CHKERRG(ierr);
542 }
543 
544 template <class S, class T0, class T1, class T2>
545 static PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, S fe_name, T0 method,
546  T1 pre_only, T2 post_only) {
547  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
549  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
550  if (pre_only) {
551  dm_field->kspCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
552  }
553  if (method) {
554  dm_field->kspCtx->get_loops_to_do_Rhs().push_back(
555  PairNameFEMethodPtr(fe_name, method));
556  }
557  if (post_only) {
558  dm_field->kspCtx->get_postProcess_to_do_Rhs().push_back(post_only);
559  }
560  CHKERR DMKSPSetComputeRHS(dm, KspRhs, dm_field->kspCtx.get());
562 }
563 
564 PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
565  MoFEM::FEMethod *method,
566  MoFEM::BasicMethod *pre_only,
567  MoFEM::BasicMethod *post_only) {
568  return DMMoFEMKSPSetComputeRHS<const char *, MoFEM::FEMethod *,
570  dm, fe_name, method, pre_only, post_only);
571 }
572 
573 PetscErrorCode
574 DMMoFEMKSPSetComputeRHS(DM dm, const std::string &fe_name,
575  boost::shared_ptr<MoFEM::FEMethod> method,
576  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
577  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
578  return DMMoFEMKSPSetComputeRHS<const std::string &,
579  boost::shared_ptr<MoFEM::FEMethod>,
580  boost::shared_ptr<MoFEM::BasicMethod>,
581  boost::shared_ptr<MoFEM::BasicMethod>>(
582  dm, fe_name, method, pre_only, post_only);
583 }
584 
585 template <class S, class T0, class T1, class T2>
586 static PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, S fe_name, T0 method,
587  T1 pre_only, T2 post_only) {
588  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
590  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
591  if (pre_only) {
592  dm_field->kspCtx->get_preProcess_to_do_Mat().push_back(pre_only);
593  }
594  if (method) {
595  dm_field->kspCtx->get_loops_to_do_Mat().push_back(
596  PairNameFEMethodPtr(fe_name, method));
597  }
598  if (post_only) {
599  dm_field->kspCtx->get_postProcess_to_do_Mat().push_back(post_only);
600  }
601  CHKERR DMKSPSetComputeOperators(dm, KspMat, dm_field->kspCtx.get());
603 }
604 
605 PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
606  MoFEM::FEMethod *method,
607  MoFEM::BasicMethod *pre_only,
608  MoFEM::BasicMethod *post_only) {
609  return DMMoFEMKSPSetComputeOperators<const char *, MoFEM::FEMethod *,
612  dm, fe_name, method, pre_only, post_only);
613 }
614 
615 PetscErrorCode
616 DMMoFEMKSPSetComputeOperators(DM dm, const std::string &fe_name,
617  boost::shared_ptr<MoFEM::FEMethod> method,
618  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
619  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
620  return DMMoFEMKSPSetComputeOperators<const std::string &,
621  boost::shared_ptr<MoFEM::FEMethod>>(
622  dm, fe_name, method, pre_only, post_only);
623 }
624 
625 template <class S, class T0, class T1, class T2>
626 static PetscErrorCode DMMoFEMSNESSetFunction(DM dm, S fe_name, T0 method,
627  T1 pre_only, T2 post_only) {
628  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
630  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
631  if (pre_only) {
632  dm_field->snesCtx->get_preProcess_to_do_Rhs().push_back(pre_only);
633  }
634  if (method) {
635  dm_field->snesCtx->get_loops_to_do_Rhs().push_back(
636  PairNameFEMethodPtr(fe_name, method));
637  }
638  if (post_only) {
639  dm_field->snesCtx->get_postProcess_to_do_Rhs().push_back(post_only);
640  }
641  CHKERR DMSNESSetFunction(dm, SnesRhs, dm_field->snesCtx.get());
643 }
644 
645 PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
646  MoFEM::FEMethod *method,
647  MoFEM::BasicMethod *pre_only,
648  MoFEM::BasicMethod *post_only) {
649  return DMMoFEMSNESSetFunction<const char *, MoFEM::FEMethod *,
651  dm, fe_name, method, pre_only, post_only);
652 }
653 
654 PetscErrorCode
655 DMMoFEMSNESSetFunction(DM dm, const std::string &fe_name,
656  boost::shared_ptr<MoFEM::FEMethod> method,
657  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
658  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
659  return DMMoFEMSNESSetFunction<const std::string &,
660  boost::shared_ptr<MoFEM::FEMethod>,
661  boost::shared_ptr<MoFEM::BasicMethod>,
662  boost::shared_ptr<MoFEM::BasicMethod>>(
663  dm, fe_name, method, pre_only, post_only);
664 }
665 
666 template <class S, class T0, class T1, class T2>
667 static PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, S fe_name, T0 method,
668  T1 pre_only, T2 post_only) {
669  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
671  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
672  if (pre_only) {
673  dm_field->snesCtx->get_preProcess_to_do_Mat().push_back(pre_only);
674  }
675  if (method) {
676  dm_field->snesCtx->get_loops_to_do_Mat().push_back(
677  PairNameFEMethodPtr(fe_name, method));
678  }
679  if (post_only) {
680  dm_field->snesCtx->get_postProcess_to_do_Mat().push_back(post_only);
681  }
682  CHKERR DMSNESSetJacobian(dm, SnesMat, dm_field->snesCtx.get());
684 }
685 
686 PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
687  MoFEM::FEMethod *method,
688  MoFEM::BasicMethod *pre_only,
689  MoFEM::BasicMethod *post_only) {
690  return DMMoFEMSNESSetJacobian<const char *, MoFEM::FEMethod *,
692  dm, fe_name, method, pre_only, post_only);
693 }
694 
695 PetscErrorCode
696 DMMoFEMSNESSetJacobian(DM dm, const std::string &fe_name,
697  boost::shared_ptr<MoFEM::FEMethod> method,
698  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
699  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
700  return DMMoFEMSNESSetJacobian<const std::string &,
701  boost::shared_ptr<MoFEM::FEMethod>,
702  boost::shared_ptr<MoFEM::BasicMethod>,
703  boost::shared_ptr<MoFEM::BasicMethod>>(
704  dm, fe_name, method, pre_only, post_only);
705 }
706 
707 template <class S, class T0, class T1, class T2>
708 static PetscErrorCode DMMoFEMTSSetIFunction(DM dm, S fe_name, T0 method,
709  T1 pre_only, T2 post_only) {
710  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
712  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
713  if (pre_only) {
714  dm_field->tsCtx->get_preProcess_to_do_IFunction().push_back(pre_only);
715  }
716  if (method) {
717  dm_field->tsCtx->get_loops_to_do_IFunction().push_back(
718  PairNameFEMethodPtr(fe_name, method));
719  }
720  if (post_only) {
721  dm_field->tsCtx->get_postProcess_to_do_IFunction().push_back(post_only);
722  }
723  CHKERR DMTSSetIFunction(dm, TsSetIFunction, dm_field->tsCtx.get());
725 }
726 
727 PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
728  MoFEM::FEMethod *method,
729  MoFEM::BasicMethod *pre_only,
730  MoFEM::BasicMethod *post_only) {
731  return DMMoFEMTSSetIFunction<const char *, MoFEM::FEMethod *,
733  dm, fe_name, method, pre_only, post_only);
735 }
736 
737 PetscErrorCode
738 DMMoFEMTSSetIFunction(DM dm, const std::string &fe_name,
739  boost::shared_ptr<MoFEM::FEMethod> method,
740  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
741  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
742  return DMMoFEMTSSetIFunction<const std::string,
743  boost::shared_ptr<MoFEM::FEMethod>,
744  boost::shared_ptr<MoFEM::BasicMethod>,
745  boost::shared_ptr<MoFEM::BasicMethod>>(
746  dm, fe_name, method, pre_only, post_only);
748 }
749 
750 template <class S, class T0, class T1, class T2>
751 static PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, S fe_name, T0 method,
752  T1 pre_only, T2 post_only) {
753  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
755  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
756  if (pre_only) {
757  dm_field->tsCtx->get_preProcess_to_do_IJacobian().push_back(pre_only);
758  }
759  if (method) {
760  dm_field->tsCtx->get_loops_to_do_IJacobian().push_back(
761  PairNameFEMethodPtr(fe_name, method));
762  }
763  if (post_only) {
764  dm_field->tsCtx->get_postProcess_to_do_IJacobian().push_back(post_only);
765  }
766  CHKERR DMTSSetIJacobian(dm, TsSetIJacobian, dm_field->tsCtx.get());
768 }
769 
770 PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
771  MoFEM::FEMethod *method,
772  MoFEM::BasicMethod *pre_only,
773  MoFEM::BasicMethod *post_only) {
774  return DMMoFEMTSSetIJacobian<const char *, FEMethod *, MoFEM::BasicMethod *,
775  MoFEM::BasicMethod *>(dm, fe_name, method,
776  pre_only, post_only);
777 }
778 
779 PetscErrorCode
780 DMMoFEMTSSetIJacobian(DM dm, const std::string &fe_name,
781  boost::shared_ptr<MoFEM::FEMethod> method,
782  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
783  boost::shared_ptr<MoFEM::BasicMethod> post_only) {
784  return DMMoFEMTSSetIJacobian<const std::string &,
785  boost::shared_ptr<MoFEM::FEMethod>,
786  boost::shared_ptr<MoFEM::BasicMethod>,
787  boost::shared_ptr<MoFEM::BasicMethod>>(
788  dm, fe_name, method, pre_only, post_only);
789 }
790 
791 PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx) {
792  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
794  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
795  *ksp_ctx = dm_field->kspCtx.get();
797 }
798 
799 PetscErrorCode
800 DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
801  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
803  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
804  const_cast<boost::shared_ptr<MoFEM::KspCtx> &>(ksp_ctx) = dm_field->kspCtx;
806 }
807 
808 PetscErrorCode DMMoFEMSetKspCtx(DM dm,
809  boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx) {
810  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
812  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
813  dm_field->kspCtx = ksp_ctx;
815 }
816 
817 PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx) {
818  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
820  DMCtx *dm_field = (DMCtx *)dm->data;
821  *snes_ctx = dm_field->snesCtx.get();
823 }
824 
825 PetscErrorCode
826 DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
827  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
829  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
830  const_cast<boost::shared_ptr<MoFEM::SnesCtx> &>(snes_ctx) = dm_field->snesCtx;
832 }
833 
834 PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
835  boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx) {
836  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
838  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
839  dm_field->snesCtx = snes_ctx;
841 }
842 
843 /** get if read mesh is partitioned
844  * \ingroup dm
845  */
846 PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned) {
847  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
849  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
850  dm_field->isPartitioned = is_partitioned;
852 }
853 
854 /** get if read mesh is partitioned
855  * \ingroup dm
856  */
857 PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned) {
858  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
860  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
861  *is_partitioned = dm_field->isPartitioned;
863 }
864 
865 PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx) {
866  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
868  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
869  *ts_ctx = dm_field->tsCtx.get();
871 }
872 
873 PetscErrorCode DMMoFEMGetTsCtx(DM dm,
874  const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
875  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
877  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
878  const_cast<boost::shared_ptr<MoFEM::TsCtx> &>(ts_ctx) = dm_field->tsCtx;
880 }
881 
882 PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> &ts_ctx) {
883  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
885  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
886  dm_field->tsCtx = ts_ctx;
888 }
889 
890 PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g) {
891  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
893  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
894  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateGhost(
895  dm_field->problemName, COL, g);
897 }
898 
899 PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l) {
900  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
902  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
903  CHKERR dm_field->mField_ptr->getInterface<VecManager>()->vecCreateSeq(
904  dm_field->problemName, COL, l);
906 }
907 
908 PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M) {
909  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
911  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
912  if (strcmp(dm->mattype, MATMPIAIJ) == 0) {
914  ->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(dm_field->problemName,
915  M);
916  } else if (strcmp(dm->mattype, MATAIJ) == 0) {
918  ->createSeqAIJWithArrays<PetscLocalIdx_mi_tag>(dm_field->problemName,
919  M);
920  } else {
921  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
922  "Matrix type not implemented");
923  }
925 }
926 
927 #if PETSC_VERSION_GE(3, 7, 0)
928 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
929  DM dm) {
930 #elif PETSC_VERSION_GE(3, 5, 3)
931 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm) {
932 #else
933 PetscErrorCode DMSetFromOptions_MoFEM(DM dm) {
934 #endif
935 
936  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
938  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
939 #if PETSC_VERSION_GE(3, 5, 3)
940  ierr = PetscOptionsHead(PetscOptionsObject, "DMMoFEM Options");
941  CHKERRG(ierr);
942 #else
943  ierr = PetscOptionsHead("DMMoFEM Options");
944  CHKERRG(ierr);
945 #endif
946  ierr = PetscOptionsBool("-dm_is_partitioned",
947  "set if mesh is partitioned (works which native MOAB "
948  "file format, i.e. h5m",
949  "DMSetUp", dm_field->isPartitioned,
950  &dm_field->isPartitioned, NULL);
951  CHKERRG(ierr);
953 }
954 
955 PetscErrorCode DMSetUp_MoFEM(DM dm) {
956  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
957  ProblemsManager *prb_mng_ptr;
959  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
960  CHKERR dm_field->mField_ptr->getInterface(prb_mng_ptr);
961 
962  if (dm_field->isCompDM) {
963  // It is composite probelm
964  CHKERR prb_mng_ptr->buildCompsedProblem(
965  dm_field->problemName, dm_field->rowCompPrb, dm_field->colCompPrb,
966  dm_field->isSquareMatrix == PETSC_TRUE, dm_field->verbosity);
967  } else {
968  if (dm_field->isPartitioned) {
970  dm_field->problemName, dm_field->isSquareMatrix == PETSC_TRUE,
971  dm_field->verbosity);
972  } else {
973  CHKERR prb_mng_ptr->buildProblem(dm_field->problemName,
974  dm_field->isSquareMatrix == PETSC_TRUE,
975  dm_field->verbosity);
976  CHKERR prb_mng_ptr->partitionProblem(dm_field->problemName,
977  dm_field->verbosity);
978  }
979  }
980 
981  // Partition finite elements
982  if (dm_field->isPartitioned) {
983  CHKERR prb_mng_ptr->partitionFiniteElements(
984  dm_field->problemName, true, 0, dm_field->sIze, dm_field->verbosity);
986  dm_field->problemName, dm_field->verbosity);
987  } else {
988  // partition finite elemnets
989  CHKERR prb_mng_ptr->partitionFiniteElements(dm_field->problemName, false,
990  -1, -1, dm_field->verbosity);
991  // Get ghost DOFs
992  CHKERR prb_mng_ptr->partitionGhostDofs(dm_field->problemName,
993  dm_field->verbosity);
994  }
995 
996  // Set flag that problem is build and partitioned
997  dm_field->isProblemBuild = PETSC_TRUE;
998 
1000 }
1001 
1002 PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm) {
1003  PetscValidHeaderSpecific(subdm, DM_CLASSID, 1);
1004  ProblemsManager *prb_mng_ptr;
1006 
1007  DMCtx *subdm_field = static_cast<DMCtx *>(subdm->data);
1008 
1009  // build sub dm problem
1010  CHKERR subdm_field->mField_ptr->getInterface(prb_mng_ptr);
1011  CHKERR prb_mng_ptr->buildSubProblem(
1012  subdm_field->problemName, subdm_field->rowFields, subdm_field->colFields,
1013  subdm_field->problemMainOfSubPtr->getName(),
1014  subdm_field->isSquareMatrix == PETSC_TRUE, subdm_field->verbosity);
1015 
1016  // partition problem
1017  subdm_field->isPartitioned = subdm_field->isPartitioned;
1018  if (subdm_field->isPartitioned) {
1019  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, true,
1020  0, subdm_field->sIze,
1021  subdm_field->verbosity);
1022  // set ghost nodes
1024  subdm_field->problemName, subdm_field->verbosity);
1025  } else {
1026  // partition finite elements
1027  CHKERR prb_mng_ptr->partitionFiniteElements(subdm_field->problemName, false,
1028  -1, -1, subdm_field->verbosity);
1029  // set ghost nodes
1030  CHKERR prb_mng_ptr->partitionGhostDofs(subdm_field->problemName,
1031  subdm_field->verbosity);
1032  }
1033 
1034  subdm_field->isProblemBuild = PETSC_TRUE;
1035 
1037 }
1038 
1039 PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec g, InsertMode mode,
1040  Vec l) {
1041  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1043  CHKERR VecGhostUpdateBegin(g, INSERT_VALUES, SCATTER_FORWARD);
1045 }
1046 
1047 PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec g, InsertMode mode, Vec l) {
1049  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1051 
1052  CHKERR VecGhostUpdateEnd(g, INSERT_VALUES, SCATTER_FORWARD);
1053 
1054  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1055  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1056  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1057 
1058  double *array_loc, *array_glob;
1059  CHKERR VecGetArray(l, &array_loc);
1060  CHKERR VecGetArray(g, &array_glob);
1061  switch (mode) {
1062  case INSERT_VALUES:
1063  cblas_dcopy(nb_dofs + nb_ghost, array_glob, 1, array_loc, 1);
1064  break;
1065  case ADD_VALUES:
1066  cblas_daxpy(nb_dofs + nb_ghost, 1, array_glob, 1, array_loc, 1);
1067  break;
1068  default:
1069  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1070  }
1071  CHKERR VecRestoreArray(l, &array_loc);
1072  CHKERR VecRestoreArray(g, &array_glob);
1074 }
1075 
1076 PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM dm, Vec l, InsertMode mode,
1077  Vec g) {
1078 
1079  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1081 
1082  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1083  int nb_dofs = dm_field->problemPtr->getNbLocalDofsRow();
1084  int nb_ghost = dm_field->problemPtr->getNbGhostDofsRow();
1085 
1086  double *array_loc, *array_glob;
1087  CHKERR VecGetArray(l, &array_loc);
1088  CHKERR VecGetArray(g, &array_glob);
1089  switch (mode) {
1090  case INSERT_VALUES:
1091  cblas_dcopy(nb_dofs + nb_ghost, array_loc, 1, array_glob, 1);
1092  break;
1093  case ADD_VALUES:
1094  cblas_daxpy(nb_dofs + nb_ghost, 1, array_loc, 1, array_glob, 1);
1095  break;
1096  default:
1097  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
1098  }
1099  CHKERR VecRestoreArray(l, &array_loc);
1100  CHKERR VecRestoreArray(g, &array_glob);
1101 
1103 }
1104 
1105 PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec l, InsertMode mode, Vec g) {
1106  //
1109 }
1110 
1111 PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
1112  char ***fieldNames, IS **fields) {
1113  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1115 
1116  if (numFields) {
1117  *numFields = 0;
1118  }
1119  if (fieldNames) {
1120  *fieldNames = NULL;
1121  }
1122  if (fields) {
1123  *fields = NULL;
1124  }
1125 
1126  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1127  const Field_multiIndex *fields_ptr;
1128  CHKERR dm_field->mField_ptr->get_fields(&fields_ptr);
1129  Field_multiIndex::iterator fit, hi_fit;
1130  fit = fields_ptr->begin();
1131  hi_fit = fields_ptr->end();
1132  *numFields = std::distance(fit, hi_fit);
1133 
1134  if (fieldNames) {
1135  CHKERR PetscMalloc1(*numFields, fieldNames);
1136  }
1137  if (fields) {
1138  CHKERR PetscMalloc1(*numFields, fields);
1139  }
1140 
1141  for (int f = 0; fit != hi_fit; fit++, f++) {
1142  if (fieldNames) {
1143  CHKERR PetscStrallocpy(fit->get()->getName().c_str(),
1144  (char **)&(*fieldNames)[f]);
1145  }
1146  if (fields) {
1147  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1148  ->isCreateProblemFieldAndRank(
1149  dm_field->problemPtr->getName(), ROW, fit->get()->getName(), 0,
1150  fit->get()->getNbOfCoeffs(), &(*fields)[f]);
1151  }
1152  }
1153 
1155 }
1156 
1157 PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
1158  IS *is) {
1159  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1161  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1162  CHKERR dm_field->mField_ptr->getInterface<ISManager>()
1163  ->isCreateProblemFieldAndRank(dm_field->problemPtr->getName(), ROW,
1164  field_name, 0, 1000, is);
1166 }
1167 
1168 PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb) {
1169  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
1171  DMCtx *dm_field = static_cast<DMCtx *>(dm->data);
1172  dm_field->verbosity = verb;
1174 }
1175 
1176 } // 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:265
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:645
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:1168
int verbosity
verbosity
Definition: DMMoFEM.hpp:770
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:564
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:241
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:1039
boost::shared_ptr< SubProblemData > & getSubData() const
Get main problem of sub-problem is.
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:283
PetscBool isSubDM
Definition: DMMoFEM.hpp:756
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:34
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
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:453
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:390
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMMoFEM.cpp:256
PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:399
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
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:433
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:142
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 ...
MoFEMErrorCode buildProblem(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
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
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:955
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
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:319
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:80
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.
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
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
Data structure to exchange data between mofem and User Loop Methods on entities.It allows to exchange...
virtual ~DMCtx()
Definition: DMMMoFEM.cpp:75
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:1076
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:908
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:605
DofIdx getNbGhostDofsRow() const
RowColData
RowColData.
Definition: definitions.h:184
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:515
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:339
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:128
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM...
Definition: DMMMoFEM.cpp:890
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[], PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMMoFEM.cpp:409
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:205
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:301
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:882
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 TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:56
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1002
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:423
useful compiler directives and definitions
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMMoFEM.cpp:817
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:465
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:92
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMMoFEM.cpp:381
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:791
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:686
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:443
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:476
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:1047
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:486
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:148
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:120
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:98
#define CHKERR
Inline error check.
Definition: definitions.h:594
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1111
Field data structure storing information about space, approximation base, coordinate systems...
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:758
MoFEMErrorCode partitionGhostDofs(const std::string name, int verb=VERBOSE)
determine ghost nodes
Interface for creating matrices and managing matricesCreating and managing matrices.
Functions to approximate hierarchical spaces.
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...
BitRefLevel getBitRefLevel() const
PetscBool isCompDM
Definition: DMMoFEM.hpp:761
MoFEM interface.
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:762
Multindex 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:865
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:1157
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:348
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.
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:737
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:372
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:727
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
Definition: DMMMoFEM.cpp:226
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:360
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMMoFEM.cpp:857
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:532
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:496
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:808
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:846
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:780
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...
Context for PETSc Time Stepping.
Matrix manager is used to build and partition problems.
Interface managing vectorsManaging problems, build and partitioning.
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1105
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 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:136
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:933
Context for PETSc SNES, i.e. nonlinear solver.
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM...
Definition: DMMMoFEM.cpp:899
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:834