subroutine add_split_source(prior)
use mod_global_parameters
use mod_ghostcells_update
use mod_physics, only: phys_req_diagonal, phys_global_source_after
use mod_comm_lib, only: mpistop
logical, intent(in) :: prior
real(dp) :: qt
real(dp) :: dr(ndim), xloc(ndim)
! Physical cell size handed to the source term. dr is the *logical*
! spacing, which is not a length in a curvilinear system and, under
! LOG_RADIUS, is not even proportional to one.
real(dp) :: dloc(ndim)
integer :: iigrid, n, ix1, ix2, ix3
integer :: ixOmin1,ixOmin2,ixOmin3,ixOmax1,ixOmax2,ixOmax3
real(dp) :: wprim(nw_phys), wnew(nw_phys), wCT(nw_phys)
logical :: src_active
src_active = .false.
if ((.not.prior).and.(sourcesplit==sourcesplit_sf .or. &
sourcesplit==sourcesplit_ssf)) return
if (prior) then
qt=global_time
else
qt=global_time+dt
end if
if(any_source_split) then
! add split source terms
ixOmin1=ixGlo1+nghostcells;ixOmin2=ixGlo2+nghostcells
ixOmin3=ixGlo3+nghostcells;ixOmax1=ixGhi1-nghostcells
ixOmax2=ixGhi2-nghostcells;ixOmax3=ixGhi3-nghostcells;
select case (sourcesplit)
case (sourcesplit_sfs)
do iigrid=1,igridstail_active
n = igrids_active(iigrid)
dr = rnode(rpdx1_:rnodehi, n)
do ix3=ixOmin3,ixOmax3
do ix2=ixOmin2,ixOmax2
do ix1=ixOmin1,ixOmax1
xloc(1:ndim) = bgeo%x(ix1, ix2, ix3, 1:ndim, n)
dloc = dr
wCT = bg(1)%w(ix1,ix2,ix3, 1:nw_phys, n)
wnew = wCT
wprim = wCT
call to_primitive( wprim )
call addsource_local(0.5d0*dt,&
0.5d0, qt, wCT,&
wprim, qt, wnew, xloc, dloc, .true. )
bg(1)%w(ix1,ix2,ix3, 1:nw_phys, n) = wnew(1:nw_phys)
! TBD
end do
end do
end do
end do
src_active = .true.
case default
write(unitterm,*)'No such type/not yet implemented sourcesplit=',sourcesplit
call mpistop("Error: Unknown type of sourcesplit!")
end select
if (.not. prior .and. associated(phys_global_source_after)) then
call phys_global_source_after(dt, qt, src_active)
end if
if (src_active) then
call getbc(qt,0.d0,ps,iwstart,nwgc,phys_req_diagonal)
end if
end if
end subroutine add_split_source