wiki:AstroBearAmr

Version 4 (modified by Brandon Shroyer, 14 years ago) ( diff )

—

  1. AMR Explained
    1. Basic Algorithm
  2. Book Keeping
  3. Dealing with MaxLevel
  4. Dealing with Lower Levels
    1. Levels 0 and below
    2. Levels -1 and below
    3. Level -2
  5. Communication
    1. Data
    2. Tree
  6. Threading

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 level n. Note that the tree structure and grid dimensions were created on the previous level n-1; InitInfos() just creates the data structures they reference.
  • ProlongateParentsData(n) — Populates the level n data structures with prolongated data from their parents on level n-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 level n with prolongated data from level n-1, AMR() copies over data from the previous generation of grids on level n. This higher-resolution data is preferable to data prolongated from the previous level, so AMR() uses it wherever it is available.
    • ApplyPhysicalBCs(n) — Apply physical boundary conditions to level n.
    • 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 level n+1 nodes, and the data referenced by the backed-up nodes will be used when ApplyOverlaps(n+1) is called (see above).
    • CreateChildrens(n)— This routine creates child nodes on level n+1 using the refinement flags set on level n. The use of "Childrens" here is not a typo; CreateChildrens() is so named because it applies the CreateChildren(Info) subroutine to each grid on level n.
    • InheritOldNodeOverlapsChildren(n) — Nested grids mean that spatial relationships (overlaps and neighbors) are inherited from parent grids. This routine passes information about the previous generation of n+1 grids to the new grids created by CreateChildrens(n). This routine is only executed on step 1 of the AMR() execution loop.
    • InheritNewNodeOverlapsChildren(n) — The children of previous level n grids will also need to send their data to the children of new level n grids.
    • InheritNeighborsChildren(n) — The children of neighboring grids on level n will likely be neighbors on level n+1. This routine passes neighbor information from level n to level n+1.
    • AdvanceGrids(n) — Performs the hyperbolic advance step on the grids of level n. This is where our numerical solvers come in.
    • AMR(n+1) — Launches AMR routine on the child level
    • ApplyChildrenData(n) — The inverse of ProlongateParentData(), this routine restricts data from the child grids onto their parent grids, providing a more accurate solution on the coarser level.
    • SyncFluxes(n) — To enforce mass conservation and the DivB constraint, the fluxes at grid boundaries need to be synchronized.
    • AccumulateFluxes(n) - Accumulates the fluxes used on level n to send back to the parent grids on level n-1.
  • 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.

Note: See TracWiki for help on using the wiki.