Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

New updates on MEGANv2.1 #2588

Merged
merged 12 commits into from
Aug 12, 2024
Merged

New updates on MEGANv2.1 #2588

Show file tree
Hide file tree
Changes from 9 commits
Commits
Show all changes
12 commits
Select commit Hold shift + click to select a range
086e218
Updates from Arctic isoprene research
HuiWangWanderInGitHub Jun 4, 2024
a54623b
Add the drought part
HuiWangWanderInGitHub Jun 5, 2024
597ccd2
Merge tag 'ctsm5.2.008' into CLM5_MEGAN_ISOP
slevis-lmwg Jul 6, 2024
d940ee0
Add ._r8 to hardwired constants
slevis-lmwg Jul 6, 2024
237f848
Revert "Add ._r8 to hardwired constants"
slevis-lmwg Jul 10, 2024
c502523
Merge tag 'ctsm5.2.012' into CLM5_MEGAN_ISOP
slevis-lmwg Jul 16, 2024
7c2f596
Revisions according to review, plus some clean-up
slevis-lmwg Jul 18, 2024
0f147f8
Name hardwired coefficients into named params
slevis-lmwg Aug 1, 2024
793678f
Merge tag 'ctsm5.2.017' into CLM5_MEGAN_ISOP
slevis-lmwg Aug 1, 2024
c5fb0bb
Documentation for two variables
slevis-lmwg Aug 1, 2024
9dfc55a
Draft ChangeLog/ChangeSum
slevis-lmwg Aug 1, 2024
cc6e364
Merge remote-tracking branch 'escomp/master' into CLM5_MEGAN_ISOP
slevis-lmwg Aug 12, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
149 changes: 64 additions & 85 deletions src/biogeochem/VOCEmissionMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,7 @@ module VOCEmissionMod
use SolarAbsorbedType , only : solarabs_type
use TemperatureType , only : temperature_type
use PatchType , only : patch
use EnergyFluxType , only : energyflux_type
!
implicit none
private
Expand Down Expand Up @@ -378,7 +379,7 @@ end subroutine InitCold
!-----------------------------------------------------------------------
subroutine VOCEmission (bounds, num_soilp, filter_soilp, &
atm2lnd_inst, canopystate_inst, photosyns_inst, temperature_inst, &
vocemis_inst)
vocemis_inst, energyflux_inst)
!
! ! NEW DESCRIPTION
! Volatile organic compound emission
Expand Down Expand Up @@ -422,6 +423,7 @@ subroutine VOCEmission (bounds, num_soilp, filter_soilp, &
type(photosyns_type) , intent(in) :: photosyns_inst
type(temperature_type) , intent(in) :: temperature_inst
type(vocemis_type) , intent(inout) :: vocemis_inst
type(energyflux_type) , intent(in) :: energyflux_inst
!
! !REVISION HISTORY:
! 4/29/11: Colette L. Heald: expand MEGAN to 20 compound classes
Expand Down Expand Up @@ -476,14 +478,7 @@ subroutine VOCEmission (bounds, num_soilp, filter_soilp, &
end if

associate( &
!dz => col%dz , & ! Input: [real(r8) (:,:) ] depth of layer (m)
!bsw => soilstate_inst%bsw_col , & ! Input: [real(r8) (:,:) ] Clapp and Hornberger "b" (nlevgrnd)
!clayfrac => soilstate_inst%clayfrac_col , & ! Input: [real(r8) (:) ] fraction of soil that is clay
!sandfrac => soilstate_inst%sandfrac_col , & ! Input: [real(r8) (:) ] fraction of soil that is sand
!watsat => soilstate_inst%watsat_col , & ! Input: [real(r8) (:,:) ] volumetric soil water at saturation (porosity) (nlevgrnd)
!sucsat => soilstate_inst%sucsat_col , & ! Input: [real(r8) (:,:) ] minimum soil suction (mm) (nlevgrnd)
!h2osoi_vol => waterstate_inst%h2osoi_vol_col , & ! Input: [real(r8) (:,:) ] volumetric soil water (m3/m3)
!h2osoi_ice => waterstate_inst%h2osoi_ice_col , & ! Input: [real(r8) (:,:) ] ice soil content (kg/m3)
btran => energyflux_inst%btran_patch , & ! Input: [real(r8) (:) ] transpiration wetness factor (0 to 1)

forc_solad => atm2lnd_inst%forc_solad_downscaled_col, & ! Input: [real(r8) (:,:) ] direct beam radiation (visible only)
forc_solai => atm2lnd_inst%forc_solai_grc , & ! Input: [real(r8) (:,:) ] diffuse radiation (visible only)
Expand Down Expand Up @@ -569,10 +564,8 @@ subroutine VOCEmission (bounds, num_soilp, filter_soilp, &
! Activity factor for LAI (Guenther et al., 2006): all species
gamma_l = get_gamma_L(fsun240(p), elai(p))

! Activity factor for soil moisture: all species (commented out for now)
! gamma_sm = get_gamma_SM(clayfrac(p), sandfrac(p), h2osoi_vol(c,:), h2osoi_ice(c,:), &
! col%dz(c,:), soilstate_inst%bsw_col(c,:), watsat(c,:), sucsat(c,:), root_depth(patch%itype(p)))
gamma_sm = 1.0_r8
! Impact of soil moisture on isoprene emission
gamma_sm = get_gamma_SM(btran(p))

! Loop through VOCs for light, temperature and leaf age activity factor & apply
! all final activity factors to baseline emission factors
Expand All @@ -599,7 +592,7 @@ subroutine VOCEmission (bounds, num_soilp, filter_soilp, &

! Activity factor for T
gamma_t = get_gamma_T(t_veg240(p), t_veg24(p),t_veg(p), ct1(class_num), ct2(class_num),&
betaT(class_num),LDF(class_num), Ceo(class_num), Eopt, topt)
betaT(class_num),LDF(class_num), Ceo(class_num), Eopt, topt, patch%itype(p))

! Activity factor for Leaf Age
gamma_a = get_gamma_A(patch%itype(p), elai240(p),elai(p),class_num)
Expand Down Expand Up @@ -818,90 +811,53 @@ function get_gamma_L(fsun240_in,elai_in)
end function get_gamma_L

!-----------------------------------------------------------------------
function get_gamma_SM(clayfrac_in, sandfrac_in, h2osoi_vol_in, h2osoi_ice_in, dz_in, &
bsw_in, watsat_in, sucsat_in, root_depth_in)
!
! Activity factor for soil moisture (Guenther et al., 2006): all species
!----------------------------------
! Calculate the mean scaling factor throughout the root depth.
! wilting point potential is in units of matric potential (mm)
! (1 J/Kg = 0.001 MPa, approx = 0.1 m)
! convert to volumetric soil water using equation 7.118 of the CLM4 Technical Note
!
! !USES:
use clm_varcon , only : denice
use clm_varpar , only : nlevsoi
!
! !ARGUMENTS:
implicit none
real(r8),intent(in) :: clayfrac_in
real(r8),intent(in) :: sandfrac_in
real(r8),intent(in) :: h2osoi_vol_in(nlevsoi)
real(r8),intent(in) :: h2osoi_ice_in(nlevsoi)
real(r8),intent(in) :: dz_in(nlevsoi)
real(r8),intent(in) :: bsw_in(nlevsoi)
real(r8),intent(in) :: watsat_in(nlevsoi)
real(r8),intent(in) :: sucsat_in(nlevsoi)
real(r8),intent(in) :: root_depth_in
!
! !LOCAL VARIABLES:
real(r8) :: get_gamma_SM
integer :: j
real(r8) :: nl ! temporary number of soil levels
real(r8) :: theta_ice ! water content in ice in m3/m3
real(r8) :: wilt ! wilting point in m3/m3
real(r8) :: theta1 ! temporary
real(r8), parameter :: deltheta1=0.06_r8 ! empirical coefficient
real(r8), parameter :: smpmax = 2.57e5_r8 ! maximum soil matrix potential
!-----------------------------------------------------------------------
function get_gamma_SM(btran_in)

if ((clayfrac_in > 0) .and. (sandfrac_in > 0)) then
get_gamma_SM = 0._r8
nl=0._r8

do j = 1,nlevsoi
if (sum(dz_in(1:j)) < root_depth_in) then
theta_ice = h2osoi_ice_in(j)/(dz_in(j)*denice)
wilt = ((smpmax/sucsat_in(j))**(-1._r8/bsw_in(j))) * (watsat_in(j) - theta_ice)
theta1 = wilt + deltheta1
if (h2osoi_vol_in(j) >= theta1) then
get_gamma_SM = get_gamma_SM + 1._r8
else if ( (h2osoi_vol_in(j) > wilt) .and. (h2osoi_vol_in(j) < theta1) ) then
get_gamma_SM = get_gamma_SM + ( h2osoi_vol_in(j) - wilt ) / deltheta1
else
get_gamma_SM = get_gamma_SM + 0._r8
end if
nl=nl+1._r8
end if
end do

if (nl > 0._r8) then
get_gamma_SM = get_gamma_SM/nl
endif
!---------------------------------------
! May 22, 2024
! Activity factor for soil moisture of Isoprene (Wang et al., 2022, JAMES)
! It is based on eq. (11) in the paper. Because the temperature response
! of isoprene has been explicitly included in CLM;

if (get_gamma_SM > 1.0_r8) then
write(iulog,*) 'healdSM > 1: gamma_SM, nl', get_gamma_SM, nl
get_gamma_SM=1.0_r8
endif
!ARGUMENTS:
implicit none
real(r8),intent(in) :: btran_in

!!!------- the drought algorithm--------
real(r8), parameter :: a1 = -7.4463_r8
real(r8), parameter :: b1 = 3.2552_r8
real(r8), parameter :: btran_threshold = 0.2_r8
real(r8) :: get_gamma_SM
!---------------------------------------

if (btran_in >= 1._r8) then
get_gamma_SM = 1._r8
else
get_gamma_SM = 1.0_r8
end if
get_gamma_SM = 1._r8 / (1._r8 + b1 * exp(a1 * (btran_in - btran_threshold)))
endif

end function get_gamma_SM

!-----------------------------------------------------------------------
function get_gamma_T(t_veg240_in, t_veg24_in,t_veg_in, ct1_in, ct2_in, betaT_in, LDF_in, Ceo_in, Eopt, topt)
function get_gamma_T(t_veg240_in, t_veg24_in,t_veg_in, ct1_in, ct2_in, betaT_in, LDF_in, Ceo_in, Eopt, topt, ivt_in)

! Activity factor for temperature
! Activity factor for temperature
!--------------------------------
! May 24, 2024 Hui updated the temperature response curves of isoprene for
! Boreal Broadleaf Deciduous Shrub and Arctic C3 grass based on
! Wang et al., 2024 (GRL) and Wang et al., 2024 (Nature Communications)
!--------------------------------
! Calculate both a light-dependent fraction as in Guenther et al., 2006 for isoprene
! of a max saturation type form. Also caculate a light-independent fraction of the
! form of an exponential. Final activity factor depends on light dependent fraction
! of compound type.
!
! !USES:
use clm_varcon, only: tfrz

! !ARGUMENTS:
implicit none
integer,intent(in) :: ivt_in
real(r8),intent(in) :: t_veg240_in
real(r8),intent(in) :: t_veg24_in
real(r8),intent(in) :: t_veg_in
Expand All @@ -917,7 +873,9 @@ function get_gamma_T(t_veg240_in, t_veg24_in,t_veg_in, ct1_in, ct2_in, betaT_in,
real(r8) :: get_gamma_T
real(r8) :: gamma_t_LDF ! activity factor for temperature
real(r8) :: gamma_t_LIF ! activity factor for temperature
real(r8) :: x ! temporary
real(r8) :: x ! temporary i
real(r8) :: bet_arc_c3
real(r8), parameter :: bet_arc_c3_max = 300._r8
real(r8), parameter :: co1 = 313._r8 ! empirical coefficient
real(r8), parameter :: co2 = 0.6_r8 ! empirical coefficient
real(r8), parameter :: co4 = 0.05_r8 ! empirical coefficient
Expand All @@ -927,19 +885,40 @@ function get_gamma_T(t_veg240_in, t_veg24_in,t_veg_in, ct1_in, ct2_in, betaT_in,
real(r8), parameter :: ct3 = 0.00831_r8 ! empirical coefficient (0.0083 in User's Guide)
real(r8), parameter :: tstd = 303.15_r8 ! std temperature [K]
real(r8), parameter :: bet = 0.09_r8 ! beta empirical coefficient [K-1]
real(r8), parameter :: std_act_energy_isopr = 95._r8 ! standard activation energy for isoprene
real(r8), parameter :: empirical_param_1 = 9.49_r8 ! empirical param for the activation energy in response to 10-day temperature change
real(r8), parameter :: empirical_param_2 = 0.53_r8 ! empirical param for the activation energy in response to 10-day temperature change
real(r8), parameter :: empirical_param_3 = 0.12_r8 ! empirical param for the emission factors of arctic C3 grass in response to 10-day temperature change
real(r8), parameter :: empirical_param_4 = 7.9_r8 ! empirical param for the emission factors of broadleaf deciduous boreal shrubs in response to 10-day temperature change
real(r8), parameter :: empirical_param_5 = 0.217_r8 ! empirical param for the emission factors of broadleaf deciduous boreal shrubs in response to 10-day temperature change
!-----------------------------------------------------------------------

! Light dependent fraction (Guenther et al., 2006)
if ( (t_veg240_in > 0.0_r8) .and. (t_veg240_in < 1.e30_r8) ) then
! topt and Eopt from eq 8 and 9:
topt = co1 + (co2 * (t_veg240_in-tstd0))
Eopt = Ceo_in * exp (co4 * (t_veg24_in-tstd0)) * exp(co4 * (t_veg240_in -tstd0))
if ( (ivt_in == nbrdlf_dcd_brl_shrub) ) then ! boreal-deciduous-shrub
! coming from BEAR-oNS campaign willows results
Eopt = empirical_param_4 * exp (empirical_param_5 * (t_veg24_in - tfrz - 24.0_r8))
else if ( (ivt_in == nc3_arctic_grass ) ) then ! boreal-grass
Eopt = exp(empirical_param_3 * (t_veg240_in - tfrz - 15._r8))
else
Eopt = Ceo_in * exp (co4 * (t_veg24_in - tstd0)) * exp(co4 * (t_veg240_in - tstd0))
endif

else
topt = topt_fix
Eopt = Eopt_fix
endif
x = ( (1._r8/topt) - (1._r8/(t_veg_in)) ) / ct3
gamma_t_LDF = Eopt * ( ct2_in * exp(ct1_in * x)/(ct2_in - ct1_in * (1._r8 - exp(ct2_in * x))) )
! for the boreal grass from BEAR-oNS campaign
if ( (ivt_in == nc3_arctic_grass ) ) then ! boreal-grass
bet_arc_c3 = std_act_energy_isopr + empirical_param_1 * exp(empirical_param_2 * (tfrz + 15._r8 - t_veg240_in))
bet_arc_c3 = min(bet_arc_c3, bet_arc_c3_max)
gamma_t_LDF = Eopt * exp(bet_arc_c3 * ((1._r8 / (tfrz + 30._r8) - 1._r8 / t_veg_in) / ct3))
else
gamma_t_LDF = Eopt * ( ct2_in * exp(ct1_in * x) / (ct2_in - ct1_in * (1._r8 - exp(ct2_in * x))) )
endif


! Light independent fraction (of exp(beta T) form)
Expand Down
4 changes: 2 additions & 2 deletions src/main/clm_driver.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -815,11 +815,11 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro
call dust_emis_inst%DustDryDep(bounds_clump, &
atm2lnd_inst, frictionvel_inst)

! VOC emission (A. Guenther's MEGAN (2006) model)
! VOC emission (A. Guenther's MEGAN (2006) model; Wang et al., 2022, 2024a, 2024b)
call VOCEmission(bounds_clump, &
filter(nc)%num_soilp, filter(nc)%soilp, &
atm2lnd_inst, canopystate_inst, photosyns_inst, temperature_inst, &
vocemis_inst)
vocemis_inst, energyflux_inst)

call t_stopf('bgc')

Expand Down
Loading