search

geodynamics / aspect

A parallel, extensible finite element code to simulate convection in both 2D and 3D models.
https://aspect.geodynamics.org/
Other
224 stars 235 forks source link

get_entropy_viscosity: Why is material model called with the average solution of the previous 2 timesteps instead of the current linearization point? #3933

Open anne-glerum opened 3 years ago

anne-glerum commented 3 years ago

As we couldn't get to a conclusion on this question during the user meeting, here a description of my question for further investigation:

Before the solving of each advection equation, the material model is called first from get_artificial_viscosity function and second from the assembly. The MaterialModelInput in both cases is different: in the get_artificial_viscosity function it is set to the average of the last 2 solutions (old_solution + old_old_solution)/2, while in the assembly it is set to the current linearization point, which is (1+dt/dt_old)old_solution - (dt/dt_old)old_old_solution. Why are different inputs used?

The get_artificial_viscosity function takes the density, specific heat and thermal conductivity from the material model output, which can depend on the temperature, pressure and composition solution, but these outputs are only used for the SUPG method in case the field is the temperature field.

tjhei commented 3 years ago

you are referring to the following, right? https://github.com/geodynamics/aspect/blob/fed300b4b4ea048f73418f26e37b6f72373261d9/source/simulator/entropy_viscosity.cc#L509-L517

This is not an extrapolation but an interpolation, which is weird. I do not remember why we do this, but maybe this is explained in the first ASPECT paper...

anne-glerum commented 3 years ago

Yes, that's what I mean. (I didn't know about linking to specific lines, very handy!) I'll have a look in the paper.

tjhei commented 3 years ago

Indeed: "Since we do not want the artificial viscosity to introduce a non-linear dependence on the current temperature T n , we make it explicit by approximating the time derivative from the two previous time levels in the BDF-2 time stepping, ∂E(T h )/∂t ≈ (E(T n−1 ) − E(T n−2 ))/k n−1 , and evaluating all other occurrences of the tempera- ture at the midpoint as (T n−1 + T n−2 )/2, including the average tem- perature."

tjhei commented 3 years ago

Maybe @bangerth remembers, but I would think this doesn't make much sense and is from the time where we made advection explicit.