advect Subroutine

  subroutine advect(idimmin,idimmax)
    use mod_global_parameters
    use mod_fix_conserve
    use mod_ghostcells_update
    use mod_physics, only: phys_req_diagonal
    use mod_comm_lib, only: mpistop
    use mod_physicaldata

    integer, intent(in) :: idimmin,idimmax
    integer             :: iigrid, igrid, ix1,ix2,ix3, iw

    call init_comm_fix_conserve(idimmin,idimmax,nwflux)
    fix_conserve_at_step = time_advance .and. levmax>levmin
    
    ! copy w instead of wold because of potential use of dimsplit or sourcesplit
     
    do iigrid=1,igridstail; igrid=igrids(iigrid);
       
       do iw = 1, nw
          do ix3 = ixGlo3, ixGhi3 
             do ix2 = ixGlo2, ixGhi2 
                do ix1 = ixGlo1, ixGhi1 
                   bg(2)%w(ix1,ix2,ix3,iw,igrid) = bg(1)%w(ix1,ix2,ix3,iw,igrid)
                end do
             end do
          end do
       end do
    end do
    
    if(stagger_grid) then
        
       do iigrid=1,igridstail; igrid=igrids(iigrid);
          
          do iw = 1, nws
             do ix3 = ps(igrid)%ixGsmin3, ps(igrid)%ixGsmax3 
                do ix2 = ps(igrid)%ixGsmin2, ps(igrid)%ixGsmax2 
                   do ix1 = ps(igrid)%ixGsmin1, ps(igrid)%ixGsmax1 
                      ps1(igrid)%ws(ix1,ix2,ix3,iw) = ps(igrid)%ws(ix1,ix2,ix3,iw)
                   end do
                end do
             end do
          end do
       end do
    end if

 istep = 0

 select case (t_stepper)
  
    case (threestep)
       select case (t_integrator)
       ! AGILE this is our integrator (default threestep)
       case (ssprk3)
          ! this is SSPRK(3,3) Gottlieb-Shu 1998 or SSP(3,2) depending on ssprk_order (3 vs 2)
          
          call advect1(flux_method,rk_beta11, idimmin,idimmax,global_time,ps,&
               bg(1),global_time,ps1,bg(2))

           
          do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
             
             do iw = 1, nw
                do ix3 = ixGlo3, ixGhi3 
                   do ix2 = ixGlo2, ixGhi2 
                      do ix1 = ixGlo1, ixGhi1 
                         bg(3)%w(ix1,ix2,ix3,iw,igrid) = rk_alfa21 * bg(1)%w(ix1,ix2,&
                              ix3,iw,igrid) + rk_alfa22 * bg(2)%w(ix1,ix2,ix3,iw,igrid)
                      end do
                   end do
                end do
             end do
             if(stagger_grid) ps2(igrid)%ws=rk_alfa21*ps(igrid)%ws+&
                  rk_alfa22*ps1(igrid)%ws
          end do

          call advect1(flux_method,rk_beta22, idimmin,idimmax,&
               global_time+rk_c2*dt,ps1,bg(2),global_time+rk_alfa22*rk_c2*dt,ps2,&
               bg(3))

           
          do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
             
             do iw = 1, nw
                do ix3 = ixGlo3, ixGhi3 
                   do ix2 = ixGlo2, ixGhi2 
                      do ix1 = ixGlo1, ixGhi1 
                         bg(1)%w(ix1,ix2,ix3,iw,igrid) = rk_alfa31 * bg(1)%w(ix1,ix2,&
                              ix3,iw,igrid) + rk_alfa33 * bg(3)%w(ix1,ix2,ix3,iw,igrid)
                      end do
                   end do
                end do
             end do
             if(stagger_grid) ps(igrid)%ws=rk_alfa31*ps(igrid)%ws+&
                  rk_alfa33*ps2(igrid)%ws
          end do

          call advect1(flux_method,rk_beta33, idimmin,idimmax,&
               global_time+rk_c3*dt,ps2,bg(3),global_time+(1.0d0-rk_beta33)*dt,&
               ps,bg(1))

       case default
          call mpistop("unkown threestep time_integrator in advect")
       end select

    case default
       call mpistop("unkown time_stepper in advect")
    end select

  end subroutine advect