EquL solver stability improvement #897
Conversation
…dd temporary fill arrays with nans inside magnetic braking
Debraheem
requested review from
annethoul,
aurimontem,
evbauer,
fxt44,
mjoyceGR,
rhdtownsend and
wmwolf
as code owners
Debraheem
changed the title
|
see https://testhub.mesastar.org/EbF%2Ffix_EquL_solver_stability/commits/5f32328 for testing info. |
|
beautiful! |
|
I'm under the impression the scaling has almost no effect and is rarely if ever toggled, rather the adjustments to the solver jump limits and making gradr_sticky provide the stability. |
|
lnTdiff_start is as far as I can tell always of order <1d-2, so the scaling is always 1d0. In that case, I case scale is always 1d0, so I'll probably scratch the residual scaling. |
| lnTdiff = dT/Tpoint ! use this in place of lnT(k-1)-lnT(k) | ||
| delm = (s% dm(k) + s% dm(k-1))/2 | ||
|
|
||
| resid = delm*dlnTdm - lnTdiff |
There was a problem hiding this comment.
I'm confused. This seems to be the only change left in hydro_temperature? So nothing changed but the inlist toggles?
There was a problem hiding this comment.
Ok, then I'd advocate of reverting even this change (its just brackets), so it's clear only inlists changed.
There was a problem hiding this comment.
ahh yes, I missed the brackets. good call.
This is either a legendary fix or a blunder. I believe I've addressed the issue causing us to have to resort to using
convergence_ignore_equL_residuals = .true..If you're curious look at the commit.
edit *** residual scaling seems to have no effect and was removed.Just light some slight tweaks to the solver tolerances for when to give up and changingconvergence_ignore_equL_residuals = .true.across the entire test_suite to.false., and changed the MESA default foruse_gold2_tolerances = .true., as we were stuck using gold_tolerances before. All the models seem fine now, except the two planets test_cases, which already resorted to this, and the pisn test case which uses a different form of the T-gradient equation anyway.This fixes issue like #563. With this pr, I can run through the entire He flash with
convergence_ignore_equL_residuals = .false.and gold2_tolerances.see also old issues: #80, #490
For example the 1M_pre_ms_to_wd test_suite now completes the he flash with LogE_err < -5, roughly 3 orders of magnitude lower than previous commits. (might just be because i switched from eps_grav to dedt form of E-equation)
(I also made some improvements to the make_co_wd test_case, substantially improving its robustness, as it was extremely sensitive)
The largest caveat is my setting of
make_gradr_sticky_in_solver_iters = .true.by default. This was already on all the massive star cc test_cases and some others. This solver hack waits until the 4th solver iteration, and then after that point, if the zone was convective and becomes radiative (not the opposite), the solver will hold gradr fixed to prevent the sign from flipping again. In practice it has almost no observable effect outside of extremely short timescales, when burning (dlnE_dt) is competing with gradT (from mlt) through the T-gradient equation. We can set this back to.false.to see where failures show up, and in practice it could remain.false.. I personally like the idea of keeping it on by default, but I will submit a test commit with it set back to false to see if it has any impact and where.please verify my changes, give the pr a try.
This pr still needs to address the photo issue in magnetic breaking, likely related to some Spruit Taylor initialization bug.