Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

brunt in rotation #372

Open
sunnywong314 opened this issue Feb 11, 2022 · 8 comments · May be fixed by #553
Open

brunt in rotation #372

sunnywong314 opened this issue Feb 11, 2022 · 8 comments · May be fixed by #553

Comments

@sunnywong314
Copy link
Collaborator

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 .

nu

@adamjermyn
Copy link
Collaborator

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
Copy link
Collaborator

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

@evbauer
Copy link
Member

evbauer commented Apr 5, 2022

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
Copy link
Collaborator

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

@orlox
Copy link
Contributor

orlox commented Apr 5, 2022

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
Copy link
Collaborator

adamjermyn commented Apr 5, 2022 via email

@warrickball
Copy link
Contributor

warrickball commented Apr 5, 2022

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
Copy link
Collaborator

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)...

@rjfarmer rjfarmer linked a pull request Jun 16, 2023 that will close this issue
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging a pull request may close this issue.

5 participants