misdc_sweep Subroutine

public subroutine misdc_sweep(this, pf, lev, t0, dt)

Uses

  • proc~~misdc_sweep~~UsesGraph proc~misdc_sweep misdc_sweep module~pf_mod_timer pf_mod_timer proc~misdc_sweep->module~pf_mod_timer module~pf_mod_dtype 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_misdc_t), intent(inout) :: this
type(pf_pfasst_t), intent(inout) :: pf
class(pf_level_t), intent(inout) :: lev
real(kind=pfdp), intent(in) :: t0
real(kind=pfdp), intent(in) :: dt

Calls

proc~~misdc_sweep~~CallsGraph proc~misdc_sweep misdc_sweep proc~start_timer start_timer proc~misdc_sweep->proc~start_timer proc~misdc_evaluate misdc_evaluate proc~misdc_sweep->proc~misdc_evaluate proc~end_timer end_timer proc~misdc_sweep->proc~end_timer

Contents

Source Code


Source Code

  subroutine misdc_sweep(this, pf, lev, t0, dt)
    use pf_mod_timer
    class(pf_misdc_t), intent(inout)    :: this
    type(pf_pfasst_t),    intent(inout) :: pf
    real(pfdp),           intent(in)    :: dt, t0
    class(pf_level_t),     intent(inout) :: lev

    integer                        :: m, n
    real(pfdp)                     :: t
    real(pfdp)                     :: dtsdc(1:lev%nnodes-1)
    class(pf_encap_t), allocatable :: S3(:)
    class(pf_encap_t), allocatable :: rhs

    call start_timer(pf, TLEVEL+lev%index-1)
    
    ! compute integrals and add fas correction
    do m = 1, lev%nnodes-1
       call lev%S(m)%setval(0.0_pfdp)
       do n = 1, lev%nnodes
          call lev%S(m)%axpy(dt*this%SdiffE(m,n), lev%F(n,1))
          call lev%S(m)%axpy(dt*this%SdiffI(m,n), lev%F(n,2))
          call lev%S(m)%axpy(dt*lev%s0mat(m,n),   lev%F(n,3))
       end do
       if (allocated(lev%tau)) then
          call lev%S(m)%axpy(1.0_pfdp, lev%tau(m))
       end if
    end do

    ! do the time-stepping
    call lev%Q(1)%copy(lev%q0)

    call misdc_evaluate(this, lev, t, 1)
    call this%f_eval(lev%Q(1), t0, lev%index, lev%F(1,1),1)
    call this%f_eval(lev%Q(1), t0, lev%index, lev%F(1,2),2)
    call this%f_eval(lev%Q(1), t0, lev%index, lev%F(1,3),3)

    call lev%ulevel%factory%create_single(rhs, lev%index,   lev%shape)

    t = t0
    dtsdc = dt * (lev%nodes(2:lev%nnodes) - lev%nodes(1:lev%nnodes-1))
    do m = 1, lev%nnodes-1
       t = t + dtsdc(m)

       call rhs%copy(lev%Q(m))
       call rhs%axpy(dtsdc(m), lev%F(m,1))
       call rhs%axpy(1.0_pfdp, lev%S(m))

       call this%f_comp(lev%Q(m+1), t, dtsdc(m), rhs, lev%index, lev%F(m+1,2),2)

       !  Now we need to do the final subtraction for the f3 piece
       call rhs%copy(Lev%Q(m+1))       
       call rhs%axpy(-1.0_pfdp*dtsdc(m), lev%F(m+1,3))

       call this%f_comp(lev%Q(m+1), t, dtsdc(m), rhs, lev%index, lev%F(m+1,3),3)
       call this%f_eval(lev%Q(m+1), t, lev%index, lev%F(m+1,1),1)
       call this%f_eval(lev%Q(m+1), t, lev%index, lev%F(m+1,2),2)
    end do
                         
    call lev%qend%copy(lev%Q(lev%nnodes))

    ! done
    call lev%ulevel%factory%destroy_single(rhs, lev%index,  lev%shape)

    call end_timer(pf, TLEVEL+lev%index-1)

  end subroutine misdc_sweep