CalculateTendency_DGModel3D Subroutine

  subroutine CalculateTendency_DGModel3D(this)
    implicit none
    class(DGModel3D),intent(inout) :: this

    call this%solution%BoundaryInterp()

    ! Post the halo exchange for the prognostic variables; static variables
    ! (nstepped+1:nvar) are carried in full by the first exchange only, and
    ! their boundary traces do not change thereafter. The MPI messages are
    ! in flight while the hooks, boundary conditions, source, and interior
    ! flux evaluations below execute.
    call this%solution%SideExchangeStart(this%mesh,this%nstepped)

    ! The hooks and methods between SideExchangeStart and SideExchangeFinish
    ! must not read extBoundary on interior (rank-shared) faces.
    ! SetBoundaryCondition writes extBoundary only on physical-boundary
    ! faces, which are disjoint from the faces the exchange fills.
    call this%PreTendencyHook() ! User-supplied
    call this%SetBoundaryCondition() ! User-supplied

    if(this%gradient_enabled) then
      ! The BR gradient consumes extBoundary (through the side averages), so
      ! the exchange must complete before the gradient is computed.
      call this%solution%SideExchangeFinish(this%mesh)
      call this%CalculateSolutionGradient()
      call this%SetGradientBoundaryCondition() ! User-supplied
      call this%solutionGradient%AverageSides()
      call this%SourceMethod() ! User supplied
      call this%BoundaryFlux() ! User supplied
      call this%FluxMethod() ! User supplied
    else
      ! Interior work overlaps with the halo exchange; BoundaryFlux is the
      ! first consumer of extBoundary and runs after the exchange completes.
      ! FluxMethod and BoundaryFlux write disjoint outputs (flux interior
      ! vs. boundarynormal), so this reordering does not change any
      ! floating-point results.
      call this%SourceMethod() ! User supplied
      call this%FluxMethod() ! User supplied
      call this%solution%SideExchangeFinish(this%mesh)
      call this%BoundaryFlux() ! User supplied
    endif

    call this%flux%MappedDGDivergence(this%fluxDivergence%interior_gpu)

    ! The tendency (source - fluxDivergence) is formed inside the fused RK/Euler
    ! update kernels (UpdateGRK_CalculateDSDt_gpu / UpdateSolution_CalculateDSDt_gpu),
    ! so there is no separate CalculateDSDt pass here.

  endsubroutine CalculateTendency_DGModel3D