v0.15.0
Loading...
Searching...
No Matches
LoopMethods.hpp
Go to the documentation of this file.
1/** \file LoopMethods.hpp
2 * \brief MoFEM loop methods interface for finite element computations
3 * \author Anonymous author(s) committing under MIT license
4 *
5 * This file defines data structures and interfaces for making loops over finite
6 * elements, entities, and degrees of freedom (DOFs) in MoFEM problems and database.
7 *
8 * \section loop_overview Loop Methods Overview
9 *
10 * The core loop methods enable efficient iteration over:
11 * - **Finite Elements**: Using FEMethod with Interface::loop_finite_elements()
12 * - **Entities**: Using EntityMethod with Interface::loop_entities()
13 * - **DOFs**: Using DofMethod with Interface::loop_dofs()
14 *
15 * \section femethod_usage FEMethod Usage in Core Interface
16 *
17 * FEMethod is the primary interface for finite element computations, designed for
18 * high-performance element-level operations. It integrates with MoFEM's
19 * core looping mechanisms:
20 *
21 * **Loop Integration:**
22 * - Called via Interface::loop_finite_elements() for each element in a mesh
23 * - User overloads FEMethod::operator() to define element-specific computations
24 * - Provides access to element geometry, DOF data, and field values
25 * - Supports assembly of global matrices/vectors from element contributions
26 *
27 * **PETSc Solver Integration:**
28 * - **KSP (Linear)**: FEMethod assembles stiffness matrices and load vectors
29 * - **SNES (Nonlinear)**: Computes residuals and Jacobians for Newton methods
30 * - **TS (Time-stepping)**: Handles time-dependent mass matrices and residuals
31 * - **TAO (Optimization)**: Evaluates objective functions and gradients
32 *
33 * **ForcesAndSources Integration:**
34 * - FEMethod works with ForcesAndSources user data operators for element computations
35 * - ForcesAndSources provides high-level element assembly and integration routines
36 * - FEMethod offers direct low-level access for performance-critical applications
37 *
38 * Each solver context provides appropriate data structures (matrices, vectors,
39 * contexts) accessible through the inheritance hierarchy (BasicMethod -> FEMethod).
40 *
41 * \section performance_design Performance Design
42 *
43 * The loop methods are optimized for computational efficiency:
44 * - Direct memory access to DOF data and field values
45 * - Minimal overhead iteration over large meshes
46 * - Cache-friendly data access patterns
47 * - Integration with PETSc's high-performance linear algebra
48 *
49 * \see Interface::loop_finite_elements, Interface::loop_entities, Interface::loop_dofs
50 *
51 */
52
53#ifndef __LOOPMETHODS_HPP__
54#define __LOOPMETHODS_HPP__
55
56namespace MoFEM {
57
58/**
59 * \brief Base data structure for PETSc-related contexts
60 * \ingroup mofem_loops
61 *
62 * This structure provides a foundation for managing PETSc data objects and contexts
63 * used throughout MoFEM finite element computations. It handles vectors, matrices,
64 * and context information needed for various PETSc solvers (KSP, SNES, TS, TAO).
65 */
66struct PetscData : public UnknownInterface {
67
68 /**
69 * \brief Query interface for type casting
70 * \param type_index Type information for interface querying
71 * \param iface Pointer to interface object
72 * \return Error code
73 */
74 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
75 UnknownInterface **iface) const;
76
77 /**
78 * \brief Default constructor
79 */
80 PetscData();
81
82 /**
83 * \brief Virtual destructor
84 */
85 virtual ~PetscData() = default;
86
87 /**
88 * \brief Enumeration for data context flags
89 *
90 * These flags indicate which PETSc data structures are currently set
91 * and available for use in computations. Multiple flags can be combined
92 * using bitwise operations.
93 */
95 CTX_SET_NONE = 0, ///< No data set
96 CTX_SET_F = 1 << 0, ///< Residual vector F is set
97 CTX_SET_A = 1 << 1, ///< Jacobian matrix A is set
98 CTX_SET_B = 1 << 2, ///< Preconditioner matrix B is set
99 CTX_SET_X = 1 << 3, ///< Solution vector X is set
100 CTX_SET_DX = 1 << 4, ///< Solution increment DX is set
101 CTX_SET_X_T = 1 << 5, ///< Time derivative X_t is set
102 CTX_SET_X_TT = 1 << 6, ///< Second time derivative X_tt is set
103 CTX_SET_TIME = 1 << 7 ///< Time value is set
104 };
105
106 using Switches = std::bitset<8>; ///< Bitset type for context switches
107
108 static constexpr Switches CtxSetNone = PetscData::Switches(CTX_SET_NONE); ///< No data switch
109 static constexpr Switches CtxSetF = PetscData::Switches(CTX_SET_F); ///< Residual vector switch
110 static constexpr Switches CtxSetA = PetscData::Switches(CTX_SET_A); ///< Jacobian matrix switch
111 static constexpr Switches CtxSetB = PetscData::Switches(CTX_SET_B); ///< Preconditioner matrix switch
112 static constexpr Switches CtxSetX = PetscData::Switches(CTX_SET_X); ///< Solution vector switch
113 static constexpr Switches CtxSetDX = PetscData::Switches(CTX_SET_DX); ///< Solution increment switch
114 static constexpr Switches CtxSetX_T = PetscData::Switches(CTX_SET_X_T); ///< First time derivative switch
115 static constexpr Switches CtxSetX_TT = PetscData::Switches(CTX_SET_X_TT); ///< Second time derivative switch
116 static constexpr Switches CtxSetTime = PetscData::Switches(CTX_SET_TIME); ///< Time value switch
117
118 Switches data_ctx; ///< Current data context switches
119
120 /**
121 * \brief Copy PETSc data from another instance
122 * \param petsc_data Source PETSc data structure
123 * \return Error code
124 */
125 MoFEMErrorCode copyPetscData(const PetscData &petsc_data);
126
127 Vec f; ///< PETSc residual vector
128 Mat A; ///< PETSc Jacobian matrix
129 Mat B; ///< PETSc preconditioner matrix
130 Vec x; ///< PETSc solution vector
131 Vec dx; ///< PETSc solution increment vector
132 Vec x_t; ///< PETSc first time derivative vector
133 Vec x_tt; ///< PETSc second time derivative vector
134};
135
136/**
137 * \brief Data structure for KSP (linear solver) context
138 * \ingroup mofem_loops
139 *
140 * This structure extends PetscData to provide context and data management
141 * specifically for PETSc KSP (Krylov Subspace) linear solvers. It stores
142 * the KSP solver instance and maintains context information about the
143 * current computation phase.
144 */
145struct KspMethod : virtual public PetscData {
146
147 /**
148 * \brief Query interface for type casting
149 * \param type_index Type information for interface querying
150 * \param iface Pointer to interface object
151 * \return Error code
152 */
153 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
154 UnknownInterface **iface) const;
155
156 /**
157 * \brief Context enumeration for KSP solver phases
158 *
159 * Indicates the current phase of KSP computation, which determines
160 * which data structures are being used and what operations are valid.
161 */
163 CTX_SETFUNCTION, ///< Setting up the linear system function
164 CTX_OPERATORS, ///< Setting up linear operators
165 CTX_KSPNONE ///< No specific KSP context
166 };
167
168 /**
169 * \brief Default constructor
170 */
171 KspMethod();
172
173 /**
174 * \brief Virtual destructor
175 */
176 virtual ~KspMethod() = default;
177
178 /**
179 * \brief Copy data from another KSP method
180 * \param ksp Source KSP method
181 * \return Error code
182 */
184
185 KSPContext ksp_ctx; ///< Current KSP computation context
186 KSP ksp; ///< PETSc KSP linear solver object
187
188 Vec &ksp_f; ///< Reference to residual vector in KSP context
189 Mat &ksp_A; ///< Reference to system matrix in KSP context
190 Mat &ksp_B; ///< Reference to preconditioner matrix in KSP context
191};
192
193/**
194 * \brief Data structure for SNES (nonlinear solver) context
195 * \ingroup mofem_loops
196 *
197 * This structure extends PetscData to provide context and data management
198 * specifically for PETSc SNES (Scalable Nonlinear Equation Solvers). It stores
199 * the SNES solver instance and maintains context information about nonlinear
200 * solution phases including function evaluation and Jacobian computation.
201 */
202struct SnesMethod : virtual protected PetscData {
203
204 /**
205 * \brief Query interface for type casting
206 * \param type_index Type information for interface querying
207 * \param iface Pointer to interface object
208 * \return Error code
209 */
210 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
211 UnknownInterface **iface) const;
212
213 /**
214 * \brief Context enumeration for SNES solver phases
215 *
216 * Indicates the current phase of SNES computation, determining which
217 * operations are being performed and what data structures are active.
218 */
220 CTX_SNESSETFUNCTION, ///< Setting up nonlinear function evaluation
221 CTX_SNESSETJACOBIAN, ///< Setting up Jacobian matrix computation
222 CTX_SNESNONE ///< No specific SNES context
223 };
224
225 /**
226 * \brief Default constructor
227 */
228 SnesMethod();
229
230 /**
231 * \brief Virtual destructor
232 */
233 virtual ~SnesMethod() = default;
234
235 /**
236 * \brief Copy SNES data from another instance
237 * \param snes Source SNES method
238 * \return Error code
239 */
241
242 SNESContext snes_ctx; ///< Current SNES computation context
243
244 SNES snes; ///< PETSc SNES nonlinear solver object
245 Vec &snes_x; ///< Reference to current solution state vector
246 Vec &snes_dx; ///< Reference to solution update/increment vector
247 Vec &snes_f; ///< Reference to residual vector
248 Mat &snes_A; ///< Reference to Jacobian matrix
249 Mat &snes_B; ///< Reference to preconditioner of Jacobian matrix
250};
251
252/**
253 * \brief Data structure for TS (time stepping) context
254 * \ingroup mofem_loops
255 *
256 * This structure extends PetscData to provide context and data management
257 * specifically for PETSc TS (Time Stepping) solvers. It manages time-dependent
258 * computations including time derivatives, time step information, and temporal
259 * integration contexts for transient finite element analyses.
260 */
261struct TSMethod : virtual protected PetscData {
262
263 /**
264 * \brief Query interface for type casting
265 * \param type_index Type information for interface querying
266 * \param iface Pointer to interface object
267 * \return Error code
268 */
269 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
270 UnknownInterface **iface) const;
271
272 /**
273 * \brief Context enumeration for TS solver phases
274 *
275 * Indicates the current phase of time stepping computation, determining
276 * which temporal operations are being performed and what data is active.
277 */
279 CTX_TSSETRHSFUNCTION, ///< Setting up right-hand side function
280 CTX_TSSETRHSJACOBIAN, ///< Setting up right-hand side Jacobian
281 CTX_TSSETIFUNCTION, ///< Setting up implicit function
282 CTX_TSSETIJACOBIAN, ///< Setting up implicit Jacobian
283 CTX_TSTSMONITORSET, ///< Setting up time step monitoring
284 CTX_TSNONE ///< No specific TS context
285 };
286
287 /**
288 * \brief Default constructor
289 */
290 TSMethod();
291
292 /**
293 * \brief Virtual destructor
294 */
295 virtual ~TSMethod() = default;
296
297 /**
298 * \brief Copy TS solver data from another instance
299 * \param ts Source TS method
300 * \return Error code
301 */
303
304 TS ts; ///< PETSc time stepping solver object
305
306 TSContext ts_ctx; ///< Current TS computation context
307
308 PetscInt ts_step; ///< Current time step number
309 PetscReal ts_a; ///< Shift parameter for U_t (see PETSc Time Solver documentation)
310 PetscReal ts_aa; ///< Shift parameter for U_tt (second time derivative)
311 PetscReal ts_t; ///< Current time value
312 PetscReal ts_dt; ///< Current time step size
313
314 Vec &ts_u; ///< Reference to current state vector U(t)
315 Vec &ts_u_t; ///< Reference to first time derivative of state vector dU/dt
316 Vec &ts_u_tt; ///< Reference to second time derivative of state vector d²U/dt²
317 Vec &ts_F; ///< Reference to residual vector F(t,U,U_t,U_tt)
318
319 Mat &ts_A; ///< Reference to Jacobian matrix: dF/dU + a*dF/dU_t + aa*dF/dU_tt
320 Mat &ts_B; ///< Reference to preconditioner matrix for ts_A
321};
322
323/**
324 * \brief Data structure for TAO (optimization) context
325 * \ingroup mofem_loops
326 *
327 * This structure extends PetscData to provide context and data management
328 * specifically for PETSc TAO (Toolkit for Advanced Optimization) solvers.
329 * It manages optimization computations including objective function evaluation,
330 * gradient computation, and Hessian matrix assembly for constrained and
331 * unconstrained optimization problems.
332 */
333struct TaoMethod : virtual protected PetscData {
334
335 /**
336 * \brief Query interface for type casting
337 * \param type_index Type information for interface querying
338 * \param iface Pointer to interface object
339 * \return Error code
340 */
341 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
342 UnknownInterface **iface) const;
343
344 /**
345 * \brief Context enumeration for TAO solver phases
346 *
347 * Indicates the current phase of optimization computation, determining
348 * which optimization operations are being performed.
349 */
351 CTX_TAO_OBJECTIVE, ///< Evaluating objective function
352 CTX_TAO_GRADIENT, ///< Computing gradient vector
353 CTX_TAO_HESSIAN, ///< Assembling Hessian matrix
354 CTX_TAO_NONE ///< No specific TAO context
355 };
356
357 /**
358 * \brief Default constructor
359 */
360 TaoMethod();
361
362 /**
363 * \brief Virtual destructor
364 */
365 virtual ~TaoMethod() = default;
366
367 /**
368 * \brief Copy TAO data from another instance
369 * \param tao Source TAO method
370 * \return Error code
371 */
373
374 TAOContext tao_ctx; ///< Current TAO computation context
375
376 Tao tao; ///< PETSc TAO optimization solver object
377 Vec &tao_x; ///< Reference to optimization variables vector
378 Vec &tao_f; ///< Reference to gradient vector
379 Mat &tao_A; ///< Reference to Hessian matrix
380 Mat &tao_B; ///< Reference to preconditioner matrix for Hessian
381};
382
383/**
384 * \brief Data structure to exchange data between MoFEM and User Loop Methods
385 * \ingroup mofem_loops
386 *
387 * This is a comprehensive base class that combines all PETSc solver contexts
388 * (KSP, SNES, TS, TAO) and provides fundamental data exchange capabilities
389 * between MoFEM and user-defined functions. It stores information about
390 * multi-indices, problem databases, and provides access to finite element
391 * data structures for loop-based computations.
392 */
394
395 /**
396 * \brief Query interface for type casting
397 * \param type_index Type information for interface querying
398 * \param iface Pointer to interface object
399 * \return Error code
400 */
401 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
402 UnknownInterface **iface) const {
404 *iface = const_cast<BasicMethod *>(this);
406 }
407
408 /**
409 * \brief Default constructor
410 */
411 BasicMethod();
412
413 /**
414 * \brief Virtual destructor
415 */
416 virtual ~BasicMethod() = default;
417
418 /**
419 * @brief Current index of processed method in the loop
420 *
421 * This counter tracks which method/element is currently being processed
422 * in the loop iteration, useful for debugging and progress monitoring.
423 */
425
426 /**
427 * @brief Total number of methods to process in the loop
428 *
429 * This value represents the total count of items (elements, entities, etc.)
430 * that will be processed in the current loop execution.
431 */
433
434 /**
435 * \brief Get current loop iteration index
436 * \return Current position in the loop (0-based)
437 */
438 inline int getNinTheLoop() const { return nInTheLoop; }
439
440 /**
441 * \brief Get total loop size
442 * \return Total number of items to process in the loop
443 */
444 inline int getLoopSize() const { return loopSize; }
445
446 /**
447 * @brief Processor rank range for distributed finite element iteration
448 *
449 * This pair stores the low and high processor ranks that define the range
450 * of processors involved in the current finite element iteration. Used for
451 * parallel processing and load balancing across MPI processes.
452 */
453 std::pair<int, int> loHiFERank;
454
455 /**
456 * @brief Get processor rank range for finite element iteration
457 *
458 * Returns the pair containing the lowest and highest processor ranks
459 * involved in the current finite element loop iteration.
460 *
461 * @return Pair<low_rank, high_rank> of processor ranges
462 */
463 inline auto getLoHiFERank() const { return loHiFERank; }
464
465 /**
466 * @brief Get lower processor rank for finite element iteration
467 *
468 * Returns the lowest processor rank in the range of processors
469 * that will process finite elements in the current loop.
470 *
471 * @return Lowest processor rank in the iteration range
472 */
473 inline auto getLoFERank() const { return loHiFERank.first; }
474
475 /**
476 * @brief Get upper processor rank for finite element iteration
477 *
478 * Returns the highest processor rank in the range of processors
479 * that will process finite elements in the current loop.
480 *
481 * @return Highest processor rank in the iteration range
482 */
483 inline auto getHiFERank() const { return loHiFERank.second; }
484
485 int rAnk; ///< Current processor rank in MPI communicator
486
487 int sIze; ///< Total number of processors in MPI communicator
488
490 *refinedEntitiesPtr; ///< Pointer to container of refined MoFEM DOF entities
491
493 *refinedFiniteElementsPtr; ///< Pointer to container of refined finite element entities
494
495 const Problem *problemPtr; ///< Raw pointer to current MoFEM problem instance
496
497 const Field_multiIndex *fieldsPtr; ///< Raw pointer to fields multi-index container
498
500 *entitiesPtr; ///< Raw pointer to container of field entities
501
502 const DofEntity_multiIndex *dofsPtr; ///< Raw pointer to container of degree of freedom entities
503
505 *finiteElementsPtr; ///< Raw pointer to container of finite elements
506
508 *finiteElementsEntitiesPtr; ///< Raw pointer to container of finite element entities
509
511 *adjacenciesPtr; ///< Raw pointer to container of adjacencies between DOFs and finite elements
512
513 /**
514 * \brief Get bit number for a specific field by name
515 *
516 * This function looks up a field by name in the fields container and returns
517 * its bit number, which is used for field identification and bit operations
518 * in MoFEM's field management system.
519 *
520 * \param field_name Name of the field to look up
521 * \return Bit number of the field, or BITFEID_SIZE if field not found
522 * \throws MoFEM exception if fieldsPtr is not set
523 */
524 inline unsigned int getFieldBitNumber(std::string field_name) const {
525 if (fieldsPtr) {
526 auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
527 if (field_it != fieldsPtr->get<FieldName_mi_tag>().end())
528 return (*field_it)->getBitNumber();
529 else
530 return BITFEID_SIZE;
531 } else {
532 THROW_MESSAGE("Pointer to fields multi-index is not set");
533 return BITFEID_SIZE;
534 }
535 }
536
537 /**
538 * @brief Copy data from another BasicMethod instance
539 *
540 * This function copies all relevant data from another BasicMethod instance
541 * including solver contexts, database pointers, and configuration settings.
542 *
543 * @param basic Source BasicMethod to copy from
544 * @return Error code indicating success or failure
545 */
547
548 /**
549 * @brief Hook function for pre-processing operations
550 *
551 * User-defined function that is executed before the main loop begins.
552 * Can be used for initialization, setup operations, or custom preprocessing.
553 */
554 boost::function<MoFEMErrorCode()> preProcessHook;
555
556 /**
557 * @brief Hook function for post-processing operations
558 *
559 * User-defined function that is executed after the main loop completes.
560 * Can be used for cleanup, finalization, or custom post-processing.
561 */
562 boost::function<MoFEMErrorCode()> postProcessHook;
563
564 /**
565 * @brief Hook function for main operator execution
566 *
567 * User-defined function that can be executed during the main operator
568 * call, providing additional customization points in the computation.
569 */
570 boost::function<MoFEMErrorCode()> operatorHook;
571
572 /**
573 * \brief Pre-processing function executed at loop initialization
574 *
575 * This virtual function is called once at the beginning of the loop.
576 * It is typically used for:
577 * - Zeroing matrices and vectors
578 * - Computing shape functions on reference elements
579 * - Preprocessing boundary conditions
580 * - Setting up initial computation state
581 *
582 * \return Error code indicating success or failure
583 */
584 virtual MoFEMErrorCode preProcess();
585
586 /**
587 * \brief Main operator function executed for each loop iteration
588 *
589 * This virtual function is called for every item (finite element, entity, etc.)
590 * in the loop. It is the core computation function typically used for:
591 * - Calculating element local matrices and vectors
592 * - Assembling contributions to global system
593 * - Element-level post-processing operations
594 *
595 * \return Error code indicating success or failure
596 */
597 virtual MoFEMErrorCode operator()();
598
599 /**
600 * \brief Post-processing function executed at loop completion
601 *
602 * This virtual function is called once at the end of the loop.
603 * It is typically used for:
604 * - Assembling matrices and vectors to global system
605 * - Computing global variables (e.g., total internal energy)
606 * - Finalizing computation results
607 * - Cleanup operations
608 *
609 * Example of iterating over DOFs:
610 * \code
611 * for(_IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP_(this,"DISPLACEMENT",it)) {
612 * // Process each DOF
613 * }
614 * \endcode
615 *
616 * \return Error code indicating success or failure
617 */
618 virtual MoFEMErrorCode postProcess();
619
620 /**
621 * @brief Get the cache weak pointer object
622 *
623 * Returns a weak pointer to the cache tuple that stores problem information
624 * on entities about DOFs. Each problem stores different information.
625 *
626 * \warning When iterating over finite elements in a TS (time-stepping) solver
627 * preprocessor, use the TS cache in the loop. Using wrong cache will cause
628 * undefined behavior or segmentation faults. This design is a necessary
629 * compromise between bug resilience and memory/performance optimization.
630 *
631 * @return Weak pointer to the cache tuple
632 */
633 inline boost::weak_ptr<CacheTuple> getCacheWeakPtr() const {
634 return cacheWeakPtr;
635 }
636
637 boost::movelib::unique_ptr<bool> vecAssembleSwitch; ///< Switch for vector assembly operations
638 boost::movelib::unique_ptr<bool> matAssembleSwitch; ///< Switch for matrix assembly operations
639
640 boost::weak_ptr<CacheTuple> cacheWeakPtr; ///< Weak pointer to cached entity data
641};
642
643/**
644 * \brief Structure for user loop methods on finite elements
645 * \ingroup mofem_loops
646 *
647 * This class provides a high-performance interface for finite element computations.
648 * It can be used to calculate stiffness matrices, residual vectors, load vectors,
649 * and other element-level quantities. While it is a low-level class, users seeking
650 * maximum speed and efficiency can use it directly.
651 *
652 * This class is designed to work with Interface::loop_finite_elements, where the
653 * user-overloaded operator FEMethod::operator() is executed for each element in
654 * the problem. The class provides additional virtual methods (preProcess and
655 * postProcess) that can be overloaded by users for custom preprocessing and
656 * postprocessing operations.
657 *
658 * Key features:
659 * - Direct access to finite element data structures
660 * - High-performance element iteration
661 * - Customizable preprocessing and postprocessing
662 * - Integration with MoFEM's data management system
663 */
664struct FEMethod : public BasicMethod {
665
666 /**
667 * \brief Query interface for type casting
668 * \param type_index Type information for interface querying
669 * \param iface Pointer to interface object
670 * \return Error code
671 */
672 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
673 UnknownInterface **iface) const {
675 *iface = const_cast<FEMethod *>(this);
677 }
678
679 /**
680 * \brief Default constructor
681 */
682 FEMethod() = default;
683
684 std::string feName; ///< Name of the finite element being processed
685
686 boost::shared_ptr<const NumeredEntFiniteElement>
687 numeredEntFiniteElementPtr; ///< Shared pointer to finite element database structure
688
689 /**
690 * @brief Get the name of the current finite element
691 *
692 * Returns the name string identifier of the finite element currently
693 * being processed in the loop.
694 *
695 * @return String containing the finite element name
696 */
697 inline auto getFEName() const { return feName; }
698
699 /**
700 * @brief Test function to determine if element should be skipped
701 *
702 * This optional hook function allows users to implement custom logic
703 * to determine whether a particular finite element should be skipped
704 * during the loop execution. If this function is set and returns false,
705 * the element will be skipped in MoFEM::Core::loop_finite_elements.
706 *
707 * \note This functionality is particularly useful for running elements
708 * only on specific bit levels or under certain conditions.
709 *
710 * @param fe_method_ptr Pointer to the current FEMethod instance
711 * @return true if element should be processed, false if it should be skipped
712 */
713 boost::function<bool(FEMethod *fe_method_ptr)> exeTestHook;
714
715 /**
716 * \brief Get pointer to DOF data for the current finite element
717 * \return Pointer to DOF data container
718 */
719 inline auto getDataDofsPtr() const {
720 return numeredEntFiniteElementPtr->getDataDofsPtr();
721 };
722
723 /**
724 * \brief Get pointer to vector DOF data for the current finite element
725 * \return Pointer to vector DOF data container
726 */
727 inline auto getDataVectorDofsPtr() const {
728 return numeredEntFiniteElementPtr->getDataVectorDofsPtr();
729 };
730
731 /**
732 * \brief Get reference to data field entities for the current finite element
733 * \return Constant reference to field entities view
734 */
736 return numeredEntFiniteElementPtr->getDataFieldEnts();
737 }
738
739 /**
740 * \brief Get shared pointer to data field entities for the current finite element
741 * \return Shared pointer to field entities view (mutable access)
742 */
743 inline boost::shared_ptr<FieldEntity_vector_view> &
745 return const_cast<NumeredEntFiniteElement *>(
748 }
749
750 /**
751 * \brief Get reference to row field entities for the current finite element
752 * \return Reference to row field entities container
753 */
754 inline auto &getRowFieldEnts() const {
755 return numeredEntFiniteElementPtr->getRowFieldEnts();
756 };
757
758 /**
759 * \brief Get shared pointer to row field entities for the current finite element
760 * \return Shared pointer to row field entities container
761 */
762 inline auto &getRowFieldEntsPtr() const {
763 return numeredEntFiniteElementPtr->getRowFieldEntsPtr();
764 };
765
766 /**
767 * \brief Get reference to column field entities for the current finite element
768 * \return Reference to column field entities container
769 */
770 inline auto &getColFieldEnts() const {
771 return numeredEntFiniteElementPtr->getColFieldEnts();
772 };
773
774 /**
775 * \brief Get shared pointer to column field entities for the current finite element
776 * \return Shared pointer to column field entities container
777 */
778 inline auto &getColFieldEntsPtr() const {
779 return numeredEntFiniteElementPtr->getColFieldEntsPtr();
780 };
781
782 /**
783 * \brief Get pointer to row DOFs for the current finite element
784 *
785 * Returns a pointer to the container of row degrees of freedom
786 * associated with the current finite element. Row DOFs are used
787 * in matrix assembly for the row indices.
788 *
789 * \return Pointer to row DOFs container
790 */
791 inline auto getRowDofsPtr() const {
792 return numeredEntFiniteElementPtr->getRowDofsPtr();
793 };
794
795 /**
796 * \brief Get pointer to column DOFs for the current finite element
797 *
798 * Returns a pointer to the container of column degrees of freedom
799 * associated with the current finite element. Column DOFs are used
800 * in matrix assembly for the column indices.
801 *
802 * \return Pointer to column DOFs container
803 */
804 inline auto getColDofsPtr() const {
805 return numeredEntFiniteElementPtr->getColDofsPtr();
806 };
807
808 /**
809 * \brief Get the number of nodes in the current finite element
810 *
811 * Returns the number of vertex nodes for the current finite element
812 * based on its entity type (e.g., 4 for tetrahedron, 8 for hexahedron).
813 *
814 * \return Number of nodes in the finite element
815 */
816 inline auto getNumberOfNodes() const;
817
818 /**
819 * \brief Get the entity handle of the current finite element
820 *
821 * Returns the MOAB entity handle that uniquely identifies the
822 * current finite element in the mesh database.
823 *
824 * \return MOAB entity handle of the finite element
825 */
826 inline EntityHandle getFEEntityHandle() const;
827
828 /**
829 * \brief Get nodal data for a specific field
830 *
831 * Extracts nodal values for the specified field from the current
832 * finite element and stores them in the provided data vector.
833 *
834 * \param field_name Name of the field to extract data for
835 * \param data Vector to store the extracted nodal data
836 * \param reset_dofs Flag to reset DOF values (default: true)
837 * \return Error code indicating success or failure
838 */
839 MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data,
840 const bool reset_dofs = true);
841
842};
843
844inline auto FEMethod::getNumberOfNodes() const {
845 return moab::CN::VerticesPerEntity(numeredEntFiniteElementPtr->getEntType());
846};
847
849 return numeredEntFiniteElementPtr->getEnt();
850}
851
852/**
853 * \brief Data structure for user loop methods on entities
854 * \ingroup mofem_loops
855 *
856 * This structure extends BasicMethod to provide entity-level loop operations.
857 * It allows exchange of data between MoFEM and user functions when iterating
858 * over mesh entities (vertices, edges, faces, volumes). It stores pointers
859 * to field and entity information for entity-based computations.
860 *
861 * Key features:
862 * - Entity-level iteration capabilities
863 * - Field data access for entities
864 * - Integration with MoFEM's entity management system
865 */
866struct EntityMethod : public BasicMethod {
867
868 /**
869 * \brief Query interface for type casting
870 * \param type_index Type information for interface querying
871 * \param iface Pointer to interface object
872 * \return Error code
873 */
874 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
875 UnknownInterface **iface) const {
877 *iface = const_cast<EntityMethod *>(this);
879 }
880
881 /**
882 * \brief Default constructor
883 */
884 EntityMethod() = default;
885
886 boost::shared_ptr<Field> fieldPtr; ///< Shared pointer to field information
887 boost::shared_ptr<FieldEntity> entPtr; ///< Shared pointer to field entity data
888};
889
890/**
891 * \brief Data structure for user loop methods on degrees of freedom (DOFs)
892 * \ingroup mofem_loops
893 *
894 * This structure extends BasicMethod to provide DOF-level loop operations.
895 * It allows exchange of data between MoFEM and user functions when iterating
896 * over degrees of freedom. It stores pointers to field, DOF entity, and
897 * numbered DOF entity information for DOF-based computations.
898 *
899 * Key features:
900 * - DOF-level iteration capabilities
901 * - Access to both raw and numbered DOF data
902 * - Field information associated with DOFs
903 * - Integration with MoFEM's DOF management system
904 */
905struct DofMethod : public BasicMethod {
906
907 /**
908 * \brief Query interface for type casting
909 * \param type_index Type information for interface querying
910 * \param iface Pointer to interface object
911 * \return Error code
912 */
913 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
914 UnknownInterface **iface) const {
916 *iface = const_cast<DofMethod *>(this);
918 }
919
920 /**
921 * \brief Default constructor
922 */
923 DofMethod() = default;
924
925 boost::shared_ptr<Field> fieldPtr; ///< Shared pointer to field information
926 boost::shared_ptr<DofEntity> dofPtr; ///< Shared pointer to DOF entity data
927 boost::shared_ptr<NumeredDofEntity> dofNumeredPtr; ///< Shared pointer to numbered DOF entity data
928};
929
930/// \deprecated name changed use DofMethod instead EntMethod
932
933} // namespace MoFEM
934
935#endif // __LOOPMETHODS_HPP__
936
937/**
938 * \defgroup mofem_loops Loops
939 * \ingroup mofem
940 */
multi_index_container< FieldEntityEntFiniteElementAdjacencyMap, indexed_by< ordered_unique< tag< Composite_Unique_mi_tag >, composite_key< FieldEntityEntFiniteElementAdjacencyMap, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getEntUniqueId >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getFeUniqueId > > >, ordered_non_unique< tag< Unique_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getEntUniqueId > >, ordered_non_unique< tag< FE_Unique_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getFeUniqueId > >, ordered_non_unique< tag< FEEnt_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, EntityHandle, &FieldEntityEntFiniteElementAdjacencyMap::getFeHandle > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, EntityHandle, &FieldEntityEntFiniteElementAdjacencyMap::getEntHandle > > > > FieldEntityEntFiniteElementAdjacencyMap_multiIndex
MultiIndex container keeps Adjacencies Element and dof entities adjacencies and vice versa.
multi_index_container< boost::shared_ptr< Field >, indexed_by< hashed_unique< tag< BitFieldId_mi_tag >, const_mem_fun< Field, const BitFieldId &, &Field::getId >, HashBit< BitFieldId >, EqBit< BitFieldId > >, ordered_unique< tag< Meshset_mi_tag >, member< Field, EntityHandle, &Field::meshSet > >, ordered_unique< tag< FieldName_mi_tag >, const_mem_fun< Field, boost::string_ref, &Field::getNameRef > >, ordered_non_unique< tag< BitFieldId_space_mi_tag >, const_mem_fun< Field, FieldSpace, &Field::getSpace > > > > Field_multiIndex
Multi-index container for field storage and retrieval.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define DEPRECATED
Definition definitions.h:17
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define BITFEID_SIZE
max number of finite elements
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< DofEntity, UId, &DofEntity::getLocalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
multi_index_container< boost::shared_ptr< EntFiniteElement >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< EntFiniteElement, UId, &EntFiniteElement::getLocalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< EntFiniteElement::interface_type_RefEntity, EntityHandle, &EntFiniteElement::getEnt > > > > EntFiniteElement_multiIndex
MultiIndex container for EntFiniteElement.
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasVector< double > VectorDouble
Definition Types.hpp:68
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
DEPRECATED typedef DofMethod EntMethod
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getEnt > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
Field_multiIndex::index< FieldName_mi_tag >::type::iterator field_it
multi_index_container< boost::shared_ptr< RefElement >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getEnt > > > > RefElement_multiIndex
multi_index_container< boost::shared_ptr< FieldEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< FieldEntity, UId, &FieldEntity::localUId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntity::interface_type_RefEntity, EntityHandle, &FieldEntity::getEnt > > > > FieldEntity_multiIndex
constexpr auto field_name
Data structure to exchange data between MoFEM and User Loop Methods.
int loopSize
Total number of methods to process in the loop.
boost::movelib::unique_ptr< bool > matAssembleSwitch
Switch for matrix assembly operations.
const FieldEntity_multiIndex * entitiesPtr
Raw pointer to container of field entities.
BasicMethod()
Default constructor.
unsigned int getFieldBitNumber(std::string field_name) const
Get bit number for a specific field by name.
virtual ~BasicMethod()=default
Virtual destructor.
MoFEMErrorCode copyBasicMethod(const BasicMethod &basic)
Copy data from another BasicMethod instance.
auto getHiFERank() const
Get upper processor rank for finite element iteration.
const FiniteElement_multiIndex * finiteElementsPtr
Raw pointer to container of finite elements.
const RefElement_multiIndex * refinedFiniteElementsPtr
Pointer to container of refined finite element entities.
int getLoopSize() const
Get total loop size.
std::pair< int, int > loHiFERank
Processor rank range for distributed finite element iteration.
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing operations.
const DofEntity_multiIndex * dofsPtr
Raw pointer to container of degree of freedom entities.
int rAnk
Current processor rank in MPI communicator.
boost::weak_ptr< CacheTuple > getCacheWeakPtr() const
Get the cache weak pointer object.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
const Field_multiIndex * fieldsPtr
Raw pointer to fields multi-index container.
int nInTheLoop
Current index of processed method in the loop.
auto getLoHiFERank() const
Get processor rank range for finite element iteration.
const RefEntity_multiIndex * refinedEntitiesPtr
Pointer to container of refined MoFEM DOF entities.
boost::movelib::unique_ptr< bool > vecAssembleSwitch
Switch for vector assembly operations.
virtual MoFEMErrorCode postProcess()
Post-processing function executed at loop completion.
const FieldEntityEntFiniteElementAdjacencyMap_multiIndex * adjacenciesPtr
Raw pointer to container of adjacencies between DOFs and finite elements.
virtual MoFEMErrorCode operator()()
Main operator function executed for each loop iteration.
virtual MoFEMErrorCode preProcess()
Pre-processing function executed at loop initialization.
const EntFiniteElement_multiIndex * finiteElementsEntitiesPtr
Raw pointer to container of finite element entities.
const Problem * problemPtr
Raw pointer to current MoFEM problem instance.
auto getLoFERank() const
Get lower processor rank for finite element iteration.
boost::function< MoFEMErrorCode()> operatorHook
Hook function for main operator execution.
int sIze
Total number of processors in MPI communicator.
boost::weak_ptr< CacheTuple > cacheWeakPtr
Weak pointer to cached entity data.
int getNinTheLoop() const
Get current loop iteration index.
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing operations.
Data structure for user loop methods on degrees of freedom (DOFs)
boost::shared_ptr< DofEntity > dofPtr
Shared pointer to DOF entity data.
DofMethod()=default
Default constructor.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
boost::shared_ptr< NumeredDofEntity > dofNumeredPtr
Shared pointer to numbered DOF entity data.
boost::shared_ptr< Field > fieldPtr
Shared pointer to field information.
Data structure for user loop methods on entities.
EntityMethod()=default
Default constructor.
boost::shared_ptr< Field > fieldPtr
Shared pointer to field information.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
boost::shared_ptr< FieldEntity > entPtr
Shared pointer to field entity data.
Structure for user loop methods on finite elements.
std::string feName
Name of the finite element being processed.
auto & getColFieldEntsPtr() const
Get shared pointer to column field entities for the current finite element.
FEMethod()=default
Default constructor.
auto & getRowFieldEnts() const
Get reference to row field entities for the current finite element.
auto getFEName() const
Get the name of the current finite element.
MoFEMErrorCode getNodeData(const std::string field_name, VectorDouble &data, const bool reset_dofs=true)
Get nodal data for a specific field.
auto getDataDofsPtr() const
Get pointer to DOF data for the current finite element.
const FieldEntity_vector_view & getDataFieldEnts() const
Get reference to data field entities for the current finite element.
auto getDataVectorDofsPtr() const
Get pointer to vector DOF data for the current finite element.
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr() const
Get shared pointer to data field entities for the current finite element.
auto getColDofsPtr() const
Get pointer to column DOFs for the current finite element.
EntityHandle getFEEntityHandle() const
Get the entity handle of the current finite element.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
auto & getRowFieldEntsPtr() const
Get shared pointer to row field entities for the current finite element.
auto getNumberOfNodes() const
Get the number of nodes in the current finite element.
auto getRowDofsPtr() const
Get pointer to row DOFs for the current finite element.
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
Shared pointer to finite element database structure.
auto & getColFieldEnts() const
Get reference to column field entities for the current finite element.
boost::function< bool(FEMethod *fe_method_ptr)> exeTestHook
Test function to determine if element should be skipped.
MultiIndex Tag for field name.
Data structure for KSP (linear solver) context.
MoFEMErrorCode copyKsp(const KspMethod &ksp)
Copy data from another KSP method.
Mat & ksp_B
Reference to preconditioner matrix in KSP context.
KSPContext ksp_ctx
Current KSP computation context.
virtual ~KspMethod()=default
Virtual destructor.
KSPContext
Context enumeration for KSP solver phases.
@ CTX_KSPNONE
No specific KSP context.
@ CTX_SETFUNCTION
Setting up the linear system function.
@ CTX_OPERATORS
Setting up linear operators.
KspMethod()
Default constructor.
KSP ksp
PETSc KSP linear solver object.
Mat & ksp_A
Reference to system matrix in KSP context.
Vec & ksp_f
Reference to residual vector in KSP context.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
Partitioned (Indexed) Finite Element in Problem.
Base data structure for PETSc-related contexts.
static constexpr Switches CtxSetA
Jacobian matrix switch.
static constexpr Switches CtxSetX
Solution vector switch.
Vec f
PETSc residual vector.
static constexpr Switches CtxSetX_TT
Second time derivative switch.
static constexpr Switches CtxSetNone
No data switch.
static constexpr Switches CtxSetF
Residual vector switch.
Vec x_t
PETSc first time derivative vector.
std::bitset< 8 > Switches
Bitset type for context switches.
virtual ~PetscData()=default
Virtual destructor.
Vec x_tt
PETSc second time derivative vector.
Vec dx
PETSc solution increment vector.
Switches data_ctx
Current data context switches.
static constexpr Switches CtxSetDX
Solution increment switch.
static constexpr Switches CtxSetX_T
First time derivative switch.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
DataContext
Enumeration for data context flags.
@ CTX_SET_NONE
No data set.
@ CTX_SET_X_T
Time derivative X_t is set.
@ CTX_SET_A
Jacobian matrix A is set.
@ CTX_SET_TIME
Time value is set.
@ CTX_SET_DX
Solution increment DX is set.
@ CTX_SET_X
Solution vector X is set.
@ CTX_SET_F
Residual vector F is set.
@ CTX_SET_X_TT
Second time derivative X_tt is set.
@ CTX_SET_B
Preconditioner matrix B is set.
Vec x
PETSc solution vector.
Mat B
PETSc preconditioner matrix.
static constexpr Switches CtxSetB
Preconditioner matrix switch.
PetscData()
Default constructor.
static constexpr Switches CtxSetTime
Time value switch.
MoFEMErrorCode copyPetscData(const PetscData &petsc_data)
Copy PETSc data from another instance.
Mat A
PETSc Jacobian matrix.
keeps basic data about problem
Data structure for SNES (nonlinear solver) context.
SNESContext
Context enumeration for SNES solver phases.
@ CTX_SNESSETJACOBIAN
Setting up Jacobian matrix computation.
@ CTX_SNESSETFUNCTION
Setting up nonlinear function evaluation.
@ CTX_SNESNONE
No specific SNES context.
MoFEMErrorCode copySnes(const SnesMethod &snes)
Copy SNES data from another instance.
Vec & snes_f
Reference to residual vector.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
Vec & snes_x
Reference to current solution state vector.
SnesMethod()
Default constructor.
virtual ~SnesMethod()=default
Virtual destructor.
Vec & snes_dx
Reference to solution update/increment vector.
Mat & snes_B
Reference to preconditioner of Jacobian matrix.
SNESContext snes_ctx
Current SNES computation context.
Mat & snes_A
Reference to Jacobian matrix.
SNES snes
PETSc SNES nonlinear solver object.
Data structure for TS (time stepping) context.
TS ts
PETSc time stepping solver object.
TSContext
Context enumeration for TS solver phases.
@ CTX_TSSETIFUNCTION
Setting up implicit function.
@ CTX_TSSETRHSFUNCTION
Setting up right-hand side function.
@ CTX_TSSETRHSJACOBIAN
Setting up right-hand side Jacobian.
@ CTX_TSTSMONITORSET
Setting up time step monitoring.
@ CTX_TSSETIJACOBIAN
Setting up implicit Jacobian.
@ CTX_TSNONE
No specific TS context.
TSMethod()
Default constructor.
PetscReal ts_t
Current time value.
Mat & ts_A
Reference to Jacobian matrix: dF/dU + a*dF/dU_t + aa*dF/dU_tt.
Vec & ts_F
Reference to residual vector F(t,U,U_t,U_tt)
PetscReal ts_dt
Current time step size.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
Vec & ts_u_tt
Reference to second time derivative of state vector d²U/dt²
TSContext ts_ctx
Current TS computation context.
PetscReal ts_a
Shift parameter for U_t (see PETSc Time Solver documentation)
Vec & ts_u_t
Reference to first time derivative of state vector dU/dt.
PetscInt ts_step
Current time step number.
Vec & ts_u
Reference to current state vector U(t)
virtual ~TSMethod()=default
Virtual destructor.
PetscReal ts_aa
Shift parameter for U_tt (second time derivative)
MoFEMErrorCode copyTs(const TSMethod &ts)
Copy TS solver data from another instance.
Mat & ts_B
Reference to preconditioner matrix for ts_A.
Data structure for TAO (optimization) context.
Mat & tao_B
Reference to preconditioner matrix for Hessian.
TAOContext tao_ctx
Current TAO computation context.
virtual ~TaoMethod()=default
Virtual destructor.
Tao tao
PETSc TAO optimization solver object.
Mat & tao_A
Reference to Hessian matrix.
TaoMethod()
Default constructor.
TAOContext
Context enumeration for TAO solver phases.
@ CTX_TAO_HESSIAN
Assembling Hessian matrix.
@ CTX_TAO_NONE
No specific TAO context.
@ CTX_TAO_GRADIENT
Computing gradient vector.
@ CTX_TAO_OBJECTIVE
Evaluating objective function.
MoFEMErrorCode copyTao(const TaoMethod &tao)
Copy TAO data from another instance.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Query interface for type casting.
Vec & tao_x
Reference to optimization variables vector.
Vec & tao_f
Reference to gradient vector.
base class for all interface classes