How can I solve the following convergence issue?

Dear users.

I’m simulating LLSVPs in a Steinberger problem, so both my mantle and my basal layer use perplex tables which are computed up to 6000 K.

The point is that, when my CMB has very high temperatures (6000 K or 5500 K) happen that in certain “special conditions” when my LLSVP reach the surface (very high Rayileigh number and low bouyancy) the simulation doesn’t converge (that doesn’t happen when my CMB has temperatures like 5000 K or lower. The behave is the same, but the timestep doesn’t go to very low as below), this is the error:

*** Timestep 11233:  t=5.36363e+08 years, dt=4731.1 years
   Dynamic core data updated.
     Tc(K)          Ri(km)         Xi             dT/dt(K/year)  dR/dt(km/year) dX/dt(1/year)
     5616.42        1211.69        0.0537547      -5.45334e-07   1.26888e-06    3.63791e-12
   Solving temperature system... 16 iterations.
   Advecting particles...  done.
   Copying properties into prescribed compositional field density_field... done.
   Rebuilding Stokes preconditioner...
   Solving Stokes system (AMG)... 794+0 iterations.

---------------------------------------------------------
TimerOutput objects finalize timed values printed to the
screen by communicating over MPI in their destructors.
Since an exception is currently uncaught, this
synchronization (and subsequent output) will be skipped
to avoid a possible deadlock.
---------------------------------------------------------


----------------------------------------------------
Exception 'ExcMessage("The time step length for the each time step needs to be positive, " "but the computed step length was: " + std::to_string(min_conduction_timestep) + ". " "Please check for non-positive material properties.")' on rank 51 on processing: 


----------------------------------------------------
Exception 'ExcMessage("The time step length for the each time step needs to be positive, " "but the computed step length was: " + std::to_string(min_conduction_timestep) + ". " "Please check for non-positive material properties.")' on rank 45 on processing: 


----------------------------------------------------
Exception 'ExcMessage("The time step length for the each time step needs to be positive, " "but the computed step length was: " + std::to_string(min_conduction_timestep) + ". " "Please check for non-positive material properties.")' on rank 3 on processing:
[19:53]
----------------------------------------------------
Exception 'ExcMessage("The time step length for the each time step needs to be positive, " "but the computed step length was: " + std::to_string(min_conduction_timestep) + ". " "Please check for non-positive material properties.")' on rank 46 on processing: 

--------------------------------------------------------
An error occurred in line <107> of file </home/francesco-radica/Documenti/aspect_dynamic_core/source/time_stepping/conduction_time_step.cc> in function
    double aspect::TimeStepping::ConductionTimeStep<dim>::execute() [with int dim = 2]
The violated condition was: 
    min_conduction_timestep > 0
Additional information: 
    The time step length for the each time step needs to be positive, but
    the computed step length was: 0.000000. Please check for non-positive
    material properties.

Stacktrace:
-----------
#0  ../../../../../../aspect-release: aspect::TimeStepping::ConductionTimeStep<2>::execute()
#1  ../../../../../../aspect-release: aspect::TimeStepping::Manager<2>::update()
#2  ../../../../../../aspect-release: aspect::Simulator<2>::run()
#3  ../../../../../../aspect-release: void run_simulator<2>(std::cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, bool, bool, bool, bool)
#4  ../../../../../../aspect-release: main
--------------------------------------------------------

Aborting!

I use the following solver:

set Use conduction timestep       = true

set CFL number                    = 0.5
set Nonlinear solver scheme       = single Advection, single Stokes

###################################################
# SUBSECTION: Solver parameters
###################################################
subsection Solver parameters
  subsection Stokes solver parameters
    set Linear solver tolerance = 1e-5
    set Number of cheap Stokes solver steps = 2000
  end
end

Is there some way I can solve? For what I understood, i’m putting my tabulated EOS to “their limits” which is 6000 K for the temperature (even if, in this case, the max temperature is 5500 K because the planet cooled). Perplex can’t code EOS up to 6000 K, I’ve already tried (get NaN values). So I hope that there is some way to fix that changing solver, bit I’ve no idea if that is possible. Reducing the tolerance doesn’t help because the problem will happen again. I had some sims which where impossible to continue even with a “0” tolerance.

Do you have some suggestions?

Thank you in advance for the help.

@Francyrad - I think a good place to start is the motivation for using a conduction time step and whether similar studies examining LLSVPs using ASPECT also used this approach?

Reducing the tolerance doesn’t help because the problem will happen again. I had some sims which where impossible to continue even with a “0” tolerance.

Indeed, this is not the right strategy to address the issue. There are numerical instabilities in the simulation, likely resulting from unphysical behavior, which is causing a myriad of issues.

Following on the note above, 5550-6000K seems far to high for a CMB temperature, but this class of simulations is not my area of expertise.

Do you have some suggestions?

I would carefully compare your model setup with the setups in recent (last 5 years) papers using ASPECT to model LLSVPs.

Dear Francesco,

I agree. The present-day CMB should not be hotter than ~4500 K (and could be much colder). For more than 5500 K, you would expect a lot of melt, so modeling the mantle as a solid is not the right approach any more.

Best,
Juliane

Hello, I’m studying the thermal evolution of the earth with the dynamic core function. So i’m assuming that the Earth’s CMB was much hotter than today, so i was doing different trials to see how much CMB can cool using different temperatures

I took the solver from here:

I’ll try to set the conduction timestep on false to see what happens. As far as I know, no other people actually simulated LLSVPs using Steinberger (exception above). I’m not using other external libraries.