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