! See copyright notice in the COPYRIGHT file.
! ****************************************************************************** !
!> author: Kannan Masilamani
!! Auxiliary functionality for musubi harvesting
!!
module mus_hvs_aux_module

  ! include musubi modules
  use tem_aux_module,                only: tem_abort
  use mus_param_module,              only: mus_param_type
  use mus_geom_module,               only: mus_geom_type
  use mus_tools_module,              only: check_streaming_layout
  use mus_interpolate_module,        only: mus_init_interpolate
  use mus_scheme_type_module,        only: mus_scheme_type
  use mus_interpolate_tools_module,  only: debug_dependencies, dump_intpLists
  use mus_time_module,               only: mus_time_homogenize
  use mus_fluid_module,              only: mus_alloc_mrt, mus_init_fluid_new
  use mus_tracking_module,           only: mus_init_tracker

  ! include treelm modules
  use tem_debug_module,         only: dbgUnit
  use tem_restart_module,       only: tem_init_restart
  use tem_tools_module,         only: tem_horizontalSpacer
  use tem_tracking_module,      only: tem_tracker
  use tem_debug_module,         only: main_debug
  use tem_logging_module,       only: logUnit
  use tem_spacetime_fun_module, only: tem_create_subTree_of_st_funList

  implicit none
  private

  public :: mus_hvs_init_aux

contains

! ****************************************************************************** !
  !> Init auxiliary features such as interpolation boundaries, restart and 
  !! the tracker
  subroutine mus_hvs_init_aux( scheme, geometry, params )
    ! ---------------------------------------------------------------------------
    !> container for the scheme
    type(mus_scheme_type), target, intent(inout)     :: scheme
    !> geometry infomation
    type(mus_geom_type), intent(inout)       :: geometry
    !> global parameters
    type(mus_param_type),intent(inout)       :: params
    ! ---------------------------------------------------------------------------
    integer :: iLevel
    integer :: minLevel, maxLevel
    ! ---------------------------------------------------------------------------
    minLevel = geometry%tree%global%minLevel
    maxLevel = geometry%tree%global%maxLevel

    ! create subTree for all spacetime function in the linked list of
    ! spacetime function
    call tem_create_subTree_of_st_funList(     &
      &       me      = scheme%st_funList,     &
      &       tree    = geometry%tree,         &
      &       bc_prop = geometry%boundary,     &
      &       stencil = scheme%layout%fStencil )

    call tem_horizontalSpacer(fUnit = logUnit(1))

    ! verify correct settings for the streaming layout
    call check_streaming_layout( minLevel, maxLevel )

    ! When the restart is read from separate restart table, we start the 
    ! simulation from the time step given in restart file. In the case,
    ! when restart is read from initial condition table, the simulation start
    ! time step is taken as the one defined in configuration file
    call mus_time_homogenize( me = params%general%simControl%now,              &
      &                       dt = params%physics%dtLvl( maxLevel ),           &
      &              readRestart = params%general%restart%controller           &
      &                             %readRestart )

    ! initialize the restart
    write(logUnit(1),*) 'Initializing restart...'
    call tem_init_restart( me           = params%general%restart, &
      &                    solver       = params%general%solver,  &
      ! &                    varSys       = scheme%varSys,          &
      &                    varMap       = scheme%stateVarMap,     &
      &                    tree         = geometry%tree           )

    ! scheme loop to initialize tracking, boundary, interpolation and output

    !> initialize fluid type which contains relaxation parameter
    !! and function pointers to get mrt paramter and nonEqScaling factor
    !! for interpolation
    select case( trim(scheme%header%kind) )
    case('lbm', 'lbm_incomp', 'lbm_nNwtn', 'lbm_incomp_nNwtn', &
      &  'isotherm_acEq')
      if (scheme%nFields > 1) then
        call tem_abort('chosen scheme kind supports only one field')
      end if
      ! initialize fluid viscosity relaxation paramters 
      call mus_init_fluid_new(                            &
        & me           = scheme%field(1)%fieldProp%fluid, &
        & physics      = params%physics,                  &
        & schemeHeader = scheme%header,                   &
        & minLevel     = minLevel,                        &
        & maxLevel     = maxLevel,                        &
        & nSolve       = scheme%pdf(:)%nElems_solve,      &
        & levelDesc    = scheme%levelDesc(:),             &
        & tNow         = params%general%simControl%now    )
    end select  

    ! \todo KM: remove old mus_alloc_mrt.
    ! also mus_setParameters called twice at initial time once here and
    ! in the beginning of mus_control_module
    if ( scheme%header%kind(1:3) == 'lbm' ) then
      do iLevel = minLevel, maxLevel
        ! allocate fluid mrt array
        call mus_alloc_mrt( mrt = scheme%field(1)%fieldProp%   &
          &                                fluid%mrt(iLevel),  &
          &                 QQ    = scheme%layout%fStencil%QQ  )
        ! set relaxation parameter for MRT
        ! \todo remove
        !KM!call mus_setParameters( scheme, iLevel, params%general%simControl%now, &
        !KM!  &                     params%physics )

      end do
    end if

    ! initialize tracking objects.
    call mus_init_tracker( scheme    = scheme,   &
      &                    geometry  = geometry, &
      &                    params    = params    )

    if( minLevel /= maxlevel ) then
      write(logUnit(1),*) 'Initializing interpolation...'
      ! initialize the interpolation
      call mus_init_interpolate(                             &
        &             intp        = scheme%intp,             &
        &             leveldesc   = scheme%leveldesc,        &
        &             scheme_kind = scheme%header%kind,      &
        &             stencil     = scheme%layout%fStencil,  &
        &             minLevel    = minLevel,                &
        &             maxLevel    = maxLevel,                &
        &             turbActive  = scheme%turbulence%active )
    end if
    if( main_debug%debugDependencies) then
      call debug_dependencies( scheme%intp, scheme%levelDesc, &
        &                      geometry%tree, params%general%proc%rank )
      call dump_intpLists( minLevel, maxLevel, scheme%intp%config%order,   &
        &                  scheme%levelDesc, params%general%proc%rank )
    end if

    call tem_horizontalSpacer( after = 1, fUnit = logUnit(1) )

  end subroutine mus_hvs_init_aux
! ****************************************************************************** !


end module mus_hvs_aux_module
! ****************************************************************************** !
