!!******************************************************************************
!!
!!  This file is part of the AMUN source code, a program to perform
!!  Newtonian or relativistic magnetohydrodynamical simulations on uniform or
!!  adaptive mesh.
!!
!!  Copyright (C) 2008-2021 Grzegorz Kowal <grzegorz@amuncode.org>
!!
!!  This program is free software: you can redistribute it and/or modify
!!  it under the terms of the GNU General Public License as published by
!!  the Free Software Foundation, either version 3 of the License, or
!!  (at your option) any later version.
!!
!!  This program is distributed in the hope that it will be useful,
!!  but WITHOUT ANY WARRANTY; without even the implied warranty of
!!  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!!  GNU General Public License for more details.
!!
!!  You should have received a copy of the GNU General Public License
!!  along with this program.  If not, see <http://www.gnu.org/licenses/>.
!!
!!******************************************************************************
!!
!! program: AMUN
!!
!!  AMUN is a code to perform numerical simulations in a fluid approximation on
!!  adaptive mesh for Newtonian and relativistic environments with or without
!!  magnetic field.
!!
!!
!!******************************************************************************
!
program amun

  use blocks         , only : get_nleafs, build_leaf_list
  use boundaries     , only : boundary_variables
  use evolution      , only : advance, initialize_time_step, new_time_step
  use evolution      , only : step, time, dt, dtp, errtol
  use helpers        , only : print_welcome, print_section, print_parameter
  use integrals      , only : store_integrals
  use io             , only : initialize_io, finalize_io
  use io             , only : read_restart_snapshot, write_restart_snapshot
  use io             , only : write_snapshot, update_dtp
  use iso_fortran_env, only : error_unit
  use mesh           , only : generate_mesh, store_mesh_stats
  use mpitools       , only : initialize_mpitools, finalize_mpitools
#ifdef MPI
  use mpitools       , only : reduce_sum
#endif /* MPI */
  use mpitools       , only : master, nprocs, nproc, check_status
  use parameters     , only : read_parameters, finalize_parameters
  use problem        , only : initialize_problem, finalize_problem
  use problem        , only : print_problem_info, execute_problem
  use problem        , only : resumed, nrun, name, rngtype
  use problem        , only : tmax, trun, nwork, nmax
  use random         , only : initialize_random, finalize_random
  use timers         , only : initialize_timers, finalize_timers
  use timers         , only : start_timer, stop_timer, set_timer, get_timer
  use timers         , only : get_timer_total, timer_enabled, timer_description
  use timers         , only : get_count, ntimers
  use user_problem   , only : user_time_statistics
  use workspace      , only : initialize_workspace, finalize_workspace

  implicit none

! the initialization and execution status flags
!
  logical :: initialization_succeeded
  logical :: verbose = .true.
  logical :: proceed = .true.
  integer :: quit    = 0
  integer :: status  = 0

! timer indices
!
  integer               :: iin, iev, itm

! iteration and time variables
!
  integer               :: i, ed, eh, em, es, ec
  integer               :: nsteps = 1
  character(len=80)     :: tmp

  real(kind=8)          :: tbeg = 0.0d+00, thrs
  real(kind=8)          :: tm_curr, tm_exec, tm_conv, tm_last = 0.0d+00

! vectors for better estimation of the remaining simulation time
!
  real(kind=8), dimension(8) :: tprv, tm_prev

#ifndef __GFORTRAN__
! the type of the function SIGNAL should be defined for Intel compiler
!
  integer(kind=4)       :: signal
#endif /* __GFORTRAN__ */

#ifdef SIGNALS
! references to functions handling signals
!
#ifdef __GFORTRAN__
  intrinsic signal
#endif /* __GFORTRAN__ */

! signal definitions
!
  integer, parameter :: SIGERR = -1
  integer, parameter :: SIGINT =  2, SIGABRT =  6, SIGTERM = 15
#endif /* SIGNALS */

! an array to store execution times
!
  real(kind=8), dimension(ntimers) :: tm

!-------------------------------------------------------------------------------
!
#ifdef SIGNALS
! signals handling
!
#ifdef __GFORTRAN__
  if (signal(SIGINT , signal_handler)     == SIGERR) status = 1
  if (signal(SIGABRT, signal_handler)     == SIGERR) status = 2
  if (signal(SIGTERM, signal_handler)     == SIGERR) status = 3
#else /* __GFORTRAN__ */
  if (signal(SIGINT , signal_handler, -1) == SIGERR) status = 1
  if (signal(SIGABRT, signal_handler, -1) == SIGERR) status = 2
  if (signal(SIGTERM, signal_handler, -1) == SIGERR) status = 3
#endif /* __GFORTRAN__ */

  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not initialize the signal handling!"
    call exit(status)
  end if
#endif /* SIGNALS */

! timers
!
  call initialize_timers()

  call set_timer('INITIALIZATION', iin)
  call set_timer('EVOLUTION'     , iev)
  call set_timer('TERMINATION'   , itm)

! parallelization
!
  call initialize_mpitools(status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not initialize module MPITOOLS!"
  end if
  if (check_status(status /= 0)) call exit(status)
  verbose = master

  call start_timer(iin)

  call print_welcome(verbose)
#ifdef MPI
  call print_section(verbose, "Parallelization")
  call print_parameter(verbose, "MPI processes", nprocs)
#endif /* MPI */

  initialization_succeeded = .false.

! read parameters
!
  call read_parameters(verbose, status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not read parameters file!"
  end if
  if (check_status(status /= 0)) go to 6000

! initialize IO early so we can recover a few parameters from
! the restart snapshots
!
  call initialize_io(verbose, status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not initialize module IO!"
  end if
  if (check_status(status /= 0)) go to 5000

! initialize module PROBLEM
!
  call initialize_problem(verbose, status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not initialize module PROBLEM!"
  end if
  if (check_status(status /= 0)) go to 4000

! initialize the random number generator
!
  call initialize_random(rngtype, 1, 0, nproc, status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not initialize module RANDOM!"
  end if
  if (check_status(status /= 0)) go to 3000

! initialize the workspace
!
  call initialize_workspace(nwork, status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Problem initializing WORKSPACE module!"
  end if
  if (check_status(status /= 0)) go to 2000

! print module information
!
  call print_problem_info(verbose)

! restore the simulation for restarted job, otherwise initiate a new one
!
  if (resumed) then

    call read_restart_snapshot(status)
    if (status /= 0) then
      write(error_unit,"('[AMUN::program]: ', a)")                             &
                       "Problem reading the restart snapshot!"
    end if
    if (check_status(status /= 0)) go to 1000

    call build_leaf_list(status)
    if (status /= 0) then
      write(error_unit,"('[AMUN::program]: ', a)")                             &
                       "Problem building the list of leafs!"
    end if
    if (check_status(status /= 0)) go to 1000

    call boundary_variables(time, 0.0d+00)

  else

! generate the initial mesh, refine that mesh to the desired level according to
! the problem setup
!
    call generate_mesh(status)
    if (status /= 0) then
      write(error_unit,"('[AMUN::program]: ', a)")                             &
                       "Problem generating the initial mesh!"
    end if
    if (check_status(status /= 0)) go to 1000

    call boundary_variables(0.0d+00, 0.0d+00)

    call initialize_time_step()

! store the initial time statistics and the snapshot
!
    call store_integrals()
    call user_time_statistics(time)
    call write_snapshot(name)

! store the initial mesh statistics
!
    call store_mesh_stats(step, time)

  end if

  call stop_timer(iin)

  initialization_succeeded = .true.

! time evolution starts here
!
  call start_timer(iev)

! print the initial progress info
!
  if (verbose) then

    tbeg       = time
    tprv(:)    = time
    tm_prev(:) = get_timer_total()

    ed = 9999
    eh =   23
    em =   59
    es =   59

    write(*,*)
    write(*,"(1x,a)"   ) "Evolving the system:"
    write(*,"(4x,'step',5x,'time',11x,'timestep',6x,'err/tol',4x," //          &
            "'blocks',7x,'ETA')")
#ifdef __INTEL_COMPILER
    write(*,"(i8,2(1x,1es14.6),1x,1es10.2,2x,i8,2x," //                        &
            "1i4.1,'d',1i2.2,'h',1i2.2,'m',1i2.2,'s',15x,a1,$)")               &
                 step, time, dt, errtol, get_nleafs(), ed, eh, em, es, char(13)
#else /* __INTEL_COMPILER */
    write(*,"(i8,2(1x,1es14.6),1x,1es10.2,2x,i8,2x," //                        &
            "1i4.1,'d',1i2.2,'h',1i2.2,'m',1i2.2,'s',15x,a1)",advance="no")    &
                 step, time, dt, errtol, get_nleafs(), ed, eh, em, es, char(13)
#endif /* __INTEL_COMPILER */

  end if

  proceed = (nsteps <= nmax) .and. (time < tmax) .and. proceed

  do while(proceed)

    call advance(status)
    if (status /= 0) then
      write(error_unit,"('[AMUN::program]: ', a)")                             &
                       "Advancing to the next step failed!"
    end if
    if (check_status(status /= 0)) go to 1000

! advance the iteration number and time
!
    time   = time   + dt
    step   = step   +  1
    nsteps = nsteps +  1

! update the time step for the precise snapshots
!
    call update_dtp()

! estimate the next time step
!
    call new_time_step()

! store mesh statistics
!
    call store_mesh_stats(step, time)

! store the time statistics and the run snapshot
!
    call store_integrals()
    call user_time_statistics(time)
    call write_snapshot(name)

! estimate the elapsed time to trigger the restart snapshot storage
!
    tm_curr = get_timer_total()
    thrs    = tm_curr / 3.6d+03

    call write_restart_snapshot(thrs, name, nrun, status)
    if (status /= 0) then
      write(error_unit,"('[AMUN::program]: ', a)")                             &
                       "Could not store the restart snapshot!"
    end if
    if (check_status(status /= 0)) go to 1000

! check the conditions for the integration interruption
!
    proceed = (nsteps <= nmax) .and. (time < tmax) .and.                       &
              (thrs <= trun) .and. .not. check_status(quit /= 0)

! print the progress info
!
    if (verbose) then

! update vectors for the remaining time estimation
!
      tm_prev(1)  = tm_curr
      tprv(1)     = time
      tm_prev     = cshift(tm_prev, 1)
      tprv        = cshift(tprv   , 1)

      if (time >= tmax .or. (tm_curr - tm_last) >= 1.0d+00) then

        ec = int((tm_curr - tm_prev(1)) * (tmax - time)                        &
                            / max(1.0d-03, time - tprv(1)), kind = 4)
        es = max(0, min(863999999, ec))
        ed = es / 86400
        es = es - 86400 * ed
        eh = es /  3600
        es = es -  3600 * eh
        em = es /    60
        es = es -    60 * em

#ifdef __INTEL_COMPILER
        write(*,"(i8,2(1x,1es14.6),1x,1es10.2,2x,i8,2x," //                    &
                "1i4.1,'d',1i2.2,'h',1i2.2,'m',1i2.2,'s',15x,a1,$)")           &
                 step, time, dt, errtol, get_nleafs(), ed, eh, em, es, char(13)
#else /* __INTEL_COMPILER */
        write(*,"(i8,2(1x,1es14.6),1x,1es10.2,2x,i8,2x," //  &
                "1i4.1,'d',1i2.2,'h',1i2.2,'m',1i2.2,'s',15x,a1)",advance="no")&
                 step, time, dt, errtol, get_nleafs(), ed, eh, em, es, char(13)
#endif /* __INTEL_COMPILER */

        tm_last = tm_curr

      end if
    end if

  end do ! main loop

  if (verbose) write(*,*)

  call stop_timer(iev)

! store the final restart snapshot
!
  call write_restart_snapshot(1.0d+16, name, nrun, status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not store the restart snapshot!"
  end if

! finalize all modules
!
  1000 continue
  call start_timer(itm)
  call finalize_workspace(status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Problem finalizing WORKSPACE module!"
  end if
  2000 continue
  call finalize_random(status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not finalize module RANDOM!"
  end if
  3000 continue
  call finalize_problem(verbose, status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Could not finalize module PROBLEM!"
  end if
  4000 continue
  call finalize_io(status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Problem finalizing IO module!"
  end if
  5000 continue
  call finalize_parameters()

  if (initialization_succeeded) then

    call stop_timer(itm)

! get the total execution time
!
    tm(1) = get_timer_total()

! get the individual timers
!
    do i = 2, ntimers
      tm(i) = get_timer(i)
    end do

#ifdef MPI
! sum up timers from all processes
!
    call reduce_sum(tm(1:ntimers))
#endif /* MPI */

    if (verbose) then

      tm_conv = 1.0d+02 / tm(1)

      write(*  ,'(a)') ''
      write(tmp,"(a)") "(2x,a32,1x,':',3x,1f16.3,' secs = ',f6.2,' %')"

      write(*,'(1x,a)') 'EXECUTION TIMINGS'
      do i = 2, ntimers
        if (timer_enabled(i)) then
          if (get_count(i) > 0) then
            write(*,tmp) timer_description(i), tm(i), tm_conv * tm(i)
          end if
        end if
      end do

! print the execution times
!
      write(tmp,"(a)") "(1x,a14,20x,':',3x,1f16.3,' secs = ',f6.2,' %')"
      write(*,tmp) 'TOTAL CPU TIME', tm(1)         , 1.0d+02
      write(*,tmp) 'TIME PER STEP ', tm(1) / nsteps, 1.0d+02 / nsteps
#ifdef MPI
      write(*,tmp) 'TIME PER CPU  ', tm(1) / nprocs, 1.0d+02 / nprocs
#endif /* MPI */

! print the total execution time
!
      tm_exec = get_timer_total()

      es = int(tm_exec)
      ec = int(1000 * (tm_exec - es))
      ed = es / 86400
      es = es - 86400 * ed
      eh = es /  3600
      es = es -  3600 * eh
      em = es /    60
      es = es -    60 * em

      write(tmp,"(a)") "(1x,a14,20x,':',3x,i14,'d'" //                         &
                                      ",i3.2,'h',i3.2,'m',i3.2,'.',i3.3,'s')"
      write(*,tmp) 'EXECUTION TIME', ed, eh, em, es, ec

    end if
  end if

  6000 continue
  call finalize_mpitools(status)
  if (status /= 0) then
    write(error_unit,"('[AMUN::program]: ', a)")                               &
                     "Problem finalizing MPITOOLS module!"
  end if

  call finalize_timers()

#ifdef SIGNALS
  contains
#ifdef __INTEL_COMPILER
!
!===============================================================================
!
! function SIGNAL_HANDLER:
! -----------------------
!
!   Function sets variable iterm after receiving a signal.
!
!   Arguments:
!
!     sig_num - the number of the signal to be handled;
!
!===============================================================================
!
integer(kind=4) function signal_handler(sig_num)

  use iso_fortran_env, only : error_unit

  implicit none

  integer(kind=4), intent(in) :: sig_num

!-------------------------------------------------------------------------------
!
  quit = sig_num
  write(error_unit,*)
  write(error_unit,*) "Received signal:", iterm
  write(error_unit,*) "Closing program..."
  signal_handler = 1

!-------------------------------------------------------------------------------
!
end function
#else /* __INTEL_COMPILER */
!
!===============================================================================
!
! subroutine SIGNAL_HANDLER:
! -------------------------
!
!   Subroutine sets variable iterm after receiving a signal.
!
!===============================================================================
!
subroutine signal_handler()

  use iso_fortran_env, only : error_unit

  implicit none

!-------------------------------------------------------------------------------
!
  quit = 15
  write(error_unit,*)
  write(error_unit,*) "Received signal: 2, 9, or 15"
  write(error_unit,*) "Closing program..."
  return

!-------------------------------------------------------------------------------
!
end subroutine
#endif /* __INTEL_COMPILER */
#endif /* SIGNALS */

!===============================================================================
!
end program