diff --git a/src/blocks.F90 b/src/blocks.F90 index a69e350..b135ce9 100644 --- a/src/blocks.F90 +++ b/src/blocks.F90 @@ -23,27 +23,67 @@ !! !! module: BLOCKS !! -!! This module allocates, deallocates, and handles blocks of the adaptive mesh -!! structures. +!! This module provides data structures, variables and subroutines to +!! construct and dynamically modify the hierarchy of blocks corresponding +!! to the simulated mesh geometry. !! !!****************************************************************************** ! module blocks +#ifdef PROFILE +! import external subroutines +! + use timers, only : set_timer, start_timer, stop_timer +#endif /* PROFILE */ + ! module variables are not implicit by default ! implicit none -! module parameters +#ifdef PROFILE +! timer indices +! + integer, save :: imi, ima, imu, imp, imq, imr, imd +#endif /* PROFILE */ + +! MODULE PARAMETERS: +! ================= +! +! ndims - the number of dimensions (2 or 3); +! nsides - the number of sides along each direction (2); +! nfaces - the number of faces at each side (2 for 2D, 4 for 3D); +! nchildren - the number of child blocks for each block (4 for 2D, 8 for 3D); ! integer(kind=4), parameter :: ndims = NDIMS integer(kind=4), parameter :: nsides = 2 integer(kind=4), parameter :: nfaces = 2**(ndims - 1) integer(kind=4), parameter :: nchildren = 2**ndims -!! BLOCK STRUCTURE POINTERS (they have to be defined before block structures) -!! -! define pointers to meta, data, and info block structures +! MODULE VARIABLES: +! ================ +! +! the identification of the last allocated block (always increases) +! + integer(kind=4), save :: last_id + +! the number of allocated meta and data blocks, and the number of leafs +! + integer(kind=4), save :: mblocks, dblocks, nleafs + +! the number of variables and fluxes stored in data blocks +! + integer(kind=4), save :: nvars, nflux + +! the spacial dimensions of allocatable data block arrays +! + integer(kind=4), save :: nx, ny, nz + +! BLOCK STRUCTURE POINTERS: +! ======================== +! +! define pointers to meta, data, and info block structures defined below; +! they have to be defined before block structures ! type pointer_meta type(block_meta), pointer :: ptr @@ -57,9 +97,12 @@ module blocks type(block_info), pointer :: ptr end type pointer_info -!! BLOCK STRUCTURES +! BLOCK STRUCTURES: +! ================ ! -! define the META block structure +! define the META block structure; each process keeps exactly same meta block +! structure all the time, so processes can know how the block structure changes +! and where to move data blocks; ! type block_meta ! pointers to the previous and next meta blocks @@ -82,48 +125,63 @@ module blocks ! type(block_data) , pointer :: data - ! the block identification + ! the identification number (unique for each + ! block) ! integer(kind=4) :: id - ! the number of associated cpu + ! the process number to which the meta block + ! is bounded ! - integer(kind=4) :: cpu + integer(kind=4) :: process ! the level of refinement ! integer(kind=4) :: level - ! the configuration flag for its children order + ! the number describing the configuration of + ! the child meta blocks ! - integer(kind=4) :: config + integer(kind=4) :: conf ! the refinement flag, -1, 0, and 1 for - ! derefinement, no change, and refinement, - ! respectively + ! the block marked to be derefined, not + ! changed, and refined, respectively ! integer(kind=4) :: refine - ! the position of the block in its parent block + ! the position of the block in its siblings + ! group ! integer(kind=4) :: pos(ndims) ! the coordinate of the lower corner of the ! block in the effective resolution units ! - integer(kind=4) :: coord(ndims) + integer(kind=4) :: coords(ndims) - ! the leaf flag + ! the leaf flag, signifying that the block is + ! the highest block in the local block + ! structure ! logical :: leaf - ! the block coordinates in the physical units + ! the flag indicates that the corresponding + ! data needs to be updated (e.g. boundaries or + ! primitive variables), therefore it is + ! usually .true. ! - real :: xmin, xmax, ymin, ymax, zmin, zmax + logical :: update + + ! the block bounds in the coordinate units + ! + real(kind=8) :: xmin, xmax, ymin, ymax, zmin, zmax end type block_meta -! define the DATA block structure +! define the DATA block structure; all data blocks are divided between +! processes, therefore the same data block cannot be associated with two +! different processes, but they can be moved from one process to another; ! type block_data ! pointers to the previous and next data blocks @@ -134,12 +192,14 @@ module blocks ! type(block_meta), pointer :: meta - ! a pointer to the array conserved variables + ! a pointer to the current conserved variable + ! array ! real, dimension(:,:,:,:) , pointer :: u ! an allocatable arrays to store all conserved - ! variables + ! variables (required two for Runge-Kutta + ! temporal integration methods) ! real, dimension(:,:,:,:) , allocatable :: u0, u1 @@ -152,13 +212,6 @@ module blocks ! real, dimension(:,:,:,:,:), allocatable :: f -#ifdef DEBUG - ! an allocatable array to store refinement - ! values - ! - real, dimension(:,:,:) , allocatable :: c -#endif /* DEBUG */ - end type block_data ! define the INFO block structure @@ -176,7 +229,9 @@ module blocks ! type(block_meta) , pointer :: neigh - ! the direction, side and face indices + ! the direction, side and face numbers + ! indicating the neighbor block orientation + ! with respect to the block ! integer(kind=4) :: direction, side, face @@ -187,87 +242,98 @@ module blocks end type block_info -!! POINTER TO THE FIST AND LAST BLOCKS IN THE LISTS -!! -! chains of meta blocks and data blocks +! POINTERS TO THE FIST AND LAST BLOCKS IN THE LISTS: +! ================================================= +! +! these pointers construct the lists of meta and data blocks; ! type(block_meta), pointer, save :: list_meta, last_meta type(block_data), pointer, save :: list_data, last_data -!! MODULE VARIABLES -!! -! the identification of the last allocated block (should always increase) -! - integer(kind=4) , save :: last_id - -! the numbers of allocated meta and data blocks, and leafs -! - integer(kind=4) , save :: mblocks, dblocks, nleafs - -! the numbers of variables and fluxes stored in data blocks -! - integer(kind=4) , save :: nvars, nflux - -! the spacial dimensions of data block allocatable arrays -! - integer(kind=4) , save :: nx, ny, nz - ! all variables and subroutines are private by default ! private -! declare public subroutines +! declare public pointers, structures, and variables ! public :: pointer_meta, pointer_info public :: block_meta, block_data, block_info public :: list_meta, list_data - public :: nchildren, ndims, nsides, nfaces - public :: initialize_blocks, finalize_blocks - public :: set_last_id, get_last_id, get_mblocks, get_dblocks, get_nleafs - public :: link_blocks, unlink_blocks - public :: append_metablock - public :: allocate_datablock, deallocate_datablock - public :: append_datablock, remove_datablock - public :: metablock_set_id, metablock_set_cpu, metablock_set_refine & - , metablock_set_config, metablock_set_level, metablock_set_position & - , metablock_set_coord, metablock_set_bounds, metablock_set_leaf - public :: datablock_set_dims - public :: refine_block, derefine_block -#ifdef DEBUG - public :: check_metablock -#endif /* DEBUG */ + public :: ndims, nsides, nfaces, nchildren +! declare public subroutines +! + public :: initialize_blocks, finalize_blocks + public :: set_block_dimensions + public :: append_metablock, remove_metablock + public :: append_datablock, remove_datablock + public :: allocate_metablock, deallocate_metablock + public :: allocate_datablock, deallocate_datablock + public :: link_blocks, unlink_blocks + public :: refine_block, derefine_block + public :: set_last_id, get_last_id, get_mblocks, get_dblocks, get_nleafs + public :: set_blocks_update + public :: metablock_set_id, metablock_set_process, metablock_set_level + public :: metablock_set_configuration, metablock_set_refinement + public :: metablock_set_position, metablock_set_coordinates + public :: metablock_set_bounds, metablock_set_leaf, metablock_unset_leaf + +!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - +! contains ! -!!============================================================================== +!=============================================================================== !! -!! INITIALIZATION/FINALIZATION SUBROUTINES +!!*** PUBLIC SUBROUTINES ***************************************************** !! +!=============================================================================== ! !=============================================================================== ! ! subroutine INITIALIZE_BLOCKS: ! ---------------------------- ! -! Subroutine initializes the variables related to the elementary block of -! the adaptive structure. +! Subroutine initializes the module structures, pointers and variables. +! +! Arguments: +! +! verbose - flag determining if the subroutine should be verbose; +! iret - return flag of the procedure execution status; ! !=============================================================================== ! - subroutine initialize_blocks() + subroutine initialize_blocks(verbose, iret) ! local variables are not implicit by default ! implicit none + +! subroutine arguments +! + logical, intent(in) :: verbose + integer, intent(inout) :: iret ! !------------------------------------------------------------------------------- ! -! nullify list pointers +#ifdef PROFILE +! set timer descriptions ! - nullify(list_meta) - nullify(list_data) - nullify(last_meta) - nullify(last_data) + call set_timer('blocks:: initialization' , imi) + call set_timer('blocks:: meta block allocation' , ima) + call set_timer('blocks:: meta block deallocation', imu) + call set_timer('blocks:: data block allocation' , imp) + call set_timer('blocks:: data block deallocation', imq) + call set_timer('blocks:: refine' , imr) + call set_timer('blocks:: derefine' , imd) + +! start accounting time for module initialization/finalization +! + call start_timer(imi) +#endif /* PROFILE */ + +! reset identification counter +! + last_id = 0 ! reset the number of meta blocks, data blocks, and leafs ! @@ -278,7 +344,7 @@ module blocks ! set the initial number of variables and fluxes ! nvars = 1 - nflux = 0 + nflux = 1 ! set the initial data block resolution ! @@ -286,9 +352,18 @@ module blocks ny = 1 nz = 1 -! reset identification counter +! nullify pointers defining the meta and data lists ! - last_id = 0 + nullify(list_meta) + nullify(list_data) + nullify(last_meta) + nullify(last_data) + +#ifdef PROFILE +! stop accounting time for module initialization/finalization +! + call stop_timer(imi) +#endif /* PROFILE */ !------------------------------------------------------------------------------- ! @@ -302,57 +377,388 @@ module blocks ! Subroutine iterates over all meta blocks and first deallocates all ! associated with them data blocks, and then their metadata structure. ! +! Arguments: +! +! iret - return flag of the procedure execution status; +! !=============================================================================== ! - subroutine finalize_blocks() + subroutine finalize_blocks(iret) ! local variables are not implicit by default ! implicit none -! a pointer to the current meta block +! subroutine arguments +! + integer, intent(inout) :: iret + +! local variables ! type(block_meta), pointer :: pmeta ! !------------------------------------------------------------------------------- ! -! assiociate pmeta pointer with the first block in the list +#ifdef PROFILE +! start accounting time for module initialization/finalization ! - pmeta => list_meta + call start_timer(imi) +#endif /* PROFILE */ +! associate the pointer with the last block on the meta block list +! + pmeta => last_meta + +! iterate until the first block on the list is reached +! do while(associated(pmeta)) -! deallocate current meta block +! deallocate the last meta block ! - call deallocate_metablock(pmeta) + call remove_metablock(pmeta) -! associate pmeta pointer with the next meta block in the list +! assign the pointer to the last block on the meta block list ! - pmeta => list_meta + pmeta => last_meta - end do + end do ! meta blocks + +! nullify pointers defining the meta and data lists +! + nullify(list_meta) + nullify(list_data) + nullify(last_meta) + nullify(last_data) + +#ifdef PROFILE +! stop accounting time for module initialization/finalization +! + call stop_timer(imi) +#endif /* PROFILE */ !------------------------------------------------------------------------------- ! end subroutine finalize_blocks ! -!!============================================================================== -!! -!! META BLOCK SUBROUTINES -!! +!=============================================================================== +! +! subroutine SET_BLOCK_DIMENSIONS: +! ------------------------------- +! +! Subroutine sets the number of variables, fluxes and block dimensions +! (without ghost cells) for arrays allocated in data blocks. +! +! Arguments: +! +! nv - the number of variables stored in %u and %q; +! nf - the number of fluxes stored in %f; +! ni - the block dimension along X; +! nj - the block dimension along Y; +! nk - the block dimension along Z; ! !=============================================================================== ! -! allocate_metablock: subroutine allocates space for one meta block and returns -! the pointer to this block + subroutine set_block_dimensions(nv, nf, ni, nj, nk) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + integer(kind=4), intent(in) :: nv, nf, ni, nj, nk +! +!------------------------------------------------------------------------------- +! +! set the number of variables and fluxes +! + nvars = nv + nflux = nf + +! set the block dimensions +! + nx = ni + ny = nj +#if NDIMS == 3 + nz = nk +#endif /* NDIMS == 3 */ + +!------------------------------------------------------------------------------- +! + end subroutine set_block_dimensions +! +!=============================================================================== +! +! subroutine APPEND_METABLOCK: +! --------------------------- +! +! Subroutine allocates memory for one meta block, appends it to the meta +! block list and returns a pointer associated with it. +! +! Arguments: +! +! pmeta - the pointer associated with the newly appended meta block; +! +!=============================================================================== +! + subroutine append_metablock(pmeta) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(out) :: pmeta +! +!------------------------------------------------------------------------------- +! +! allocate memory for the new meta block +! + call allocate_metablock(pmeta) + +! check if there are any blocks in the meta block list +! + if (associated(last_meta)) then + +! add the new block to the end of the list +! + pmeta%prev => last_meta + last_meta%next => pmeta + + else + +! there are no blocks in the list, so add this one as the first block +! + list_meta => pmeta + + end if + +! update the pointer to the last block on the list +! + last_meta => pmeta + +!------------------------------------------------------------------------------- +! + end subroutine append_metablock +! +!=============================================================================== +! +! subroutine REMOVE_METABLOCK: +! --------------------------- +! +! Subroutine removes a meta block associated with the input pointer from +! the meta block list, and deallocates space used by it. +! +! Arguments: +! +! pmeta - the pointer pointing to the meta block which will be removed; +! +!=============================================================================== +! + subroutine remove_metablock(pmeta) + +! import external procedures +! + use error , only : print_error + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta +! +!------------------------------------------------------------------------------- +! +! check if the pointer is actually associated with any block +! + if (associated(pmeta)) then + +! if this is the first block in the list, update the list_meta pointer +! + if (pmeta%id == list_meta%id) list_meta => pmeta%next + +! if this is the last block in the list, update the last_meta pointer +! + if (pmeta%id == last_meta%id) last_meta => pmeta%prev + +! update the %next and %prev pointers of the previous and next blocks, +! respectively +! + if (associated(pmeta%prev)) pmeta%prev%next => pmeta%next + if (associated(pmeta%next)) pmeta%next%prev => pmeta%prev + +! deallocate memory used by the meta block +! + call deallocate_metablock(pmeta) + + else + +! the argument contains a null pointer, so print an error +! + call print_error("blocks::remove_metablock" & + , "Null pointer argument to meta block!") + end if + +!------------------------------------------------------------------------------- +! + end subroutine remove_metablock +! +!=============================================================================== +! +! subroutine APPEND_DATABLOCK: +! --------------------------- +! +! Subroutine allocates memory for one data block, appends it to the data +! block list and returns a pointer associated with it. +! +! Arguments: +! +! pdata - the pointer associated with the newly appended data block; +! +!=============================================================================== +! + subroutine append_datablock(pdata) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_data), pointer, intent(out) :: pdata +! +!------------------------------------------------------------------------------- +! +! allocate memory for the new data block +! + call allocate_datablock(pdata) + +! check if there are any blocks in the data block list +! + if (associated(last_data)) then + +! add the new block to the end of the list +! + pdata%prev => last_data + last_data%next => pdata + + else + +! there are no blocks in the list, so add this one as the first block +! + list_data => pdata + + end if + +! update the pointer to the last block on the list +! + last_data => pdata + +!------------------------------------------------------------------------------- +! + end subroutine append_datablock +! +!=============================================================================== +! +! subroutine REMOVE_DATABLOCK: +! --------------------------- +! +! Subroutine removes a data block associated with the input pointer from +! the data block list, and deallocates space used by it. +! +! Arguments: +! +! pdata - the pointer pointing to the data block which will be removed; +! +!=============================================================================== +! + subroutine remove_datablock(pdata) + +! import external procedures +! + use error , only : print_error + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_data), pointer, intent(inout) :: pdata +! +!------------------------------------------------------------------------------- +! +! check if the pointer is actually associated with any block +! + if (associated(pdata)) then + +! check if the data block has associated meta block +! + if (associated(pdata%meta)) then + +! if this is the first block in the list, update the list_data pointer +! + if (pdata%meta%id == list_data%meta%id) list_data => pdata%next + +! if this is the last block in the list, update the last_data pointer +! + if (pdata%meta%id == last_data%meta%id) last_data => pdata%prev + +! update the %next and %prev pointers of the previous and next blocks, +! respectively +! + if (associated(pdata%prev)) pdata%prev%next => pdata%next + if (associated(pdata%next)) pdata%next%prev => pdata%prev + + else ! %meta associated + +! there is no meta block associated, so print an error +! + call print_error("blocks::remove_datablock" & + , "No meta block associated with the data block!") + + end if ! %meta associated + +! deallocate the associated data block +! + call deallocate_datablock(pdata) + + else + +! the argument contains a null pointer, so print an error +! + call print_error("blocks::remove_datablock" & + , "Null pointer argument to data block!") + end if + +!------------------------------------------------------------------------------- +! + end subroutine remove_datablock +! +!=============================================================================== +! +! subroutine ALLOCATE_METABLOCK: +! ----------------------------- +! +! Subroutine allocates memory for one meta block, initializes its fields +! and returns a pointer associated with it. +! +! Arguments: +! +! pmeta - the pointer associated with the newly allocated meta block; ! !=============================================================================== ! subroutine allocate_metablock(pmeta) +! local variables are not implicit by default +! implicit none -! output arguments +! subroutine arguments ! type(block_meta), pointer, intent(out) :: pmeta @@ -362,60 +768,85 @@ module blocks ! !------------------------------------------------------------------------------- ! -! allocate block structure +#ifdef PROFILE +! start accounting time for the meta block allocation +! + call start_timer(ima) +#endif /* PROFILE */ + +! allocate the meta block structure for one object ! allocate(pmeta) -! nullify pointers +! nullify fields pointing to previous and next block on the meta block list ! nullify(pmeta%prev) nullify(pmeta%next) + +! nullify the field pointing to the parent +! nullify(pmeta%parent) - nullify(pmeta%data) + +! nullify fields pointing to children +! do i = 1, nchildren nullify(pmeta%child(i)%ptr) end do - do k = 1, nfaces + +! nullify fields pointing to neighbors +! + do i = 1, ndims do j = 1, nsides - do i = 1, ndims + do k = 1, nfaces nullify(pmeta%neigh(i,j,k)%ptr) end do end do end do +! nullify the field pointing to the associated data block +! + nullify(pmeta%data) + ! set unique ID ! - pmeta%id = increase_id() + pmeta%id = increase_id() -! unset the CPU number of current block, level, the configuration, refine and -! leaf flags +! unset the process number, level, the children configuration, refine, leaf, +! and update flags ! - pmeta%cpu = -1 - pmeta%level = -1 - pmeta%config = -1 - pmeta%refine = 0 - pmeta%leaf = .false. + pmeta%process = -1 + pmeta%level = -1 + pmeta%conf = -1 + pmeta%refine = 0 + pmeta%leaf = .false. + pmeta%update = .true. ! initialize the position in the parent block ! - pmeta%pos(:) = -1 + pmeta%pos(:) = -1 -! initialize the coordinate +! initialize the effective coordinates ! - pmeta%coord(:) = 0 + pmeta%coords(:) = 0 -! initialize bounds of the block +! initialize coordinate bounds of the block ! - pmeta%xmin = 0.0 - pmeta%xmax = 1.0 - pmeta%ymin = 0.0 - pmeta%ymax = 1.0 - pmeta%zmin = 0.0 - pmeta%zmax = 1.0 + pmeta%xmin = 0.0d+00 + pmeta%xmax = 1.0d+00 + pmeta%ymin = 0.0d+00 + pmeta%ymax = 1.0d+00 + pmeta%zmin = 0.0d+00 + pmeta%zmax = 1.0d+00 ! increase the number of allocated meta blocks ! - mblocks = mblocks + 1 + mblocks = mblocks + 1 + +#ifdef PROFILE +! stop accounting time for the meta block allocation +! + call stop_timer(ima) +#endif /* PROFILE */ !------------------------------------------------------------------------------- ! @@ -423,16 +854,30 @@ module blocks ! !=============================================================================== ! -! deallocate_metablock: subroutine deallocates space occupied by a given meta -! block +! subroutine DEALLOCATE_METABLOCK: +! ------------------------------- +! +! Subroutine releases memory used by the meta block associated with +! the pointer argument. +! +! Arguments: +! +! pmeta - the pointer associated with the meta block which will be +! deallocated; ! !=============================================================================== ! subroutine deallocate_metablock(pmeta) +! import external procedures +! + use error , only : print_error + +! local variables are not implicit by default +! implicit none -! input arguments +! subroutine arguments ! type(block_meta), pointer, intent(inout) :: pmeta @@ -441,275 +886,263 @@ module blocks integer :: i, j, k ! !------------------------------------------------------------------------------- +! +#ifdef PROFILE +! start accounting time for the meta block deallocation +! + call start_timer(imu) +#endif /* PROFILE */ + +! check if the pointer is actually associated with any block ! if (associated(pmeta)) then -! if this is the first block in the list, update the list_meta pointer +! decrease the number of leafs ! - if (pmeta%id .eq. list_meta%id) & - list_meta => pmeta%next + if (pmeta%leaf) nleafs = nleafs - 1 -! if this is the last block in the list, update the last_meta pointer +! decrease the number of allocated meta blocks ! - if (pmeta%id .eq. last_meta%id) & - last_meta => pmeta%prev + mblocks = mblocks - 1 -! update the pointer of previous and next blocks +! nullify fields pointing to previous and next block on the meta block list ! - if (associated(pmeta%prev)) & - pmeta%prev%next => pmeta%next + nullify(pmeta%prev) + nullify(pmeta%next) - if (associated(pmeta%next)) & - pmeta%next%prev => pmeta%prev +! nullify the field pointing to the parent +! + nullify(pmeta%parent) -! nullify children +! nullify fields pointing to children ! do i = 1, nchildren nullify(pmeta%child(i)%ptr) end do -! nullify neighbors +! nullify fields pointing to neighbors ! - do k = 1, nfaces + do i = 1, ndims do j = 1, nsides - do i = 1, ndims + do k = 1, nfaces nullify(pmeta%neigh(i,j,k)%ptr) end do end do end do -! if corresponding data block is allocated, deallocate it too +! if there is a data block is associated, remove it ! - if (associated(pmeta%data)) & - call remove_datablock(pmeta%data) + if (associated(pmeta%data)) call remove_datablock(pmeta%data) -! nullify pointers +! nullify the field pointing to the associated data block ! - nullify(pmeta%next) - nullify(pmeta%prev) nullify(pmeta%data) - nullify(pmeta%parent) -! free and nullify the block +! release the memory occupied by the block ! deallocate(pmeta) + +! nullify the pointer to the deallocated meta block +! nullify(pmeta) -! decrease the number of allocated blocks -! - mblocks = mblocks - 1 + else +! the argument contains a null pointer, so print an error +! + call print_error("blocks::deallocate_metablock" & + , "Null pointer argument to meta block!") end if +#ifdef PROFILE +! stop accounting time for the meta block deallocation +! + call stop_timer(imu) +#endif /* PROFILE */ + !------------------------------------------------------------------------------- ! end subroutine deallocate_metablock ! !=============================================================================== ! -! append_metablock: subroutine allocates space for one meta block and appends it -! to the meta block list +! subroutine ALLOCATE_DATABLOCK: +! ----------------------------- +! +! Subroutine allocates memory for one data block, initializes its fields +! and returns a pointer associated with it. +! +! Arguments: +! +! pdata - the pointer associated with the newly allocated data block; ! !=============================================================================== ! - subroutine append_metablock(pmeta) + subroutine allocate_datablock(pdata) +! local variables are not implicit by default +! implicit none -! output arguments +! subroutine arguments ! - type(block_meta), pointer, intent(out) :: pmeta + type(block_data), pointer, intent(out) :: pdata ! !------------------------------------------------------------------------------- ! -! allocate block +#ifdef PROFILE +! start accounting time for the data block allocation ! - call allocate_metablock(pmeta) + call start_timer(imp) +#endif /* PROFILE */ -! add to the list +! allocate the block structure ! - if (associated(last_meta)) then - pmeta%prev => last_meta - last_meta%next => pmeta + allocate(pdata) + +! nullify field pointing to the previous and next blocks on the data block list +! + nullify(pdata%prev) + nullify(pdata%next) + +! nullify the field pointing to the associate meta block list +! + nullify(pdata%meta) + +! allocate space for conserved variables +! + allocate(pdata%u0(nvars,nx,ny,nz), pdata%u1(nvars,nx,ny,nz)) + +! allocate space for primitive variables +! + allocate(pdata%q(nvars,nx,ny,nz)) + +! allocate space for numerical fluxes +! + allocate(pdata%f(ndims,nflux,nx,ny,nz)) + +! initiate the conserved variable pointer +! + pdata%u => pdata%u0 + +! increase the number of allocated data blocks +! + dblocks = dblocks + 1 + +#ifdef PROFILE +! stop accounting time for the data block allocation +! + call stop_timer(imp) +#endif /* PROFILE */ + +!------------------------------------------------------------------------------- +! + end subroutine allocate_datablock +! +!=============================================================================== +! +! subroutine DEALLOCATE_DATABLOCK: +! ------------------------------- +! +! Subroutine releases memory used by the data block associated with +! the pointer argument. +! +! Arguments: +! +! pdata - the pointer associated with the data block which will be +! deallocated; +! +!=============================================================================== +! + subroutine deallocate_datablock(pdata) + +! import external procedures +! + use error , only : print_error + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_data), pointer, intent(inout) :: pdata +! +!------------------------------------------------------------------------------- +! +#ifdef PROFILE +! start accounting time for the data block deallocation +! + call start_timer(imq) +#endif /* PROFILE */ + +! check if the pointer is actually associated with any block +! + if (associated(pdata)) then + +! decrease the number of allocated data blocks +! + dblocks = dblocks - 1 + +! nullify field pointing to the previous and next blocks on the data block list +! + nullify(pdata%prev) + nullify(pdata%next) + +! nullify the field pointing to the associate meta block list +! + nullify(pdata%meta) + +! nullify pointer to the current conserved variable array +! + nullify(pdata%u) + +! deallocate conserved variables +! + if (allocated(pdata%u0)) deallocate(pdata%u0) + if (allocated(pdata%u1)) deallocate(pdata%u1) + +! deallocate primitive variables +! + if (allocated(pdata%q )) deallocate(pdata%q ) + +! deallocate numerical fluxes +! + if (allocated(pdata%f )) deallocate(pdata%f ) + +! release the memory occupied by the block +! + deallocate(pdata) + +! nullify the pointer to the deallocated meta block +! + nullify(pdata) + else - list_meta => pmeta - end if -! set the pointer to the last block in the list +! the argument contains a null pointer, so print an error ! - last_meta => pmeta + call print_error("blocks::deallocate_datablock" & + , "Null pointer argument to data block!") + + end if ! pdata associated with a data block + +#ifdef PROFILE +! stop accounting time for the data block deallocation +! + call stop_timer(imq) +#endif /* PROFILE */ !------------------------------------------------------------------------------- ! - end subroutine append_metablock -! -!!============================================================================== -!! -!! TOOL SUBROUTINES -!! -! -!=============================================================================== -! -! increase_id: function increases the last identification by 1 and returns its -! value -! -!=============================================================================== -! - function increase_id() - - implicit none - -! return variable -! - integer(kind=4) :: increase_id -! -!------------------------------------------------------------------------------- -! -! increase ID by 1 -! - last_id = last_id + 1 - -! return ID -! - increase_id = last_id - - return - -!------------------------------------------------------------------------------- -! - end function increase_id -! -!=============================================================================== -! -! set_last_id: subroutine sets the last identification value -! -!=============================================================================== -! - subroutine set_last_id(id) - - use error, only : print_error - - implicit none - -! input argument -! - integer(kind=4), intent(in) :: id -! -!------------------------------------------------------------------------------- -! - if (last_id .gt. id) then - call print_error("blocks::set_last_id" & - , "New last_id must be larger than old one!") - else - last_id = id - end if - -!------------------------------------------------------------------------------- -! - end subroutine set_last_id -! -!=============================================================================== -! -! get_last_id: function returns the last identification value -! -!=============================================================================== -! - function get_last_id() - - implicit none - -! return variable -! - integer(kind=4) :: get_last_id -! -!------------------------------------------------------------------------------- -! - get_last_id = last_id - - return - -!------------------------------------------------------------------------------- -! - end function get_last_id -! -!=============================================================================== -! -! get_mblocks: function returns the number of meta blocks -! -!=============================================================================== -! - function get_mblocks() - - implicit none - -! return variable -! - integer(kind=4) :: get_mblocks -! -!------------------------------------------------------------------------------- -! - get_mblocks = mblocks - - return - -!------------------------------------------------------------------------------- -! - end function get_mblocks -! -!=============================================================================== -! -! get_dblocks: function returns the number of data blocks -! -!=============================================================================== -! - function get_dblocks() - - implicit none - -! return variable -! - integer(kind=4) :: get_dblocks -! -!------------------------------------------------------------------------------- -! - get_dblocks = dblocks - - return - -!------------------------------------------------------------------------------- -! - end function get_dblocks -! -!=============================================================================== -! -! get_nleafs: function returns the number of leafs -! -!=============================================================================== -! - function get_nleafs() - - implicit none - -! return variable -! - integer(kind=4) :: get_nleafs -! -!------------------------------------------------------------------------------- -! - get_nleafs = nleafs - - return - -!------------------------------------------------------------------------------- -! - end function get_nleafs + end subroutine deallocate_datablock ! !=============================================================================== ! ! subroutine LINK_BLOCKS: ! ---------------------- ! -! Subroutine links a data block to meta block. +! Subroutine links meta and data blocks. ! ! Arguments: ! @@ -731,7 +1164,7 @@ module blocks ! !------------------------------------------------------------------------------- ! -! set the pointers +! associate the corresponging pointers ! pmeta%data => pdata pdata%meta => pmeta @@ -745,7 +1178,7 @@ module blocks ! subroutine UNLINK_BLOCKS: ! ------------------------ ! -! Subroutine unlinks meta block with a data block linked to it. +! Subroutine unlinks meta and data blocks. ! ! Arguments: ! @@ -767,7 +1200,7 @@ module blocks ! !------------------------------------------------------------------------------- ! -! clean up the block pointers +! nullify the corresponging pointers ! nullify(pmeta%data) nullify(pdata%meta) @@ -778,322 +1211,25 @@ module blocks ! !=============================================================================== ! -! metablock_set_id: subroutine sets the identification value +! subroutine REFINE_BLOCK: +! ----------------------- ! -!=============================================================================== -! - subroutine metablock_set_id(pmeta, id) - - implicit none - -! input arguments -! - type(block_meta), pointer, intent(inout) :: pmeta - integer(kind=4) , intent(in) :: id -! -!------------------------------------------------------------------------------- -! -! set the id field -! - pmeta%id = id - -! check if the id is larger then last_id, if so reset last_id to id -! - if (last_id .lt. id) last_id = id - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_id -! -!=============================================================================== -! -! metablock_set_cpu: subroutine sets the cpu number -! -!=============================================================================== -! - subroutine metablock_set_cpu(pmeta, cpu) - - implicit none - -! input arguments -! - type(block_meta), pointer, intent(inout) :: pmeta - integer(kind=4) , intent(in) :: cpu -! -!------------------------------------------------------------------------------- -! -! set the cpu field -! - pmeta%cpu = cpu - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_cpu -! -!=============================================================================== -! -! metablock_set_refine: subroutine sets the refine flag -! -!=============================================================================== -! - subroutine metablock_set_refine(pmeta, refine) - - use error, only : print_error - - implicit none - -! input arguments -! - type(block_meta), pointer, intent(inout) :: pmeta - integer(kind=4) , intent(in) :: refine -! -!------------------------------------------------------------------------------- -! -! check if the refine value is correct -! - if (abs(refine) .gt. 1) then - -! print error about wrong refine flag -! - call print_error("blocks::metablock_set_refine" & - , "New refine flag is incorrect!") - - else - -! set the refine field -! - pmeta%refine = refine - - end if - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_refine -! -!=============================================================================== -! -! metablock_set_leaf: subroutine marks the block as a leaf -! -!=============================================================================== -! - subroutine metablock_set_leaf(pmeta) - - implicit none - -! input/output arguments -! - type(block_meta), pointer, intent(inout) :: pmeta -! -!------------------------------------------------------------------------------- -! -! set the leaf flag -! - pmeta%leaf = .true. - -! increase the number of leafs -! - nleafs = nleafs + 1 - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_leaf -! -!=============================================================================== -! -! metablock_unset_leaf: subroutine unmarks the block as a leaf -! -!=============================================================================== -! - subroutine metablock_unset_leaf(pmeta) - - implicit none - -! input/output arguments -! - type(block_meta), pointer, intent(inout) :: pmeta -! -!------------------------------------------------------------------------------- -! -! set the leaf flag -! - pmeta%leaf = .false. - -! decrease the number of leafs -! - nleafs = nleafs - 1 - -!------------------------------------------------------------------------------- -! - end subroutine metablock_unset_leaf -! -!=============================================================================== -! -! metablock_set_config: subroutine sets the configuration flag -! -!=============================================================================== -! - subroutine metablock_set_config(pmeta, config) - - implicit none - -! input/output arguments -! - type(block_meta), pointer, intent(inout) :: pmeta - integer(kind=4) , intent(in) :: config -! -!------------------------------------------------------------------------------- -! -! set the config flag -! - pmeta%config = config - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_config -! -!=============================================================================== -! -! metablock_set_level: subroutine sets the level of data block -! -!=============================================================================== -! - subroutine metablock_set_level(pmeta, level) - - implicit none - -! input/output arguments -! - type(block_meta), pointer, intent(inout) :: pmeta - integer(kind=4) , intent(in) :: level -! -!------------------------------------------------------------------------------- -! -! set the refinement level -! - pmeta%level = level - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_level -! -!=============================================================================== -! -! metablock_set_position: subroutine sets the position of the meta block in the -! parent block -! -!=============================================================================== -! - subroutine metablock_set_position(pmeta, px, py, pz) - - implicit none - -! input/output arguments -! - type(block_meta), pointer, intent(inout) :: pmeta - integer(kind=4) , intent(in) :: px, py, pz -! -!------------------------------------------------------------------------------- -! -! set the position in the parent block -! - pmeta%pos(1) = px - pmeta%pos(2) = py -#if NDIMS == 3 - pmeta%pos(3) = pz -#endif /* NDIMS == 3 */ - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_position -! -!=============================================================================== -! -! metablock_set_coord: subroutine sets the coordinates of the meta block -! -!=============================================================================== -! - subroutine metablock_set_coord(pmeta, px, py, pz) - - implicit none - -! input/output arguments -! - type(block_meta), pointer, intent(inout) :: pmeta - integer(kind=4) , intent(in) :: px, py, pz -! -!------------------------------------------------------------------------------- -! -! set the coordinates -! - pmeta%coord(1) = px - pmeta%coord(2) = py -#if NDIMS == 3 - pmeta%coord(3) = pz -#endif /* NDIMS == 3 */ - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_coord -! -!=============================================================================== -! -! metablock_set_bounds: subroutine sets the bounds of data block -! -!=============================================================================== -! - subroutine metablock_set_bounds(pmeta, xmin, xmax, ymin, ymax, zmin, zmax) - - implicit none - -! input/output arguments -! - type(block_meta), pointer, intent(inout) :: pmeta - real , intent(in) :: xmin, xmax - real , intent(in) :: ymin, ymax - real , intent(in) :: zmin, zmax -! -!------------------------------------------------------------------------------- -! -! set bounds of the block -! - pmeta%xmin = xmin - pmeta%xmax = xmax - pmeta%ymin = ymin - pmeta%ymax = ymax - pmeta%zmin = zmin - pmeta%zmax = zmax - -!------------------------------------------------------------------------------- -! - end subroutine metablock_set_bounds -! -!!============================================================================== -!! -!! DATA BLOCK SUBROUTINES -!! -! -!=============================================================================== -! -! allocate_datablock: subroutine allocates space for one data block and returns -! the pointer to this block -! -!=============================================================================== -! -!=============================================================================== -! -! subroutine ALLOCATE_DATABLOCK: -! ----------------------------- -! -! Subroutine allocates space for one data block and returns a pointer -! associated with it. +! Subroutine creates children of the current block and initializes their +! configuration, pointers and fields. ! ! Arguments: ! -! pdata - the pointer associated with the created data block; +! pmeta - a pointer to meta block for which children will be created; +! res - the resolution of the block; +! fdata - a flag indicating if data blocks for children should be allocated; ! !=============================================================================== ! - subroutine allocate_datablock(pdata) + subroutine refine_block(pmeta, res, fdata) + +! import external procedures +! + use error , only : print_error ! local variables are not implicit by default ! @@ -1101,526 +1237,46 @@ module blocks ! subroutine arguments ! - type(block_data), pointer, intent(out) :: pdata -! -!------------------------------------------------------------------------------- -! -! allocate the block structure -! - allocate(pdata) - -! nullify all pointers -! - nullify(pdata%prev) - nullify(pdata%next) - nullify(pdata%meta) - -! allocate the space for conserved variables -! - allocate(pdata%u0(nvars,nx,ny,nz)) - allocate(pdata%u1(nvars,nx,ny,nz)) - -! allocate the space for primitive variables -! - allocate(pdata%q(nvars,nx,ny,nz)) - -! initiate the conserved variable pointer -! - pdata%u => pdata%u0 - -! allocate the space for numerical fluxes -! - if (nflux > 0) allocate(pdata%f(ndims,nflux,nx,ny,nz)) - -#ifdef DEBUG -! allocate the space for the refinement criterion array -! - allocate(pdata%c(nx,ny,nz)) -#endif /* DEBUG */ - -! increase the number of allocated meta blocks -! - dblocks = dblocks + 1 - -!------------------------------------------------------------------------------- -! - end subroutine allocate_datablock -! -!=============================================================================== -! -! subroutine DEALLOCATE_DATABLOCK: -! ------------------------------- -! -! Subroutine deallocates space of the data block associated with the input -! pointer. -! -! Arguments: -! -! pdata - the pointer pointing to the data block for deallocating; -! -!=============================================================================== -! - subroutine deallocate_datablock(pdata) - -! local variables are not implicit by default -! - implicit none - -! subroutine arguments -! - type(block_data), pointer, intent(inout) :: pdata -! -!------------------------------------------------------------------------------- -! -! check if the input pointer is associated with a data block -! - if (associated(pdata)) then - -! deallocate conservative variables -! - if (allocated(pdata%u0)) deallocate(pdata%u0) - if (allocated(pdata%u1)) deallocate(pdata%u1) - -! deallocate primitive variables -! - if (allocated(pdata%q)) deallocate(pdata%q) - -! deallocate numerical fluxes -! - if (allocated(pdata%f)) deallocate(pdata%f) - -#ifdef DEBUG -! deallocate the refinement critarion array -! - if (allocated(pdata%c)) deallocate(pdata%c) -#endif /* DEBUG */ - -! nullify pointers -! - nullify(pdata%u) - nullify(pdata%next) - nullify(pdata%prev) - nullify(pdata%meta) - -! free and nullify the block -! - deallocate(pdata) - nullify(pdata) - -! decrease the number of allocated blocks -! - dblocks = dblocks - 1 - - end if ! pdata associated with a data block - -!------------------------------------------------------------------------------- -! - end subroutine deallocate_datablock -! -!=============================================================================== -! -! subroutine APPEND_DATABLOCK: -! --------------------------- -! -! Subroutine allocates space for one data block and appends it to the data -! block list returning a pointer associated with it. -! -! Arguments: -! -! pdata - the pointer associated with the created data block; -! -!=============================================================================== -! - subroutine append_datablock(pdata) - -! local variables are not implicit by default -! - implicit none - -! subroutine arguments -! - type(block_data), pointer, intent(out) :: pdata -! -!------------------------------------------------------------------------------- -! -! allocate the data block -! - call allocate_datablock(pdata) - -! add the allocated block to the data block list -! - if (associated(last_data)) then - pdata%prev => last_data - last_data%next => pdata - else - list_data => pdata - end if - -! set the pointer to the last block in the list -! - last_data => pdata - -!------------------------------------------------------------------------------- -! - end subroutine append_datablock -! -!=============================================================================== -! -! subroutine REMOVE_DATABLOCK: -! --------------------------- -! -! Subroutine removes a data block associated with the input pointer from -! the data block list, and deallocates space used by this block. -! -! Arguments: -! -! pdata - the pointer pointing to the data block for removing; -! -!=============================================================================== -! - subroutine remove_datablock(pdata) - -! local variables are not implicit by default -! - implicit none - -! subroutine arguments -! - type(block_data), pointer, intent(inout) :: pdata -! -!------------------------------------------------------------------------------- -! -! check if the input pointer is associated with a data block -! - if (associated(pdata)) then - -! remove from the meta block list if the meta pointer is set -! - if (associated(pdata%meta)) then - -! if this is the first block in the list, update the list_data pointer -! - if (pdata%meta%id == list_data%meta%id) list_data => pdata%next - -! if this is the last block in the list, update the last_data pointer -! - if (pdata%meta%id == last_data%meta%id) last_data => pdata%prev - -! update the pointer of previous and next blocks -! - if (associated(pdata%prev)) pdata%prev%next => pdata%next - - if (associated(pdata%next)) pdata%next%prev => pdata%prev - - end if ! %meta associated - -! deallocate the associated data block -! - call deallocate_datablock(pdata) - - end if ! pdata associated with a data block - -!------------------------------------------------------------------------------- -! - end subroutine remove_datablock -! -!=============================================================================== -! -! datablock_set_dims: subroutine sets the number of variables and dimensions -! for arrays allocated in data blocks -! -!=============================================================================== -! - subroutine datablock_set_dims(nv, nf, ni, nj, nk) - - implicit none - -! input arguments -! - integer(kind=4), intent(in) :: nv, nf, ni, nj, nk -! -!------------------------------------------------------------------------------- -! - nvars = nv - nflux = nf - nx = ni - ny = nj -#if NDIMS == 3 - nz = nk -#endif /* NDIMS == 3 */ - -!------------------------------------------------------------------------------- -! - end subroutine datablock_set_dims -! -!!============================================================================== -!! -!! REFINEMENT/DEREFINEMENT SUBROUTINES -!! -! -!=============================================================================== -! -! refine_block: subroutine refines selected block -! -!=============================================================================== -! - subroutine refine_block(pblock, res, falloc_data) - - use error, only : print_error - - implicit none - -! input parameters -! - type(block_meta), pointer , intent(inout) :: pblock + type(block_meta), pointer , intent(inout) :: pmeta integer(kind=4), dimension(3), intent(in) :: res - logical , intent(in) :: falloc_data + logical , intent(in) :: fdata ! pointers ! - type(block_meta), pointer :: pneigh, pchild, pfirst, plast + type(block_meta), pointer :: pnext, pneigh, pchild type(block_data), pointer :: pdata +! local variables +! + integer :: p, q, i, j, k, ic, jc, kc + real :: xln, yln, zln, xmn, xmx, ymn, ymx, zmn, zmx + ! local arrays ! integer, dimension(nchildren) :: config, order integer, dimension(ndims,nsides,nfaces) :: set - -! local variables -! - integer :: p, i, j, k, ic, jc, kc - real :: xln, yln, zln, xmn, xmx, ymn, ymx, zmn, zmx ! !------------------------------------------------------------------------------- ! +#ifdef PROFILE +! start accounting time for the block refinement +! + call start_timer(imr) +#endif /* PROFILE */ + ! check if pointer is associated ! - if (associated(pblock)) then + if (associated(pmeta)) then -! unset block leaf flag +! store the pointer to the next block on the list ! - call metablock_unset_leaf(pblock) - -! reset refinement flag -! - pblock%refine = 0 - -! iterate over all child blocks -! - do p = 1, nchildren - -! create child meta and data blocks -! - call allocate_metablock(pblock%child(p)%ptr) - -! set it as a leaf -! - call metablock_set_leaf(pblock%child(p)%ptr) - -! assign pointer to the parent block -! - pblock%child(p)%ptr%parent => pblock - -! increase the refinement level -! - pblock%child(p)%ptr%level = pblock%level + 1 - -! copy the parent cpu number to each child -! - pblock%child(p)%ptr%cpu = pblock%cpu - - end do - -! assign neighbors of the child blocks -! -! interior of the block -! - do p = 1, nfaces - -! X direction (left side) -! - pblock%child(2)%ptr%neigh(1,1,p)%ptr => pblock%child(1)%ptr - pblock%child(4)%ptr%neigh(1,1,p)%ptr => pblock%child(3)%ptr -#if NDIMS == 3 - pblock%child(6)%ptr%neigh(1,1,p)%ptr => pblock%child(5)%ptr - pblock%child(8)%ptr%neigh(1,1,p)%ptr => pblock%child(7)%ptr -#endif /* NDIMS == 3 */ - pneigh => pblock%neigh(1,1,1)%ptr - if (associated(pneigh)) then - if (pneigh%id .eq. pblock%id) then - pblock%child(1)%ptr%neigh(1,1,p)%ptr => pblock%child(2)%ptr - pblock%child(3)%ptr%neigh(1,1,p)%ptr => pblock%child(4)%ptr -#if NDIMS == 3 - pblock%child(5)%ptr%neigh(1,1,p)%ptr => pblock%child(6)%ptr - pblock%child(7)%ptr%neigh(1,1,p)%ptr => pblock%child(8)%ptr -#endif /* NDIMS == 3 */ - end if - end if - -! X direction (right side) -! - pblock%child(1)%ptr%neigh(1,2,p)%ptr => pblock%child(2)%ptr - pblock%child(3)%ptr%neigh(1,2,p)%ptr => pblock%child(4)%ptr -#if NDIMS == 3 - pblock%child(5)%ptr%neigh(1,2,p)%ptr => pblock%child(6)%ptr - pblock%child(7)%ptr%neigh(1,2,p)%ptr => pblock%child(8)%ptr -#endif /* NDIMS == 3 */ - pneigh => pblock%neigh(1,2,1)%ptr - if (associated(pneigh)) then - if (pneigh%id .eq. pblock%id) then - pblock%child(2)%ptr%neigh(1,2,p)%ptr => pblock%child(1)%ptr - pblock%child(4)%ptr%neigh(1,2,p)%ptr => pblock%child(3)%ptr -#if NDIMS == 3 - pblock%child(6)%ptr%neigh(1,2,p)%ptr => pblock%child(5)%ptr - pblock%child(8)%ptr%neigh(1,2,p)%ptr => pblock%child(7)%ptr -#endif /* NDIMS == 3 */ - end if - end if - -! Y direction (left side) -! - pblock%child(3)%ptr%neigh(2,1,p)%ptr => pblock%child(1)%ptr - pblock%child(4)%ptr%neigh(2,1,p)%ptr => pblock%child(2)%ptr -#if NDIMS == 3 - pblock%child(7)%ptr%neigh(2,1,p)%ptr => pblock%child(5)%ptr - pblock%child(8)%ptr%neigh(2,1,p)%ptr => pblock%child(6)%ptr -#endif /* NDIMS == 3 */ - pneigh => pblock%neigh(2,1,1)%ptr - if (associated(pneigh)) then - if (pneigh%id .eq. pblock%id) then - pblock%child(1)%ptr%neigh(2,1,p)%ptr => pblock%child(3)%ptr - pblock%child(2)%ptr%neigh(2,1,p)%ptr => pblock%child(4)%ptr -#if NDIMS == 3 - pblock%child(5)%ptr%neigh(2,1,p)%ptr => pblock%child(7)%ptr - pblock%child(6)%ptr%neigh(2,1,p)%ptr => pblock%child(8)%ptr -#endif /* NDIMS == 3 */ - end if - end if - -! Y direction (right side) -! - pblock%child(1)%ptr%neigh(2,2,p)%ptr => pblock%child(3)%ptr - pblock%child(2)%ptr%neigh(2,2,p)%ptr => pblock%child(4)%ptr -#if NDIMS == 3 - pblock%child(5)%ptr%neigh(2,2,p)%ptr => pblock%child(7)%ptr - pblock%child(6)%ptr%neigh(2,2,p)%ptr => pblock%child(8)%ptr -#endif /* NDIMS == 3 */ - pneigh => pblock%neigh(2,2,1)%ptr - if (associated(pneigh)) then - if (pneigh%id .eq. pblock%id) then - pblock%child(3)%ptr%neigh(2,2,p)%ptr => pblock%child(1)%ptr - pblock%child(4)%ptr%neigh(2,2,p)%ptr => pblock%child(2)%ptr -#if NDIMS == 3 - pblock%child(7)%ptr%neigh(2,2,p)%ptr => pblock%child(5)%ptr - pblock%child(8)%ptr%neigh(2,2,p)%ptr => pblock%child(6)%ptr -#endif /* NDIMS == 3 */ - end if - end if - -#if NDIMS == 3 -! Z direction (left side) -! - pblock%child(5)%ptr%neigh(3,1,p)%ptr => pblock%child(1)%ptr - pblock%child(6)%ptr%neigh(3,1,p)%ptr => pblock%child(2)%ptr - pblock%child(7)%ptr%neigh(3,1,p)%ptr => pblock%child(3)%ptr - pblock%child(8)%ptr%neigh(3,1,p)%ptr => pblock%child(4)%ptr - pneigh => pblock%neigh(3,1,1)%ptr - if (associated(pneigh)) then - if (pneigh%id .eq. pblock%id) then - pblock%child(1)%ptr%neigh(3,1,p)%ptr => pblock%child(5)%ptr - pblock%child(2)%ptr%neigh(3,1,p)%ptr => pblock%child(6)%ptr - pblock%child(3)%ptr%neigh(3,1,p)%ptr => pblock%child(7)%ptr - pblock%child(4)%ptr%neigh(3,1,p)%ptr => pblock%child(8)%ptr - end if - end if - -! Z direction (right side) -! - pblock%child(1)%ptr%neigh(3,2,p)%ptr => pblock%child(5)%ptr - pblock%child(2)%ptr%neigh(3,2,p)%ptr => pblock%child(6)%ptr - pblock%child(3)%ptr%neigh(3,2,p)%ptr => pblock%child(7)%ptr - pblock%child(4)%ptr%neigh(3,2,p)%ptr => pblock%child(8)%ptr - pneigh => pblock%neigh(3,2,1)%ptr - if (associated(pneigh)) then - if (pneigh%id .eq. pblock%id) then - pblock%child(5)%ptr%neigh(3,2,p)%ptr => pblock%child(1)%ptr - pblock%child(6)%ptr%neigh(3,2,p)%ptr => pblock%child(2)%ptr - pblock%child(7)%ptr%neigh(3,2,p)%ptr => pblock%child(3)%ptr - pblock%child(8)%ptr%neigh(3,2,p)%ptr => pblock%child(4)%ptr - end if - end if -#endif /* NDIMS == 3 */ - end do - -! prepare set array -! -#if NDIMS == 2 - set(1,1,:) = (/ 1, 3 /) - set(1,2,:) = (/ 2, 4 /) - set(2,1,:) = (/ 1, 2 /) - set(2,2,:) = (/ 3, 4 /) -#endif /* NDIMS == 2 */ -#if NDIMS == 3 - set(1,1,:) = (/ 1, 3, 5, 7 /) - set(1,2,:) = (/ 2, 4, 6, 8 /) - set(2,1,:) = (/ 1, 2, 5, 6 /) - set(2,2,:) = (/ 3, 4, 7, 8 /) - set(3,1,:) = (/ 1, 2, 3, 4 /) - set(3,2,:) = (/ 5, 6, 7, 8 /) -#endif /* NDIMS == 3 */ - -! set pointers to neighbors and update neighbors pointers -! - do i = 1, ndims - do j = 1, nsides - do k = 1, nfaces - pneigh => pblock%neigh(i,j,k)%ptr - pchild => pblock%child(set(i,j,k))%ptr - - if (associated(pneigh)) then - if (pneigh%id .ne. pblock%id) then - -! point to the right neighbor -! - do p = 1, nfaces - pchild%neigh(i,j,p)%ptr => pneigh - end do - -! neighbor level is the same as the refined block -! - if (pneigh%level .eq. pblock%level) then - pneigh%neigh(i,3-j,k)%ptr => pchild - end if - -! neighbor level is the same as the child block -! - if (pneigh%level .gt. pblock%level) then - do p = 1, nfaces - pneigh%neigh(i,3-j,p)%ptr => pchild - end do - end if - - end if - - end if - - end do - end do - end do - -! reset neighbor pointers of the parent block -! - do i = 1, ndims - do j = 1, nsides - do k = 1, nfaces - nullify(pblock%neigh(i,j,k)%ptr) - end do - end do - end do + pnext => pmeta%next +!! PREPARE CHILD CONFIGURATION PARAMETERS +!! ! set corresponding configuration of the new blocks ! - select case(pblock%config) + select case(pmeta%conf) case(0) #if NDIMS == 2 @@ -1720,366 +1376,1972 @@ module blocks end select -! set blocks configurations +! calculate sizes of the child blocks ! - do p = 1, nchildren - pblock%child(order(p))%ptr%config = config(p) - end do - -! connect blocks in chain -! - do p = 2, nchildren - pblock%child(order(p ))%ptr%prev => pblock%child(order(p-1))%ptr - pblock%child(order(p-1))%ptr%next => pblock%child(order(p ))%ptr - end do - -! insert this chain after the parent block -! - pneigh => pblock%next - pfirst => pblock%child(order( 1))%ptr - plast => pblock%child(order(nchildren))%ptr - if (associated(pneigh)) then - pneigh%prev => plast - plast%next => pneigh - else - last_meta => plast - end if - - pblock%next => pfirst - pfirst%prev => pblock - -! calculate the size of new blocks -! - xln = 0.5 * (pblock%xmax - pblock%xmin) - yln = 0.5 * (pblock%ymax - pblock%ymin) + xln = 0.5d+00 * (pmeta%xmax - pmeta%xmin) + yln = 0.5d+00 * (pmeta%ymax - pmeta%ymin) #if NDIMS == 3 - zln = 0.5 * (pblock%zmax - pblock%zmin) + zln = 0.5d+00 * (pmeta%zmax - pmeta%zmin) #else /* NDIMS == 3 */ - zln = (pblock%zmax - pblock%zmin) + zln = (pmeta%zmax - pmeta%zmin) #endif /* NDIMS == 3 */ -! iterate over all children and allocate data blocks +!! ALLOCATE CHILDREN AND APPEND THEM TO THE META LIST +!! +! iterate over the number of children in the reverse configuration order, i.e. +! the allocated blocks are inserted after the parent block following +! the reversed Hilbert curve +! + do p = nchildren, 1, -1 + +! insert a new meta block after pmeta and associate it with pchild +! + call insert_metablock_after(pmeta, pchild) + +! set the child configuration number +! + call metablock_set_configuration(pchild, config(p)) + +! associate the parent's children array element with the freshly created +! meta block +! + pmeta%child(order(p))%ptr => pchild + + end do ! nchildren + +! iterate over all children ! do p = 1, nchildren -! assign a pointer to the current child +! associate a pointer with the current child ! - pchild => pblock%child(p)%ptr + pchild => pmeta%child(p)%ptr + +! associate the parent field with pmeta +! + pchild%parent => pmeta + +! mark the child as the leaf +! + call metablock_set_leaf(pchild) + +! mark the child to be updated +! + call metablock_set_update(pchild) + +! set the child refinement level +! + call metablock_set_level(pchild, pmeta%level + 1) + +! set the child process number +! + call metablock_set_process(pchild, pmeta%process) ! calculate the block position indices ! - i = mod((p - 1) ,2) - j = mod((p - 1) / 2,2) - k = mod((p - 1) / 4,2) + q = p - 1 + i = mod(q ,2) + j = mod(q / 2,2) + k = mod(q / 4,2) + +! calculate the block coordinates in effective resolution units +! + ic = pmeta%coords(1) + i * res(1) + jc = pmeta%coords(2) + j * res(2) +#if NDIMS == 3 + kc = pmeta%coords(3) + k * res(3) +#endif /* NDIMS == 3 */ + +! calculate block bounds +! + xmn = pmeta%xmin + xln * i + ymn = pmeta%ymin + yln * j + zmn = pmeta%zmin + zln * k + + xmx = xmn + xln + ymx = ymn + yln + zmx = zmn + zln ! set the block position ! call metablock_set_position(pchild, i, j, k) -! set the block coordinates +! set the effective resolution coordinates ! - ic = pblock%coord(1) + i * res(1) - jc = pblock%coord(2) + j * res(2) -#if NDIMS == 3 - kc = pblock%coord(3) + k * res(3) -#endif /* NDIMS == 3 */ - call metablock_set_coord(pchild, ic, jc, kc) + call metablock_set_coordinates(pchild, ic, jc, kc) -! calculate block bounds -! - xmn = pblock%xmin + xln * i - ymn = pblock%ymin + yln * j - zmn = pblock%zmin + zln * k - - xmx = xmn + xln - ymx = ymn + yln - zmx = zmn + zln - -! set block bounds +! set the child block bounds ! call metablock_set_bounds(pchild, xmn, xmx, ymn, ymx, zmn, zmx) - end do + end do ! nchildren -! allocate data blocks if necessary +!! ASSIGN PROPER NEIGHBORS FOR THE CHILDREN IN THE INTERIOR OF THE PARENT BLOCK +!! +! iterate over faces and update the interior of the block ! - if (falloc_data) then + do p = 1, nfaces -! iterate over all children and allocate data blocks +! X direction (left side) +! + pmeta%child(2)%ptr%neigh(1,1,p)%ptr => pmeta%child(1)%ptr + pmeta%child(4)%ptr%neigh(1,1,p)%ptr => pmeta%child(3)%ptr +#if NDIMS == 3 + pmeta%child(6)%ptr%neigh(1,1,p)%ptr => pmeta%child(5)%ptr + pmeta%child(8)%ptr%neigh(1,1,p)%ptr => pmeta%child(7)%ptr +#endif /* NDIMS == 3 */ + +! associate pneigh with a neighbor +! + pneigh => pmeta%neigh(1,1,1)%ptr + +! if neighbor and associated and points to parent block, it corresponds to +! periodic boundaries at the lowest level +! + if (associated(pneigh)) then + + if (pneigh%id == pmeta%id) then + + pmeta%child(1)%ptr%neigh(1,1,p)%ptr => pmeta%child(2)%ptr + pmeta%child(3)%ptr%neigh(1,1,p)%ptr => pmeta%child(4)%ptr +#if NDIMS == 3 + pmeta%child(5)%ptr%neigh(1,1,p)%ptr => pmeta%child(6)%ptr + pmeta%child(7)%ptr%neigh(1,1,p)%ptr => pmeta%child(8)%ptr +#endif /* NDIMS == 3 */ + + end if + + end if + +! X direction (right side) +! + pmeta%child(1)%ptr%neigh(1,2,p)%ptr => pmeta%child(2)%ptr + pmeta%child(3)%ptr%neigh(1,2,p)%ptr => pmeta%child(4)%ptr +#if NDIMS == 3 + pmeta%child(5)%ptr%neigh(1,2,p)%ptr => pmeta%child(6)%ptr + pmeta%child(7)%ptr%neigh(1,2,p)%ptr => pmeta%child(8)%ptr +#endif /* NDIMS == 3 */ + +! associate pneigh with a neighbor +! + pneigh => pmeta%neigh(1,2,1)%ptr + +! if neighbor and associated and points to parent block, it corresponds to +! periodic boundaries at the lowest level +! + if (associated(pneigh)) then + + if (pneigh%id == pmeta%id) then + + pmeta%child(2)%ptr%neigh(1,2,p)%ptr => pmeta%child(1)%ptr + pmeta%child(4)%ptr%neigh(1,2,p)%ptr => pmeta%child(3)%ptr +#if NDIMS == 3 + pmeta%child(6)%ptr%neigh(1,2,p)%ptr => pmeta%child(5)%ptr + pmeta%child(8)%ptr%neigh(1,2,p)%ptr => pmeta%child(7)%ptr +#endif /* NDIMS == 3 */ + + end if + + end if + +! Y direction (left side) +! + pmeta%child(3)%ptr%neigh(2,1,p)%ptr => pmeta%child(1)%ptr + pmeta%child(4)%ptr%neigh(2,1,p)%ptr => pmeta%child(2)%ptr +#if NDIMS == 3 + pmeta%child(7)%ptr%neigh(2,1,p)%ptr => pmeta%child(5)%ptr + pmeta%child(8)%ptr%neigh(2,1,p)%ptr => pmeta%child(6)%ptr +#endif /* NDIMS == 3 */ + +! associate pneigh with a neighbor +! + pneigh => pmeta%neigh(2,1,1)%ptr + +! if neighbor and associated and points to parent block, it corresponds to +! periodic boundaries at the lowest level +! + if (associated(pneigh)) then + + if (pneigh%id == pmeta%id) then + + pmeta%child(1)%ptr%neigh(2,1,p)%ptr => pmeta%child(3)%ptr + pmeta%child(2)%ptr%neigh(2,1,p)%ptr => pmeta%child(4)%ptr +#if NDIMS == 3 + pmeta%child(5)%ptr%neigh(2,1,p)%ptr => pmeta%child(7)%ptr + pmeta%child(6)%ptr%neigh(2,1,p)%ptr => pmeta%child(8)%ptr +#endif /* NDIMS == 3 */ + + end if + + end if + +! Y direction (right side) +! + pmeta%child(1)%ptr%neigh(2,2,p)%ptr => pmeta%child(3)%ptr + pmeta%child(2)%ptr%neigh(2,2,p)%ptr => pmeta%child(4)%ptr +#if NDIMS == 3 + pmeta%child(5)%ptr%neigh(2,2,p)%ptr => pmeta%child(7)%ptr + pmeta%child(6)%ptr%neigh(2,2,p)%ptr => pmeta%child(8)%ptr +#endif /* NDIMS == 3 */ + +! associate pneigh with a neighbor +! + pneigh => pmeta%neigh(2,2,1)%ptr + +! if neighbor and associated and points to parent block, it corresponds to +! periodic boundaries at the lowest level +! + if (associated(pneigh)) then + + if (pneigh%id == pmeta%id) then + + pmeta%child(3)%ptr%neigh(2,2,p)%ptr => pmeta%child(1)%ptr + pmeta%child(4)%ptr%neigh(2,2,p)%ptr => pmeta%child(2)%ptr +#if NDIMS == 3 + pmeta%child(7)%ptr%neigh(2,2,p)%ptr => pmeta%child(5)%ptr + pmeta%child(8)%ptr%neigh(2,2,p)%ptr => pmeta%child(6)%ptr +#endif /* NDIMS == 3 */ + + end if + + end if + +#if NDIMS == 3 +! Z direction (left side) +! + pmeta%child(5)%ptr%neigh(3,1,p)%ptr => pmeta%child(1)%ptr + pmeta%child(6)%ptr%neigh(3,1,p)%ptr => pmeta%child(2)%ptr + pmeta%child(7)%ptr%neigh(3,1,p)%ptr => pmeta%child(3)%ptr + pmeta%child(8)%ptr%neigh(3,1,p)%ptr => pmeta%child(4)%ptr + +! associate pneigh with a neighbor +! + pneigh => pmeta%neigh(3,1,1)%ptr + +! if neighbor and associated and points to parent block, it corresponds to +! periodic boundaries at the lowest level +! + if (associated(pneigh)) then + + if (pneigh%id == pmeta%id) then + + pmeta%child(1)%ptr%neigh(3,1,p)%ptr => pmeta%child(5)%ptr + pmeta%child(2)%ptr%neigh(3,1,p)%ptr => pmeta%child(6)%ptr + pmeta%child(3)%ptr%neigh(3,1,p)%ptr => pmeta%child(7)%ptr + pmeta%child(4)%ptr%neigh(3,1,p)%ptr => pmeta%child(8)%ptr + + end if + + end if + +! Z direction (right side) +! + pmeta%child(1)%ptr%neigh(3,2,p)%ptr => pmeta%child(5)%ptr + pmeta%child(2)%ptr%neigh(3,2,p)%ptr => pmeta%child(6)%ptr + pmeta%child(3)%ptr%neigh(3,2,p)%ptr => pmeta%child(7)%ptr + pmeta%child(4)%ptr%neigh(3,2,p)%ptr => pmeta%child(8)%ptr + +! associate pneigh with a neighbor +! + pneigh => pmeta%neigh(3,2,1)%ptr + +! if neighbor and associated and points to parent block, it corresponds to +! periodic boundaries at the lowest level +! + if (associated(pneigh)) then + + if (pneigh%id == pmeta%id) then + + pmeta%child(5)%ptr%neigh(3,2,p)%ptr => pmeta%child(1)%ptr + pmeta%child(6)%ptr%neigh(3,2,p)%ptr => pmeta%child(2)%ptr + pmeta%child(7)%ptr%neigh(3,2,p)%ptr => pmeta%child(3)%ptr + pmeta%child(8)%ptr%neigh(3,2,p)%ptr => pmeta%child(4)%ptr + + end if + + end if +#endif /* NDIMS == 3 */ + + end do ! nfaces + +!! UPDATE NEIGHBORS AND EXTERNAL NEIGHBORS OF CHILDREN +!! +! prepare set array +! +#if NDIMS == 2 + set(1,1,:) = (/ 1, 3 /) + set(1,2,:) = (/ 2, 4 /) + set(2,1,:) = (/ 1, 2 /) + set(2,2,:) = (/ 3, 4 /) +#endif /* NDIMS == 2 */ +#if NDIMS == 3 + set(1,1,:) = (/ 1, 3, 5, 7 /) + set(1,2,:) = (/ 2, 4, 6, 8 /) + set(2,1,:) = (/ 1, 2, 5, 6 /) + set(2,2,:) = (/ 3, 4, 7, 8 /) + set(3,1,:) = (/ 1, 2, 3, 4 /) + set(3,2,:) = (/ 5, 6, 7, 8 /) +#endif /* NDIMS == 3 */ + +! set pointers to neighbors and update neighbors' pointers +! + do i = 1, ndims + do j = 1, nsides + +! prepare reverse side index +! + q = 3 - j + +! iterate over all faces +! + do k = 1, nfaces + +! associate pointers with the neighbor and child +! + pneigh => pmeta%neigh(i,j,k)%ptr + pchild => pmeta%child(set(i,j,k))%ptr + +! check if neighbor is associated +! + if (associated(pneigh)) then + +! check if the parent block does not point to itself (periodic boundaries) +! + if (pneigh%id /= pmeta%id) then + +! point the child neigh field to the right neighbor +! + do p = 1, nfaces + pchild%neigh(i,j,p)%ptr => pneigh + end do + +! update neighbor pointer if it is at the same level +! + if (pneigh%level == pmeta%level) then + pneigh%neigh(i,q,k)%ptr => pchild + end if + +! update neighbor pointer if it is at higher level +! + if (pneigh%level > pmeta%level) then + do p = 1, nfaces + pneigh%neigh(i,q,p)%ptr => pchild + end do + end if + +! if neighbor has lower level than parent, something is wrong, since lower +! levels should be already refined +! + if (pneigh%level < pmeta%level) then + call print_error("blocks::refine_block" & + , "Neighbor found at lower level!") + end if + + end if ! pmeta and pneigh point to different blocks + + end if ! pneigh is associated + + end do ! nfaces + end do ! nsides + end do ! ndims + +! mark all neighbors to be updated as well +! + call neighbors_set_update(pmeta) + +!! ASSOCIATE DATA BLOCKS IF NECESSARY +!! +! allocate data blocks if requested +! + if (fdata) then + +! iterate over all children ! do p = 1, nchildren ! assign a pointer to the current child ! - pchild => pblock%child(p)%ptr + pchild => pmeta%child(p)%ptr -! allocate data block +! allocate new data block and append it to the data block list ! - call allocate_datablock(pdata) + call append_datablock(pdata) -! associate with the meta block +! associate the new data block with the current child ! call link_blocks(pchild, pdata) - end do + end do ! nchildren -! connect blocks in chain + end if ! allocate data blocks for children + +!! RESET PARENT'S FIELDS +!! +! unset the block leaf flag ! - do p = 2, nchildren - pblock%child(order(p ))%ptr%data%prev => & - pblock%child(order(p-1))%ptr%data - pblock%child(order(p-1))%ptr%data%next => & - pblock%child(order(p ))%ptr%data - end do + call metablock_unset_leaf(pmeta) -! insert this chain after the parent block +! reset the refinement flag ! - pdata => pblock%data%next + call metablock_set_refinement(pmeta, 0) - pfirst => pblock%child(order( 1))%ptr - plast => pblock%child(order(nchildren))%ptr - - if (associated(pdata)) then - pdata%prev => plast%data - plast%data%next => pdata - else - last_data => plast%data - end if - - pblock%data%next => pfirst%data - pfirst%data%prev => pblock%data +! nullify the parent's neighbor pointers +! + do i = 1, ndims + do j = 1, nsides + do k = 1, nfaces + nullify(pmeta%neigh(i,j,k)%ptr) + end do + end do + end do +! restore the pointer to the current block +! + if (associated(pnext)) then + pmeta => pnext%prev + else + pmeta => last_meta end if -! point the current block to the last created one -! - pblock => plast + else ! pmeta is not associated - else - -! terminate program if the pointer passed by argument is not associated +! it's impossible to refine since there is not block associated with +! the argument pointer ! call print_error("blocks::refine_block" & - , "Input pointer is not associated! Terminating!") + , "No block associated with the argument pointer!") end if +#ifdef PROFILE +! stop accounting time for the block refinement +! + call stop_timer(imr) +#endif /* PROFILE */ + !------------------------------------------------------------------------------- ! end subroutine refine_block ! !=============================================================================== ! -! derefine_block: subroutine derefines selected block +! subroutine DEREFINE_BLOCK: +! ------------------------- +! +! Subroutine derefines the current block by distrying all its children and +! restoring the block configuration, pointers and fields. +! +! Arguments: +! +! pmeta - a pointer to derefined meta block; ! !=============================================================================== ! - subroutine derefine_block(pblock) + subroutine derefine_block(pmeta) +! local variables are not implicit by default +! implicit none -! input parameters +! subroutine arguments ! - type(block_meta), pointer, intent(inout) :: pblock - -! local variables -! - integer :: i, j, k, l, p - -! local arrays -! - integer, dimension(ndims, nsides, nfaces) :: arr + type(block_meta), pointer, intent(inout) :: pmeta ! local pointers ! type(block_meta), pointer :: pchild, pneigh + +! local variables +! + integer :: i, j, k, l, p + +! local saved variables +! + logical, save :: first = .true. + +! local arrays +! + integer, dimension(ndims, nsides, nfaces), save :: arr ! !------------------------------------------------------------------------------- ! +#ifdef PROFILE +! start accounting time for the block derefinement +! + call start_timer(imd) +#endif /* PROFILE */ + +! prepare saved variables at the first execution +! + if (first) then + ! prepare reference array ! #if NDIMS == 3 - arr(1,1,:) = (/ 1, 3, 5, 7 /) - arr(1,2,:) = (/ 2, 4, 6, 8 /) - arr(2,1,:) = (/ 1, 2, 5, 6 /) - arr(2,2,:) = (/ 3, 4, 7, 8 /) - arr(3,1,:) = (/ 1, 2, 3, 4 /) - arr(3,2,:) = (/ 5, 6, 7, 8 /) + arr(1,1,:) = (/ 1, 3, 5, 7 /) + arr(1,2,:) = (/ 2, 4, 6, 8 /) + arr(2,1,:) = (/ 1, 2, 5, 6 /) + arr(2,2,:) = (/ 3, 4, 7, 8 /) + arr(3,1,:) = (/ 1, 2, 3, 4 /) + arr(3,2,:) = (/ 5, 6, 7, 8 /) #else /* NDIMS == 3 */ - arr(1,1,:) = (/ 1, 3 /) - arr(1,2,:) = (/ 2, 4 /) - arr(2,1,:) = (/ 1, 2 /) - arr(2,2,:) = (/ 3, 4 /) + arr(1,1,:) = (/ 1, 3 /) + arr(1,2,:) = (/ 2, 4 /) + arr(2,1,:) = (/ 1, 2 /) + arr(2,2,:) = (/ 3, 4 /) #endif /* NDIMS == 3 */ -! iterate over all boundaries of the parent block +! reset the first execution flag +! + first = .false. + + end if + +! iterate over dimensions, sides, and faces ! do i = 1, ndims do j = 1, nsides do k = 1, nfaces -! calculate the right child number +! get the current child index ! p = arr(i,j,k) -! assign the pointer to the current neighbor +! associate a pointer with the neighbor ! - pneigh => pblock%child(p)%ptr%neigh(i,j,k)%ptr + pneigh => pmeta%child(p)%ptr%neigh(i,j,k)%ptr -! assign the right neighbor to the current neighbor pointer +! update the parent neighbor field ! - pblock%neigh(i,j,k)%ptr => pneigh + pmeta%neigh(i,j,k)%ptr => pneigh -! update the neighbor fields of neighbors +! update the neigh field of the neighbor ! if (associated(pneigh)) then + l = 3 - j do p = 1, nfaces - pneigh%neigh(i,l,p)%ptr => pblock + pneigh%neigh(i,l,p)%ptr => pmeta end do - end if - end do - end do - end do + end if ! pneigh is associated -! deallocate child blocks + end do ! nfaces + end do ! nsides + end do ! ndims + +! iterate over children ! do p = 1, nchildren - call metablock_unset_leaf(pblock%child(p)%ptr) - call deallocate_metablock(pblock%child(p)%ptr) - end do -! set the leaf flag of parent block +! remove the child from the meta block list ! - call metablock_set_leaf(pblock) + call remove_metablock(pmeta%child(p)%ptr) + + end do ! nchild + +! update the parent leaf flag +! + call metablock_set_leaf(pmeta) ! reset the refinement flag of the parent block ! - pblock%refine = 0 + call metablock_set_refinement(pmeta, 0) + +! mark the parent to be updated +! + call metablock_set_update(pmeta) + +! mark all neighbors to be updated as well +! + call neighbors_set_update(pmeta) + +#ifdef PROFILE +! stop accounting time for the block derefinement +! + call stop_timer(imd) +#endif /* PROFILE */ !------------------------------------------------------------------------------- ! end subroutine derefine_block -#ifdef DEBUG -!!============================================================================== -!! -!! DEBUG SUBROUTINES -!! ! !=============================================================================== ! -! check_metablock: subroutine checks if the meta block has proper structure +! subroutine SET_LAST_ID: +! ---------------------- +! +! Subroutine sets the last identification number. This subroutine should +! be only used when the job is resumed. +! +! Arguments: +! +! id - the identification number to set; ! !=============================================================================== ! - subroutine check_metablock(pblock, string) + subroutine set_last_id(id) +! import external procedures +! + use error , only : print_error + +! local variables are not implicit by default +! implicit none -! input parameters +! subroutine arguments ! - type(block_meta), pointer, intent(in) :: pblock - character(len=*) , intent(in) :: string + integer(kind=4), intent(in) :: id +! +!------------------------------------------------------------------------------- +! +! check if the new last_id is larger than the already existing +! + if (last_id > id) then -! local variables + call print_error("blocks::set_last_id" & + , "New last_id must be larger than the old one!") + else + +! set the last identification number ! - integer :: p, i, j, k + last_id = id + + end if + +!------------------------------------------------------------------------------- +! + end subroutine set_last_id +! +!=============================================================================== +! +! function GET_LAST_ID: +! -------------------- +! +! Function returns the last identification number. +! +! +!=============================================================================== +! + function get_last_id() result(id) + +! local variables are not implicit by default +! + implicit none + +! return variable +! + integer(kind=4) :: id +! +!------------------------------------------------------------------------------- +! +! set the return value +! + id = last_id + +!------------------------------------------------------------------------------- +! + end function get_last_id +! +!=============================================================================== +! +! function GET_MBLOCKS: +! -------------------- +! +! Function returns the number of meta blocks. +! +! +!=============================================================================== +! + function get_mblocks() result(nr) + +! local variables are not implicit by default +! + implicit none + +! return variable +! + integer(kind=4) :: nr +! +!------------------------------------------------------------------------------- +! +! set the return value +! + nr = mblocks + +!------------------------------------------------------------------------------- +! + end function get_mblocks +! +!=============================================================================== +! +! function GET_DBLOCKS: +! -------------------- +! +! Function returns the number of data blocks. +! +! +!=============================================================================== +! + function get_dblocks() result(nr) + +! local variables are not implicit by default +! + implicit none + +! return variable +! + integer(kind=4) :: nr +! +!------------------------------------------------------------------------------- +! +! set the return value +! + nr = dblocks + +!------------------------------------------------------------------------------- +! + end function get_dblocks +! +!=============================================================================== +! +! function GET_NLEAFS: +! ------------------- +! +! Function returns the number of leafs. +! +! +!=============================================================================== +! + function get_nleafs() result(nr) + +! local variables are not implicit by default +! + implicit none + +! return variable +! + integer(kind=4) :: nr +! +!------------------------------------------------------------------------------- +! +! set the return value +! + nr = nleafs + +!------------------------------------------------------------------------------- +! + end function get_nleafs +! +!=============================================================================== +! +! subroutine SET_BLOCKS_UPDATE: +! ---------------------------- +! +! Subroutine sets the update flag of all meta block in the list. +! +! Arguments: +! +! flag - the flag to be set; +! +!=============================================================================== +! + subroutine set_blocks_update(flag) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + logical, intent(in) :: flag ! local pointers ! - type(block_meta), pointer :: ptemp + type(block_meta), pointer :: pmeta ! !------------------------------------------------------------------------------- ! -! check block ID +! associate the pointer with the first block on the meta block list ! - ptemp => pblock - if (ptemp%id .le. 0 .or. ptemp%id .gt. last_id) then - print *, '' - print *, '' - print *, trim(string) - print *, 'wrong meta block id = ', ptemp%id - stop + pmeta => list_meta + +! iterate over all blocks in the list +! + do while(associated(pmeta)) + +! mark the block for update +! + pmeta%update = flag + +! associate the pointer with the next block on the list +! + pmeta => pmeta%next + + end do ! meta blocks + +!------------------------------------------------------------------------------- +! + end subroutine set_blocks_update +! +!=============================================================================== +! +! subroutine METABLOCK_SET_ID: +! --------------------------- +! +! Subroutine sets the identification number of the meta block pointed by +! the input argument. This subroutine should be used only when resuming jobs. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! id - the identification number to set; +! +!=============================================================================== +! + subroutine metablock_set_id(pmeta, id) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta + integer(kind=4) , intent(in) :: id +! +!------------------------------------------------------------------------------- +! +! set the meta block %id field +! + pmeta%id = id + +! check if the last identification number is smaller than id, if so set +! the value of last_id to id +! + if (last_id < id) last_id = id + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_id +! +!=============================================================================== +! +! subroutine METABLOCK_SET_PROCESS: +! -------------------------------- +! +! Subroutine sets the process number of the meta block pointed by +! the input argument. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! np - the process number; +! +!=============================================================================== +! + subroutine metablock_set_process(pmeta, np) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta + integer(kind=4) , intent(in) :: np +! +!------------------------------------------------------------------------------- +! +! set the block's %process field +! + pmeta%process = np + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_process +! +!=============================================================================== +! +! subroutine METABLOCK_SET_LEVEL: +! ------------------------------ +! +! Subroutine sets the refinement level number of the meta block pointed +! by the input argument. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! lv - the refinement level number; +! +!=============================================================================== +! + subroutine metablock_set_level(pmeta, lv) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta + integer(kind=4) , intent(in) :: lv +! +!------------------------------------------------------------------------------- +! +! set the block's refinement level +! + pmeta%level = lv + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_level +! +!=============================================================================== +! +! subroutine METABLOCK_SET_CONFIGURATION: +! -------------------------------------- +! +! Subroutine sets the children block configuration number of the meta block +! pointed by the input argument. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! cf - the configuration number; +! +!=============================================================================== +! + subroutine metablock_set_configuration(pmeta, cf) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta + integer(kind=4) , intent(in) :: cf +! +!------------------------------------------------------------------------------- +! +! set the block's children configuration number +! + pmeta%conf = cf + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_configuration +! +!=============================================================================== +! +! subroutine METABLOCK_SET_REFINEMENT: +! ----------------------------------- +! +! Subroutine sets the refinement flag of the meta block pointed by +! the input argument. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! rf - the refinement flag; +! +!=============================================================================== +! + subroutine metablock_set_refinement(pmeta, rf) + +! import external procedures +! + use error , only : print_error + +! local variables are not implicit by default +! + implicit none + + type(block_meta), pointer, intent(inout) :: pmeta + integer(kind=4) , intent(in) :: rf +! +!------------------------------------------------------------------------------- +! +! check if the refinement value is correct +! + if (abs(rf) > 1) then + +! print error about wrong refine flag +! + call print_error("blocks::metablock_set_refinement" & + , "The refinement value is wrong!") + + else + +! set the block's refinement field +! + pmeta%refine = rf + end if -! check prev ID +!------------------------------------------------------------------------------- ! - ptemp => pblock%prev - if (associated(ptemp)) then - if (ptemp%id .le. 0 .or. ptemp%id .gt. last_id) then - print *, '' - print *, '' - print *, trim(string) - print *, 'wrong previous block id = ', ptemp%id, pblock%id - stop + end subroutine metablock_set_refinement +! +!=============================================================================== +! +! subroutine METABLOCK_SET_POSITION: +! --------------------------------- +! +! Subroutine sets the position coordinates in the parent block of +! the meta block pointed by the input argument. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! px, py, pz - the block position coordinates; +! +!=============================================================================== +! + subroutine metablock_set_position(pmeta, px, py, pz) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta + integer(kind=4) , intent(in) :: px, py, pz +! +!------------------------------------------------------------------------------- +! +! set the block's position in the parent block +! + pmeta%pos(1) = px + pmeta%pos(2) = py +#if NDIMS == 3 + pmeta%pos(3) = pz +#endif /* NDIMS == 3 */ + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_position +! +!=============================================================================== +! +! subroutine METABLOCK_SET_COORDINATES: +! ------------------------------------ +! +! Subroutine sets the effective resolution coordinates of the meta block +! pointed by the input argument. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! px, py, pz - the effective resolution coordinates; +! +!=============================================================================== +! + subroutine metablock_set_coordinates(pmeta, px, py, pz) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta + integer(kind=4) , intent(in) :: px, py, pz +! +!------------------------------------------------------------------------------- +! +! set the block's effective resolution coordinates +! + pmeta%coords(1) = px + pmeta%coords(2) = py +#if NDIMS == 3 + pmeta%coords(3) = pz +#endif /* NDIMS == 3 */ + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_coordinates +! +!=============================================================================== +! +! subroutine METABLOCK_SET_BOUNDS: +! ------------------------------- +! +! Subroutine sets the physical coordinate bounds of the meta block pointed +! by the input argument. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! xmn, xmx - the coordinate bounds along X; +! ymn, ymx - the coordinate bounds along Y; +! zmn, zmx - the coordinate bounds along Z; +! +!=============================================================================== +! + subroutine metablock_set_bounds(pmeta, xmn, xmx, ymn, ymx, zmn, zmx) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta + real(kind=8) , intent(in) :: xmn, xmx + real(kind=8) , intent(in) :: ymn, ymx + real(kind=8) , intent(in) :: zmn, zmx +! +!------------------------------------------------------------------------------- +! +! set the block's coordinate bounds +! + pmeta%xmin = xmn + pmeta%xmax = xmx + pmeta%ymin = ymn + pmeta%ymax = ymx + pmeta%zmin = zmn + pmeta%zmax = zmx + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_bounds +! +!=============================================================================== +! +! subroutine METABLOCK_SET_LEAF: +! ----------------------------- +! +! Subroutine marks the meta block pointed by the input argument as the leaf. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! +!=============================================================================== +! + subroutine metablock_set_leaf(pmeta) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta +! +!------------------------------------------------------------------------------- +! +! set the block's leaf flag +! + pmeta%leaf = .true. + +! increase the number of leafs +! + nleafs = nleafs + 1 + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_leaf +! +!=============================================================================== +! +! subroutine METABLOCK_UNSET_LEAF: +! ------------------------------- +! +! Subroutine marks the meta block pointed by the input argument as non-leaf. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! +!=============================================================================== +! + subroutine metablock_unset_leaf(pmeta) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta +! +!------------------------------------------------------------------------------- +! +! unset the block's leaf flag +! + pmeta%leaf = .false. + +! decrease the number of leafs +! + nleafs = nleafs - 1 + +!------------------------------------------------------------------------------- +! + end subroutine metablock_unset_leaf +! +!=============================================================================== +! +! subroutine METABLOCK_SET_UPDATE: +! ------------------------------- +! +! Subroutine marks the meta block pointed by the input argument to be updated. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! +!=============================================================================== +! + subroutine metablock_set_update(pmeta) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta +! +!------------------------------------------------------------------------------- +! +! set the block's update flag +! + pmeta%update = .true. + +!------------------------------------------------------------------------------- +! + end subroutine metablock_set_update +! +!=============================================================================== +! +! subroutine METABLOCK_UNSET_UPDATE: +! --------------------------------- +! +! Subroutine marks the meta block pointed by the input argument to not +! be updated. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! +!=============================================================================== +! + subroutine metablock_unset_update(pmeta) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta +! +!------------------------------------------------------------------------------- +! +! unset the block's update flag +! + pmeta%update = .false. + +!------------------------------------------------------------------------------- +! + end subroutine metablock_unset_update +! +!=============================================================================== +!! +!!*** PRIVATE SUBROUTINES **************************************************** +!! +!=============================================================================== +! +!=============================================================================== +! +! subroutine INSERT_METABLOCK_AFTER: +! --------------------------------- +! +! Subroutine allocates memory for one meta block, inserts it to the meta +! block list after the provided pointer and returns a pointer associated +! with it. +! +! Arguments: +! +! pmeta - the pointer associated with the newly appended meta block; +! pprev - the pointer after which the new block has to be inserted; +! +!=============================================================================== +! + subroutine insert_metablock_after(pprev, pmeta) + +! import external procedures +! + use error , only : print_error + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(in) :: pprev + type(block_meta), pointer, intent(out) :: pmeta +! +!------------------------------------------------------------------------------- +! +! allocate memory for the new meta block +! + call allocate_metablock(pmeta) + +! if pprev is associated, insert the new block after it +! + if (associated(pprev)) then + +! associate the %prev and %next pointers +! + pmeta%prev => pprev + pmeta%next => pprev%next + +! update the pointer of the next and previous blocks +! + if (associated(pprev%next)) pprev%next%prev => pmeta + pprev%next => pmeta + +! check if last_meta is associated +! + if (associated(last_meta)) then + +! update the last_meta pointer if necessary +! + if (pprev%id == last_meta%id) last_meta => pmeta + + else + +! strange situation, pprev is associated, but last_meta not +! + call print_error("blocks::intert_metablock_after" & + , "Argument pprev is associated but last_meta is not!") + + end if + + else + +! if pprev is null and list_meta is associated, there is something wrong +! + if (associated(list_meta)) then + +! strange situation, pprev is null but list_meta is associated +! + call print_error("blocks::intert_metablock_after" & + , "Argument pprev is null but list_meta is associated!") + + else + +! pprev and list_meta are nulls, so add the first block to the list by +! associating list_meta and last_meta +! + list_meta => pmeta + last_meta => pmeta + + end if + + end if + +!------------------------------------------------------------------------------- +! + end subroutine insert_metablock_after +! +!=============================================================================== +! +! subroutine INSERT_METABLOCK_BEFORE: +! ---------------------------------- +! +! Subroutine allocates memory for one meta block, inserts it to the meta +! block list before the provided pointer and returns a pointer associated +! with it. +! +! Arguments: +! +! pmeta - the pointer associated with the newly appended meta block; +! pnext - the pointer before which the new block has to be inserted; +! +!=============================================================================== +! + subroutine insert_metablock_before(pnext, pmeta) + +! import external procedures +! + use error , only : print_error + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(in) :: pnext + type(block_meta), pointer, intent(out) :: pmeta +! +!------------------------------------------------------------------------------- +! +! allocate memory for the new meta block +! + call allocate_metablock(pmeta) + +! if pnext is associated, insert the new block before it +! + if (associated(pnext)) then + +! associate the %prev and %next pointers +! + pmeta%prev => pnext%prev + pmeta%next => pnext + +! update the pointer of the next and previous blocks +! + if (associated(pnext%prev)) pnext%prev%next => pmeta + pnext%prev => pmeta + +! check if list_meta is associated +! + if (associated(list_meta)) then + +! update the list_meta pointer if necessary +! + if (pnext%id == list_meta%id) list_meta => pmeta + + else + +! strange situation, pnext is associated, but list_meta not +! + call print_error("blocks::intert_metablock_before" & + , "Argument pnext is associated but list_meta is not!") + + end if + + else + +! if pnext is null and last_meta is associated, there is something wrong +! + if (associated(last_meta)) then + +! strange situation, pnext is null but last_meta is associated +! + call print_error("blocks::intert_metablock_before" & + , "Argument pnext is null but last_meta is associated!") + + else + +! pnext and last_meta are nulls, so add the first block to the list by +! associating list_meta and last_meta +! + list_meta => pmeta + last_meta => pmeta + + end if + + end if + +!------------------------------------------------------------------------------- +! + end subroutine insert_metablock_before +! +!=============================================================================== +! +! function INCREASE_ID: +! -------------------- +! +! Function increases the last identification number by 1 and returns its +! value. +! +! +!=============================================================================== +! + function increase_id() result(id) + +! local variables are not implicit by default +! + implicit none + +! return variable +! + integer(kind=4) :: id +! +!------------------------------------------------------------------------------- +! +! increase the last identification number by 1 +! + last_id = last_id + 1 + +! return its value +! + id = last_id + +!------------------------------------------------------------------------------- +! + end function increase_id +! +!=============================================================================== +! +! subroutine NEIGHBORDS_SET_UPDATE: +! -------------------------------- +! +! Subroutine marks all the neighbors (including edge and corner ones) of +! the meta block pointed by the input argument to be updated. +! +! Arguments: +! +! pmeta - a pointer to the updated meta block; +! +!=============================================================================== +! + subroutine neighbors_set_update(pmeta) + +! local variables are not implicit by default +! + implicit none + +! subroutine arguments +! + type(block_meta), pointer, intent(inout) :: pmeta + +! local pointers +! + type(block_meta), pointer :: pneigh +! +!------------------------------------------------------------------------------- +! +#if NDIMS == 3 +!! 3D CASE: WALK AROUND THE CORNERS +!! +! around (0,0,0) corner +! +! (1,1,1) - X = (1,2,1) - Y = (2,1,2) +! + pneigh => pmeta%neigh(1,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) end if end if -! check next ID +! (3,1,1) - Z = (3,2,1) - X = (1,1,3) ! - ptemp => pblock%next - if (associated(ptemp)) then - if (ptemp%id .le. 0 .or. ptemp%id .gt. last_id) then - print *, '' - print *, '' - print *, trim(string) - print *, 'wrong next block id = ', ptemp%id, pblock%id - stop + pneigh => pmeta%neigh(3,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,1,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) end if end if -! check parent ID +! (2,1,1) - Y = (2,2,1) - Z = (3,1,3) - X = [1,1,4] ! - ptemp => pblock%parent - if (associated(ptemp)) then - if (ptemp%id .le. 0 .or. ptemp%id .gt. last_id) then - print *, '' - print *, '' - print *, trim(string) - print *, 'wrong parent block id = ', ptemp%id, pblock%id - stop - end if - end if + pneigh => pmeta%neigh(2,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) -! check children IDs -! - do p = 1, nchildren - ptemp => pblock%child(p)%ptr - if (associated(ptemp)) then - if (ptemp%id .le. 0 .or. ptemp%id .gt. last_id) then - print *, '' - print *, '' - print *, trim(string) - print *, 'wrong child block id = ', ptemp%id, pblock%id, p - stop + pneigh => pneigh%neigh(3,1,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,1,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) end if end if - end do + end if -! check neighbors IDs +! around (0,1,0) corner ! - do i = 1, ndims - do j = 1, nsides - do k = 1, nfaces - ptemp => pblock%neigh(i,j,k)%ptr - if (associated(ptemp)) then - if (ptemp%id .le. 0 .or. ptemp%id .gt. last_id) then - print *, '' - print *, '' - print *, trim(string) - print *, 'wrong neighbor id = ', ptemp%id, pblock%id, i, j, k - stop - end if - end if - end do - end do - end do +! (2,2,1) + Y = (2,1,1) - X = (1,1,1) +! + pneigh => pmeta%neigh(2,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (3,1,3) - Z = (3,2,3) + Y = (2,2,3) +! + pneigh => pmeta%neigh(3,1,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,2,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (1,1,2) - X = (1,2,2) - Z = (3,1,4) + Y = [2,2,4] +! + pneigh => pmeta%neigh(1,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(3,1,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,2,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + end if + +! around (1,1,0) corner +! +! (1,2,2) + X = (1,1,2) + Y = (2,2,1) +! + pneigh => pmeta%neigh(1,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (3,1,4) - Z = (3,2,4) + X = (1,2,4) +! + pneigh => pmeta%neigh(3,1,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,2,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (2,2,2) + Y = (2,1,2) - Z = (3,1,2) + X = [1,2,3] +! + pneigh => pmeta%neigh(2,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(3,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,2,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + end if + +! around (1,0,0) corner +! +! (2,1,2) - Y = (2,2,2) + X = (1,2,2) +! + pneigh => pmeta%neigh(2,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (3,1,2) - Z = (3,2,2) - Y = (2,1,4) +! + pneigh => pmeta%neigh(3,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,1,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (1,2,1) + X = (1,1,1) - Z = (3,1,1) - Y = [2,1,3] +! + pneigh => pmeta%neigh(1,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(3,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,1,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + end if + +! around (0,0,1) corner +! +! (1,1,3) - X = (1,2,3) + Z = (3,2,2) +! + pneigh => pmeta%neigh(1,1,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(3,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (2,1,3) - Y = (2,2,3) - X = (1,1,4) +! + pneigh => pmeta%neigh(2,1,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,1,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (3,2,1) + Z = (3,1,1) - Y = (2,1,1) - X = [1,1,2] +! + pneigh => pmeta%neigh(3,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + end if + +! around (0,1,1) corner +! +! (2,2,3) + Y = (2,1,3) + Z = (3,2,1) +! + pneigh => pmeta%neigh(2,2,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(3,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (1,1,4) - X = (1,2,4) + Y = (2,2,4) +! + pneigh => pmeta%neigh(1,1,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,2,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (3,2,3) + Z = (3,1,3) - X = (1,1,2) + Z = [2,2,2] +! + pneigh => pmeta%neigh(3,2,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + end if + +! around (1,1,1) corner +! +! (1,2,4) + X = (1,1,4) + Z = (3,2,3) +! + pneigh => pmeta%neigh(1,2,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(3,2,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (2,2,4) + Y = (2,1,4) + X = (1,2,3) +! + pneigh => pmeta%neigh(2,2,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,2,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (3,2,4) + Z = (3,1,4) + Y = (2,2,2) + X = [1,2,1] +! + pneigh => pmeta%neigh(3,2,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + end if + +! around (1,0,1) corner +! +! (2,1,4) - Y = (2,2,4) + Z = (3,2,4) +! + pneigh => pmeta%neigh(2,1,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(3,2,4)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (3,2,2) + Z = (3,1,2) + X = (1,2,1) +! + pneigh => pmeta%neigh(3,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(1,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (1,2,3) + X = (1,1,3) - Y = (2,1,3) + Z = [3,2,3] +! + pneigh => pmeta%neigh(1,2,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,1,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(3,2,3)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + end if + +#else /* NDIMS == 3 */ +!! 2D CASE +!! +! (0,0)-(0,½) edge and (0,0) corner +! + pneigh => pmeta%neigh(1,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (0,½)-(0,1) edge and (0,1) corner +! + pneigh => pmeta%neigh(1,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (1,0)-(1,½) edge and (1,0) corner +! + pneigh => pmeta%neigh(1,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (1,½)-(1,1) edge and (1,1) corner +! + pneigh => pmeta%neigh(1,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + + pneigh => pneigh%neigh(2,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + end if + +! (0,0)-(½,0) edge +! + pneigh => pmeta%neigh(2,1,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + +! (½,0)-(1,0) edge +! + pneigh => pmeta%neigh(2,1,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + +! (0,1)-(½,1) edge +! + pneigh => pmeta%neigh(2,2,1)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if + +! (½,1)-(1,1) edge +! + pneigh => pmeta%neigh(2,2,2)%ptr + if (associated(pneigh)) then + call metablock_set_update(pneigh) + end if +#endif /* NDIMS == 3 */ !------------------------------------------------------------------------------- ! - end subroutine check_metablock -#endif /* DEBUG */ + end subroutine neighbors_set_update !=============================================================================== ! diff --git a/src/boundaries.F90 b/src/boundaries.F90 index 3269373..e1e830a 100644 --- a/src/boundaries.F90 +++ b/src/boundaries.F90 @@ -460,7 +460,7 @@ module boundaries ! check if the block and neighbor belong to the same process, if so, update ! fluxes directly ! - if (pmeta%cpu == nproc .and. pneigh%cpu == nproc) then + if (pmeta%process == nproc .and. pneigh%process == nproc) then #endif /* MPI */ ! update directional flux from the neighbor @@ -523,8 +523,8 @@ module boundaries ! increase the counter for the number of blocks to exchange ! - block_counter(idir,pmeta%cpu,pneigh%cpu) = & - block_counter(idir,pmeta%cpu,pneigh%cpu) + 1 + block_counter(idir,pmeta%process,pneigh%process) = & + block_counter(idir,pmeta%process,pneigh%process) + 1 ! allocate a new info object ! @@ -546,17 +546,17 @@ module boundaries ! check if the list is empty ! - if (associated(block_array(idir,pmeta%cpu,pneigh%cpu)%ptr))& - then + if (associated(block_array(idir,pmeta%process,pneigh%process)%ptr)) then ! if it is, associate the newly created block with it ! - pinfo%prev => block_array(idir,pmeta%cpu,pneigh%cpu)%ptr + pinfo%prev => & + block_array(idir,pmeta%process,pneigh%process)%ptr end if ! %ptr associated ! point the list to the newly created block ! - block_array(idir,pmeta%cpu,pneigh%cpu)%ptr => pinfo + block_array(idir,pmeta%process,pneigh%process)%ptr => pinfo end if ! pmeta and pneigh on local process #endif /* MPI */ @@ -970,30 +970,37 @@ module boundaries ! if (pmeta%leaf .and. pmeta%level == ilev) then +! process only if this block is marked for update +! + if (pmeta%update) then + #ifdef MPI ! check if the current block belongs to the local process ! - if (pmeta%cpu == nproc) then + if (pmeta%process == nproc) then #endif /* MPI */ ! iterate over all neighbors ! - do iside = 1, nsides + do iside = 1, nsides ! assign a neighbor pointer to the current neighbor ! - pneigh => pmeta%neigh(idir,iside,1)%ptr + pneigh => pmeta%neigh(idir,iside,1)%ptr ! make sure that the neighbor is not associated, then apply specific boundaries ! - if (.not. associated(pneigh)) & + if (.not. associated(pneigh)) & call boundary_specific(pmeta%data, idir, iside) - end do ! sides + end do ! sides #ifdef MPI - end if ! block belong to the local process + end if ! block belong to the local process #endif /* MPI */ + + end if ! pmeta is marked for update + end if ! leaf ! assign the pointer to the next block on the list @@ -1137,96 +1144,104 @@ module boundaries ! if (pneigh%level == pmeta%level) then +! process only if this block and its neighbor are marked for update +! + if (pmeta%update .and. pneigh%update) then + ! copy blocks only for the first face ! - if (iface == 1) then + if (iface == 1) then #ifdef MPI ! check if the current meta block and its neighbor belong to the same process ! - if (pmeta%cpu == pneigh%cpu) then + if (pmeta%process == pneigh%process) then ! check if the current meta block belongs to the current process ! - if (pmeta%cpu == nproc) then + if (pmeta%process == nproc) then #endif /* MPI */ ! assign a pointer to the data structure of the current block ! - pdata => pmeta%data + pdata => pmeta%data ! update boundaries of the current block ! - select case(idir) - case(1) - if (iside == 1) then - call boundary_copy(pdata & + select case(idir) + case(1) + if (iside == 1) then + call boundary_copy(pdata & , pneigh%data%u(:,iel:ie,:,:), idir, iside) - else - call boundary_copy(pdata & + else + call boundary_copy(pdata & , pneigh%data%u(:,ib:ibu,:,:), idir, iside) - end if - case(2) - if (iside == 1) then - call boundary_copy(pdata & + end if + case(2) + if (iside == 1) then + call boundary_copy(pdata & , pneigh%data%u(:,:,jel:je,:), idir, iside) - else - call boundary_copy(pdata & + else + call boundary_copy(pdata & , pneigh%data%u(:,:,jb:jbu,:), idir, iside) - end if + end if #if NDIMS == 3 - case(3) - if (iside == 1) then - call boundary_copy(pdata & + case(3) + if (iside == 1) then + call boundary_copy(pdata & , pneigh%data%u(:,:,:,kel:ke), idir, iside) - else - call boundary_copy(pdata & + else + call boundary_copy(pdata & , pneigh%data%u(:,:,:,kb:kbu), idir, iside) - end if + end if #endif /* NDIMS == 3 */ - end select + end select #ifdef MPI - end if ! pmeta on the current process + end if ! pmeta on the current process - else ! block and neighbor belong to different processes + else ! block and neighbor belong to different processes ! increase the counter for number of blocks to exchange ! - block_counter(pmeta%cpu,pneigh%cpu) = & - block_counter(pmeta%cpu,pneigh%cpu) + 1 + block_counter(pmeta%process,pneigh%process) = & + block_counter(pmeta%process,pneigh%process) + 1 ! allocate a new info object ! - allocate(pinfo) + allocate(pinfo) ! fill out its fields ! - pinfo%block => pmeta - pinfo%neigh => pneigh - pinfo%direction = idir - pinfo%side = iside - pinfo%face = iface - pinfo%level_difference = pmeta%level - pneigh%level + pinfo%block => pmeta + pinfo%neigh => pneigh + pinfo%direction = idir + pinfo%side = iside + pinfo%face = iface + pinfo%level_difference = pmeta%level - pneigh%level ! nullify pointer fields ! - nullify(pinfo%prev) - nullify(pinfo%next) + nullify(pinfo%prev) + nullify(pinfo%next) ! if the list is not empty append the newly created block ! - if (associated(block_array(pmeta%cpu,pneigh%cpu)%ptr)) & - pinfo%prev => block_array(pmeta%cpu,pneigh%cpu)%ptr + if (associated(block_array(pmeta%process & + ,pneigh%process)%ptr)) & + pinfo%prev => block_array(pmeta%process & + ,pneigh%process)%ptr ! point the list to the newly created block ! - block_array(pmeta%cpu,pneigh%cpu)%ptr => pinfo + block_array(pmeta%process,pneigh%process)%ptr => pinfo - end if ! block and neighbor belong to different processes + end if ! block and neighbor belong to different processes #endif /* MPI */ - end if ! iface = 1 + end if ! iface = 1 + + end if ! pmeta and pneigh marked for update end if ! neighbor at the same level @@ -1605,121 +1620,129 @@ module boundaries ! if (pmeta%level < pneigh%level) then +! process only if this block and its neighbor are marked for update +! + if (pmeta%update .and. pneigh%update) then + #ifdef MPI ! check if the current meta block and its neighbor belong to the same process ! - if (pmeta%cpu == pneigh%cpu) then + if (pmeta%process == pneigh%process) then ! check if the current meta block belongs to the current process ! - if (pmeta%cpu == nproc) then + if (pmeta%process == nproc) then #endif /* MPI */ ! process each direction separatelly ! - select case(idir) + select case(idir) - case(1) + case(1) ! prepare indices of the neighbor slices used for the boundary update ! - if (iside == 1) then - il = ie - nd + 1 - iu = ie - else - il = ib - iu = ib + nd - 1 - end if - jl = 1 - ju = jm - kl = 1 - ku = km + if (iside == 1) then + il = ie - nd + 1 + iu = ie + else + il = ib + iu = ib + nd - 1 + end if + jl = 1 + ju = jm + kl = 1 + ku = km - case(2) + case(2) ! prepare indices of the neighbor slices used for the boundary update ! - if (iside == 1) then - jl = je - nd + 1 - ju = je - else - jl = jb - ju = jb + nd - 1 - end if - il = 1 - iu = im - kl = 1 - ku = km + if (iside == 1) then + jl = je - nd + 1 + ju = je + else + jl = jb + ju = jb + nd - 1 + end if + il = 1 + iu = im + kl = 1 + ku = km #if NDIMS == 3 - case(3) + case(3) ! prepare indices of the neighbor slices used for the boundary update ! - if (iside == 1) then - kl = ke - nd + 1 - ku = ke - else - kl = kb - ku = kb + nd - 1 - end if - il = 1 - iu = im - jl = 1 - ju = jm + if (iside == 1) then + kl = ke - nd + 1 + ku = ke + else + kl = kb + ku = kb + nd - 1 + end if + il = 1 + iu = im + jl = 1 + ju = jm #endif /* NDIMS == 3 */ - end select + end select ! assign a pointer to the associate data block ! - pdata => pmeta%data + pdata => pmeta%data ! update boundaries of the current block ! - call boundary_restrict(pdata & + call boundary_restrict(pdata & , pneigh%data%u(:,il:iu,jl:ju,kl:ku) & , idir, iside, iface) #ifdef MPI - end if ! block on the current processor + end if ! block on the current processor - else ! block and neighbor on different processors + else ! block and neighbor on different processors ! increase the counter for number of blocks to exchange ! - block_counter(pmeta%cpu,pneigh%cpu) = & - block_counter(pmeta%cpu,pneigh%cpu) + 1 + block_counter(pmeta%process,pneigh%process) = & + block_counter(pmeta%process,pneigh%process) + 1 ! allocate a new info object ! - allocate(pinfo) + allocate(pinfo) ! fill out its fields ! - pinfo%block => pmeta - pinfo%neigh => pneigh - pinfo%direction = idir - pinfo%side = iside - pinfo%face = iface - pinfo%level_difference = pmeta%level - pneigh%level + pinfo%block => pmeta + pinfo%neigh => pneigh + pinfo%direction = idir + pinfo%side = iside + pinfo%face = iface + pinfo%level_difference = pmeta%level - pneigh%level ! nullify pointers ! - nullify(pinfo%prev) - nullify(pinfo%next) + nullify(pinfo%prev) + nullify(pinfo%next) ! if the list is not empty append the created block ! - if (associated(block_array(pmeta%cpu,pneigh%cpu)%ptr)) & - pinfo%prev => block_array(pmeta%cpu,pneigh%cpu)%ptr + if (associated(block_array(pmeta%process & + ,pneigh%process)%ptr)) & + pinfo%prev => block_array(pmeta%process & + ,pneigh%process)%ptr ! point the list to the last created block ! - block_array(pmeta%cpu,pneigh%cpu)%ptr => pinfo + block_array(pmeta%process,pneigh%process)%ptr => pinfo - end if ! block and neighbor on different processors + end if ! block and neighbor on different processors #endif /* MPI */ + end if ! pmeta and pneigh marked for update + end if ! block at lower level than neighbor end if ! neighbor associated @@ -2112,116 +2135,124 @@ module boundaries ! if (pneigh%level < pmeta%level) then +! process only if this block and its neighbor are marked for update +! + if (pmeta%update .and. pneigh%update) then + ! perform update only for the first face, since all faces point the same block ! - if (iface == 1) then + if (iface == 1) then #ifdef MPI ! check if the current meta block and its neighbor belong to the same process ! - if (pmeta%cpu == pneigh%cpu) then + if (pmeta%process == pneigh%process) then ! check if the current meta block belong to the current process ! - if (pmeta%cpu == nproc) then + if (pmeta%process == nproc) then #endif /* MPI */ ! find the neighbor side and face pointing to the current block ! - nside = 3 - iside - nface = 1 - do while(pmeta%id /= & + nside = 3 - iside + nface = 1 + do while(pmeta%id /= & pneigh%neigh(idir,nside,nface)%ptr%id) - nface = nface + 1 - end do + nface = nface + 1 + end do ! prepare indices of the neighbor slices used for the boundary update ! - il = 1 - iu = im - jl = 1 - ju = jm - kl = 1 - ku = km + il = 1 + iu = im + jl = 1 + ju = jm + kl = 1 + ku = km - select case(idir) - case(1) - if (iside == 1) then - il = ie - nh - iu = ie + 1 - else - il = ib - 1 - iu = ib + nh - end if - case(2) - if (iside == 1) then - jl = je - nh - ju = je + 1 - else - jl = jb - 1 - ju = jb + nh - end if - case(3) - if (iside == 1) then - kl = ke - nh - ku = ke + 1 - else - kl = kb - 1 - ku = kb + nh - end if - end select + select case(idir) + case(1) + if (iside == 1) then + il = ie - nh + iu = ie + 1 + else + il = ib - 1 + iu = ib + nh + end if + case(2) + if (iside == 1) then + jl = je - nh + ju = je + 1 + else + jl = jb - 1 + ju = jb + nh + end if + case(3) + if (iside == 1) then + kl = ke - nh + ku = ke + 1 + else + kl = kb - 1 + ku = kb + nh + end if + end select ! assign a pointer to the associated data block ! - pdata => pmeta%data + pdata => pmeta%data ! update boundaries of the current block from its neighbor ! - call boundary_prolong(pdata & + call boundary_prolong(pdata & , pneigh%data%u(:,il:iu,jl:ju,kl:ku) & , idir, iside, nface) #ifdef MPI - end if ! pmeta on the current process + end if ! pmeta on the current process - else ! block and neighbor belong to different processes + else ! block and neighbor belong to different processes ! increase the counter for the number of blocks to exchange ! - block_counter(pmeta%cpu,pneigh%cpu) = & - block_counter(pmeta%cpu,pneigh%cpu) + 1 + block_counter(pmeta%process,pneigh%process) = & + block_counter(pmeta%process,pneigh%process) + 1 ! allocate a new info object ! - allocate(pinfo) + allocate(pinfo) ! fill out its fields ! - pinfo%block => pmeta - pinfo%neigh => pneigh - pinfo%direction = idir - pinfo%side = iside - pinfo%face = iface - pinfo%level_difference = pmeta%level - pneigh%level + pinfo%block => pmeta + pinfo%neigh => pneigh + pinfo%direction = idir + pinfo%side = iside + pinfo%face = iface + pinfo%level_difference = pmeta%level - pneigh%level ! nullify pointers ! - nullify(pinfo%prev) - nullify(pinfo%next) + nullify(pinfo%prev) + nullify(pinfo%next) ! if the list is not empty append the newly created info object to it ! - if (associated(block_array(pmeta%cpu,pneigh%cpu)%ptr)) & - pinfo%prev => block_array(pmeta%cpu,pneigh%cpu)%ptr + if (associated(block_array(pmeta%process & + ,pneigh%process)%ptr)) & + pinfo%prev => block_array(pmeta%process & + ,pneigh%process)%ptr ! point the list to the newly created info object ! - block_array(pmeta%cpu,pneigh%cpu)%ptr => pinfo + block_array(pmeta%process,pneigh%process)%ptr => pinfo - end if ! block and neighbor belong to different processes + end if ! block and neighbor belong to different processes #endif /* MPI */ - end if ! iface = 1 + end if ! iface = 1 + + end if ! pmeta and pneigh marked for update end if ! neighbor belongs to lower level diff --git a/src/domains.F90 b/src/domains.F90 index 5184673..4a3d6c4 100644 --- a/src/domains.F90 +++ b/src/domains.F90 @@ -136,8 +136,8 @@ module domains use blocks , only : pointer_meta, block_meta, block_data & , append_metablock, append_datablock & , link_blocks, metablock_set_leaf & - , metablock_set_config, metablock_set_level & - , metablock_set_coord, metablock_set_bounds + , metablock_set_configuration, metablock_set_level & + , metablock_set_coordinates, metablock_set_bounds use blocks , only : nsides, nfaces use boundaries , only : xlbndry, xubndry, ylbndry, yubndry, zlbndry, zubndry use coordinates, only : xmin, xmax, ymin, ymax, zmin, zmax @@ -228,8 +228,9 @@ module domains ! set the configuration type ! - call metablock_set_config(block_array(loc(1),loc(2),loc(3))%ptr & - , cfg(loc(1),loc(2),loc(3))) + call metablock_set_configuration( & + block_array(loc(1),loc(2),loc(3))%ptr & + , cfg(loc(1),loc(2),loc(3))) ! increase the block number ! @@ -305,7 +306,7 @@ module domains ! set block coordinates ! - call metablock_set_coord(pmeta, il, jl, kl) + call metablock_set_coordinates(pmeta, il, jl, kl) ! set the bounds ! diff --git a/src/driver.F90 b/src/driver.F90 index 6d15a3c..5ad23d7 100644 --- a/src/driver.F90 +++ b/src/driver.F90 @@ -364,7 +364,7 @@ program amun ! initialize block module ! - call initialize_blocks() + call initialize_blocks(master, iret) ! initialize boundaries module and print info ! @@ -633,7 +633,7 @@ program amun ! deallocate block structure ! - call finalize_blocks() + call finalize_blocks(iret) ! finalize the random number generator ! diff --git a/src/evolution.F90 b/src/evolution.F90 index 2ae7439..341c4ac 100644 --- a/src/evolution.F90 +++ b/src/evolution.F90 @@ -211,6 +211,7 @@ module evolution ! include external procedures ! + use blocks , only : set_blocks_update use boundaries , only : boundary_variables use mesh , only : update_mesh @@ -240,6 +241,10 @@ module evolution ! if (toplev > 1) then +! set all meta blocks to not be updated +! + call set_blocks_update(.false.) + ! check refinement and refine ! call update_mesh() @@ -248,11 +253,15 @@ module evolution ! call boundary_variables() - end if ! toplev > 1 - ! update primitive variables ! - call update_variables() + call update_variables() + +! set all meta blocks to be updated +! + call set_blocks_update(.true.) + + end if ! toplev > 1 !------------------------------------------------------------------------------- ! @@ -694,6 +703,7 @@ module evolution ! include external variables ! + use blocks , only : block_meta, list_meta use blocks , only : block_data, list_data ! local variables are not implicit by default @@ -702,22 +712,30 @@ module evolution ! local pointers ! - type(block_data), pointer :: pblock + type(block_meta), pointer :: pmeta + type(block_data), pointer :: pdata ! !------------------------------------------------------------------------------- ! +! associate the pointer with the first block on the data block list +! + pdata => list_data + ! iterate over all data blocks ! - pblock => list_data - do while (associated(pblock)) + do while (associated(pdata)) + +! associate pmeta with the corresponding meta block +! + pmeta => pdata%meta ! convert conserved variables to primitive ones for the current block ! - call update_primitive_variables(pblock%u, pblock%q) + if (pmeta%update) call update_primitive_variables(pdata%u, pdata%q) ! assign pointer to the next block ! - pblock => pblock%next + pdata => pdata%next end do diff --git a/src/io.F90 b/src/io.F90 index 62852a2..30e33fb 100644 --- a/src/io.F90 +++ b/src/io.F90 @@ -2294,12 +2294,12 @@ module io if (associated(pmeta%data) ) dat(l) = 1 id (l) = pmeta%id - cpu(l) = pmeta%cpu + cpu(l) = pmeta%process lev(l) = pmeta%level - cfg(l) = pmeta%config + cfg(l) = pmeta%conf ref(l) = pmeta%refine pos(l,:) = pmeta%pos(:) - cor(l,:) = pmeta%coord(:) + cor(l,:) = pmeta%coords(:) if (pmeta%leaf) lea(l) = 1 @@ -2417,10 +2417,10 @@ module io use blocks , only : block_meta, list_meta use blocks , only : nchildren, nsides, nfaces use blocks , only : get_mblocks - use blocks , only : metablock_set_id, metablock_set_cpu & - , metablock_set_refine, metablock_set_config & + use blocks , only : metablock_set_id, metablock_set_process & + , metablock_set_refinement, metablock_set_configuration & , metablock_set_level, metablock_set_position & - , metablock_set_coord, metablock_set_bounds & + , metablock_set_coordinates, metablock_set_bounds & , metablock_set_leaf use error , only : print_error use hdf5 , only : hid_t, hsize_t @@ -2549,14 +2549,14 @@ module io block_array(id(l))%ptr => pmeta - call metablock_set_id (pmeta, id (l)) - call metablock_set_cpu (pmeta, min(lcpu, cpu(l))) - call metablock_set_refine (pmeta, ref(l)) - call metablock_set_config (pmeta, cfg(l)) - call metablock_set_level (pmeta, lev(l)) - call metablock_set_position(pmeta, pos(l,1), pos(l,2), pos(l,3)) - call metablock_set_coord (pmeta, cor(l,1), cor(l,2), cor(l,3)) - call metablock_set_bounds (pmeta, xmn(l), xmx(l), ymn(l), ymx(l) & + call metablock_set_id (pmeta, id (l)) + call metablock_set_process (pmeta, min(lcpu, cpu(l))) + call metablock_set_refinement (pmeta, ref(l)) + call metablock_set_configuration(pmeta, cfg(l)) + call metablock_set_level (pmeta, lev(l)) + call metablock_set_position (pmeta, pos(l,1), pos(l,2), pos(l,3)) + call metablock_set_coordinates (pmeta, cor(l,1), cor(l,2), cor(l,3)) + call metablock_set_bounds (pmeta, xmn(l), xmx(l), ymn(l), ymx(l) & , zmn(l), zmx(l)) if (lea(l) .eq. 1) call metablock_set_leaf(pmeta) @@ -3065,7 +3065,7 @@ module io ! fill in the coordinate array ! - cor(l,:) = pdata%meta%coord(:) + cor(l,:) = pdata%meta%coords(:) ! fill in the bounds array ! diff --git a/src/makefile b/src/makefile index fb219ca..0b47b66 100644 --- a/src/makefile +++ b/src/makefile @@ -122,7 +122,7 @@ clean-all: clean-bak clean-data clean-exec clean-logs clean-modules \ #------------------------------------------------------------------------------- -blocks.o : blocks.F90 error.o +blocks.o : blocks.F90 error.o timers.o boundaries.o : boundaries.F90 blocks.o coordinates.o equations.o error.o \ interpolations.o mpitools.o timers.o constants.o : constants.F90 diff --git a/src/mesh.F90 b/src/mesh.F90 index 8f3e484..066fca1 100644 --- a/src/mesh.F90 +++ b/src/mesh.F90 @@ -86,7 +86,7 @@ module mesh ! import external procedures and variables ! - use blocks , only : datablock_set_dims + use blocks , only : set_block_dimensions use coordinates , only : xmin, xmax, ymin, ymax, zmin, zmax use coordinates , only : toplev, im, jm, km use equations , only : nv @@ -131,7 +131,7 @@ module mesh ! set data block dimensions ! - call datablock_set_dims(nv, nv, im, jm, km) + call set_block_dimensions(nv, nv, im, jm, km) ! only master prepares the mesh statistics file ! @@ -387,9 +387,9 @@ module mesh ! increase the block level and process counts ! - ldist(pmeta%level) = ldist(pmeta%level) + 1 + ldist(pmeta%level) = ldist(pmeta%level) + 1 #ifdef MPI - cdist(pmeta%cpu+1) = cdist(pmeta%cpu+1) + 1 + cdist(pmeta%process+1) = cdist(pmeta%process+1) + 1 #endif /* MPI */ end if ! the leaf @@ -701,7 +701,7 @@ module mesh ! assign the process number to the current block ! - pmeta%cpu = np + pmeta%process = np ! check if the current block is the leaf ! @@ -714,7 +714,7 @@ module mesh ! if the block belongs to the current process, append a new data block, link it ! to the current meta block and initialize the problem ! - if (pmeta%cpu == nproc) then + if (pmeta%process == nproc) then ! append new data block ! @@ -1175,13 +1175,14 @@ module mesh ! check if the parent blocks is on the same processor as the next block, if not ! move it to the same processor ! - if (pmeta%cpu /= pmeta%next%cpu) pmeta%cpu = pmeta%next%cpu + if (pmeta%process /= pmeta%next%process) & + pmeta%process = pmeta%next%process ! find the case when child blocks are spread across at least 2 processors ! flag = .false. do p = 1, nchildren - flag = flag .or. (pmeta%child(p)%ptr%cpu /= pmeta%cpu) + flag = flag .or. (pmeta%child(p)%ptr%process /= pmeta%process) end do if (flag) then @@ -1192,16 +1193,16 @@ module mesh ! generate the tag for communication ! - itag = pmeta%child(p)%ptr%cpu * nprocs + pmeta%cpu & + itag = pmeta%child(p)%ptr%process * nprocs + pmeta%process & + nprocs + p + 1 ! if the current children is not on the same processor, then ... ! - if (pmeta%child(p)%ptr%cpu /= pmeta%cpu) then + if (pmeta%child(p)%ptr%process /= pmeta%process) then ! if the meta block is on the same process ! - if (pmeta%cpu == nproc) then + if (pmeta%process == nproc) then ! allocate data blocks for children on the processor which will receive data ! @@ -1210,8 +1211,8 @@ module mesh ! receive the data ! - call receive_real_array(size(rbuf), pmeta%child(p)%ptr%cpu & - , itag, rbuf, iret) + call receive_real_array(size(rbuf) & + , pmeta%child(p)%ptr%process, itag, rbuf, iret) ! coppy buffer to data ! @@ -1221,7 +1222,7 @@ module mesh ! send data to the right processor and deallocate data block ! - if (pmeta%child(p)%ptr%cpu == nproc) then + if (pmeta%child(p)%ptr%process == nproc) then ! copy data to buffer ! @@ -1229,7 +1230,8 @@ module mesh ! send data ! - call send_real_array(size(rbuf), pmeta%cpu, itag, rbuf, iret) + call send_real_array(size(rbuf), pmeta%process & + , itag, rbuf, iret) ! deallocate data block ! @@ -1239,7 +1241,7 @@ module mesh ! set the current processor of the block ! - pmeta%child(p)%ptr%cpu = pmeta%cpu + pmeta%child(p)%ptr%process = pmeta%process end if ! if child is are on different processes @@ -1274,7 +1276,7 @@ module mesh if (associated(pparent)) then #ifdef MPI - if (pmeta%cpu .eq. nproc) then + if (pmeta%process .eq. nproc) then #endif /* MPI */ if (.not. associated(pparent%data)) then call append_datablock(pdata) @@ -1313,7 +1315,7 @@ module mesh if (pmeta%refine .eq. 1) then pparent => pmeta #ifdef MPI - if (pmeta%cpu .eq. nproc) then + if (pmeta%process .eq. nproc) then #endif /* MPI */ call refine_block(pmeta, res(pmeta%level + 1,:), .true.) call prolong_block(pparent) @@ -1435,7 +1437,7 @@ module mesh ! check if the block belongs to another process ! - if (pmeta%cpu /= np) then + if (pmeta%process /= np) then ! check if the block is the leaf ! @@ -1443,11 +1445,11 @@ module mesh ! generate a tag for communication ! - itag = pmeta%cpu * nprocs + np + nprocs + 1 + itag = pmeta%process * nprocs + np + nprocs + 1 ! sends the block to the right process ! - if (nproc == pmeta%cpu) then + if (nproc == pmeta%process) then ! copy data to buffer ! @@ -1464,7 +1466,7 @@ module mesh ! send data block ! - end if ! nproc == pmeta%cpu + end if ! nproc == pmeta%process ! receive the block from another process ! @@ -1477,7 +1479,7 @@ module mesh ! receive the data ! - call receive_real_array(size(rbuf), pmeta%cpu, itag, rbuf, iret) + call receive_real_array(size(rbuf), pmeta%process, itag, rbuf, iret) ! coppy the buffer to data block ! @@ -1490,9 +1492,9 @@ module mesh ! set new processor number ! - pmeta%cpu = np + pmeta%process = np - end if ! pmeta%cpu /= np + end if ! pmeta%process /= np ! increase the number of blocks on the current process; if it exceeds the ! allowed number reset the counter and increase the processor number @@ -1882,65 +1884,6 @@ module mesh !------------------------------------------------------------------------------- ! end subroutine restrict_block -#ifdef DEBUG -! -!=============================================================================== -! -! check_mesh: subroutine checks if the block structure is correct -! subroutine CHECK_MESH: -! --------------------- -! -! Subroutine checks if the meta block structure is correct. -! -! Arguments: -! -! string - the identification string; -! -!=============================================================================== -! - subroutine check_mesh(string) - -! import external procedures and variables -! - use blocks , only : block_meta, list_meta - use blocks , only : check_metablock - -! local variables are not implicit by default -! - implicit none - -! input arguments -! - character(len=*), intent(in) :: string - -! local pointers -! - type(block_meta), pointer :: pmeta - -!------------------------------------------------------------------------------- -! -! assign the pointer with the first block on the list -! - pmeta => list_meta - -! iterate over all meta blocks -! - do while(associated(pmeta)) - -! check the current block -! - call check_metablock(pmeta, string) - -! assign the pointer with the next block on the meta block list -! - pmeta => pmeta%next - - end do ! over meta blocks - -!------------------------------------------------------------------------------- -! - end subroutine check_mesh -#endif /* DEBUG */ !=============================================================================== ! diff --git a/src/schemes.F90 b/src/schemes.F90 index b058411..47904f3 100644 --- a/src/schemes.F90 +++ b/src/schemes.F90 @@ -1118,7 +1118,7 @@ module schemes ! calculate the flux along the Z direction ! - do j = jbl, ieu + do j = jbl, jeu do i = ibl, ieu ! copy directional variable vectors to pass to the one dimensional solver @@ -2901,37 +2901,10 @@ module schemes else ! sm = 0 -! conservative variables for the left intermediate state +! when Sₘ = 0 all variables are continuous, therefore the flux reduces +! to the HLL one ! - ui(idn) = wl(idn) / sl - ui(imx) = 0.0d+00 - ui(imy) = ui(idn) * vy - ui(imz) = ui(idn) * vz - ui(ibx) = bx - ui(iby) = by - ui(ibz) = bz - ui(ibp) = ql(ibp,i) - ui(ien) = (wl(ien) - bx * vb) / sl - -! the left intermediate flux -! - f(:,i) = sl * ui(:) - wl(:) - -! conservative variables for the right intermediate state -! - ui(idn) = wr(idn) / sr - ui(imx) = 0.0d+00 - ui(imy) = ui(idn) * vy - ui(imz) = ui(idn) * vz - ui(ibx) = bx - ui(iby) = by - ui(ibz) = bz - ui(ibp) = qr(ibp,i) - ui(ien) = (wr(ien) - bx * vb) / sr - -! the right intermediate flux -! - f(:,i) = 0.5d+00 * (f(:,i) + (sr * ui(:) - wr(:))) + f(:,i) = (sl * wr(:) - sr * wl(:)) / srml end if ! sm = 0 @@ -3348,37 +3321,10 @@ module schemes else ! sm = 0 -! conservative variables for the inmost left intermediate state +! in the case when Sₘ = 0 and Bₓ² > 0, all variables are continuous, therefore +! the flux can be averaged from the Alfvén waves using the simple HLL formula ! - ui(idn) = dnl - ui(imx) = 0.0d+00 - ui(imy) = dnl * vy - ui(imz) = dnl * vz - ui(ibx) = bx - ui(iby) = by - ui(ibz) = bz - ui(ibp) = ql(ibp,i) - ui(ien) = (wcl(ien) - bx * vb) / cal - -! the inmost left intermediate flux -! - f(:,i) = cal * ui(:) - wcl(:) - -! conservative variables for the inmost right intermediate state -! - ui(idn) = dnr - ui(imx) = 0.0d+00 - ui(imy) = dnr * vy - ui(imz) = dnr * vz - ui(ibx) = bx - ui(iby) = by - ui(ibz) = bz - ui(ibp) = qr(ibp,i) - ui(ien) = (wcr(ien) - bx * vb) / car - -! the inmost right intermediate flux -! - f(:,i) = 0.5d+00 * (f(:,i) + (car * ui(:) - wcr(:))) + f(:,i) = (sml * wcr(:) - smr * wcl(:)) / (smr - sml) end if ! sm = 0