diff --git a/FV3GFS/FV3GFS_io.F90 b/FV3GFS/FV3GFS_io.F90 index b2d558ec..f01b4603 100644 --- a/FV3GFS/FV3GFS_io.F90 +++ b/FV3GFS/FV3GFS_io.F90 @@ -4791,6 +4791,30 @@ subroutine gfdl_diag_register(Time, Sfcprop, Gfs_diag, Model, Cldprop, Atm_block enddo !Sofar added: end + idx = idx + 1 + Diag(idx)%axes = 2 + Diag(idx)%name = 'Rb' + Diag(idx)%desc = 'Surface layer bulk Richardson number' + Diag(idx)%unit = '1' + Diag(idx)%mod_name = 'gfs_phys' + Diag(idx)%intpl_method = 'bilinear' + allocate (Diag(idx)%data(nblks)) + do nb = 1,nblks + Diag(idx)%data(nb)%var2 => Gfs_diag(nb)%rb(:) + enddo + + idx = idx + 1 + Diag(idx)%axes = 2 + Diag(idx)%name = 'zol' + Diag(idx)%desc = 'Dimensionless surface layer stability' + Diag(idx)%unit = '1' + Diag(idx)%mod_name = 'gfs_phys' + Diag(idx)%intpl_method = 'bilinear' + allocate (Diag(idx)%data(nblks)) + do nb = 1,nblks + Diag(idx)%data(nb)%var2 => Gfs_diag(nb)%zol(:) + enddo + idx = idx + 1 Diag(idx)%axes = 2 Diag(idx)%name = 'dpt2m' diff --git a/GFS_layer/GFS_physics_driver.F90 b/GFS_layer/GFS_physics_driver.F90 index 90ead827..8f7c7604 100644 --- a/GFS_layer/GFS_physics_driver.F90 +++ b/GFS_layer/GFS_physics_driver.F90 @@ -1202,7 +1202,7 @@ subroutine GFS_physics_driver & Sfcprop%ffmm, Sfcprop%ffhh, Sfcprop%uustar, & wind, Tbd%phy_f2d(1,Model%num_p2d), fm10, fh2, & sigmaf, vegtype, Sfcprop%shdmax, Model%ivegsrc, & - tsurf, flag_iter) !, Model%redrag, Model%z0s_max, & + tsurf, flag_iter, Diag%zol) !, Model%redrag, Model%z0s_max, & !Model%do_z0_moon, Model%do_z0_hwrf15, & !Model%do_z0_hwrf17, Model%do_z0_hwrf17_hwonly, & !Model%wind_th_hwrf) @@ -1225,7 +1225,7 @@ subroutine GFS_physics_driver & tsurf, flag_iter, Model%redrag, Model%z0s_max, & Model%do_z0_moon, Model%do_z0_hwrf15, & Model%do_z0_hwrf17, Model%do_z0_hwrf17_hwonly, & - Model%wind_th_hwrf) + Model%wind_th_hwrf, Diag%zol) else ! GFS original sfc_diff modified by kgao call sfc_diff (im,Statein%pgr, Statein%ugrs, Statein%vgrs,& @@ -1239,8 +1239,9 @@ subroutine GFS_physics_driver & Model%z0s_max, & Model%do_z0_moon, Model%do_z0_hwrf15, & Model%do_z0_hwrf17, Model%do_z0_hwrf17_hwonly, & - Model%wind_th_hwrf) + Model%wind_th_hwrf, Diag%zol) endif + Diag%rb = rb !endif !endif diff --git a/GFS_layer/GFS_typedefs.F90 b/GFS_layer/GFS_typedefs.F90 index c75ed8b5..eeaedfb7 100644 --- a/GFS_layer/GFS_typedefs.F90 +++ b/GFS_layer/GFS_typedefs.F90 @@ -1330,6 +1330,8 @@ module GFS_typedefs real (kind=kind_phys), pointer :: evap (:) => null() !< sfc moisture flux real (kind=kind_phys), pointer :: u10n (:) => null() !< 10 meter u/v neutral wind speed !Sofar added 10/19/23 real (kind=kind_phys), pointer :: v10n (:) => null() !< 10 meter u/v neutral wind speed !Sofar added 10/19/23 + real (kind=kind_phys), pointer :: zol (:) => null() !< Dimensionless surface layer stability + real (kind=kind_phys), pointer :: rb (:) => null() !< Surface layer bulk richardson number real (kind=kind_phys), pointer :: dpt2m (:) => null() !< 2 meter dew point temperature real (kind=kind_phys), pointer :: zlvl (:) => null() !< layer 1 height (m) real (kind=kind_phys), pointer :: psurf (:) => null() !< surface pressure (Pa) @@ -4078,6 +4080,8 @@ subroutine diag_create (Diag, IM, Model) allocate (Diag%evap (IM)) allocate (Diag%u10n (IM)) ! Sofar added: 10/19/13 allocate (Diag%v10n (IM)) ! Sofar added: 10/19/13 + allocate (Diag%zol (IM)) + allocate (Diag%rb (IM)) allocate (Diag%dpt2m (IM)) allocate (Diag%zlvl (IM)) allocate (Diag%psurf (IM)) @@ -4376,6 +4380,8 @@ subroutine diag_phys_zero (Diag, Model, linit, iauwindow_center) Diag%evap = zero Diag%u10n = zero ! Sofar added: 10/19/23 Diag%v10n = zero ! Sofar added: 10/19/23 + Diag%zol = zero + Diag%rb = zero Diag%dpt2m = zero Diag%zlvl = zero Diag%psurf = zero diff --git a/gsmphys/sfc_diff_coupled.f b/gsmphys/sfc_diff_coupled.f index dfe08bcb..1c1b7f18 100644 --- a/gsmphys/sfc_diff_coupled.f +++ b/gsmphys/sfc_diff_coupled.f @@ -4,7 +4,7 @@ subroutine sfc_diff_coupled(im,ps,u1,v1,t1,q1,z1, & stress,fm,fh, & ustar,wind,ddvel,fm10,fh2, & sigmaf,vegtype,shdmax,ivegsrc, - & tsurf,flag_iter) !,redrag, + & tsurf,flag_iter, zol) !,redrag, ! & z0s_max) ! & do_z0_moon, do_z0_hwrf15, do_z0_hwrf17, ! & do_z0_hwrf17_hwonly, wind_th_hwrf) @@ -26,7 +26,7 @@ subroutine sfc_diff_coupled(im,ps,u1,v1,t1,q1,z1, &, prsl1, prslki, stress &, fm, fh, ustar, wind, ddvel &, fm10, fh2, sigmaf, shdmax - &, tsurf, snwdph + &, tsurf, snwdph, zol integer, dimension(im) ::vegtype, islimsk logical flag_iter(im) @@ -145,7 +145,7 @@ subroutine sfc_diff_coupled(im,ps,u1,v1,t1,q1,z1, call monin_obukhov_similarity & (z1(i), snwdph(i), thv1, wind(i), z0max, ztmax, tvs, & rb(i), fm(i), fh(i), fm10(i), fh2(i), - & cm(i), ch(i), stress(i), ustar(i)) + & cm(i), ch(i), stress(i), ustar(i), zol(i)) elseif (islimsk(i) == 0) then @@ -167,7 +167,7 @@ subroutine sfc_diff_coupled(im,ps,u1,v1,t1,q1,z1, call monin_obukhov_similarity & (z1(i), snwdph(i), thv1, wind(i), z0max, ztmax, tvs, & rb(i), fm(i), fh(i), fm10(i), fh2(i), - & cm(i), ch(i), tem1, tem2) !stress(i), ustar(i)) + & cm(i), ch(i), tem1, tem2, zol(i)) !stress(i), ustar(i)) ! kgao: use ustar from coupler to get stress stress(i) = ustar(i) * ustar(i) diff --git a/gsmphys/sfc_diff_gfdl.f b/gsmphys/sfc_diff_gfdl.f index 87aeaf7f..bb2bc69b 100644 --- a/gsmphys/sfc_diff_gfdl.f +++ b/gsmphys/sfc_diff_gfdl.f @@ -22,7 +22,7 @@ subroutine sfc_diff_gfdl(im,ps,u1,v1,t1,q1,z1, & tsurf,flag_iter,redrag, & z0s_max, & do_z0_moon, do_z0_hwrf15, do_z0_hwrf17, - & do_z0_hwrf17_hwonly, wind_th_hwrf) + & do_z0_hwrf17_hwonly, wind_th_hwrf, zol) ! oct 2019 - a clean and updated version by Kun Gao at GFDL (Kun.Gao@noaa.gov) @@ -49,8 +49,9 @@ subroutine sfc_diff_gfdl(im,ps,u1,v1,t1,q1,z1, &, ddvel &, fm10, fh2, sigmaf, shdmax &, tsurf, snwdph - &, fm_neutral, fm10_neutral ! Sofar added Spring 2023 - real(kind=kind_phys) :: ws1, ws10n ! Sofar added 10/19/23 + &, fm_neutral, fm10_neutral, + & zol + real(kind=kind_phys) :: ws1, ws10n integer, dimension(im) ::vegtype, islimsk logical flag_iter(im) @@ -174,7 +175,7 @@ subroutine sfc_diff_gfdl(im,ps,u1,v1,t1,q1,z1, & (z1(i), snwdph(i), thv1, wind(i), z0max, ztmax, tvs, & rb(i), fm(i), fh(i), fm10(i), fh2(i), & fm_neutral(i), fm10_neutral(i), !(ADDED by Sofar) - & cm(i), ch(i), stress(i), ustar(i)) + & cm(i), ch(i), stress(i), ustar(i), zol(i)) elseif (islimsk(i) == 0) then ! over water @@ -205,7 +206,7 @@ subroutine sfc_diff_gfdl(im,ps,u1,v1,t1,q1,z1, & (z1(i), snwdph(i), thv1, wind(i), z0max, ztmax, tvs, & rb(i), fm(i), fh(i), fm10(i), fh2(i), & fm_neutral(i), fm10_neutral(i), !(ADDED by Sofar) - & cm(i), ch(i), stress(i), ustar(i)) + & cm(i), ch(i), stress(i), ustar(i), zol(i)) ! === iteration 2 @@ -305,7 +306,7 @@ subroutine sfc_diff_gfdl(im,ps,u1,v1,t1,q1,z1, & (z1(i), snwdph(i), thv1, wind(i), z0max, ztmax, tvs, & rb(i), fm(i), fh(i), fm10(i), fh2(i), & fm_neutral(i), fm10_neutral(i), !(ADDED by Sofar) - & cm(i), ch(i), stress(i), ustar(i)) + & cm(i), ch(i), stress(i), ustar(i), zol(i)) z0rl(i) = 100.0 * z0max ztrl(i) = 100.0 * ztmax @@ -511,7 +512,7 @@ subroutine monin_obukhov_similarity & ( z1, snwdph, thv1, wind, z0max, ztmax, tvs, & rb, fm, fh, fm10, fh2, & fm_neutral, fm10_neutral, !(ADDED by Sofar) - & cm, ch, stress, ustar) + & cm, ch, stress, ustar, zol) ! --- input ! z1 - lowest model level height @@ -530,6 +531,7 @@ subroutine monin_obukhov_similarity ! cm, ch - surface exchange coefficients for momentum and heat ! stress - surface wind stress ! ustar - surface frictional velocity +! zol - dimensionless stability use machine , only : kind_phys use physcons, grav => con_g @@ -585,6 +587,7 @@ subroutine monin_obukhov_similarity fh2 = log((ztmax+2.) * tem2) hlinf = rb * fm * fm / fh hlinf = min(max(hlinf,ztmin1),ztmax1) + zol = hlinf fm_neutral = fm !(ADDED by Sofar) fm10_neutral = fm10 !(ADDED by Sofar)