Bring AQM model into UFS

.gitmodules

@@ -8,8 +8,10 @@
   branch = main
 [submodule "ccpp/physics"]
   path = ccpp/physics
-  url = https://github.com/NCAR/ccpp-physics
-  branch = main
+  #url = https://github.com/NCAR/ccpp-physics
+  #branch = main
+  url = https://github.com/BrianCurtis-NOAA/ccpp-physics
+  branch = feature/aqm
 [submodule "upp"]
   path = upp
   url = https://github.com/NOAA-EMC/UPP

ccpp/data/GFS_typedefs.F90

@@ -295,6 +295,7 @@ module GFS_typedefs
     real (kind=kind_phys), pointer :: ffmm   (:)   => null()  !< fm parameter from PBL scheme
     real (kind=kind_phys), pointer :: ffhh   (:)   => null()  !< fh parameter from PBL scheme
     real (kind=kind_phys), pointer :: f10m   (:)   => null()  !< fm at 10m - Ratio of sigma level 1 wind and 10m wind
+    real (kind=kind_phys), pointer :: rca     (:)  => null()  !< canopy resistance
     real (kind=kind_phys), pointer :: tprcp  (:)   => null()  !< sfc_fld%tprcp - total precipitation
     real (kind=kind_phys), pointer :: srflag (:)   => null()  !< sfc_fld%srflag - snow/rain flag for precipitation
     real (kind=kind_phys), pointer :: slc    (:,:) => null()  !< liquid soil moisture
@@ -647,6 +648,7 @@ module GFS_typedefs
     logical              :: cplocn2atm      !< default yes ocn->atm coupling
     logical              :: cplwav          !< default no cplwav collection
     logical              :: cplwav2atm      !< default no wav->atm coupling
+    logical              :: cplaqm          !< default no cplaqm collection
     logical              :: cplchm          !< default no cplchm collection
     logical              :: use_cice_alb    !< default .false. - i.e. don't use albedo imported from the ice model
     logical              :: cpl_imp_mrg     !< default no merge import with internal forcings
@@ -1840,6 +1842,9 @@ module GFS_typedefs
     ! Extended output diagnostics for Thompson MP
     real (kind=kind_phys), pointer :: thompson_ext_diag3d (:,:,:) => null() ! extended diagnostic 3d output arrays from Thompson MP
 
+    ! Diagnostics for coupled air quality model
+    real (kind=kind_phys), pointer :: aod   (:)   => null()    !< instantaneous aerosol optical depth ( n/a )
+
     ! Auxiliary output arrays for debugging
     real (kind=kind_phys), pointer :: aux2d(:,:)  => null()    !< auxiliary 2d arrays in output (for debugging)
     real (kind=kind_phys), pointer :: aux3d(:,:,:)=> null()    !< auxiliary 2d arrays in output (for debugging)
@@ -2583,6 +2588,12 @@ subroutine sfcprop_create (Sfcprop, IM, Model)
       Sfcprop%dt_cool = zero
       Sfcprop%qrain   = zero
     endif
+    if (Model%lsm == Model%lsm_noah) then
+      allocate (Sfcprop%xlaixy   (IM))
+      allocate (Sfcprop%rca      (IM))
+      Sfcprop%xlaixy     = clear_val
+      Sfcprop%rca        = clear_val
+    end if
     if (Model%lsm == Model%lsm_ruc .or. Model%lsm == Model%lsm_noahmp) then
       allocate(Sfcprop%raincprv  (IM))
       allocate(Sfcprop%rainncprv (IM))
@@ -3002,6 +3013,23 @@ subroutine coupling_create (Coupling, IM, Model)
       Coupling%pfl_lsan  = clear_val
     endif
 
+    ! -- additional coupling options for air quality
+    if (Model%cplaqm .and. .not.Model%cplflx) then
+      !--- outgoing instantaneous quantities
+      allocate (Coupling%dtsfci_cpl  (IM))
+      allocate (Coupling%dqsfci_cpl  (IM))
+      allocate (Coupling%nswsfci_cpl (IM))
+      allocate (Coupling%t2mi_cpl    (IM))
+      allocate (Coupling%q2mi_cpl    (IM))
+      allocate (Coupling%psurfi_cpl  (IM))
+      Coupling%dtsfci_cpl  = clear_val
+      Coupling%dqsfci_cpl  = clear_val
+      Coupling%nswsfci_cpl = clear_val
+      Coupling%t2mi_cpl    = clear_val
+      Coupling%q2mi_cpl    = clear_val
+      Coupling%psurfi_cpl  = clear_val
+    endif
+
     !--- stochastic physics option
     if (Model%do_sppt .or. Model%ca_global) then
       allocate (Coupling%sppt_wts  (IM,Model%levs))
@@ -3145,6 +3173,7 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
     logical              :: cplocn2atm     = .true.          !< default yes cplocn2atm coupling (turn on the feedback from ocn to atm)
     logical              :: cplwav         = .false.         !< default no cplwav collection
     logical              :: cplwav2atm     = .false.         !< default no cplwav2atm coupling
+    logical              :: cplaqm         = .false.         !< default no cplaqm collection
     logical              :: cplchm         = .false.         !< default no cplchm collection
     logical              :: use_cice_alb   = .false.         !< default no cice albedo
     logical              :: cpl_imp_mrg    = .false.         !< default no merge import with internal forcings
@@ -3621,7 +3650,7 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
     logical :: do_spp       = .false.
 
 !--- aerosol scavenging factors
-    integer, parameter :: max_scav_factors = 25
+    integer, parameter :: max_scav_factors = 183
     character(len=40)  :: fscav_aero(max_scav_factors)
 
     real(kind=kind_phys) :: radar_tten_limits(2) = (/ limit_unspecified, limit_unspecified /)
@@ -3635,8 +3664,8 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
                                naux3d, aux2d_time_avg, aux3d_time_avg, fhcyc,               &
                                thermodyn_id, sfcpress_id,                                   &
                           !--- coupling parameters
-                               cplflx, cplice, cplocn2atm, cplwav, cplwav2atm, cplchm,      &
-                               cpl_imp_mrg, cpl_imp_dbg,                                    &
+                               cplflx, cplice, cplocn2atm, cplwav, cplwav2atm, cplaqm,      &
+                               cplchm, cpl_imp_mrg, cpl_imp_dbg,                            &
                                use_cice_alb,                                                &
 #ifdef IDEA_PHYS
                                lsidea, weimer_model, f107_kp_size, f107_kp_interval,        &
@@ -3948,7 +3977,8 @@ subroutine control_initialize (Model, nlunit, fn_nml, me, master, &
     Model%cplocn2atm       = cplocn2atm
     Model%cplwav           = cplwav
     Model%cplwav2atm       = cplwav2atm
-    Model%cplchm           = cplchm
+    Model%cplaqm           = cplaqm
+    Model%cplchm           = cplchm .or. cplaqm
     Model%use_cice_alb     = use_cice_alb
     Model%cpl_imp_mrg      = cpl_imp_mrg
     Model%cpl_imp_dbg      = cpl_imp_dbg
@@ -5777,6 +5807,7 @@ subroutine control_print(Model)
       print *, ' cplocn2atm        : ', Model%cplocn2atm
       print *, ' cplwav            : ', Model%cplwav
       print *, ' cplwav2atm        : ', Model%cplwav2atm
+      print *, ' cplaqm            : ', Model%cplaqm
       print *, ' cplchm            : ', Model%cplchm
       print *, ' use_cice_alb      : ', Model%use_cice_alb
       print *, ' cpl_imp_mrg       : ', Model%cpl_imp_mrg
@@ -7028,6 +7059,13 @@ subroutine diag_create (Diag, IM, Model)
       Diag%thompson_ext_diag3d = clear_val
     endif
 
+    ! Air quality diagnostics
+    ! -- initialize diagnostic variables
+    if (Model%cplaqm) then
+      allocate (Diag%aod(IM))
+      Diag%aod = zero
+    end if
+
     ! Auxiliary arrays in output for debugging
     if (Model%naux2d>0) then
       allocate (Diag%aux2d(IM,Model%naux2d))

ccpp/data/GFS_typedefs.meta

@@ -955,6 +955,13 @@
   dimensions = (horizontal_loop_extent)
   type = real
   kind = kind_phys
+[rca]
+  standard_name = aerodynamic_resistance_in_canopy
+  long_name = canopy resistance
+  units = s m-1
+  dimensions = (horizontal_loop_extent)
+  type = real
+  kind = kind_phys
 [tprcp]
   standard_name = nonnegative_lwe_thickness_of_precipitation_amount_on_dynamics_timestep
   long_name = total precipitation amount in each time step
@@ -1378,7 +1385,7 @@
   dimensions = (horizontal_loop_extent)
   type = real
   kind = kind_phys
-  active = (control_for_land_surface_scheme == identifier_for_noahmp_land_surface_scheme .or. (control_for_land_surface_scheme == identifier_for_ruc_land_surface_scheme .and. flag_for_reading_leaf_area_index_from_input))
+  active = (control_for_land_surface_scheme == identifier_for_noah_land_surface_scheme .or. control_for_land_surface_scheme == identifier_for_noahmp_land_surface_scheme .or. (control_for_land_surface_scheme == identifier_for_ruc_land_surface_scheme .and. flag_for_reading_leaf_area_index_from_input))
 [xsaixy]
   standard_name = stem_area_index
   long_name = stem area index
@@ -2077,15 +2084,15 @@
   dimensions = (horizontal_loop_extent)
   type = real
   kind = kind_phys
-  active = (flag_for_surface_flux_coupling)
+  active = (flag_for_surface_flux_coupling .or. flag_for_air_quality_coupling)
 [dqsfci_cpl]
   standard_name = surface_upward_latent_heat_flux_for_coupling
   long_name = instantaneous sfc latent heat flux
   units = W m-2
   dimensions = (horizontal_loop_extent)
   type = real
   kind = kind_phys
-  active = (flag_for_surface_flux_coupling)
+  active = (flag_for_surface_flux_coupling .or. flag_for_air_quality_coupling)
 [dlwsfci_cpl]
   standard_name = surface_downwelling_longwave_flux_for_coupling
   long_name = instantaneous sfc downward lw flux
@@ -2149,7 +2156,7 @@
   dimensions = (horizontal_loop_extent)
   type = real
   kind = kind_phys
-  active = (flag_for_surface_flux_coupling)
+  active = (flag_for_surface_flux_coupling .or. flag_for_air_quality_coupling)
 [nnirbmi_cpl]
   standard_name = surface_net_downwelling_direct_nir_shortwave_flux_for_coupling
   long_name = instantaneous net nir beam sfc downward sw flux
@@ -2189,15 +2196,15 @@
   dimensions = (horizontal_loop_extent)
   type = real
   kind = kind_phys
-  active = (flag_for_surface_flux_coupling)
+  active = (flag_for_surface_flux_coupling .or. flag_for_air_quality_coupling)
 [q2mi_cpl]
   standard_name = specific_humidity_at_2m_for_coupling
   long_name = instantaneous Q2m
   units = kg kg-1
   dimensions = (horizontal_loop_extent)
   type = real
   kind = kind_phys
-  active = (flag_for_surface_flux_coupling)
+  active = (flag_for_surface_flux_coupling .or. flag_for_air_quality_coupling)
 [u10mi_cpl]
   standard_name = x_wind_at_10m_for_coupling
   long_name = instantaneous U10m
@@ -2229,7 +2236,7 @@
   dimensions = (horizontal_loop_extent)
   type = real
   kind = kind_phys
-  active = (flag_for_surface_flux_coupling)
+  active = (flag_for_surface_flux_coupling .or. flag_for_air_quality_coupling)
 [ulwsfcin_cpl]
   standard_name = surface_upwelling_longwave_flux_from_coupled_process
   long_name = surface upwelling LW flux for coupling
@@ -2747,6 +2754,12 @@
   units = flag
   dimensions = ()
   type = logical
+[cplaqm]
+  standard_name = flag_for_air_quality_coupling
+  long_name = flag controlling cplaqm collection (default off)
+  units = flag
+  dimensions = ()
+  type = logical
 [cplchm]
   standard_name = flag_for_chemistry_coupling
   long_name = flag controlling cplchm collection (default off)

cpl/module_cplfields.F90

@@ -26,7 +26,7 @@ module module_cplfields
   !  l : model levels (3D)
   !  s : surface (2D)
   !  t : tracers (4D)
-  integer,          public, parameter :: NexportFields = 105
+  integer,          public, parameter :: NexportFields = 111
   type(ESMF_Field), target, public    :: exportFields(NexportFields)
 
   type(FieldInfo), dimension(NexportFields), public, parameter :: exportFieldsInfo = [ &
@@ -110,6 +110,12 @@ module module_cplfields
     FieldInfo("lake_fraction                            ", "s"), &
     FieldInfo("ocean_fraction                           ", "s"), &
     FieldInfo("surface_snow_area_fraction               ", "s"), &
+    FieldInfo("canopy_moisture_storage                  ", "s"), &
+    FieldInfo("inst_aerodynamic_conductance             ", "s"), &
+    FieldInfo("inst_canopy_resistance                   ", "s"), &
+    FieldInfo("leaf_area_index                          ", "s"), &
+    FieldInfo("temperature_of_soil_layer                ", "g"), &
+    FieldInfo("height                                   ", "s"), &
 
 
     !  For JEDI
@@ -142,7 +148,7 @@ module module_cplfields
     FieldInfo("t2m                                      ", "s") ]
 
 ! Import Fields ----------------------------------------
-  integer,          public, parameter :: NimportFields = 42
+  integer,          public, parameter :: NimportFields = 43
   logical,          public            :: importFieldsValid(NimportFields)
   type(ESMF_Field), target, public    :: importFields(NimportFields)
 
@@ -164,6 +170,7 @@ module module_cplfields
     FieldInfo("inst_ice_vis_dif_albedo                  ", "s"), &
     FieldInfo("inst_ice_vis_dir_albedo                  ", "s"), &
     FieldInfo("wave_z0_roughness_length                 ", "s"), &
+    FieldInfo("inst_tracer_diag_aod                     ", "s"), &
 
     !  For JEDI
     ! dynamics
@@ -219,7 +226,15 @@ module module_cplfields
     "ice_fraction_in_atm             ", &
     "lake_fraction                   ", &
     "ocean_fraction                  ", &
-    "surface_snow_area_fraction      "  &
+    "surface_snow_area_fraction      ", &
+    "inst_vegetation_area_frac       ", &
+    "canopy_moisture_storage         ", &
+    "inst_aerodynamic_conductance    ", &
+    "inst_canopy_resistance          ", &
+    "leaf_area_index                 ", &
+    "soil_type                       ", &
+    "temperature_of_soil_layer       ", &
+    "height                          "  &
     ]
 
   ! Methods