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