v0.9.1
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 Duplicate internal data struture
43  *
44  * @param dm
45  * @param dm_duplicate
46  * @return PetscErrorCode
47  */
48 PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate);
49 
50 /**
51  * \brief Must be called by user to set Sub DM MoFEM data structures
52  * \ingroup dm
53  */
54 PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[]);
55 
56 /**
57  * \brief Get pointer to MoFEM::Interface
58  * @param dm Distributed mesh manager
59  * @param m_field_ptr Pointer to pointer of field interface
60  * @return Error code
61  * \ingroup dm
62  */
63 PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr);
64 
65 /**
66  * \brief Get pointer to problem data structure
67  * \ingroup dm
68  */
69 PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr);
70 
71 /**
72  * If this is set to PETSC_TRUE problem is deleted with DM
73  * @param dm the DM object
74  * @param destroy if PETSC_TRUE problem is destroyed
75  * @return error code
76  */
77 PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem);
78 
79 /**
80  * Get if problem will be destroyed with DM
81  * @param dm the DM object
82  * @param destroy return if PETSC_TRUE problem is destroyed
83  * @return error code
84  */
85 PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem);
86 
87 /**
88  * \brief set squared problem
89  * \ingroup dm
90 
91  It if true is assumed that matrix has the same indexing on rows and
92  columns. This reduces interprocessor communication.
93 
94  */
95 PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem);
96 
97 /**
98  * \brief get squared problem
99  * \ingroup dm
100 
101  It if true is assumed that matrix has the same indexing on rows and
102  columns. This reduces interprocessor communication.
103 
104  */
105 PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem);
106 
107 /**
108  * \brief Resolve shared entities
109  *
110  * @param dm dm
111  * @param fe_name finite element for which shared entities are resolved
112  * @return error code
113  *
114  * \note This function is valid for parallel algebra and serial mesh. It
115  * should be run collectively, i.e. on all processors.
116  *
117  * This allows for tag reduction or tag exchange, f.e.
118  *
119  * \code
120  * CHKERR DMMoFEMResolveSharedFiniteElements(dm,"SHELL_ELEMENT");
121  * Tag th;
122  * CHKERR mField.get_moab().tag_get_handle("ADAPT_ORDER",th);
123  * ParallelComm* pcomm =
124  * ParallelComm::get_pcomm(&mField.get_moab(),MYPCOMM_INDEX);
125  * // CHKERR pcomm->reduce_tags(th,MPI_SUM,prisms);
126  * CHKERR pcomm->exchange_tags(th,prisms);
127  * \endcode
128  *
129  * \ingroup dm
130  */
131 PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[]);
132 
133 /**
134  * @deprecated Use DMMoFEMResolveSharedFiniteElements
135  */
136 DEPRECATED PetscErrorCode DMMoFEMResolveSharedEntities(DM dm,
137  const char fe_name[]);
138 
139 /**
140  * \brief Get finite elements layout in the problem
141  *
142  * In layout is stored information how many elements is on each processor, for
143  * more information look int petsc documentation
144  * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/IS/PetscLayoutCreate.html#PetscLayoutCreate>
145  *
146  * @param dm discrete manager for this problem
147  * @param fe_name finite element name
148  * @param layout pointer to layout, for created layout user takes
149  * responsibility for destroying it.
150  * @return error code
151  *
152  * \ingroup dm
153  */
154 PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[],
155  PetscLayout *layout);
156 
157 /**
158  * \brief add element to dm
159  * \ingroup dm
160  *
161  * \note add_file is a collective, should be executed on all processors.
162  * Otherwise could lead to deadlock.
163  *
164  */
165 PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[]);
166 
167 /**
168  * \brief unset element from dm
169  * \ingroup dm
170  */
171 PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[]);
172 
173 /**
174  * \brief set local (or ghosted) vector values on mesh for partition only
175  * \ingroup dm
176 
177  * \param l vector
178  * \param mode see petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
179  * \param scatter_mode see petsc manual for ScatterMode (The available modes
180  are: SCATTER_FORWARD or SCATTER_REVERSE)
181  *
182  * SCATTER_REVERSE set data to field entities from V vector.
183  *
184  * SCATTER_FORWARD set vector V from data field entities
185 
186  */
187 PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
188  ScatterMode scatter_mode);
189 
190 /**
191  * \brief set ghosted vector values on all existing mesh entities
192  * \ingroup dm
193 
194  * \param g vector
195  * \param mode see petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
196  * \param scatter_mode see petsc manual for ScatterMode (The available modes
197  are: SCATTER_FORWARD or SCATTER_REVERSE)
198  *
199  * SCATTER_REVERSE set data to field entities from V vector.
200  *
201  * SCATTER_FORWARD set vector V from data field entities
202 
203  */
204 PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
205  ScatterMode scatter_mode);
206 
207 /**
208  * \brief execute finite element method for each element in dm (problem)
209  * \ingroup dm
210  */
211 PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method);
212 
213 /**
214  * \brief execute finite element method for each element in dm (problem)
215  * \ingroup dm
216  */
217 PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method);
218 
219 /**
220  * \brief Executes FEMethod for finite elements in DM
221  * @param dm MoFEM discrete manager
222  * @param fe_name name of finite element
223  * @param method pointer to \ref MoFEM::FEMethod
224  * @param low_rank lowest rank of processor
225  * @param up_rank upper run of processor
226  * @return Error code
227  * \ingroup dm
228  */
229 PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[],
230  MoFEM::FEMethod *method,
231  int low_rank, int up_rank);
232 
233 /**
234  * \brief Executes FEMethod for finite elements in DM
235  * @param dm MoFEM discrete manager
236  * @param fe_name name of finite element
237  * @param method pointer to \ref MoFEM::FEMethod
238  * @param low_rank lowest rank of processor
239  * @param up_rank upper run of processor
240  * @return Error code
241  * \ingroup dm
242  */
243 PetscErrorCode
244 DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const std::string fe_name,
245  boost::shared_ptr<MoFEM::FEMethod> method,
246  int low_rank, int up_rank);
247 
248 /**
249  * \brief Executes FEMethod for finite elements in DM
250  * @param dm MoFEM discrete manager
251  * @param fe_name name of element
252  * @param method pointer to \ref MOFEM::FEMethod
253  * @return Error code
254  * \ingroup dm
255  */
256 PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[],
257  MoFEM::FEMethod *method);
258 
259 /**
260  * \brief Executes FEMethod for finite elements in DM
261  * @param dm MoFEM discrete manager
262  * @param fe_name name of element
263  * @param method pointer to \ref MOFEM::FEMethod
264  * @return Error code
265  * \ingroup dm
266  */
267 PetscErrorCode
268 DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
269  boost::shared_ptr<MoFEM::FEMethod> method);
270 
271 /**
272  * \brief execute method for dofs on field in problem
273  * \ingroup dm
274  */
275 PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
276  MoFEM::DofMethod *method);
277 
278 // /**
279 // * \brief Set compute operator for KSP solver via sub-matrix and IS
280 // *
281 // * @param dm DM
282 // * @return error code
283 // *
284 // * \ingroup dm
285 // */
286 // PetscErrorCode DMMoFEMKSPSetComputeOperatorsViaSubMatrixbByIs(DM dm);
287 
288 /**
289  * \brief set KSP right hand side evaluation function
290  * \ingroup dm
291  */
292 PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
293  MoFEM::FEMethod *method,
294  MoFEM::BasicMethod *pre_only,
295  MoFEM::BasicMethod *post_only);
296 
297 /**
298  * \brief set KSP right hand side evaluation function
299  * \ingroup dm
300  */
301 PetscErrorCode
302 DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
303  boost::shared_ptr<MoFEM::FEMethod> method,
304  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
305  boost::shared_ptr<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 DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
321  MoFEM::FEMethod *method,
322  MoFEM::BasicMethod *pre_only,
323  MoFEM::BasicMethod *post_only);
324 
325 /**
326  * \brief Set KSP operators and push mofem finite element methods
327  *
328  * @param dm DM
329  * @param fe_name finite element name
330  * @param method method on the element (executed for each element in the
331  * problem which given name)
332  * @param pre_only method for pre-process before element method
333  * @param post_only method for post-process after element method
334  * @return error code
335  *
336  * \ingroup dm
337  */
338 PetscErrorCode
339 DMMoFEMKSPSetComputeOperators(DM dm, const std::string fe_name,
340  boost::shared_ptr<MoFEM::FEMethod> method,
341  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
342  boost::shared_ptr<MoFEM::BasicMethod> post_only);
343 
344 /**
345  * \brief set SNES residual evaluation function
346  * \ingroup dm
347  */
348 PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
349  MoFEM::FEMethod *method,
350  MoFEM::BasicMethod *pre_only,
351  MoFEM::BasicMethod *post_only);
352 
353 /**
354  * \brief set SNES residual evaluation function
355  * \ingroup dm
356  */
357 PetscErrorCode
358 DMMoFEMSNESSetFunction(DM dm, const std::string fe_name,
359  boost::shared_ptr<MoFEM::FEMethod> method,
360  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
361  boost::shared_ptr<MoFEM::BasicMethod> post_only);
362 
363 /**
364  * \brief set SNES Jacobian evaluation function
365  * \ingroup dm
366  */
367 PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
368  MoFEM::FEMethod *method,
369  MoFEM::BasicMethod *pre_only,
370  MoFEM::BasicMethod *post_only);
371 
372 /**
373  * \brief set SNES Jacobian evaluation function
374  * \ingroup dm
375  */
376 PetscErrorCode
377 DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
378  boost::shared_ptr<MoFEM::FEMethod> method,
379  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
380  boost::shared_ptr<MoFEM::BasicMethod> post_only);
381 
382 /**
383  * \brief set TS implicit function evaluation function
384  * \ingroup dm
385  */
386 PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
387  MoFEM::FEMethod *method,
388  MoFEM::BasicMethod *pre_only,
389  MoFEM::BasicMethod *post_only);
390 
391 /**
392  * \brief set TS implicit function evaluation function
393  * \ingroup dm
394  */
395 PetscErrorCode
396 DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
397  boost::shared_ptr<MoFEM::FEMethod> method,
398  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
399  boost::shared_ptr<MoFEM::BasicMethod> post_only);
400 
401 /**
402  * \brief set TS Jacobian evaluation function
403  * \ingroup dm
404  */
405 PetscErrorCode
406 DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
407  boost::shared_ptr<MoFEM::FEMethod> method,
408  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
409  boost::shared_ptr<MoFEM::BasicMethod> post_only);
410 
411 /**
412  * \brief set TS Jacobian evaluation function
413  * \ingroup dm
414  */
415 PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
416  MoFEM::FEMethod *method,
417  MoFEM::BasicMethod *pre_only,
418  MoFEM::BasicMethod *post_only);
419 
420 /**
421  * @brief set TS the right hand side function
422  *
423  * <a
424  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetRHSFunction.html#TSSetRHSFunction>See
425  * petsc documentation</a>
426  *
427  * @param dm
428  * @param fe_name
429  * @param method
430  * @param pre_only
431  * @param post_only
432  * @return PetscErrorCode
433  */
434 PetscErrorCode
435 DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name,
436  boost::shared_ptr<MoFEM::FEMethod> method,
437  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
438  boost::shared_ptr<MoFEM::BasicMethod> post_only);
439 
440 /**
441  * @brief set TS the right hand side jacobian
442  *
443  * <a
444  * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetRHSJacobian.html>See
445  * petsc documentation</a>
446  *
447  * @param dm
448  * @param fe_name
449  * @param method
450  * @param pre_only
451  * @param post_only
452  * @return PetscErrorCode
453  */
454 PetscErrorCode
455 DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
456  boost::shared_ptr<MoFEM::FEMethod> method,
457  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
458  boost::shared_ptr<MoFEM::BasicMethod> post_only);
459 
460 /**
461  * \brief set TS implicit function evaluation function
462  * \ingroup dm
463  */
464 PetscErrorCode
465 DMMoFEMTSSetI2Function(DM dm, const std::string fe_name,
466  boost::shared_ptr<MoFEM::FEMethod> method,
467  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
468  boost::shared_ptr<MoFEM::BasicMethod> post_only);
469 
470 /**
471  * \brief set TS Jacobian evaluation function
472  * \ingroup dm
473  */
474 PetscErrorCode
475 DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name,
476  boost::shared_ptr<MoFEM::FEMethod> method,
477  boost::shared_ptr<MoFEM::BasicMethod> pre_only,
478  boost::shared_ptr<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() = default;
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 
922 /**
923  * @brief Get smart matrix from DM
924  * \ingroup dm
925  *
926  */
927 auto smartCreateDMMatrix = [](DM dm) {
929  ierr = DMCreateMatrix_MoFEM(dm, a);
930  CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
931  return a;
932 };
933 
934 /**
935  * @brief Get smart vector from DM
936  * \ingroup dm
937  *
938  */
939 auto smartCreateDMVector = [](DM dm) {
942  CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
943  return v;
944 };
945 
946 /**
947  * @deprecated Use smartCreateDMVector
948  *
949  * @param dm
950  * @return DEPRECATED smartCreateDMVector
951  */
952 DEPRECATED inline auto smartCreateDMDVector(DM dm) {
953  return smartCreateDMVector(dm);
954 }
955 
956 } // namespace MoFEM
957 
958 #endif //__DMMMOFEM_H
959 
960 /**
961  * \defgroup dm Distributed mesh manager
962  * \brief Implementation of PETSc DM, managing interactions between mesh data
963  *structures and vectors and matrices
964  *
965  * DM objects are used to manage communication between the algebraic structures
966  *in PETSc (Vec and Mat) and mesh data structures in PDE-based (or other)
967  * simulations.
968  *
969  * DM is abstract interface, here is it particular implementation for MoFEM
970  *code.
971  *
972  * \ingroup mofem
973  **/
std::string problemName
Problem name.
Definition: DMMoFEM.hpp:884
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:253
structure for User Loop Methods on finite elementsIt can be used to calculate stiffness matrices,...
MPI_Comm getCommFromPetscObject(PetscObject obj)
Get the Comm From Petsc Object object.
Definition: AuxPETSc.hpp:189
PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES residual evaluation function
Definition: DMMMoFEM.cpp:637
Deprecated interface functions.
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:1342
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:556
DEPRECATED PetscErrorCode DMMoFEMResolveSharedEntities(DM dm, const char fe_name[])
Definition: DMMMoFEM.cpp:397
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1213
boost::shared_ptr< std::map< std::string, std::pair< EntityType, EntityType > > > mapTypeRow
Definition: DMMoFEM.hpp:905
auto smartCreateDMMatrix
Get smart matrix from DM.
Definition: DMMoFEM.hpp:927
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
Definition: DMMMoFEM.cpp:271
PetscBool isSubDM
Definition: DMMoFEM.hpp:896
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, const char fe_name[])
Resolve shared entities.
Definition: DMMMoFEM.cpp:387
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:445
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:378
Interface for Time Stepping (TS) solver.
Definition: TsCtx.hpp:27
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition: DMMMoFEM.cpp:244
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[], EntityType lo_type=MBVERTEX, EntityType hi_type=MBMAXTYPE)
Definition: DMMMoFEM.cpp:221
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:425
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:1119
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:307
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: DMMMoFEM.cpp:36
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
Definition: DMMMoFEM.cpp:160
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...
keeps basic data about problemThis is low level structure with information about problem,...
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1250
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition: DMMMoFEM.cpp:1053
PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set KSP operators and push mofem finite element methods.
Definition: DMMMoFEM.cpp:597
PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:772
RowColData
RowColData.
Definition: definitions.h:192
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:507
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:55
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:327
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:84
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM.
Definition: DMMMoFEM.cpp:1026
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, const char fe_name[], PetscLayout *layout)
Get finite elements layout in the problem.
Definition: DMMMoFEM.cpp:401
PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:844
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:801
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
Definition: DMMMoFEM.cpp:177
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:289
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:198
auto smartCreateDMVector
Get smart vector from DM.
Definition: DMMoFEM.hpp:939
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:1018
PetscErrorCode DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS Jacobian evaluation function
Definition: DMMMoFEM.cpp:866
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition: DMMMoFEM.cpp:1166
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:415
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition: DMMMoFEM.cpp:953
virtual ~DMCtx()=default
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:457
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:48
PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
set TS the right hand side jacobian
Definition: DMMMoFEM.cpp:831
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition: DMMMoFEM.cpp:369
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:927
PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set SNES Jacobian evaluation function
Definition: DMMMoFEM.cpp:678
PetscErrorCode DMMoFEMUnSetElement(DM dm, const char fe_name[])
unset element from dm
Definition: DMMMoFEM.cpp:435
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:468
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1221
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition: DMMMoFEM.cpp:478
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:105
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition: DMMMoFEM.cpp:76
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:54
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition: DMMMoFEM.cpp:1285
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
Definition: DMMMoFEM.cpp:915
std::vector< std::string > colFields
Definition: DMMoFEM.hpp:898
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:72
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:50
std::vector< std::string > rowCompPrb
Definition: DMMoFEM.hpp:902
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:1001
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition: DMMMoFEM.cpp:1331
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:336
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:877
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition: DMMMoFEM.cpp:360
PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
set TS implicit function evaluation function
Definition: DMMMoFEM.cpp:719
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:348
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition: DMMMoFEM.cpp:993
static const MOFEMuuid IDD_DMCTX
Definition: DMMoFEM.hpp:862
DEPRECATED auto smartCreateDMDVector(DM dm)
Definition: DMMoFEM.hpp:952
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition: DMMMoFEM.cpp:524
PetscErrorCode DMoFEMLoopFiniteElementsUpAndLowRank(DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank, int up_rank)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:488
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:944
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition: DMMMoFEM.cpp:982
Data structure to exchange data between mofem and User Loop Methods.It allows to exchange data betwee...
FTensor::Index< 'l', 2 > l
Definition: ContactOps.hpp:29
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition: DMMMoFEM.cpp:1279
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition: DMMMoFEM.cpp:1097
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVectorsets the routine to create a local vector associated with the shell DM.
Definition: DMMMoFEM.cpp:1044
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > &snes_ctx)
Set MoFEM::SnesCtx data structure.
Definition: DMMMoFEM.cpp:970