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