MODULE trcdms_medusa !!====================================================================== !! *** MODULE trcdms_medusa *** !! TOP : MEDUSA !!====================================================================== !! History : !! - ! 2014-08 (J. Palmieri - A. Yool) added for UKESM1 project !! - ! 2017-05 (A. Yool) add extra Anderson scheme !! - ! 2018-10 (A. Yool) Add air-sea DMS flux !!---------------------------------------------------------------------- #if defined key_medusa && defined key_roam !!---------------------------------------------------------------------- !! MEDUSA DMS surface concentration !!---------------------------------------------------------------------- !! trc_dms_medusa : !!---------------------------------------------------------------------- USE oce_trc USE trc USE sms_medusa USE lbclnk USE prtctl_trc ! Print control for debugging USE in_out_manager ! I/O manager USE yomhook, ONLY: lhook, dr_hook USE parkind1, ONLY: jprb, jpim IMPLICIT NONE PRIVATE PUBLIC trc_dms_medusa ! called in trc_bio_medusa PUBLIC dms_flux_ocn ! called in air_sea !!* Substitution # include "domzgr_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/TOP 2.0 , LOCEAN-IPSL (2007) !! $Id$ !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) !!---------------------------------------------------------------------- CONTAINS !======================================================================= ! SUBROUTINE trc_dms_medusa( chn, chd, mld, xqsr, xdin, xlim, & !! inputs & dms_andr, dms_simo, dms_aran, dms_hall, dms_andm) !! outputs ! !======================================================================= !! !! Title : Calculates DMS ocean surface concentration !! Author : Julien Palmieri and Andrew Yool !! Date : 08/08/14 !! !! DMS module is called in trc_bio's huge jk,jj,ji loop !! --> DMS concentration is calculated in a specific cell !! (no need of ji,jj,jk) !! !! AXY (13/03/15): amend to include all four schemes tested !! during winter/spring 2015; these are: !! !! 1. Anderson et al. (2001); this uses fields !! of surface chl, irradiance and nutrients !! to empirically estimate DMS via a broken !! stick approach !! !! 2. Simo & Dachs (2002); this uses fields of !! surface chl and mixed layer depth !! !! 3. Aranami & Tsunogai (2004); this is an !! embellishment of Simo & Dachs !! !! 4. Halloran et al. (2010); this is an !! alternative embellishment of Sim & Dachs !! and is included because it is formally !! published (and different from the above) !! !! AXY (25/05/17): add extra "corrected" Anderson scheme !! !! 5. As Anderson et al. (2001) but modified to !! more accurately reflect nutrient limitation !! status of phytoplankton community !! !! AXY (08/07/15): amend to remove Julien's original calculation !! as this is now superfluous; the four schemes !! are calculated and one is chosen to be passed !! to the atmosphere in trc_bio_medusa !! !======================================================================= USE yomhook, ONLY: lhook, dr_hook USE parkind1, ONLY: jprb, jpim IMPLICIT NONE ! REAL(wp), INTENT( in ) :: chn !! non-diatom chlorophyll (mg/m3) REAL(wp), INTENT( in ) :: chd !! diatom chlorophyll (mg/m3) REAL(wp), INTENT( in ) :: mld !! mix layer depth (m) REAL(wp), INTENT( in ) :: xqsr !! surface irradiance (W/m2) REAL(wp), INTENT( in ) :: xdin !! surface DIN (mmol N/m3) REAL(wp), INTENT( in ) :: xlim !! surface DIN limitation (mmol N/m3) REAL(wp), INTENT( inout ) :: dms_andr !! DMS surface concentration (nmol/L) REAL(wp), INTENT( inout ) :: dms_simo !! DMS surface concentration (nmol/L) REAL(wp), INTENT( inout ) :: dms_aran !! DMS surface concentration (nmol/L) REAL(wp), INTENT( inout ) :: dms_hall !! DMS surface concentration (nmol/L) REAL(wp), INTENT( inout ) :: dms_andm !! DMS surface concentration (nmol/L) ! REAL(wp) :: CHL, cmr, sw_dms REAL(wp) :: Jterm, Qterm !! temporary variables REAL(wp) :: fq1,fq2,fq3 INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 REAL(KIND=jprb) :: zhook_handle CHARACTER(LEN=*), PARAMETER :: RoutineName='TRC_DMS_MEDUSA' IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle) ! !======================================================================= ! ! AXY (13/03/15): per remarks above, the following calculations estimate ! DMS using all of the schemes examined for UKESM1 ! CHL = 0.0 CHL = chn+chd !! mg/m3 cmr = CHL / mld ! ! AXY (13/03/15): Anderson et al. (2001) !! JPALM --19-12-2017-- Tunable through the namelist !! within dmsmin - dmscut - dmsslp Jterm = xqsr + 1.0e-6 !! this next line makes a hard-coded assumption about the !! half-saturation constant of MEDUSA (which should be !! done properly; perhaps even scaled with the proportion !! of diatoms and non-diatoms) Qterm = xdin / (xdin + 0.5) fq1 = log10(CHL * Jterm * Qterm) if (fq1 > dmscut) then dms_andr = (dmsslp * (fq1 - dmscut)) + dmsmin else dms_andr = dmsmin endif ! ! AXY (13/03/15): Simo & Dachs (2002) fq1 = (-1.0 * log(mld)) + 5.7 fq2 = (55.8 * cmr) + 0.6 if (cmr < 0.02) then dms_simo = fq1 else dms_simo = fq2 endif ! ! AXY (13/03/15): Aranami & Tsunogai (2004) fq1 = 60.0 / mld fq2 = (55.8 * cmr) + 0.6 if (cmr < 0.02) then dms_aran = fq1 else dms_aran = fq2 endif ! ! AXY (13/03/15): Halloran et al. (2010) fq1 = (-1.0 * log(mld)) + 5.7 fq2 = (55.8 * cmr) + 0.6 fq3 = (90.0 / mld) if (cmr < 0.02) then dms_hall = fq1 else dms_hall = fq2 endif if (mld > 182.5) then dms_hall = fq3 endif ! ! AXY (25/05/17): modified Anderson et al. (2001) Jterm = xqsr + 1.0e-6 !! this version fixes the hard-coded assumption above Qterm = xlim fq1 = log10(CHL * Jterm * Qterm) if (fq1 > 1.72) then dms_andm = (8.24 * (fq1 - 1.72)) + 2.29 else dms_andm = 2.29 endif IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle) END SUBROUTINE trc_dms_medusa !======================================================================= !======================================================================= !======================================================================= !======================================================================= ! SUBROUTINE dms_flux_ocn( wind_10m, tstar, dms_conc, i_dms_flux, & !! inputs & f_dms ) !! outputs ! !======================================================================= !! !! Title : Calculates DMS air-sea exchange !! Author : Andrew Yool, based on UKMO original code !! Date : 11/10/18 !! !! Air-sea DMS flux is normally calculated by the UM atmosphere, !! as part of its aerosols; however, the OMIP simulation for CMIP6 !! is ocean-only and does not include the UM; consequently, this !! code has been added to permit ocean-only UKESM1 to produce an !! air-sea DMS flux in addition to surface DMS which, hitherto, !! was all it would produce; code is largely copy-pasted from the !! UKMO original (UM code block is dms_flux_4A.F90); the code !! here is hard-wired to use single input values (i.e. not a 2D !! area) and make use of the Liss & Merlivat (1986) function !! !! This DMS function is called from air_sea.F90 !--------------------------------------------------------------------- ! Purpose: To calculate the flux of DMS (as kg m-2 s-1 of sulphur) ! from the ocean surface as a function of its concentration ! in seawater and of windspeed. The sea-air exchange can ! be determined according to one of three commonly-used ! parametrization schemes, those of Liss & Merlivat (1986), ! Wanninkhof (1992) or Nightingale et al. (2000). The routine ! is called by Aero_Ctl. ! ! Method: The Schmidt number ! for DMS is calculated as in Saltzman et al. (1993), and ! used with the windspeed to determine the mass transfer (or ! "piston") velocity according to the desired parametrization. ! This is then used to determine the sea-air mass flux of DMS ! as a function of sea-water DMS concentration. High surface ! temperatures (caused by the land portion of a gridbox when ! coastal tiling is not active) cause negative Sc values which ! would give a floating-point error in the k_DMS calculation, ! so the Tstar values are capped. This shouldn't be a problem ! when coastal tiling is on as then the Tstar values passed in ! are those for sea only. ! ! Code Owner: Please refer to the UM file CodeOwners.txt ! This file belongs in section: Aerosols ! ! Code Description: ! Language: Fortran 90 ! This code is written to UMDP3 v8 programming standards ! !--------------------------------------------------------------------- USE yomhook, ONLY: lhook, dr_hook USE parkind1, ONLY: jprb, jpim IMPLICIT NONE ! REAL(wp), INTENT( in ) :: wind_10m !! 10m wind (m/s) REAL(wp), INTENT( in ) :: tstar !! SST (degrees C) REAL(wp), INTENT( in ) :: dms_conc !! surface DMS (nmol / l) INTEGER, INTENT(in) :: i_dms_flux !! gas transfer choice ! REAL(wp), INTENT( inout ) :: f_dms !! DMS flux (kg S m-2 s-1) ! Local variables: REAL :: sc ! Schmidt number REAL :: k_dms ! Piston velocity of DMS (cm h-1) REAL :: t_c ! Surface temperature in degrees Celsius ! Piston velocities for gases with Schmidt numbers of 600 & 660 resp. (cm h-1) REAL :: k_600 REAL :: k_660 REAL :: n ! Schmidt number exponent REAL, PARAMETER :: t_max = 47.0 !! Max T to avoid breaking the Sc fit (C) INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 REAL(KIND=jprb) :: zhook_handle CHARACTER(LEN=*), PARAMETER :: RoutineName='DMS_FLUX_OCN' IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle) ! Calculate the Schmidt number (Sc): t_c = MIN(tstar, t_max) sc = 2674.0 - (147.12*t_c) + (3.726*t_c**2) & - (0.038*t_c**3) ! Determine the mass transfer (or "piston") velocity (k_DMS) over sea ! according to the specified parametrization scheme: if (i_dms_flux .eq. 1) then ! ---------------------------------------------------------------------- ! Liss & Merlivat (1986) IF (wind_10m .le. 3.6) THEN k_600 = 0.17 * wind_10m n = -2.0/3.0 ELSEIF ( wind_10m .gt. 3.6 .AND. wind_10m .le. 13.0 ) THEN k_600 = (2.85 * wind_10m) - 9.65 n = -0.5 ELSE k_600 = (5.90 * wind_10m) - 49.3 n = -0.5 END IF k_dms = k_600 * (sc / 600.0)**n elseif (i_dms_flux .eq. 2) then ! ---------------------------------------------------------------------- ! Wanninkhof (1992) k_660 = 0.31 * wind_10m**2 n = -0.5 k_dms = k_660 * (sc / 660.0)**n elseif (i_dms_flux .eq. 3) then ! ---------------------------------------------------------------------- ! Nightingale et al. (2000) k_600 = (0.222 * wind_10m**2) + (0.333 * wind_10m) n = -0.5 k_dms = k_600 * (sc / 600.0)**n else ! ---------------------------------------------------------------------- ! You shouldn't be here k_dms = 0.0 endif ! Finally, calculate the sea-air flux of DMS as a function of k_DMS ! and dissolved DMS concentration. The former requires a conversion ! from cm hour-1 to ms-1, and the latter from nanomoles per litre to ! kg[S] m-3, to return the flux in kg[S] m-2 sec-1. f_dms = (k_dms / 3.6e5) * (dms_conc * 32.0e-9) IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle) END SUBROUTINE dms_flux_ocn !======================================================================= !======================================================================= !======================================================================= #else !!====================================================================== !! Dummy module : No MEDUSA bio-model !!====================================================================== CONTAINS !======================================================================= ! SUBROUTINE trc_dms_medusa( kt ) !! EMPTY Routine ! ! INTEGER, INTENT( in ) :: kt INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0 INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1 REAL(KIND=jprb) :: zhook_handle CHARACTER(LEN=*), PARAMETER :: RoutineName='TRC_DMS_MEDUSA' IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle) ! WRITE(*,*) 'trc_dms_medusa: You should not have seen this print! error?' IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle) END SUBROUTINE trc_dms_medusa #endif !!====================================================================== END MODULE trcdms_medusa