 |
deal.II version GIT relicensing-7040-g409d9ea789 2026-11-06 04:40:01+00:00
|
|
deal.II version GIT relicensing-7040-g409d9ea789 2026-11-06 04:40:01+00:00
|
Finite element operators are built from contributions on cells and, for some formulations, on faces. deal.II offers several ways to organize this computation. The choice concerns both the traversal of the mesh and whether the global matrix is stored explicitly.
In step-3, an explicit loop visits each active cell. FEValues evaluates shape functions and quadrature data, nested loops compute a local matrix and right hand side, and their entries are added to the global objects. This exposes all parts of assembly and is a useful starting point for a new formulation. Applications with hanging nodes or other constraints can use AffineConstraints::distribute_local_to_global() to transfer local contributions, as shown in step-6.
For example, the cell loop for a scalar Laplace matrix looks as follows. Here fe is the finite element, constraints is closed, and system_matrix has the required size and sparsity pattern. FEValues requests gradients and quadrature weights; the local objects are reused between cells.
WorkStream::run() separates the computation of local contributions from their insertion into global objects. A worker uses scratch data to compute each cell's contribution and stores the result in copy data. A copier then inserts that result into the global matrix and right hand side. This separation allows local work to run on multiple threads. For the iterator-range overload, insertion into shared objects is serialized. The colored overload can also run copiers concurrently: its coloring must ensure that copiers within one color do not write to the same memory. The local finite element calculations can remain the same as in an explicit cell loop. See step-9 and the shared-memory parallelism documentation for examples and details of the WorkStream model.
The iterator-range overload replaces the explicit cell loop by a worker and a copier. In this sketch, ScratchData owns an FEValues object and CopyData owns the local matrix and DoF indices. The application function assemble_cell resets and fills the copy data using the local calculation above; it only writes to its scratch and copy arguments. Both data types must satisfy WorkStream's copying requirements (see step-9).
In the scheme outlined above, one has a loop over all cells (either explicitly, or by letting WorkStream::run() iterate over them) and on each cell one assembles the necessary terms that pertain to this cell. Sometimes, as in the example above, these terms really only contain a contribution from each cell. But sometimes, one also has contributions from boundary faces (say, when dealing with Neumann or Robin boundary conditions) or interior faces (say, when dealing with discontinuous Galerkin schemes). In these cases, computing all of these contributions in one place might not be the best software design – it might be better to split different contributions into separate functions.
MeshWorker::mesh_loop() builds on WorkStream and provides separate callbacks for cells, boundary faces, and interior faces. Flags select which contributions to compute and how interior faces are visited. This is useful for formulations with face terms, such as discontinuous Galerkin methods: the application supplies the local integrals and the loop handles traversal. Scratch and copy data still separate local computation from insertion into global objects. See step-12 and the The MeshWorker interface module.
The traversal can include all three kinds of contributions in one call:
Here the application has also defined boundary_worker and face_worker with the signatures documented by MeshWorker::mesh_loop(). They append face contributions to CopyData; copy_cell_and_face_data inserts both cell and face contributions into global objects. This requires additional face storage in the scratch and copy data, as illustrated in step-12.
When an iterative solver only needs the action of an operator on a vector, it is also possible to avoid assembling and storing the global matrix. MatrixFree organizes cell and face loops and the data needed for them; FEEvaluation and FEFaceEvaluation evaluate finite element functions and integrate local contributions. These contributions are accumulated into a destination vector in each operator application. For tensor-product elements, sum factorization can reduce the cost of local evaluation, particularly at higher polynomial degrees. See step-37 and the Matrix-free infrastructure module.
This approach changes the operator representation, not just the assembly loop. It saves matrix storage but requires an operator implementation and a compatible preconditioner. An assembled matrix remains useful when a solver or preconditioner requires its entries. Applications can also combine matrix-free operator application with assembled auxiliary matrices. For the same constant-coefficient Laplace operator, the following member functions of an application class LaplaceOperator replace matrix assembly and multiplication. The class has an initialized MatrixFree<dim> member data, a compile-time polynomial degree fe_degree, and the alias VectorType for LinearAlgebra::distributed::Vector<double>. Its vectors use the layout supplied by MatrixFree::initialize_dof_vector().
The loop runs over cell batches. The final true requests zeroing dst before accumulating contributions; cell_loop() also handles the required communication for these distributed vectors. Boundary conditions and constraint treatment must match the operator used by the solver, as in step-37.
1.9.8