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