The amr package uses 5 principal datastructures. These are declared within modules in the header files TREE.F, PHYSICALDATA.F, and WORKSPACE.F. The modules are included inside the various routines in the package through USE statements. Each processor has their own private version of these datastructures.
Definition of these datastructures is accomplished by user editing of the pre-processor include file PARAMESH_PREPROCESSOR.FH. Some pre-processor variables are then used to set fortran parameters in the module paramesh_dimensions.F, which in turn are used to dimension arrays within these datastructure.
The `tree datastructure' stores all the
information associated with the tree which manages the grid blocks and
the relationships between the grid blocks. It is defined in the module
TREE.F90
as
!----------------------------------------------------------------------The tree organizes a set of up to MAXBLOCKS sub-grid blocks on each processor. All the grid blocks are assumed to be logically cartesian with a uniform logical size. Each block has a level of refinement associated with it. The set of level 0 grid blocks cover the computational domain without overlap. Each grid can be the parent of 2**NDIM offspring grids which completely cover the sub-domain of their parents, where NDIM is the physical dimension of the simulation. The linear resolution varies by a factor of 2 between successive levels of refinement. At no point in the computational domain do we allow the resolution to jump by more than one level of refinement.
! PARAMESH - an adaptive mesh library.
! Copyright (C) 2003
!
! Use of the PARAMESH software is governed by the terms of the
! usage agreement which can be found in the file
! 'PARAMESH_USERS_AGREEMENT' in the main paramesh directory.
!----------------------------------------------------------------------
!
!!****h* headers/tree
!!
!! NAME
!!
!! physicaldata
!!
!! SYNOPSIS
!!
!! module tree
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!!
!! USES
!!
!! paramesh_dimensions
!!
!! DESCRIPTION
!!
!! Fortran 90 Module which 'holds' the tree data which PARAMESH uses
!! for a quad or oct-tree data structure, implemented on a parallel computer.
!! This data includes information such which blocks are neighbors of other
!! blocks, which block is another block's parent, and which blocks are child
!! blocks of a particular block.
!!
!! A convention is established for numbering the neighbors (or faces
!! of a block. The first neighbor is at lower x coordinate, the
!! second at higher x, the third at lower y, fourth at higher y, fifth
!! at lower z and the sixth at higher z.
!!
!! The convention by which the children of a block are numbered is the
!! same as the fortran array ordering, so that the first child is
!! at lower x, y and z coordinate, the second child is at upper x
!! but lower y and z, the third is at lower x, upper y and lower z,
!! and so on.
!!
!!
!! When a block has a refined neighbor we will need to know which children
!! of this neighbor are to provide guard cell information. The id's of the
!! correct children are stored in kchild using the conventions described
!! above. For example, if we are working on the 3rd neighbor of the
!! current block and it is at finer refinement level, then we must access
!! the children designated by kchild(:,3), in this case children 1, 2, 5
!! and 6.
!!
!! The tree organizes a set of up to maxblocks_tr grids on each processor.
!! All the grids are assumed to be cartesian with a uniform size. Each
!! grid has a level of refinement associated with it. The set of level 0
!! grids cover the computational domain without overlap. Each grid
!! can be the parent of 2**d offspring grids which completely cover
!! the sub-domain of their parents, where d is the physical dimension
!! of the simulation. The linear resolution varies by a factor of 2
!! between successive levels of refinement. At no point do we allow the
!! resolution to jump by more than one level of refinement.
!!
!!
!!
!! AUTHORS
!!
!! Peter MacNeice and Kevin Olson
!!
!!***
!!REORDER(5): unk, facevar[xyz], tfacevar[xyz]
!!REORDER(4): recvar[xyz]f
#include "paramesh_preprocessor.fh"
!-----------------------------------------------------------------
! tree module
module tree
use paramesh_dimensions
private
public :: maxblocks_tr
public :: nchild, nfaces, mchild, mfaces, mdim, mflags
public :: mboundaries,nboundaries
! Block limit to be used in manipulating tree when modifying the grid.
#ifndef LIBRARY
integer, parameter :: maxblocks_tr=4*maxblocks
#else
integer :: maxblocks_tr
#endif
! Number of children of a node
#ifndef LIBRARY
integer, parameter :: nchild=2**ndim
#else
integer :: nchild
#endif
! Number of faces on a grid block
#ifndef LIBRARY
integer, parameter :: nfaces=2*ndim
#else
integer :: nfaces
#endif
! Parameters used to define array sizes
integer, parameter :: mdim=3,mchild=2**mdim,mfaces=2*mdim
! Parameters used to declare the number of block marker flags needed
#ifndef LIBRARY
! --<< USER EDIT >>--
integer, parameter :: mflags=MFLAGS
#else
integer, save :: mflags
#endif
! Parameters used to declare the number of boundary regions where boundary
! conditions are to be applied. Typically: 2*ndim
#ifndef LIBRARY
! --<< USER EDIT >>--
integer, parameter :: nboundaries=NBOUNDARIES
#else
integer,save :: nboundaries
#endif
integer, parameter :: mboundaries=100
public :: neigh,child,which_child
public :: parent,lrefine,lnblocks,new_lnblocks
public :: nodetype,empty,bflags,newchild,derefine,refine
public :: stay,work_block,coord,bsize,bnd_box
public :: grid_xmin,grid_xmax,grid_ymin,grid_ymax
public :: grid_zmin,grid_zmax
public :: lrefine_max,lrefine_min
public :: level_cell_sizes
! Variables for storing tree datastructure
#ifndef LIBRARY
integer, save :: neigh(2,mfaces,maxblocks_tr)
integer, save :: child(2,mchild,maxblocks_tr)
integer, save :: which_child(maxblocks_tr)
integer, save :: parent(2,maxblocks_tr),lrefine(maxblocks_tr)
#else
integer, allocatable, save :: neigh(:,:,:)
integer, allocatable, save :: child(:,:,:)
integer, allocatable, save :: which_child(:)
integer, allocatable, save :: parent(:,:),lrefine(:)
#endif
integer, save :: lnblocks,new_lnblocks
#ifndef LIBRARY
integer, save :: nodetype(maxblocks_tr)
integer, save :: empty(maxblocks_tr),bflags(mflags,maxblocks_tr)
logical, save :: newchild(maxblocks_tr)
logical, save :: derefine(maxblocks_tr),refine(maxblocks_tr)
logical, save :: stay(maxblocks_tr)
#else
integer, allocatable, save :: nodetype(:)
integer, allocatable, save :: empty(:),bflags(:,:)
logical, allocatable, save :: newchild(:)
logical, allocatable, save :: derefine(:),refine(:)
logical, allocatable, save :: stay(:)
#endif
#ifndef LIBRARY
real, save :: work_block(maxblocks_tr)
real, save :: coord(mdim,maxblocks_tr)
real, save :: bsize(mdim,maxblocks_tr)
real, save :: bnd_box(2,mdim,maxblocks_tr)
#else
real, allocatable, save :: work_block(:)
real, allocatable, save :: coord(:,:)
real, allocatable, save :: bsize(:,:)
real, allocatable, save :: bnd_box(:,:,:)
#endif
real,save :: grid_xmin,grid_xmax
real,save :: grid_ymin,grid_ymax
real,save :: grid_zmin,grid_zmax
#ifndef LIBRARY
real,save :: level_cell_sizes(mdim,maxlevels)
#else
real, allocatable, save :: level_cell_sizes(:,:)
#endif
integer, save :: lrefine_max,lrefine_min
! flag to record grid change
public :: grid_changed, grid_analysed_mpi
integer, save :: grid_changed, grid_analysed_mpi
! added for surrblks calculation with mpi
public :: boundary_box,boundary_index
#ifndef LIBRARY
real, save :: boundary_box(2,mdim,mboundaries)
integer, save :: boundary_index(mboundaries)
#else
real, allocatable,save :: boundary_box(:,:,:)
integer, allocatable, save :: boundary_index(:)
#endif
! added for use with mpi block buffering
public :: strt_buffer,last_buffer
public :: strt_buffer_tree,last_buffer_tree
public :: laddress,surr_blks
#ifdef SAVE_MORTS
public :: surr_morts
#endif
integer, save :: strt_buffer,last_buffer
integer, save :: strt_buffer_tree,last_buffer_tree
#ifndef LIBRARY
integer, save :: surr_blks(3,3,1+2*k2d,1+2*k3d,maxblocks_alloc)
#ifdef SAVE_MORTS
integer, save :: surr_morts(6,3,1+2*k2d,1+2*k3d,maxblocks_alloc)
#endif
integer, save :: laddress(1:2,1:maxblocks_alloc)
#else
integer, allocatable, save :: surr_blks(:,:,:,:,:)
#ifdef SAVE_MORTS
integer, allocatable, save :: surr_morts(:,:,:,:,:)
#endif
integer, allocatable, save :: laddress(:,:)
#endif
! arrays to store info about block neighbors which are boundaries
public :: bc_block_neighs,bc_block_neighs_send
public :: bc_block_neighs_length
public :: bc_block_neighs_status
integer,save,allocatable :: bc_block_neighs(:,:)
integer,save,allocatable :: bc_block_neighs_send(:,:)
integer,save :: bc_block_neighs_length
integer,save :: bc_block_neighs_status
! DECLARE variables which are targets
target refine, derefine, newchild, empty
target lrefine, nodetype, work_block
target parent, coord, bsize, neigh
target child, bnd_box, stay
end module tree
A convention is established for numbering the neighbors (or faces of a block. The first neighbor is at lower x coordinate, the second at higher x, the third at lower y, fourth at higher y, fifth at lower z and the sixth at higher z.
The convention by which the children of a block are numbered is the same as the fortran array ordering, so that the first child is at lower x, y and z coordinate, the second child is at upper x but lower y and z, the third is at lower x, upper y and lower z, and so on.
The datastructures which store the model solution are declared in PHYSICALDATA.F90.
| Data Location | Array Names |
|---|---|
| Cell Center | UNK |
| Cell x Face Center | FACEVARX |
| Cell y Face Center | FACEVARY |
| Cell z Face Center | FACEVARZ |
| Cell x Edge Center | UNK_E_X |
| Cell y Edge Center | UNK_E_Y |
| Cell z Edge Center | UNK_E_Z |
| Cell Corner | UNK_N |
The location in a grid cell of each type of data in the solution
is shown in this diagram 
The declaration of UNK inside PHYSICALDATA.F90 is
! the solution for cell-centered quantities.
public :: unk
#ifndef LIBRARY
real,save :: unk(nvar, &
& il_bnd:iu_bnd,jl_bnd:ju_bnd,kl_bnd:ku_bnd, &
& maxblocks)
#else
real,allocatable,save :: unk(:,:,:,:,:)
#endif
target :: unk
The first dimension of UNK spans the data words required in the
solution
at each grid point. The second dimension is the x-coordinate index, the
third dimension is the y-coordinate index and the fourth dimension is
the z-coordinate index. The final dimension is the grid block number.
UNK can
be declared as a POINTER target so that a user can use pointers
to connect their variable names with UNK (to do this uncomment the
last line above and comment out the second last line).
NVAR is a user defined parameter. It should be set to the number of cell-centered data words per grid cell. It is set to the value defined for N_VAR in PARAMESH_PREPROCESSOR.FH or in the file amr_runtime_parameters if you are using PARAMESH in Library mode. If you are using a multi-step timestep integration algorithm (eg predictor-corrector) then the recommended value for NVAR is nvarp*(nphase + 1) where nvarp denotes the number of physical variables (ie. 1D hydro would have nvarp=3, for mass, momentum and energy), nphase is the number of stages in an integration timestep (ie predictor-corrector would have nphase=2). Similar considerations apply for nfacevar.
The array UNK may or may not have permanent memory allocated to store guardcell data. This choice is made by defining or not defining the preprocessor variable NO_PERMANENT_GUARDCELLS in PARAMESH_PREPROCESSOR.FH or in equivalent varirable in the file amr_runtime_parameters. If NO_PERMANENT_GUARDCELLS is defined then NPGS is set to 0, otherwise NPGS=1. The coordinate index ranges for UNK run from 1 to NXB+2*NGUARD*NPGS in the x-direction, 1 to NYB+2*NGUARD*NPGS in the y-direction and 1 to NZB+2*NGUARD*NPGS in the z-direction. NXB, NYB, NZB and NGUARD are parameters also defined by the user by setting the equivalent pre-processor variables in PARAMESH_PREPROCESSOR.FH.
The number of data words needed on a cell face is set by NFACEVAR. It is set to the value defined for N_FACEVAR in PARAMESH_PREPROCESSOR.FH or in amr_runtime_parameters. If there is no face centered data in your calculation set this to 0. The face-centered variables are also declared inside PHYSICALDATA.F90.
! the solution for cell-face-centered quantities.
public :: facevarx,facevary,facevarz
#ifndef LIBRARY
real,save :: facevarx(nbndvar, &
& il_bnd:iu_bnd+1,jl_bnd:ju_bnd, &
& kl_bnd:ku_bnd, &
& maxblocksf)
real,save :: facevary(nbndvar, &
& il_bnd:iu_bnd,jl_bnd:ju_bnd+k2d, &
& kl_bnd:ku_bnd, &
& maxblocksf)
real,save :: facevarz(nbndvar, &
& il_bnd:iu_bnd,jl_bnd:ju_bnd, &
& kl_bnd:ku_bnd+k3d, &
& maxblocksf)
#else
real,allocatable,save :: facevarx(:,:,:,:,:)
real,allocatable,save :: facevary(:,:,:,:,:)
real,allocatable,save :: facevarz(:,:,:,:,:)
#endif
target :: facevarx,facevary,facevarz
FACEVARX stores data located at the center of grid cell faces which are
perpendicular to the x direction, FACEVARY on faces perpendicular to
the y direction and FACEVARZ on faces perpendicular to
the z direction.
The parameter IFACE_OFF defines the indexing convention which relates the indeces associated with face defined variables to the indeces of the cells which they bound. For example, if IFACE_OFF=0 then the data in FACEVARX(:,i,:,:,:) is located on the faces of grid cells at x(i-1/2), and FACEVARX(:,i+1,:,:,:) on the faces of grid cells at x(i+1/2). On the other hand if IFACE_OFF=-1 then the data in FACEVARX(:,i-1,:,:,:) is located on the faces of grid cells at x(i-1/2), and FACEVARX(:,i,:,:,:) on the faces of grid cells at x(i+1/2).
The number of data words needed on a cell edge is set by NVAREDGE. It is set to the value defined for N_VAR_EDGE in PARAMESH_PREPROCESSOR.FH or in amr_runtime_parameters. If there is no edge centered data in your calculation set this to 0. The edge-centered variables are declared inside PHYSICALDATA.F90.
public :: unk_e_x,unk_e_y,unk_e_z
#ifndef LIBRARY
real,save :: unk_e_x(nbndvare, &
& il_bnd:iu_bnd,jl_bnd:ju_bnd+k2d, &
& kl_bnd:ku_bnd+k3d, &
& maxblocksue)
real,save :: unk_e_y(nbndvare, &
& il_bnd:iu_bnd+1,jl_bnd:ju_bnd, &
& kl_bnd:ku_bnd+k3d, &
& maxblocksue)
real,save :: unk_e_z(nbndvare, &
& il_bnd:iu_bnd+1,jl_bnd:ju_bnd+k2d, &
& kl_bnd:ku_bnd, &
& maxblocksue)
#else
real,allocatable,save :: unk_e_x(:,:,:,:,:)
real,allocatable,save :: unk_e_y(:,:,:,:,:)
real,allocatable,save :: unk_e_z(:,:,:,:,:)
#endif
target :: unk_e_x,unk_e_y,unk_e_z
The array UNK_E_X stores data at centers of cell edges aligned in the x coordinate direction, UNK_E_Y stores data at centers of cell edges aligned in the y coordinate direction and UNK_E_Z for z aligned edges.
The number of data words needed on a cell corner is set by NVARCORN. It is set to the value defined for N_VAR_CORN in PARAMESH_PREPROCESSOR.FH If there is no cell corner data in your calculation set this to 0. The cell corner variables are declared inside PHYSICALDATA.F.
! the solution for cell-corner based quantities.
public :: unk_n
#ifndef LIBRARY
real,save :: unk_n(nbndvarc, &
& il_bnd:iu_bnd+1,jl_bnd:ju_bnd+k2d, &
& kl_bnd:ku_bnd+k3d, &
& maxblocksn)
#else
real,allocatable,save :: unk_n(:,:,:,:,:)
#endif
target :: unk_n
There is a datastructure defined for use in enforcement of conservation laws and circulation integral constraints. This is used when conservation laws within the users model require consistent use of data at the common boundaries between grid blocks at different refinement levels.
The first set of arrays (FLUX_X, FLUX_Y, FLUX_Z) have the same
structure as the FACEVAR arrays, except that they only exist at block
boundaries.
Fluxes computed at block boundaries are stored by the user in these
arrays
at the time they are applied to update the solution.
These arrays are used in the AMR_FLUX_CONSERVE routines which guarantee
that
consistent fluxes are used at the common boundary between grid blocks
at
different refinement levels.
In addition we declare a second set of arrays (TFLUX_X, TFLUX_Y, TFLUX_Z) with the same structure as (FLUX_X, FLUX_Y, FLUX_Z). When the conservation enforcing routine amr_flux_conserve is called, the FLUX_ arrays are copied to their TFLUX_ equivalents, and then the FLUX_ arrays are modified at locations of refinement jumps to ensure conservation. The final step requires the user to modify their updated solution by the difference between the modified FLUX_ and the TFLUX_ arrays.
! storage used for fluxes at block boundaries. This is used when conservation ! constraints need to be imposed.The parameter NFLUXVAR defines the number of data words required on a typical block boundary grid cell face. If none are needed set this to 0. In the case of FLUX_X for example, the x index runs from 1 to 2, with 1 designating the left grid block boundary perpendicular to the x direction and 2 designating the right grid block boundary perpendicular to the x direction.
public :: nfluxvar,nfluxes,maxblocksfl
public :: flux_x,flux_y,flux_z
public :: tflux_x,tflux_y,tflux_z
#ifndef LIBRARY
integer, parameter :: nfluxvar = N_FLUX_VAR
integer, parameter :: nfluxes=max(1,nfluxvar)
integer, parameter :: maxblocksfl= &
& 1+(maxblocks-1)*min(1,nfluxvar)
real, save :: flux_x(nfluxes,1:2, &
& jl_bndi:ju_bndi,kl_bndi:ku_bndi,maxblocksfl)
real, save :: flux_y(nfluxes,il_bndi:iu_bndi, &
& 1:2,kl_bndi:ku_bndi,maxblocksfl)
real, save :: flux_z(nfluxes,il_bndi:iu_bndi, &
& jl_bndi:ju_bndi,1:2,maxblocksfl)
real, save :: tflux_x(nfluxes,1:2, &
& jl_bndi:ju_bndi, &
& kl_bndi:ku_bndi,maxblocksfl)
real, save :: tflux_y(nfluxes,il_bndi:iu_bndi, &
& 1:2,kl_bndi:ku_bndi,maxblocksfl)
real, save :: tflux_z(nfluxes,il_bndi:iu_bndi, &
& jl_bndi:ju_bndi,1:2,maxblocksfl)
#else
integer,save :: nfluxvar
integer,save :: nfluxes
integer,save :: maxblocksfl
real, allocatable, save :: flux_x(:,:,:,:,:)
real, allocatable, save :: flux_y(:,:,:,:,:)
real, allocatable, save :: flux_z(:,:,:,:,:)
real, allocatable, save :: tflux_x(:,:,:,:,:)
real, allocatable, save :: tflux_y(:,:,:,:,:)
real, allocatable, save :: tflux_z(:,:,:,:,:)
#endif
target :: flux_x, flux_y, flux_z
The second set of arrays provided for use in maintaining consistency across block boundaries at different refinement levels, are defined on the edges of grid cells at block boundaries.
These are typically used in MHD codes to store v x B type terms. The naming convention is that BEDGE_FACEX_Y stores data on a cell edge on the block boundary perpendicular to the x direction, with the edge pointing in the y direction. Again, if no edge values are needed, set N_EDGE_VAR to 0.
! storage used for cell edges at block boundaries.
! This is used when quantities located at cell edge centers need to
! be used consistently at the boundaries between blocks at different
! refinement levels.
public :: nedgevar1,nedgevar,nedges,maxblockse
public :: bedge_facex_y,bedge_facex_z,bedge_facey_x
public :: bedge_facey_z,bedge_facez_x,bedge_facez_y
public :: recvarx1e,recvary1e,recvarz1e
public :: recvarx2e,recvary2e,recvarz2e
#ifndef LIBRARY
integer, parameter :: nedgevar1=N_EDGE_VAR
integer, parameter :: nedgevar=max(nedgevar1,nvaredge)
integer, parameter :: nedges=max(1,nedgevar)
integer, parameter :: maxblockse= &
& 1+(maxblocks-1)*min(1,nedgevar)
real, save :: bedge_facex_y(nedges,1:2,jl_bnd:ju_bnd+1, &
& kl_bnd:ku_bnd+1,maxblockse)
real, save :: bedge_facex_z(nedges,1:2,jl_bnd:ju_bnd+1, &
& kl_bnd:ku_bnd+1,maxblockse)
real, save :: bedge_facey_x(nedges,il_bnd:iu_bnd+1,1:2, &
& kl_bnd:ku_bnd+1,maxblockse)
real, save :: bedge_facey_z(nedges,il_bnd:iu_bnd+1,1:2, &
& kl_bnd:ku_bnd+1,maxblockse)
real, save :: bedge_facez_x(nedges,il_bnd:iu_bnd+1, &
& jl_bnd:ju_bnd+1,1:2,maxblockse)
real, save :: bedge_facez_y(nedges,il_bnd:iu_bnd+1, &
& jl_bnd:ju_bnd+1, &
& 1:2,maxblockse)
#else
integer, save :: nedgevar1
integer, save :: nedgevar
integer, save :: nedges
integer, save :: maxblockse
real, allocatable, save :: bedge_facex_y(:,:,:,:,:)
real, allocatable, save :: bedge_facex_z(:,:,:,:,:)
real, allocatable, save :: bedge_facey_x(:,:,:,:,:)
real, allocatable, save :: bedge_facey_z(:,:,:,:,:)
real, allocatable, save :: bedge_facez_x(:,:,:,:,:)
real, allocatable, save :: bedge_facez_y(:,:,:,:,:)
#endif
public :: tbedge_facex_y,tbedge_facex_z,tbedge_facey_x
public :: tbedge_facey_z,tbedge_facez_x,tbedge_facez_y
#ifndef LIBRARY
real, save :: tbedge_facex_y(nedges,1:2,jl_bnd:ju_bnd+1, &
& kl_bnd:ku_bnd+1,maxblockse)
real, save :: tbedge_facex_z(nedges,1:2,jl_bnd:ju_bnd+1, &
& kl_bnd:ku_bnd+1,maxblockse)
real, save :: tbedge_facey_x(nedges,il_bnd:iu_bnd+1,1:2, &
& kl_bnd:ku_bnd+1,maxblockse)
real, save :: tbedge_facey_z(nedges,il_bnd:iu_bnd+1,1:2, &
& kl_bnd:ku_bnd+1,maxblockse)
real, save :: tbedge_facez_x(nedges,il_bnd:iu_bnd+1, &
& jl_bnd:ju_bnd+1, &
& 1:2,maxblockse)
real, save :: tbedge_facez_y(nedges,il_bnd:iu_bnd+1, &
& jl_bnd:ju_bnd+1, &
& 1:2,maxblockse)
#else
real, allocatable, save :: tbedge_facex_y(:,:,:,:,:)
real, allocatable, save :: tbedge_facex_z(:,:,:,:,:)
real, allocatable, save :: tbedge_facey_x(:,:,:,:,:)
real, allocatable, save :: tbedge_facey_z(:,:,:,:,:)
real, allocatable, save :: tbedge_facez_x(:,:,:,:,:)
real, allocatable, save :: tbedge_facez_y(:,:,:,:,:)
#endif
An `working block' was added in version 2.0. This is to enable PARAMESH to operate without allocating permanent memory to store guardcell information on every block.
The working block does have memory allocated for guardcell storage. The idea is that any modifications to the solution are made in the working block. To advance the solution on any block, block A say, we first copy block A into the working block on the local processor and fill the guardcells of the working block with data appropriate for A. We can then advance the solution in the working block, and when we are done, copy the solution back to the A's permanent storage, without its guardcells.
To support this each array in the solution datastructure has an equivalent array associated with the working block. These arrays have the same names but with a suffix 1, ie UNK1, FACEVARX1, FACEVARY1, FACEVARZ1, UNK_E_X1, UNK_E_Y1, UNK_E_Z1 and UNK_N1. They are declared in PHYSICALDATA.F. These arrays have 5 dimensions. The fifth dimension can take the value 1 or 2. When using these arrays it is expected that a user should only need to use then with the fifth index set to 1. (PARAMESH sometimes needs to store data for the working block's parent. This it by setting the fifth index to 2. It is not expected that a user should need to refer directly to the working block parent's data.)
! This file declares the storage space used to handle the `current
! working block' when the user decides not to reserve permanent
! storage space for guardcells for all blocks, but instead to
! fill guardcells as needed. This strategy requires 2 working blocks,
! one for the leaf node and one for its parent.
public :: unk1,facevarx1,facevary1,facevarz1
public :: unk_e_x1,unk_e_y1,unk_e_z1
public :: unk_n1
#ifndef LIBRARY
! the solution for cell-centered quantities.
real, save :: unk1(nvar,il_bnd1:iu_bnd1,jl_bnd1:ju_bnd1, &
& kl_bnd1:ku_bnd1, &
& npblks)
! the solution for cell-face-centered quantities.
real, save :: facevarx1(nbndvar,il_bnd1:iu_bnd1+1, &
& jl_bnd1:ju_bnd1, &
& kl_bnd1:ku_bnd1,npblks)
real, save :: facevary1(nbndvar,il_bnd1:iu_bnd1, &
& jl_bnd1:ju_bnd1+k2d, &
& kl_bnd1:ku_bnd1,npblks)
real, save :: facevarz1(nbndvar,il_bnd1:iu_bnd1, &
& jl_bnd1:ju_bnd1, &
& kl_bnd1:ku_bnd1+k3d,npblks)
! the solution for cell-edge-centered quantities.
real,save :: unk_e_x1(nbndvare, &
& il_bnd1:iu_bnd1,jl_bnd1:ju_bnd1+k2d, &
& kl_bnd1:ku_bnd1+k3d, &
& npblks)
real,save :: unk_e_y1(nbndvare, &
& il_bnd1:iu_bnd1+1,jl_bnd1:ju_bnd1, &
& kl_bnd1:ku_bnd1+k3d, &
& npblks)
real,save :: unk_e_z1(nbndvare, &
& il_bnd1:iu_bnd1+1,jl_bnd1:ju_bnd1+k2d, &
& kl_bnd1:ku_bnd1, &
& npblks)
! the solution for cell-corner based quantities.
real,save :: unk_n1(nbndvarc, &
& il_bnd1:iu_bnd1+1,jl_bnd1:ju_bnd1+k2d, &
& kl_bnd1:ku_bnd1+k3d, &
& npblks)
#else
! the solution for cell-centered quantities.
real, save, allocatable :: unk1(:,:,:,:,:)
! the solution for cell-face-centered quantities.
real, save, allocatable :: facevarx1(:,:,:,:,:)
real, save, allocatable :: facevary1(:,:,:,:,:)
real, save, allocatable :: facevarz1(:,:,:,:,:)
! the solution for cell-edge-centered quantities.
real, save, allocatable :: unk_e_x1(:,:,:,:,:)
real, save, allocatable :: unk_e_y1(:,:,:,:,:)
real, save, allocatable :: unk_e_z1(:,:,:,:,:)
! the solution for cell-corner based quantities.
real, save, allocatable :: unk_n1(:,:,:,:,:)
#endif
target :: unk1
target :: facevarx1, facevary1, facevarz1
target :: unk_e_x1, unk_e_y1, unk_e_z1
target :: unk_n1
Not to be confused with the working block!
The final datastructure is a workspace array called WORK, which is
declared inside WORKSPACE.F, and whose structure is
similar to that of UNK. There are 2 differences. As with UNK, you can
specify multiple cell-centered variables per grid cell.
However the index which discriminates between variables is the last
dimension (ie the fifth). The second difference is that the number of
guard cells associated with WORK is set by the parameter NGUARD_WORK,
which
can differ from NGUARD, the number of guard cells provided for UNK.
WORK is provided to enable the user to perform computations which need only one data word from each grid cell. For example, consider a routine which tests variation in the mass density between grid cells in order to decide if grid blocks need to be refined or derefined. Let us further suppose that the test wants to look 4 grid cells (say) in any direction to detect features in the solution which will arrive shortly and will require refinement. This means that 4 guard cells need to be filled with correct data from neighboring blocks. It is inefficient to use UNK to pass this information because we do not want to communicate all the other unneeded data words at each grid cell during the guard cell filling operation, and we do not want to set NGUARD to 4 just for this purpose, if it is not required to be this large in the rest of the code. Instead we simply copy the density from UNK into WORK(:,:,:,:,IOPTW), where IOPTW is greater than 1 and less than or equal to NVAR_WORK+1, and then instruct the guardcell filling routine (AMR_GUARDCELL or AMR_1BLK_GUARDCELL) to operate on WORK (IOPT=IOPTW) instead of UNK (IOPT =1).
! workspace arrays<>
public :: work
#ifndef LIBRARY
real, save :: work(ilw:iuw,jlw:juw,klw:kuw,maxblocksw, &
& nvar_work)
#else
real, allocatable, save :: work(:,:,:,:,:)
#endif
! common block storing the solution for cell-centered quantities.
public :: work1
#ifndef LIBRARY
real, save :: work1(ilw1:iuw1,jlw1:juw1,klw1:kuw1,npblks)
#else
real, allocatable, save :: work1(:,:,:,:)
#endif
The subroutines in the package form three groups. The first is a set of higher level routines which for most applications should be the only routines that an application developer needs to call. The second group of routines are provided by the user who can craft them from templates provided in the templates subdirectory. The remaining routines are lower level and are called by the higher level routines. Because we use modules to make the various datastructures available to routines, most routines have very simple argument lists.
This routine copies a global solution update from the time synchronized global solution arrays, into the arrays used during the solution update, as is required when using NO_PERMANENT_GUARDCELLS and the amr_1blk_guardcell routines.
subroutine amr_1blk_copy_soln(level)
Arguments:
level input integer If -1 then blocks at all refinement
levels are copied, otherwise only blocks
at level are copied.
Subroutines called by this routine:
amr_gsurrounding_blks
mpi_amr_global_domain_limits
mpi_amr_morton_limits
mpi_morton_bnd
mpi_amr_gsurr_blks
This routine manages the transfer of guard cell data for a specific single block.
!! NAME
!!
!! amr_1blk_guardcell
!!
!! SYNOPSIS
!!
!! call amr_1blk_guarcell(mype, iopt, nlayers, lb, pe,
!! lcc, lfc, lec, lnc, l_srl_only,
!! icoord, ldiag)
!! call amr_1blk_guarcell(mype, iopt, nlayers, lb, pe,
!! lcc, lfc, lec, lnc, l_srl_only,
!! icoord, ldiag,
!! nlayersx, nlayersy, nlayersz)
!!
!! call amr_1blk_guarcell(integer, integer, integer, integer, integer,
!! logical, logical, logical, logical, logical,
!! integer, logical,
!! optional integer, optional integer, optional integer)
!!
!! ARGUMENTS
!!
!! integer, intent(in) :: mype
!! The local processor number.
!!
!! integer, intent(in) :: iopt
!! A switch to control which data source is to be used:
!! iopt=1 will use 'unk', 'facevarx', 'facevary', 'facevarz',
!! 'unk_e_x', 'unk_e_y', 'unk_e_z', and 'unk_n'
!! iopt>=2 will use 'work'
!!
!! integer, intent(in) :: nlayers
!! The number of guard cell layers at each boundary
!! (Note: this argument is deprecated not longer functions !).
!!
!! integer, intent(in) :: lb
!! The block number on processor pe selected for guardcell filling.
!!
!! integer, intent(in) :: pe
!! The processor storing the block selected for guarcell filling.
!!
!! logical, intent(in) :: lcc
!! A logical switch controlling whether unk or work data is filled.
!!
!! logical, intent(in) :: lfc
!! A logical switch controlling whether facevarx(y)(z) data is filled.
!!
!! logical, intent(in) :: lec
!! A logical switch controlling whether unk_e_x(y)(z) data is filled.
!!
!! logical, intent(in) :: lnc
!! A logical switch controlling whether unk_n data is filled.
!!
!! logical, intent(in) :: l_srl_only
!! A logical switch which, if true, switches off the filling from
!! coarse neighbors. This is used during restriction when odd
!! block sizes are in use.
!!
!! integer, intent(in) :: icoord
!! An integer switch used to select which faces of the block are to be
!! considered. If icoord=0 all faces are considered. If icoord=1 only
!! faces perp. to the x-axis are considered, if icoord=2 only faces
!! perp. to the y-axis are considered, and if icoord=3 only faces perp.
!! to the z-axis are considered.
!!
!! logical, intent(in) :: ldiag
!! A logical switch which controls whether guardcells corresponding to
!! neighbor blocks diagonally opposite block edges and corners are filled.
!!
!! optional, integer, intent(in) :: nlayersx, nlayersy, nlayersz
!! Optional arguments to specify the number of guardcells to fill in each
!! coordninate direction. If not specified, then the default is to fill all
!! available guardcells.
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!! mpif.h
!!
!! USES
!!
!! paramesh_dimensions
!! physicaldata
!! tree
!! timings
!! workspace
!! mpi_morton
!! paramesh_mpi_interfaces
!! paramesh_interfaces
!!
!! CALLS
!!
!! mpi_amr_local_surr_blks_lkup
!! mpi_amr_1blk_guardcell_c_to_f
!! mpi_amr_get_remote_block
!! amr_perm_to_1blk
!! amr_1blk_guardcell_srl
!! amr_1blk_guardcell_f_to_c
!!
!! RETURNS
!!
!! Does not return anything. Upon exit, guardcells are filled with data for
!! the block secified in the call to this routine.
!!
!! DESCRIPTION
!!
!! This routine manages the transfer of guard cell data for a
!! specific single block. It uses the morton numbering scheme for the
!! determination of the neighbor relations and differs in that from
!! the older routine amr_1blk_guardcell. The user can choose how many layers
!! of guardcells are filled on each face of a block by specifying the optional
!! arguments 'nlayersx', 'nlayersy', and 'nlayersz'.
!!
!! NOTES
!!
!! This routine must be used with care !
!! This routine was written to be used in a code as illustrated
!! in the following snippet of pseudo-code
!!
!! .
!! .
!! .
!! synchronization pt
!!
!! loop over grid blocks
!! (copy current block into working block and fill its guardcells)
!! (perform some set of operations on block)
!! (store result from working block)
!! end loop
!!
!! synchronization pt
!! .
!! .
!! .
!!
!! Caveat 1:
!! If you are using this routine, you must remember to save the solution
!! at the first synchronization point (ie call amr_1blk_copy_soln), so
!! that each block uses the same time synchronized solution during its
!! guardcell filling.
!!
!! Caveat 2:
!! It is implicitly assumed that the parent blocks on all leaf nodes
!! have valid data. (This is necessary to ensure that a general restriction
!! operator can be supported in the neighborhood of a jump in refinement.)
!! If ADVANCE_ALL_LEVELS is defined, then this will generally be true.
!! However, if the solution is being time-advanced on leaf blocks only
!! this may not be true. In this case you should call amr_restrict
!! before the first synchronization pt in the example above.
!! If you are using blocks with an even number of grid points and the
!! default restriction operators this is not necessary.
!!
!! AUTHORS
!!
!! Michael Gehmeyr & Peter MacNeice (November 1999) with modifications by
!! Kevin Olson for layered guardcell filling.
!!
!!***
This routine copies data from the 1-block working arrays with guardcells to the permanent data arrays, which may or may not have permanent guardcells, depending on whether NO_PERMANENT_GUARDCELLS is defined in paramesh_preprocessor.fh.
subroutine amr_1blk_to_perm(lcc,lfc,lec,lnc,lb,iopt,idest)
Arguments:
lcc input logical Routine will act on cell centered data
if this is .true., otherwise set to
.false.
lfc input logical Routine will act on cell face-centered
data if this is .true., otherwise set to
.false.
lec input logical Routine will act on cell edge-centered
data if this is .true., otherwise set to
.false.
lnc input logical Routine will act on cell corner
data if this is .true., otherwise set to
.false.
lb input integer The block number for which guardcell
data is required.
iopt input integer If iopt=1 then AMR_GUARDCELL fills the
guard cell elements of the array UNK.
If iopt is greater than 2 then the
routine acts on the layer of the array
WORK with 5th index set to IOPT-1.
idest input integer The working block data structure has
two layers, one for the current block
and one for a copy of its parent.
These two copies are accesses by setting
idest=1 for the current block and
idest=2 for its parent.
Subroutines called by this routine:
None
This routine copies data from the 1-block working arrays with guardcells to the temporary copy T_UNK, T_FACEVAR*, of the permanent data arrays, which may or may not have permanent guardcells, depending on whether NO_PERMANENT_GUARDCELLS is defined in paramesh_preprocessor.fh.
subroutine amr_1blk_t_to_perm(lcc,lfc,lec,lnc,lb,idest)
Arguments:
lcc input logical Routine will act on cell centered data
if this is .true., otherwise set to
.false.
lfc input logical Routine will act on cell face-centered
data if this is .true., otherwise set to
.false.
lec input logical Routine will act on cell edge-centered
data if this is .true., otherwise set to
.false.
lnc input logical Routine will act on cell corner
data if this is .true., otherwise set to
.false.
lb input integer The block number for which guardcell
data is required.
idest input integer The working block data structure has
two layers, one for the current block
and one for a copy of its parent.
These two copies are accesses by setting
idest=1 for the current block and
idest=2 for its parent.
Subroutines called by this routine:
None
This routine computes length, area and volume properties of grid cells for the selected block, lb on processor pe. It is required if you are using curvilinear coordinates. Cartesian, polar(2D), cylindrical, and spherical coordinates are currently included. The cell geometry info for block (lb,pe) is stored in a number of arrays declared in the module PHYSICALDATA.F, cell volume in cell_vol, cell face areas in cell_area1, cell_area2, and cell_area3, and cell side lengths in cell_leng1, cell_leng2, and cell_leng3.
subroutine amr_block_geometry(lb,pe)
Arguments:
lb input integer The block number for which geometry
data will be computed.
pe input integer The processor address of the block
for which geometry data will be computed.
Subroutines called by this routine: None.
Subroutine to read checkpoint files written using amr_checkpoint_wr.
!!****f* mpi_source/amr_checkpoint_re
!! NAME
!!
!! amr_checkpoint_re
!!
!! SYNOPSIS
!!
!! call amr_checkpoint_re(file_num)
!! call amr_checkpoint_re(file_num, l_with_guardcells, check_format,
!! user_attr1, user_attr2, user_attr3,
!! user_attr4, user_attr5)
!!
!! call amr_checkpoint_re(integer, optional logical, optional char*80,
!! optional real, optional real, optional real,
!! optional real, optional real)
!!
!! ARGUMENTS
!!
!! integer, intent(in) :: file_num
!! An integer number which be appended to the end of the file name.
!!
!! optional, logical, intent(in) :: l_with_guardcells
!! If true, then guardcells are included in the checkpoint file. Otherwise
!! (the default) they are not included.
!!
!! optional, character(len=80), intent(in) :: check_format
!! Argument describing what type of output has been used. Currently only supports
!! fortran binary(the default), parallel hdf5 output, or parallel mpiio in
!! in native binary.
!! To read and hdf5 checkpoint file:
!! character(len=80) :: checkf
!! checkf = 'hdf5'
!! call amr_checkpoint_re(...., check_format=checkf, ....)
!!
!! optional, real, intent(out) :: user_attr1(2,3,4,5)
!! Arguments which allow the user to add some extra information to the file.
!! Currently only 5 real numbers can be added. Note also that the intent for these
!! arguments is different than for the amr_checkpoint_wr routine.
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!! mpif.h
!!
!! USES
!!
!! paramesh_interfaces
!!
!! CALLS
!!
!! amr_checkpoint_re_default
!! amr_checkpoint_re_hdf5
!! amr_checkpoint_re_mpiio
!! amr_abort
!!
!! RETURNS
!!
!! Does not return anything. Upon exit a checkpoint file has been read in.
!!
!! DESCRIPTION
!!
!! Subroutine to read checkpoint files written by amr_checkpoint_wr routine.
!! It read in the tree data structure and data stored in the PARAMESH blocks.
!! Optionally, a user may read a small amout of attribute data fime the file.
!! If the default bevaviour is selected, reads are done serially by processor 0.
!! The data is sent to the other processors from processor 0.
!! The default behaviour USES UNFORMATTED DIRECT I/O.
!! Optionally, the user can specify if the file will be read in using the HDF5
!! portable data format. If supported on the system, selecting hdf5 will
!! also result in the file will be read in parallel.
!!
!! The files read in must have names of the form 'paramesh_chk_######' or
!! 'paramesh_chk_######.hdf5'. where '######' is the file_num argument passed into
!! this routine.
!!
!! AUTHORS
!!
!! Peter MacNeice (1997) with modifications by Kevin Olson for parallel hdf5 and
!! mpiio (2004-2005).
!!
!!***
Subroutine to write checkpoint files to be read using amr_checkpoint_wr.
!! NAME
!!
!! amr_checkpoint_wr
!!
!! SYNOPSIS
!!
!! call amr_checkpoint_wr(file_num)
!! call amr_checkpoint_wr(file_num, l_with_guardcells, check_format,
!! user_attr1, user_attr2, user_attr3,
!! user_attr4, user_attr5)
!!
!! call amr_checkpoint_wr(integer, optional logical, optional char*80,
!! optional real, optional real, optional real,
!! optional real, optional real)
!!
!! ARGUMENTS
!!
!! integer, intent(in) :: file_num
!! An integer number which be appended to the end of the file name.
!!
!! optional, logical, intent(in) :: l_with_guardcells
!! If true, then guardcells are included in the checkpoint file. Otherwise
!! (the default) they are not included.
!!
!! optional, character(len=80), intent(in) :: check_format
!! Argument describing what type of output to use. Currently supports
!! fortran binary(the default), parallel hdf5 output, or parallel mpiio
!! in native binary.
!! To produce parallel hdf5 output:
!! character(len=80) :: checkf
!! checkf = 'hdf5'
!! call amr_checkpoint_wr(...., check_format=checkf, ....)
!!
!! optional, real, intent(in) :: user_attr1(2,3,4,5)
!! Arguments which allow the user to add some extra information to the file.
!! Currently only 5 real numbers can be added.
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!! mpif.h
!!
!! USES
!!
!! paramesh_interfaces
!!
!! CALLS
!!
!! amr_checkpoint_wr_default
!! amr_checkpoint_wr_hdf5
!! amr_checkpoint_wr_mpiio
!!
!! RETURNS
!!
!! Does not return anything. Upon exit a checkpoint file has been written.
!!
!! DESCRIPTION
!!
!! Subroutine to checkpoint runs using PARAMESH.
!! It writes out the tree data structure and data stored in PARAMESH blocks.
!! Optionally, a user may add a small amout of attribute data to the files written.
!! If the default bevaviour is selected, writes are done serially by processor 0.
!! I.e. data is collected from other processors and then written out.
!! The default behaviour USES UNFORMATTED DIRECT I/O.
!! Optionally, the user can specify if the file will be written using the HDF5
!! portable data format. If supported on the system, selecting hdf5 output will
!! also result in the file will be written in parallel.
!!
!! The files produced will have names of the form 'paramesh_chk_######' or
!! 'paramesh_chk_######.hdf5'. where '######' is the file_num argument passed into
!! this routine.
!!
!! AUTHORS
!!
!! Peter MacNeice (1997) with modifications by Kevin Olson for parallel hdf5 and
!! mpiio (2004-2005).
!!
!!***
This subroutine closes the amr package. This does any tidying up which is required. The most notable task performed is closing the mpi package if you are using the mpi version.
subroutine amr_close()
Arguments: None.
Subroutines called by this routine:
comm_finish
This routine is responsible for ensuring that variables defined on the edges of cells are used in a consistent way at the boundaries between grid blocks at different refinement levels. The edge values at the boundary with a neighboring block at higher refinement level are replaced with the average of the equivalent edge values computed on the more refined block. This will often be required to guarantee conservation, and in some MHD codes may be required to ensure div B remains 0. This routine is really only a wrapper which calls one or other of the routines AMR_EDGE_AVERAGE_UDT or AMR_EDGE_AVERAGE_VDT. If the preprocessor variable VAR_DT has not been defined, then the routine AMR_EDGE_AVERAGE_UDT is used. This assumes that a uniform timestep is being used everywhere. If VAR_DT is defined (inside paramesh_preprocessor.fh) then different timesteps may be used on different blocks. Details of the variations that are permitted are discussed in the Users Guide. In this case AMR_EDGE_AVERAGE_VDT would be used.
!! NAME
!!
!! amr_edge_average
!!
!! SYNOPSIS
!!
!! call amr_edge_average(mype, lfullblocks, nsub)
!!
!! call amr_edge_average(integer, logical, integer)
!!
!! ARGUMENTS
!!
!! integer, intent(in) :: mype
!! The calling processors id.
!!
!! integer, intent(in) :: nsub
!! Controls something for variable timestep algorithms ???.
!!
!! logical, intent(in) :: lfullblock
!! Controls which arrays this routine works on. If .true. then the unk_e_x,y,z
!! arrays are used. Otherwise, the bedge_* arrays are used.
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!! mpif.h
!!
!! USES
!!
!! paramesh_dimensions
!! physicaldata
!! tree
!! paramesh_interfaces
!!
!! CALLS
!!
!! amr_edge_average_udt
!! amr_edge_diagonal_check
!! amr_edge_average_vdt
!!
!! RETURNS
!!
!! Nothing returned.
!!
!! DESCRIPTION
!!
!! This is a wrapper routine which makes the appropriate calls to the
!! routines which manage edge data consistency at the boundaries between
!! grid blocks of different refinement level. This routine is called by the user
!! when they wish to enforce consistency of edge centered variables (either stored
!! in the bedge_* arrays or in the unk_e_* arrays) at jumps in refinement. This
!! routine is analogous to the subroutine 'amr_flux_conserve'.
!!
!! If lfullblock is true, then this routine will update the main edge-centered
!! datastructure arrays, unk_e_x[y][z], otherwise it simply operates on
!! the temporary edge-centered data computed on block boundary faces only.
!!
!! AUTHORS
!!
!! Written by Peter MacNeice (July 1997) and modified February 2001.
!!
!!***
This routine is responsible for ensuring that fluxes are used in a consistent way at the boundaries between grid blocks at different refinement levels. The fluxes at the boundary with a neighboring block at higher refinement level are replaced with the equivalent fluxes computed on the more refined block. This will often be required to guarantee conservation. This routine is really only a wrapper which calls one or other of the routines AMR_FLUX_CONSERVE_UDT or AMR_FLUX_CONSERVE_VDT. If the preprocessor variable VAR_DT has not been defined, then the routine AMR_FLUX_CONSERVE_UDT is used. This assumes that a uniform timestep is being used everywhere. If VAR_DT is defined (inside paramesh_preprocessor.fh) then different timesteps may be used on different blocks. Details of the variations that are permitted are discussed in the Users Guide. In this case AMR_FLUX_CONSERVE_VDT would be used.
!! NAME
!!
!! amr_flux_conserve
!!
!! SYNOPSIS
!!
!! call amr_flux_conserve (mype, nsub)
!! call amr_flux_conserve (mype, nsub, flux_dir)
!!
!! call amr_flux_conserve (integer, integer, optional integer)
!!
!! ARGUMENTS
!!
!! integer, intent(in) :: mype
!! The calling processor.
!!
!! integer, intent(in) :: nsub
!! The current time subcycle. If this is 1 then this info is used to
!! reset the temporary boundary flux arrays to 0. This argument only has
!! an effect if variable time steps are being used.
!!
!! optional, integer, intent(in) :: flux_dir
!! Option integer which selects which coordinate direction to apply
!! the flux conservation operation to:
!! If flux_dir = 1 -> x direction
!! flux_dir = 2 -> y direction
!! flux_dir = 3 -> z direction
!! If this argument is not specified, then the default behaviour is
!! to operate on all directions. Using this argument can be useful
!! for Strang-split schemes to improve performance.
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!!
!! USES
!!
!! paramesh_dimensions
!! physicaldata
!! tree
!!
!! CALLS
!!
!! amr_flux_conserve_udt
!! amr_flux_conserve_vdt
!!
!! RETURNS
!!
!! Does not return anything. Upon exit the fluxes stored in the arrays
!! flux_x, flux_y, or flux_z are corrected at jumps in refinement. This
!! is either an averaging proceedure or a sum as selected by the user
!! by adjusting the preprocessor variables which control this behaviour
!! in 'paramesh_preprocessor.fh'. This can also be controled through the
!! equivalent logical variables read at runtime from the file
!! 'amr_runtime_parameters' if the 'LIBRARY' preprocessor flag is defined.
!!
!! DESCRIPTION
!!
!! This is a wrapper routine which makes the appropriate call to the
!! routines which manage flux conservation at the boundaries between
!! grid blocks of different refinement level.
!!
!! These routines get block boundary data from neighbors who are
!! parents of leaf blocks. This is required in flux conserving schemes
!! where the coarser block needs to use the same fluxes and mean pressures
!! as will be used on the finer blocks across their shared boundary.
!!
!! AUTHORS
!!
!! Peter MacNeice (1997) with modifications by Kevin Olson for
!! directional guardcell filling.
!!
!!***
This routine globally fills guard cells with the data which they require from neighboring grid blocks. It only operates if NO_PERMANENT_GUARDCELLS is undefined. Basically, the routine calls AMR_RESTRICT, and then loops over grid blocks calling AMR_1BLK_GUARDCELL for each in turn.
AMR_GUARDCELL can operate on two different datastructures, UNK etc, or WORK, depending on the value assigned to IOPT. Occassionally it may be useful to be able to compute with one particular data word at each grid cell, such as density, for example. By copying the density into the array WORK, you can fill WORK's guard cells and perform this computation without needlessly filling the guard cell values for all the other data words stored in UNK etc.
Note, the routine AMR_INITIALIZE must be called once, before any calls to AMR_GUARDCELL can be made.
!! NAME
!!
!! amr_guardcell
!!
!! SYNOPSIS
!!
!! call amr_guarcell(mype, iopt, nlayers)
!! call amr_guarcell(mype, iopt, nlayers, nlayersx, nlayersy, nlayersz)
!!
!! call amr_guarcell(integer, integer, integer,
!! optional integer, optional integer, optional integer)
!!
!! ARGUMENTS
!!
!! integer, intent(in) :: mype
!! The calling processor.
!!
!! integer, intent(in) :: iopt
!! Selects whether to fill the guardcells for the arrays unk,
!! facevarx, facevary, facevarz, unk_e_x, unk_e_y, unk_e_z, and unk_n
!! (if iopt = 1) or work (if iopt 2).
!!
!! integer, intent(in) :: nlayers
!! Dummy variable which does nothing. Included for consistency with older
!! PARAMESH versions.
!!
!! optional, integer, intent(in) :: nlayersx, nlayersy, nlayersz
!! Optional integers which select how many guardcell layers to fill in each
!! direction. If these varaibles are not passed in, then the default is to
!! fill all allocated guardcell space.
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!!
!! USES
!!
!! paramesh_dimensions
!! physicaldata
!! workspace
!! tree
!! paramesh_interfaces
!! paramesh_mpi_interfaces
!!
!! CALLS
!!
!! amr_1blk_guardcell_reset
!! amr_restrict
!! amr_1blk_guardcell
!! mpi_amr_comm_setup
!!
!! RETURNS
!!
!! Does not return anything. Upon exit, guardcells are filled with data for
!! all blocks.
!!
!! DESCRIPTION
!!
!! Routine to perform guardcell filling in the case that permanent guardcells
!! are user. Calling this routine will result in the guardcells being filled
!! for all child blocks if 'advance_all_levels' is turned OFF or in ALL the blocks
!! having their guardcells filled if 'advance_all_levels' is turned ON.
!! The user can choose how many layers of guardcells are filled on each face
!! of a block by specifying the optional arguments 'nlayersx', 'nlayersy', and
!! 'nlayersz'.
!!
!! AUTHORS
!!
!! Peter MacNeice (1997) with modifications by Kevin Olson
!!
!!***
This subroutine performs application independent initializations of the amr package. Amongst other things, this will pre-compute and store some of the coefficients used in the interpolation operations performed inside AMR_PROLONG. This routine must be called before any calls to AMR_PROLONG or AMR_GUARDCELL are issued.
During the initialization of the coefficients used in prolongation, the user must choose whether to implement the special `conservative' treatment for the grid cells which touch a block boundary. The pre-processor variable CONSERVE controls this choice, and should be defined in PARAMESH_PREPROCESSOR.FH if this option is required. This issue is explained in detail in the listing for the routine AMR_PROLONG_FUN_INIT below.
subroutine amr_initialize()
Arguments: None.
Subroutines called by this routine:
amr_1blk_guardcell_reset
amr_prolong_fun_init
amr_bcset_init
comm_start
This subroutine copies data from the permanent solution datastructure into the working block.
subroutine amr_perm_to_1blk( lcc,lfc,lec,lnc,lb,pe,iopt,idest)
Arguments:
lcc input logical Routine will act on cell centered data
if this is .true., otherwise set to
.false.
lfc input logical Routine will act on cell face-centered
data if this is .true., otherwise set to
.false.
lec input logical Routine will act on cell edge data
if this is .true., otherwise set to
.false.
lnc input logical Routine will act on cell corner
data if this is .true., otherwise set to
.false.
lb input integer The block number for which guardcell
data is required.
pe input integer The processor on which block lb is stored.
iopt input integer If iopt=1 then AMR_GUARDCELL fills the
guard cell elements of the array UNK.
If iopt is greater than 2 then the
routine acts on the layer of the array
WORK with 5th index set to IOPT-1.
idest input integer The working block data structure has
two layers, one for the current block
and one for a copy of its parent.
These two copies are accesses by setting
idest=1 for the current block and
idest=2 for its parent.
Subroutines called by this routine:
amr_mpi_find_blk_in_buffer
mpi_set_message_limits
This routine initializes some or all of the arrays UNK, FACEVAR, UNK_E and UNK_N, or WORK on newly created leaf blocks. It does this by interpolating from the same arrays on the parent of the newly created block. AMR_PROLONG can operate on two different datastructures, UNK etc, or WORK, depending on the value assigned to IOPT. The interpolation involved in this process is actually performed through calls to the set of routines AMR_1BLK_XX_PROL_GEN...., when XX is one of CC, FC, EC NC. These routines are wrapper routines which in turn can call a selection of different interpolation functions. The user can select a particular function for each component of the solution by editing these wrapper routines and following the detailed comments in the routine listings. This structure enables a user to easily add their own interpolation function if they wish.
Note, the routine AMR_INITIALIZE must be called once, before any calls to AMR_PROLONG can be made.
!! NAME
!!
!! amr_prolong
!!
!! SYNOPSIS
!!
!! call amr_prolong (mype, iopt, nlayers)
!!
!! call amr_prolong (integer, integer, integer)
!!
!! ARGUMENTS
!!
!! integer, intent(in) :: mype
!! Current processor number
!!
!! integer, intent(in) :: iopt
!! Switch to select which datastructures are updated. If iopt=1
!! then this routine acts on 'unk', 'facevarx(y,z)', 'unk_e_x(y,z)',
!! and 'unk_n'.
!! If iopt=2 only 'work' is updated.
!!
!! integer, intent(in) :: nlayers
!! Number of layers of guard cells at a block boundary.
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!!
!! USES
!!
!! paramesh_dimensions
!! physicaldata
!! tree
!! workspace
!! Mpiv_morton
!! paramesh_mpi_interfaces
!! paramesh_interfaces
!!
!! CALLS
!!
!! amr_1blk_cc_prol_gen_work_fun
!! amr_1blk_fc_prol_gen_fun
!! amr_1blk_ec_prol_gen_fun
!! amr_1blk_nc_prol_gen_fun
!! amr_1blk_copy_soln
!! amr_1blk_guardcell_reset
!! amr_1blk_guardcell
!! amr_1blk_cc_prol_gen_unk_fun
!! comm_int_max_to_all
!! comm_int_min_to_all
!! amr_1blk_fc_prol_dbz
!! mpi_amr_comm_setup
!!
!! RETURNS
!!
!! Upon exit, prolongation (i.e. interpolation) from parent blocks to
!! their newly created child (marked by the 'newchild' flag) has been
!! performed.
!!
!! DESCRIPTION
!!
!! This routine interpolates data from parent blocks to any newly created
!! child blocks which are marked with the 'newchild' flag set to true.
!!
!! AUTHORS
!!
!! Peter MacNeice (July 1997)
!!
!! Modified by Michael L. Rilee, November 2002, *dbz*
!! Initial support for divergenceless prolongation
!! Modified by Michael L. Rilee, December 2002, *clean_divb*
!! Support for projecting field onto divergenceless field
!!
!!***
This routine performs any refinement or derefinement that has been requested. It rebuilds the block tree to accomodate these changes, and can distribute the blocks amongst the processors using a Peano-Hilbert space filling curve which maximizes the locality associated with inter-block communications.
This routine guarantees that the refinement level will not change by more than one level across any block boundary in the computational domain.
!! NAME
!!
!! amr_refine_derefine
!!
!! SYNOPSIS
!!
!! call amr_refine_derefine()
!!
!!
!! ARGUMENTS
!!
!! NO ARGUMENTS
!!
!! INCLUDES
!!
!! paramesh_preprocessor.fh
!! mpif.h
!!
!! USES
!!
!! paramesh_dimensions
!! physicaldata
!! tree
!! mpi_morton
!! timings
!! io
!! paramesh_interfaces
!!
!! CALLS
!!
!! amr_check_refine
!! amr_check_derefine
!! amr_refine_blocks
!! amr_derefine_blocks
!! amr_morton_order
!!
!! RETURNS
!!
!! Does not return anything. Upon exit blocks marked for refinement have been
!! created, blocks marked for derefinement have been eliminated, and the blocks
!! have been reordered to acheive load balance.
!!
!! DESCRIPTION
!!
!! Subroutine to refine or derefine blocks. This routine is called by the user
!! who sets the refine or derefine flags to be true or false. These flags are
!! logical variables called 'refine' and 'derefine' and are stored for each block.
!! If a block is marked for refinement, amr_refine_derefine will create its new
!! child blocks. Also, tests will be executed to see if any other blocks need to
!! also refine (by calling amr_check_refine) to ensure that the a jump in refinement of
!! more than one level is not created.
!!
!! If a block is marked for derefinement, amr_refine_derefine first checks to make
!! that the block can derefine and not create a jump in refinement of more than one
!! level. A check is also run to check that all the siblings of the derefining
!! block's siblings are also marked for derefinement. If these tests succeed, the
!! block is removed from the list of blocks.
!!
!! Once these operations are completed, the routine 'amr_morton_order' is called
!! and the tree data structure is reorganized to acheive load balance using a
!! morton space filling curve. After this routine is called, the routine
!! 'amr_redist_blk' is called, which actually moves the block data into the correct
!! positions in the morton order list of blocks.
!!
!! Finally, the routine 'amr_morton_process' is called. This routine computes the
!! communications patterns needed for guardcell filling, restriction, and
!! prologation and stores them for later use.
!!
!! AUTHORS
!!
!! Kevin Olson (1997)
!!
!!***
This routine reorders an existing set of blocks across the processors using a Peano-Hilbert space filling curve which maximizes the locality associated with inter-block communications. It does not move solution data. Its primary use is in ordering the blocks when creating an initial grid with multiple root blocks by hand.
subroutine AMR_REORDER_GRID()
Arguments: None.
Subroutines called by this routine:
amr_morton_order
amr_morton_process
This routine updates the arrays UNK or WORK on the parents of leaf blocks by averaging the data provided by the leaf blocks. Each leaf block provides its data to its parent. The parent block then computes an averaged representation of this data to fill the appropriate half(in 1D) quadrant (in 2D) or octant (in 3D). The averaging is performed through calls to a number of restriction routines listed immediately below.
!! NAME
!!
!! amr_restrict
!!