Test FOODIE with the integration of Burgers equation.
| Type | Attributes | Name | Initial | |||
|---|---|---|---|---|---|---|
| real(kind=R_P), | parameter | :: | CFL | = | 0.1_R_P | CFL value. |
| integer(kind=I_P), | parameter | :: | Ni | = | 100 | Number of grid nodes. |
| type(command_line_interface) | :: | cli | Command line interface handler. |
|||
| type(burgers) | :: | domain | Burgers field domain. |
|||
| integer(kind=I_P) | :: | error | Error handler. |
|||
| real(kind=R_P), | allocatable | :: | final_state(:) | Final state. |
||
| real(kind=R_P) | :: | h | Space step discretization. |
|||
| real(kind=R_P) | :: | initial_state(1:Ni) | Initial state. |
|||
| real(kind=R_P), | parameter | :: | nu | = | 1._R_P | Viscosity. |
| logical | :: | plots | Flag for activating plots saving. |
|||
| logical | :: | results | Flag for activating results saving. |
|||
| character(len=99) | :: | solver | Solver used. |
|||
| real(kind=R_P), | parameter | :: | t_final | = | 0.6_R_P | Final time. |
| real(kind=R_P) | :: | x(1:Ni) | Nodes values. |
Initialize the field.
Save plots of results.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| character(len=*), | intent(in) | :: | title | Plot title. |
||
| character(len=*), | intent(in) | :: | filename | Output filename. |
Test explicit Adams-Bashforth class of ODE solvers.
Test explicit forward Euler ODE solver.
Test explicit leapfrog class of ODE solvers.
Test explicit low storage Runge-Kutta class of ODE solvers.
Test explicit TVD/SSP Runge-Kutta class of ODE solvers.
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