add ReadInterpolate thorn to Einstein Toolkit

Issue #2416 resolved
Roland Haas created an issue

ReadInterpolate is A FileReader like thorn that uses InterolateLocalUniform to interpolate the data read in onto the new grid.

itis available from and has been used in publications (

Comments (22)

  1. Yosef Zlochower

    The code seems to require that cctk_lsh[…] = cctk_ash[…] How restrictive is this in practice?

  2. Roland Haas reporter

    Not a whole lot of restriction given that other old thorns do the same. Given that this is a new thorn thorn, it should not do so. I will push a fix.

  3. Roland Haas reporter

    @Yosef Zlochower it would be good to add the thorn to the master thornlist soon so that it can mature a bit and be used / compiled on multiple clusters. This does not signify acceptance yet, just that is is harmless enough to be compiled by default.

    In your function as the reviewer, do you judge the code in the thorn to be well written enough to enable it in the master thorn list?

  4. Roland Haas reporter

    @Yosef Zlochower please add any comments about the thorn. Any comments are fine, from places where to code is unclear, to dead code, to missing / incorrect comments, to missing documentation or tests or demos.

  5. Yosef Zlochower

    Yes, please add it the the master thornlist. I think the only outstanding issues is that there needs to be documentation added and a testcase constructed that shows the output values are as expected. For the latter, I was thinking of proposing filling fields with low-order polynomials and then determining if the interpolated values match the polynomials to roundoff precision

  6. Roland Haas reporter

    I added documentation and a synthetic test case that shows that interpolation is indeed working correctly.

  7. Roland Haas reporter

    I removed interpolation in time, so the testsuite does not test it.

    The original code was (very) wrong (and cannot have ever worked).

    I coded up a version that wold stand a chance of working in a branch "rhaas/timerinterp" (just pushed to GitHub):

    but that one will not work either right now due to the CarpetIOHDF5 checkpoint files (the only one with multiple timelevels which would make interpolation in time possible) lacking information what time the actual time levels actually correspond to. There is not enough information present in the files to deduce what time it is based on the cctk_time attribute that is stored, the delta_time attribute in global attributes and the refinement level information, at least not without making assumptions about time refinement factors. So in the end I
    shelved time interpolation. If someone needs it and is willing to spend some time working with me testing it the branch can go in as a feature update sometime after the release once it has been actually used in production.

    The code can still read in multiple time levels but will only do interpolation in space, requiring that the time step is the same in the source and target simulation (which more or less limits its usefulness to going from a Cartesian to a Llama simulation or so).

  8. Roland Haas reporter

    The synthetic test data tests reading on multiple time levels but does not use time interpolation.

  9. Roland Haas reporter

    @Yosef Zlochower any verdict on the code yet? Once it is reviewed ok, it is put up for a final vote in the next ET call and the formally included in the ET.

  10. Yosef Zlochower

    My only concern was that the synthetic data function is an explicit function of time, but my understanding is that time interpolation doesn’t work.

  11. Roland Haas reporter

    Right, the synthetic data is a function of time so that when i read in three time levels (without interpolation in time) I can actually check that the code did not incorrectly read in data from the dataset for timelevel 0 into all three levels. Reading this now, I of course wonder if the option ReadInterpolate::read_only_timelevel_0 shouldn’t be removed since (I think) it is already included in the options Carpet::init_each_timelevel (equivalent to read_only_timelevel_0 being false) and Carpet::init_fill_timelevels (equivalent to read_only_timelevel_0 being true).

    I will put up ReadInterpolate for vote then. This week with an announcement that there will be a vote next week.

  12. Log in to comment