Version 2 (modified by 14 years ago) ( diff ) | ,
---|
AMR Explained
At first glance the AMR routine in amr_control.f90 is a little intimidating, but it can be understood by examining the most basic features first.
Please note that this is an explanation of AstroBEAR's AMR implementation, not adaptive mesh refinement itself. Readers who are unfamiliar with AMR might wish to look over some the AMR resources we have on the How AstroBEAR Works page.
Basic Algorithm
First we'll start with a simplified version of the AMR algorithm. This version focuses only on the most essential features:
RECURSIVE SUBROUTINE AMR(n) INTEGER :: n, nSteps, step nSteps = 2 CALL InitInfos(n) CALL ProlongateParentsData(n) DO step=1,nSteps levels(n)%step=step IF (step == 2) CALL UpdateOverlaps(n) CALL ApplyOverlaps(n,step) CALL ApplyPhysicalBCs(n) CALL SetErrFlags(n) IF (step == 2) CALL AgeNodesChildren(n) CALL BackupNodes(n+1) CALL CreateChildrens(n) IF (step == 1) THEN CALL InheritOldNodeOverlapsChildren(n) CALL InheritNewNodeOverlapsChildren(n) ELSE CALL InheritOverlapsOldChildren(n) CALL InheritOverlapsNewChildren(n) END IF CALL InheritNeighborsChildren(n) CALL AdvanceGrids(n) CALL AMR(n+1) CALL ApplyChildrenData(n) CALL SyncFluxes(n) CALL AccumulateFluxes(n) IF (step == 2) CALL NullifyNeighbors(n) END DO CALL CoarsenDataForParents(n) END SUBROUTINE AMR
In this example, the parameter n
represents the current level of the operation.
Remember as we step through the AMR()
subroutine that it is recursive. At each step on level n
, AMR(n)
within the DO nSteps
loop). calls itself on the next level up (see the calls to AMR(n+1)
during each step on level n
. This is because for each step on level n >= 0
, there are two steps on level n+1
, and the levels above the base level are regridded at each step.
IMPORTANT: Because AMR()
is recursive, each call to AMR()
at level n
assumes that certain steps were carried out by the call at level -1
. This can make the structure of AMR()
a little confusing, especially once parallelization is included.
AMR(n)
calls two subroutines before stepping through the simulation to finish initializing level n
's data:
InitInfos(n)
— Allocates grid data (InfoDef
) structures for the grids on leveln
. Note that the tree structure and grid dimensions were created on the previous leveln-1
;InitInfos()
just creates the data structures they reference.
ProlongateParentsData(n)
— Populates the leveln
data structures with prolongated data from their parents on leveln-1
.
Once the data has been constructed, we can begin the work of advancing the simulation. Each level n
will take two steps for every single step taken by n-1
, but the process is slightly different on each step for levels above the base level.
- Step 1
ApplyOverlaps(n)
— After initializing leveln
with prolongated data from leveln-1
,AMR()
copies over data from the previous generation of grids on leveln
. This higher-resolution data is preferable to data prolongated from the previous level, soAMR()
uses it wherever it is available.ApplyPhysicalBCs(n)
— Apply physical boundary conditions to leveln
.SetErrFlags(n)
— Determine which regions to refine. Refinement regions are determined by the physical processes involved, as well as specific conditions imposed by the problem modules.BackupNodes(n+1)
— Caches the nodes on the child level (n+1). We are about to create the new leveln+1
nodes, and the data referenced by the backed-up nodes will be used whenApplyOverlaps(n+1)
is called (see above).CreateChildrens(n)
— A favorite of hillbillies everywhere, this routine creates child nodes on leveln+1
using the refinement flags set on leveln
.InheritOldNodeOverlapsChildren(n)
— Nested grids mean that spatial relationships (overlaps and neighbors) are inherited from parent grids. This routine passes information about the previous generation ofn+1
grids to the new grids created byCreateChildrens(n)
. This routine is only executed on step 1 of theAMR()
execution loop.InheritNewNodeOverlapsChildren(n)
- This inherits the relationships going the other way. The children of previous level n grids will need to send their data to the children of new level n grids.- InheritNeighborsChildren - Neighbors of children will be children of neighbors
- AdvanceGrids - Performs hyperbolic advance of data structures
- AMR(n+1) - Launches AMR routine on the child level
- ApplyChildrenData - Inverse of ProlongateParentData - Applies restricted data from child grids to improve solution on coarser parent grids.
- SyncFluxes - To enforce mass conservation and DivB constraint - common fluxes used at grid boundaries need to be synchronized.
- AccumulateFluxes - We need to accumulate the used fluxes to send back to our parent grids.
- Step 2
- UpdateOverlaps - On the second step we don't need to receive overlap data from the previous generation of grids, however we do need to 'ghost' data with our current overlap grids or neighbors. So we think of our neighbors as our current overlaps. (Note we use the same nodelist but point to it twice with node%neighbors and node%overlaps - this is the reason we have to just NullifyNeighbors later)
- ApplyOverlaps - Need to ghost data from current overlaps (neighbor grids)
- ApplyPhysicalBCs - Apply physical boundary conditions.
- SetErrFlags - Determine which regions to refine
- AgeNodesChildren - Because of the nested grids giving us inherited relationships, we need to backup the relationships connecting us to the previous child grids on level n+1, as well as backing up the nodes themselves
- BackupNodes - Since we are about to create a new set of child nodes, we need to backup the previous nodes on the child level (n+1)
- CreateChildrens - Create child nodes on level n+1
- InheritOverlapsOldChildren - Nested grids mean that spatial relationships (overlaps/neighbors) are inherited from parent grids. On the second step the previous generation of level n+1 grids are the old children of the current generation of level n grids.
- InheritOverlapsNewChildren - This inherits the relationships going the other way. The old children of level n grids will need to send their data to the new children of level n grids.
- InheritNeighborsChildren - Neighbors of children will be children of neighbors
- AdvanceGrids - Performs hyperbolic advance of data structures
- AMR(n+1) - Launches AMR routine on the child level
- ApplyChildrenData - Inverse of ProlongateParentData - Applies restricted data from child grids to improve solution on coarser parent grids.
- SyncFluxes - To enforce mass conservation and DivB constraint - common fluxes used at grid boundaries need to be synchronized.
- AccumulateFluxes - We need to accumulate the used fluxes to send back to our parent grids.
- NullifyNeighbors - Since we called UpdateOverlaps we linked to the same nodelist with two sets of pointers (overlaps and neighbors) so we need to just nullify one of them - and destroy the other. (The destruction happens later in BackupNodes)
- CoarsenDataForParent - Finally we need to coarsen our data for parent grids.
Book Keeping
Now we will add in eight routines that add functionality and perform some basic book keeping
RECURSIVE SUBROUTINE AMR(n) INTEGER :: n, nSteps, step nSteps = 2 CALL InitInfos(n) CALL ProlongateParentsData(n) CALL ChildMaskOverlaps(n) DO step=1,nSteps levels(n)%step=step IF (step == 2) CALL UpdateOverlaps(n) CALL ApplyOverlaps(n,step) CALL ProlongationFixups(n) IF (lParticles) CALL ParticleUpdate(n) CALL ApplyPhysicalBCs(n) CALL SetErrFlags(n) IF (step == 2) CALL AgeNodesChildren(n) CALL BackupNodes(n+1) CALL CreateChildrens(n) IF (step == 1) THEN CALL InheritOldNodeOverlapsChildren(n) CALL InheritNewNodeOverlapsChildren(n) ELSE CALL InheritOverlapsOldChildren(n) CALL InheritOverlapsNewChildren(n) END IF CALL InheritNeighborsChildren(n) CALL AdvanceGrids(n) IF (lElliptic) CALL Elliptic(n) CALL PrintAdvance(n) CALL AMR(n+1) CALL ApplyChildrenData(n) CALL RestrictionFixups(n) CALL AfterFixups(n) CALL UpdateChildMasks(n) CALL SyncFluxes(n) CALL AccumulateFluxes(n) IF (step == 2) CALL NullifyNeighbors(n) END DO CALL CoarsenDataForParents(n) END SUBROUTINE AMR
- ParticleUpdate - If there are sink particles then we update the particles here
- Elliptic - If there are elliptic equations to solve we solve them here
- PrintAdvance - Just prints the 'Advancing level n …' line to standard out
- ProlongationFixups - It is better to complete the prolongation of the aux fields after receiving overlaps. This guarantees child grids have divergence free fields consistent with both their neighbors and their parents.
- ChildMaskOverlaps - This sets ChildMask to 0 for ghost cells that are refined by neighbors
- UpdateChildMask - This sets ChildMask to 1 for cells that are refined by children and sets ChildMask to NEIGHBOR_CHILD for cells that are refined by neighbor children
- RestrictionFixups - This updates cell centered representations of aux fields after receiving restricted data from children
- AfterFixups - This allows for user defined routines to be applied after a grid has been fully updated.
Dealing with MaxLevel
Up to this point we've assumed we are on an intermediate level of the AMR tree. What is different if we are on the highest level MaxLevel?
RECURSIVE SUBROUTINE AMR(n) INTEGER :: n, nSteps, step nSteps = 2 CALL InitInfos(n) CALL ProlongateParentsData(n) CALL ChildMaskOverlaps(n) DO step=1,nSteps levels(n)%step=step IF (step == 2) CALL UpdateOverlaps(n) CALL ApplyOverlaps(n,step) CALL ProlongationFixups(n) IF (lParticles) CALL ParticleUpdate(n) CALL ApplyPhysicalBCs(n) IF (n < MaxLevel) THEN CALL SetErrFlags(n) IF (step == 2 CALL AgeNodesChildren(n) CALL BackupNodes(n+1) CALL CreateChildrens(n) IF (step == 1) THEN CALL InheritOldNodeOverlapsChildren(n) CALL InheritNewNodeOverlapsChildren(n) ELSE CALL InheritOverlapsOldChildren(n) CALL InheritOverlapsNewChildren(n) END IF CALL InheritNeighborsChildren(n) END IF CALL AdvanceGrids(n) IF (lElliptic) CALL Elliptic(n) CALL PrintAdvance(n) IF (n < MaxLevel) CALL AMR(n+1) IF (n < MaxLevel) CALL ApplyChildrenData(n) CALL RestrictionFixups(n) CALL AfterFixups(n) IF (n < MaxLevel) CALL UpdateChildMasks(n) CALL SyncFluxes(n) CALL AccumulateFluxes(n) IF (step == 2) CALL NullifyNeighbors(n) END DO CALL CoarsenDataForParents(n) END SUBROUTINE AMR
We've basically wrapped the following routines that deal with child nodes inside of conditionals that prevent their execution on the MaxLevel
- SetErrFlags - We won't be refining any regions so no need to set error flags
- AgeNodesChildren - We have no children to age
- BackupNodes - There are no level MaxLevel+1 nodes to backup
- CreateChildrens - We don't create level MaxLevel+1 grids
- InheritOverlaps - Since level MaxLevel nodes have no children there is no inheritting that needs to be done
- AMR - We are on the MaxLevel
- ApplyChildrenData - No child data to apply
- UpdateChildMask - No children to modify the childmask
Dealing with Lower Levels
RECURSIVE SUBROUTINE AMR(n) INTEGER :: n, nSteps, step IF (n <= 0) nSteps=1 IF (n > 0) nSteps = 2 IF (n > -2) THEN CALL InitInfos(n) CALL ProlongateParentsData(n) IF (n > -1) CALL ChildMaskOverlaps(n) END IF DO step=1,nSteps levels(n)%step=step IF (step == 2) CALL UpdateOverlaps(n) IF (n > -2) CALL ApplyOverlaps(n,step) IF (n > 0) CALL ProlongationFixups(n) IF (n > -1 .AND. lParticles) CALL ParticleUpdate(n) IF (n > -1) CALL ApplyPhysicalBCs(n) END IF IF (n < MaxLevel) THEN IF (n > -1) THEN CALL SetErrFlags(n) END IF IF (step == 2 .OR. n == -2) THEN CALL AgeNodesChildren(n) END IF CALL BackupNodes(n+1) CALL CreateChildrens(n) IF (n == -2) THEN CALL InheritOverlapsOldChildren(n) CALL InheritNeighborsChildren(n) CALL InheritOverlapsNewChildren(n) ELSE IF (step == 1) THEN CALL InheritOldNodeOverlapsChildren(n) CALL InheritNewNodeOverlapsChildren(n) CALL InheritNeighborsChildren(n) ELSE CALL InheritOverlapsOldChildren(n) CALL InheritNeighborsChildren(n) CALL InheritOverlapsNewChildren(n) END IF END IF END IF IF (n > -1) THEN CALL AdvanceGrids(n) IF (lElliptic) CALL Elliptic(n) CALL PrintAdvance(n) END IF IF (n < MaxLevel) CALL AMR(n+1) IF (n < MaxLevel) CALL ApplyChildrenData(n) IF (n > -1) THEN CALL RestrictionFixups(n) CALL AfterFixups(n) END IF IF (n > -1) THEN IF (n < MaxLevel) CALL UpdateChildMasks(n) CALL SyncFluxes(n) END IF IF (n > 0) CALL AccumulateFluxes(n) IF (step == 2) CALL NullifyNeighbors(n) END DO IF (n > -2) CALL CoarsenDataForParents(n) END SUBROUTINE AMR
Levels 0 and below
The root Level (level 0) represents the lowest level of hydrodynamic data (although there are lower levels of nodes). As such all levels from the root level down do not need to call:
- ProlongationFixups
- AccumulateFluxes
Were in not for costmap data, these levels would neither need to call
- ProlongateParentsData
- CoarsenDataForParents
Levels -1 and below
Levels -1 and below do not need to call any routines related solely to hydrodynamic variables. This includes in addition to the routines above:
- ParticleUpdate
- ApplyPhysicalBCs
- SetErrFlags
- AdvanceGrids
- Elliptic
- PrintAdvance
- RestrictionFixups
- AfterFixups
- SyncFluxes
Additionally since the entire domain is refined at the root level, levels < 0 do not need to maintain the childmask array. So these levels do not need to call:
- ChildMaskOverlaps
- UpdateChildMasks
Level -2
The Level 2 grid is persistent so it does not need to be initialized or overlapped. So it does not need to call
- InitInfos
- ApplyOverlaps
Additionally the level 2 grid has no parent nodes so there is no need to call parent-related routines
- ProlongateParentsData
- CoarsenDataForParents
Finally since the level 2 is persistent, it behaves like a higher level grid in between steps so it always calls
- AgeNodesChildren
- InheritOverlapsOldChildren
- InheritOverlapsNewChildren
- InheritNeighborsNewChildren
Communication
Data
There are essentially four basic data routines that involve sharing of data between grids
- ProlongateParentsData - Parent to Child (Inter-Level)
- ApplyChildrenData - Child to Parent (Inter-Level)
- ApplyOverlaps - Old Grids to Current Grids (Intra-Level)
- SyncFluxes - Current Grids to Current Grids (Intra-Level)
For parallel applications this requires some degree of communication. In order to overlap computation with communication, it is good to post the sends as soon as the data is available - and to do as much computation as possible until having to wait for the receives to complete. When the sends are checked for completion and when the receives are first posted is somewhat arbitrary. It is reasonable to post the receives before you expect the sends to post and to complete the sends sometime after you expect the receives to have finished.
For each operation there is likely to be a degree of local sharing between grids. The basic approach therefore is to post the receives followed by the sends. Then perform the local sharing before waiting on the receives to complete, and then the sends. Sometimes the posting of the receives is shifted earlier, and the completion of the sends is put off until later. For example the parallel version of ApplyOverlaps is
CALL PostRecvOverlaps ... CALL PostSendOverlaps CALL ApplyOverlaps CALL CompRecvOverlaps ... CALL CompSendOverlaps
Tree
In a similar manner there are five tree operations that require some communication between nodes
- CreateChildren
- InheritNeighborsChildren
- InheritOldNodeOverlapsChildren
- InheritNewNodeOverlapsChildren
- InheritOverlapsOldChildren
- InheritOverlapsNewChildren
As in the case with data operations, each of these requires four communication calls in order to overlap the computation with communication. In all of these cases, it is node's children that are being communicated - since this is the only tree data that is created locally.
Threading
There are several threading options for parallelizing the hydro advance across levels. There are currently three basic approaches to address this
- Threading the Advances - The advancing of each level can be done independently although higher level threads should have higher priorities
- Threading the AMR levels - Each AMR level can also be thought of as an independent thread. Unfortunately this approach requires threads to communicate with other threads on different processors. This requires MPI to be extremely thread safe
- PseudoThreading - This is essentially careful scheduling of the advances to try and mimic the switching that would occur under a threaded implementation. This has the advance of not requiring any external libraries.
For more information on threading see the Scrambler Threading page.