Long term simulation of 2-D model with gravity

[KejingHuang]

Hi all

I constructed a simple 2-D model with an elastic crust, an elastic loading body and a viscoelastic mantle as the picture below shows.

My target is to see the subsidence rate of the surface under the locale loading body after a half million years. However, when I tried to insert gravity following section 7.17.4 of the manual, the velocity output
of the domain turns out to be asymmetrical and strange under a symmetrical loading after 20 000 years. The figures below show the velocity results after 15 000 year’s simulation and 500 000 year’s simulation, respectively.

Do you know what’s wrong with the simulation or how can I fulfill my target reasonably?

Kejing

[baagaard]

It would be helpful if you could supply more information about your problem and what you are trying to do. What is the Maxwell relaxation time and time step in your simulation? What are you trying to resolve on the 500k year time scale?

As for the asymmetry, is the mesh perfectly symmetric? Is the asymmetry a function of the time step size and solver tolerances? If you are taking many thousands of time steps, then roundoff errors may be accumulating on the order of the solver tolerances.

[KejingHuang]

Thank you for your quick reply.

I’m trying to resolve the surface defamation rate under the loading of a volcano which formed ~500k years ago (simplified as a rectangular body at the top) and make a comparison to the current observations. The Maxwell relaxation time and time step are 100 yr and 20yr, respectively. And since the deformation rate continues to decrease with time, I tried to apply this long running time to obtain a stable solution.

And the mesh is constructed under perfectly symmetric from my perspective. I tried 1e-12 and 1e-25 of the solver tolerance and the asymmetry both gradually take place after ~20k yr’s simulation. And the KSP Residual converged to less than 1e-27 after ~30 iterations and doesn’t go smaller if I give more strict solver tolerance.

Maybe I misunderstand something to estimate the surface deformation like this. Do you have further suggestions?

Sincerely
Kejing

[baagaard]

Solver tolerances of 1.0e-20 or below are close to the numerical precision available with 64 bit floating point values. With only 5 time steps per Maxwell time, you are not resolving the Maxwell viscoelastic relaxation very well. With a Maxwell relaxation time of 100 years and a total duration of 500k years, you are trying to resolve the deformation that is 5000 times the Maxwell relaxation time. If this is indeed the appropriate Maxwell time, then you should be able to find the stable value after about 10k years (before any asymmetry appears).

I suggest reducing your time step and reducing the duration of your simulation to find the state of the system. In other words, the relaxation time is short compared to 500k years, so it doesn’t make sense to try to resolve the rate of deformation after 500k years in a numerical model.

[KejingHuang]

Thanks for your suggestions.

I just tried 5 years for time step and 10k years for total simulation time without other modifications. The results are symmetric as I expected. But while the the surface loading weight is comparable to the real case, the resulting surface subsidence rate is about 10cm/yr (if the default unit of velocity output is m/s), which is far more than the observations (2~3 mm/yr). Is there any simulation step that I missed?

[baagaard]

How well constrained are all of the other parameters in the simulation, such as shear modulus? When comparing against observations, you need to consider all sources of uncertainty.

[seob]

hi
i saw your model then i have a question to you.

When my model was lowercrust and uppermantle of viscoelastic material in gravity, the cell collapsed.

To prevent this, I used initial stress as in step 2 of the 7.17.5 chapter in manual.

I wonder if you have used this method too.