Block-structured AMR primer

Introduction

Understanding the ins and outs of adaptive mesh refinement is by no means trivial. In SOMAR, we’ve tried to automate things to the point that you can just change some input file parameters and watch the grids adapt as your simulation runs. But are you sure you are maximizing efficiency? Are you sure you are refining where you need to refine? Do you even need to use AMR or are you just making things complicated for yourself? Perhaps static mesh refinement is good enough. With these issues in mind, let’s look at the basics of SOMAR’s adaptive algorithm – when to use it, where to use it, and how the different levels of refinement work with each other.

What does “block-structured” mean?

“Block-structured” means that we cannot refine a single point. If amr.blockFactor = 16, then the fine level’s grids will be made of boxes that have at least 16 cells on a side. This leads to an inefficiency since more cells will be refined than were requested. This is not always a bad thing. First of all, underresolved cells are usually clustered, so refining an entire box that encloses the cluster can be quite practical. Second, iterating over a set of disjoint cells is rarely as efficient as iterating over a box of cells because rectangular boxes have simple dimensions that your processor’s prefetch cache can use to its advantage. Third, if you find a feature that needs refinement, it’s usually a good idea to refine additional surrounding cells to prevent the feature from advecting back into unrefined territory.

Block-structured AMR. — Fig. 1 2D block-structured AMR with two levels of refinement (three levels total).

Having said this, there are variants of AMR on structured grids that can refine individual cells, such as quadtrees in 2D or octrees in 3D. But the capability is often used to restrict the refined regions to a very tight region around complex obstacles, not to refine individual cells.

What gets refined?

Part of the difficulty with AMR is correctly identifying underresolved patches of the domain. There are several solutions that are commonly used. These include, but are not limited to, calculating

undivided differences of state variables,
vorticity, and
error via Richardson extrapolation.

Of these options, undivided differences are the simplest. Given any discretized state variable, \(q\), we compute the jumps between adjacent cells. If the magnitude of the jump is above some threshold, we tag both cells for refinement. A nice detail of this tagging method is that \(|\Delta q|\) will decrease with resolution. The downside of the method is that it is often difficult to identify proper thresholds for each state variable, as they depend on the physics of the flow.

As an example, consider the following input file snippet.

amr.velTagTol     = 1.0e-2
# amr.bTagTol       =
# amr.TTagTol       =
# amr.STagTol       =
# amr.bpertTagTol   =
# amr.TpertTagTol   =
# amr.SpertTagTol   =
# amr.scalarsTagTol =

This will refine any cell whose \(|\Delta u|\), \(|\Delta v|\), or \(|\Delta w|\) exceeds \(1.0 \times 10^{-2}\). If you prefer tagging on density differences, then uncomment and set amr.bTagTol. If you have a background stratification and would like to tag on temperature perturbations (pert = total - background), then use amr.TpertTagTol. Any number of these can be set at once, but don’t go crazy. Only set what makes sense for your simulation.

Todo

Tagging on vorticity

Tagging on the error computed via Richardson extrapolation is detailed in (CITE Dan Martin), but is not current implemented in SOMAR.

There are also methods that are specialized to the physics being studied. For example, turbulence in a stratified flow is expected when the gradient Richardson number, defined as the ratio of the local stratification to the horizontal shear squared, becomes less than 0.25 (Miles, 1961; Howard, 1961). In practice, we tag those cells that satisfy \(\text{Ri} \leq \mathcal{O}(1)\) in anticipation of turbulence. Unlike other tagging strategies, this strategy is robust in the sense that the tolerance is flow-independent. We can use the same threshold, \(\text{Ri} \leq \mathcal{O}(1)\), for all stratified flows and expect satisfactory results.

Some disturbances effect the entire water column. An internal solitary wave is a good example. In this case, an option is to project the buoyancy deviation, \(b'(x,y,z,t) = b(x,y,z,t) - \bar{b}(z)\), onto the vertical structure function

\[A(x,y,t) = \int b'(x,y,z,t) \, \varphi_0(z) \, dz,\]

where the integral is evaluated over the vertical domain extent. We can then increase the resolution over entire vertical sections where \(|A(x,y,t)|\) exceeds some threshold. This method was used successfully in (CITE upcoming JAMES paper).

Should I use AMR?

Whether or not you should use AMR in your project is, of course, dependent on what you are trying to simulate. We’ve tried to make SOMAR’s adaptive capabilities easy to turn on and tweak. This way, once you create a single-level simulation, you can easily experiment with the adaptive capabilities. So, my advice is to follow the guidelines below, and if you think you may benefit from AMR, start with a slightly unresolved simulation on a single level, then try adding adaptive grids to achieve your goals. If you end up deciding against adaptive grids in favor of global refinement, you can make it so right from your input file.

Do not use AMR if…

you need to resolve the entire bottom boundary layer. Try grid stretching instead.
the features that require high resolution are not localized. The overhead of AMR is only efficient when a small fraction of the domain is refined.

Use AMR if…

providing global resolution would require too much memory.
providing global resolution would make the simulation intractible.
the fine-scale features are highly localized in space and/or time.
you need to provide some extra resolution for a localized LES or non-linear interaction.
to increase resolution near small-scale features of an immersed boundary.
to increase resolution near a small portion of a boundary layer that is known to generate small-scale features, such as the generation site of an internal wave.

Glossary

level: A collection of grids and data at a specific resolution. Figure Fig. 1 shows three distinct levels. Note that only the level 0 spans the entire computational domain.
composite: The collection of valid grids or data over the entire AMR hierarchy.
refinement ratio: The amount of refinement between two levels. In figure Fig. 1, the refinement ratio between levels 0 and 1 is (2,1) and the refinement ratio between levels 1 and 2 is (1, 4).
block factor: The minimum size of a single grid. Small values allow you to keep refinement tight around an obstacle or feature. Large values allow the multigrid Poisson solver to achieve deeper V-cycles, potentially speeding up the simulation.
grid buffer: Block structured AMR does not allow us to jump multiple levels at a single coarse-fine interface. Notice in figure Fig. 1, levels 0 and 2 to not directly touch – level 1 provides a buffer region between them.
composite timestep: A complete timestep taken by a level and all of uts overlying finer levels.
subcycling: A level’s maximum allowable timestep is linked to its resolution – finer resolutions require smaller timesteps. When subcycling, each level takes the largest timestep it can, without being restricted by the timestep of its overlying finer levels. Subcycling allows an AMR simulation to be refned in time as well as space.
coarse-fine interface (CFI): In an N-dimensional simulation, the coarse-fine interfaces is the (N-1)-dimensional regions where two different resolutions meet.
ghost layer: A halo of cells around a grid that are used to help enforce boundary conditions or perform data exchanges among neighboring MPI ranks. Second-order derivative and interpolation stencils typically only require one ghost layer, but higher-order stencils can use more.
valid: The collection of a level’s cells/data that is not in a ghost layer and is not covered by finer cells/data.
invalid: The collection of cells that are not in valid regions. These cells are either covered by finer grids or exist in ghost layers.
layout: An arrangement of boxes that composes a single level.
patch/box/grid: These terms are interchangeable. They represent a single, rectangular collection of cells. In Fig. 1, levels 0 and 1 may be made of one or more grids, but level 2 must me made of at least two grids.
effective resolution: If we construct our adaptive grids appropriately, applying adequate resolution to all key features of a simulation, then the results should be equivalent to a single-level simulation at the finest level’s resolution. This is the effective resolution. It is equivalent to the base resolution times all of its refinement ratios. The effective resolution of a simulation using the grids in figure Fig. 1 would be 32-by-64.
average down: Fine-level data is typically assumed more accurate than coarse-level data. We often need to correct the coarse-level data by interpolating the fine data down to the coarse-level’s grids wherever it is available. The simplest and most useful interpolation scheme is to replace each coarse-level cell with an average of the data in the overlying fine level cells. This is often needed at the end of a composite timestep when synchronizing levels, or during the coarsening phase of a V-cycle.
interpolate up: Given a coarse-level of data, we can use a constant, linear, or quadratic interpolation scheme to fill the overlying fine level where it exists. This is often needed when a new fine level grid appears and needs to be initialized with data, or during the prolongation phase of a V-cycle.
refluxing: The AMR timestepping scheme can be viewed as a predictor-corrector method. First, the coarse level is timestepped to produce a prediction of the final state. Then, the fine levels are timestepped to produce more accurate corrections. These corrections usually push more or less fluid through the CFI than was predicted by the coarse level. Refluxing corrects these inconsistencies by computing the flux mismatches at the CFI and adding or removing heat/salt/mass/scalar concentrations from adjacent coarse cells.
tagging: The process of identifying which of a level’s cells require refinement. If no cells are tagged, then the finer level is either destroyed or not created. If cells are tagged, then SOMAR will ensure that the resulting fine level’s layout will be created or changed to include every tagged cell.
initialization: Once a fine level is created or changed, its new regions will require accurate data – the level needs to be initialized. This is a three-step process. First, we interpolate up from the coarse level to the new, empty regions. Second, we project the new velocity field to enforce the incompressibility constraint. Third, we run a very small timestep on the new level to generate a smooth pressure estimate. This last step can be made fast by increasing the tolerance of the pressure Poisson solver and using a forward Euler scheme.
synchronization: Once two or more levels are timestepped and reach the same time, they need to be synchronized. First, we reflux all cell-centered scalar fields to ensure global conservation. Second, we project the velocity on all synchronizing levels to enforce incopressibility at the CFI. Third, we correct all invalid data by averaging down from the finest to the coarsest synchronizing levels.