reconstruct eos/kap on faces for mlt/tdc

[Debraheem]

This branch was birthed out of doing some testing with tdc pulsations, where I was benchmarking where TDC decoupled versus RSP.

At the moment, TDC and mlt used P and kap and other quantites which are averaged onto cell faces. Probably more thermodynamically consistent is to use the face averaged T and rho to recalculate the opacity or eos, which is what use_face_values_eos_and_kap_mlt_tdc does in this pr.

Claim: using use_face_values_eos_and_kap_mlt_tdc allows mesa to converge significantly better on courser mesh, when the opacity or other eos quantity is not well resolved at the cell face.

Example:

Below i show a 9M cepheid which would be fine in RSP, because RSP is cell centered, but in MESA-TDC features in the mesh from a poorly interpolated opacity cause ringing in kappa -> gradr -> L -> Lc -> (unphysical oscillations in the pulsations). By using use_face_values_eos_and_kap_mlt_tdc, these features are eliminated. The only other solution I could find was to double or triple the mesh resolution in MESA to get similar smoothness to RSP at this particular resolution, and sometimes that wasn't enough.

---

With use_face_values_eos_and_kap_mlt_tdc the face_averaged primitives T and rho etc, are used to computes the face eos and kap quantities using an additional eos call, instead of averaging the Pressure and opacity and other thermodynamic quantities on cell faces. This generates a more accurate reconstruction at the cost of an extra eos call. I refactored things such that we minimize the number of times eos calls are made though. This is done by storing the face values in global star pointers.

This is a dev control for now, but if it is useful for low resolution meshes, it could be potentially useful as a global default, pending testing and time.

The TDC_Cepheid_9M test case can be removed before this branch is merged, but it serves as a test case for now. (We can squash this pr as well, if we merge.)

[Debraheem]

Here is a graphic of what this pr is attempting to do:

[Debraheem]

A few more extreme examples of what undersampling an opacity bump might look like. The point is that, this undersampling can happen with any reconstructed quantity be it from the opacity or any of the many eos quantities including the face pressure. Having a mesh that accurately resolves all of these is much more complicated, and potentially expensive than reconstructing them self consistently. However, In the limit of infinite mesh density, the old averaging procedure should reproduce the results of this new approach.

[Debraheem]

This is all made worse when the model is not hydrostatic, and a fixed lagrangian mass shell moves back and forth across one of these bumps, in which case the model can switch back and forth between under to over predicting the opacity depending on the location of the cells with respect to a feature in the opacity or any other quantity. My example I shared above of a pulsating Cepheid is one where the opacity being underesolved leads to oscillations kap_face -> gradr -> oscillations in the Lc and Lr coming out mlt/tdc.

Obviously remesh can help avoid this, but if you're remeshing scheme is not helping conserve a particular eos/kap quantity of interest, you might still run into this issue and not know it.

[fxt44]

oohh, i like this idea of interpolating the (T, rho) primitives, then calling the microphysics.

[warrickball]

Fundamentally I like this change too. I appreciate that there must be some cost in additional EoS, opacity, etc. calls but I expect the self-consistency of the models (incl. the equilibrium models for pulsation codes) is well worth it. We might even be able to achieve the same fidelity with coarser meshes, longer timesteps and fewer solver iterations, and therefore even net some run time back.

I'm trying to digest the code in what little time I have but I'll echo someone else's point (maybe Mattias?) that it isn't clear why a lot of new face variables have the mlt_tdc_ prefix. E.g. The value of $c_P$ on a face would reasonably be called Cp_face for any purpose. Whether or not it's used in TDC/MLT seems irrelevant?

[Debraheem]

Fundamentally I like this change too. I appreciate that there must be some cost in additional EoS, opacity, etc. calls but I expect the self-consistency of the models (incl. the equilibrium models for pulsation codes) is well worth it. We might even be able to achieve the same fidelity with coarser meshes, longer timesteps and fewer solver iterations, and therefore even net some run time back.

I'm trying to digest the code in what little time I have but I'll echo someone else's point (maybe Mattias?) that it isn't clear why a lot of new face variables have the mlt_tdc_ prefix. E.g. The value of c P on a face would reasonably be called Cp_face for any purpose. Whether or not it's used in TDC/MLT seems irrelevant?

Thanks for the feedback! I agree on the naming. This was intentionally done, since there are still many places in mesa where face means wrapping rather than reconstructing and I wanted to be clear about which ones were reconstructed while building this branch. The prefix mlt_tdc indicates that it is a variable which is only really used in turb. We can definitely change the naming, but it cannot just be _face until I do another full sweep across mesa and ensure it always means the same thing. Otherwise _face means different things in different routines and it would be confusing to implement. at the moment some routines like rotation still build their own scale height and brunt term, so there are still some inconsistencies in MESA between routines (predating this branch). Perhaps this pr is an opportunity for me to ensure mesa is using the same quantities everywhere. There is still some work to be done!