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 * SCATTER_FORWARD set vector V from data field entities
190 *
191 * This code will update ghosted entities on all processors
192 * \code
193 * auto d_vector = createDMVector(dm);
194 * // fill d_vector on owned entities
195 CHKERR DMoFEMMeshToLocalVector(dm, d_vector, INSERT_VALUES,
196 SCATTER_FORWARD);
197 // update ghost values in d_vector
198 CHKERR VecGhostUpdateBegin(d_vector, INSERT_VALUES, SCATTER_FORWARD);
199 CHKERR VecGhostUpdateEnd(d_vector, INSERT_VALUES, SCATTER_FORWARD);
200 // now d_vector has updated ghost values, and following save values of
201 // vector to fields entities
202 CHKERR DMoFEMMeshToLocalVector(dm, d_vector, INSERT_VALUES,
203 SCATTER_REVERSE);
204 * \endcode
205 *
206 * \note you use DMoFEMMeshToGlobalVector as scatter reverse
207 *
208 * For implementations see VecManager::setLocalGhostVector.
209 */
210PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode,
211 ScatterMode scatter_mode);
212
213/**
214 * \brief set ghosted vector values on all existing mesh entities
215 * \ingroup dm
216
217 * \param g vector
218 * \param mode see petsc manual for VecSetValue (ADD_VALUES or INSERT_VALUES)
219 * \param scatter_mode see petsc manual for ScatterMode (The available modes
220 are: SCATTER_FORWARD or SCATTER_REVERSE)
221 *
222 * SCATTER_REVERSE set data to field entities from V vector.
223 *
224 * SCATTER_FORWARD set vector V from data field entities
225
226 */
227PetscErrorCode DMoFEMMeshToGlobalVector(DM dm, Vec g, InsertMode mode,
228 ScatterMode scatter_mode);
229
230/**
231 * \brief execute finite element method for each element in dm (problem)
232 * \ingroup dm
233 */
234PetscErrorCode DMoFEMPreProcessFiniteElements(DM dm, MoFEM::FEMethod *method);
235
236/**
237 * \brief execute finite element method for each element in dm (problem)
238 * \ingroup dm
239 */
240PetscErrorCode DMoFEMPostProcessFiniteElements(DM dm, MoFEM::FEMethod *method);
241
242/**
243 * \brief Executes FEMethod for finite elements in DM
244 * @param dm MoFEM discrete manager
245 * @param fe_name name of finite element
246 * @param method pointer to \ref MoFEM::FEMethod
247 * @param low_rank lowest rank of processor
248 * @param up_rank upper run of processor
249 * @return Error code
250 * \ingroup dm
251 */
253 DM dm, const char fe_name[], MoFEM::FEMethod *method, int low_rank,
254 int up_rank, CacheTupleWeakPtr cache_ptr = CacheTupleSharedPtr());
255
256/**
257 * \brief Executes FEMethod for finite elements in DM
258 * @param dm MoFEM discrete manager
259 * @param fe_name name of finite element
260 * @param method pointer to \ref MoFEM::FEMethod
261 * @param low_rank lowest rank of processor
262 * @param up_rank upper run of processor
263 * @return Error code
264 * \ingroup dm
265 */
267 DM dm, const std::string fe_name, boost::shared_ptr<MoFEM::FEMethod> method,
268 int low_rank, int up_rank,
270
271/**
272 * \brief Executes FEMethod for finite elements in DM
273 * @param dm MoFEM discrete manager
274 * @param fe_name name of element
275 * @param method pointer to \ref MOFEM::FEMethod
276 * @return Error code
277 * \ingroup dm
278 */
279PetscErrorCode
280DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method,
282
283/**
284 * \brief Executes FEMethod for finite elements in DM
285 * @param dm MoFEM discrete manager
286 * @param fe_name name of element
287 * @param method pointer to \ref MOFEM::FEMethod
288 * @return Error code
289 * \ingroup dm
290 */
291PetscErrorCode
292DMoFEMLoopFiniteElements(DM dm, const std::string fe_name,
293 boost::shared_ptr<MoFEM::FEMethod> method,
295
296/**
297 * \brief execute method for dofs on field in problem
298 * \ingroup dm
299 */
300PetscErrorCode DMoFEMLoopDofs(DM dm, const char field_name[],
301 MoFEM::DofMethod *method);
302
303// /**
304// * \brief Set compute operator for KSP solver via sub-matrix and IS
305// *
306// * @param dm DM
307// * @return error code
308// *
309// * \ingroup dm
310// */
311// PetscErrorCode DMMoFEMKSPSetComputeOperatorsViaSubMatrixbByIs(DM dm);
312
313/**
314 * \brief set KSP right hand side evaluation function
315 * \ingroup dm
316 */
317PetscErrorCode DMMoFEMKSPSetComputeRHS(DM dm, const char fe_name[],
318 MoFEM::FEMethod *method,
319 MoFEM::BasicMethod *pre_only,
320 MoFEM::BasicMethod *post_only);
321
322/**
323 * \brief set KSP right hand side evaluation function
324 * \ingroup dm
325 */
326PetscErrorCode
327DMMoFEMKSPSetComputeRHS(DM dm, const std::string fe_name,
328 boost::shared_ptr<MoFEM::FEMethod> method,
329 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
330 boost::shared_ptr<MoFEM::BasicMethod> post_only);
331
332/**
333 * \brief Set KSP operators and push mofem finite element methods
334 *
335 * @param dm DM
336 * @param fe_name finite element name
337 * @param method method on the element (executed for each element in the
338 * problem which given name)
339 * @param pre_only method for pre-process before element method
340 * @param post_only method for post-process after element method
341 * @return error code
342 *
343 * \ingroup dm
344 */
345PetscErrorCode DMMoFEMKSPSetComputeOperators(DM dm, const char fe_name[],
346 MoFEM::FEMethod *method,
347 MoFEM::BasicMethod *pre_only,
348 MoFEM::BasicMethod *post_only);
349
350/**
351 * \brief Set KSP operators and push mofem finite element methods
352 *
353 * @param dm DM
354 * @param fe_name finite element name
355 * @param method method on the element (executed for each element in the
356 * problem which given name)
357 * @param pre_only method for pre-process before element method
358 * @param post_only method for post-process after element method
359 * @return error code
360 *
361 * \ingroup dm
362 */
363PetscErrorCode
364DMMoFEMKSPSetComputeOperators(DM dm, const std::string fe_name,
365 boost::shared_ptr<MoFEM::FEMethod> method,
366 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
367 boost::shared_ptr<MoFEM::BasicMethod> post_only);
368
369/**
370 * \brief set SNES residual evaluation function
371 * \ingroup dm
372 */
373PetscErrorCode DMMoFEMSNESSetFunction(DM dm, const char fe_name[],
374 MoFEM::FEMethod *method,
375 MoFEM::BasicMethod *pre_only,
376 MoFEM::BasicMethod *post_only);
377
378/**
379 * \brief set SNES residual evaluation function
380 * \ingroup dm
381 */
382PetscErrorCode
383DMMoFEMSNESSetFunction(DM dm, const std::string fe_name,
384 boost::shared_ptr<MoFEM::FEMethod> method,
385 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
386 boost::shared_ptr<MoFEM::BasicMethod> post_only);
387
388/**
389 * \brief set SNES Jacobian evaluation function
390 * \ingroup dm
391 */
392PetscErrorCode DMMoFEMSNESSetJacobian(DM dm, const char fe_name[],
393 MoFEM::FEMethod *method,
394 MoFEM::BasicMethod *pre_only,
395 MoFEM::BasicMethod *post_only);
396
397/**
398 * \brief set SNES Jacobian evaluation function
399 * \ingroup dm
400 */
401PetscErrorCode
402DMMoFEMSNESSetJacobian(DM dm, const std::string fe_name,
403 boost::shared_ptr<MoFEM::FEMethod> method,
404 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
405 boost::shared_ptr<MoFEM::BasicMethod> post_only);
406
407/**
408 * \brief set TS implicit function evaluation function
409 * \ingroup dm
410 */
411PetscErrorCode DMMoFEMTSSetIFunction(DM dm, const char fe_name[],
412 MoFEM::FEMethod *method,
413 MoFEM::BasicMethod *pre_only,
414 MoFEM::BasicMethod *post_only);
415
416/**
417 * \brief set TS implicit function evaluation function
418 * \ingroup dm
419 */
420PetscErrorCode
421DMMoFEMTSSetIFunction(DM dm, const std::string fe_name,
422 boost::shared_ptr<MoFEM::FEMethod> method,
423 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
424 boost::shared_ptr<MoFEM::BasicMethod> post_only);
425
426/**
427 * \brief set TS Jacobian evaluation function
428 * \ingroup dm
429 */
430PetscErrorCode
431DMMoFEMTSSetIJacobian(DM dm, const std::string fe_name,
432 boost::shared_ptr<MoFEM::FEMethod> method,
433 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
434 boost::shared_ptr<MoFEM::BasicMethod> post_only);
435
436/**
437 * \brief set TS Jacobian evaluation function
438 * \ingroup dm
439 */
440PetscErrorCode DMMoFEMTSSetIJacobian(DM dm, const char fe_name[],
441 MoFEM::FEMethod *method,
442 MoFEM::BasicMethod *pre_only,
443 MoFEM::BasicMethod *post_only);
444
445/**
446 * @brief set TS the right hand side function
447 *
448 * <a
449 * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetRHSFunction.html#TSSetRHSFunction>See
450 * petsc documentation</a>
451 *
452 * @param dm
453 * @param fe_name
454 * @param method
455 * @param pre_only
456 * @param post_only
457 * @return PetscErrorCode
458 */
459PetscErrorCode
460DMMoFEMTSSetRHSFunction(DM dm, const std::string fe_name,
461 boost::shared_ptr<MoFEM::FEMethod> method,
462 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
463 boost::shared_ptr<MoFEM::BasicMethod> post_only);
464
465PetscErrorCode DMMoFEMTSSetRHSFunction(DM dm, const char fe_name[],
466 MoFEM::FEMethod *method,
467 MoFEM::BasicMethod *pre_only,
468 MoFEM::BasicMethod *post_only);
469
470/**
471 * @brief set TS the right hand side jacobian
472 *
473 * <a
474 * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSSetRHSJacobian.html>See
475 * petsc documentation</a>
476 *
477 * @param dm
478 * @param fe_name
479 * @param method
480 * @param pre_only
481 * @param post_only
482 * @return PetscErrorCode
483 */
484PetscErrorCode
485DMMoFEMTSSetRHSJacobian(DM dm, const std::string fe_name,
486 boost::shared_ptr<MoFEM::FEMethod> method,
487 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
488 boost::shared_ptr<MoFEM::BasicMethod> post_only);
489
490PetscErrorCode DMMoFEMTSSetRHSJacobian(DM dm, const char fe_name[],
491 MoFEM::FEMethod *method,
492 MoFEM::BasicMethod *pre_only,
493 MoFEM::BasicMethod *post_only);
494
495/**
496 * \brief set TS implicit function evaluation function
497 * \ingroup dm
498 */
499PetscErrorCode
500DMMoFEMTSSetI2Function(DM dm, const std::string fe_name,
501 boost::shared_ptr<MoFEM::FEMethod> method,
502 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
503 boost::shared_ptr<MoFEM::BasicMethod> post_only);
504/**
505 * \brief set TS implicit function evaluation function
506 * \ingroup dm
507 */
508PetscErrorCode DMMoFEMTSSetI2Function(DM dm, const char fe_name[],
509 MoFEM::FEMethod *method,
510 MoFEM::BasicMethod *pre_only,
511 MoFEM::BasicMethod *post_only);
512
513/**
514 * \brief set TS Jacobian evaluation function
515 * \ingroup dm
516 */
517PetscErrorCode
518DMMoFEMTSSetI2Jacobian(DM dm, const std::string fe_name,
519 boost::shared_ptr<MoFEM::FEMethod> method,
520 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
521 boost::shared_ptr<MoFEM::BasicMethod> post_only);
522/**
523 * \brief set TS Jacobian evaluation function
524 * \ingroup dm
525 */
526PetscErrorCode
527DMMoFEMTSSetI2Jacobian(DM dm, const char fe_name[],
528 MoFEM::FEMethod *method,
529 MoFEM::BasicMethod *pre_only,
530 MoFEM::BasicMethod *post_only);
531
532/**
533 * @brief Set Monitor To TS solver
534 *
535 * <a
536 * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSMonitorSet.html>See
537 * PETSc documentaton here</a>
538 *
539 * @param dm
540 * @param ts time solver
541 * @param fe_name
542 * @param method
543 * @param pre_only
544 * @param post_only
545 * @return PetscErrorCod
546 */
547PetscErrorCode
548DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name,
549 boost::shared_ptr<MoFEM::FEMethod> method,
550 boost::shared_ptr<MoFEM::BasicMethod> pre_only,
551 boost::shared_ptr<MoFEM::BasicMethod> post_only);
552
553/**
554 * @brief Set Monitor To TS solver
555 *
556 * <a
557 * href=https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/TSMonitorSet.html>See
558 * PETSc documentaton here</a>
559 *
560 * @param dm
561 * @param ts time solver
562 * @param fe_name
563 * @param method
564 * @param pre_only
565 * @param post_only
566 * @return PetscErrorCod
567 */
568PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const char fe_name[],
569 MoFEM::FEMethod *method,
570 MoFEM::BasicMethod *pre_only,
571 MoFEM::BasicMethod *post_only);
572
573/**
574 * \brief get MoFEM::KspCtx data structure
575 * \ingroup dm
576 */
577PetscErrorCode DMMoFEMGetKspCtx(DM dm, MoFEM::KspCtx **ksp_ctx);
578
579/**
580 * \brief get MoFEM::KspCtx data structure
581 * \ingroup dm
582 */
583PetscErrorCode
584DMMoFEMGetKspCtx(DM dm, const boost::shared_ptr<MoFEM::KspCtx> &ksp_ctx);
585
586/**
587 * \brief set MoFEM::KspCtx data structure
588 * \ingroup dm
589 */
590PetscErrorCode DMMoFEMSetKspCtx(DM dm,
591 boost::shared_ptr<MoFEM::KspCtx> ksp_ctx);
592
593/**
594 * \brief get MoFEM::SnesCtx data structure
595 * \ingroup dm
596 */
597PetscErrorCode DMMoFEMGetSnesCtx(DM dm, MoFEM::SnesCtx **snes_ctx);
598
599/**
600 * \brief get MoFEM::SnesCtx data structure
601 * \ingroup dm
602 */
603PetscErrorCode
604DMMoFEMGetSnesCtx(DM dm, const boost::shared_ptr<MoFEM::SnesCtx> &snes_ctx);
605
606/**
607 * \brief Set MoFEM::SnesCtx data structure
608 * \ingroup dm
609
610 */
611PetscErrorCode DMMoFEMSetSnesCtx(DM dm,
612 boost::shared_ptr<MoFEM::SnesCtx> snes_ctx);
613
614/**
615 * \brief get MoFEM::TsCtx data structure
616 * \ingroup dm
617 */
618PetscErrorCode DMMoFEMGetTsCtx(DM dm, MoFEM::TsCtx **ts_ctx);
619
620/**
621 * \brief get MoFEM::TsCtx data structure
622 * \ingroup dm
623 */
624PetscErrorCode DMMoFEMGetTsCtx(DM dm,
625 const boost::shared_ptr<MoFEM::TsCtx> &ts_ctx);
626
627/**
628 * \brief Set MoFEM::TsCtx data structure
629 * \ingroup dm
630
631 It take over pointer, do not delete it, DM will destroy pointer
632 when is destroyed.
633
634 */
635PetscErrorCode DMMoFEMSetTsCtx(DM dm, boost::shared_ptr<MoFEM::TsCtx> ts_ctx);
636
637/** sets if read mesh is partitioned
638 * \ingroup dm
639 */
640PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned);
641
642/** get if read mesh is partitioned
643 * \ingroup dm
644 */
645PetscErrorCode DMMoFEMGetIsPartitioned(DM dm, PetscBool *is_partitioned);
646
647/**
648 * \brief Set operators for MoFEM dm
649 * @param dm
650 * @return error code
651 * \ingroup dm
652 */
653PetscErrorCode DMSetOperators_MoFEM(DM dm);
654
655/**
656 * \brief Create dm data structure with MoFEM data structure
657 * \ingroup dm
658 */
659PetscErrorCode DMCreate_MoFEM(DM dm);
660
661/**
662 * \brief Destroys dm with MoFEM data structure
663 * \ingroup dm
664 */
665PetscErrorCode DMDestroy_MoFEM(DM dm);
666
667/**
668 * \brief DMShellSetCreateGlobalVector
669 * \ingroup dm
670 *
671 * sets the routine to create a global vector
672 * associated with the shell DM
673 */
674PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g);
675
676/**
677 * \brief DMShellSetCreateGlobalVector
678 * \ingroup dm
679 *
680 * sets the routine to create a global vector
681 * associated with the shell DM
682 */
683PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, SmartPetscObj<Vec> &g_ptr);
684
685/**
686 * \brief DMShellSetCreateLocalVector
687 * \ingroup dm
688 *
689 * sets the routine to create a local vector
690 * associated with the shell DM
691 */
692PetscErrorCode DMCreateLocalVector_MoFEM(DM dm, Vec *l);
693
694/**
695 * DMShellSetCreateMatrix
696 * \ingroup dm
697 *
698 * sets the routine to create a matrix associated with the shell DM
699 */
700PetscErrorCode DMCreateMatrix_MoFEM(DM dm, Mat *M);
701
702/**
703 * DMShellSetCreateMatrix
704 * \ingroup dm
705 *
706 * sets the routine to create a matrix associated with the shell DM
707 */
708PetscErrorCode DMCreateMatrix_MoFEM(DM dm, SmartPetscObj<Mat> &M);
709
710/**
711 * Set options for MoFEM DM
712 * \ingroup dm
713 */
714#if PETSC_VERSION_GE(3, 17, 0)
715PetscErrorCode DMSetFromOptions_MoFEM(DM dm,
716 PetscOptionItems *PetscOptionsObject);
717#elif PETSC_VERSION_GE(3, 7, 0)
718PetscErrorCode DMSetFromOptions_MoFEM(PetscOptionItems *PetscOptionsObject,
719 DM dm);
720#elif PETSC_VERSION_GE(3, 5, 3)
721PetscErrorCode DMSetFromOptions_MoFEM(PetscOptions *PetscOptionsObject, DM dm);
722#else
723PetscErrorCode DMSetFromOptions_MoFEM(DM dm);
724#endif
725
726/**
727 * sets up the MoFEM structures inside a DM object
728 * \ingroup dm
729 */
730PetscErrorCode DMSetUp_MoFEM(DM dm);
731
732/**
733 * Sets up the MoFEM structures inside a DM object for sub dm
734 * \ingroup dm
735 */
736PetscErrorCode DMSubDMSetUp_MoFEM(DM subdm);
737
738/**
739 * Add field to sub dm problem on rows
740 * \ingroup dm
741 */
742PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[]);
743
744/** @copydoc DMMoFEMAddSubFieldRow(DM dm, const char field_name[]) */
745PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, std::string field_name);
746
747/**
748 * @brief Add field to sub dm problem on rows
749 * \ingroup dm
750 *
751 * @param dm
752 * @param field_name
753 * @param m
754 * @return PetscErrorCode
755 */
756PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[],
757 boost::shared_ptr<Range> r_ptr);
758
759/** @copydoc DMMoFEMAddSubFieldRow(DM dm, const char field_name[], boost::shared_ptr<Range> r_ptr) */
760PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, std::string field_name,
761 boost::shared_ptr<Range> r_ptr);
762
763/**
764 * Add field to sub dm problem on columns
765 * \ingroup dm
766 */
767PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[]);
768
769/** @copydoc DMMoFEMAddSubFieldCol(DM dm, const char field_name[]) */
770PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, std::string field_name);
771
772/**
773 * @brief Add field to sub dm problem on columns
774 *
775 * @param dm
776 * @param field_name
777 * @param range of entities
778 * @return PetscErrorCode
779 */
780PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[],
781 boost::shared_ptr<Range> r_ptr);
782
783
784/** @copydoc DMMoFEMAddSubFieldCol(DM dm, const char field_name[], boost::shared_ptr<Range> r_ptr); */
785PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, std::string field_name,
786 boost::shared_ptr<Range> r_ptr);
787
788/**
789 * Return true if this DM is sub problem
790 \ingroup dm
791 * @param dm the DM object
792 * @param is_subproblem true if subproblem
793 * @return error code
794 */
795PetscErrorCode DMMoFEMGetIsSubDM(DM dm, PetscBool *is_sub_dm);
796
797/**
798 * \brief get sub problem is
799 * @param dm has to be created with DMMoFEMSetSquareProblem
800 * @param is return is on the row
801 * @return error code
802 *
803 * Returns IS with global indices of the DM used to create SubDM
804 *
805 */
806PetscErrorCode DMMoFEMGetSubRowIS(DM dm, IS *is);
807
808/**
809 * \brief get sub problem is
810 * @param dm has to be created with DMMoFEMSetSquareProblem
811 * @param is return is on the row
812 * @return error code
813 *
814 * Returns IS with global indices of the DM used to create SubDM
815 *
816 */
817PetscErrorCode DMMoFEMGetSubColIS(DM dm, IS *is);
818
819/**
820 * \brief Add problem to composite DM on row
821 \ingroup dm
822 *
823 * This create block on row with DOFs from problem of given name
824 *
825 * @param dm the DM object
826 * @param prb_name add problem name
827 * @return error code
828 */
829PetscErrorCode DMMoFEMAddRowCompositeProblem(DM dm, const char prb_name[]);
830
831/**
832 * \brief Add problem to composite DM on col
833 * \ingroup dm
834 *
835 * This create block on col with DOFs from problem of given name
836 *
837 * @param dm the DM object
838 * @param prb_name add problem name
839 * @return error code
840 */
841PetscErrorCode DMMoFEMAddColCompositeProblem(DM dm, const char prb_name[]);
842
843/**
844 * \brief Get if this DM is composite DM
845 * \ingroup dm
846 *
847 * @param dm the DM object
848 * @param is_comp_dm return true if composite problem here
849 * @return error code
850 */
851PetscErrorCode DMMoFEMGetIsCompDM(DM dm, PetscBool *is_comp_dm);
852
853/**
854 * destroy the MoFEM structure
855 * \ingroup dm
856 */
857PetscErrorCode DMDestroy_MoFEM(DM dm);
858
859/**
860 * DMShellSetGlobalToLocal
861 * \ingroup dm
862 *
863 * the routine that begins the global to local scatter
864 */
865PetscErrorCode DMGlobalToLocalBegin_MoFEM(DM dm, Vec, InsertMode, Vec);
866
867/**
868 * DMShellSetGlobalToLocal
869 * \ingroup dm
870 *
871 * the routine that begins the global to local scatter
872 */
873PetscErrorCode DMGlobalToLocalEnd_MoFEM(DM dm, Vec, InsertMode, Vec);
874
875/**
876 * DMShellSetLocalToGlobal
877 * \ingroup dm
878 *
879 * the routine that begins the local to global scatter
880 */
881PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec);
882
883/**
884 * DMShellSetLocalToGlobal
885 * \ingroup dm
886 *
887 * the routine that ends the local to global scatter
888 */
889PetscErrorCode DMLocalToGlobalBegin_MoFEM(DM, Vec, InsertMode, Vec);
890
891/**
892 * DMShellSetLocalToGlobal
893 * \ingroup dm
894 *
895 * the routine that ends the local to global scatter
896 */
897PetscErrorCode DMLocalToGlobalEnd_MoFEM(DM, Vec, InsertMode, Vec);
898
899/**
900 * Creates a set of IS objects with the global indices of dofs for each field
901 * @param dm The number of fields (or NULL if not requested)
902 *
903 * Output:
904 * @param numFields The number of fields (or NULL if not requested)
905 * @param fieldNames The name for each field (or NULL if not requested)
906 * @param fields The global indices for each field (or NULL if not
907 requested)
908 *
909 * @return error code
910
911 * \note The user is responsible for freeing all requested arrays. In
912 particular,
913 * every entry of names should be freed with PetscFree(), every entry of fields
914 * should be destroyed with ISDestroy(), and both arrays should be freed with
915 * PetscFree().
916
917 \ingroup dm
918
919 */
920PetscErrorCode DMCreateFieldIS_MoFEM(DM dm, PetscInt *numFields,
921 char ***fieldNames, IS **fields);
922
923/**
924 * \brief get field is in the problem
925 * @param dm the DM objects
926 * @param rc ROW or COL (no difference is problem is squared)
927 * @param field_name name of the field
928 * @param is returned the IS object
929 * @return error code
930 *
931 * \code
932 * IS is;
933 * ierr = DMMoFEMGetFieldIS(dm,ROW,"DISP",&is_disp); CHKERRG(ierr);
934 * \endcode
935 *
936 *
937 \ingroup dm
938 */
939PetscErrorCode DMMoFEMGetFieldIS(DM dm, RowColData rc, const char field_name[],
940 IS *is);
941
942/**
943 * @brief Set verbosity level
944 *
945 * @param dm
946 * @param verb see VERBOSITY_LEVELS for list of the levels
947 * @return PetscErrorCode
948 */
949PetscErrorCode DMMoFEMSetVerbosity(DM dm, const int verb);
950
951struct BlockStructure;
952
953/** @brief Set data for block mat
954 *
955 * \note You can reset data by setting nullptr
956 *
957 */
959 boost::shared_ptr<BlockStructure>);
960
961/**
962 * @brief Get data for block mat
963 *
964 * @param dm
965 * @return MoFEMErrorCode
966 */
967MoFEMErrorCode DMMoFEMGetBlocMatData(DM dm, boost::shared_ptr<BlockStructure> &);
968
969/**
970 * @brief Create block matrix
971 \ingroup dm
972 *
973 * @param dm
974 * @param mat
975 * @return MoFEMErrorCode
976 */
978
979/**
980 * @brief Create block matrix
981 \ingroup dm
982 *
983 * @param dm
984 * @param mat smart pointer
985 * @return MoFEMErrorCode
986 */
988
989/**
990 * @brief Set data for nest schur (see specialisation in Schur.hpp)
991 *
992 * \note You can reset data by setting nullptr
993 *
994 * @tparam T = NestSchurData
995 * @param dm
996 * @return MoFEMErrorCode
997 */
998template <typename T>
999MoFEMErrorCode DMMoFEMSetNestSchurData(DM dm, boost::shared_ptr<T>);
1000
1001/**
1002 * @brief Create nest schur matrix
1003 *
1004 * @param dm
1005 * @param mat
1006 * @return MoFEMErrorCode
1007 */
1009
1010/**
1011 * @brief Create matrix for hybridised system
1012 *
1013 * @param dm
1014 * @param mat
1015 * @return MoFEMErrorCode
1016 */
1018
1019/**
1020 * \brief PETSc Discrete Manager data structure
1021 *
1022 * This structure should not be accessed or modified by user. Is not
1023 * available from outside MoFEM DM manager. However user can inherit dat
1024 * class and add data for additional functionality.
1025 *
1026 * This is part of implementation for PETSc interface, see more details in
1027 * <http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DM/index.html>
1028 *
1029 * \ingroup dm
1030 *
1031 */
1032struct DMCtx : public UnknownInterface {
1033
1034 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
1035 UnknownInterface **iface) const;
1036
1037 virtual ~DMCtx() = default;
1038
1039 boost::shared_ptr<KspCtx> kspCtx; ///< data structure KSP
1040 boost::shared_ptr<SnesCtx> snesCtx; ///< data structure SNES
1041 boost::shared_ptr<TsCtx> tsCtx; ///< data structure for TS solver
1042
1043protected:
1044 DMCtx() = default;
1045};
1046
1047/**
1048 * @brief Get the Interface Ptr object
1049 *
1050 * @param dm
1051 * @return auto
1052 */
1053inline auto getInterfacePtr(DM dm) {
1054 MoFEM::Interface *m_field_ptr;
1055 CHK_THROW_MESSAGE(DMoFEMGetInterfacePtr(dm, &m_field_ptr),
1056 "Can not get interface ptr from DM");
1057 return m_field_ptr;
1058};
1059
1060/**
1061 * @brief get problem pointer from DM
1062 *
1063 */
1064inline auto getProblemPtr(DM dm) {
1065 const MoFEM::Problem *problem_ptr;
1066 CHK_THROW_MESSAGE(DMMoFEMGetProblemPtr(dm, &problem_ptr),
1067 "Get cot get problem ptr from DM");
1068 return problem_ptr;
1069};
1070
1071/**
1072 * @brief Get smart matrix from DM
1073 * \ingroup dm
1074 *
1075 */
1076inline auto createDMMatrix(DM dm) {
1079 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1080 return a;
1081};
1082
1083/**
1084 * @brief Get smart hybridised L2 matrix from DM
1085 *
1086 * @param dm
1087 * @return auto
1088 */
1089inline auto createDMHybridisedL2Matrix(DM dm) {
1092 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1093 return a;
1094};
1095
1096inline auto createDMBlockMat(DM dm) {
1097 Mat raw_a;
1098 ierr = DMMoFEMCreateBlockMat(dm, &raw_a);
1099 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1100 return SmartPetscObj<Mat>(raw_a);
1101};
1102
1103inline auto createDMNestSchurMat(DM dm) {
1104 Mat raw_a;
1105 ierr = DMMoFEMCreateNestSchurMat(dm, &raw_a);
1106 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1107 return SmartPetscObj<Mat>(raw_a);
1108};
1109
1110
1111/** @deprecated use createDMMatrix */
1112DEPRECATED inline auto smartCreateDMMatrix(DM dm) { return createDMMatrix(dm); }
1113
1114/**
1115 * @brief Get smart vector from DM
1116 * \ingroup dm
1117 *
1118 */
1119inline auto createDMVector(DM dm) {
1122 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1123 return v;
1124};
1125
1126/** @deprecated use createDMVector*/
1127DEPRECATED inline auto smartCreateDMVector(DM dm) { return createDMVector(dm); }
1128
1129/**
1130 * @brief Get KSP context data structure used by DM
1131 *
1132 */
1133inline auto getDMKspCtx(DM dm) {
1134 boost::shared_ptr<MoFEM::KspCtx> ksp_ctx;
1135 ierr = DMMoFEMGetKspCtx(dm, ksp_ctx);
1136 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1137 return ksp_ctx;
1138};
1139
1140/** @deprecated getDMKspCtx */
1141DEPRECATED inline auto smartGetDMKspCtx(DM dm) { return getDMKspCtx(dm); }
1142
1143/**
1144 * @brief Get SNES context data structure used by DM
1145 *
1146 */
1147inline auto getDMSnesCtx(DM dm) {
1148 boost::shared_ptr<MoFEM::SnesCtx> snes_ctx;
1149 ierr = DMMoFEMGetSnesCtx(dm, snes_ctx);
1150 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1151 return snes_ctx;
1152};
1153
1154/** @deprecated use smartGetDMSnesCtx*/
1155DEPRECATED inline auto smartGetDMSnesCtx(DM dm) { return getDMSnesCtx(dm); }
1156
1157/**
1158 * @brief Get TS context data structure used by DM
1159 *
1160 */
1161inline auto getDMTsCtx(DM dm) {
1162 boost::shared_ptr<MoFEM::TsCtx> ts_ctx;
1164 CHKERRABORT(getCommFromPetscObject(reinterpret_cast<PetscObject>(dm)), ierr);
1165 return ts_ctx;
1166};
1167
1168/** @deprecated use getDMTsCtx */
1169DEPRECATED inline auto smartGetDMTsCtx(DM dm) { return getDMTsCtx(dm); }
1170
1171/**
1172 * @brief Get sub problem data structure
1173 *
1174 * @param dm
1175 * @return auto
1176 */
1177inline auto getDMSubData(DM dm) {
1178 auto prb_ptr = getProblemPtr(dm);
1179 return prb_ptr->getSubData();
1180};
1181
1182
1183} // namespace MoFEM
1184
1185#endif //__DMMOFEM_H
1186
1187/**
1188 * \defgroup dm Distributed mesh manager
1189 * \brief Implementation of PETSc DM, managing interactions between mesh data
1190 *structures and vectors and matrices
1191 *
1192 * DM objects are used to manage communication between the algebraic structures
1193 *in PETSc (Vec and Mat) and mesh data structures in PDE-based (or other)
1194 * simulations.
1195 *
1196 * DM is abstract interface, here is it particular implementation for MoFEM
1197 *code.
1198 *
1199 * \ingroup mofem
1200 **/
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 KSP right hand side evaluation function
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:1119
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:1076
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:1112
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:1046
auto getDMTsCtx(DM dm)
Get TS context data structure used by DM.
Definition DMMoFEM.hpp:1161
DEPRECATED auto smartGetDMKspCtx(DM dm)
Definition DMMoFEM.hpp:1141
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:1127
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:1133
auto createDMHybridisedL2Matrix(DM dm)
Get smart hybridised L2 matrix from DM.
Definition DMMoFEM.hpp:1089
auto getDMSubData(DM dm)
Get sub problem data structure.
Definition DMMoFEM.hpp:1177
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:1169
auto getInterfacePtr(DM dm)
Get the Interface Ptr object.
Definition DMMoFEM.hpp:1053
MoFEMErrorCode DMMoFEMSetNestSchurData(DM dm, boost::shared_ptr< NestSchurData >)
Definition DMMoFEM.cpp:1554
DEPRECATED auto smartGetDMSnesCtx(DM dm)
Definition DMMoFEM.hpp:1155
auto getDMSnesCtx(DM dm)
Get SNES context data structure used by DM.
Definition DMMoFEM.hpp:1147
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:1103
PetscErrorCode DMMoFEMSwapDMCtx(DM dm, DM dm_swap)
Swap internal data struture.
Definition DMMoFEM.cpp:196
boost::shared_ptr< CacheTuple > CacheTupleSharedPtr
auto getProblemPtr(DM dm)
get problem pointer from DM
Definition DMMoFEM.hpp:1064
PetscErrorCode DMMoFEMDuplicateDMCtx(DM dm, DM dm_duplicate)
Duplicate internal data struture.
Definition DMMoFEM.cpp:180
boost::weak_ptr< CacheTuple > CacheTupleWeakPtr
auto createDMBlockMat(DM dm)
Definition DMMoFEM.hpp:1096
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:1032
virtual ~DMCtx()=default
boost::shared_ptr< SnesCtx > snesCtx
data structure SNES
Definition DMMoFEM.hpp:1040
boost::shared_ptr< TsCtx > tsCtx
data structure for TS solver
Definition DMMoFEM.hpp:1041
boost::shared_ptr< KspCtx > kspCtx
data structure KSP
Definition DMMoFEM.hpp:1039
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