rk_step Subroutine

  1. pf_imk.f90
  2. pf_mod_imk
  3. rk_step

public subroutine rk_step(this, pf, t0, dt)

Uses

  • proc~~rk_step~~UsesGraph proc~rk_step rk_step module~pf_mod_hooks pf_mod_hooks proc~rk_step->module~pf_mod_hooks module~pf_mod_timer pf_mod_timer proc~rk_step->module~pf_mod_timer module~pf_mod_dtype pf_mod_dtype module~pf_mod_hooks->module~pf_mod_dtype module~pf_mod_timer->module~pf_mod_dtype iso_c_binding iso_c_binding module~pf_mod_dtype->iso_c_binding

Arguments

Type IntentOptional AttributesName
class(pf_imk_t), intent(inout) :: this

Inputs
type(pf_pfasst_t), intent(inout), target:: pf

PFASST structure
real(kind=pfdp), intent(in) :: t0

Time at beginning of time step
real(kind=pfdp), intent(in) :: dt

time step size

Calls

proc~~rk_step~~CallsGraph proc~rk_step rk_step proc~call_hooks call_hooks proc~rk_step->proc~call_hooks proc~pf_residual pf_residual proc~rk_step->proc~pf_residual proc~start_timer start_timer proc~call_hooks->proc~start_timer proc~end_timer end_timer proc~call_hooks->proc~end_timer proc~pf_residual->proc~start_timer proc~pf_residual->proc~end_timer

Called by

proc~~rk_step~~CalledByGraph proc~rk_step rk_step proc~imk_sweep imk_sweep proc~imk_sweep->proc~rk_step

Contents

Source Code

rk_step

Source Code

  subroutine rk_step(this, pf, t0, dt)
    use pf_mod_timer
    use pf_mod_hooks

    !>  Inputs
    class(pf_imk_t), intent(inout) :: this
    type(pf_pfasst_t), intent(inout),target :: pf      !!  PFASST structure
    real(pfdp),        intent(in   ) :: t0             !!  Time at beginning of time step
    real(pfdp),        intent(in   ) :: dt             !!  time step size

    !>  Local variables
    class(pf_level_t), pointer :: lev    !!  points to current level

    integer     :: m        !!  Loop variables
    real(pfdp)  :: t        !!  Time at nodes

    t = t0 + dt

    lev => pf%levels(1)
    do m = 1, 5
       call lev%Q(m)%setval(0.0_pfdp)
    end do

    call call_hooks(pf, 1, PF_PRE_SWEEP)

    call lev%Q(1)%copy(lev%q0, flags=1)
    call this%f_eval(lev%Q(1), t, lev%index, this%A(1))
    ! commutator_p flags=1 hack copies Q(1)%y -> Q(1)%array
    ! all subsequent RK stages are done on Q(m)%array
    call this%commutator_p(this%A(1), lev%Q(1), lev%F(1,1), flags=1)

    call lev%Q(2)%axpy(1.0_pfdp, lev%Q(1))
    call lev%Q(2)%axpy(0.5_pfdp*dt, lev%F(1,1))
    call this%f_eval(lev%Q(2), t, lev%index, this%A(2))
    call this%commutator_p(this%A(2), lev%Q(2), lev%F(2,1))

    call lev%Q(3)%axpy(1.0_pfdp, lev%Q(1))
    call lev%Q(3)%axpy(0.5_pfdp*dt, lev%F(2,1))
    call this%f_eval(lev%Q(3), t, lev%index, this%A(3))
    call this%commutator_p(this%A(3), lev%Q(3), lev%F(3,1))

    call lev%Q(4)%axpy(1.0_pfdp, lev%Q(1))
    call lev%Q(4)%axpy(dt, lev%F(3,1))
    call this%f_eval(lev%Q(4), t, lev%index, this%A(4))
    call this%commutator_p(this%A(4), lev%Q(4), lev%F(4,1))

    call lev%Q(5)%axpy(1.0_pfdp, lev%Q(1))
    call lev%Q(5)%axpy(dt/6.0_pfdp, lev%F(1,1))
    call lev%Q(5)%axpy(dt/3.0_pfdp, lev%F(2,1))
    call lev%Q(5)%axpy(dt/3.0_pfdp, lev%F(3,1))
    call lev%Q(5)%axpy(dt/6.0_pfdp, lev%F(4,1))

    ! commutator_p flags=2 hack copies Q(m)%array -> Q(m)%y
    ! only the Q argument in this case matters
    ! in order to get solution back to qend
    call this%commutator_p(this%A(1), lev%Q(1), lev%F(1,1), flags=2)
    call this%commutator_p(this%A(2), lev%Q(2), lev%F(1,1), flags=2)
    call this%commutator_p(this%A(3), lev%Q(3), lev%F(1,1), flags=2)
    call this%commutator_p(this%A(4), lev%Q(4), lev%F(1,1), flags=2)
    call this%commutator_p(this%A(5), lev%Q(5), lev%F(1,1), flags=2)

    call pf_residual(pf, lev, dt)
    call lev%qend%copy(lev%Q(lev%nnodes), flags=1)

    call call_hooks(pf, 1, PF_POST_SWEEP)
  end subroutine rk_step