Open gassmoeller opened 4 years ago
Specifically, the way we currently see the mesh deformation is as the time integral of the mesh velocity, and we compute updates to the deformation via the box rule as timestep * velocity
in the code above. We should compute this integral with a better scheme, for example the trapezoidal rule (=Crank-Nicolson) or using a 3-point interpolation (similar to the BDF-2 scheme).
The free surface is currently advected using a forward euler scheme in https://github.com/geodynamics/aspect/blob/95a4d4567e4dc8f5eb3faf856d260bd752cfd9f6/source/mesh_deformation/interface.cc#L582
We should investigate if we can use better time-stepping schemes like Crank-Nicholson. We already store the old Stokes solution, but we may also need to store old mesh velocities.