brunt in rotation

[sunnywong314]

In $MESA_DIR/star/private/rotation_mixing_info.f90, lines 595-596, the total Brunt and Brunt composition term are calculated again. I was wondering if it would be equivalent to simply replace those lines with s% brunt_N2 and s% brunt_N2_composition_term ?

The attached image shows a MESA profile comparing these, where I let a WD accrete at constant Mdot, turn Tayler-Spruit on only, and run with small timesteps (0.002 yr) until the model settles into a quasi-steady state. I extracted the necessary lines from rotation_mixing_info and put them as extra_profile_columns_data in my run_star_extras. In the middle panel, the N2_mu calculated (orange solid) is much noisier than the usual s% brunt_N2_composition_term (green solid). For this WD I averaged the inner 0.85 Msun abundances, so there should be no N2_mu in the inner 0.85 Msun.

This may have an impact on the Tayler-Spruit viscosity because it is calculated as the harmonic sum between the structure term and composition term. In the final panel, I take the shear profile from MESA (q in top panel), and then calculate the Tayler-Spruit viscosity, one using lines 595-596 (orange line) and one using brunt_N2_composition term from profile_columns.list (green). The former is agrees with the viscosity that the MESA profile gives (blue line), but the latter is what the viscosity should be as the inner 0.85 Msun has no composition gradient .

[adamjermyn]

Hi Sunny,

My initial reaction is that what you're proposing is a clear improvement and we should do it. It's also a straightforward change, which is nice.

My hesitation is that I know there's funny business in MESA very specifically in terms of when we calculate/update the Brunt and when we calculate/update mixing coefficients, so the two disagreeing could be quite bad.

At minimum it seems like we should replace the code in those lines with a call to the routines in brunt.f90 that calculate the brunt (and we should, along the way, factor those so that they expose a 'pure' version that doesn't try to write the brunt into the star pointer).

I'm tagging @orlox, who knows more about the rotation guts than I do, and @evbauer, who was very recently working on stuff in brunt.f90.

[adamjermyn]

If no one objects I'd like to go ahead and make this change... @evbauer? @orlox?

[evbauer]

Yes, I think this makes sense. I think the difference mostly comes down to the inclusion of smooth_brunt_B in brunt.f90, so just be aware that this has had some smoothing applied rather than being the raw Brunt from the structure of the model. But that's probably a reasonable thing to use here?

[adamjermyn]

Ok, I've implemented this in the branch brunt_in_rotation. If that passes I'll merge it in.

[orlox]

This ok for me as well, though it would be ideal to do a run with fpe checks. In the way rotation is split up, with various things computed at specific stages, it can be easy to introduce fpe issues.

[adamjermyn]

Do we currently pass on main with FPE checks turned on? If so very happy to trigger a run… otherwise I'm not quite sure what I'd be looking for.

-Adam

On Tue, Apr 05, 2022 at 1:45 PM, Pablo Marchant @.***> wrote:

This ok for me as well, though it would be ideal to do a run with fpe checks. In the way rotation is split up, with various things computed at specific stages, it can be easy to introduce fpe issues.

— Reply to this email directly, view it on GitHub https://github.com/MESAHub/mesa/issues/372#issuecomment-1089107507, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABPR5H2KMN7PKH55CA6WPR3VDR347ANCNFSM5OFGODNQ . You are receiving this because you commented.Message ID: @.***>

[warrickball]

Do we currently pass on main with FPE checks turned on?

No. I wrote about some results in #bugs-and-problems on April 1.

[adamjermyn]

Ok. In that case I don't think we should wait for FPE checks to pass.

However in the meantime I see some failures on testhub. Looks like 1M_thermohaline and ppisn. On 1M_thermohaline the failure is:

                             log total_angular_momentum    4.9034800969549750D+01    4.5000000000000000D+01    5.5000000000000000D+01
                                       log center_omega   -3.2919413454258080D+00   -4.0000000000000000D+00   -2.0000000000000000D+00
                                      log he_core_omega   -4.2901355684847662D+00   -5.0000000000000000D+00   -2.0000000000000000D+00
                                          surface j_rot    1.6466216760682634D+01    5.0000000000000000D+00    2.5000000000000000D+01
                                          surface v_rot    8.1786942180389999D-01    0.0000000000000000D+00    1.0000000000000000D+00

                                    avg near 0.245 Msun
                                       *** BAD *** logT    6.4685621317362161D+00    7.0000000000000000D+00    7.5000000000000000D+00
                                     *** BAD *** logRho   -1.2418988442735306D+00    5.0000000000000000D-01    2.0000000000000000D+00
                                              log j_rot    1.4353261037447641D+01    1.0000000000000000D+01    2.0000000000000000D+01
                                       *** BAD *** D_ST   -9.8999999999999986D+01    1.0000000000000000D+00    8.0000000000000000D+00
                                      *** BAD *** nu_ST   -9.8999999999999986D+01    4.0000000000000000D+00    9.0000000000000000D+00

which to me looks an awful lot like these fields aren't filled and give NaN's (which get written out as -99)...