!>\file seas_ngac_mod.F90 !! This file contains the ngac sea-salt module. !------------------------------------------------------------------------- ! NASA/GSFC, Data Assimilation Office, Code 910.3, GEOS/DAS ! ! Adapted by NOAA/GSD/ESRL ! !------------------------------------------------------------------------- !BOP ! ! !MODULE: seas_ngac_mod.F90 --- Calculate the Seasalt Emissions ! ! !INTERFACE: ! module seas_ngac_mod ! !USES: ! use chem_comm_mod, only : chem_comm_isroot use machine , only : kind_phys implicit none ! !PUBLIC TYPES: ! PRIVATE ! ! !PUBLIC MEMBER FUNCTIONS: ! PUBLIC SeasaltEmission ! !CONSTANTS real(kind=kind_phys), parameter :: r80fac = 1.65 ! ratio of radius(RH=0.8)/radius(RH=0.) [Gerber] real(kind=kind_phys), parameter :: rhop = 2200. ! dry seasalt density [kg m-3] ! ! !DESCRIPTION: ! ! This module implements the sea salt aerosol emission parameterizations. ! For all variants, emissions are some function of wind speed (and possibly ! other dynamical parameters) and the sea salt particle radius. Here, ! we assume the model passes in dry radius (or dry radius of size bin edges). ! Output is the mass emission flux (kg m-2 s-1) into that radius bin. ! ! !REVISION HISTORY: ! ! 30Mar2010 Colarco First crack! ! !EOP !------------------------------------------------------------------------- CONTAINS ! ! !IROUTINE: SeasaltEmission - Master driver to compute the sea salt emissions ! ! !INTERFACE: ! subroutine SeasaltEmission ( rLow, rUp, method, w10m, ustar, pi, & memissions, nemissions, rc ) ! !DESCRIPTION: Calculates the seasalt mass emission flux every timestep. ! The particular method (algorithm) used for the calculation is based ! on the value of "method" passed on input. Mostly these algorithms are ! a function of wind speed and particle size (nominally at 80% RH). ! Routine is called once for each size bin, passing in the edge radii ! "rLow" and "rUp" (in dry radius, units of um). Returned in the emission ! mass flux [kg m-2 s-1]. A sub-bin assumption is made to break (possibly) ! large size bins into a smaller space. ! ! !USES: implicit NONE ! !INPUT PARAMETERS: real(kind=kind_phys), intent(in) :: rLow, rUp ! Dry particle bin edge radii [um] real(kind=kind_phys), intent(in) :: w10m ! 10-m wind speed [m s-1] real(kind=kind_phys), intent(in) :: ustar ! friction velocity [m s-1] real(kind=kind_phys), intent(in) :: pi ! ratio of a circle's circumference to its diameter integer, intent(in) :: method ! Algorithm to use ! !OUTPUT PARAMETERS: real(kind=kind_phys), intent(inout) :: memissions ! Mass Emissions Flux [kg m-2 s-1] real(kind=kind_phys), intent(inout) :: nemissions ! Number Emissions Flux [# m-2 s-1] integer, intent(out) :: rc ! Error return code: ! 0 - all is well ! 1 - ! !Local Variables integer :: ir real(kind=kind_phys) :: w ! Intermediary wind speed [m s-1] real(kind=kind_phys) :: r, dr ! sub-bin radius spacing (dry, um) real(kind=kind_phys) :: rwet, drwet ! sub-bin radius spacing (rh=80%, um) real(kind=kind_phys) :: aFac, bFac, scalefac, rpow, exppow, wpow integer, parameter :: nr = 10 ! Number of (linear) sub-size bins character(len=*), parameter :: myname = 'SeasaltEmission' ! Define the sub-bins (still in dry radius) dr = (rUp - rLow)/nr r = rLow + 0.5*dr ! Loop over size bins nemissions = 0. memissions = 0. do ir = 1, nr rwet = r80fac * r drwet = r80fac * dr select case(method) case(1) ! Gong 2003 aFac = 4.7*(1.+30.*rwet)**(-0.017*rwet**(-1.44)) bFac = (0.433-log10(rwet))/0.433 scalefac = 1. rpow = 3.45 exppow = 1.607 wpow = 3.41 w = w10m case(2) ! Gong 1997 aFac = 3. bFac = (0.380-log10(rwet))/0.650 scalefac = 1. rpow = 1.05 exppow = 1.19 wpow = 3.41 w = w10m case(3) ! GEOS5 2012 aFac = 4.7*(1.+30.*rwet)**(-0.017*rwet**(-1.44)) bFac = (0.433-log10(rwet))/0.433 scalefac = 33.0e3 rpow = 3.45 exppow = 1.607 wpow = 3.41 - 1. w = ustar case default ! if(chem_comm_isroot()) print *, 'SeasaltEmission missing algorithm method' rc = 1 return end select ! Number emissions flux (# m-2 s-1) nemissions = nemissions + SeasaltEmissionGong( rwet, drwet, w, scalefac, aFac, bFac, rpow, exppow, wpow ) ! Mass emissions flux (kg m-2 s-1) scalefac = scalefac * 4./3.*pi*rhop*r**3.*1.e-18 memissions = memissions + SeasaltEmissionGong( rwet, drwet, w, scalefac, aFac, bFac, rpow, exppow, wpow ) r = r + dr end do rc = 0 end subroutine SeasaltEmission ! Function to compute sea salt emissions following the Gong style ! parameterization. Functional form is from Gong 2003: ! dN/dr = scalefac * 1.373 * (w^wpow) * (r^-aFac) * (1+0.057*r^rpow) * 10^(exppow*exp(-bFac^2)) ! where r is the particle radius at 80% RH, dr is the size bin width at 80% RH, and w is the wind speed function SeasaltEmissionGong ( r, dr, w, scalefac, aFac, bFac, rpow, exppow, wpow ) real(kind=kind_phys), intent(in) :: r, dr ! Wet particle radius, bin width [um] real(kind=kind_phys), intent(in) :: w ! Grid box mean wind speed [m s-1] (10-m or ustar wind) real(kind=kind_phys), intent(in) :: scalefac, aFac, bFac, rpow, exppow, wpow real(kind=kind_phys) :: SeasaltEmissionGong ! Initialize SeasaltEmissionGong = 0. ! Particle size distribution function SeasaltEmissionGong = scalefac * 1.373*r**(-aFac)*(1.+0.057*r**rpow) & *10**(exppow*exp(-bFac**2.))*dr ! Apply wind speed function SeasaltEmissionGong = w**wpow * SeasaltEmissionGong end function SeasaltEmissionGong end module seas_ngac_mod