Unexpected behavior when using multilevel MustRefineParticles

[tumlinson]

tl;dr - if you want accurate results with multilevel zooms using MustRefineParticles, you need to make a small code change to Grid_SetFlaggingField.C

Our group has encountered unexpected, and possibly unintended, behavior when Enzo is given multi-level zoom ICs using MustRefineParticles. We discovered this issue in the course of creating new ICs for low density environments. We began by running DM-only ICs, which behaved more or less as expected. When these same ICs were coverted to gas (using MUSIC and the enzo-mrp-music script to do all this), we found that enzo did not refine the same small DM halos as it had in the DM run. Further diagnosis led to the conclusion that cells/grids that "should" refine based on DM are not being refined unless they also meet the gas refinement criterion inside the MRP region, which for these small halos only happens late if at all. This basically makes it impossible to simulate small halos with the same setups and refinement logic we use in, e.g. the main FOGGIE runs.

Further investigation with the help of @bwoshea and @brittonsmith led to us to focus on Grid_SetFlaggingField, which contains the relevant logic (but not all of it, it's a bit convoluted). We eventually figured out that by flipping the sense of the "RestrictFlag" in the final call to FlagCellsToBeRefinedByMass in Grid_SetFlaggingField.C from TRUE to FALSE (line 324 of the main enzo-dev branch) removes the behavior. In effect, this is allowing cells/grids to refine on DM only, even if the gas refinement criterion is not met, leading to the physically expected behavior. Our hypothesis is that the default logic here is attempting to "unflag" cells/grids that are outside the region occupied by MRPs and that this logic goes awry somewhere. This is consistent with the observed pattern, especially the fact that the grids that get added back with the fix in place lie outside the MRP region.

We checked that this unexpected behavior was not unique to these ICs by restarting one of the main FOGGIE runs with the 'workaround' included, and found that indeed this fix brings out more refined structure in the outskirts of the Tempest halo, where the same circumstances hold (second figure).

I am filing this as an issue rather than a PR with the workaround fix because (1) I don't totally understand the behavior well enough to guarantee we found the root cause and (2) the workaround gives us results that we like, but at a cost of many more high resolution grids that others might not want.

I've attached plots showing the key behavior with and without the fix and can supply more information to anyone interested in either applying this fix or understanding it.

enzo issue MRPFix Workaround.pdf

[jwise77]

I'll look into this problem early next week. I've run several multi-level zooms with MRPs and haven't seen this behavior, but I'll see if I can reproduce this issue with one of my test zooms. You're correct in that the routine should unflag the cells that don't overlap with the MRPs. The workaround fix isn't ideal since the point of these particles was to restrict the AMR grids to only the region with the MRPs. I'll see what I can do since I helped with the original implementation.

[jwise77]

Hi @tumlinson, I can't reproduce the behavior you're seeing. Can you post the parameter file? Originally I thought the problem was that you were only using the DM mass and MRPs for refinement, whereas I usually include the gas mass for refinement, too. But I ran my test simulation (see below) with refinement methods = "4 8", and I still get the expected behavior. I'm also including my parameter file, so you can inspect it, too. I will continue to look into this issue.

[jwise77]

Upon closer inspection at later times, both the RefineMethod = "4 8" and "2 4 6 8" are producing grids where the MRPs don't exist, where the DM+MRP refinement produces many more grids for some reason. I'm using the latest main branch for these tests, btw.

Frames in the movie: MRPs, "4 8" DM+MRP refinement, "2 4 6 8" gas+DM+MRP refinement

cycle.mp4