Merge pull request NCAR#73 from SamuelTrahanNOAA/feature/pert-scale

SamuelTrahanNOAA · web-flow · commit 6c2d4e81fb7f · 2021-01-15T12:17:15.000-05:00
gsd-chem stochastic emissions: scale factor, remove CA support, remove emis_multiplier
diff --git a/gsdchem/gsd_chem_anthropogenic_wrapper.F90 b/gsdchem/gsd_chem_anthropogenic_wrapper.F90
@@ -41,9 +41,8 @@ end subroutine gsd_chem_anthropogenic_wrapper_finalize
     subroutine gsd_chem_anthropogenic_wrapper_run(im, kte, kme, ktau, dt,               &
                    pr3d, ph3d,phl3d, prl3d, tk3d, spechum,emi_in,                       &
                    ntrac,ntso2,ntsulf,ntpp25,ntbc1,ntoc1,ntpp10,                        &
-                   gq0,qgrs,abem,chem_opt_in,kemit_in,                                  &
-                   emis_amp_anthro, emis_multiplier, do_sppt_emis, ca_global_emis,      &
-                   errmsg,errflg)
+                   gq0,qgrs,abem,chem_opt_in,kemit_in,pert_scale_anthro,                &
+                   emis_amp_anthro,do_sppt_emis,sppt_wts,errmsg,errflg)
 
     implicit none
 
@@ -52,10 +51,9 @@ subroutine gsd_chem_anthropogenic_wrapper_run(im, kte, kme, ktau, dt,
     integer,        intent(in) :: ntrac
     integer,        intent(in) :: ntso2,ntpp25,ntbc1,ntoc1,ntpp10
     integer,        intent(in) :: ntsulf
-    real(kind_phys),intent(in) :: dt, emis_amp_anthro
-
-    logical,        intent(in) :: ca_global_emis, do_sppt_emis
-    real, optional, intent(in) :: emis_multiplier(:)
+    real(kind_phys),intent(in) :: dt, emis_amp_anthro, pert_scale_anthro
+    real(kind_phys), optional, intent(in) :: sppt_wts(:,:)
+    logical,        intent(in) :: do_sppt_emis
 
     integer, parameter :: ids=1,jds=1,jde=1, kds=1
     integer, parameter :: ims=1,jms=1,jme=1, kms=1
@@ -102,10 +100,8 @@ subroutine gsd_chem_anthropogenic_wrapper_run(im, kte, kme, ktau, dt,
    !ppm2ugkg(p_so2 ) = 1.e+03_kind_phys * mw_so2_aer / mwdry
     ppm2ugkg(p_sulf) = 1.e+03_kind_phys * mw_so4_aer / mwdry
 
-    if(do_sppt_emis .or. ca_global_emis) then
-      do i = ims, im
-        random_factor(i,jms) = min(10.0,max(0.0,(emis_multiplier(i)-1.0)*emis_amp_anthro + 1.0))
-      enddo
+    if(do_sppt_emis) then
+      random_factor(:,jms) = pert_scale_anthro*max(min(1+(sppt_wts(:,kme/2)-1)*emis_amp_anthro,2.0),0.0)
     else
       random_factor = 1.0
     endif
diff --git a/gsdchem/gsd_chem_anthropogenic_wrapper.meta b/gsdchem/gsd_chem_anthropogenic_wrapper.meta
@@ -220,32 +220,24 @@
   type = integer
   intent = in
   optional = F
-[emis_amp_anthro]
-  standard_name = anthropogenic_emissions_perturbation_amplitude
-  long_name = multiplier of emissions random perturbation of anthropogenic emissions
+[pert_scale_anthro]
+  standard_name = anthropogenic_scaling_factor
+  long_name = Scaling factor for emissions of anthropogenic emissions
   units = none
   dimensions = ()
   type = real
   kind = kind_phys
   intent = in
   optional = F
-[emis_multiplier]
-  standard_name = gsd_chem_ca_global_emis_multiplier
-  long_name = fraction of emissions to generate based on cellular automata
-  units = frac
-  dimensions = (horizontal_dimension)
+[emis_amp_anthro]
+  standard_name = anthropogenic_emissions_perturbation_amplitude
+  long_name = multiplier of emissions random perturbation of anthropogenic emissions
+  units = none
+  dimensions = ()
   type = real
   kind = kind_phys
   intent = in
   optional = F
-[ca_global_emis]
-  standard_name = flag_for_tracer_emissions_global_cellular_automata
-  long_name = switch for global ca applied to chemical tracer emissions
-  units = flag
-  dimensions = ()
-  type = logical
-  intent = in
-  optional = F
 [do_sppt_emis]
   standard_name = flag_for_stochastic_emissions_perturbations
   long_name = flag for stochastic emissions perturbations
@@ -254,6 +246,15 @@
   type = logical
   intent = in
   optional = F
+[sppt_wts]
+  standard_name = weights_for_stochastic_sppt_perturbation
+  long_name = weights for stochastic sppt perturbation
+  units = none
+  dimensions = (horizontal_dimension,vertical_dimension)
+  type = real
+  kind = kind_phys
+  intent = in
+  optional = F
 [errmsg]
   standard_name = ccpp_error_message
   long_name = error message for error handling in CCPP
diff --git a/gsdchem/gsd_chem_dust_wrapper.F90 b/gsdchem/gsd_chem_dust_wrapper.F90
@@ -50,20 +50,19 @@ subroutine gsd_chem_dust_wrapper_run(im, kte, kme, ktau, dt, garea, land,   &
                    ntdust1,ntdust2,ntdust3,ntdust4,ntdust5,ndust,               &
                    gq0,qgrs,duem,                                               &
                    chem_opt_in,dust_opt_in,dust_calcdrag_in,                    &
-                   dust_alpha_in,dust_gamma_in,                                 &
-                   emis_amp_dust, emis_multiplier, do_sppt_emis, ca_global_emis,&
-                   errmsg,errflg)
+                   dust_alpha_in,dust_gamma_in,pert_scale_dust,                 &
+                   emis_amp_dust, do_sppt_emis, sppt_wts, errmsg,errflg)
 
     implicit none
 
 
     integer,        intent(in) :: im,kte,kme,ktau,nsoil
     integer,        intent(in) :: nseasalt,ntrac
     integer,        intent(in) :: ntdust1,ntdust2,ntdust3,ntdust4,ntdust5,ndust
-    real(kind_phys),intent(in) :: dt, emis_amp_dust
+    real(kind_phys),intent(in) :: dt, emis_amp_dust, pert_scale_dust
 
-    logical,        intent(in) :: ca_global_emis, do_sppt_emis
-    real, optional, intent(in) :: emis_multiplier(:)
+    logical,        intent(in) :: do_sppt_emis
+    real(kind_phys), optional, intent(in) :: sppt_wts(:,:)
 
     integer, parameter :: ids=1,jds=1,jde=1, kds=1
     integer, parameter :: ims=1,jms=1,jme=1, kms=1
@@ -135,10 +134,8 @@ subroutine gsd_chem_dust_wrapper_run(im, kte, kme, ktau, dt, garea, land,   &
     ite=im
     kde=kte
 
-    if(do_sppt_emis .or. ca_global_emis) then
-      do i = ims, im
-        random_factor(i,jms) = min(10.0,max(0.0,(emis_multiplier(i)-1.0)*emis_amp_dust + 1.0))
-      enddo
+    if(do_sppt_emis) then
+      random_factor(:,jms) = pert_scale_dust*max(min(1+(sppt_wts(:,kme/2)-1)*emis_amp_dust,2.0),0.0)
     else
       random_factor = 1.0
     endif
diff --git a/gsdchem/gsd_chem_dust_wrapper.meta b/gsdchem/gsd_chem_dust_wrapper.meta
@@ -439,32 +439,24 @@
   kind = kind_phys
   intent = in
   optional = F
-[emis_amp_dust]
-  standard_name = dust_emissions_perturbation_amplitude
-  long_name = multiplier of emissions random perturbation of dust emissions
+[pert_scale_dust]
+  standard_name = dust_emissions_scaling_factor
+  long_name = Scaling factor for emissions of dust emissions
   units = none
   dimensions = ()
   type = real
   kind = kind_phys
   intent = in
   optional = F
-[emis_multiplier]
-  standard_name = gsd_chem_ca_global_emis_multiplier
-  long_name = fraction of emissions to generate based on cellular automata
-  units = frac
-  dimensions = (horizontal_dimension)
+[emis_amp_dust]
+  standard_name = dust_emissions_perturbation_amplitude
+  long_name = multiplier of emissions random perturbation of dust emissions
+  units = none
+  dimensions = ()
   type = real
   kind = kind_phys
   intent = in
   optional = F
-[ca_global_emis]
-  standard_name = flag_for_tracer_emissions_global_cellular_automata
-  long_name = switch for global ca applied to chemical tracer emissions
-  units = flag
-  dimensions = ()
-  type = logical
-  intent = in
-  optional = F
 [do_sppt_emis]
   standard_name = flag_for_stochastic_emissions_perturbations
   long_name = flag for stochastic emissions perturbations
@@ -473,6 +465,15 @@
   type = logical
   intent = in
   optional = F
+[sppt_wts]
+  standard_name = weights_for_stochastic_sppt_perturbation
+  long_name = weights for stochastic sppt perturbation
+  units = none
+  dimensions = (horizontal_dimension,vertical_dimension)
+  type = real
+  kind = kind_phys
+  intent = in
+  optional = F
 [errmsg]
   standard_name = ccpp_error_message
   long_name = error message for error handling in CCPP
diff --git a/gsdchem/gsd_chem_plume_wrapper.F90 b/gsdchem/gsd_chem_plume_wrapper.F90
@@ -44,24 +44,23 @@ subroutine gsd_chem_plume_wrapper_run(im, kte, kme, ktau, dt,                &
                    pr3d, ph3d,phl3d, prl3d, tk3d, us3d, vs3d, spechum,           &
                    w,vegtype,fire_GBBEPx,fire_MODIS,                             &
                    ntrac,ntso2,ntpp25,ntbc1,ntoc1,ntpp10,                        &
-                   gq0,qgrs,ebu,abem,                                            &
-                   biomass_burn_opt_in,plumerise_flag_in,plumerisefire_frq_in,   &
-                   emis_amp_plume, emis_multiplier, do_sppt_emis,                &
-                   ca_global_emis, errmsg,errflg)
+                   gq0,qgrs,ebu,abem,biomass_burn_opt_in,plumerise_flag_in,      &
+                   plumerisefire_frq_in,pert_scale_plume,                        &
+                   emis_amp_plume, do_sppt_emis, sppt_wts, errmsg,errflg)
 
     implicit none
 
 
     integer,        intent(in) :: im,kte,kme,ktau
     integer,        intent(in) :: ntrac,ntso2,ntpp25,ntbc1,ntoc1,ntpp10
-    real(kind_phys),intent(in) :: dt, emis_amp_plume
+    real(kind_phys),intent(in) :: dt, emis_amp_plume, pert_scale_plume
 
     integer, parameter :: ids=1,jds=1,jde=1, kds=1
     integer, parameter :: ims=1,jms=1,jme=1, kms=1
     integer, parameter :: its=1,jts=1,jte=1, kts=1
 
-    logical,        intent(in) :: ca_global_emis, do_sppt_emis
-    real, optional, intent(in) :: emis_multiplier(:)
+    logical,        intent(in) :: do_sppt_emis
+    real(kind_phys), optional, intent(in) :: sppt_wts(:,:)
     integer, dimension(im), intent(in) :: vegtype    
     real(kind_phys), dimension(im,    5), intent(in) :: fire_GBBEPx
     real(kind_phys), dimension(im,   13), intent(in) :: fire_MODIS
@@ -175,10 +174,8 @@ subroutine gsd_chem_plume_wrapper_run(im, kte, kme, ktau, dt,                &
     if (biomass_burn_opt == BURN_OPT_ENABLE) then
       jp = jte
 
-      if(plumerise_flag == FIRE_OPT_GBBEPx .and. (do_sppt_emis .or. ca_global_emis)) then
-        do i = ims, im
-          random_factor(i) = min(10.0,max(0.0,(emis_multiplier(i)-1.0)*emis_amp_plume + 1.0))
-        enddo
+      if(plumerise_flag == FIRE_OPT_GBBEPx .and. do_sppt_emis) then
+        random_factor(:) = pert_scale_plume*max(min(1+(sppt_wts(:,kme/2)-1)*emis_amp_plume,2.0),0.0)
       endif
 
       factor3 = 0._kind_phys
diff --git a/gsdchem/gsd_chem_plume_wrapper.meta b/gsdchem/gsd_chem_plume_wrapper.meta
@@ -249,6 +249,15 @@
   kind = kind_phys
   intent = inout
   optional = F
+[pert_scale_plume]
+  standard_name = plume_emissions_scaling_factor
+  long_name = Scaling factor for emissions of plume rising
+  units = none
+  dimensions = ()
+  type = real
+  kind = kind_phys
+  intent = in
+  optional = F
 [biomass_burn_opt_in]
   standard_name = gsd_chem_biomass_burn_opt
   long_name = gsd chem biomass burning option
@@ -280,23 +289,6 @@
   dimensions = ()
   type = real
   kind = kind_phys
-[emis_multiplier]
-  standard_name = gsd_chem_ca_global_emis_multiplier
-  long_name = fraction of emissions to generate based on cellular automata
-  units = frac
-  dimensions = (horizontal_dimension)
-  type = real
-  kind = kind_phys
-  intent = in
-  optional = F
-[ca_global_emis]
-  standard_name = flag_for_tracer_emissions_global_cellular_automata
-  long_name = switch for global ca applied to chemical tracer emissions
-  units = flag
-  dimensions = ()
-  type = logical
-  intent = in
-  optional = F
 [do_sppt_emis]
   standard_name = flag_for_stochastic_emissions_perturbations
   long_name = flag for stochastic emissions perturbations
@@ -305,6 +297,15 @@
   type = logical
   intent = in
   optional = F
+[sppt_wts]
+  standard_name = weights_for_stochastic_sppt_perturbation
+  long_name = weights for stochastic sppt perturbation
+  units = none
+  dimensions = (horizontal_dimension,vertical_dimension)
+  type = real
+  kind = kind_phys
+  intent = in
+  optional = F
 [errmsg]
   standard_name = ccpp_error_message
   long_name = error message for error handling in CCPP
diff --git a/gsdchem/gsd_chem_seas_wrapper.F90 b/gsdchem/gsd_chem_seas_wrapper.F90
@@ -19,19 +19,7 @@ module gsd_chem_seas_wrapper
 
 !> \brief Brief description of the subroutine
 !!
-      subroutine gsd_chem_seas_wrapper_init(ca_global_emis,do_sppt_emis,im,emis_multiplier,errmsg,errflg)
-        implicit none
-        logical, intent(in) :: ca_global_emis, do_sppt_emis
-        real, intent(out) :: emis_multiplier(:)
-        character(len=*), intent(out) :: errmsg
-        integer, intent(out) :: errflg, im
-
-        errmsg=''
-        errflg=0
-
-        if(ca_global_emis .or. do_sppt_emis) then
-          emis_multiplier=1.0
-        endif
+      subroutine gsd_chem_seas_wrapper_init
       end subroutine gsd_chem_seas_wrapper_init
 
 !> \brief Brief description of the subroutine
@@ -51,23 +39,21 @@ end subroutine gsd_chem_seas_wrapper_finalize
 !!
 !>\section gsd_chem_seas_wrapper GSD Chemistry Scheme General Algorithm
 !> @{
-    subroutine gsd_chem_seas_wrapper_run(im, kte, kme, ktau, dt, garea,       &
-                   land, u10m, v10m, ustar, rlat, rlon, tskin,                &
-                   pr3d, ph3d,prl3d, tk3d, us3d, vs3d, spechum,               &
-                   nseasalt,ntrac,ntss1,ntss2,ntss3,ntss4,ntss5,              &
-                   gq0,qgrs,ssem,seas_opt_in, emis_amp_seas, emis_multiplier, &
-                   ca1, ca_global_emis, do_sppt_emis, sppt_wts, errmsg, errflg)
+    subroutine gsd_chem_seas_wrapper_run(im, kte, kme, ktau, dt, garea,          &
+                   land, u10m, v10m, ustar, rlat, rlon, tskin,                   &
+                   pr3d, ph3d,prl3d, tk3d, us3d, vs3d, spechum,                  &
+                   nseasalt,ntrac,ntss1,ntss2,ntss3,ntss4,ntss5,                 &
+                   gq0,qgrs,ssem,seas_opt_in, pert_scale_seas, emis_amp_seas,    &
+                   do_sppt_emis, sppt_wts, errmsg, errflg)
 
     implicit none
 
-
     integer,        intent(in) :: im,kte,kme,ktau
     integer,        intent(in) :: nseasalt,ntrac,ntss1,ntss2,ntss3,ntss4,ntss5
     real(kind_phys),intent(in) :: dt
 
-    logical,        intent(in) :: ca_global_emis, do_sppt_emis
-    real, optional, intent(inout) :: emis_multiplier(:)
-    real, intent(in)    :: ca1(im), emis_amp_seas
+    logical,        intent(in) :: do_sppt_emis
+    real(kind=kind_phys), intent(in)    :: emis_amp_seas, pert_scale_seas
     real(kind_phys), optional, intent(in) :: sppt_wts(:,:)
 
 
@@ -100,7 +86,7 @@ subroutine gsd_chem_seas_wrapper_run(im, kte, kme, ktau, dt, garea,       &
 
     integer :: ide, ime, ite, kde
     real(kind_phys), dimension(1:num_chem) :: ppm2ugkg
-    real(kind_phys) :: ca1_scaled, random_factor(ims:im,jms:jme)
+    real(kind_phys) :: random_factor(ims:im,jms:jme)
 
 !>-- local variables
     integer :: i, j, jp, k, kp, n
@@ -122,28 +108,10 @@ subroutine gsd_chem_seas_wrapper_run(im, kte, kme, ktau, dt, garea,       &
    !ppm2ugkg(p_so2 ) = 1.e+03_kind_phys * mw_so2_aer / mwdry
     ppm2ugkg(p_sulf) = 1.e+03_kind_phys * mw_so4_aer / mwdry
 
-    random_factor = 1
-
-    if (do_sppt_emis) then
-      do i = ims, im
-        emis_multiplier(i) = sppt_wts(i,kme/2)
-      enddo
-    elseif (ca_global_emis) then
-      do i = ims, im
-        ! ca1(i) is always precisely 0 or 2
-        if(ca1(i)<1.0) then
-          ca1_scaled=0.9
-        else
-          ca1_scaled=1.0/0.9
-        endif
-        emis_multiplier(i) = max(0.8,min(1.2,emis_multiplier(i)*0.95 + ca1_scaled*0.05))
-      enddo
-    endif
-
-    if(do_sppt_emis .or. ca_global_emis) then
-      do i=ims,im
-        random_factor(i,jms) = min(10.0,max(0.0,(emis_multiplier(i)-1.0)*emis_amp_seas + 1.0))
-      enddo
+    if(do_sppt_emis) then
+      random_factor(:,jms) = pert_scale_seas*max(min(1+(sppt_wts(:,kme/2)-1)*emis_amp_seas,2.0),0.0)
+    else
+      random_factor = 1.0
     endif
 
 !>- get ready for chemistry run
diff --git a/gsdchem/gsd_chem_seas_wrapper.meta b/gsdchem/gsd_chem_seas_wrapper.meta
