Numerical Quadrature and Assembly

Finite element matrices are computed cell by cell and then assembled into a global sparse matrix.

Reference element mapping. Each cell ${\Omega }_e$ is the image of a reference cell $\hat{K}$ under a map $F_e$ with Jacobian $J_e$. Shape functions and gradients transform as

$${\phi }_i={\hat{\phi }}_i\circ F_e^{-1},\nabla {\phi }_i=J_e^{-\top }\hat{\nabla }{\hat{\phi }}_i.$$

Quadrature. Integrals over the reference cell are replaced by weighted sums at quadrature points ${\hat{x}}_q$ with weights $w_q$:

$${\int }_{{\Omega }_e}f\mathrm{d}x\approx \sum _{q}f(F_e({\hat{x}}_q))w_q|\det J_e({\hat{x}}_q)|.$$

Gauss rules with $m$ points per direction integrate polynomials of degree $2m-1$ exactly. For $Q_1$ stiffness matrices with a constant or $Q_1$ conductivity, $2\times 2$ Gauss points are standard.

Assembly loop.

A ← 0
for each cell e:
    A_e ← local matrix from quadrature
    A[dofs(e), dofs(e)] += A_e

The same pattern assembles the Mass Matrix, the Stiffness Matrix, the Weighted Stiffness Matrix, load vectors, and the gradient of the data misfit (see Functional Derivative of the Data Misfit). Boundary terms use a loop over boundary facets with facet quadrature (see Boundary Mass and Stiffness Matrices).

The sparsity pattern depends only on the mesh and the degrees of freedom. It can be computed once and the values overwritten in each iteration.

References

1. S. C. Brenner, L. R. Scott (2008). The Mathematical Theory of Finite Element Methods, 3rd ed. Springer. doi:10.1007/978-0-387-75934-0
2. K. Carlsson, F. Ekre, and contributors. Ferrite.jl (software). github.com/Ferrite-FEM/Ferrite.jl