Skip to content

Eos MESA updated tables + entropy tables for GRTDE#850

Open
AliPourmand wants to merge 45 commits into
danieljprice:mainfrom
AliPourmand:eos_mesa_entropy_tde
Open

Eos MESA updated tables + entropy tables for GRTDE#850
AliPourmand wants to merge 45 commits into
danieljprice:mainfrom
AliPourmand:eos_mesa_entropy_tde

Conversation

@AliPourmand

@AliPourmand AliPourmand commented Jul 8, 2026

Copy link
Copy Markdown
Contributor

Description:
I have updated EOS MESA tables to be made with version: mesa-24.08.1
keeping EOS MESA tables of Reichardt 2020 as a hidden option that can be modified.

Additionally, I have changed cons2primsolver so that now the gr code can use eos_mesa by adding new tables tabulated with entropy. This has been tested and it works.

Along with bug fixes and changes in datafile among others.
Components modified:

  • Setup (src/setup)
  • [ x] Main code (src/main)
  • Moddump utilities (src/utils/moddump)
  • Analysis utilities (src/utils/analysis)
  • Test suite (src/tests)
  • Documentation (docs/)
  • Build/CI (build/ or github actions)

Type of change:

  • [ x] Bug fix
  • Physics improvements
  • Better initial conditions
  • Performance improvements
  • Documentation update
  • Better testing
  • Code cleanup / refactor
  • Other (please describe)

Testing:
By running sample eos mesa stars and comparing with old runs.
By running eos_mesa runs with grtde

Did you run the bots? no

Did you update relevant documentation in the docs directory? no

Did you add comments such that the purpose of the code is understandable? yes

I

Summary by CodeRabbit

  • New Features
    • Added support for an additional GR-specific MESA EOS table set, enabling pressure/temperature and internal-energy lookups from density+entropy.
    • Extended entropy-based EOS recovery to support ieos=10.
    • Introduced eosmesa_version to select the appropriate Zenodo-hosted MESA EOS table set.
  • Bug Fixes
    • Improved GR compatibility by accepting eos_type=10 and broadening unsupported-EOS error messaging.
    • Fixed MESA entropy handling to use the table’s entropy value with correct cgs unit conversion.
    • Refined GR EOS interpolation/extrapolation and updated opacity table selection to use the correct dataset key.

…sity for GR, I want to change as little as possible though and preferably get away with GR users setting the same EOS number too
…d now, in eos, eos_mesa and eos_mesa_microphysics, now I need to add the entropy tables here, and also see if they get allocated correctly, I am only adding 4 variables in these tables so they occupy less size, other than that the compositions will be the same, but the minimum density should go as low as 1e-25 g/cm as fitz' simulations show. also this is not tested so I should do that
…adly I don't have fitz to double check it rn, but I think it should work. Also, I still need to figure out why my tables aren't being read.
… there is an issue probably with the tables themselves which I'll need to test
…sa tables, because they are specific entropy divided by avo*kboltz
… in log, so you need to specify which are and areny
…ed in eosvars but now has to somehow go to be used in relaxstar
…es for eosmesa, this well affect classical star relaxing as well as grstar and make them dependent on entropy
…also to see whehter Daniel is happy with this at all
…axation because the entropy at that stage shouldn't be tampered with
…o false by default to use thomas reichardt's old tables
…anged opacs table to be read separately for eos mesa as I have not updated them so they would just be downloaded from Thomas old link
…osmesa_version. also the gr tables will be uploaded in the same place that the classic tables will be
@coderabbitai

coderabbitai Bot commented Jul 8, 2026

Copy link
Copy Markdown

Review Change Stack

Note

Reviews paused

It looks like this branch is under active development. To avoid overwhelming you with review comments due to an influx of new commits, CodeRabbit has automatically paused this review. You can configure this behavior by changing the reviews.auto_review.auto_pause_after_reviewed_commits setting.

Use the following commands to manage reviews:

  • @coderabbitai resume to resume automatic reviews.
  • @coderabbitai review to trigger a single review.

Use the checkboxes below for quick actions:

  • ▶️ Resume reviews
  • 🔍 Trigger review
📝 Walkthrough

Walkthrough

This PR adds GR MESA EOS support for ieos=10, including new table accessors, interpolation paths, solver wiring, and version-based selection of MESA EOS data URLs.

Changes

MESA GR EOS Support

Layer / File(s) Summary
GR MESA table storage and interpolation routines
src/main/eos_mesa_microphysics.f90
Adds GR grid parameters and storage, GR table loading, GR value lookup, GR spline interpolation, deallocation, and updated opacity datafile keys.
GR initialization and accessor subroutines in eos_mesa
src/main/eos_mesa.f90
Imports the GR runtime flag, conditionally initializes GR EOS tables in init_eos_mesa, and adds GR pressure/temperature and internal-energy accessors.
EOS module public API and inversion routines
src/main/eos.f90
Extends public exports, accepts eos_type=10, updates MESA entropy handling, extends pressure inversion, and adds internal-energy recovery helpers for ieos=10.
cons2primsolver ieos=10 wiring
src/main/cons2primsolver.f90
Imports get_u_from_rho_s and adds ieos==10 branches in the entropy-based Newton iteration setup and post-convergence handling, plus the updated fatal message.
MESA EOS table version selection
src/main/datafiles.f90
Adds eosmesa_version and updates map_dir_to_web to choose the MESA EOS table URL by version, with a stop on unrecognized values.

Estimated code review effort: 4 (Complex) | ~60 minutes

🚥 Pre-merge checks | ✅ 5
✅ Passed checks (5 passed)
Check name Status Explanation
Title check ✅ Passed The title clearly summarizes the main change: updated MESA tables and added entropy tables for GR use.
Description check ✅ Passed The description covers the main sections and testing, though the unit-test and related-issues prompts are left incomplete.
Docstring Coverage ✅ Passed No functions found in the changed files to evaluate docstring coverage. Skipping docstring coverage check.
Linked Issues check ✅ Passed Check skipped because no linked issues were found for this pull request.
Out of Scope Changes check ✅ Passed Check skipped because no linked issues were found for this pull request.
✨ Finishing Touches
🧪 Generate unit tests (beta)
  • Create PR with unit tests

Thanks for using CodeRabbit! It's free for OSS, and your support helps us grow. If you like it, consider giving us a shout-out.

❤️ Share

Comment @coderabbitai help to get the list of available commands.

@gemini-code-assist

Copy link
Copy Markdown
Contributor

Warning

Gemini encountered an error creating the review. You can try again by commenting /gemini review.

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 7

Caution

Some comments are outside the diff and can’t be posted inline due to platform limitations.

⚠️ Outside diff range comments (2)
src/main/eos.f90 (1)

1223-1235: 🎯 Functional Correctness | 🟡 Minor | ⚡ Quick win

Stale fatal message: ieos=10 is now handled but not mentioned.

The case default fatal in get_p_from_rho_s still reads "only implemented for eos 2 and 12", even though case(10) was just added above. The analogous message in cons2primsolver.f90 (Line 215) was correctly updated to mention eos 10; this one was missed.

🐛 Proposed fix
- call fatal('eos','[get_p_from_rho_s] only implemented for eos 2 and 12')
+ call fatal('eos','[get_p_from_rho_s] only implemented for eos 2, 10 and 12')
🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/eos.f90` around lines 1223 - 1235, Update the stale fatal message in
get_p_from_rho_s so it matches the actual supported EOS cases in the select case
on ieos. Since case(10) now routes to get_eos_ptemp_from_rhos_mesa_gr, change
the case default call to fatal to mention eos 10 as supported as well, mirroring
the updated message used in cons2primsolver.f90.
src/main/cons2primsolver.f90 (1)

194-219: 🚀 Performance & Scalability | 🟠 Major | ⚡ Quick win

Redundant per-iteration MESA GR table lookup for ieos==10.

Inside the Newton loop, case(10) calls get_u_from_rho_s every iteration, but its result (u) is never used before being overwritten next iteration — for the entropy_s branch, enth is updated purely from gamfac*P/dens (Line 236), and gamfac is only ever refreshed inside case(12), never case(10). Contrast this with case(12), where u is genuinely needed to update gamma/gamfac for the next iteration.

The final u for ieos==10 is already (correctly) recomputed once after convergence (Lines 264-267), so the in-loop call just adds an unnecessary table interpolation on every Newton iteration for every GR-MESA particle.

⚡ Proposed fix
        select case(ieos)
        case (10)
-         ! inputs and outputs are all in code units
-         call get_u_from_rho_s(ieos,en,dens,u)
-
+         ! u is not needed until after convergence; computed post-loop
+         ! to avoid redundant per-iteration table lookups
        case (12)
🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/cons2primsolver.f90` around lines 194 - 219, The Newton-loop branch
in cons2primsolver for ien_entropy_s / case (10) is doing a redundant
get_u_from_rho_s lookup every iteration even though u is not used to update enth
or gamfac there; only case (12) needs u to refresh gamma and gamfac. Remove the
per-iteration get_u_from_rho_s call from the ieos==10 path in the select case
block and keep the single post-convergence recomputation of u for ieos==10 so
the final value is still set correctly.
🧹 Nitpick comments (1)
src/main/eos_mesa.f90 (1)

150-171: 🩺 Stability & Availability | 🔵 Trivial | ⚡ Quick win

Missing extrapolation/error signalling for the new GR accessors.

get_eos_ptemp_from_rhos_mesa_gr and get_eos_u_from_rhos_mesa_gr call getvalue_mesa_gr without an ierr output, unlike get_eos_pressure_temp_gamma1_mesa, whose ierr is used in eos.f90 (case 10) to warn when "extrapolating off tables". Silent off-table extrapolation in the GR hot loop (called every Newton iteration in cons2primsolver.f90) could produce unphysical results without any diagnostic.

If getvalue_mesa_gr supports an optional ierr argument (mirroring getvalue_mesa), consider threading it through here and up through get_p_from_rho_s/get_u_from_rho_s in eos.f90.

🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/eos_mesa.f90` around lines 150 - 171, The new GR accessors
`get_eos_ptemp_from_rhos_mesa_gr` and `get_eos_u_from_rhos_mesa_gr` currently
ignore table-boundary failures by calling `getvalue_mesa_gr` without surfacing
any `ierr` status. Thread an optional error/status output through these
routines, mirroring the existing `get_eos_pressure_temp_gamma1_mesa` handling,
and propagate it upward through `get_p_from_rho_s` and `get_u_from_rho_s` in
`eos.f90` so off-table extrapolation can be warned about instead of happening
silently in the GR solve path.
🤖 Prompt for all review comments with AI agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

Inline comments:
In `@src/main/datafiles.f90`:
- Around line 70-77: The invalid eosmesa_version branch in datafiles.f90
currently uses a plain stop, which can strand other MPI ranks at barrier_mpi().
Replace the stop path in the select case around eosmesa_version with the
project’s MPI-aware error/abort handling so all ranks exit cleanly; keep the
change localized to the version switch and use the existing MPI
error/termination helper already used elsewhere in the codebase.

In `@src/main/eos_mesa_microphysics.f90`:
- Around line 850-856: `deallocate_arrays_mesa` is missing cleanup for the GR
composition arrays, so update the deallocation block to also release
`mesa_eos_gr_z` and `mesa_eos_gr_h` alongside the existing GR arrays. Use the
same allocation checks already used for `mesa_eos_gr0`, `mesa_eos_gr_logss`, and
the other `mesa_eos_gr_*` symbols, and keep the GR teardown consistent with the
non-GR cleanup in `deallocate_arrays_mesa`.
- Around line 683-693: In getvalue_mesa_gr, the GR fallback branch is
exponentiating all columns, which incorrectly transforms ivout=4 (Gamma1)
outside the table; update the fallback logic to keep column 4 in linear form,
matching the existing spline path where only the first three log10 columns are
converted back. Also update deallocate_arrays_mesa to release the GR grid arrays
allocated by get_eos_constants_mesa_gr, specifically mesa_eos_gr_z and
mesa_eos_gr_h, alongside the other mesa_eos_gr_* arrays.

In `@src/main/eos_mesa.f90`:
- Around line 144-171: The header comments for get_eos_ptemp_from_rhos_mesa_gr
and get_eos_u_from_rhos_mesa_gr still say "NOT TESTED YET!", which is now stale.
Update or remove those banner comments so they match the current tested status
described for the new GR-table EOS routines, keeping the documentation in
eos_mesa accurate for future readers.
- Around line 61-67: The GR initialization path in eos_mesa.f90 can hide a
failure from get_eos_constants_mesa_gr because read_eos_mesa_gr may reset ierr
back to zero. Update the gr branch in eos_mesa so it checks ierr immediately
after calling get_eos_constants_mesa_gr and returns at once on failure, before
invoking read_eos_mesa_gr; keep the existing final ierr guard only for the read
step.

In `@src/main/eos.f90`:
- Around line 167-171: The fatal error message in the GR compatibility check is
out of sync with the allowed `eos_type` values in `eos` because it omits
`ieos=4`. Update the `fatal` call in the `eos` logic so the message matches the
expanded `any((/2,4,10,11,12/)==eos_type)` check and explicitly includes `4`
alongside the other supported values.
- Around line 1280-1290: get_u_from_rho_s_vec is missing the same caller-size
guard that get_entropy_vec uses, so a mismatched n can cause out-of-bounds
access. Update get_u_from_rho_s_vec to validate that n matches the actual sizes
of S, rho, and u before the loop, and call fatal with a clear message if they
differ. Use the existing get_entropy_vec size-check pattern and keep the fix
localized to get_u_from_rho_s_vec/get_u_from_rho_s.

---

Outside diff comments:
In `@src/main/cons2primsolver.f90`:
- Around line 194-219: The Newton-loop branch in cons2primsolver for
ien_entropy_s / case (10) is doing a redundant get_u_from_rho_s lookup every
iteration even though u is not used to update enth or gamfac there; only case
(12) needs u to refresh gamma and gamfac. Remove the per-iteration
get_u_from_rho_s call from the ieos==10 path in the select case block and keep
the single post-convergence recomputation of u for ieos==10 so the final value
is still set correctly.

In `@src/main/eos.f90`:
- Around line 1223-1235: Update the stale fatal message in get_p_from_rho_s so
it matches the actual supported EOS cases in the select case on ieos. Since
case(10) now routes to get_eos_ptemp_from_rhos_mesa_gr, change the case default
call to fatal to mention eos 10 as supported as well, mirroring the updated
message used in cons2primsolver.f90.

---

Nitpick comments:
In `@src/main/eos_mesa.f90`:
- Around line 150-171: The new GR accessors `get_eos_ptemp_from_rhos_mesa_gr`
and `get_eos_u_from_rhos_mesa_gr` currently ignore table-boundary failures by
calling `getvalue_mesa_gr` without surfacing any `ierr` status. Thread an
optional error/status output through these routines, mirroring the existing
`get_eos_pressure_temp_gamma1_mesa` handling, and propagate it upward through
`get_p_from_rho_s` and `get_u_from_rho_s` in `eos.f90` so off-table
extrapolation can be warned about instead of happening silently in the GR solve
path.
🪄 Autofix (Beta)

Fix all unresolved CodeRabbit comments on this PR:

  • Push a commit to this branch (recommended)
  • Create a new PR with the fixes

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro Plus

Run ID: 11d48c7b-d454-4506-af1e-6fff3d5f9582

📥 Commits

Reviewing files that changed from the base of the PR and between 6a853b3 and b288c9c.

📒 Files selected for processing (5)
  • src/main/cons2primsolver.f90
  • src/main/datafiles.f90
  • src/main/eos.f90
  • src/main/eos_mesa.f90
  • src/main/eos_mesa_microphysics.f90

Comment thread src/main/datafiles.f90
Comment thread src/main/eos_mesa_microphysics.f90
Comment thread src/main/eos_mesa_microphysics.f90 Outdated
Comment thread src/main/eos_mesa.f90 Outdated
Comment thread src/main/eos_mesa.f90 Outdated
Comment thread src/main/eos.f90
Comment thread src/main/eos.f90 Outdated

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Caution

Some comments are outside the diff and can’t be posted inline due to platform limitations.

⚠️ Outside diff range comments (1)
src/main/eos.f90 (1)

1222-1233: 🎯 Functional Correctness | 🟡 Minor | ⚡ Quick win

Update the unsupported-EOS fatal message.

get_p_from_rho_s now implements 2, 5, 10, and 12, but the fatal text still says only 2 and 12.

🐛 Proposed fix
-    call fatal('eos','[get_p_from_rho_s] only implemented for eos 2 and 12')
+    call fatal('eos','[get_p_from_rho_s] only implemented for eos 2, 5, 10 and 12')
🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/eos.f90` around lines 1222 - 1233, The fatal message in
get_p_from_rho_s is outdated and still claims only EOS 2 and 12 are supported,
even though the select case now handles 2, 5, 10, and 12. Update the
default-case fatal call in eos.f90 so the error text accurately lists the
supported EOS values, keeping the message aligned with the select case branches
in get_p_from_rho_s.
🤖 Prompt for all review comments with AI agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

Outside diff comments:
In `@src/main/eos.f90`:
- Around line 1222-1233: The fatal message in get_p_from_rho_s is outdated and
still claims only EOS 2 and 12 are supported, even though the select case now
handles 2, 5, 10, and 12. Update the default-case fatal call in eos.f90 so the
error text accurately lists the supported EOS values, keeping the message
aligned with the select case branches in get_p_from_rho_s.

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro Plus

Run ID: 8dd0020c-cf90-4f33-ae6f-75f1bb733bd9

📥 Commits

Reviewing files that changed from the base of the PR and between b288c9c and a55be3e.

📒 Files selected for processing (1)
  • src/main/eos.f90

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Caution

Some comments are outside the diff and can’t be posted inline due to platform limitations.

⚠️ Outside diff range comments (2)
src/main/eos.f90 (2)

1210-1213: 🎯 Functional Correctness | 🟡 Minor | ⚡ Quick win

Stale fatal message: case default now omits ieos=5 and ieos=10.

The select case now handles 2, 5, 10, and 12, but the fallthrough message still reads "only implemented for eos 2 and 12".

🐛 Proposed fix
-    call fatal('eos','[get_p_from_rho_s] only implemented for eos 2 and 12')
+    call fatal('eos','[get_p_from_rho_s] only implemented for eos 2, 5, 10 and 12')
🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/eos.f90` around lines 1210 - 1213, The fallback fatal message in
get_p_from_rho_s is stale because the select case now supports more EOS values
than the message says. Update the case default branch in the get_p_from_rho_s
logic to reflect the current supported set handled by the surrounding select
case, namely eos 2, 5, 10, and 12, so the fatal text is accurate and consistent.

1196-1213: 🎯 Functional Correctness | 🟠 Major | ⚡ Quick win

Restore the missing kboltz factor for GR entropy
get_entropy already stores S as entropy divided by kboltz, so cgss = S*unit_ergg passes entropy_cgs/kboltz into the GR MESA wrappers. getvalue_mesa_gr expects cgs entropy and converts it with s/(kboltz*avogadro), so this path is off by another factor of kboltz. Apply the same fix in get_u_from_rho_s too.

🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/eos.f90` around lines 1196 - 1213, The GR entropy conversion is
missing the kboltz factor in the EOS-to-MESA path, so the entropy passed from
get_p_from_rho_s (and the matching get_u_from_rho_s helper) is still scaled
incorrectly before calling the GR wrappers. Update the cgs entropy calculation
near cgss in both routines to restore the entropy units expected by
get_eos_ptemp_from_rhos_mesa_gr / getvalue_mesa_gr, and keep the fix aligned
with how get_entropy stores S and how the ieos = 10 branch handles GR
conversions.
🧹 Nitpick comments (1)
src/main/eos.f90 (1)

1224-1255: 📐 Maintainability & Code Quality | 🔵 Trivial | 💤 Low value

Doc header describes temperature, but the routine returns internal energy u.

The banner says "Calculate temperature given density and entropy", yet get_u_from_rho_s outputs u. Correct the comment to avoid confusion. (Units concern already raised at the get_p_from_rho_s root cause.)

🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/eos.f90` around lines 1224 - 1255, The header comment for
get_u_from_rho_s is misleading because it says the routine calculates
temperature, but the subroutine actually returns internal energy u. Update the
banner/docstring near get_u_from_rho_s to describe internal energy from density
and entropy, and keep the implementation and unit conversion logic unchanged.
🤖 Prompt for all review comments with AI agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

Outside diff comments:
In `@src/main/eos.f90`:
- Around line 1210-1213: The fallback fatal message in get_p_from_rho_s is stale
because the select case now supports more EOS values than the message says.
Update the case default branch in the get_p_from_rho_s logic to reflect the
current supported set handled by the surrounding select case, namely eos 2, 5,
10, and 12, so the fatal text is accurate and consistent.
- Around line 1196-1213: The GR entropy conversion is missing the kboltz factor
in the EOS-to-MESA path, so the entropy passed from get_p_from_rho_s (and the
matching get_u_from_rho_s helper) is still scaled incorrectly before calling the
GR wrappers. Update the cgs entropy calculation near cgss in both routines to
restore the entropy units expected by get_eos_ptemp_from_rhos_mesa_gr /
getvalue_mesa_gr, and keep the fix aligned with how get_entropy stores S and how
the ieos = 10 branch handles GR conversions.

---

Nitpick comments:
In `@src/main/eos.f90`:
- Around line 1224-1255: The header comment for get_u_from_rho_s is misleading
because it says the routine calculates temperature, but the subroutine actually
returns internal energy u. Update the banner/docstring near get_u_from_rho_s to
describe internal energy from density and entropy, and keep the implementation
and unit conversion logic unchanged.

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro Plus

Run ID: cf09ade9-bbf1-4d5b-a0db-cabad4f04136

📥 Commits

Reviewing files that changed from the base of the PR and between 1365a7a and 4149bea.

📒 Files selected for processing (1)
  • src/main/eos.f90

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Caution

Some comments are outside the diff and can’t be posted inline due to platform limitations.

⚠️ Outside diff range comments (2)
src/main/eos_mesa_microphysics.f90 (1)

559-569: 🩺 Stability & Availability | 🟡 Minor | ⚡ Quick win

Fix the GR X-selection indices
mesa_eos_gr_data_exists(i, ...) uses the exhausted file-loop index, so this reads past the first dimension; j-2 is also 0 when j==2. Use the selected Z index (nz2) and a 1-based column index here.

🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/eos_mesa_microphysics.f90` around lines 559 - 569, The GR
X-selection logic in mesa_eos_microphysics should stop indexing
mesa_eos_gr_data_exists with the exhausted file-loop variable and an invalid
zero column. Update the nx2 selection block to use the chosen Z index (nz2)
instead of i, and ensure the column lookup is 1-based when checking adjacent
entries in the mesa_eos_gr_data_exists array. Keep the fix localized to the do j
loop that computes nx2 so the boundary cases at j==2 and j==mesa_eos_gr_nh
remain valid.
src/main/eos.f90 (1)

1205-1213: 🎯 Functional Correctness | 🟡 Minor | ⚡ Quick win

Update stale fatal message in get_p_from_rho_s.

The routine now handles ieos 2, 5, 10 and 12, but the default-branch message still reads "only implemented for eos 2 and 12", which will mislead users hitting an unsupported ieos.

🐛 Proposed fix
- call fatal('eos','[get_p_from_rho_s] only implemented for eos 2 and 12')
+ call fatal('eos','[get_p_from_rho_s] only implemented for eos 2,5,10 and 12')
🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

In `@src/main/eos.f90` around lines 1205 - 1213, Update the stale fatal message in
get_p_from_rho_s so it matches the actual ieos cases handled by the select case.
The default branch currently says only eos 2 and 12 are implemented, but the
routine also supports 5 and 10; revise the fatal('eos', ...) text in that
function to accurately list all supported values so users get the correct
unsupported-model message.
🤖 Prompt for all review comments with AI agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

Outside diff comments:
In `@src/main/eos_mesa_microphysics.f90`:
- Around line 559-569: The GR X-selection logic in mesa_eos_microphysics should
stop indexing mesa_eos_gr_data_exists with the exhausted file-loop variable and
an invalid zero column. Update the nx2 selection block to use the chosen Z index
(nz2) instead of i, and ensure the column lookup is 1-based when checking
adjacent entries in the mesa_eos_gr_data_exists array. Keep the fix localized to
the do j loop that computes nx2 so the boundary cases at j==2 and
j==mesa_eos_gr_nh remain valid.

In `@src/main/eos.f90`:
- Around line 1205-1213: Update the stale fatal message in get_p_from_rho_s so
it matches the actual ieos cases handled by the select case. The default branch
currently says only eos 2 and 12 are implemented, but the routine also supports
5 and 10; revise the fatal('eos', ...) text in that function to accurately list
all supported values so users get the correct unsupported-model message.

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro Plus

Run ID: ce80e499-8e30-4af3-a11b-5b07056dd231

📥 Commits

Reviewing files that changed from the base of the PR and between 4149bea and b50ffd7.

📒 Files selected for processing (3)
  • src/main/eos.f90
  • src/main/eos_mesa.f90
  • src/main/eos_mesa_microphysics.f90

@danieljprice

Copy link
Copy Markdown
Owner

@AliPourmand the test is failing because the binary test now takes more than 6 hours to run when it has to download the new MESA tables. Also the opacity table repeatedly downloads to the run directory because the mesa_opac directory does not exist. Need to put a README in an empty eos/mesa_opac directory to avoid this

Are the data files a lot larger than the previous tables? (and if so, do they need to be this large?)

@AliPourmand

AliPourmand commented Jul 17, 2026

Copy link
Copy Markdown
Contributor Author

@danieljprice The old link on zenodo for the tables are 29.4 MB and the new one is 33.2MB so the size shouldn't be the problem.
I added a README file, I hope that fixes it.

@AliPourmand

Copy link
Copy Markdown
Contributor Author

@danieljprice I tested running the binary workflow. When doing the final test giantbd.setup, the particle number of the companion is not specified in the .setup file which results in it assigning only 92 particles, which is why it gets stuck during relaxation:

WARNING: particle mass is already set, using np = 92 placed 92 particles in random-but-symmetric sphere

Unfortunately me manually trying to set np2=5000 didn't fix this. I have never seen this issue before in my previous pull requests. Has anything in setup_binary or this online file giantbd.setup changed recently?

@AliPourmand

Copy link
Copy Markdown
Contributor Author

@danieljprice I get a similar problem on my main branch too. I wasn't getting this error the last time I pulled, so I think something else changed in the past 1-2 weeks?

@danieljprice

Copy link
Copy Markdown
Owner

this test is not failing on main

@AliPourmand

Copy link
Copy Markdown
Contributor Author

@danieljprice I have pinpointed the problem, in the last test # now also run a test with a giant + brown dwarf binary (https://github.com/danieljprice/phantom/issues/777)
it works with the old tables and my code, but fails with my new tables.
The issue appears to be specific to the brown dwarf component — the giant star does not show the same problem. Looking at BDProfile.data, the electron degeneracy parameter (eta) becomes very large in the interior, indicating that the brown dwarf is entering a strongly degenerate regime.

image

Unfortunately, the specific internal energy is not included among the default columns in BDProfile.data, so I cannot directly compare the MESA profile internal energy against the EOSMESA tables. However, my suspicion is that the problematic particles are reaching the region above the red line in the plot, where the EOS tables transition towards the zero-temperature degeneracy limit. Above the zero-temperature branch was filled with an arbitrary maximum value because the tables were not intended to cover this region.

This is similar to the issue we encountered for white dwarfs, and we suggest using the Helmholtz EOS instead for any highly degenerate object such as Brown Dwarfs too.

The old tables did contain values in this region, but these were not physically accurate — they were simply continuous extrapolated values. This likely allowed the simulation to continue, although particles could still drift into this regime.

I suspect the new tables are failing because they correctly expose this unsupported region rather than providing continuous but unphysical values.

patchedzerotempmaps_EV_panel_1_z0 02_x0 80

@danieljprice

Copy link
Copy Markdown
Owner

This makes sense. Please discuss the correct approach with Mike Lau who contributed the brown dwarf example. It's not possible to use the zero T eos here as there is both a brown dwarf and a giant in the simulation. Ideally the tables should still work for this situation.

@AliPourmand

Copy link
Copy Markdown
Contributor Author

My concern is that this may not be possible with the current MESA EOS implementation. During Thomas Reichardt's time the MESA's EOS module may have behaved differently, but with the current version there are regions above the zero-temperature boundary where the MESA EOS solver simply does not return valid values:

eosfullmaps_EV_panel_1_z0 02_x0 80

Just to remind you, this is what the profile would look like for a white dwarf with a highly degenerate core. Even before Phantom starts evolving the system, some particles already lie slightly outside the tabulated EOS domain, probably due to effects such as interpolation between neighbouring metallicity tables:

twdgrid

@themikelau @danieljprice what do you think about this? I would strongly suggest we not "force" to make this work, because if it was an easy fix, we would have been using eos_mesa instead of eos_helmholtz for our white dwarf simulations.

@themikelau

Copy link
Copy Markdown
Collaborator

I did not actually contribute the brown dwarf example, and so I don't think you need my confirmation to determine the outcome of this discussion. My opinion is that we should continue the previous approach of throwing an extrapolation warning when falling outside the limits of the EoS table. But bothphantomsetup and phantom should still run even if they are producing unphysical results.

@AliPourmand

AliPourmand commented Jul 20, 2026

Copy link
Copy Markdown
Contributor Author

@themikelau I think there are two distinct issues here.

The extrapolation issue is the case where the EOS table have not been tabulated far enough into a physically valid region (e.g., below a certain internal energy), and in that case an extrapolating makes sense.

The issue I am referring to above the zero-temperature line is different. It corresponds to combinations of density and internal energy that are thermodynamically impossible for a given composition. In other words, there is no physical EOS solution for those states.

In that situation, I don't think Phantom should continue evolving with an extrapolated value, because the particle has entered a state that cannot exist according.

@danieljprice I would rather change this test to be an MS star instead of a brown dwarf.
If you think it is a good idea, I am happy to add a warning inside EOS-MESA that if someone uses a degenerate star, it stops the code and tells the user:
eos_mesa not valid for degenerate stars: use eos=15 (helmholtz), or 25 (zeroTemp) instead

There are two catches though: This would force the users to add the MESA internal energy as a column in their profiles too, which Phantom does not require right now.

Also, to do this I'd have to add the equation for zerotemperature in the eos-mesa subroutines too. (They exist in eos_zeroTemp but I'd rather not make eos_mesa and eos_zeroTemp dependent on each other)

The two catches are why I hadn't suggested it earlier.

@danieljprice

Copy link
Copy Markdown
Owner

we should discuss , we should 100% find a way make the test pass rather than changing the test

@AliPourmand

Copy link
Copy Markdown
Contributor Author

@danieljprice Sure, we can discuss.
Just to add some context, the figure below is the result of running the current main branch with the old tables, i.e. the version where this test passes. As you can see, essentially all of the brown dwarf particles lie above the zero-temperature line, meaning they occupy thermodynamically impossible states according to the EOS. So although the test passed, it was doing so with particles already outside the physically allowed region of the EOS.

P.S. This particular setup also represents the brown dwarf with only 92 SPH particles, compared to about 5000 particles for the giant. At that resolution, the brown dwarf is not really being modelled in a physically meaningful way anyway.

all_particles_phase_z0 04_x0 80

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 this pull request may close these issues.

3 participants