diff --git a/physics/dcyc2.f b/physics/dcyc2.f index ad9365851..dfa9f02ed 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -178,7 +178,7 @@ subroutine dcyc2t3_run & & sfcnirbmu,sfcnirdfu,sfcvisbmu,sfcvisdfu, & & sfcnirbmd,sfcnirdfd,sfcvisbmd,sfcvisdfd, & & im, levs, deltim, fhswr, & - & dry, icy, wet, & + & dry, icy, wet, damp_LW_fluxadj, lfnc_k, lfnc_p0, & & minGPpres, use_LW_jacobian, sfculw, fluxlwUP_jac, & & t_lay, t_lev, p_lay, p_lev, flux2D_lwUP, flux2D_lwDOWN, & & pert_radtend, do_sppt,ca_global, & @@ -213,10 +213,11 @@ subroutine dcyc2t3_run & ! integer, intent(in) :: ipr ! logical lprnt logical, dimension(:), intent(in) :: dry, icy, wet - logical, intent(in) :: use_LW_jacobian, pert_radtend + logical, intent(in) :: use_LW_jacobian, damp_LW_fluxadj, & + & pert_radtend logical, intent(in) :: do_sppt,ca_global real(kind=kind_phys), intent(in) :: solhr, slag, cdec, sdec, & - & deltim, fhswr, minGPpres + & deltim, fhswr, minGPpres, lfnc_k, lfnc_p0 real(kind=kind_phys), dimension(:), intent(in) :: & & sinlat, coslat, xlon, coszen, tf, tsflw, sfcdlw, & @@ -253,11 +254,19 @@ subroutine dcyc2t3_run & integer, intent(out) :: errflg ! --- locals: - integer :: i, k, nstp, nstl, it, istsun(im),iSFC + integer :: i, k, nstp, nstl, it, istsun(im),iSFC,iTOA real(kind=kind_phys) :: cns, coszn, tem1, tem2, anginc, & & rstl, solang, dT real(kind=kind_phys), dimension(im,levs+1) :: flxlwup_adj, & & flxlwdn_adj, t_lev2 + real(kind=kind_phys) :: fluxlwnet_adj,fluxlwnet,dT_sfc, & + &fluxlwDOWN_jac,lfnc,c1 + ! Length scale for flux-adjustment scaling + real(kind=kind_phys), parameter :: & + & L = 1. + ! Scaling factor for downwelling LW Jacobian profile. + real(kind=kind_phys), parameter :: & + & gamma = 0.2 ! !===> ... begin here ! @@ -267,9 +276,11 @@ subroutine dcyc2t3_run & ! Vertical ordering? if (p_lev(1,1) .lt. p_lev(1, levs)) then - iSFC = levs + iSFC = levs + 1 + iTOA = 1 else iSFC = 1 + iTOA = levs + 1 endif tem1 = fhswr / deltim @@ -375,32 +386,47 @@ subroutine dcyc2t3_run & call cmp_tlev(im, levs, minGPpres, p_lay, t_lay, p_lev, tsfc, & & t_lev2) + ! Compute adjusted net LW flux foillowing Hogan and Bozzo 2015 (10.1002/2015MS000455) + ! Here we assume that the profile of the downwelling LW Jaconiam has the same shape + ! as the upwelling, but scaled and offset. + ! The scaling factor is 0.2 + ! The profile of the downwelling Jacobian (J) is offset so that + ! J_dn_sfc / J_up_sfc = scaling_factor + ! J_dn_toa / J_up_sfc = 0 ! - ! Adjust up/downward fluxes (at layer interfaces). - ! - do k = 1, levs+1 - do i = 1, im - dT = t_lev2(i,k) - t_lev(i,k) - flxlwup_adj(i,k) = flux2D_lwUP(i,k) + & - & fluxlwUP_jac(i,k)*dT - enddo - enddo - ! - ! Compute new heating rate (within each layer). - ! - do k = 1, levs - htrlw(1:im,k) = & - & (flxlwup_adj(1:im,k+1) - flxlwup_adj(1:im,k) - & - & flux2D_lwDOWN(1:im,k+1) + flux2D_lwDOWN(1:im,k)) * & - & con_g / (con_cp * (p_lev(1:im,k+1) - p_lev(1:im,k))) - enddo - - ! - ! Add radiative heating rates to physics heating rate - ! - do k = 1, levs - do i = 1, im - dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + htrlw(i,k) + ! Optionally, the flux adjustment can be damped with height using a logistic function + ! fx ~ L / (1 + exp(-k*dp)), where dp = p - p0 + ! L = 1, fix scale between 0-1. - Fixed + ! k = 1 / pressure decay length (Pa) - Controlled by namelist + ! p0 = Transition pressure (Pa) - Controlled by namelsit + do i = 1, im + c1 = fluxlwUP_jac(i,iTOA) / fluxlwUP_jac(i,iSFC) + dT_sfc = t_lev2(i,iSFC) - t_lev(i,iSFC) + do k = 1, levs + ! LW net flux + fluxlwnet = (flux2D_lwUP(i, k+1) - flux2D_lwUP(i, k) - & + & flux2D_lwDOWN(i,k+1) + flux2D_lwDOWN(i,k)) + ! Downward LW Jacobian (Eq. 9) + fluxlwDOWN_jac = gamma * & + & (fluxlwUP_jac(i,k)/fluxlwUP_jac(i,iSFC) - c1) / & + & (1 - c1) + ! Adjusted LW net flux(Eq. 10) + fluxlwnet_adj = fluxlwnet + dT_sfc* & + & (fluxlwUP_jac(i,k)/fluxlwUP_jac(i,iSFC) - & + & fluxlwDOWN_jac) + ! Adjusted LW heating rate + htrlw(i,k) = fluxlwnet_adj * con_g / & + & (con_cp * (p_lev(i,k+1) - p_lev(i,k))) + + ! Add radiative heating rates to physics heating rate. Optionally, scaled w/ height + ! using a logistic function + if (damp_LW_fluxadj) then + lfnc = L / (1+exp(-(p_lev(i,k) - lfnc_p0)/lfnc_k)) + else + lfnc = 1. + endif + dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + & + & htrlw(i,k)*lfnc + (1.-lfnc)*hlw(i,k) enddo enddo else diff --git a/physics/dcyc2.meta b/physics/dcyc2.meta index a460db7ab..91e01a2d2 100644 --- a/physics/dcyc2.meta +++ b/physics/dcyc2.meta @@ -370,6 +370,32 @@ type = logical intent = in optional = F +[damp_LW_fluxadj] + standard_name = flag_to_damp_RRTMGP_LW_jacobian_flux_adjustment + long_name = logical flag to control RRTMGP LW calculation + units = flag + dimensions = () + type = logical + intent = in + optional = F +[lfnc_k] + standard_name = transition_pressure_length_scale_for_flux_damping + long_name = depth of transition layer in logistic function for LW flux adjustment damping + units = Pa + dimensions = () + type = real + kind = kind_phys + intent = in + optional = F +[lfnc_p0] + standard_name = transition_pressure_for_flux_damping + long_name = transition pressure for LW flux adjustment damping + units = Pa + dimensions = () + type = real + kind = kind_phys + intent = in + optional = F [sfculw] standard_name = surface_upwelling_longwave_flux_on_radiation_time_step long_name = total sky sfc upward lw flux