Skip to content

RePrimAnd Con2Prim#126

Merged
MChabanov merged 13 commits intoEinsteinToolkit:developmentfrom
jaykalinani:rpa
Apr 7, 2026
Merged

RePrimAnd Con2Prim#126
MChabanov merged 13 commits intoEinsteinToolkit:developmentfrom
jaykalinani:rpa

Conversation

@jaykalinani
Copy link
Copy Markdown
Collaborator

@jaykalinani jaykalinani commented Jan 30, 2026

  • Add RPA support to Con2PrimFactory.

  • Code tested on Vista:

  1. Magnetic rotor test in AMR. Left RPA, right Noble/Palenzuela
Screenshot 2026-01-29 at 8 03 51 PM
  1. magTOV + Z4c + AMR
Noble_RPA
  1. Spherical shock test with Bz=1, vw_lim=10
RPA_Bz1_tol1e-8_iter100_vwlim10 Screenshot 2026-03-04 at 2 59 42 AM

@jaykalinani jaykalinani requested review from MChabanov and lwJi January 30, 2026 02:06
Copy link
Copy Markdown
Collaborator

@MChabanov MChabanov left a comment

Choose a reason for hiding this comment

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

Some comments from my side! Code looks very good!

Copy link
Copy Markdown
Collaborator

@MChabanov MChabanov left a comment

Choose a reason for hiding this comment

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

Great effort! I'm very impressed. Most comments are minor and are mainly related to consistency with the old code. Thanks!

atmo.set(pv, cv, glo);
return;
}

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

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

Note that this functionality already exists in line 226 of con2prim.cxx

if (cv.dens <= sqrt_detg * rho_atmo_cut) {

cv.DYe);
set_to_nan(pv, cv);
return;
}
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

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

isfinite(cv.dens) is checked below in line 116. cv.dens <= 0 should be captured before

if (cv.dens <= sqrt_detg * rho_atmo_cut) {

Hence, I'm not sure if we need this check here?

if (std::isfinite(a) && std::isfinite(b) && b >= a) {
const CCTK_REAL width = b - a;
const CCTK_REAL scale =
fmax(CCTK_REAL(1.0), fmax(std::abs(a), std::abs(b)));
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

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

Same comment as before for Palenzuela and Entropy:
To my understanding RPA solves for 1/(hW) which has the limits 0 < mu < 1/hmin ...
That means you will likely choose 1 for the scale ...
I would be good to double-check that, in case that's true the better lower limit would be again 0.0, so, replace line 181 with

fmax(CCTK_REAL(0.0), fmax(std::abs(a), std::abs(b)));

Copy link
Copy Markdown
Collaborator

@MChabanov MChabanov left a comment

Choose a reason for hiding this comment

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

All good!

@MChabanov MChabanov merged commit 452f1e9 into EinsteinToolkit:development Apr 7, 2026
3 checks passed
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.

2 participants