v0.15.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, 17, 0)
698PetscErrorCode DMSetFromOptions_MoFEM(DM dm,
699 PetscOptionItems *PetscOptionsObject);
700#elif PETSC_VERSION_GE(3, 7, 0)
701PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
702 DM dm);
703#elif PETSC_VERSION_GE(3, 5, 3)
704PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm);
705#else
706PetscErrorCode DMSetFromOptions_MoFEM(DM dm);
707#endif
708
709/**
710 * sets up the MoFEM structures inside a DM object
711 * \ingroup dm
712 */
713PetscErrorCode DMSetUp_MoFEM(DM dm);
714
715/**
716 * Sets up the MoFEM structures inside a DM object for sub dm
717 * \ingroup dm
718 */
719PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm);
720
721/**
722 * Add field to sub dm problem on rows
723 * \ingroup dm
724 */
725PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[]);
726
727/** @copydoc DMMoFEMAddSubFieldRow(DM dm, const char field_name[]) */
728PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, std::string field_name);
729
730/**
731 * @brief Add field to sub dm problem on rows
732 * \ingroup dm
733 *
734 * @param dm
735 * @param field_name
736 * @param m
737 * @return PetscErrorCode
738 */
739PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
740 boost::shared_ptr<Range> r_ptr);
741
742/** @copydoc DMMoFEMAddSubFieldRow(DM dm, const char field_name[], boost::shared_ptr<Range> r_ptr) */
743PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, std::string field_name,
744 boost::shared_ptr<Range> r_ptr);
745
746/**
747 * Add field to sub dm problem on columns
748 * \ingroup dm
749 */
750PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[]);
751
752/** @copydoc DMMoFEMAddSubFieldCol(DM dm, const char field_name[]) */
753PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, std::string field_name);
754
755/**
756 * @brief Add field to sub dm problem on columns
757 *
758 * @param dm
759 * @param field_name
760 * @param range of entities
761 * @return PetscErrorCode
762 */
763PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
764 boost::shared_ptr<Range> r_ptr);
765
766
767/** @copydoc DMMoFEMAddSubFieldCol(DM dm, const char field_name[], boost::shared_ptr<Range> r_ptr); */
768PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, std::string field_name,
769 boost::shared_ptr<Range> r_ptr);
770
771/**
772 * Return true if this DM is sub problem
773 \ingroup dm
774 * @param dm the DM object
775 * @param is_subproblem true if subproblem
776 * @return error code
777 */
778PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm);
779
780/**
781 * \brief get sub problem is
782 * @param dm has to be created with DMMoFEMSetSquareProblem
783 * @param is return is on the row
784 * @return error code
785 *
786 * Returns IS with global indices of the DM used to create SubDM
787 *
788 */
789PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is);
790
791/**
792 * \brief get sub problem is
793 * @param dm has to be created with DMMoFEMSetSquareProblem
794 * @param is return is on the row
795 * @return error code
796 *
797 * Returns IS with global indices of the DM used to create SubDM
798 *
799 */
800PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is);
801
802/**
803 * \brief Add problem to composite DM on row
804 \ingroup dm
805 *
806 * This create block on row with DOFs from problem of given name
807 *
808 * @param dm the DM object
809 * @param prb_name add problem name
810 * @return error code
811 */
812PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]);
813
814/**
815 * \brief Add problem to composite DM on col
816 * \ingroup dm
817 *
818 * This create block on col with DOFs from problem of given name
819 *
820 * @param dm the DM object
821 * @param prb_name add problem name
822 * @return error code
823 */
824PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]);
825
826/**
827 * \brief Get if this DM is composite DM
828 * \ingroup dm
829 *
830 * @param dm the DM object
831 * @param is_comp_dm return true if composite problem here
832 * @return error code
833 */
834PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm);
835
836/**
837 * destroy the MoFEM structure
838 * \ingroup dm
839 */
840PetscErrorCode DMDestroy_MoFEM(DM dm);
841
842/**
843 * DMShellSetGlobalToLocal
844 * \ingroup dm
845 *
846 * the routine that begins the global to local scatter
847 */
848PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec);
849
850/**
851 * DMShellSetGlobalToLocal
852 * \ingroup dm
853 *
854 * the routine that begins the global to local scatter
855 */
856PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec);
857
858/**
859 * DMShellSetLocalToGlobal
860 * \ingroup dm
861 *
862 * the routine that begins the local to global scatter
863 */
864PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec);
865
866/**
867 * DMShellSetLocalToGlobal
868 * \ingroup dm
869 *
870 * the routine that ends the local to global scatter
871 */
872PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec);
873
874/**
875 * DMShellSetLocalToGlobal
876 * \ingroup dm
877 *
878 * the routine that ends the local to global scatter
879 */
880PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec);
881
882/**
883 * Creates a set of IS objects with the global indices of dofs for each field
884 * @param dm The number of fields (or NULL if not requested)
885 *
886 * Output:
887 * @param numFields The number of fields (or NULL if not requested)
888 * @param fieldNames The name for each field (or NULL if not requested)
889 * @param fields The global indices for each field (or NULL if not
890 requested)
891 *
892 * @return error code
893
894 * \note The user is responsible for freeing all requested arrays. In
895 particular,
896 * every entry of names should be freed with PetscFree(), every entry of fields
897 * should be destroyed with ISDestroy(), and both arrays should be freed with
898 * PetscFree().
899
900 \ingroup dm
901
902 */
903PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
904 char ***fieldNames, IS **fields);
905
906/**
907 * \brief get field is in the problem
908 * @param dm the DM objects
909 * @param rc ROW or COL (no difference is problem is squared)
910 * @param field_name name of the field
911 * @param is returned the IS object
912 * @return error code
913 *
914 * \code
915 * IS is;
916 * ierr = DMMoFEMGetFieldIS(dm,ROW,"DISP",&is_disp); CHKERRG(ierr);
917 * \endcode
918 *
919 *
920 \ingroup dm
921 */
922PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
923 IS *is);
924
925/**
926 * @brief Set verbosity level
927 *
928 * @param dm
929 * @param verb see VERBOSITY_LEVELS for list of the levels
930 * @return PetscErrorCode
931 */
932PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb);
933
934struct BlockStructure;
935
936/** @brief Set data for block mat
937 *
938 * \note You can reset data by setting nullptr
939 *
940 */
942 boost::shared_ptr<BlockStructure>);
943
944/**
945 * @brief Get data for block mat
946 *
947 * @param dm
948 * @return MoFEMErrorCode
949 */
950MoFEMErrorCode DMMoFEMGetBlocMatData(DM dm, boost::shared_ptr<BlockStructure> &);
951
952/**
953 * @brief Create block matrix
954 \ingroup dm
955 *
956 * @param dm
957 * @param mat
958 * @return MoFEMErrorCode
959 */
961
962/**
963 * @brief Create block matrix
964 \ingroup dm
965 *
966 * @param dm
967 * @param mat smart pointer
968 * @return MoFEMErrorCode
969 */
971
972/**
973 * @brief Set data for nest schur (see specialisation in Schur.hpp)
974 *
975 * \note You can reset data by setting nullptr
976 *
977 * @tparam T = NestSchurData
978 * @param dm
979 * @return MoFEMErrorCode
980 */
981template <typename T>
982MoFEMErrorCode DMMoFEMSetNestSchurData(DM dm, boost::shared_ptr<T>);
983
984/**
985 * @brief Create nest schur matrix
986 *
987 * @param dm
988 * @param mat
989 * @return MoFEMErrorCode
990 */
992
993/**
994 * @brief Create matrix for hybridised system
995 *
996 * @param dm
997 * @param mat
998 * @return MoFEMErrorCode
999 */
1001
1002/**
1003 * \brief PETSc Discrete Manager data structure
1004 *
1005 * This structure should not be accessed or modified by user. Is not
1006 * available from outside MoFEM DM manager. However user can inherit dat
1007 * class and add data for additional functionality.
1008 *
1009 * This is part of implementation for PETSc interface, see more details in
1010 * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DM/index.html>
1011 *
1012 * \ingroup dm
1013 *
1014 */
1015struct DMCtx : public UnknownInterface {
1016
1017 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
1018 UnknownInterface **iface) const;
1019
1020 virtual ~DMCtx() = default;
1021
1022 boost::shared_ptr<KspCtx> kspCtx; ///< data structure KSP
1023 boost::shared_ptr<SnesCtx> snesCtx; ///< data structure SNES
1024 boost::shared_ptr<TsCtx> tsCtx; ///< data structure for TS solver
1025
1026protected:
1027 DMCtx() = default;
1028};
1029
1030/**
1031 * @brief Get the Interface Ptr object
1032 *
1033 * @param dm
1034 * @return auto
1035 */
1036inline auto getInterfacePtr(DM dm) {
1037 MoFEM::Interface *m_field_ptr;
1038 CHK_THROW_MESSAGE(DMoFEMGetInterfacePtr(dm, &m_field_ptr),
1039 "Can not get interface ptr from DM");
1040 return m_field_ptr;
1041};
1042
1043/**
1044 * @brief get problem pointer from DM
1045 *
1046 */
1047inline auto getProblemPtr(DM dm) {
1048 const MoFEM::Problem *problem_ptr;
1049 CHK_THROW_MESSAGE(DMMoFEMGetProblemPtr(dm, &problem_ptr),
1050 "Get cot get problem ptr from DM");
1051 return problem_ptr;
1052};
1053
1054/**
1055 * @brief Get smart matrix from DM
1056 * \ingroup dm
1057 *
1058 */
1059inline auto createDMMatrix(DM dm) {
1062 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1063 return a;
1064};
1065
1066/**
1067 * @brief Get smart hybridised L2 matrix from DM
1068 *
1069 * @param dm
1070 * @return auto
1071 */
1072inline auto createDMHybridisedL2Matrix(DM dm) {
1075 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1076 return a;
1077};
1078
1079inline auto createDMBlockMat(DM dm) {
1080 Mat raw_a;
1081 ierr = DMMoFEMCreateBlockMat(dm, &raw_a);
1082 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1083 return SmartPetscObj<Mat>(raw_a);
1084};
1085
1086inline auto createDMNestSchurMat(DM dm) {
1087 Mat raw_a;
1088 ierr = DMMoFEMCreateNestSchurMat(dm, &raw_a);
1089 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1090 return SmartPetscObj<Mat>(raw_a);
1091};
1092
1093
1094/** @deprecated use createDMMatrix */
1095DEPRECATED inline auto smartCreateDMMatrix(DM dm) { return createDMMatrix(dm); }
1096
1097/**
1098 * @brief Get smart vector from DM
1099 * \ingroup dm
1100 *
1101 */
1102inline auto createDMVector(DM dm) {
1105 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1106 return v;
1107};
1108
1109/** @deprecated use createDMVector*/
1110DEPRECATED inline auto smartCreateDMVector(DM dm) { return createDMVector(dm); }
1111
1112/**
1113 * @brief Get KSP context data structure used by DM
1114 *
1115 */
1116inline auto getDMKspCtx(DM dm) {
1117 boost::shared_ptr<MoFEM::KspCtx> ksp_ctx;
1118 ierr = DMMoFEMGetKspCtx(dm, ksp_ctx);
1119 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1120 return ksp_ctx;
1121};
1122
1123/** @deprecated getDMKspCtx */
1124DEPRECATED inline auto smartGetDMKspCtx(DM dm) { return getDMKspCtx(dm); }
1125
1126/**
1127 * @brief Get SNES context data structure used by DM
1128 *
1129 */
1130inline auto getDMSnesCtx(DM dm) {
1131 boost::shared_ptr<MoFEM::SnesCtx> snes_ctx;
1132 ierr = DMMoFEMGetSnesCtx(dm, snes_ctx);
1133 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1134 return snes_ctx;
1135};
1136
1137/** @deprecated use smartGetDMSnesCtx*/
1138DEPRECATED inline auto smartGetDMSnesCtx(DM dm) { return getDMSnesCtx(dm); }
1139
1140/**
1141 * @brief Get TS context data structure used by DM
1142 *
1143 */
1144inline auto getDMTsCtx(DM dm) {
1145 boost::shared_ptr<MoFEM::TsCtx> ts_ctx;
1147 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1148 return ts_ctx;
1149};
1150
1151/** @deprecated use getDMTsCtx */
1152DEPRECATED inline auto smartGetDMTsCtx(DM dm) { return getDMTsCtx(dm); }
1153
1154/**
1155 * @brief Get sub problem data structure
1156 *
1157 * @param dm
1158 * @return auto
1159 */
1160inline auto getDMSubData(DM dm) {
1161 auto prb_ptr = getProblemPtr(dm);
1162 return prb_ptr->getSubData();
1163};
1164
1165
1166} // namespace MoFEM
1167
1168#endif //__DMMOFEM_H
1169
1170/**
1171 * \defgroup dm Distributed mesh manager
1172 * \brief Implementation of PETSc DM, managing interactions between mesh data
1173 *structures and vectors and matrices
1174 *
1175 * DM objects are used to manage communication between the algebraic structures
1176 *in PETSc (Vec and Mat) and mesh data structures in PDE-based (or other)
1177 * simulations.
1178 *
1179 * DM is abstract interface, here is it particular implementation for MoFEM
1180 *code.
1181 *
1182 * \ingroup mofem
1183 **/
constexpr double a
RowColData
RowColData.
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define DEPRECATED
Definition definitions.h:17
PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition DMMoFEM.cpp:1391
PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec)
Definition DMMoFEM.cpp:1399
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition DMMoFEM.cpp:1113
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:215
PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned)
Definition DMMoFEM.cpp:1124
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition DMMoFEM.cpp:488
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:708
PetscErrorCode DMMoFEMGetProblemFiniteElementLayout(DM dm, std::string fe_name, PetscLayout *layout)
Get finite elements layout in the problem.
Definition DMMoFEM.cpp:467
PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l)
DMShellSetCreateLocalVector.
Definition DMMoFEM.cpp:1177
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition DMMoFEM.cpp:450
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:790
PetscErrorCode DMDestroy_MoFEM(DM dm)
Destroys dm with MoFEM data structure.
Definition DMMoFEM.cpp:79
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:114
PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition DMMoFEM.cpp:546
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:514
PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on col.
Definition DMMoFEM.cpp:382
PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec)
Definition DMMoFEM.cpp:1427
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition DMMoFEM.cpp:422
PetscErrorCode DMCreate_MoFEM(DM dm)
Create dm data structure with MoFEM data structure.
Definition DMMoFEM.cpp:71
PetscErrorCode DMMoFEMSetSnesCtx(DM dm, boost::shared_ptr< MoFEM::SnesCtx > snes_ctx)
Set MoFEM::SnesCtx data structure.
Definition DMMoFEM.cpp:1101
PetscErrorCode DMMoFEMSetKspCtx(DM dm, boost::shared_ptr< MoFEM::KspCtx > ksp_ctx)
set MoFEM::KspCtx data structure
Definition DMMoFEM.cpp:1075
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
Definition DMMoFEM.cpp:238
PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm)
Definition DMMoFEM.cpp:322
PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[])
Add problem to composite DM on row.
Definition DMMoFEM.cpp:364
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:749
PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M)
Definition DMMoFEM.cpp:1187
PetscErrorCode DMSetFromOptions_MoFEM(DM dm)
Definition DMMoFEM.cpp:1275
PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx)
get MoFEM::TsCtx data structure
Definition DMMoFEM.cpp:1132
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition DMMoFEM.cpp:43
PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx)
get MoFEM::KspCtx data structure
Definition DMMoFEM.cpp:1058
PetscErrorCode DMMoFEMGetSquareProblem(DM dm, PetscBool *square_problem)
get squared problem
Definition DMMoFEM.cpp:480
PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[], MoFEM::DofMethod *method)
execute method for dofs on field in problem
Definition DMMoFEM.cpp:595
PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[], MoFEM::FEMethod *method, MoFEM::BasicMethod *pre_only, MoFEM::BasicMethod *post_only)
Set compute operator for KSP solver via sub-matrix and IS.
Definition DMMoFEM.cpp:627
PetscErrorCode DMMoFEMResolveSharedFiniteElements(DM dm, std::string fe_name)
Resolve shared entities.
Definition DMMoFEM.cpp:458
PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr< MoFEM::TsCtx > ts_ctx)
Set MoFEM::TsCtx data structure.
Definition DMMoFEM.cpp:1149
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:576
PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[], IS *is)
get field is in the problem
Definition DMMoFEM.cpp:1506
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:843
auto createDMVector(DM dm)
Get smart vector from DM.
Definition DMMoFEM.hpp:1102
MoFEMErrorCode DMMoFEMCreateBlockMat(DM dm, Mat *mat)
Create block matrix.
Definition DMMoFEM.cpp:1543
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition DMMoFEM.cpp:410
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
Definition DMMoFEM.cpp:280
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:1007
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:965
auto createDMMatrix(DM dm)
Get smart matrix from DM.
Definition DMMoFEM.hpp:1059
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition DMMoFEM.cpp:1157
PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx)
get MoFEM::SnesCtx data structure
Definition DMMoFEM.cpp:1084
PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode, ScatterMode scatter_mode)
set ghosted vector values on all existing mesh entities
Definition DMMoFEM.cpp:525
PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm)
Get if this DM is composite DM.
Definition DMMoFEM.cpp:402
PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm)
Definition DMMoFEM.cpp:1344
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:668
PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
Definition DMMoFEM.cpp:1461
PetscErrorCode DMMoFEMUnSetElement(DM dm, std::string fe_name)
unset element from dm
Definition DMMoFEM.cpp:505
PetscErrorCode DMSetUp_MoFEM(DM dm)
Definition DMMoFEM.cpp:1297
PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec)
Definition DMMoFEM.cpp:1456
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:557
PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method)
execute finite element method for each element in dm (problem)
Definition DMMoFEM.cpp:536
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition DMMoFEM.cpp:49
const double v
phase velocity of light in medium (cm/ns)
MoFEM::TsCtx * ts_ctx
FTensor::Index< 'l', 3 > l
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
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:1095
MPI_Comm getCommFromPetscObject(PetscObject obj)
Get the Comm From Petsc Object object.
boost::weak_ptr< CacheTuple > CacheTupleWeakPtr
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:1046
auto getDMTsCtx(DM dm)
Get TS context data structure used by DM.
Definition DMMoFEM.hpp:1144
DEPRECATED auto smartGetDMKspCtx(DM dm)
Definition DMMoFEM.hpp:1124
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:872
PetscErrorCode DMMoFEMSetDestroyProblem(DM dm, PetscBool destroy_problem)
Definition DMMoFEM.cpp:434
MoFEMErrorCode DMMoFEMCreateHybridL2Mat(DM dm, SmartPetscObj< Mat > &mat)
Create matrix for hybridised system.
Definition DMMoFEM.cpp:1574
PetscErrorCode DMMoFEMGetDestroyProblem(DM dm, PetscBool *destroy_problem)
Definition DMMoFEM.cpp:442
MoFEMErrorCode DMMoFEMGetBlocMatData(DM dm, boost::shared_ptr< BlockStructure > &)
Get data for block mat.
Definition DMMoFEM.cpp:1534
PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is)
get sub problem is
Definition DMMoFEM.cpp:330
DEPRECATED auto smartCreateDMVector(DM dm)
Definition DMMoFEM.hpp:1110
MoFEMErrorCode DMMoFEMCreateNestSchurMat(DM dm, Mat *mat)
Create nest schur matrix.
Definition DMMoFEM.cpp:1564
PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is)
get sub problem is
Definition DMMoFEM.cpp:347
auto getDMKspCtx(DM dm)
Get KSP context data structure used by DM.
Definition DMMoFEM.hpp:1116
auto createDMHybridisedL2Matrix(DM dm)
Get smart hybridised L2 matrix from DM.
Definition DMMoFEM.hpp:1072
auto getDMSubData(DM dm)
Get sub problem data structure.
Definition DMMoFEM.hpp:1160
MoFEMErrorCode DMMoFEMSetBlocMatData(DM dm, boost::shared_ptr< BlockStructure >)
Set data for block mat.
Definition DMMoFEM.cpp:1525
DEPRECATED auto smartGetDMTsCtx(DM dm)
Definition DMMoFEM.hpp:1152
boost::shared_ptr< CacheTuple > CacheTupleSharedPtr
auto getInterfacePtr(DM dm)
Get the Interface Ptr object.
Definition DMMoFEM.hpp:1036
MoFEMErrorCode DMMoFEMSetNestSchurData(DM dm, boost::shared_ptr< NestSchurData >)
Definition DMMoFEM.cpp:1554
DEPRECATED auto smartGetDMSnesCtx(DM dm)
Definition DMMoFEM.hpp:1138
auto getDMSnesCtx(DM dm)
Get SNES context data structure used by DM.
Definition DMMoFEM.hpp:1130
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:912
PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb)
Set verbosity level.
Definition DMMoFEM.cpp:1517
auto createDMNestSchurMat(DM dm)
Definition DMMoFEM.hpp:1086
PetscErrorCode DMMoFEMSwapDMCtx(DM dm, DM dm_swap)
Swap internal data struture.
Definition DMMoFEM.cpp:196
auto getProblemPtr(DM dm)
get problem pointer from DM
Definition DMMoFEM.hpp:1047
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
Definition DMMoFEM.cpp:180
auto createDMBlockMat(DM dm)
Definition DMMoFEM.hpp:1079
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:1015
virtual ~DMCtx()=default
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition DMMoFEM.hpp:1023
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition DMMoFEM.hpp:1024
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition DMMoFEM.hpp:1022
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition DMMoFEM.cpp:19
DMCtx()=default
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:14
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:17
base class for all interface classes