fdrmrc / Polydeal

C++ implementation of Polygonal Discontinuous Galerkin method within the deal.II Finite Element library.
https://fdrmrc.github.io/Polydeal/
Other
0 stars 0 forks source link

Add new Mapping based on agglomerates #133

Closed fdrmrc closed 2 months ago

fdrmrc commented 2 months ago

So far, we have been using a MappingFEField machinery which required the definition of an auxiliary DoFHandler and an Eulerian vector describing the vertices of boxes associated with agglomerates. After this setup, this mapping was given to FEValues-like objects with properly scaled quadrature weights and normals since our weights and normals were already in real space. This was rather unfortunate and was explained in https://github.com/fdrmrc/Polydeal/issues/58#issuecomment-1871155392


This PR depends on #132. First relevant commit is 1869db1b65aa32d9637befef55d9ab40adc665a6.

It provides a new MappingBox class responsible for mappings associated with the bounding box of agglomerates and aware of the fact that we are in real space. It describes the affine transformation $F(\hat{x}) = J\hat{x} + c$ from reference space to the (real) bounding box. In this case $J$ is diagonal with entries $h_x,h_y,h_z$. It is very close to what MappingCartesian does but it considers the fact that quantities such as weights and normals must not be scaled in our case, as they are collected from real space.

Its construction requires only the vector of local BoundingBox<dim> objects and a map which translates from the (master) cell index to the associated box index (which coincides with the polytope index). We can then construct FEValues-like objects with such a MappingBox object.

This approach completely removes the need for the mentioned MappingFEField and allows removing all the tricks with scaling factors that we had to do before. This is clear by looking at changes in source/agglomeration_handler.cc.

@luca-heltai What do you think? The performance gain with this new design is quite large:

Wall-clock time to assemble DG stiffness matrix

$$ \begin{aligned} & \text {Wall-clock time [s] for assemble system()}\ &\begin{array}{cccc} \hline \hline \text { Degree } & \text { Before this PR } & \text { This PR } \ \hline 1 & 11.23 & 0.618 \ 2 & 17.68 & 2.98 \ 3 & 40.04 & 22.06 \ \hline \end{array} \end{aligned} $$

fdrmrc commented 2 months ago

The failing test will be fixed in #134