Anomalies at mesh refinement boundaries in viscoelastic model

[bobmyhill]

Related to #6248, but this is a solver-agnostic problem, and in 2D.

While working on #6829, I noticed that one of the tests that changed had peculiar anomalies at mesh refinement boundaries (boundary_traction_function_cartesian_free_surface.prm). The anomalies themselves do not originate in that PR, and are not fixed by it.

I very slightly tweaked the prm file so that I could visualize the output:
boundary_traction_function_cartesian_free_surface.txt

Removing the boundary traction part of the prm and fixing the elastic timestep to 2.5 yrs changes the result, but does not remove the anomalies:
free_surface.txt

The anomalies are not present when I use a purely viscous material model ("Simpler"):
boundary_traction_function_cartesian_free_surface_simpler.txt

[tiannh7]

Hi @bobmyhill,

As you may have seen, I’ve actually encountered some tricky issues myself when working with traction + free surface, such as in Velocity and Strain Rate Anomalies at Mesh Refinement Boundaries in 3D Spherical Models with Solver Differences. Most of these have been resolved, although the fixes haven’t been pushed promptly.

The core of this issue is that for non-box geometries, usually Q2 or even higher-order elements are required for mesh deformation. @gassmoeller has done most of the work on this in Improve free surface accuracy, but for GMG there is still a mismatch between the order of the Laplace equation and the mesh deformation, for example:

Assert(mesh_deformation_fe.degree == 1, ExcNotImplemented());
// To be efficient, the operations performed in the matrix-free implementation require
// knowledge of loop lengths at compile time, which are given by the degree of the finite element.
const unsigned int mesh_deformation_fe_degree = 1;

using SystemOperatorType = dealii::MatrixFreeOperators::
                           LaplaceOperator<dim, mesh_deformation_fe_degree, mesh_deformation_fe_degree + 1, dim>;

I have extended and fixed Improve free surface accuracy on my side, and testing shows that it resolves this problem. This work is reflected in my latest commit in Fix mesh deformation for GMG by adding support for higher-order FE degrees.

I really appreciate everyone’s help, which has made this progress possible.

Thanks again for the ongoing support!

[bobmyhill]

Dear @tiannh7 ,

Thanks for your message! So am I right in thinking that the fix is in #6284 for AMG, and #6836 for GMG? I notice that #6284 has been [WIP] for a while - do you or @gassmoeller have a clear picture of what needs to be done to get that merged?

[tiannh7]

Hi Bob,

So am I right in thinking that the fix is in #6284 for AMG, and #6836 for GMG?

Yes, for #6284 AMG already achieves near-correct results. In #6836 I made some minor changes by adding support for the degree parameter to make it more flexible, and also ensured the Amg_data parameter is properly passed to the solver.

do you or @gassmoeller have a clear picture of what needs to be done to get that merged?

If @gassmoeller doesn’t mind, it might be best to merge #6836 as a separate PR to resolve this issue directly.

[bobmyhill]

We can certainly do that. It would also be good to get #6284 moving again, especially if it is nearly ready for review.