AMR refinement criteria¶
Several parameters controlling refinement can be set in the setrun function. See Specifying AMRClaw run-time parameters in setrun.py for further description of these. Many of the parameters discussed below are attributes of rundata.amrdata in setrun.py.
Every regrid_interval time steps on each level, the error is estimated in all cells on grids at this level. Cells where some refinement criteria are satisfied are flagged for refinement. Default options for flagging are described below. Additional cells surrounding the flagged cells are also flagged to insure that moving features of the solution (e.g. shock waves) do not escape from the region of refinement before the next regridding time. The number of buffer cells flagged is specified by regrid_buffer_width and the number of steps between regridding on each level is specified by regrid_interval. Typically these are equal (assuming the Courant number is close to 1) and taken to be some small integer such as 2 or 3.
In addition to flagging individual cells based on the behavior of the solution, it is also possible to specify that certain regions of the domain should always be refined to a certain level (and/or never refined above some level). This is described further in Specifying AMR regions. These regions are used in conjunction with the methods described below to determine whether or not a given cell should be flagged for refinement.
The cells that have been flagged are then clustered into rectangular regions to form grids at the next finer level. The clustering is done in light of the tradeoffs between a few large grids (which usually means refinement of many additional cells that were not flagged) or many small grids (which typically results in fewer fine grid cells but more grids and hence more overhead and less efficient looping over shorter rows of cells). The parameter clustering_cutoff in amrNez.data is used to control this tradeoff. At least this fraction of the fine grid cells should result from coarse cells that were flagged as needing refinement. The value clustering_cutoff = 0.7 is usually reasonable.
Flagging criteria¶
Two possible approaches to flagging individual cells for refinement (based on the behavior of the solution) are built into AMRClaw. (A different default approach is used in GeoClaw, see Flagging criteria in GeoClaw).
Note: Starting in v5.6.0, a new approach is also available, see Guiding AMR with adjoint flagging.
flag2refine¶
One approach to flagging cells for refinement (the default used in most examples) is to set flag2refine == True and specify a tolerance flag2refine_tol. This indicates that the library subroutine $CLAW/amrclaw/src/Nd/flag2refine.f90 should be used to flag cells for refinement. This routine computes the maximum max-norm of the undivided difference of \(q_{i,j}\) based its four neighbors in two space dimensions (or 6 neighbors in 3d). If this is greater than the specified tolerance, then the cell is flagged for refinement (subject to limitations imposed by “regions”). The undivided difference (not divided by the mesh width) is used, e.g. \(|q(m,i+1,j) - q(m,i-1,j)|\) for each component \(m\).
Note that the user can change the criterion used for flagging cells by modifying this routine – best done by copying the library routine to your application directory and modifying the Makefile to point to the modified version.
Richardson extrapolation¶
The second approach to flagging individual cells is based on using Richardson extrapolation to estimate the error in each cell. This is used if flag_richardson == True. In this case a cell is flagged if the error estimate exceeds the value flag_richardson_tol. Richardson estimation requires taking two time steps on the current grid and comparing the result with what’s obtained by taking one step on a coarsened grid. One time step on the fine grid is re-used, so only one additional time step on the fine grid and one on a coarsened grid are required. It is somewhat more expensive than the flag2refine approach, but may be more useful for cases where the solution is smooth and undivided differences do not identify the regions of greatest error.
Note: Both approaches can be used together: if flag2refine == True and flag_richardson == True then a cell will be flagged if either of the corresponding specified tolerances is exceeded.
Specifying AMR regions¶
In addition to specifying a tolerance or other criteria for flagging individual cells as described above, it is possible to specify regions of the domain so that all points in the region, over some time interval also specified, will be refined to at least some level minlevel and at most some level maxlevel. These are specified through the parameter rundata.regiondata.regions in setrun.py. This is a list of lists, each of which specifies a region. A new region can be added via:
rundata.regiondata.regions.append([minlevel,maxlevel,t1,t2,x1,x2,y1,y2])
This indicates that over the time period from t1 to t2, cells in the rectangle x1 <= x <= x2 and y1 <= y <= y2 should be refined to at least minlevel and at most maxlevel.
To determine whether a grid cell lies in one of the regions specified, the center of the grid cell is used. If a mapped grid is being used, the limits for the regions should be in terms of the computational grid coordinates, not the physical coordinates.
If a cell center lies in more than one specified region, then the cell will definitely be flagged for refinement at level L (meaning it should be covered by a Level L+1 grid) if L+1 <= minlevel for any of the regions, regardless of whether the general flagging criteria hold or not. This means the smallest of the various minlevel parameters for any region covering this point will take effect. Conversely it will not be flagged for refinement if L+1 > maxlevel for all regions that cover this point. This means the largest of the various maxlevel parameters for any region covering this point will take effect. (However, note that since flagged cells are buffered as described above by flagging some adjacent cells, a cell may still end up flagged for refinement even if the above tests say it should not be.)
For example, suppose that amr_levels_max = 6 has been specified along with these two regions:
rundata.regiondata.regions.append([2, 5, 10.0, 30.0, 0.0, 0.5, 0.0, 0.5])
rundata.regiondata.regions.append([3, 4, 20.0, 40.0, 0.2, 1.0, 0.2, 1.0])
The first region specifies that from time 10 to 30 there should be at least 2 levels and at most 5 levels of refinement for points in the spatial domain 0 < x < 0.5 and 0 < y < 0.5.
The second region specifies that from time 20 to 40 there should be at least 3 level and at most 4 levels of refinement for points in the spatial domain 0.2 < x < 1.0 and 0.2 < y < 1.0.
Note that these regions overlap in both space and time, and in regions of overlap the maximum of the minlevel and also the maximum of the maxlevel parameters applies. So in the above example, from time 20 to 30 there will be at least 3 levels and at most 5 levels in the region of overlap, 0.2 < x < 0.5 and 0.2 < y < 0.5.
Within these regions, how many levels are chosen at each point will be determined by the error flagging criteria, i.e. by the default or user-supplied routine flag2refine, or as determined by Richardson extrapolation, as described above.
Points that are not covered by either region are not constrained by the regions at all. With amr_levels_max = 6 then they might be refined to any level from 1 to 6 depending on the error flagging criteria.
Implementation¶
It is perhaps easiest to understand how this works by summarizing the implementation. Note the order in which different flagging criteria are checked was modified in Version 5.5.0 in order to avoid the more expensive options for grid cells that are either forbidden from being refined or forced to be refined based on regions they lie in.
The regridding algorithm from level L to L+1 loops over all grid patches (in parallel when OpenMP is used with more than one thread). The cells on each patch are initially marked with amrflags(i,j) = UNSET, a special value (set in amr_module.f90).
In flagging based on regions:
If the current level is less than the maximum of all minlevel values for regions that contain the cell, then it is marked with amrflags(i,j) = DOFLAG.
If the current level is greater than or equal to the maximum of all maxlevel values for regions that contain the cell, then it is marked with amrflags(i,j) = DONTFLAG.
If there are any cells in the patch that are still marked as UNSET after checking all the regions, then if the setrun parameter flag2refine is True, then the routine flag2refine is called. The user may wish to create their own version of this routine. The default library version was modified with the addition of the adjoint option in Version 5.6.0 (see Guiding AMR with adjoint flagging), and does one of two things:
If adjoint_flagging then it checks the inner product of the forward solution with all adjoints over the specified time period and if the magnitude is greater than tolsp the cell is marked DOFLAG.
Otherwise, the undivided difference of all components of q in each coordinate direction is computed, e.g. abs(q(:,i+1,j) - q(:,i-1,j)) and abs(q(:,i,j+1) - q(:,i,j-1)) in 2d, and if the maximum of these is greater than tolsp the cell is marked DOFLAG.
If there are still any cells in the patch that are still marked as UNSET, then if the setrun parameter flag_richardson is True, then the routine errest is called. This does flagging based on estimates of the one-step error produced by Richardson extrapolation using the solution on the current grid and on a coarsened grid. If adjoint_flagging then these estimates are applied to the inner product of the error estimate with the adjoint solutions over the relevant time period. In either case, the setrun parameter flag_richardson_tol is used as the tolerance.
Flagging criteria in GeoClaw¶
In GeoClaw, a special flag2refine subroutine is defined. A standard approach for modeling tsunamis propagating across the ocean is to refine anywhere that the surface elevation of the ocean \(\eta = h+B\) is greater in absolute value than some specified wave_tolerance, e.g. 0.1 meter as set, for example, in the setrun.py file of $CLAW/geoclaw/examples/tsunami/chile2010. This wave_tolerance parameter can be set for any GeoClaw application.
Often this ends up refining the entire ocean when in fact only some waves are of interest. In this case one can use regions as described in Specifying AMR regions to limit refinement to certain space-time regions.
Alternatively, starting in Version 5.6.0 one can use adjoint flagging (see Guiding AMR with adjoint flagging) to better select the waves that will reach a particular location over a specified time range, including those that might reflect off distant shores.
Generally one must also use regions to allow (or force) much higher levels of refinement over small regions along the coast if one is doing detailed inundation modeling of a particular community.
