module sulemis 1,3
!-----------------------------------------------------------------------
! Purpose:
! Add the sulfur emissions to the tracer array. This is done here rather
! than in the vertical diffusions routines because the emissions are added
! to more than just the surface layer (stack emissions above 100 m).
!
! Author: module coded by B. Eaton.
!-----------------------------------------------------------------------
use shr_kind_mod
, only: r8 => shr_kind_r8
use pmgrid
use ppgrid
implicit none
save
private
public :: &
addsulemis ! Add emissions for SO2, SO4, and DMS to the tracer array.
!##############################################################################
contains
!##############################################################################
subroutine addsulemis( lchnk, ncol, lat, lon, dtime, gravit, rpdel, zi, chtr) 1,8
!-----------------------------------------------------------------------
! Purpose:
! Add emissions for SO2, SO4, and DMS to the tracer array.
!
! The emissions are specified from data provided by Benkovitz.
! SO2 and SO4 are emitted at 2 levels (below and above 100m).
! DMS is emitted from the surface only.
!
! Author: M. Barth
!-----------------------------------------------------------------------
use scyc
, only: useGEIA
use dmsbnd
, only: dmsbndget
use soxbnd
, only: soxbndget
use sulbnd
, only: sulbndget
#ifdef MATCH
use histout, only: outfld
#else
use history
, only: outfld
#endif
#ifdef SCYC_MASSBGT
use massbgt, only: aal
#endif
implicit none
integer, intent(in) :: lchnk ! chunk
integer, intent(in) :: ncol ! number of columns
integer, intent(in) :: lat(pcols) ! model latitude index
integer, intent(in) :: lon(pcols) ! model longitude index
real(r8), intent(in) ::&
dtime, &
gravit, &
rpdel(pcols,pver), &! reciprocal of pdel
zi(pcols,pver)
real(r8), intent(inout) ::&
chtr(pcols,pver,3) ! chemical tracer array (SO2,SO4,DMS)
! obuf(*) ! history output buffer
! Local variables:
integer :: i, m, il, kp
real(r8) ::&
scl, &! scale factor to put flux in correct units
sflx(pcols,6), &! sulfur emissions for a given latitude
! 1,4 SO2 emissions at 0m, >100m
! 2,5 SO4 emissions at 0m, >100m
! 3,6 DMS emissions at 0m, >100m
tmp(pcols,2)
!-----------------------------------------------------------------------
! Get latitude slice of emissions data.
if ( useGEIA() ) then
call sulbndget
( pcols, ncol, lat, lon, sflx )
else
call soxbndget
( ncol, lchnk, tmp )
call dmsbndget
( ncol, lchnk, sflx(1,3) )
do i = 1, ncol
! Split SOx emissions, 98% --> SO2, 2% --> SO4
sflx(i,1) = .98*tmp(i,1) * 2. ! *2. converts kg S --> kg SO2
sflx(i,4) = .98*tmp(i,2) * 2.
sflx(i,2) = .02*tmp(i,1) * 3. ! *3. converts kg S --> kg SO4
sflx(i,5) = .02*tmp(i,2) * 3.
end do
end if
! Put fluxes on history file.
do i = 1, ncol
tmp(i,1) = sflx(i,1) + sflx(i,4)
tmp(i,2) = sflx(i,2) + sflx(i,5)
end do
#ifdef SCYC_MASSBGT
call aal( 'so2sf', lat, tmp(:,1) )
call aal( 'so4sf', lat, tmp(:,2) )
call aal( 'dmssf', lat, sflx(:,3) )
#endif
! Add surface flux to lowest tracer level.
do il = 1, 3
if( il .eq. 1 ) m = 1
if( il .eq. 2 ) m = 2
if( il .eq. 3 ) m = 3
do i = 1, ncol
scl = gravit*rpdel(i,pver) ! scales emissions from kg/m2-s to kg/kg-s
chtr(i,pver,m) = chtr(i,pver,m) + dtime*sflx(i,il)*scl
end do
end do
! Add "above 100 m emissions" to bottom and next to bottom levels as
! appropriate. (SO2 and SO4 only)
do il = 4, 5
if( il .eq. 4 ) m = 1
if( il .eq. 5 ) m = 2
do i = 1, ncol
if( 100. .lt. zi(i,pver) ) then
scl = gravit*rpdel(i,pver)
chtr(i,pver,m) = chtr(i,pver,m) + &
dtime*sflx(i,il)*scl * (zi(i,pver)-100.)/ &
(zi(i,pver-1)-100.)
scl = gravit*rpdel(i,pver-1)
chtr(i,pver-1,m) = chtr(i,pver-1,m) + &
dtime*sflx(i,il)*scl * (zi(i,pver-1)-zi(i,pver))/ &
(zi(i,pver-1)-100.)
else
scl = gravit*rpdel(i,pver-1)
chtr(i,pver-1,m) = chtr(i,pver-1,m) +dtime*sflx(i,il)*scl
endif
end do
end do
end subroutine addsulemis
!#######################################################################
end module sulemis