Skip to content

Fix MUSCL THINC right-state using already-overwritten left-state values#1181

Draft
sbryngelson wants to merge 8 commits intoMFlowCode:masterfrom
sbryngelson:fix/muscl-thinc-overwrite
Draft

Fix MUSCL THINC right-state using already-overwritten left-state values#1181
sbryngelson wants to merge 8 commits intoMFlowCode:masterfrom
sbryngelson:fix/muscl-thinc-overwrite

Conversation

@sbryngelson
Copy link
Member

@sbryngelson sbryngelson commented Feb 21, 2026

User description

Summary

Severity: HIGH — right-state reconstruction uses overwritten left-state values.

File: src/simulation/m_muscl.fpp, lines 285-301

In the MUSCL-THINC interface reconstruction, the left-state values (vL_rs_vf for contxb, contxe, and advxb) are overwritten during left reconstruction. The right-state reconstruction then divides by these already-overwritten values instead of the original cell-center values. While the math accidentally cancels (left THINC value / left THINC value = original ratio), this relies on exact cancellation and is fragile.

Before

! Left reconstruction overwrites vL_rs_vf
vL_rs_vf(j,k,l,contxb) = vL_rs_vf(j,k,l,contxb) / vL_rs_vf(j,k,l,advxb) * aTHINC
vL_rs_vf(j,k,l,advxb) = aTHINC
! Right reconstruction uses the overwritten values
vR_rs_vf(j,k,l,contxb) = vL_rs_vf(j,k,l,contxb) / vL_rs_vf(j,k,l,advxb) * aTHINC
!                         ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ overwritten!

After

! Save original density ratios before any overwrites
rho_b = vL_rs_vf(j,k,l,contxb) / vL_rs_vf(j,k,l,advxb)
rho_e = vL_rs_vf(j,k,l,contxe) / (1._wp - vL_rs_vf(j,k,l,advxb))
! Both reconstructions use the saved ratios
vL_rs_vf(j,k,l,contxb) = rho_b * aTHINC_L
vR_rs_vf(j,k,l,contxb) = rho_b * aTHINC_R

Why this went undetected

The algebraic cancellation makes this produce correct results in most cases. The refactor makes the intent explicit and eliminates dependence on accidental algebraic cancellation.

Test plan

  • Run MUSCL-THINC interface compression test
  • Verify left and right reconstructed states are physically consistent

🤖 Generated with Claude Code

Fixes #1200

Summary by CodeRabbit

  • Refactor
    • Revised interface-compression density update: precompute density-ratio factors and apply multiplicative updates consistently across left/right reconstruction branches in all three spatial directions. Enhances numerical consistency, stability, and parallel execution suitability.

CodeAnt-AI Description

Fix MUSCL-THINC using original density ratios for both left and right reconstructions

What Changed

  • Right-side THINC reconstruction now uses saved original density ratios instead of values overwritten by the left reconstruction, so left and right interface states are computed consistently.
  • Added small guards when computing density ratios to avoid divide-by-zero when the interface fraction approaches 0 or 1.
  • Marked per-thread THINC scalars as private in the GPU parallel loop to prevent data races on GPU backends.

Impact

✅ Fewer incorrect interface reconstructions
✅ Fewer simulation crashes from divide-by-zero at interfaces
✅ Fewer GPU data races during parallel reconstruction

💡 Usage Guide

Checking Your Pull Request

Every time you make a pull request, our system automatically looks through it. We check for security issues, mistakes in how you're setting up your infrastructure, and common code problems. We do this to make sure your changes are solid and won't cause any trouble later.

Talking to CodeAnt AI

Got a question or need a hand with something in your pull request? You can easily get in touch with CodeAnt AI right here. Just type the following in a comment on your pull request, and replace "Your question here" with whatever you want to ask:

@codeant-ai ask: Your question here

This lets you have a chat with CodeAnt AI about your pull request, making it easier to understand and improve your code.

Example

@codeant-ai ask: Can you suggest a safer alternative to storing this secret?

Preserve Org Learnings with CodeAnt

You can record team preferences so CodeAnt AI applies them in future reviews. Reply directly to the specific CodeAnt AI suggestion (in the same thread) and replace "Your feedback here" with your input:

@codeant-ai: Your feedback here

This helps CodeAnt AI learn and adapt to your team's coding style and standards.

Example

@codeant-ai: Do not flag unused imports.

Retrigger review

Ask CodeAnt AI to review the PR again, by typing:

@codeant-ai: review

Check Your Repository Health

To analyze the health of your code repository, visit our dashboard at https://app.codeant.ai. This tool helps you identify potential issues and areas for improvement in your codebase, ensuring your repository maintains high standards of code health.

Copilot AI review requested due to automatic review settings February 21, 2026 03:23
@codeant-ai codeant-ai bot added the size:S This PR changes 10-29 lines, ignoring generated files label Feb 21, 2026
@codeant-ai
Copy link
Contributor

codeant-ai bot commented Feb 21, 2026

Nitpicks 🔍

🔒 No security issues identified
⚡ Recommended areas for review

  • Division by zero
    The new code computes rho_b and rho_e by dividing by vL_rs_vf_${XYZ}$(...,advxb) and by (1 - vL_rs_vf_${XYZ}$(...,advxb)). If advxb equals 0 or 1 (or is extremely close), these divisions can produce infinities or very large values, causing NaNs/instabilities downstream. The denominators need a safe-guard (clamping or small epsilon) before division.

  • Unbounded ratios
    rho_b and rho_e are used directly to scale contxb/contxe after being computed. If the denominator is very small the computed ratios may be nonphysical (negative, >1, or huge), which will produce unphysical volume fractions when later multiplied by aTHINC. Consider clamping or validating rho_b/rho_e to a physically consistent range before applying.

  • GPU private list change
    The GPU parallel-loop private list was extended to include rho_b and rho_e. Verify these variables are fully local and initialized per iteration on the device and that their use introduces no unintended dependencies or extra register pressure. Also confirm the added temporaries do not exceed device resource limits for large collapse values.

Copy link
Contributor

Copilot AI left a comment

Choose a reason for hiding this comment

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

Pull request overview

Fixes MUSCL-THINC right-state reconstruction by preserving pre-overwrite density ratios from the left state and reusing them for both left and right reconstructions.

Changes:

  • Save pre-reconstruction density ratios (rho_b, rho_e) before THINC overwrites vL_rs_vf_* values.
  • Update left/right reconstructions to use the saved ratios instead of dividing by potentially overwritten left-state fields.

cubic-dev-ai[bot]

This comment was marked as off-topic.

@codecov
Copy link

codecov bot commented Feb 21, 2026

Codecov Report

❌ Patch coverage is 0% with 6 lines in your changes missing coverage. Please review.
✅ Project coverage is 44.06%. Comparing base (df28255) to head (fb5cb32).
⚠️ Report is 29 commits behind head on master.

Files with missing lines Patch % Lines
src/simulation/m_muscl.fpp 0.00% 6 Missing ⚠️
Additional details and impacted files
@@           Coverage Diff           @@
##           master    #1181   +/-   ##
=======================================
  Coverage   44.05%   44.06%           
=======================================
  Files          70       70           
  Lines       20496    20494    -2     
  Branches     1989     1989           
=======================================
  Hits         9030     9030           
+ Misses      10328    10326    -2     
  Partials     1138     1138           

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@sbryngelson sbryngelson requested review from anandrdbz and wilfonba and removed request for anandrdbz February 21, 2026 14:11
@codeant-ai codeant-ai bot added size:S This PR changes 10-29 lines, ignoring generated files and removed size:S This PR changes 10-29 lines, ignoring generated files labels Feb 22, 2026
Copy link
Contributor

@coderabbitai coderabbitai bot left a comment

Choose a reason for hiding this comment

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

Actionable comments posted: 2

🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.

Inline comments:
In `@src/simulation/m_muscl.fpp`:
- Around line 250-256: The GPU parallel loop's private clause for the MUSCL
reconstruction is missing scalars A, B, and C which leads to cross-thread data
races; update the private list in the GPU_PARALLEL_LOOP invocation (the clause
currently listing j,k,l,aCL,aC,aCR,aTHINC,moncon,sign,qmin,qmax,rho_b,rho_e) to
also include A, B, and C so that the loop-local temporaries A, B, C used in the
MUSCL/THINC computation are private per thread (look for the muscl_dir
conditional and the GPU_PARALLEL_LOOP macro surrounding the loop body where
A/B/C are assigned and used).
- Around line 282-284: Guard against near-zero denominators when computing rho_b
and rho_e by clamping or early-checking the MUSCL-reconstructed state
vL_rs_vf_${XYZ}$ (j,k,l,advxb) using the existing sgm_eps (or ic_eps) tolerance:
ensure advxb is not below sgm_eps before dividing for rho_b and ensure (1 -
advxb) is not below sgm_eps before dividing for rho_e; if either check fails,
use a safe fallback (e.g., set rho_b/rho_e to a clipped value or the
cell-average aC) so that rho_b, rho_e computations in this block cannot produce
NaN/Inf.

@codeant-ai codeant-ai bot added size:S This PR changes 10-29 lines, ignoring generated files and removed size:S This PR changes 10-29 lines, ignoring generated files labels Feb 23, 2026
@sbryngelson sbryngelson force-pushed the fix/muscl-thinc-overwrite branch from 9992699 to ac18720 Compare February 23, 2026 14:48
@codeant-ai codeant-ai bot added size:S This PR changes 10-29 lines, ignoring generated files and removed size:S This PR changes 10-29 lines, ignoring generated files labels Feb 23, 2026
@MFlowCode MFlowCode deleted a comment from github-actions bot Feb 23, 2026
@MFlowCode MFlowCode deleted a comment from github-actions bot Feb 23, 2026
@sbryngelson sbryngelson force-pushed the fix/muscl-thinc-overwrite branch from f6e0b83 to 34c2a78 Compare February 24, 2026 02:35
Copy link
Contributor

@coderabbitai coderabbitai bot left a comment

Choose a reason for hiding this comment

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

♻️ Duplicate comments (1)
src/simulation/m_muscl.fpp (1)

282-301: ⚠️ Potential issue | 🟡 Minor

Guard rho_b/rho_e denominators against near-zero reconstructed advxb.

Even with the aC gate, the reconstructed advxb can be pushed to 0 or 1, yielding divide-by-zero/Inf in rho_b/rho_e. Clamp with sgm_eps before division.

🛡️ Proposed fix
-                                rho_b = vL_rs_vf_${XYZ}$ (j, k, l, contxb)/vL_rs_vf_${XYZ}$ (j, k, l, advxb)
-                                rho_e = vL_rs_vf_${XYZ}$ (j, k, l, contxe)/(1._wp - vL_rs_vf_${XYZ}$ (j, k, l, advxb))
+                                rho_b = vL_rs_vf_${XYZ}$ (j, k, l, contxb)/ &
+                                        max(vL_rs_vf_${XYZ}$ (j, k, l, advxb), sgm_eps)
+                                rho_e = vL_rs_vf_${XYZ}$ (j, k, l, contxe)/ &
+                                        max(1._wp - vL_rs_vf_${XYZ}$ (j, k, l, advxb), sgm_eps)
🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@src/simulation/m_muscl.fpp` around lines 282 - 301, The division that
computes rho_b and rho_e can hit divide-by-zero when the reconstructed
advxb/advxe are exactly 0 or 1; before computing rho_b =
vL_rs_vf_${XYZ}(...,contxb)/vL_rs_vf_${XYZ}(...,advxb) and rho_e =
vL_rs_vf_${XYZ}(...,contxe)/(1._wp - vL_rs_vf_${XYZ}(...,advxb)), clamp the
denominators with sgm_eps (or a similar small threshold) to ensure they are at
least sgm_eps (e.g. denom_b = max(vL_rs_vf_${XYZ}(...,advxb), sgm_eps); denom_e
= max(1._wp - vL_rs_vf_${XYZ}(...,advxb), sgm_eps)) and then use denom_b/denom_e
in place of the raw denominator so subsequent assignments to vL_rs_vf_${XYZ} and
vR_rs_vf_${XYZ} (contxb, contxe, advxb, advxe) never produce Inf or NaN.
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.

Duplicate comments:
In `@src/simulation/m_muscl.fpp`:
- Around line 282-301: The division that computes rho_b and rho_e can hit
divide-by-zero when the reconstructed advxb/advxe are exactly 0 or 1; before
computing rho_b = vL_rs_vf_${XYZ}(...,contxb)/vL_rs_vf_${XYZ}(...,advxb) and
rho_e = vL_rs_vf_${XYZ}(...,contxe)/(1._wp - vL_rs_vf_${XYZ}(...,advxb)), clamp
the denominators with sgm_eps (or a similar small threshold) to ensure they are
at least sgm_eps (e.g. denom_b = max(vL_rs_vf_${XYZ}(...,advxb), sgm_eps);
denom_e = max(1._wp - vL_rs_vf_${XYZ}(...,advxb), sgm_eps)) and then use
denom_b/denom_e in place of the raw denominator so subsequent assignments to
vL_rs_vf_${XYZ} and vR_rs_vf_${XYZ} (contxb, contxe, advxb, advxe) never produce
Inf or NaN.

ℹ️ Review info

Configuration used: Path: .coderabbit.yaml

Review profile: CHILL

Plan: Pro

📥 Commits

Reviewing files that changed from the base of the PR and between 1d0265f and 34c2a78.

📒 Files selected for processing (1)
  • src/simulation/m_muscl.fpp

@sbryngelson sbryngelson force-pushed the fix/muscl-thinc-overwrite branch from 34c2a78 to 8cada11 Compare February 24, 2026 16:03
@codeant-ai codeant-ai bot added size:S This PR changes 10-29 lines, ignoring generated files and removed size:S This PR changes 10-29 lines, ignoring generated files labels Feb 24, 2026
Copy link
Contributor

@coderabbitai coderabbitai bot left a comment

Choose a reason for hiding this comment

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

♻️ Duplicate comments (1)
src/simulation/m_muscl.fpp (1)

282-284: ⚠️ Potential issue | 🟡 Minor

Near-zero denominator guard still missing on rho_b / rho_e

The outer guard checks the cell-average aC (from v_rs_ws_*_muscl), but the denominators on lines 283–284 use the MUSCL-reconstructed left state vL_rs_vf_${XYZ}$(j,k,l,advxb), which can legitimately be smaller than ic_eps (or 1 - advxb < ic_eps) when a neighbouring cell is nearly pure (TVD limiters allow the face value to reach the neighbour average). The existing sgm_eps in line 278 covers the same risk pattern; applying it here is consistent.

🛡️ Proposed fix
-                                rho_b = vL_rs_vf_${XYZ}$ (j, k, l, contxb)/vL_rs_vf_${XYZ}$ (j, k, l, advxb)
-                                rho_e = vL_rs_vf_${XYZ}$ (j, k, l, contxe)/(1._wp - vL_rs_vf_${XYZ}$ (j, k, l, advxb))
+                                rho_b = vL_rs_vf_${XYZ}$ (j, k, l, contxb)/ &
+                                        max(vL_rs_vf_${XYZ}$ (j, k, l, advxb), sgm_eps)
+                                rho_e = vL_rs_vf_${XYZ}$ (j, k, l, contxe)/ &
+                                        max(1._wp - vL_rs_vf_${XYZ}$ (j, k, l, advxb), sgm_eps)
🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@src/simulation/m_muscl.fpp` around lines 282 - 284, The division for rho_b
and rho_e uses the MUSCL-reconstructed face value vL_rs_vf_${XYZ}(j,k,l,advxb)
which can be near zero or near one; add the same sgm_eps safeguard used earlier
(e.g., max(vL_rs_vf_${XYZ}(...,advxb), sgm_eps) and max(1._wp -
vL_rs_vf_${XYZ}(...,advxb), sgm_eps)) so rho_b =
vL_rs_vf_${XYZ}(j,k,l,contxb)/max(vL_rs_vf_${XYZ}(j,k,l,advxb), sgm_eps) and
rho_e = vL_rs_vf_${XYZ}(j,k,l,contxe)/max(1._wp - vL_rs_vf_${XYZ}(j,k,l,advxb),
sgm_eps); keep other logic and naming unchanged (references: rho_b, rho_e,
sgm_eps, vL_rs_vf_${XYZ}).
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.

Duplicate comments:
In `@src/simulation/m_muscl.fpp`:
- Around line 282-284: The division for rho_b and rho_e uses the
MUSCL-reconstructed face value vL_rs_vf_${XYZ}(j,k,l,advxb) which can be near
zero or near one; add the same sgm_eps safeguard used earlier (e.g.,
max(vL_rs_vf_${XYZ}(...,advxb), sgm_eps) and max(1._wp -
vL_rs_vf_${XYZ}(...,advxb), sgm_eps)) so rho_b =
vL_rs_vf_${XYZ}(j,k,l,contxb)/max(vL_rs_vf_${XYZ}(j,k,l,advxb), sgm_eps) and
rho_e = vL_rs_vf_${XYZ}(j,k,l,contxe)/max(1._wp - vL_rs_vf_${XYZ}(j,k,l,advxb),
sgm_eps); keep other logic and naming unchanged (references: rho_b, rho_e,
sgm_eps, vL_rs_vf_${XYZ}).

ℹ️ Review info

Configuration used: Path: .coderabbit.yaml

Review profile: CHILL

Plan: Pro

📥 Commits

Reviewing files that changed from the base of the PR and between 34c2a78 and 8cada11.

📒 Files selected for processing (1)
  • src/simulation/m_muscl.fpp

The right reconstruction divides by left-state values that were already
overwritten during left reconstruction. Save original density ratios
(contxb/advxb and contxe/(1-advxb)) before left reconstruction and
reuse them for both left and right states.

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
@sbryngelson sbryngelson marked this pull request as draft February 25, 2026 03:26
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from coderabbitai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from Copilot AI Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from Copilot AI Feb 26, 2026
Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>
@github-actions
Copy link

Claude Code Review

Head SHA: fb5cb32850842ee09d8265cea514ab2db8b87d32
Files changed: 1 — src/simulation/m_muscl.fpp (+11 / -10)


Summary

  • Refactors MUSCL-THINC density reconstruction to pre-compute rho_b/rho_e ratios before any overwrites, eliminating dependence on the accidental cancellation pattern in the right-state reconstruction.
  • Adds A, B, C, rho_b, rho_e to the GPU parallel loop's private list — the omission of A, B, C was a pre-existing data race on GPU.
  • Fixes a comment typo (indicalindex).
  • The mathematics is unchanged in all non-degenerate cases; the code is cleaner and the GPU race is genuinely fixed.

Findings

1. A, B, C missing from original private list — pre-existing GPU data race (fixed here)
m_muscl.fpp (old line ~255, new line ~256)

A, B, and C are computed inside the GPU_PARALLEL_LOOP body but were absent from the private list in the original code. On GPU this causes a race condition — multiple threads write to the same stack-allocated scalars. The PR correctly adds all three. This is the most important correctness fix.


2. PR description claims divide-by-zero guards were added, but none appear in the diff
m_muscl.fpp, new lines ~283–284

The description states:

"Added small guards when computing density ratios to avoid divide-by-zero when the interface fraction approaches 0 or 1."

The actual diff has no such guards:

rho_b = vL_rs_vf_${XYZ}$ (j, k, l, contxb)/vL_rs_vf_${XYZ}$ (j, k, l, advxb)
rho_e = vL_rs_vf_${XYZ}$ (j, k, l, contxe)/(1._wp - vL_rs_vf_${XYZ}$ (j, k, l, advxb))

If the cell-center advxb is exactly 0 (or 1), rho_b (or rho_e) is 0/0. The pre-existing ic_eps clamps only protect aTHINC after reconstruction, not the input denominator here. The original code shared this vulnerability (same division was done inline), so this PR does not regress it, but the description is misleading. Consider either adding the advertised guards (e.g., max(advxb, ic_eps)) or removing the claim from the description.


3. Mathematical equivalence confirmation (informational)

The original right-state code accidentally cancelled correctly:

vR(contxb) = vL[overwritten](contxb) / vL[overwritten](advxb) * aTHINC_R
           = (rho_b * aTHINC_L) / aTHINC_L * aTHINC_R
           = rho_b * aTHINC_R   ✓

The cancellation holds for contxe too. The PR is therefore a correctness-neutral refactor for CPU paths (same results, clearer intent) and a genuine race-condition fix for GPU paths. The severity described ("HIGH — right-state uses overwritten values") is overstated for CPU; the GPU data race on A/B/C is the higher-severity issue.


Verdict

The core refactor is sound and improves readability and GPU safety. The GPU private-list fix for A, B, C is a real bug fix. The only concern is the misleading description about divide-by-zero guards — either add the guards or correct the description before merging.

@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from coderabbitai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from coderabbitai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from codeant-ai bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from github-actions bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from github-actions bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from github-actions bot Feb 26, 2026
@MFlowCode MFlowCode deleted a comment from github-actions bot Feb 26, 2026
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

size:S This PR changes 10-29 lines, ignoring generated files

Development

Successfully merging this pull request may close these issues.

MUSCL THINC right-state reconstruction uses already-overwritten left-state values

3 participants