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