tem_IdOfCoord Function

public pure function tem_IdOfCoord(coord, offset) result(nTreeID)

This function calculates the tree ID for a given x,y and z on the given refinement level. This ID in the complete tree does not have to be in the list of leafs (treeIDlist) An example of this procedure is following: 1. Convert coordinates into binary representation: (x,y,z) = (5,9,1) = (0101,1001,0001) 2. Interleaving the bits by 3 digits for each direction: x(0101): 0 1 0 1 y(1001): 1 0 0 1 z(0001): 0 0 0 1 Combining these bits results in a binary number: 010 001 000 111, which is 1095 when seen as a 10-base number. 3. This number is the cell position in a single level sorted element list. Its treeID can be obtained by adding the correspoding level offset.

Periodic displacement for requested elements outside the cube

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: coord(4)

3D Coordinates to translate
integer(kind=long_k), intent(in), optional :: offset

Return Value integer(kind=long_k)

resulting treeID


Source Code

  pure function tem_IdOfCoord(coord, offset) result(nTreeID)
    ! ---------------------------------------------------------------------------
    !> 3D Coordinates to translate
    integer, intent(in) :: coord(4)
    !>
    integer(kind=long_k), intent(in), optional :: offset
    !> resulting treeID
    integer(kind=long_k) :: nTreeID
    ! ---------------------------------------------------------------------------
    integer :: iLevel    ! level iterator
    ! this variable defines the upper threshold for x, y and z coordinates
    ! for a given refinement level e.g. for level = 1 only x, y and z
    ! from 0 to 1 are allowed, otherwise this neighbor is outside of the
    ! geometry
    integer :: upperBoundary
    integer :: myCoord(3)
    integer :: fak2
    integer(kind=long_k) :: fak8
    ! ---------------------------------------------------------------------------

    ! the offset has to be added to the level summation
    ! @todo: level offset is hard coded in array firstIdAtLevel,
    !        thus no need to pass offset as an argument.
    if (present(offset)) then
      nTreeID = offset
    else
      nTreeID = firstIdAtLevel(coord(4))
    end if

    !! Periodic displacement for requested elements outside the cube
    upperBoundary = 2**coord(4)
    myCoord = mod(coord(1:3) + upperBoundary, upperBoundary)

    ! To compute the tree ID of a neighbor on the same refinement level the
    ! following formula is used:
    ! newTreeID = sum(i=0, level) {8**i * (mod( (x / (2**level) ) , 2)
    !                                      + 2 * mod( ( y / (2**level) ) , 2)
    !                                      + 4 * mod( (z / (2**level) ) , 2) ) }
    fak8 = 1
    fak2 = 1
    do iLevel = 0, coord(4) - 1
      nTreeID = nTreeID + fak8 * (      mod( (myCoord(1) / fak2), 2) &
        &                         + 2 * mod( (myCoord(2) / fak2), 2) &
        &                         + 4 * mod( (myCoord(3) / fak2), 2) )
      fak2 = fak2 * 2
      fak8 = fak8 * 8
    end do

  end function tem_IdOfCoord