v0.15.0
Loading...
Searching...
No Matches
BiLinearFormsIntegrators.hpp
Go to the documentation of this file.
1/** \file BiLinearFormsIntegrators.hpp
2 * \brief Bilinear forms integrators for finite element assembly
3 * \ingroup mofem_form
4 *
5 * This file provides template aliases for bilinear form integrators used in finite
6 * element assembly. Bilinear forms represent the system matrices (left-hand side) in
7 * finite element formulations through the mathematical expression \f$(u, v)_\Omega\f$,
8 * where \f$u\f$ and \f$v\f$ are trial and test functions respectively.
9 *
10 * The integrators encompass fundamental finite element operations:
11 * - **Gradient operators**: \f$(v_{,i}, u_{,j})_\Omega\f$ for diffusion and elasticity
12 * - **Mass matrices**: \f$(v, u)_\Omega\f$ for dynamic and parabolic problems
13 * - **Mixed formulations**: Coupling different field types for constrained problems
14 * - **Tensor operations**: Anisotropic material behavior with full tensor coupling
15 * - **Hybridization**: Interface constraints for discontinuous Galerkin methods
16 *
17 * All template aliases reference concrete implementations in BiLinearFormsIntegratorsImpl.hpp,
18 * providing a clean interface while maintaining computational efficiency. The operators
19 * support various coordinate systems, field dimensions, and specializations for different
20 * problem types in computational mechanics and multi-physics simulations.
21
22 \todo Some operators could be optimized. To do that, we need to write tests
23 and use Valgrind to profile code, checking cache misses. For example, some
24 operators should have iteration first over columns, then rows for optimal
25 memory access patterns. Since these operators are used in many problems,
26 implementation efficiency is critical for overall performance.
27
28*/
29
30#ifndef __BILINEAR_FORMS_INTEGRATORS_HPP__
31#define __BILINEAR_FORMS_INTEGRATORS_HPP__
32
33namespace MoFEM {
34
35/**
36 * @brief Bilinear integrator form
37 * @ingroup mofem_forms
38 *
39 * @tparam EleOp
40 * @tparam A
41 * @tparam I
42 */
43template <typename EleOp>
44template <AssemblyType A>
45template <IntegrationType I>
47
48 /**
49 * @brief Gradient-gradient bilinear form: \f$(v_{,i},\beta(\mathbf{x}) u_{,j})_\Omega\f$
50 * @ingroup mofem_forms
51 *
52 * This operator integrates the product of test function gradients with trial
53 * function gradients, scaled by a coefficient \f$\beta(\mathbf{x})\f$. This is
54 * fundamental for diffusion, heat conduction, and elasticity problems. Applications include:
55 * - Laplacian operators in diffusion equations
56 * - Stiffness matrices in structural mechanics
57 * - Gradient-based regularization terms
58 *
59 * @note \f$\beta(\mathbf{x})\f$ is a spatially varying coefficient
60 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
61 *
62 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
63 * @tparam FIELD_DIM Dimension of the field
64 * @tparam SPACE_DIM Spatial dimension of the problem
65 */
66 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
68
69 /**
70 * @brief Mass matrix bilinear form: \f$(v_i,\beta(\mathbf{x}) u_j)_\Omega\f$
71 * @ingroup mofem_forms
72 *
73 * This operator integrates the product of test and trial functions, scaled by
74 * a coefficient \f$\beta(\mathbf{x})\f$. This forms the mass matrix in dynamic
75 * problems and identity-type operators in various formulations. Applications include:
76 * - Mass matrices in dynamic finite element analysis
77 * - Identity operators in parabolic equations
78 * - Penalty terms and stabilization methods
79 * - Time-dependent term discretization
80 *
81 * @note FIELD_DIM = 4 or 9 assumes OpMass is for tensorial field 2x2 or 3x3
82 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
83 *
84 * @todo Consider another template parameter RANK to distinguish between
85 * scalar, vectorial and tensorial fields
86 *
87 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
88 * @tparam FIELD_DIM Dimension of the field
89 */
90 template <int BASE_DIM, int FIELD_DIM>
92
93 /**
94 * @brief Cached mass matrix bilinear form (same as OpMass with caching)
95 * @ingroup mofem_forms
96 *
97 * This is a cached version of OpMass that stores intermediate computations
98 * to improve performance in repeated evaluations. Useful when the same
99 * mass matrix structure is evaluated multiple times.
100 *
101 * @copydetails OpMass
102 */
103 template <int BASE_DIM, int FIELD_DIM>
105
106 /**
107 * @brief Symmetric tensor gradient bilinear form: \f$(v_k,D_{ijkl} u_{,l})_\Omega\f$
108 * @ingroup mofem_forms
109 *
110 * This operator integrates test functions with trial function gradients through
111 * a symmetric fourth-order tensor \f$D_{ijkl}\f$ with minor and major symmetry.
112 * Essential for anisotropic diffusion and linear elasticity. Applications include:
113 * - Elasticity stiffness matrices with material symmetry
114 * - Anisotropic diffusion operators
115 * - Constitutive relations in continuum mechanics
116 * - Stress-strain coupling in solid mechanics
117 *
118 * @note \f$D_{ijkl}\f$ is a tensor with minor & major symmetry
119 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
120 *
121 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
122 * @tparam FIELD_DIM Dimension of the field
123 * @tparam SPACE_DIM Spatial dimension of the problem
124 * @tparam S Additional parameter for specialization (default: 0)
125 */
126 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
129
130 /**
131 * @brief Second-order gradient symmetric tensor bilinear form: \f$(v_{,ij},D_{ijkl} u_{,kl})_\Omega\f$
132 * @ingroup mofem_forms
133 *
134 * This operator integrates second-order gradients (Hessians) of test and trial
135 * functions through a symmetric fourth-order tensor \f$D_{ijkl}\f$. Used in
136 * higher-order partial differential equations and plate/shell theories. Applications include:
137 * - Biharmonic operators in plate bending problems
138 * - Fourth-order diffusion equations
139 * - Higher-order regularization in optimization
140 * - Gradient damage models in mechanics
141 *
142 * @note \f$D_{ijkl}\f$ is a tensor with no symmetries
143 * @note Requires C1 continuous finite elements
144 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
145 *
146 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
147 * @tparam FIELD_DIM Dimension of the field
148 * @tparam SPACE_DIM Spatial dimension of the problem
149 * @tparam S Additional parameter for specialization (default: 0)
150 */
151 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
154 OpBase>;
155
156 /**
157 * @brief General tensor gradient bilinear form: \f$(v_{,j},D_{ijkl} u_{,l})_\Omega\f$
158 * @ingroup mofem_forms
159 *
160 * This operator integrates test function gradients with trial function gradients
161 * through a general fourth-order tensor \f$D_{ijkl}\f$ with no symmetry assumptions.
162 * This is the most general form for anisotropic material behavior. Applications include:
163 * - General anisotropic elasticity without symmetry
164 * - Non-symmetric material tensors in advanced materials
165 * - Coupled multi-physics problems with complex constitutive relations
166 * - Plasticity with kinematic hardening
167 *
168 * @note \f$D_{ijkl}\f$ is a tensor with no symmetries
169 * @note More computationally expensive than symmetric versions
170 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
171 *
172 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
173 * @tparam FIELD_DIM Dimension of the field
174 * @tparam SPACE_DIM Spatial dimension of the problem
175 * @tparam S Additional parameter for specialization (default: 0)
176 */
177 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 0>
180
181 /**
182 * @brief Mixed divergence-scalar bilinear form: \f$(\lambda_{i,i},u)_\Omega\f$
183 * @ingroup mofem_forms
184 *
185 * This operator integrates the divergence of vector test functions with scalar
186 * trial functions. Essential for mixed formulations enforcing divergence constraints
187 * such as incompressibility. Applications include:
188 * - Pressure-velocity coupling in incompressible flow
189 * - Mixed formulations for Darcy flow
190 * - Constraint enforcement in optimization problems
191 * - Divergence-free condition implementation
192 *
193 * @note \f$\lambda_{i,i}\f$ is the divergence of vector test function
194 * @note \f$u\f$ is a scalar trial function
195 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
196 *
197 * @tparam SPACE_DIM Spatial dimension of the problem
198 */
199 template <int SPACE_DIM>
201
202 /**
203 * @brief Mixed tensor divergence-vector bilinear form: \f$(\lambda_{ij,j},u_{i})_\Omega\f$
204 * @ingroup mofem_forms
205 *
206 * This operator integrates the divergence of tensor test functions with vector
207 * trial functions. Used in mixed formulations with stress-displacement coupling
208 * and equilibrium enforcement. Applications include:
209 * - Stress-displacement mixed formulations in elasticity
210 * - Mixed finite element methods for mechanics
211 * - Equilibrium constraint enforcement
212 * - Multi-physics coupling with tensor fields
213 *
214 * @note \f$\lambda_{ij,j}\f$ is the divergence of tensor test function
215 * @note \f$u_i\f$ is a vector trial function
216 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
217 *
218 * @tparam SPACE_DIM Spatial dimension of the problem
219 * @tparam CoordSys Coordinate system type (default: CARTESIAN)
220 */
221 template <int SPACE_DIM, CoordinateTypes CoordSys = CARTESIAN>
223
224 /**
225 * @brief Mixed scalar-divergence bilinear form: \f$(\lambda,u_{i,i})_\Omega\f$
226 * @ingroup mofem_forms
227 *
228 * This operator integrates scalar test functions with the divergence of vector
229 * trial functions. This is the transpose of OpMixDivTimesScalar and enforces
230 * the same constraints but with different test/trial function roles. Applications include:
231 * - Symmetric mixed formulations
232 * - Lagrange multiplier methods for divergence constraints
233 * - Saddle point problems in fluid mechanics
234 * - Mixed Poisson formulations
235 *
236 * @note \f$\lambda\f$ is a scalar test function
237 * @note \f$u_{i,i}\f$ is the divergence of vector trial function
238 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
239 *
240 * @tparam SPACE_DIM Spatial dimension of the problem
241 * @tparam COORDINATE_SYSTEM Coordinate system type (default: CARTESIAN)
242 */
243 template <int SPACE_DIM, CoordinateTypes COORDINATE_SYSTEM = CARTESIAN>
246
247 /**
248 * @brief Mixed vector-gradient bilinear form: \f$(\lambda_{i},u_{,j})_\Omega\f$
249 * @ingroup mofem_forms
250 *
251 * This operator integrates vector test functions with the gradient of trial
252 * functions. Used in mixed formulations where vector fields couple with
253 * gradient terms. Applications include:
254 * - Mixed formulations with vector Lagrange multipliers
255 * - Gradient constraint enforcement
256 * - Vector-scalar coupling in multi-physics problems
257 * - Stabilization terms in finite element methods
258 *
259 * @note \f$\lambda_i\f$ is a vector test function
260 * @note \f$u_{,j}\f$ is the gradient of trial function
261 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
262 *
263 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
264 * @tparam FIELD_DIM Dimension of the field
265 * @tparam SPACE_DIM Spatial dimension of the problem
266 */
267 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
270
271 /**
272 * @brief Mixed tensor-gradient bilinear form: \f$(\lambda_{ij},u_{i,j})_\Omega\f$
273 * @ingroup mofem_forms
274 *
275 * This operator integrates tensor test functions with the gradient of vector
276 * trial functions. Essential for mixed formulations in solid mechanics where
277 * stress tensors couple with displacement gradients. Applications include:
278 * - Stress-strain mixed formulations
279 * - Mixed finite element methods in elasticity
280 * - Constitutive relation enforcement
281 * - Coupled stress-displacement problems
282 *
283 * @note \f$\lambda_{ij}\f$ is a tensor test function
284 * @note \f$u_{i,j}\f$ is the gradient of vector trial function
285 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
286 *
287 * @tparam SPACE_DIM Spatial dimension of the problem
288 */
289 template <int SPACE_DIM>
291
292 /**
293 * @brief Hybridization constraint bilinear form for broken spaces
294 * @ingroup mofem_forms
295 *
296 * This operator assembles constraint matrices during hybridization of broken
297 * finite element spaces. Hybridization introduces Lagrange multipliers on
298 * element interfaces to enforce continuity in discontinuous Galerkin methods.
299 * The operator handles the coupling between interior element solutions and
300 * interface unknowns. Applications include:
301 * - Discontinuous Galerkin (DG) methods
302 * - Hybridizable discontinuous Galerkin (HDG) methods
303 * - Mixed finite element methods with broken spaces
304 * - Static condensation procedures
305 *
306 * @note Used specifically for broken space hybridization procedures
307 * @note Handles interface coupling in DG methods
308 * @note Implementation details are in BiLinearFormsIntegratorsImpl.hpp
309 *
310 * @tparam FIELD_DIM Dimension of the field
311 */
312 template <int FIELD_DIM>
315};
316
317} // namespace MoFEM
318
319#endif //__BILINEAR_FORMS_INTEGRATORS_HPP__
constexpr int SPACE_DIM
constexpr int FIELD_DIM
constexpr int BASE_DIM
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
constexpr IntegrationType I