burgers.f90 Source File

Test FOODIE with the integration of Burgers equation.

This File Depends On

sourcefile~~burgers.f90~~EfferentGraph sourcefile~burgers.f90 burgers.f90 sourcefile~foodie.f90 foodie.f90 sourcefile~foodie.f90->sourcefile~burgers.f90 sourcefile~type_burgers.f90 type_burgers.f90 sourcefile~foodie.f90->sourcefile~type_burgers.f90 sourcefile~type_burgers.f90->sourcefile~burgers.f90 sourcefile~foodie_integrator_adams_bashforth.f90 foodie_integrator_adams_bashforth.f90 sourcefile~foodie_integrator_adams_bashforth.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_adams_bashforth_moulton.f90 foodie_integrator_adams_bashforth_moulton.f90 sourcefile~foodie_integrator_adams_bashforth.f90->sourcefile~foodie_integrator_adams_bashforth_moulton.f90 sourcefile~foodie_integrator_multistage_object.f90 foodie_integrator_multistage_object.f90 sourcefile~foodie_integrator_multistage_object.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_runge_kutta_low_storage.f90 foodie_integrator_runge_kutta_low_storage.f90 sourcefile~foodie_integrator_multistage_object.f90->sourcefile~foodie_integrator_runge_kutta_low_storage.f90 sourcefile~foodie_integrator_euler_explicit.f90 foodie_integrator_euler_explicit.f90 sourcefile~foodie_integrator_multistage_object.f90->sourcefile~foodie_integrator_euler_explicit.f90 sourcefile~foodie_integrator_runge_kutta_lssp.f90 foodie_integrator_runge_kutta_lssp.f90 sourcefile~foodie_integrator_multistage_object.f90->sourcefile~foodie_integrator_runge_kutta_lssp.f90 sourcefile~foodie_integrator_runge_kutta_ssp.f90 foodie_integrator_runge_kutta_ssp.f90 sourcefile~foodie_integrator_multistage_object.f90->sourcefile~foodie_integrator_runge_kutta_ssp.f90 sourcefile~foodie_integrator_runge_kutta_embedded.f90 foodie_integrator_runge_kutta_embedded.f90 sourcefile~foodie_integrator_multistage_object.f90->sourcefile~foodie_integrator_runge_kutta_embedded.f90 sourcefile~foodie_integrator_lmm_ssp.f90 foodie_integrator_lmm_ssp.f90 sourcefile~foodie_integrator_lmm_ssp.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_adams_bashforth_moulton.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_lmm_ssp_vss.f90 foodie_integrator_lmm_ssp_vss.f90 sourcefile~foodie_integrator_lmm_ssp_vss.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_runge_kutta_low_storage.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_backward_differentiation_formula.f90 foodie_integrator_backward_differentiation_formula.f90 sourcefile~foodie_integrator_backward_differentiation_formula.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_multistage_multistep_object.f90 foodie_integrator_multistage_multistep_object.f90 sourcefile~foodie_integrator_multistage_multistep_object.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_ms_runge_kutta_ssp.f90 foodie_integrator_ms_runge_kutta_ssp.f90 sourcefile~foodie_integrator_multistage_multistep_object.f90->sourcefile~foodie_integrator_ms_runge_kutta_ssp.f90 sourcefile~foodie_error_codes.f90 foodie_error_codes.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_adams_bashforth.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_lmm_ssp.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_adams_bashforth_moulton.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_lmm_ssp_vss.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_runge_kutta_low_storage.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_backward_differentiation_formula.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_ms_runge_kutta_ssp.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_euler_explicit.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_runge_kutta_lssp.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_runge_kutta_ssp.f90 sourcefile~foodie_integrator_adams_moulton.f90 foodie_integrator_adams_moulton.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_adams_moulton.f90 sourcefile~foodie_integrator_leapfrog.f90 foodie_integrator_leapfrog.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_leapfrog.f90 sourcefile~foodie_error_codes.f90->sourcefile~foodie_integrator_runge_kutta_embedded.f90 sourcefile~foodie_integrator_ms_runge_kutta_ssp.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_euler_explicit.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_runge_kutta_lssp.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_runge_kutta_ssp.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_adams_moulton.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_adams_moulton.f90->sourcefile~foodie_integrator_adams_bashforth_moulton.f90 sourcefile~foodie_integrator_leapfrog.f90->sourcefile~foodie.f90 sourcefile~foodie_integrand_object.f90 foodie_integrand_object.F90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_adams_bashforth.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_multistage_object.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_lmm_ssp.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_adams_bashforth_moulton.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_lmm_ssp_vss.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_runge_kutta_low_storage.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_backward_differentiation_formula.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_multistage_multistep_object.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_ms_runge_kutta_ssp.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_euler_explicit.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_runge_kutta_lssp.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_runge_kutta_ssp.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_adams_moulton.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_leapfrog.f90 sourcefile~foodie_integrator_multistep_object.f90 foodie_integrator_multistep_object.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_multistep_object.f90 sourcefile~foodie_integrand_object.f90->sourcefile~foodie_integrator_runge_kutta_embedded.f90 sourcefile~foodie_integrator_multistep_object.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_multistep_object.f90->sourcefile~foodie_integrator_adams_bashforth.f90 sourcefile~foodie_integrator_multistep_object.f90->sourcefile~foodie_integrator_lmm_ssp.f90 sourcefile~foodie_integrator_multistep_object.f90->sourcefile~foodie_integrator_adams_bashforth_moulton.f90 sourcefile~foodie_integrator_multistep_object.f90->sourcefile~foodie_integrator_lmm_ssp_vss.f90 sourcefile~foodie_integrator_multistep_object.f90->sourcefile~foodie_integrator_backward_differentiation_formula.f90 sourcefile~foodie_integrator_multistep_object.f90->sourcefile~foodie_integrator_adams_moulton.f90 sourcefile~foodie_integrator_multistep_object.f90->sourcefile~foodie_integrator_leapfrog.f90 sourcefile~foodie_integrator_object.f90 foodie_integrator_object.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_adams_bashforth.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_multistage_object.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_lmm_ssp.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_adams_bashforth_moulton.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_lmm_ssp_vss.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_runge_kutta_low_storage.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_backward_differentiation_formula.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_multistage_multistep_object.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_ms_runge_kutta_ssp.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_euler_explicit.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_runge_kutta_lssp.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_runge_kutta_ssp.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_adams_moulton.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_leapfrog.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_multistep_object.f90 sourcefile~foodie_integrator_object.f90->sourcefile~foodie_integrator_runge_kutta_embedded.f90 sourcefile~foodie_integrator_runge_kutta_embedded.f90->sourcefile~foodie.f90
Help

Source Code


Source Code

!< Test FOODIE with the integration of Burgers equation.
program integrate_burgers
!-----------------------------------------------------------------------------------------------------------------------------------
!< Test FOODIE with the integration of Burgers equation.
!-----------------------------------------------------------------------------------------------------------------------------------

!-----------------------------------------------------------------------------------------------------------------------------------
use flap, only : command_line_interface
use foodie, only : integrator_adams_bashforth, &
                   integrator_euler_explicit, &
                   integrator_leapfrog, &
                   integrator_runge_kutta_ls, &
                   integrator_runge_kutta_tvd
use IR_Precision, only : R_P, I_P, FR_P, str
use pyplot_module, only :  pyplot
use type_burgers, only : burgers
!-----------------------------------------------------------------------------------------------------------------------------------

!-----------------------------------------------------------------------------------------------------------------------------------
implicit none
type(command_line_interface) :: cli                 !< Command line interface handler.
type(burgers)                :: domain              !< Burgers field domain.
real(R_P),    parameter      :: CFL=0.1_R_P         !< CFL value.
real(R_P),    parameter      :: t_final=0.6_R_P     !< Final time.
real(R_P),    parameter      :: nu=1._R_P           !< Viscosity.
integer(I_P), parameter      :: Ni=100              !< Number of grid nodes.
real(R_P)                    :: h                   !< Space step discretization.
real(R_P)                    :: initial_state(1:Ni) !< Initial state.
real(R_P)                    :: x(1:Ni)             !< Nodes values.
real(R_P), allocatable       :: final_state(:)      !< Final state.
integer(I_P)                 :: error               !< Error handler.
character(99)                :: solver              !< Solver used.
logical                      :: plots               !< Flag for activating plots saving.
logical                      :: results             !< Flag for activating results saving.
!-----------------------------------------------------------------------------------------------------------------------------------

!-----------------------------------------------------------------------------------------------------------------------------------
! setting Command Line Interface
call cli%init(progname    = 'burgers',                                             &
              authors     = 'Fortran-FOSS-Programmers',                            &
              license     = 'GNU GPLv3',                                           &
              description = 'Test FOODIE library on Burgers equation integration', &
              examples    = ["burgers --solver euler --results  ",                 &
                             "burgers --solver ls-runge-kutta -r",                 &
                             "burgers --solver adams-bashforth  ",                 &
                             "burgers --solver all --plots -r   "])
call cli%add(switch='--solver', switch_ab='-s', help='ODE solver used', required=.true., act='store')
call cli%add(switch='--results', switch_ab='-r', help='Save results', required=.false., act='store_true', def='.false.')
call cli%add(switch='--plots', switch_ab='-p', help='Save plots of results', required=.false., act='store_true', def='.false.')
! parsing Command Line Interface
call cli%parse(error=error)
call cli%get(switch='-s', val=solver, error=error) ; if (error/=0) stop
call cli%get(switch='-r', val=results, error=error) ; if (error/=0) stop
call cli%get(switch='-p', val=plots, error=error) ; if (error/=0) stop
! create Burgers field initial state
call init()
! integrate Burgers equation
select case(trim(adjustl(solver)))
case('adams-bashforth')
  call test_ab()
case('euler')
  call test_euler()
case('leapfrog')
  call test_leapfrog()
case('ls-runge-kutta')
  call test_ls_rk()
case('tvd-runge-kutta')
  call test_tvd_rk()
case('all')
  call test_ab()
  call test_euler()
  call test_leapfrog()
  call test_ls_rk()
  call test_tvd_rk()
case default
  print "(A)", 'Error: unknown solver "'//trim(adjustl(solver))//'"'
  print "(A)", 'Valid solver names are:'
  print "(A)", '  + adams-bashforth'
  print "(A)", '  + euler'
  print "(A)", '  + leapfrog'
  print "(A)", '  + ls-runge-kutta'
  print "(A)", '  + tvd-runge-kutta'
  print "(A)", '  + all'
endselect
stop
!-----------------------------------------------------------------------------------------------------------------------------------
contains
  subroutine init()
  !---------------------------------------------------------------------------------------------------------------------------------
  !< Initialize the field.
  !---------------------------------------------------------------------------------------------------------------------------------
  real(R_P), parameter :: pi=4._R_P * atan(1._R_P)  !< Pi greek.
  integer(I_P)         :: i                         !< Space counter.
  !---------------------------------------------------------------------------------------------------------------------------------

  !---------------------------------------------------------------------------------------------------------------------------------
  h = 2._R_P * pi / Ni
  do i=1, Ni
    x(i) = h * (i - 1)
    initial_state(i) = 10._R_P * sin(x(i))
  enddo
  return
  !---------------------------------------------------------------------------------------------------------------------------------
  endsubroutine init

  subroutine save_results(title, filename)
  !---------------------------------------------------------------------------------------------------------------------------------
  !< Save plots of results.
  !---------------------------------------------------------------------------------------------------------------------------------
  character(*), intent(IN) :: title    !< Plot title.
  character(*), intent(IN) :: filename !< Output filename.
  integer(I_P)             :: rawfile  !< Raw file unit for saving results.
  type(pyplot)             :: plt      !< Plot file handler.
  integer(I_P)             :: i        !< Counter.
  !---------------------------------------------------------------------------------------------------------------------------------

  !---------------------------------------------------------------------------------------------------------------------------------
  if (results) then
    open(newunit=rawfile, file=filename//'.dat')
    write(rawfile, '(A)')'# '//title
    write(rawfile, '(A)')'# VARIABLES: "x" "U"'
    do i=1, Ni
      write(rawfile, '(2('//FR_P//',1X))')x(i), final_state(i)
    enddo
    close(rawfile)
  endif
  if (plots) then
    call plt%initialize(grid=.true., xlabel='x', title=title)
    call plt%add_plot(x=x, y=final_state, label='U', linestyle='b-', linewidth=1)
    call plt%savefig(filename//'.png')
  endif
  return
  !---------------------------------------------------------------------------------------------------------------------------------
  endsubroutine save_results

  subroutine test_ab()
  !---------------------------------------------------------------------------------------------------------------------------------
  !< Test explicit Adams-Bashforth class of ODE solvers.
  !---------------------------------------------------------------------------------------------------------------------------------
  type(integrator_runge_kutta_tvd) :: rk_integrator         !< Runge-Kutta integrator.
  integer, parameter               :: rk_stages=5           !< Runge-Kutta stages number.
  type(burgers)                    :: rk_stage(1:rk_stages) !< Runge-Kutta stages.
  type(integrator_adams_bashforth) :: ab_integrator         !< Adams-Bashforth integrator.
  integer, parameter               :: ab_steps=4            !< Adams-Bashforth steps number.
  type(burgers)                    :: previous(1:ab_steps)  !< Previous time steps solutions.
  integer                          :: step                  !< Time steps counter.
  real(R_P)                        :: dt                    !< Time step.
  real(R_P)                        :: t(1:ab_steps)         !< Times.
  integer(I_P)                     :: s                     !< AB steps counter.
  integer(I_P)                     :: ss                    !< AB substeps counter.
  !---------------------------------------------------------------------------------------------------------------------------------

  !---------------------------------------------------------------------------------------------------------------------------------
  print "(A)", 'Integrating Burgers equation by means of Adams-Bashforth class of solvers'
  do s=1, ab_steps
    print "(A)", ' AB-'//trim(str(.true.,s))
    call ab_integrator%init(steps=s)
    select case(s)
    case(1, 2, 3)
      call rk_integrator%init(stages=s)
    case(4)
      call rk_integrator%init(stages=5)
    endselect
    call domain%init(initial_state=initial_state, Ni=Ni, h=h, nu=nu, steps=s)
    dt = domain%dt(CFL=CFL)
    t = 0._R_P
    step = 1
    do while(t(s)<t_final)
      if (s>=step) then
        ! the time steps from 1 to s - 1 must be computed with other scheme...
        call rk_integrator%integrate(U=domain, stage=rk_stage, dt=dt, t=t(s))
        previous(step) = domain
        if (step>1) then
          t(step) = t(step-1) + dt
        else
          t(step) = dt
        endif
      else
        call ab_integrator%integrate(U=domain, previous=previous(1:s), Dt=Dt, t=t)
        do ss=1, s-1
          t(ss) = t(ss + 1)
        enddo
        t(s) = t(s) + dt
      endif
      step = step + 1
    enddo
    final_state = domain%output()
    call save_results(title='FOODIE test: Burgers equation integration, t='//str(n=t_final)//' explicit '//&
                            'Adams-Bashforth '//trim(str(.true., s))//' steps', &
                      filename='burgers_integration-ab-'//trim(str(.true., s)))
  enddo
  print "(A)", 'Finish!'
  return
  !---------------------------------------------------------------------------------------------------------------------------------
  endsubroutine test_ab

  subroutine test_euler()
  !---------------------------------------------------------------------------------------------------------------------------------
  !< Test explicit forward Euler ODE solver.
  !---------------------------------------------------------------------------------------------------------------------------------
  type(integrator_euler_explicit) :: euler_integrator !< Euler integrator.
  real(R_P)                       :: dt               !< Time step.
  real(R_P)                       :: t                !< Time.
  !---------------------------------------------------------------------------------------------------------------------------------

  !---------------------------------------------------------------------------------------------------------------------------------
  print "(A)", 'Integrating Burgers equation by means of explicit Euler solver'
  call domain%init(initial_state=initial_state, Ni=Ni, h=h, nu=nu)
  dt = domain%dt(CFL=CFL)
  t = 0._R_P
  do while(t<t_final)
    call euler_integrator%integrate(U=domain, dt=dt, t=t)
    t = t + dt
  enddo
  final_state = domain%output()
  call save_results(title='FOODIE test: Burgers equation integration, t='//str(n=t_final)//' explicit Euler', &
                    filename='burgers_integration-euler')
  print "(A)", 'Finish!'
  return
  !---------------------------------------------------------------------------------------------------------------------------------
  endsubroutine test_euler

  subroutine test_leapfrog()
  !---------------------------------------------------------------------------------------------------------------------------------
  !< Test explicit leapfrog class of ODE solvers.
  !---------------------------------------------------------------------------------------------------------------------------------
  type(integrator_runge_kutta_tvd) :: rk_integrator         !< Runge-Kutta integrator.
  integer, parameter               :: rk_stages=2           !< Runge-Kutta stages number.
  type(burgers)                    :: rk_stage(1:rk_stages) !< Runge-Kutta stages.
  type(burgers)                    :: filter                !< Filter displacement.
  type(integrator_leapfrog)        :: lf_integrator         !< Leapfrog integrator.
  type(burgers)                    :: previous(1:2)         !< Previous time steps solutions.
  real(R_P)                        :: dt                    !< Time step.
  real(R_P)                        :: t                     !< Time.
  integer                          :: step                  !< Time steps counter.
  !---------------------------------------------------------------------------------------------------------------------------------

  !---------------------------------------------------------------------------------------------------------------------------------
  print "(A)", 'Integrating Burgers equation by means of leapfrog (RAW filtered) class of solvers'
  call lf_integrator%init(nu=1.0_R_P, alpha=0._R_P)
  call rk_integrator%init(stages=rk_stages)
  call domain%init(initial_state=initial_state, Ni=Ni, h=h, nu=nu, steps=2)
  dt = domain%dt(CFL=CFL)
  t = 0._R_P
  step = 1
  do while(t<t_final)
    if (2>=step) then
      ! the time steps from 1 to 2 must be computed with other scheme...
      call rk_integrator%integrate(U=domain, stage=rk_stage, dt=dt, t=t)
      previous(step) = domain
    else
      call lf_integrator%integrate(U=domain, previous=previous, dt=dt, t=t, filter=filter)
    endif
    t = t + dt
    step = step + 1
  enddo
  final_state = domain%output()
  call save_results(title='FOODIE test: Burgers equation integration, t='//str(n=t_final)//' explicit leapfrog scheme',&
                    filename='burgers_integration-lf-RAW-filter')
  print "(A)", 'Finish!'
  return
  !---------------------------------------------------------------------------------------------------------------------------------
  endsubroutine test_leapfrog

  subroutine test_ls_rk()
  !---------------------------------------------------------------------------------------------------------------------------------
  !< Test explicit low storage Runge-Kutta class of ODE solvers.
  !---------------------------------------------------------------------------------------------------------------------------------
  type(integrator_runge_kutta_ls) :: rk_integrator         !< Runge-Kutta integrator.
  integer, parameter              :: rk_stages=5           !< Runge-Kutta stages number.
  integer, parameter              :: registers=2           !< Runge-Kutta stages number.
  type(burgers)                   :: rk_stage(1:registers) !< Runge-Kutta stages.
  real(R_P)                       :: dt                    !< Time step.
  real(R_P)                       :: t                     !< Time.
  integer(I_P)                    :: s                     !< RK stages counter.
  !---------------------------------------------------------------------------------------------------------------------------------

  !---------------------------------------------------------------------------------------------------------------------------------
  print "(A)", 'Integrating Burgers equation by means of low storage (2N) Runge-Kutta class of solvers'
  do s=1, rk_stages
    if (s==2) cycle ! 2 stages not yet implemented
    if (s==3) cycle ! 3 stages not yet implemented
    if (s==4) cycle ! 4 stages not yet implemented
    print "(A)", ' RK-'//trim(str(.true.,s))
    call rk_integrator%init(stages=s)
    call domain%init(initial_state=initial_state, Ni=Ni, h=h, nu=nu)
    dt = domain%dt(CFL=CFL)
    t = 0._R_P
    do while(t<t_final)
      call rk_integrator%integrate(U=domain, stage=rk_stage, dt=dt, t=t)
      t = t + dt
    enddo
    final_state = domain%output()
    call save_results(title='FOODIE test: Burgers equation integration, t='//str(n=t_final)//' explicit low storage Runge-Kutta '//&
                            trim(str(.true., s))//' stages', &
                      filename='burgers_integration-lsrk-'//trim(str(.true., s)))
  enddo
  print "(A)", 'Finish!'
  return
  !---------------------------------------------------------------------------------------------------------------------------------
  endsubroutine test_ls_rk

  subroutine test_tvd_rk()
  !---------------------------------------------------------------------------------------------------------------------------------
  !< Test explicit TVD/SSP Runge-Kutta class of ODE solvers.
  !---------------------------------------------------------------------------------------------------------------------------------
  type(integrator_runge_kutta_tvd) :: rk_integrator         !< Runge-Kutta integrator.
  integer, parameter               :: rk_stages=5           !< Runge-Kutta stages number.
  type(burgers)                    :: rk_stage(1:rk_stages) !< Runge-Kutta stages.
  real(R_P)                        :: dt                    !< Time step.
  real(R_P)                        :: t                     !< Time.
  integer(I_P)                     :: s                     !< RK stages counter.
  !---------------------------------------------------------------------------------------------------------------------------------

  !---------------------------------------------------------------------------------------------------------------------------------
  print "(A)", 'Integrating Burgers equation by means of TVD/SSP Runge-Kutta class of solvers'
  do s=1, rk_stages
    if (s==4) cycle ! 4 stages not yet implemented
    print "(A)", ' RK-'//trim(str(.true.,s))
    call rk_integrator%init(stages=s)
    call domain%init(initial_state=initial_state, Ni=Ni, h=h, nu=nu)
    dt = domain%dt(CFL=CFL)
    t = 0._R_P
    do while(t<t_final)
      call rk_integrator%integrate(U=domain, stage=rk_stage(1:s), dt=dt, t=t)
      t = t + dt
    enddo
    final_state = domain%output()
    call save_results(title='FOODIE test: Burgers equation integration, t='//str(n=t_final)//' explicit TVD Runge-Kutta '//&
                            trim(str(.true., s))//' stages', &
                      filename='burgers_integration-tvdrk-'//trim(str(.true., s)))
  enddo
  print "(A)", 'Finish!'
  return
  !---------------------------------------------------------------------------------------------------------------------------------
  endsubroutine test_tvd_rk
endprogram integrate_burgers