- removed comment
GRHydro atmosphere bug
Issue #889
closed
I have fixed a bug in GRHydro's atmosphere handling after reconstruction. After reconstruction, the reconstructed plus and minus face values are tested for whether they drop below atmosphere level and are then reset. In particular, this means that plus and minus face values for cell i don't generally coincide with the minus/plus face values of the corresponding cells i-1 and i+1. If the reconstruction previously found that the Riemann problem was trivial between, say rhominus(i+1) and rhoplus(i), when changing rhoplus(i), it is not anymore! However, the mask "trivial_rp" testing for a trivial Riemann problem is not changed! This leads to inconsistent behavior at the boundary of the atmosphere.
The symptom was a drift in the center of mass of a static and perfectly symmetric TOV star.
Keyword: GRHydro,
Keyword: atmosphere
Comments (9)
-
-
- removed comment
If I see this correctly the patch removes the loop vrom 1 to cctk_lsh and replaces it with a loop over only the interior - so the loop boundaries should be ok. However, this means that we should check whether a sync is necessary now. A possible way to circumvent that is to still loop over the entire grid, but add an explicit check for the lower bound when we would otherwise run outside of the GF. We don't have to reset something there anyway.
-
- removed comment
If it is necessary to compare values at different grid points (or different faces), then a sync is necessary. If you e.g. skip the test at the first point, then you are not testing whether this point is consistent with the respective neighbour points (that doesn't exist on this process). However, the neighbouring process will still perform this test.
-
- removed comment
GRHydro contains 'non-nice' stuff like functions that explicitly set 2..lsh-1 (in Fortran index notation) without sync because it is 'known' that the 'following' function only needs and reads that portion which makes a sync after the first function not necessary. That is why I said we have to carefully check.
-
- removed comment
Christian Reisswig and I discussed the patch and came up with a new version. This one extents the region of the grid it affects by one cell at the upper edge which means it matches the MON regions in GRHydro_PPM.F90 lines 770ff. It now also declares both the (i,j,k) and the (i-xoffset,j-yoffset,k-zoffset) Riemann Problems to be potentially non-trivial since rhominus affects the problem at (i-xoffset,j-yoffset,k-zoffset) and rhoplus affects the one at (i,j,k).
I am also noting that this loop does not take MHD into account and that another option would be to look at a single interface and its left/right values then reset only that Riemann Problem (to the average of the cell values left and right of the interface). That would be more of a new feature than an obvious bugfix though.
Ok to apply in this form?
-
- changed status to open
- removed comment
-
- changed status to open
- removed comment
-
- changed status to resolved
- removed comment
applied.
-
- edited description
- changed status to closed
- Log in to comment
- Assignee
- –
- Type
- Priority
- Status
- closed
- Component
- EinsteinToolkit thorn
- Milestone
- –
- Version
- –
- Votes
- 0
- Watchers
- 0
The reasoning given seems correct though I am not sure if the patch is correct. First you seem to loop over everything (from 1 to cctk_lsh) but then use the [xyz]offset to access elements other than i,j,k which seems as if you might access memory outside of the array. Second (and of this I am not sure since it depends on how the plus/minus arrays line up with the ordinary ones), I would expect that if a plus or minus value is reset that the Riemann problem affected differs since plus is the value of the field at the right side of a cell but minus is the value at the left side (I believe to remember). Otherwise: very good catch :-)
The diff issues of $890 also apply :-) this time to GRHydro_Reconstruct.F90.