v0.15.0
Loading...
Searching...
No Matches
LinearFormsIntegrators.hpp
Go to the documentation of this file.
1/** \file LinearFormsIntegrators.hpp
2 * \brief Linear forms integrators
3 * \ingroup mofem_form
4 *
5 * This file provides template aliases for linear form integrators used in finite
6 * element assembly. Linear forms represent the right-hand side vectors in finite
7 * element formulations \f$(v, f)_\Omega\f$ where \f$v\f$ is a test function and
8 * \f$f\f$ is a source term or field. The integrators handle various mathematical
9 * operations including scalar/vector/tensor field integration, mixed formulations,
10 * boundary conditions, and convective terms. All aliases reference implementations
11 * in LinearFormsIntegratorsImpl.hpp for actual computational details.
12 *
13*/
14
15#ifndef __LINEAR_FORMS_INTEGRATORS_HPP__
16#define __LINEAR_FORMS_INTEGRATORS_HPP__
17
18
19
20namespace MoFEM {
21
22/**
23 * @brief Linear integrator form
24 * @ingroup mofem_forms
25 *
26 * @tparam EleOp
27 * @tparam A
28 * @tparam I
29 */
30template <typename EleOp>
31template <AssemblyType A>
32template <IntegrationType I>
33struct FormsIntegrators<EleOp>::Assembly<A>::LinearForm {
34
35 /**
36 * @brief Integrate \f$(v,f(\mathbf{x}))_\Omega\f$ where f is a scalar function
37 * @ingroup mofem_forms
38 *
39 * This operator integrates the product of test functions with a scalar source
40 * function defined at integration points. Commonly used for:
41 * - Body forces in mechanics (gravitational, thermal expansion)
42 * - Source terms in heat transfer and diffusion problems
43 * - Prescribed loads and boundary conditions
44 *
45 * @note \f$f(\mathbf{x})\f$ is a scalar function of coordinates
46 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
47 *
48 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
49 * @tparam FIELD_DIM Dimension of the field being integrated
50 * @tparam S Source function specialization type (default: SourceFunctionSpecialization)
51 */
52 template <int BASE_DIM, int FIELD_DIM,
53 typename S = SourceFunctionSpecialization>
54 using OpSource =
56
57 /**
58 * @brief Scalar field integrator \f$(v_i,f)_\Omega\f$ where f is a scalar
59 * @ingroup mofem_forms
60 *
61 * This operator integrates the product of vector test functions with a scalar
62 * field value. Each component of the test function is multiplied by the same
63 * scalar value. Useful for:
64 * - Uniform pressure loads in structural mechanics
65 * - Constant temperature fields in thermal problems
66 * - Scalar material properties affecting vector fields
67 *
68 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
69 *
70 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
71 * @tparam S Stride or scaling factor (default: 1)
72 */
73 template <int BASE_DIM, int S = 1>
75
76 /**
77 * @brief Deprecated alias for OpBaseTimesScalar
78 * @deprecated Use OpBaseTimesScalar instead for consistency
79 *
80 * @tparam BASE_DIM Dimension of the basis functions
81 * @tparam S Stride or scaling factor (default: 1)
82 */
83 template <int BASE_DIM, int S = 1>
85
86 /**
87 * @brief Vector field integrator \f$(v,f_i)_\Omega\f$ where f is a vector
88 * @ingroup mofem_forms
89 *
90 * This operator integrates the product of test functions with vector field
91 * values. Each test function component is multiplied by the corresponding
92 * component of the vector field. Applications include:
93 * - Vector body forces (wind loads, electromagnetic forces)
94 * - Distributed loads with directional components
95 * - Multi-component source terms in coupled problems
96 *
97 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
98 *
99 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
100 * @tparam FIELD_DIM Dimension of the vector field
101 * @tparam S Stride or scaling factor
102 */
103 template <int BASE_DIM, int FIELD_DIM, int S>
106 //! [Source operator]
107
108 //! [Grad times tensor]
109
110 /**
111 * @brief Integrate \f$(v_{,i},f_{ij})_\Omega\f$ where f is a tensor
112 * @ingroup mofem_forms
113 *
114 * This operator integrates the product of test function gradients with tensor
115 * field values. The gradient of test functions contracts with the tensor field.
116 * Common applications include:
117 * - Stress-based formulations in mechanics
118 * - Flux calculations in transport problems
119 * - Gradient-dependent source terms
120 *
121 * @note \f$f_{ij}\f$ is a tensor field at integration points
122 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
123 *
124 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
125 * @tparam FIELD_DIM Dimension of the field
126 * @tparam SPACE_DIM Spatial dimension of the problem
127 * @tparam S Stride or scaling factor (default: 1)
128 */
129 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1>
132
133 /**
134 * @brief Integrate \f$(v_{,i},f_{ij})_\Omega\f$ where f is a symmetric tensor
135 * @ingroup mofem_forms
136 *
137 * This operator is specialized for symmetric tensors, providing computational
138 * efficiency when the tensor field has symmetry properties. Applications include:
139 * - Symmetric stress tensors in solid mechanics
140 * - Symmetric diffusion tensors in transport problems
141 * - Symmetric material property tensors
142 *
143 * @note \f$f_{ij}\f$ is a symmetric tensor field at integration points
144 * @note Implementation details are in LinearFormsIntegratorsImpl.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 Stride or scaling factor (default: 1)
150 */
151 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1>
154
155 /**
156 * @brief Mixed formulation integrator: \f$(\lambda_{ij,j},u_{i})_\Omega\f$
157 * @ingroup mofem_forms
158 *
159 * This operator integrates the divergence of a tensor field \f$\lambda_{ij}\f$
160 * with vector test functions \f$u_i\f$. Essential for mixed formulations where
161 * stress equilibrium or vector field divergence is enforced. Applications include:
162 * - Stress-displacement coupling in solid mechanics
163 * - Mixed formulations with tensor constraints
164 * - Multi-physics problems with vector-tensor coupling
165 *
166 * @note \f$\lambda_{ij,j}\f$ is the divergence of tensor field at integration points
167 * @note \f$u_i\f$ is the vector test function
168 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
169 *
170 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
171 * @tparam FIELD_DIM Dimension of the field
172 * @tparam SPACE_DIM Spatial dimension of the problem
173 * @tparam CoordSys Coordinate system type (default: CARTESIAN)
174 */
175 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM,
176 CoordinateTypes CoordSys = CARTESIAN>
179
180 /**
181 * @brief Mixed formulation integrator: \f$(\lambda_{ij},u_{i,j})_\Omega\f$
182 * @ingroup mofem_forms
183 *
184 * This operator integrates tensor field values with the gradient of vector
185 * test functions. The tensor field contracts with the gradient of test functions.
186 * Applications include:
187 * - Stress-strain coupling in mechanics
188 * - Mixed formulations with gradient constraints
189 * - Tensor-gradient field interactions
190 *
191 * @note \f$\lambda_{ij}\f$ is a tensor field at integration points
192 * @note \f$u_{i,j}\f$ is the gradient of vector test function
193 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
194 *
195 * @tparam SPACE_DIM Spatial dimension of the problem
196 */
197 template <int SPACE_DIM>
199
200 /**
201 * @brief Mixed formulation integrator: \f$(u_{i},\lambda_{ij,j})_\Omega\f$
202 * @ingroup mofem_forms
203 *
204 * This operator integrates vector test functions with the divergence of a
205 * tensor field. This is the adjoint form of OpMixDivTimesU, commonly used
206 * in mixed formulations for equilibrium enforcement. Applications include:
207 * - Adjoint stress equilibrium in mechanics
208 * - Dual formulations with tensor divergence
209 * - Vector-tensor coupling in mixed problems
210 *
211 * @note \f$u_i\f$ is a vector test function
212 * @note \f$\lambda_{ij,j}\f$ is the divergence of tensor field at integration points
213 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
214 *
215 * @tparam SPACE_DIM Spatial dimension of the problem
216 */
217 template <int SPACE_DIM>
220
221 /**
222 * @brief Boundary integrator: normal vector times scalar function
223 * @ingroup mofem_forms
224 *
225 * This operator multiplies vector test functions by the normal vector on
226 * a face and integrates with a scalar field. Essential for natural boundary
227 * conditions in mixed formulations. Applications include:
228 * - Natural boundary conditions in fluid mechanics
229 * - Surface flux conditions in transport problems
230 * - Mixed boundary value problems
231 *
232 * @note This operator is used on domain boundaries (faces)
233 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
234 *
235 * @tparam SPACE_DIM Spatial dimension of the problem
236 */
237 template <int SPACE_DIM>
240
241 /**
242 * @brief Boundary integrator: normal vector times vector field
243 * @ingroup mofem_forms
244 *
245 * This operator multiplies vector test functions by the normal vector on
246 * a face and integrates with a vector field. Used for natural boundary
247 * conditions involving vector quantities. Applications include:
248 * - Vector-valued boundary conditions in mechanics
249 * - Surface traction conditions in solid mechanics
250 * - Vector flux conditions in multi-physics problems
251 *
252 * @note This operator is used on domain boundaries (faces)
253 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
254 *
255 * @tparam SPACE_DIM Spatial dimension of the problem
256 */
257 template <int SPACE_DIM>
260
261 /**
262 * @brief Convective term integrator: \f$(v, u_i \mathbf{y}_{,i})\f$
263 * @ingroup mofem_forms
264 *
265 * This operator integrates convective terms where a velocity field \f$u_i\f$
266 * transports a scalar or vector field \f$\mathbf{y}\f$. The convective operator
267 * \f$u_i \mathbf{y}_{,i}\f$ represents advection. Applications include:
268 * - Advection in fluid mechanics
269 * - Convective transport in heat transfer
270 * - Material transport in multi-physics problems
271 *
272 * @note \f$u_i\f$ is the velocity field at integration points
273 * @note \f$\mathbf{y}_{,i}\f$ is the gradient of the transported field
274 * @note \f$\mathbf{y}\f$ can be scalar or vector field
275 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
276 *
277 * @tparam BASE_DIM Dimension of the basis functions (1D, 2D, 3D)
278 * @tparam FIELD_DIM Dimension of the field
279 * @tparam SPACE_DIM Spatial dimension of the problem
280 */
281 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
284
285 /**
286 * @brief Hybridization constraint integrator for broken spaces
287 * @ingroup mofem_forms
288 *
289 * This operator assembles the right-hand side for constraint matrix C
290 * during hybridization of broken finite element spaces. Hybridization
291 * introduces Lagrange multipliers to enforce continuity across element
292 * boundaries in discontinuous Galerkin methods. Applications include:
293 * - Discontinuous Galerkin methods
294 * - Mixed finite element methods with broken spaces
295 * - Hybridizable discontinuous Galerkin (HDG) methods
296 *
297 * @note Used specifically for broken space hybridization
298 * @note Implementation details are in LinearFormsIntegratorsImpl.hpp
299 *
300 * @tparam FIELD_DIM Dimension of the field
301 */
302 template <int FIELD_DIM>
305
306 /**
307 * @brief Assemble Rhs for contraint matrix C^T whole hybridisation
308 *
309 * @tparam FIELD_DIM
310 */
311 template <int FIELD_DIM>
314};
315
316} // namespace MoFEM
317
318#endif // __LINEAR_FORMS_INTEGRATORS_HPP__
constexpr int SPACE_DIM
constexpr int FIELD_DIM
CoordinateTypes
Coodinate system.
@ CARTESIAN
constexpr int BASE_DIM
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
Tensor field time gradient of vector field.
Vector field time divergence of tensor.
Multiply vector times normal on the face times scalar function.
Multiply vector times normal on the face times vector field.