tem_tracking_getData Subroutine

  subroutine tem_tracking_getData(varMap, subTree, varSys, mesh, time, nDofs, &
    &                             res)
    !> varMap from tem_tracking_instance_type
    type(tem_varMap_type) :: varMap

    !> subTree from tem_tracking_instance_type
    type(tem_subTree_type) :: subTree

    !> Variable system describing available data.
    type(tem_varsys_type), intent(in) :: varsys

    !> Mesh definition of the input data.
    type(treelmesh_type), intent(in) :: mesh

    !> Time information for the current data.
    type(tem_time_type), intent(in) :: time

    !> Number of degrees of freedom.
    integer, intent(in) :: nDofs

    !> Tracked data, has to match the subtree definition.
    !!
    !! The memory layout is like this:
    !!  1. All variable components
    !!  2. nDofs
    !!  3. nElems (subtree%nElems)
    real(kind=rk), intent(out) :: res(:)
    ! -------------------------------------------------------------------- !
    integer :: maxComponents
    integer :: nComponents
    integer :: compOff
    integer :: elemOff
    integer :: nElems
    integer :: nScalars, nVars
    integer :: elemSize
    integer :: nChunks
    integer :: chunksize
    integer :: nChunkElems
    integer :: res_size
    integer :: buf_start, buf_end
    integer :: e_start, d_start, t_start
    integer :: iElem, iChunk, iDoF, iVar
    integer :: varpos
    real(kind=rk), allocatable :: tmpdat(:)
    ! -------------------------------------------------------------------- !

    nElems = subTree%nElems
    nScalars = varMap%nScalars
    nVars = varMap%varpos%nVals

    ! Need to obtain the data variable for variable, and store it in an
    ! intermediate array, because all components should be put together in the
    ! res array.
    ! The temporary array therefore needs to be sufficiently large to store the
    ! maximal number of components.
    maxComponents = maxval(varSys%method%val(varMap%varPos &
      &                          %val(:nVars))%nComponents )

    ! Number of elements to fit into a single chunk.
    chunkSize = min( io_buffer_size / (maxComponents*nDofs), nElems )

    ! Size of a single element
    elemsize = nScalars*nDofs

    if ( (nElems > 0) .and. (chunkSize == 0) ) then
      write(logUnit(0),*) 'The chosen io_buffer_size of ', io_buffer_size
      write(logUnit(0),*) 'is too small for outputting ', maxComponents
      write(logUnit(0),*) 'scalar values with ', nDofs
      write(logUnit(0),*) 'degrees of freedom!'
      write(logUnit(0),*) 'Please increase the io_buffer_size to at least ', &
        &                 real(maxComponents*nDofs) / real(131072), ' MB!'
      call tem_abort()
    end if

    ! Using a temporary array to store the variables and transfer them to res
    ! in the correct ordering afterwards.
    allocate(tmpdat(chunkSize*maxComponents*nDofs))

    nChunks = 0
    if (chunkSize > 0) then
      nChunks = ceiling( real(nElems, kind=rk)      &
        &                / real(chunkSize, kind=rk) )
    end if

    chunks: do iChunk=1,nChunks
      elemOff = ((iChunk-1)*chunkSize)
      nChunkElems = min(chunkSize, nElems - elemOff)
      buf_start = elemOff + 1
      buf_end = elemOff + nChunkElems

      compOff = 0
      vars: do iVar=1,varMap%varPos%nVals
        varpos = varMap%varPos%val(iVar)
        nComponents = varSys%method%val(varPos)%nComponents
        res_size = nChunkElems * nDofs * nComponents
        ! derive the quantities for all the elements in the current chunk
        call varSys%method%val(varpos)%get_element(                          &
          &                                varSys  = varSys,                 &
          &                                elemPos = subtree%map2global(     &
          &                                              buf_start:buf_end), &
          &                                time    = time,                   &
          &                                tree    = mesh,                   &
          &                                nElems  = nChunkElems,            &
          &                                nDofs   = nDofs,                  &
          &                                res     = tmpdat(:res_size)       )
        do iElem=1,nChunkElems
          e_start = (elemOff+iElem-1)*elemsize
          t_start = (iElem-1)*nComponents*nDofs
          do iDof=1,nDofs
            d_start = (iDof-1)*nScalars + compOff
            res( (e_start+d_start+1) : (e_start+d_start+nComponents) ) &
              &  = tmpdat( t_start + (iDof-1)*nComponents + 1 &
              &            :t_start + iDof*nComponents        )
          end do
        end do
        ! Increase the component offset for the next variables.
        compOff = compOff + nComponents
      end do vars
    end do chunks

  end subroutine tem_tracking_getData