GRHydroToyGPU: Possible bug in flux_dens

[Federico G Lopez Armengol]

Possible bug in branch rhass/GRHydroToyGPU , file GRHydroToyGPU/src/fluxes.cxx , in lines 195 and 196:

Comparing the current expressions of flux_dens_m and flux_dens_p with equation (7.216) of Rezzolla & Zanotti (2013), I see a possible misplacement of alp_r[0-1] :

const CCTK_REAL flux_dens_m = dens_r[0] * (vel_r[0] - beta_r[0] / alp_r[0]);
const CCTK_REAL flux_dens_p = dens_r[1] * (vel_r[1] - beta_r[1] / alp_r[1]);

const CCTK_REAL flux_dens_m = dens_r[0] * (alp_r[0] * vel_r[0] - beta_r[0]);
const CCTK_REAL flux_dens_p = dens_r[1] * (alp_r[1] * vel_r[1] - beta_r[1]);

Actually, the current expression is valid only if sqrt_detg_r is the determinant of the four metric, but that would conflict with the definition of dens_r in line 191.

‌

[Roland Haas]

Is this still relevant? The equations

        //computing fluxes of conserved variabes
        const array<CCTK_REAL, 2> flux_dens = {
                dens_rc[0] * (vel_rc[0] - beta_avg / alp_avg),
                dens_rc[1] * (vel_rc[1] - beta_avg / alp_avg)
        };

in the code follow those in GRHydro code.

The current code is from c55ab39 "Few more modifications to fluxes.cxx" of cactusamrex

The GRHydro docs http://einsteintoolkit.org/thornguide/EinsteinEvolve/GRHydro/documentation.html#x1-80003 and (for the sources and actual ODE) eq. (16) and (17) of https://arxiv.org/abs/1304.5544 (note the “extra” lapse in (16) and (17)) have an extra lapse like the non-crossed out expression in the description. GRHydro’s code adds that one back when it compute the actual RHS:

            densrhs(i,j,k) = densrhs(i,j,k) + &
                (alp_l * densflux(i-xoffset,j-yoffset,k-zoffset) - &
                 alp_r * densflux(i,j,k)) * idx

‌

[Erik Schnetter]

All other things being equal, accessing the lapse in the RHS calculations increases memory bandwidth and reduces performance somewhat.

[Federico G Lopez Armengol]

If the lapse is added when computing the RHS, then I agree the code is correct. However, such implementation was confusing for me at least, besides non-optimal as Erik mentioned.

I see Jay is using the same expression as GRHydro in his version of fluxes.cxx, but I do not see any extra alpha in his rhs.cxx, so there could be a missing alpha in such version.

[Roland Haas]

It may be interesting to ask one of the original Whisky authors (which is not me) if they remember if there was a particular reason for what GRHydro (nee Whisky) does. These original authors may be Bruno (2nd gen admittedly), or maybe someone else at the AEI a that time (with the primary original author probably being Luca Baiotti, now in Japan).

[Erik Schnetter]

Or Ian Hawke, now in Southampton.

[Roland Haas]

Email thread on this is: http://lists.einsteintoolkit.org/pipermail/carpetx-developers/2021-November/000015.html

Basically either one is fine, though Ian says:

I would not do this now - keeping the fully densitized terms together is cleaner and a better way to extend to high order.

which I take to mean to not do what GRHydro does:

 densrhs(i,j,k) = densrhs(i,j,k) + &
                (alp_l * densflux(i-xoffset,j-yoffset,k-zoffset) - &
                 alp_r * densflux(i,j,k)) * idx

So we should check that we are consistent internally and if not make it consistent such that:

densrhs(i,j,k) = densrhs(i,j,k) + &
                (densflux(i-xoffset,j-yoffset,k-zoffset) - &
                 densflux(i,j,k)) * idx

‌

[Roland Haas]

Ping. HydroToy will likely go away or be archived since we now have “real” codes. However the issue, if it exists in other codes, may be interesting good to resolve.

[Jay Kalinani]

There are still a few bugs/design choices in GRHydroToyGPU that one should be careful of, since the code was still partially written before moving to AsterX. For instance, now in AsterX, we have a different definition of vtilde (vtilde_asterx = alp*vel - beta) which takes into account the multiplication of alp while computing densflux. While in Spritz and GRHydro, the multiplication by alp is done later explicitly as mentioned above, since their definition of vtilde is different (vtilde_spritz = vel - beta/alp ). I believe the bug in GRHydroToyGPU thorn regarding the multiplication of alp was not fixed after this thorn was abandoned.

‌