HAFS related modifications for GFS sa-TKE EDMF PBL scheme and for coupling without a sea ice component

physics/GFS_surface_composites.F90

@@ -28,7 +28,7 @@ end subroutine GFS_surface_composites_pre_finalize
 !! \htmlinclude GFS_surface_composites_pre_run.html
 !!
    subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, lsm, lsm_noahmp, lsm_ruc, frac_grid,      &
-                                 flag_cice, cplflx, cplwav2atm, landfrac, lakefrac, lakedepth, oceanfrac, frland,         &
+                                 flag_cice, cplflx, cplice, cplwav2atm, landfrac, lakefrac, lakedepth, oceanfrac, frland, &
                                  dry, icy, lake, use_flake, ocean, wet, hice, cice, zorlo, zorll, zorli,                  &
                                  snowd,            snowd_lnd, snowd_ice, tprcp, tprcp_wat,                                &
                                  tprcp_lnd, tprcp_ice, uustar, uustar_wat, uustar_lnd, uustar_ice,                        &
@@ -43,7 +43,7 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, lsm
       ! Interface variables
       integer,                             intent(in   ) :: im, lkm, kdt
       integer,                             intent(in   ) :: lsm, lsm_noahmp, lsm_ruc
-      logical,                             intent(in   ) :: flag_init, flag_restart, frac_grid, cplflx, cplwav2atm
+      logical,                             intent(in   ) :: flag_init, flag_restart, frac_grid, cplflx, cplice, cplwav2atm
       logical, dimension(:),              intent(inout)  :: flag_cice
       logical,              dimension(:), intent(inout)  :: dry, icy, lake, use_flake, ocean, wet
       real(kind=kind_phys), dimension(:), intent(in   )  :: landfrac, lakefrac, lakedepth, oceanfrac
@@ -155,7 +155,14 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, lsm
               if (cice(i) >= min_seaice) then
                 icy(i)   = .true.
                 tisfc(i) = max(timin, min(tisfc(i), tgice))
-                if (cplflx)  then
+                ! This cplice namelist option was added to deal with the
+                ! situation of the FV3ATM-HYCOM coupling without an active sea
+                ! ice (e.g., CICE6) component. By default, the cplice is true
+                ! when cplflx is .true. (e.g., for the S2S application).
+                ! Whereas, for the HAFS FV3ATM-HYCOM coupling, cplice is set as
+                ! .false.. In the future HAFS FV3ATM-MOM6 coupling, the cplflx
+                ! could be .true., while cplice being .false.. 
+                if (cplice .and. cplflx)  then
                   islmsk_cice(i) = 4
                   flag_cice(i)   = .true.
                 else
@@ -173,7 +180,11 @@ subroutine GFS_surface_composites_pre_run (im, flag_init, flag_restart, lkm, lsm
               endif
               if (cice(i) < one) then
                 wet(i) = .true. ! some open ocean
-                if (.not. cplflx .and. icy(i)) tsfco(i) = max(tisfc(i), tgice)
+                if (cplice) then
+                  if (.not. cplflx .and. icy(i)) tsfco(i) = max(tisfc(i), tgice)
+                else
+                  if (icy(i)) tsfco(i) = max(tisfc(i), tgice)
+                endif
               endif
             else
               if (cice(i) >= min_lakeice) then

physics/GFS_surface_composites.meta

@@ -87,6 +87,14 @@
   type = logical
   intent = in
   optional = F
+[cplice]
+  standard_name = flag_for_sea_ice_coupling
+  long_name = flag controlling cplice collection (default on)
+  units = flag
+  dimensions = ()
+  type = logical
+  intent = in
+  optional = F
 [cplwav2atm]
   standard_name = flag_for_one_way_ocean_wave_coupling_to_atmosphere
   long_name = flag controlling ocean wave coupling to the atmosphere (default off)

physics/satmedmfvdifq.F

@@ -75,6 +75,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntiw,ntke,          &
      &     prsi,del,prsl,prslk,phii,phil,delt,                          &
      &     dspheat,dusfc,dvsfc,dtsfc,dqsfc,hpbl,dkt,dku,                &
      &     kinver,xkzm_m,xkzm_h,xkzm_s,dspfac,bl_upfr,bl_dnfr,          &
+     &     rlmx,elmx,sfc_rlm,                                           &
      &     ntoz,ntqv,dtend,dtidx,index_of_temperature,index_of_x_wind,  &
      &     index_of_y_wind,index_of_process_pbl,gen_tend,ldiag3d,       &
      &     errmsg,errflg)
@@ -86,6 +87,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntiw,ntke,          &
 !
 !----------------------------------------------------------------------
       integer, intent(in)  :: im, km, ntrac, ntcw, ntiw, ntke, ntoz,ntqv
+      integer, intent(in)  :: sfc_rlm
       integer, intent(in)  :: kinver(:)
       integer, intent(out) :: kpbl(:)
       logical, intent(in)  :: gen_tend,ldiag3d
@@ -94,6 +96,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntiw,ntke,          &
      &                                    eps,epsm1
       real(kind=kind_phys), intent(in) :: delt, xkzm_m, xkzm_h, xkzm_s
       real(kind=kind_phys), intent(in) :: dspfac, bl_upfr, bl_dnfr
+      real(kind=kind_phys), intent(in) :: rlmx, elmx
       real(kind=kind_phys), intent(inout) :: dv(:,:),     du(:,:),      &
      &                     tdt(:,:),    rtg(:,:,:)
       real(kind=kind_phys), intent(in) ::                               &
@@ -212,7 +215,6 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntiw,ntke,          &
      &                     rdt,     rdz,    qmin,   qlmin,
      &                     rimin,   rbcr,   rbint,  tdzmin,
      &                     rlmn,    rlmn0,  rlmn1,  rlmn2,
-     &                     rlmx,    elmx,
      &                     ttend,   utend,  vtend,  qtend,
      &                     zfac,    zfmin,  vk,     spdk2,
      &                     tkmin,   tkbmx,  xkgdx,
@@ -232,7 +234,6 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntiw,ntke,          &
       parameter(vk=0.4,rimin=-100.)
       parameter(rbcr=0.25,zolcru=-0.02,tdzmin=1.e-3)
       parameter(rlmn=30.,rlmn0=5.,rlmn1=5.,rlmn2=10.)
-      parameter(rlmx=300.,elmx=300.)
       parameter(prmin=0.25,prmax=4.0)
       parameter(pr0=1.0,prtke=1.0,prscu=0.67)
       parameter(f0=1.e-4,crbmin=0.15,crbmax=0.35)
@@ -1029,6 +1030,12 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntiw,ntke,          &
             zk = tem / ptem
           endif 
           elm(i,k) = zk*rlam(i,k)/(rlam(i,k)+zk)
+
+!>  - If sfc_rlm=1, use zk for elm within surface layer
+          if ( sfc_rlm == 1 ) then
+            if ( sfcflg(i) .and.
+     &           zl(i,k) < min(100.0,hpbl(i)*0.05) ) elm(i,k)=zk
+          endif
 !
           dz = zi(i,k+1) - zi(i,k)
           tem = max(gdx(i),dz)

physics/satmedmfvdifq.meta

@@ -611,6 +611,32 @@
   kind = kind_phys
   intent = in
   optional = F
+[rlmx]
+  standard_name = maximum_allowed_mixing_length_in_boundary_layer_mass_flux_scheme
+  long_name = maximum allowed mixing length in boundary layer mass flux scheme
+  units = m
+  dimensions = ()
+  type = real
+  kind = kind_phys
+  intent = in
+  optional = F
+[elmx]
+  standard_name = maximum_allowed_dissipation_mixing_length_in_boundary_layer_mass_flux_scheme
+  long_name = maximum allowed dissipation mixing length in boundary layer mass flux scheme
+  units = m
+  dimensions = ()
+  type = real
+  kind = kind_phys
+  intent = in
+  optional = F
+[sfc_rlm]
+  standard_name = choice_of_near_surface_mixing_length_in_boundary_layer_mass_flux_scheme
+  long_name = choice of near surface mixing length in boundary layer mass flux scheme
+  units = none
+  dimensions = ()
+  type = integer
+  intent = in
+  optional = F
 [ntoz]
   standard_name = index_of_ozone_mixing_ratio_in_tracer_concentration_array
   long_name = tracer index for ozone mixing ratio

physics/sfc_diff.f

@@ -426,7 +426,8 @@ subroutine sfc_diff_run (im,rvrdm1,eps,epsm1,grav,                &  !intent(in)
                  z0rl_wat(i) = 1.0e-4_kp
               endif
 
-            elseif (z0rl_wav(i) <= 1.0e-7_kp) then
+            elseif (z0rl_wav(i) <= 1.0e-7_kp .or.                       &
+     &              z0rl_wav(i) > 1.0_kp) then
 !             z0 = (charnock / grav) * ustar_wat(i) * ustar_wat(i)
               tem1 = 0.11 * vis / ustar_wat(i)
               z0 = tem1 + (charnock/grav)*ustar_wat(i)*ustar_wat(i)