tem_CoordOfId Function

public pure function tem_CoordOfId(TreeID, offset) result(coord)

The following function provides the spatial index triple for a given ID in the complete mesh, on its refinement level. The returned coordinates start at 0. The fourth entry is the level on which the coordinates are located in the tree. The steps are following: 1. calculate the refinement level, store to coord(4). 2. Calcualte the level offset. 3. the Morton index is obtained by subtracting the offset from treeID. 4. Apply bit interleaving rule to generate coordinates.

Arguments

Type IntentOptional Attributes Name
integer(kind=long_k), intent(in) :: TreeID

input element ID

integer(kind=long_k), intent(in), optional :: offset

First Elem of current level, if known, to avoid recomputations

Return Value integer, (4)

coordinate of element return value


Source Code

  pure function tem_CoordOfId(TreeID, offset) result(coord)
    ! ---------------------------------------------------------------------------
    !> input element ID
    integer(kind=long_k), intent(in) :: TreeID
    !> First Elem of current level, if known, to avoid recomputations
    integer(kind=long_k), intent(in), optional :: offset
    !> coordinate of element return value
    integer :: coord(4)
    ! ---------------------------------------------------------------------------
    integer(kind=long_k) :: tElem ! temporary element
    integer(kind=long_k) :: fak(3)
    integer :: bitlevel
    integer :: i
    ! ---------------------------------------------------------------------------

    coord(:) = 0
    coord(4) = tem_LevelOf(treeID)
    fak(1) = 1
    fak(2) = 2
    fak(3) = 4
    bitlevel = 1

    ! offset has to be subtracted to get the "real" TreeID
    if (present(offset)) then
      tElem = TreeID - offset
    else
      tElem = TreeID - FirstIdAtLevel(coord(4))
    end if

    ! get x coordinate from
    !
    ! x = sum(iLevel=0, inf) { 2**iLevel * mod(tElem / (8**iLevel) ),2)   }
    ! y = sum(iLevel=0, inf) { 2**iLevel * mod(tElem /(2 * (8**iLevel) ),2)   }
    ! z = sum(iLevel=0, inf) { 2**iLevel * mod(tElem /(4 * (8**iLevel) ),2)   }
    !
    do while ((tElem / fak(1)) > 0)
      do i=1,3
        coord(i) = coord(i) + bitlevel * int(mod(tElem / fak(i), 2_long_k))
      end do
      bitlevel = bitlevel * 2
      fak = fak * 8
    end do
  end function tem_CoordOfId