Hi all,

I have a problem with numerical instabilities causing spurious strainrates. The model setup is as follows:

- isoviscous, closed, free-slip box, 660x2640 km
- cold lithosphere at surface
- 250-K thermal instability at the bottom boundary at x=1000 km (to develop mantle plume)

I’ve stripped the model setup down to a small model with very few ingredients, but the problem still remains. The following images illustrate the problem. The strainrate field at t=0 looks like:

and a corresponding temperature field that looks smooth:

The corresponding input file can be found here:

The mesh is 128x32 elements, without mesh refinement, and solver tolerance is set low (1e-9). The visco-plastic material model is used (it is isoviscous now, nbut the original problem was rheologically more complex). Local mesh refinement seems to make it worse.

The problem occurs (or can be noticed) only when overall flow velocities are low. Velocities at t=0 are small, but not beyond the point of resolving, I would have thought (up to 0.4 mm/yr, although this is probably at least in part due to the numerical instability). The problem remains for for later model times, until velocities are large enough to overshadow this problem.

Has anyone encountered this problem, and if so, is there a solution?

Many thanks,

Jeroen van Hunen