Which novel data structures are used in adaptive FEM?

Question

A lot of adaptive FEM libraries use more advanced mesh data structures to handle adding/removing nodes, edges, triangles, tetrahedra, etc. For example, the p4est library uses octree data structures for adaptive mesh refinement; you wouldn't often find octrees used for computations on a static mesh.

What changes on the linear algebra side for adaptive FEM?

The most blunt way I can conceive of would be to completely rebuild all the system matrices whenever the mesh is refined or coarsened. If mesh adaptation is a sufficiently infrequent operation, then the expense of doing so is ultimately amortized over the rest of the computation. One could easily leverage existing sparse linear algebra software (PETSc, Trilinos, etc.) with this approach.

Is this blunt method the most commonly used, or are there libraries that manage to reuse or modify the old matrix during refinement? After all, most of the mesh and the corresponding matrices are unchanged during a mesh adaptation.

Answer 1

Yes, the most common approach is to rebuild. Data structures that are modifiable in-place tend to be less efficient once set up, and reallocation is actually quite cheap compared to reassembly (e.g., due to nonlinearity) so it's really a fine solution. Outside of relatively rare niches with very easy solves, attempts to use dynamic data structures in the solvers will only make your application slower. It is a common perception among those who forget to measure or model performance, however.

Answer 2

Like Jed already said, reuse of linear algebra components such as matrices and vectors is not commonly done. It is also not necessary: setting up these components is comparably very cheap relative to the cost of solving linear systems.

If you look for things that change when going from static to adaptively refined meshes, then the biggest obstacle is dealing with hanging nodes. In deal.II, this is handled by the ConstraintMatrix class, which takes several 1000 lines of code. You can find some of the description of what this class does in the paper by myself and Oliver Kayser-Herold (linked to from my publications page). No other component (other than handling the mesh, of course) has required so much adaptation when going from fixed to adaptive meshes.