Colors Bugfix: use Fe/H instead of zbase

colors/README.rst

@@ -22,7 +22,7 @@ How does the MESA colors module work?
 
 The module operates by coupling the stellar structure model with pre-computed grids of stellar atmospheres.
 
-1.  **Interpolation**: At each timestep, the module takes the star's current surface parameters—Effective Temperature (:math:`T_{\rm eff}`), Surface Gravity (:math:`\log g`), and Metallicity ([M/H])—and queries a user-specified library of stellar atmospheres (defined in ``stellar_atm``). It interpolates within this grid to construct a specific Spectral Energy Distribution (SED) for the stars current features.
+1.  **Interpolation**: At each timestep, the module takes the star's current surface parameters, Effective Temperature (:math:`T_{\rm eff}`), Surface Gravity (:math:`\log g`), and Metallicity ([Fe/H]), and queries a user-specified library of stellar atmospheres (defined in ``stellar_atm``). It interpolates within this grid to construct a specific Spectral Energy Distribution (SED) for the stars current features.
 
 2.  **Convolution**: This specific SED is then convolved with filter transmission curves (defined in ``instrument``) to calculate the flux passing through each filter.
 
@@ -57,7 +57,7 @@ stellar_atm
 
 Specifies the path to the directory containing the grid of stellar atmosphere models. This directory must contain:
 
-1.  **lookup_table.csv**: A map linking filenames to physical parameters (:math:`T_{\rm eff}`, :math:`\log g`, [M/H]).
+1.  **lookup_table.csv**: A map linking filenames to physical parameters (:math:`T_{\rm eff}`, :math:`\log g`, [Fe/H]).
 2.  **SED files**: The actual spectra (text or binary format).
 3.  **flux_cube.bin**: (Optional but recommended) A binary cube for rapid interpolation.
 
@@ -87,6 +87,51 @@ The distance to the star in centimeters.
 
     distance = 5.1839d20
 
+
+z_over_x_ref
+------------
+
+**Default:** ``0.0230057173972`` (the GS98 solar :math:`Z/X` ratio)
+
+Reference metal-to-hydrogen ratio used to convert the model photospheric
+composition into the metallicity coordinate passed to the atmosphere grid:
+
+.. math::
+
+   [Fe/H] = \log_{10}\left(\frac{(Z/X)_\mathrm{phot}}{(Z/X)_\mathrm{ref}}\right)
+
+The default is chosen to match the default ``stellar_atm =
+'/data/colors_data/stellar_models/Kurucz2003all/'`` grid, which is based on
+ATLAS9 models using the Grevesse & Sauval (1998) solar abundance scale. The
+numeric value ``0.0230057173972`` is computed from the MESA GS98 constants
+``Z_solar = 0.0169`` and ``Y_solar = 0.2485``, so
+``(Z/X)_ref = 0.0169 / (1 - 0.0169 - 0.2485)``.
+
+If you switch to a different atmosphere library, you must set ``z_over_x_ref``
+so that its reference composition matches the metallicity labels in that
+library's ``lookup_table.csv``. For custom abundance patterns, making this
+choice is the responsibility of the user.
+
+If the photospheric metallicity, expressed in terms of ``Z/X``, falls outside
+the metallicity range available in the tabulated atmosphere lookup table, then
+MESA uses the nearest metallicity in the table for the atmosphere
+interpolation. If the photospheric hydrogen mass fraction or metal mass
+fraction is not positive, then MESA cannot form
+``[Fe/H] = log10((Z/X)/(Z/X)_ref)``, and so it uses the lowest metallicity in
+the table instead. This is a fallback safeguard and not a dedicated treatment
+for H-free atmospheres.
+
+For the shipped ``Kurucz2003all`` table in this checkout, the tabulated
+metallicity range is ``[Fe/H] = -2.5`` to ``[Fe/H] = 4.0``. With the default
+``z_over_x_ref = 2.30057173972d-2``, the lower edge corresponds to an adopted
+``Z/X`` of about ``7.275d-5``.
+
+**Example:**
+
+.. code-block:: fortran
+
+      z_over_x_ref = 2.30057173972d-2
+
 make_csv
 --------
 
@@ -216,10 +261,10 @@ Below are the default values for the colors module parameters as defined in ``co
       vega_sed = '/data/colors_data/stellar_models/vega_flam.csv'
       stellar_atm = '/data/colors_data/stellar_models/Kurucz2003all/'
       distance = 3.0857d19  ! 10 parsecs in cm (Absolute Magnitude)
+      z_over_x_ref = 2.30057173972d-2  ! GS98 solar Z/X for Kurucz2003all
       make_csv = .false.
       colors_results_directory = 'SED'
       mag_system = 'Vega'
-      vega_sed = '/data/colors_data/stellar_models/vega_flam.csv'
 
 Visual Summary of Data Flow
 ===========================

colors/defaults/colors.defaults

@@ -13,8 +13,6 @@
    ! ~~~~~~~~~~~~~~
    ! ``vega_sed``
    ! ~~~~~~~~~~~~
-   ! ``stellar_atm``
-   ! ~~~~~~~~~~~~~~~
    ! ``distance``
    ! ~~~~~~~~~~~~
    ! ``make_csv``
@@ -24,36 +22,83 @@
    ! ``mag_system``
    ! ~~~~~~~~~~~~~~
 
-   ! If ``use_colors`` is true, the colors module is turned on, which will calculate
-   ! bolometric and synthetic magnitudes by interpolating stellar atmosphere model grids and convolving with photometric filter transmission curves.
-   ! Vega SED for Vega photometric system is used for photometric zero points.
-   ! Stellar distance is given in cm.
+   ! If ``use_colors`` is true, then MESA will calculate bolometric and
+   ! synthetic magnitudes by interpolating stellar atmosphere model grids and
+   ! convolving the resulting spectra with the selected photometric filters.
+   ! The Vega SED is used for photometric zero points when
+   ! ``mag_system = 'Vega'``. The stellar distance is given in cm.
+
    ! ::
 
       use_colors = .false.
       instrument = '/data/colors_data/filters/Generic/Johnson'
-      stellar_atm = '/data/colors_data/stellar_models/Kurucz2003all/'
       distance = 3.0857d19
       make_csv = .false.
       colors_results_directory = 'SED'
       mag_system = 'Vega'
       vega_sed = '/data/colors_data/stellar_models/vega_flam.csv'
 
 
+   ! ``stellar_atm``
+   ! ~~~~~~~~~~~~~~~
+
+   ! Path to the directory containing the stellar atmosphere lookup table and
+   ! spectra used by the colors module. This directory should contain a
+   ! ``lookup_table.csv`` file and the corresponding SED files. The default
+   ! uses the ``Kurucz2003all`` ATLAS9 atmosphere grid shipped with MESA.
+
+   ! ::
+
+      stellar_atm = '/data/colors_data/stellar_models/Kurucz2003all/'
+
+
+   ! ``z_over_x_ref``
+   ! ~~~~~~~~~~~~~~~~
+
+   ! ``z_over_x_ref`` is the reference metal-to-hydrogen ratio used to map
+   ! the photospheric composition onto the atmosphere-table metallicity axis:
+   ! ``[Fe/H] = log10((Z/X)/z_over_x_ref)``.
+
+   ! The default matches the GS98 solar mixture used by the default
+   ! ``Kurucz2003all`` ATLAS9 atmosphere grid and equals
+   ! ``0.0169/(1 - 0.0169 - 0.2485) = 2.30057173972d-2``.
+
+   ! If the photospheric metallicity, expressed as ``Z/X``, falls outside the
+   ! metallicity range available in the tabulated atmosphere lookup table,
+   ! then MESA uses the nearest metallicity in the table for the atmosphere
+   ! interpolation. If the photospheric hydrogen mass fraction or metal mass
+   ! fraction is not positive, then MESA cannot form
+   ! ``log10((Z/X)/z_over_x_ref)``, so it uses the lowest metallicity in the
+   ! table instead. This is a fallback safeguard and not a dedicated treatment
+   ! for H-free atmospheres.
+
+   ! For the shipped ``Kurucz2003all`` table in this checkout, the
+   ! metallicity range is ``[Fe/H] = -2.5`` to ``[Fe/H] = 4.0``. With the
+   ! default ``z_over_x_ref``, the lower edge corresponds to
+   ! Z/X ~ 7.275d-5.
+
+   ! ::
+
+      z_over_x_ref = 2.30057173972d-2
+
+
 ! Extra inlist controls
 ! ---------------------
 
    ! One can split a colors inlist into pieces using the following parameters.
    ! It works recursively, so the extras can read extras too.
 
 
-         ! ``read_extra_colors_inlist(1..5)``
-         ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-         ! ``extra_colors_inlist_name(1..5)``
-         ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+   ! ``read_extra_colors_inlist(1..5)``
+   ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+   ! ``extra_colors_inlist_name(1..5)``
+   ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
-         ! If ``read_extra_colors_inlist(i)`` is true, then read ``&colors`` from the file ``extra_colors_inlist_name(i)``.
-         ! ::
+   ! If ``read_extra_colors_inlist(i)`` is true, then read ``&colors``
+   ! from the file ``extra_colors_inlist_name(i)``.
+   !
+
+   ! ::
 
-      read_extra_colors_inlist(:) = .false.
-      extra_colors_inlist_name(:) = 'undefined'
+   read_extra_colors_inlist(:) = .false.
+   extra_colors_inlist_name(:) = 'undefined'

colors/private/colors_ctrls_io.f90

@@ -38,6 +38,7 @@ module colors_ctrls_io
    character(len=32) :: mag_system
 
    real(dp) :: distance
+   real(dp) :: z_over_x_ref
    logical :: make_csv
    logical :: use_colors
 
@@ -46,6 +47,7 @@ module colors_ctrls_io
       vega_sed, &
       stellar_atm, &
       distance, &
+      z_over_x_ref, &
       make_csv, &
       mag_system, &
       colors_results_directory, &
@@ -155,6 +157,7 @@ subroutine store_controls(rq, ierr)
       rq%vega_sed = vega_sed
       rq%stellar_atm = stellar_atm
       rq%distance = distance
+      rq%z_over_x_ref = z_over_x_ref
       rq%make_csv = make_csv
       rq%colors_results_directory = colors_results_directory
       rq%use_colors = use_colors
@@ -191,6 +194,7 @@ subroutine set_controls_for_writing(rq)
       vega_sed = rq%vega_sed
       stellar_atm = rq%stellar_atm
       distance = rq%distance
+      z_over_x_ref = rq%z_over_x_ref
       make_csv = rq%make_csv
       colors_results_directory = rq%colors_results_directory
       use_colors = rq%use_colors

colors/private/colors_history.f90

@@ -112,7 +112,8 @@ subroutine data_for_colors_history_columns( &
             filter_name = trim(remove_dat(color_filter_names(i)))
             names(i + filter_offset) = filter_name
 
-            if (t_eff >= 0 .and. metallicity >= 0) then
+            ! Negative [M/H] values are valid for metal-poor atmosphere grids.
+            if (t_eff >= 0) then
                ! Select precomputed zero-point based on magnitude system
                select case (trim(cs%mag_system))
                case ('VEGA', 'Vega', 'vega')
@@ -151,4 +152,4 @@ subroutine data_for_colors_history_columns( &
 
    end subroutine data_for_colors_history_columns
 
-end module colors_history
+end module colors_history

colors/public/colors_def.f90

@@ -46,6 +46,7 @@ module colors_def
       character(len=256) :: mag_system
       real(dp) :: metallicity
       real(dp) :: distance
+      real(dp) :: z_over_x_ref
       logical :: make_csv
       logical :: use_colors
       integer :: handle
@@ -206,4 +207,4 @@ subroutine do_free_colors_tables
 
    end subroutine do_free_colors_tables
 
-end module colors_def
+end module colors_def

star/private/history.f90

@@ -70,7 +70,8 @@ subroutine do_history_info(s, write_flag, ierr)
          num_extra_cols, num_binary_cols, num_extra_binary_cols, num_colors_cols,num_extra_header_items, n
       integer, parameter :: num_epsnuc_out = 12
       real(dp) :: &
-         epsnuc_out(num_epsnuc_out), csound_surf, v_surf, envelope_fraction_left
+         epsnuc_out(num_epsnuc_out), csound_surf, v_surf, envelope_fraction_left, m_div_h, &
+         min_m_div_h, max_m_div_h
       integer :: mixing_regions, mix_relr_regions, burning_regions, burn_relr_regions
       integer, pointer :: mixing_type(:), mix_relr_type(:), burning_type(:), burn_relr_type(:)
       character (len = maxlen_history_column_name), pointer, dimension(:) :: &
@@ -267,7 +268,27 @@ subroutine do_history_info(s, write_flag, ierr)
          end if
          colors_col_names(1:num_colors_cols) = 'unknown'
          colors_col_vals(1:num_colors_cols) = -1d99
-         call data_for_colors_history_columns(s%T(1), log10(s%grav(1)), s%R(1), s%kap_rq%Zbase, &
+
+         if (colors_settings% z_over_x_ref <= 0d0) then
+            write(*, *) 'colors error: z_over_x_ref must be positive'
+            ierr = -1
+            call dealloc
+            return
+         end if
+
+         min_m_div_h = minval(colors_settings% lu_meta)
+         max_m_div_h = maxval(colors_settings% lu_meta)
+
+         ! Map the current photospheric Z/X onto the atmosphere-table metallicity axis.
+         if (s% X(s% photosphere_cell_k) > 0d0 .and. s% Z(s% photosphere_cell_k) > 0d0) then
+            m_div_h = log10((s% Z(s% photosphere_cell_k)/s% X(s% photosphere_cell_k)) / &
+               colors_settings% z_over_x_ref)
+            m_div_h = max(min_m_div_h, min(max_m_div_h, m_div_h))
+         else
+            m_div_h = min_m_div_h
+         end if
+
+         call data_for_colors_history_columns(s%T(1), log10(s%grav(1)), s%R(1), m_div_h, &
             s% colors_handle, num_colors_cols, colors_col_names, colors_col_vals, ierr)
          if (ierr /= 0) then
             call dealloc