v0.8.23
DMMoFEM.hpp
Go to the documentation of this file.
1 /** \file DMMoFEM.hpp
2  \brief Discrete manager interface for MoFEM
3  */
4 
5 /* This file is part of MoFEM.
6  * MoFEM is free software: you can redistribute it and/or modify it under
7  * the terms of the GNU Lesser General Public License as published by the
8  * Free Software Foundation, either version 3 of the License, or (at your
9  * option) any later version.
10  *
11  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
12  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
14  * License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public
17  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
18 
19 #ifndef __DMMMOFEM_H
20 #define __DMMMOFEM_H
21 
22 #define DM_NO_ELEMENT "DMNONEFE"
23 
24 namespace MoFEM {
25 
26 /**
27  * \brief Register MoFEM problem
28  * \ingroup dm
29  */
30 PetscErrorCode DMRegister_MoFEM(const char sname[]);
31 
32 /**
33  * \brief Must be called by user to set MoFEM data structures
34  * \ingroup dm
35  */
36 PetscErrorCode DMMoFEMCreateMoFEM(
37  DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[],
38  const MoFEM::BitRefLevel bit_level,
39  const MoFEM::BitRefLevel bit_mask = MoFEM::BitRefLevel().set());
40 
41 /**
42  * \brief Must be called by user to set Sub DM MoFEM data structures
43  * \ingroup dm
44  */
45 PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]);
46 
47 /**
48  * \brief Get pointer to MoFEM::Interface
49  * @param dm Distributed mesh manager
50  * @param m_field_ptr Pointer to pointer of field interface
51  * @return Error code
52  * \ingroup dm
53  */
54 PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr);
55 
56 /**
57  * \brief Get pointer to problem data structure
58  * \ingroup dm
59  */
60 PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr);
61 
62 /**
63  * If this is set to PETSC_TRUE problem is deleted with DM
64  * @param dm the DM object
65  * @param destroy if PETSC_TRUE problem is destroyed
66  * @return error code
67  */
68 PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem);
69 
70 /**
71  * Get if problem will be destroyed with DM
72  * @param dm the DM object
73  * @param destroy return if PETSC_TRUE problem is destroyed
74  * @return error code
75  */
76 PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem);
77 
78 /**
79  * \brief set squared problem
80  * \ingroup dm
81 
82  It if true is assumed that matrix has the same indexing on rows and
83  columns. This reduces interprocessor communication.
84 
85  */
86 PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem);
87 
88 /**
89  * \brief get squared problem
90  * \ingroup dm
91 
92  It if true is assumed that matrix has the same indexing on rows and
93  columns. This reduces interprocessor communication.
94 
95  */
96 PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem);
97 
98 /**
99  * \brief Resolve shared entities
100  * @param dm dm
101  * @param fe_name finite element for which shared entities are resolved
102  * @return error code
103  *
104 
105  * This allows for tag reduction or tag exchange, f.e.
106 
107  \code
108  ierr = DMMoFEMGetSquareProblem(dm,"SHELL_ELEMENT"); CHKERRG(ierr);
109  Tag th;
110  rval = mField.get_moab().tag_get_handle("ADAPT_ORDER",th); CHKERRQ_MOAB(rval);
111  ParallelComm* pcomm =
112  ParallelComm::get_pcomm(&mField.get_moab(),MYPCOMM_INDEX);
113  // rval = pcomm->reduce_tags(th,MPI_SUM,prisms);
114  rval = pcomm->exchange_tags(th,prisms);
115  \endcode
116 
117  * \ingroup dm
118  */
119 PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[]);
120 
121 /**
122  * \brief Get finite elements layout in the problem
123  *
124  * In layout is stored information how many elements is on each processor, for
125  * more information look int petsc documentation
126  * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/IS/PetscLayoutCreate.html#PetscLayoutCreate>
127  *
128  * @param dm discrete manager for this problem
129  * @param fe_name finite element name
130  * @param layout pointer to layout, for created layout user takes
131  * responsibility for destroying it.
132  * @return error code
133  *
134  * \ingroup dm
135  */
136 PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[],
137  PetscLayout *layout);
138 
139 /**
140  * \brief add element to dm
141  * \ingroup dm
142  *
143  * \note add_file is a collective, should be executed on all processors.
144  * Otherwise could lead to deadlock.
145  *
146  */
147 PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[]);
148 
149 /**
150  * \brief unset element from dm
151  * \ingroup dm
152  */
153 PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[]);
154 
155 /**
156  * \brief set local (or ghosted) vector values on mesh for partition only
157  * \ingroup dm
158 
159  * \param l vector
160  * \param mode see petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
161  * \param scatter_mode see petsc manual for ScatterMode (The available modes
162  are: SCATTER_FORWARD or SCATTER_REVERSE)
163  *
164  * SCATTER_REVERSE set data to field entities from V vector.
165  *
166  * SCATTER_FORWARD set vector V from data field entities
167 
168  */
169 PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
170  ScatterMode scatter_mode);
171 
172 /**
173  * \brief set ghosted vector values on all existing mesh entities
174  * \ingroup dm
175 
176  * \param g vector
177  * \param mode see petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
178  * \param scatter_mode see petsc manual for ScatterMode (The available modes
179  are: SCATTER_FORWARD or SCATTER_REVERSE)
180  *
181  * SCATTER_REVERSE set data to field entities from V vector.
182  *
183  * SCATTER_FORWARD set vector V from data field entities
184 
185  */
186 PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
187  ScatterMode scatter_mode);
188 
189 /**
190  * \brief execute finite element method for each element in dm (problem)
191  * \ingroup dm
192  */
193 PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method);
194 
195 /**
196  * \brief execute finite element method for each element in dm (problem)
197  * \ingroup dm
198  */
199 PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method);
200 
201 /**
202  * \brief Executes FEMethod for finite elements in DM
203  * @param dm MoFEM discrete manager
204  * @param fe_name name of finite element
205  * @param method pointer to \ref MoFEM::FEMethod
206  * @param low_rank lowest rank of processor
207  * @param up_rank upper run of processor
208  * @return Error code
209  * \ingroup dm
210  */
211 PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
212  MoFEM::FEMethod *method,
213  int low_rank, int up_rank);
214 
215 /**
216  * \brief Executes FEMethod for finite elements in DM
217  * @param dm MoFEM discrete manager
218  * @param fe_name name of finite element
219  * @param method pointer to \ref MoFEM::FEMethod
220  * @param low_rank lowest rank of processor
221  * @param up_rank upper run of processor
222  * @return Error code
223  * \ingroup dm
224  */
225 PetscErrorCode
226 DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const std::string fe_name,
227  boost::shared_ptr<MoFEM::FEMethod> method,
228  int low_rank, int up_rank);
229 
230 /**
231  * \brief Executes FEMethod for finite elements in DM
232  * @param dm MoFEM discrete manager
233  * @param fe_name name of element
234  * @param method pointer to \ref MOFEM::FEMethod
235  * @return Error code
236  * \ingroup dm
237  */
238 PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
239  MoFEM::FEMethod *method);
240 
241 /**
242  * \brief Executes FEMethod for finite elements in DM
243  * @param dm MoFEM discrete manager
244  * @param fe_name name of element
245  * @param method pointer to \ref MOFEM::FEMethod
246  * @return Error code
247  * \ingroup dm
248  */
249 PetscErrorCode
250 DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
251  boost::shared_ptr<MoFEM::FEMethod> method);
252 
253 /**
254  * \brief execute method for dofs on field in problem
255  * \ingroup dm
256  */
257 PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
258  MoFEM::DofMethod *method);
259 
260 // /**
261 // * \brief Set compute operator for KSP solver via sub-matrix and IS
262 // *
263 // * @param dm DM
264 // * @return error code
265 // *
266 // * \ingroup dm
267 // */
268 // PetscErrorCode DMMoFEMKSPSetComputeOperatorsViaSubMatrixbByIs(DM dm);
269 
270 /**
271  * \brief set KSP right hand side evaluation function
272  * \ingroup dm
273  */
274 PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
275  MoFEM::FEMethod *method,
276  MoFEM::BasicMethod *pre_only,
277  MoFEM::BasicMethod *post_only);
278 
279 /**
280  * \brief set KSP right hand side evaluation function
281  * \ingroup dm
282  */
283 PetscErrorCode
284 DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
285  boost::shared_ptr<MoFEM::FEMethod> method,
286  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
287  boost::shared_ptr<MoFEM::BasicMethod> post_only);
288 
289 /**
290  * \brief Set KSP operators and push mofem finite element methods
291  *
292  * @param dm DM
293  * @param fe_name finite element name
294  * @param method method on the element (executed for each element in the
295  * problem which given name)
296  * @param pre_only method for pre-process before element method
297  * @param post_only method for post-process after element method
298  * @return error code
299  *
300  * \ingroup dm
301  */
302 PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
303  MoFEM::FEMethod *method,
304  MoFEM::BasicMethod *pre_only,
305  MoFEM::BasicMethod *post_only);
306 
307 /**
308  * \brief Set KSP operators and push mofem finite element methods
309  *
310  * @param dm DM
311  * @param fe_name finite element name
312  * @param method method on the element (executed for each element in the
313  * problem which given name)
314  * @param pre_only method for pre-process before element method
315  * @param post_only method for post-process after element method
316  * @return error code
317  *
318  * \ingroup dm
319  */
320 PetscErrorCode
321 DMMoFEMKSPSetComputeOperators(DM dm, const std::string fe_name,
322  boost::shared_ptr<MoFEM::FEMethod> method,
323  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
324  boost::shared_ptr<MoFEM::BasicMethod> post_only);
325 
326 /**
327  * \brief set SNES residual evaluation function
328  * \ingroup dm
329  */
330 PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
331  MoFEM::FEMethod *method,
332  MoFEM::BasicMethod *pre_only,
333  MoFEM::BasicMethod *post_only);
334 
335 /**
336  * \brief set SNES residual evaluation function
337  * \ingroup dm
338  */
339 PetscErrorCode
340 DMMoFEMSNESSetFunction(DM dm, const std::string fe_name,
341  boost::shared_ptr<MoFEM::FEMethod> method,
342  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
343  boost::shared_ptr<MoFEM::BasicMethod> post_only);
344 
345 /**
346  * \brief set SNES Jacobian evaluation function
347  * \ingroup dm
348  */
349 PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
350  MoFEM::FEMethod *method,
351  MoFEM::BasicMethod *pre_only,
352  MoFEM::BasicMethod *post_only);
353 
354 /**
355  * \brief set SNES Jacobian evaluation function
356  * \ingroup dm
357  */
358 PetscErrorCode
359 DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
360  boost::shared_ptr<MoFEM::FEMethod> method,
361  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
362  boost::shared_ptr<MoFEM::BasicMethod> post_only);
363 
364 /**
365  * \brief set TS implicit function evaluation function
366  * \ingroup dm
367  */
368 PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
369  MoFEM::FEMethod *method,
370  MoFEM::BasicMethod *pre_only,
371  MoFEM::BasicMethod *post_only);
372 
373 /**
374  * \brief set TS implicit function evaluation function
375  * \ingroup dm
376  */
377 PetscErrorCode
378 DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
379  boost::shared_ptr<MoFEM::FEMethod> method,
380  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
381  boost::shared_ptr<MoFEM::BasicMethod> post_only);
382 
383 /**
384  * \brief set TS Jacobian evaluation function
385  * \ingroup dm
386  */
387 PetscErrorCode
388 DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
389  boost::shared_ptr<MoFEM::FEMethod> method,
390  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
391  boost::shared_ptr<MoFEM::BasicMethod> post_only);
392 
393 /**
394  * \brief set TS Jacobian evaluation function
395  * \ingroup dm
396  */
397 PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
398  MoFEM::FEMethod *method,
399  MoFEM::BasicMethod *pre_only,
400  MoFEM::BasicMethod *post_only);
401 
402 /**
403  * @brief set TS the right hand side function
404  *
405  * <a
406  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetRHSFunction.html#TSSetRHSFunction>See
407  * petsc documentation</a>
408  *
409  * @param dm
410  * @param fe_name
411  * @param method
412  * @param pre_only
413  * @param post_only
414  * @return PetscErrorCode
415  */
416 PetscErrorCode
417 DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name,
418  boost::shared_ptr<MoFEM::FEMethod> method,
419  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
420  boost::shared_ptr<MoFEM::BasicMethod> post_only);
421 
422 /**
423  * @brief set TS the right hand side function
424  *
425  * <a
426  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetRHSFunction.html#TSSetRHSFunction>See
427  * petsc documentation</a>
428  *
429  * @param dm
430  * @param fe_name
431  * @param method
432  * @param pre_only
433  * @param post_only
434  * @return PetscErrorCode
435  */
436 PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const char fe_name[],
437  MoFEM::FEMethod *method,
438  MoFEM::BasicMethod *pre_only,
439  MoFEM::BasicMethod *post_only);
440 
441 /**
442  * @brief set TS the right hand side jacobian
443  *
444  * <a
445  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetRHSJacobian.html>See
446  * petsc documentation</a>
447  *
448  * @param dm
449  * @param fe_name
450  * @param method
451  * @param pre_only
452  * @param post_only
453  * @return PetscErrorCode
454  */
455 PetscErrorCode
456 DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
457  boost::shared_ptr<MoFEM::FEMethod> method,
458  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
459  boost::shared_ptr<MoFEM::BasicMethod> post_only);
460 
461 /**
462  * @brief set TS the right hand side jacobian
463  *
464  * <a
465  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetRHSJacobian.html>See
466  * petsc documentation</a>
467  *
468  * @param dm
469  * @param fe_name
470  * @param method
471  * @param pre_only
472  * @param post_only
473  * @return PetscErrorCode
474  */
475 PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const char fe_name[],
476  MoFEM::FEMethod *method,
477  MoFEM::BasicMethod *pre_only,
478  MoFEM::BasicMethod *post_only);
479 
480 /**
481  * @brief Set Monitor To TS solver
482  *
483  * <a
484  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSMonitorSet.html>See
485  * PETSc documentaton here</a>
486  *
487  * @param dm
488  * @param ts time solver
489  * @param fe_name
490  * @param method
491  * @param pre_only
492  * @param post_only
493  * @return PetscErrorCod
494  */
495 PetscErrorCode
496 DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name,
497  boost::shared_ptr<MoFEM::FEMethod> method,
498  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
499  boost::shared_ptr<MoFEM::BasicMethod> post_only);
500 
501 /**
502  * @brief Set Monitor To TS solver
503  *
504  * <a
505  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSMonitorSet.html>See
506  * PETSc documentaton here</a>
507  *
508  * @param dm
509  * @param ts time solver
510  * @param fe_name
511  * @param method
512  * @param pre_only
513  * @param post_only
514  * @return PetscErrorCod
515  */
516 PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const char fe_name[],
517  MoFEM::FEMethod *method,
518  MoFEM::BasicMethod *pre_only,
519  MoFEM::BasicMethod *post_only);
520 
521 /**
522  * \brief get MoFEM::KspCtx data structure
523  * \ingroup dm
524  */
525 PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx);
526 
527 /**
528  * \brief get MoFEM::KspCtx data structure
529  * \ingroup dm
530  */
531 PetscErrorCode
532 DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx);
533 
534 /**
535  * \brief set MoFEM::KspCtx data structure
536  * \ingroup dm
537  */
538 PetscErrorCode DMMoFEMSetKspCtx(DM dm,
539  boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx);
540 
541 /**
542  * \brief get MoFEM::SnesCtx data structure
543  * \ingroup dm
544  */
545 PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx);
546 
547 /**
548  * \brief get MoFEM::SnesCtx data structure
549  * \ingroup dm
550  */
551 PetscErrorCode
552 DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx);
553 
554 /**
555  * \brief Set MoFEM::SnesCtx data structure
556  * \ingroup dm
557 
558  */
559 PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
560  boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx);
561 
562 /**
563  * \brief get MoFEM::TsCtx data structure
564  * \ingroup dm
565  */
566 PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx);
567 
568 /**
569  * \brief get MoFEM::TsCtx data structure
570  * \ingroup dm
571  */
572 PetscErrorCode DMMoFEMGetTsCtx(DM dm,
573  const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx);
574 
575 /**
576  * \brief Set MoFEM::TsCtx data structure
577  * \ingroup dm
578 
579  It take over pointer, do not delete it, DM will destroy pointer
580  when is destroyed.
581 
582  */
583 PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> &ts_ctx);
584 
585 /** sets if read mesh is partitioned
586  * \ingroup dm
587  */
588 PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned);
589 
590 /** get if read mesh is partitioned
591  * \ingroup dm
592  */
593 PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned);
594 
595 /**
596  * \brief Set operators for MoFEM dm
597  * @param dm
598  * @return error code
599  * \ingroup dm
600  */
601 PetscErrorCode DMSetOperators_MoFEM(DM dm);
602 
603 /**
604  * \brief Create dm data structure with MoFEM data structure
605  * \ingroup dm
606  */
607 PetscErrorCode DMCreate_MoFEM(DM dm);
608 
609 /**
610  * \brief Destroys dm with MoFEM data structure
611  * \ingroup dm
612  */
613 PetscErrorCode DMDestroy_MoFEM(DM dm);
614 
615 /**
616  * \brief DMShellSetCreateGlobalVector
617  * \ingroup dm
618  *
619  * sets the routine to create a global vector
620  * associated with the shell DM
621  */
622 PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g);
623 
624 /**
625  * \brief DMShellSetCreateGlobalVector
626  * \ingroup dm
627  *
628  * sets the routine to create a global vector
629  * associated with the shell DM
630  */
631 PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj<Vec> &g_ptr);
632 
633 /**
634  * \brief DMShellSetCreateLocalVector
635  * \ingroup dm
636  *
637  * sets the routine to create a local vector
638  * associated with the shell DM
639  */
640 PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l);
641 
642 /**
643  * DMShellSetCreateMatrix
644  * \ingroup dm
645  *
646  * sets the routine to create a matrix associated with the shell DM
647  */
648 PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M);
649 
650 /**
651  * DMShellSetCreateMatrix
652  * \ingroup dm
653  *
654  * sets the routine to create a matrix associated with the shell DM
655  */
656 PetscErrorCode DMCreateMatrix_MoFEM(DM dm, SmartPetscObj<Mat> &M);
657 
658 /**
659  * Set options for MoFEM DM
660  * \ingroup dm
661  */
662 #if PETSC_VERSION_GE(3, 7, 0)
663 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
664  DM dm);
665 #elif PETSC_VERSION_GE(3, 5, 3)
666 PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm);
667 #else
668 PetscErrorCode DMSetFromOptions_MoFEM(DM dm);
669 #endif
670 
671 /**
672  * sets up the MoFEM structures inside a DM object
673  * \ingroup dm
674  */
675 PetscErrorCode DMSetUp_MoFEM(DM dm);
676 
677 /**
678  * Sets up the MoFEM structures inside a DM object for sub dm
679  * \ingroup dm
680  */
681 PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm);
682 
683 /**
684  * Add field to sub dm problem on rows
685  * \ingroup dm
686  */
687 PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
688  EntityType lo_type = MBVERTEX,
689  EntityType hi_type = MBMAXTYPE);
690 
691 /**
692  * Add field to sub dm problem on columns
693  * \ingroup dm
694  */
695 PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
696  EntityType lo_type = MBVERTEX,
697  EntityType hi_type = MBMAXTYPE);
698 
699 /**
700  * Return true if this DM is sub problem
701  \ingroup dm
702  * @param dm the DM object
703  * @param is_subproblem true if subproblem
704  * @return error code
705  */
706 PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm);
707 
708 /**
709  * \brief get sub problem is
710  * @param dm has to be created with DMMoFEMSetSquareProblem
711  * @param is return is on the row
712  * @return error code
713  *
714  * Returns IS with global indices of the DM used to create SubDM
715  *
716  */
717 PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is);
718 
719 /**
720  * \brief get sub problem is
721  * @param dm has to be created with DMMoFEMSetSquareProblem
722  * @param is return is on the row
723  * @return error code
724  *
725  * Returns IS with global indices of the DM used to create SubDM
726  *
727  */
728 PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is);
729 
730 /**
731  * \brief Add problem to composite DM on row
732  \ingroup dm
733  *
734  * This create block on row with DOFs from problem of given name
735  *
736  * @param dm the DM object
737  * @param prb_name add problem name
738  * @return error code
739  */
740 PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]);
741 
742 /**
743  * \brief Add problem to composite DM on col
744  * \ingroup dm
745  *
746  * This create block on col with DOFs from problem of given name
747  *
748  * @param dm the DM object
749  * @param prb_name add problem name
750  * @return error code
751  */
752 PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]);
753 
754 /**
755  * \brief Get if this DM is composite DM
756  * \ingroup dm
757  *
758  * @param dm the DM object
759  * @param is_comp_dm return true if composite problem here
760  * @return error code
761  */
762 PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm);
763 
764 /**
765  * destroy the MoFEM structure
766  * \ingroup dm
767  */
768 PetscErrorCode DMDestroy_MoFEM(DM dm);
769 
770 /**
771  * DMShellSetGlobalToLocal
772  * \ingroup dm
773  *
774  * the routine that begins the global to local scatter
775  */
776 PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec);
777 
778 /**
779  * DMShellSetGlobalToLocal
780  * \ingroup dm
781  *
782  * the routine that begins the global to local scatter
783  */
784 PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec);
785 
786 /**
787  * DMShellSetLocalToGlobal
788  * \ingroup dm
789  *
790  * the routine that begins the local to global scatter
791  */
792 PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec);
793 
794 /**
795  * DMShellSetLocalToGlobal
796  * \ingroup dm
797  *
798  * the routine that ends the local to global scatter
799  */
800 PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec);
801 
802 /**
803  * DMShellSetLocalToGlobal
804  * \ingroup dm
805  *
806  * the routine that ends the local to global scatter
807  */
808 PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec);
809 
810 /**
811  * Creates a set of IS objects with the global indices of dofs for each field
812  * @param dm The number of fields (or NULL if not requested)
813  *
814  * Output:
815  * @param numFields The number of fields (or NULL if not requested)
816  * @param fieldNames The name for each field (or NULL if not requested)
817  * @param fields The global indices for each field (or NULL if not
818  requested)
819  *
820  * @return error code
821 
822  * \note The user is responsible for freeing all requested arrays. In
823  particular,
824  * every entry of names should be freed with PetscFree(), every entry of fields
825  * should be destroyed with ISDestroy(), and both arrays should be freed with
826  * PetscFree().
827 
828  \ingroup dm
829 
830  */
831 PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
832  char ***fieldNames, IS **fields);
833 
834 /**
835  * \brief get field is in the problem
836  * @param dm the DM objects
837  * @param rc ROW or COL (no difference is problem is squared)
838  * @param field_name name of the field
839  * @param is returned the IS object
840  * @return error code
841  *
842  * \code
843  * IS is;
844  * ierr = DMMoFEMGetFieldIS(dm,ROW,"DISP",&is_disp); CHKERRG(ierr);
845  * \endcode
846  *
847  *
848  \ingroup dm
849  */
850 PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
851  IS *is);
852 
853 /**
854  * @brief Set verbosity level
855  *
856  * @param dm
857  * @param verb see VERBOSITY_LEVELS for list of the levels
858  * @return PetscErrorCode
859  */
860 PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb);
861 
863 
864 /**
865  * \brief PETSc Discrete Manager data structure
866  *
867  * This structure should not be accessed or modified by user. Is not available
868  * from outside MoFEM DM manager. However user can inherit dat class and
869  * add data for additional functionality.
870  *
871  * This is part of implementation for PETSc interface, see more details in
872  * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DM/index.html>
873  *
874  * \ingroup dm
875  *
876  */
877 struct DMCtx : public UnknownInterface {
878 
880  UnknownInterface **iface) const;
881 
882  Interface *mField_ptr; ///< MoFEM interface
883  PetscBool isProblemBuild; ///< True if problem is build
884  std::string problemName; ///< Problem name
885 
886  // Options
887  PetscBool isPartitioned; ///< true if read mesh is on parts
888  PetscBool isSquareMatrix; ///< true if rows equals to cols
889 
890  int rAnk, sIze;
891 
892  // pointer to data structures
893  const Problem *problemPtr; ///< pinter to problem data structure
894 
895  // sub problem
896  PetscBool isSubDM;
897  std::vector<std::string> rowFields;
898  std::vector<std::string> colFields;
899  const Problem *problemMainOfSubPtr; ///< pinter to main problem to sub-problem
900 
901  PetscBool isCompDM;
902  std::vector<std::string> rowCompPrb;
903  std::vector<std::string> colCompPrb;
904  boost::shared_ptr<std::map<std::string, std::pair<EntityType, EntityType>>>
906  boost::shared_ptr<std::map<std::string, std::pair<EntityType, EntityType>>>
908 
909  PetscBool destroyProblem; ///< If true destroy problem with DM
910 
911  DMCtx();
912  virtual ~DMCtx();
913 
914  int verbosity; ///< verbosity
916 
917  boost::shared_ptr<KspCtx> kspCtx; ///< data structure KSP
918  boost::shared_ptr<SnesCtx> snesCtx; ///< data structure SNES
919  boost::shared_ptr<TsCtx> tsCtx; ///< data structure for TS solver
920 };
921 } // namespace MoFEM
922 
923 #endif //__DMMMOFEM_H
924 
925 /***************************************************************************/ /**
926  * \defgroup dm Distributed mesh manager
927  * \brief Implementation of PETSc DM, managing interactions between mesh data structures and vectors and matrices
928  *
929  * DM objects are used to manage communication between the algebraic structures in
930  * PETSc (Vec and Mat) and mesh data structures in PDE-based (or other)
931  * simulations.
932  *
933  * DM is abstract interface, here is it particular implementation for MoFEM code.
934  *
935  * \ingroup mofem
936  ******************************************************************************/
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:884
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:242
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:622
MoFEM interface unique ID.
PetscBool isSquareMatrix
true if rows equals to cols
Definition: DMMoFEM.hpp:888
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition: DMMMoFEM.cpp:1303
int verbosity
verbosity
Definition: DMMoFEM.hpp:914
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:541
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1174
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeRow
Definition: DMMoFEM.hpp:905
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:260
PetscBool isSubDM
Definition: DMMoFEM.hpp:896
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:430
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:367
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMMoFEM.cpp:233
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:210
PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:376
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:410
Interface * mField_ptr
MoFEM interface.
Definition: DMMoFEM.hpp:882
base class for all interface classes
PetscBool isPartitioned
true if read mesh is on parts
Definition: DMMoFEM.hpp:887
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1080
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:296
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: DMMMoFEM.cpp:41
const Problem * problemPtr
pinter to problem data structure
Definition: DMMoFEM.hpp:893
std::vector< std::string > colCompPrb
Definition: DMMoFEM.hpp:903
Data structure to exchange data between mofem and User Loop Methods on entities.It allows to exchange...
virtual ~DMCtx()
Definition: DMMMoFEM.cpp:36
keeps basic data about problemThis is low level structure with information about problem,...
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1211
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1014
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:582
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:757
RowColData
RowColData.
Definition: definitions.h:185
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:492
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition: DMMMoFEM.cpp:316
PetscBool destroyProblem
If true destroy problem with DM.
Definition: DMMoFEM.hpp:909
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition: DMMMoFEM.cpp:89
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM.
Definition: DMMMoFEM.cpp:987
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[], PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMMoFEM.cpp:386
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition: DMMoFEM.hpp:919
PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS the right hand side function
Definition: DMMMoFEM.cpp:796
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:166
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:278
std::vector< std::string > rowFields
Definition: DMMoFEM.hpp:897
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:187
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:979
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1127
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:400
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMMoFEM.cpp:914
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition: DMMMoFEM.cpp:442
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition: DMMoFEM.hpp:917
PetscBool isProblemBuild
True if problem is build.
Definition: DMMoFEM.hpp:883
int referenceNumber
Definition: DMMoFEM.hpp:915
Interface for nonlinear (SNES) solver.
Definition: SnesCtx.hpp:27
const Problem * problemMainOfSubPtr
pinter to main problem to sub-problem
Definition: DMMoFEM.hpp:899
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:53
PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS the right hand side jacobian
Definition: DMMMoFEM.cpp:836
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMMoFEM.cpp:358
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition: DMMoFEM.hpp:918
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:888
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:663
PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[])
unset element from dm
Definition: DMMMoFEM.cpp:420
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:453
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1182
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:463
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:109
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:81
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:59
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1246
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
Definition: DMMMoFEM.cpp:876
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:898
PetscBool isCompDM
Definition: DMMoFEM.hpp:901
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:902
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:962
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMMoFEM.cpp:1292
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:325
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:877
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:349
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:704
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:337
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMMoFEM.cpp:954
static const MOFEMuuid IDD_DMCTX
Definition: DMMoFEM.hpp:862
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMMoFEM.cpp:509
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:473
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeCol
Definition: DMMoFEM.hpp:907
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > &ksp_ctx)
set MoFEM::KspCtx data structure
Definition: DMMMoFEM.cpp:905
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:943
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1240
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1058
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM.
Definition: DMMMoFEM.cpp:1005
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > &snes_ctx)
Set MoFEM::SnesCtx data structure.
Definition: DMMMoFEM.cpp:931