From df30cf7c01d3fcd5f5a97ea4569220d1b3da7daf Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Sun, 9 May 2021 19:32:06 +0000 Subject: [PATCH 01/14] Limit full LW flux profile adjustment to below 100hPa. --- physics/dcyc2.f | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index ad9365851..a3d7cf193 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -380,9 +380,12 @@ subroutine dcyc2t3_run & ! 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 + flxlwup_adj(i,k) = flux2D_lwUP(i,k) + if (p_lev(i,k) .gt. 10000.) then + dT = t_lev2(i,k) - t_lev(i,k) + flxlwup_adj(i,k) = flux2D_lwUP(i,k) + & + & fluxlwUP_jac(i,k)*dT + endif enddo enddo ! From 6961b10546035b55021132e62ce139333eed9cc0 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Sun, 9 May 2021 23:50:22 +0000 Subject: [PATCH 02/14] Added more safeguards against out-of-bounds temperature to GP inputs. --- physics/GFS_rrtmgp_pre.F90 | 12 ++++++++---- physics/GFS_rrtmgp_pre.meta | 9 +++++++++ physics/dcyc2.f | 9 +++++---- physics/dcyc2.meta | 18 ++++++++++++++++++ physics/radiation_tools.F90 | 20 ++++++++++++++------ physics/rrtmgp_lw_gas_optics.F90 | 4 +++- physics/rrtmgp_lw_gas_optics.meta | 9 +++++++++ 7 files changed, 66 insertions(+), 15 deletions(-) diff --git a/physics/GFS_rrtmgp_pre.F90 b/physics/GFS_rrtmgp_pre.F90 index 88e534595..af7e5f1a0 100644 --- a/physics/GFS_rrtmgp_pre.F90 +++ b/physics/GFS_rrtmgp_pre.F90 @@ -98,8 +98,8 @@ end subroutine GFS_rrtmgp_pre_init !! subroutine GFS_rrtmgp_pre_run(nCol, nLev, nTracers, i_o3, lsswr, lslwr, fhswr, fhlwr, & xlat, xlon, prsl, tgrs, prslk, prsi, qgrs, tsfc, con_eps, con_epsm1, con_fvirt, & - con_epsqs, minGPpres, minGPtemp, raddt, p_lay, t_lay, p_lev, t_lev, tsfg, tsfa, & - qs_lay, q_lay, tv_lay, relhum, tracer, gas_concentrations, errmsg, errflg) + con_epsqs, minGPpres, minGPtemp, maxGPtemp, raddt, p_lay, t_lay, p_lev, t_lev, tsfg, & + tsfa, qs_lay, q_lay, tv_lay, relhum, tracer, gas_concentrations, errmsg, errflg) ! Inputs integer, intent(in) :: & @@ -112,6 +112,7 @@ subroutine GFS_rrtmgp_pre_run(nCol, nLev, nTracers, i_o3, lsswr, lslwr, fhswr, f lslwr ! Call LW radiation real(kind_phys), intent(in) :: & minGPtemp, & ! Minimum temperature allowed in RRTMGP. + maxGPtemp, & ! Maximum temperature allowed in RRTMGP. minGPpres, & ! Minimum pressure allowed in RRTMGP. fhswr, & ! Frequency of SW radiation call. fhlwr ! Frequency of LW radiation call. @@ -208,11 +209,14 @@ subroutine GFS_rrtmgp_pre_run(nCol, nLev, nTracers, i_o3, lsswr, lslwr, fhswr, f if (t_lay(iCol,iLay) .le. minGPtemp) then t_lay(iCol,iLay) = minGPtemp + epsilon(minGPtemp) endif + if (t_lay(iCol,iLay) .ge. maxGPtemp) then + t_lay(iCol,iLay) = maxGPtemp - epsilon(maxGPtemp) + endif enddo enddo ! Temperature at layer-interfaces - call cmp_tlev(nCol,nLev,minGPpres,p_lay,t_lay,p_lev,tsfc,t_lev) + call cmp_tlev(nCol,nLev,minGPpres,minGPtemp,maxGPtemp,p_lay,t_lay,p_lev,tsfc,t_lev) ! Compute a bunch of thermodynamic fields needed by the cloud microphysics schemes. ! Relative humidity, saturation mixing-ratio, vapor mixing-ratio, virtual temperature, @@ -273,7 +277,7 @@ subroutine GFS_rrtmgp_pre_run(nCol, nLev, nTracers, i_o3, lsswr, lslwr, fhswr, f ! ####################################################################################### ! Setup surface ground temperature and ground/air skin temperature if required. ! ####################################################################################### - tsfg(1:NCOL) = tsfc(1:NCOL) + tsfg(1:NCOL) = t_lev(1:NCOL,iSFC) tsfa(1:NCOL) = t_lay(1:NCOL,iSFC) end subroutine GFS_rrtmgp_pre_run diff --git a/physics/GFS_rrtmgp_pre.meta b/physics/GFS_rrtmgp_pre.meta index 8096aef2a..895bbc630 100644 --- a/physics/GFS_rrtmgp_pre.meta +++ b/physics/GFS_rrtmgp_pre.meta @@ -239,6 +239,15 @@ kind = kind_phys intent = in optional = F +[maxGPtemp] + standard_name = maximum_temperature_in_RRTMGP + long_name = maximum temperature allowed in RRTMGP + units = K + dimensions = () + type = real + kind = kind_phys + intent = in + optional = F [raddt] standard_name = time_step_for_radiation long_name = radiation time step diff --git a/physics/dcyc2.f b/physics/dcyc2.f index a3d7cf193..4678efa0b 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, minGPtemp, maxGPtemp, & & 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, & @@ -216,7 +216,8 @@ subroutine dcyc2t3_run & logical, intent(in) :: use_LW_jacobian, 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, & + & minGPtemp, maxGPtemp real(kind=kind_phys), dimension(:), intent(in) :: & & sinlat, coslat, xlon, coszen, tf, tsflw, sfcdlw, & @@ -372,8 +373,8 @@ subroutine dcyc2t3_run & ! ! Compute temperatute at level interfaces. ! - call cmp_tlev(im, levs, minGPpres, p_lay, t_lay, p_lev, tsfc, & - & t_lev2) + call cmp_tlev(im, levs, minGPpres, minGPtemp, maxGPtemp, p_lay,& + & t_lay, p_lev, tsfc, t_lev2) ! ! Adjust up/downward fluxes (at layer interfaces). diff --git a/physics/dcyc2.meta b/physics/dcyc2.meta index a460db7ab..25b06cc83 100644 --- a/physics/dcyc2.meta +++ b/physics/dcyc2.meta @@ -353,6 +353,24 @@ type = logical intent = in optional = F +[minGPtemp] + standard_name = minimum_temperature_in_RRTMGP + long_name = minimum temperature allowed in RRTMGP + units = K + dimensions = () + type = real + kind = kind_phys + intent = in + optional = F +[maxGPtemp] + standard_name = maximum_temperature_in_RRTMGP + long_name = maximum temperature allowed in RRTMGP + units = K + dimensions = () + type = real + kind = kind_phys + intent = in + optional = F [minGPpres] standard_name = minimum_pressure_in_RRTMGP long_name = minimum pressure allowed in RRTMGP diff --git a/physics/radiation_tools.F90 b/physics/radiation_tools.F90 index c6524aab6..a8d3f5457 100644 --- a/physics/radiation_tools.F90 +++ b/physics/radiation_tools.F90 @@ -2,20 +2,16 @@ module radiation_tools use machine, only: & kind_phys ! Working type implicit none - - real(kind_phys) :: & - rrtmgp_minP, & ! Minimum pressure allowed in RRTMGP - rrtmgp_minT ! Minimum temperature allowed in RRTMGP contains ! ######################################################################################### ! ######################################################################################### - subroutine cmp_tlev(nCol,nLev,minP,p_lay,t_lay,p_lev,tsfc,t_lev) + subroutine cmp_tlev(nCol,nLev,minP,minT,maxT,p_lay,t_lay,p_lev,tsfc,t_lev) ! Inputs integer, intent(in) :: & nCol,nLev real(kind_phys),intent(in) :: & - minP + minP,minT,maxT real(kind_phys),dimension(nCol),intent(in) :: & tsfc real(kind_phys),dimension(nCol,nLev),intent(in) :: & @@ -78,6 +74,18 @@ subroutine cmp_tlev(nCol,nLev,minP,p_lay,t_lay,p_lev,tsfc,t_lev) t_lev(1:NCOL,iTOA+1) = t_lay(1:NCOL,iTOA) endif + ! Bound temperature at layer interfaces + do iCol=1,NCOL + do iLay=1,nLev+1 + if (t_lev(iCol,iLay) .le. minT) then + t_lev(iCol,iLay) = minT + epsilon(minT) + endif + if (t_lev(iCol,iLay) .ge. maxT) then + t_lev(iCol,iLay) = maxT - epsilon(maxT) + endif + enddo + enddo + end subroutine cmp_tlev ! ######################################################################################### diff --git a/physics/rrtmgp_lw_gas_optics.F90 b/physics/rrtmgp_lw_gas_optics.F90 index a116ad772..d7201e026 100644 --- a/physics/rrtmgp_lw_gas_optics.F90 +++ b/physics/rrtmgp_lw_gas_optics.F90 @@ -76,7 +76,7 @@ module rrtmgp_lw_gas_optics !! \htmlinclude rrtmgp_lw_gas_optics_init.html !! subroutine rrtmgp_lw_gas_optics_init(rrtmgp_root_dir, rrtmgp_lw_file_gas, mpicomm, & - mpirank, mpiroot, gas_concentrations, minGPpres, minGPtemp, errmsg, errflg) + mpirank, mpiroot, gas_concentrations, minGPpres, minGPtemp, maxGPtemp, errmsg, errflg) ! Inputs type(ty_gas_concs), intent(inout) :: & @@ -96,6 +96,7 @@ subroutine rrtmgp_lw_gas_optics_init(rrtmgp_root_dir, rrtmgp_lw_file_gas, mpicom errflg ! CCPP error code real(kind_phys), intent(out) :: & minGPtemp, & ! Minimum temperature allowed by RRTMGP. + maxGPtemp, & ! Maximum temperature allowed by RRTMG. minGPpres ! Minimum pressure allowed by RRTMGP. ! Local variables @@ -450,6 +451,7 @@ subroutine rrtmgp_lw_gas_optics_init(rrtmgp_root_dir, rrtmgp_lw_file_gas, mpicom ! temperature (GFS_rrtmgp_pre.F90) minGPpres = lw_gas_props%get_press_min() minGPtemp = lw_gas_props%get_temp_min() + maxGPtemp = lw_gas_props%get_temp_max() end subroutine rrtmgp_lw_gas_optics_init diff --git a/physics/rrtmgp_lw_gas_optics.meta b/physics/rrtmgp_lw_gas_optics.meta index c92567e14..823501cfa 100644 --- a/physics/rrtmgp_lw_gas_optics.meta +++ b/physics/rrtmgp_lw_gas_optics.meta @@ -92,6 +92,15 @@ kind = kind_phys intent = out optional = F +[maxGPtemp] + standard_name = maximum_temperature_in_RRTMGP + long_name = maximum temperature allowed in RRTMGP + units = K + dimensions = () + type = real + kind = kind_phys + intent = out + optional = F ######################################################################## [ccpp-arg-table] From 6ebe85ef4a4d51a8e8bf7921e8079d395f90d89a Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Tue, 11 May 2021 18:40:02 +0000 Subject: [PATCH 03/14] Apply constant heating-rate adjustment above 100hPa --- physics/dcyc2.f | 71 ++++++++++++++++++++++++++++++------------------- 1 file changed, 44 insertions(+), 27 deletions(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index 4678efa0b..368272ff1 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -254,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,c1,c2,c3,c4,c5,ku + ! Pressure limit for LW flux adjustment + real(kind=kind_phys), parameter :: & + & plim_fluxAdj_upper = 10000. + ! Scaling factor for downwelling LW Jacobian profile. + real(kind=kind_phys), parameter :: & + & c0 = 0.2 + logical :: init_lev ! !===> ... begin here ! @@ -268,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 @@ -376,34 +386,41 @@ subroutine dcyc2t3_run & call cmp_tlev(im, levs, minGPpres, minGPtemp, maxGPtemp, 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 - flxlwup_adj(i,k) = flux2D_lwUP(i,k) - if (p_lev(i,k) .gt. 10000.) then - dT = t_lev2(i,k) - t_lev(i,k) - flxlwup_adj(i,k) = flux2D_lwUP(i,k) + & - & fluxlwUP_jac(i,k)*dT + do i = 1, im + c1 = fluxlwUP_jac(i,iSFC) + c2 = fluxlwUP_jac(i,iTOA) / c1 + c3 = t_lev2(i,iSFC) - t_lev(i,iSFC) + init_lev = .true. + do k = 1, levs + ! Only apply the Jacobian adjustment below plim_fluxAdj_upper + if (p_lev(i,k) .gt. plim_fluxAdj_upper) then + c4 = fluxlwUP_jac(i,k)/c1 + fluxlwnet = (flux2D_lwUP(i,k+1) - flux2D_lwUP(i,k) -& + & flux2D_lwDOWN(i,k+1) + flux2D_lwDOWN(i,k)) + ! (Eq. 9) + c5 = c0 * (c4 - c2) / (1 - c2) + ! (Eq. 10) + fluxlwnet_adj = fluxlwnet + c3*(c4-c5) + ! Compute adjusted heating rate + htrlw(i,k) = fluxlwnet_adj * con_g / & + & (con_cp * (p_lev(i,k+1) - p_lev(i,k))) + + ! Store vertical index for plim_fluxAdj_upper + ku = k + ! Above, offset the heating rate by he same amount as in plim_fluxAdj_upper + else + htrlw(i,k) = hlw(i,k) + (htrlw(i,ku)-hlw(i,ku)) endif - 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 + ! Add radiative heating rates to physics heating rate dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + htrlw(i,k) enddo enddo From 230d479e58cb380de9d24e8b67778db468250ba3 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Tue, 11 May 2021 19:47:23 +0000 Subject: [PATCH 04/14] Add vetical decay to impact of flux adjustment above threshold. --- physics/dcyc2.f | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index 368272ff1..f671cf1f2 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -266,7 +266,6 @@ subroutine dcyc2t3_run & ! Scaling factor for downwelling LW Jacobian profile. real(kind=kind_phys), parameter :: & & c0 = 0.2 - logical :: init_lev ! !===> ... begin here ! @@ -398,7 +397,6 @@ subroutine dcyc2t3_run & c1 = fluxlwUP_jac(i,iSFC) c2 = fluxlwUP_jac(i,iTOA) / c1 c3 = t_lev2(i,iSFC) - t_lev(i,iSFC) - init_lev = .true. do k = 1, levs ! Only apply the Jacobian adjustment below plim_fluxAdj_upper if (p_lev(i,k) .gt. plim_fluxAdj_upper) then @@ -417,7 +415,8 @@ subroutine dcyc2t3_run & ku = k ! Above, offset the heating rate by he same amount as in plim_fluxAdj_upper else - htrlw(i,k) = hlw(i,k) + (htrlw(i,ku)-hlw(i,ku)) + htrlw(i,k) = hlw(i,k)+(p_lev(i,k)/plim_fluxAdj_upper)*& + & (htrlw(i,ku)-hlw(i,ku)) endif ! Add radiative heating rates to physics heating rate From 5f7d6970b7d420601dc86c2db91d953feebe1a7b Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Mon, 24 May 2021 16:56:23 +0000 Subject: [PATCH 05/14] Added logistic function to damp the LW flux adjustment with height --- physics/dcyc2.f | 64 ++++++++++++++++++++++++---------------------- physics/dcyc2.meta | 26 +++++++++++++++++++ 2 files changed, 60 insertions(+), 30 deletions(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index f671cf1f2..09c80a97e 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -178,10 +178,10 @@ subroutine dcyc2t3_run & & sfcnirbmu,sfcnirdfu,sfcvisbmu,sfcvisdfu, & & sfcnirbmd,sfcnirdfd,sfcvisbmd,sfcvisdfd, & & im, levs, deltim, fhswr, & - & dry, icy, wet, minGPtemp, maxGPtemp, & - & 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, & + & dry, icy, wet, minGPtemp, maxGPtemp, damp_LW_fluxadj, & + & lfnc_k_grad, 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, & ! & dry, icy, wet, lprnt, ipr, & ! --- input/output: & dtdt,dtdtnp,htrlw, & @@ -213,11 +213,12 @@ 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, & - & minGPtemp, maxGPtemp + & deltim, fhswr, minGPpres, minGPtemp, maxGPtemp, lfnc_k_grad, & + & lfnc_p0 real(kind=kind_phys), dimension(:), intent(in) :: & & sinlat, coslat, xlon, coszen, tf, tsflw, sfcdlw, & @@ -259,10 +260,11 @@ subroutine dcyc2t3_run & & rstl, solang, dT real(kind=kind_phys), dimension(im,levs+1) :: flxlwup_adj, & & flxlwdn_adj, t_lev2 - real(kind=kind_phys) :: fluxlwnet_adj,fluxlwnet,c1,c2,c3,c4,c5,ku - ! Pressure limit for LW flux adjustment + real(kind=kind_phys) :: fluxlwnet_adj,fluxlwnet,c1,c2,c3,c4,c5, & + & dP,lfnc + ! Length scale for flux-adjustment scaling real(kind=kind_phys), parameter :: & - & plim_fluxAdj_upper = 10000. + & L = 1. ! Scaling factor for downwelling LW Jacobian profile. real(kind=kind_phys), parameter :: & & c0 = 0.2 @@ -393,34 +395,36 @@ subroutine dcyc2t3_run & ! J_dn_sfc / J_up_sfc = scaling_factor ! J_dn_toa / J_up_sfc = 0 ! + ! 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. + ! k (steepness) and p0 (midpoint) are controlled via namelist do i = 1, im c1 = fluxlwUP_jac(i,iSFC) c2 = fluxlwUP_jac(i,iTOA) / c1 c3 = t_lev2(i,iSFC) - t_lev(i,iSFC) do k = 1, levs ! Only apply the Jacobian adjustment below plim_fluxAdj_upper - if (p_lev(i,k) .gt. plim_fluxAdj_upper) then - c4 = fluxlwUP_jac(i,k)/c1 - fluxlwnet = (flux2D_lwUP(i,k+1) - flux2D_lwUP(i,k) -& - & flux2D_lwDOWN(i,k+1) + flux2D_lwDOWN(i,k)) - ! (Eq. 9) - c5 = c0 * (c4 - c2) / (1 - c2) - ! (Eq. 10) - fluxlwnet_adj = fluxlwnet + c3*(c4-c5) - ! Compute adjusted heating rate - htrlw(i,k) = fluxlwnet_adj * con_g / & - & (con_cp * (p_lev(i,k+1) - p_lev(i,k))) - - ! Store vertical index for plim_fluxAdj_upper - ku = k - ! Above, offset the heating rate by he same amount as in plim_fluxAdj_upper + c4 = fluxlwUP_jac(i,k)/c1 + fluxlwnet = (flux2D_lwUP(i,k+1) - flux2D_lwUP(i,k) - & + & flux2D_lwDOWN(i,k+1) + flux2D_lwDOWN(i,k)) + ! (Eq. 9) + c5 = c0 * (c4 - c2) / (1 - c2) + ! (Eq. 10) + fluxlwnet_adj = fluxlwnet + c3*(c4-c5) + ! Compute adjusted 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 + dp = p_lev(i,k) - lfnc_p0 + lfnc = L / (1+exp(-lfnc_k_grad*exp(1.)*dp/lfnc_p0)) else - htrlw(i,k) = hlw(i,k)+(p_lev(i,k)/plim_fluxAdj_upper)*& - & (htrlw(i,ku)-hlw(i,ku)) + lfnc = 1. endif - - ! Add radiative heating rates to physics heating rate - dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + htrlw(i,k) + dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + htrlw(i,k)*lfnc enddo enddo else diff --git a/physics/dcyc2.meta b/physics/dcyc2.meta index 25b06cc83..dceb9ce77 100644 --- a/physics/dcyc2.meta +++ b/physics/dcyc2.meta @@ -388,6 +388,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_grad] + standard_name = steepness_of_flux_damping + long_name = steepness of logistic function for damping the LW flux adjustment + units = none + dimensions = () + type = real + kind = kind_phys + intent = in + optional = F +[lfnc_p0] + standard_name = midpoint_used_for_flux_damping + long_name = midpoint for damping the LW flux adjustment + 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 From 2f6e70814fdaae4ba33b00c7d5c7ec421a8e69e8 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Mon, 24 May 2021 17:37:55 +0000 Subject: [PATCH 06/14] Reorganized RRTMGP aerosol optics. --- ...tics.F90 => GFS_rrtmgp_aerosol_optics.F90} | 66 ++++--- ...cs.meta => GFS_rrtmgp_aerosol_optics.meta} | 50 ++--- physics/rrtmgp_lw_aerosol_optics.meta | 173 ------------------ physics/rrtmgp_sw_aerosol_optics.F90 | 119 ------------ 4 files changed, 73 insertions(+), 335 deletions(-) rename physics/{rrtmgp_lw_aerosol_optics.F90 => GFS_rrtmgp_aerosol_optics.F90} (59%) rename physics/{rrtmgp_sw_aerosol_optics.meta => GFS_rrtmgp_aerosol_optics.meta} (91%) delete mode 100644 physics/rrtmgp_lw_aerosol_optics.meta delete mode 100644 physics/rrtmgp_sw_aerosol_optics.F90 diff --git a/physics/rrtmgp_lw_aerosol_optics.F90 b/physics/GFS_rrtmgp_aerosol_optics.F90 similarity index 59% rename from physics/rrtmgp_lw_aerosol_optics.F90 rename to physics/GFS_rrtmgp_aerosol_optics.F90 index df0e77163..194c33ba1 100644 --- a/physics/rrtmgp_lw_aerosol_optics.F90 +++ b/physics/GFS_rrtmgp_aerosol_optics.F90 @@ -1,8 +1,8 @@ -module rrtmgp_lw_aerosol_optics +module GFS_rrtmgp_aerosol_optics use machine, only: kind_phys use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp - use mo_optical_props, only: ty_optical_props_1scl - use radiation_tools, only: check_error_msg + use mo_optical_props, only: ty_optical_props_1scl,ty_optical_props_2str + use radiation_tools, only: check_error_msg use rrtmgp_sw_gas_optics, only: sw_gas_props use rrtmgp_lw_gas_optics, only: lw_gas_props use module_radiation_aerosols, only: & @@ -14,34 +14,37 @@ module rrtmgp_lw_aerosol_optics implicit none - public rrtmgp_lw_aerosol_optics_init, rrtmgp_lw_aerosol_optics_run, rrtmgp_lw_aerosol_optics_finalize + public GFS_rrtmgp_aerosol_optics_init, GFS_rrtmgp_aerosol_optics_run, GFS_rrtmgp_aerosol_optics_finalize contains ! ######################################################################################### - ! SUBROUTINE rrtmgp_lw_aerosol_optics_init() + ! SUBROUTINE GFS_rrtmgp_aerosol_optics_init() ! ######################################################################################### - subroutine rrtmgp_lw_aerosol_optics_init() - end subroutine rrtmgp_lw_aerosol_optics_init + subroutine GFS_rrtmgp_aerosol_optics_init() + end subroutine GFS_rrtmgp_aerosol_optics_init ! ######################################################################################### - ! SUBROUTINE rrtmgp_lw_aerosol_optics_run() + ! SUBROUTINE GFS_rrtmgp_aerosol_optics_run() ! ######################################################################################### -!! \section arg_table_rrtmgp_lw_aerosol_optics_run -!! \htmlinclude rrtmgp_lw_aerosol_optics.html +!! \section arg_table_GFS_rrtmgp_aerosol_optics_run +!! \htmlinclude GFS_rrtmgp_aerosol_optics.html !! - subroutine rrtmgp_lw_aerosol_optics_run(doLWrad, nCol, nLev, nTracer, nTracerAer,& - p_lev, p_lay, p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, & - aerodp, lw_optical_props_aerosol, errmsg, errflg) + subroutine GFS_rrtmgp_aerosol_optics_run(doLWrad, nCol, nLev, nDay, nTracer, nTracerAer, & + idxday, p_lev, p_lay, p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, & + aerodp, lw_optical_props_aerosol, sw_optical_props_aerosol, errmsg, errflg) ! Inputs logical, intent(in) :: & doLWrad ! Logical flag for longwave radiation call integer, intent(in) :: & nCol, & ! Number of horizontal grid points + nDay, & ! Number of daylit points nLev, & ! Number of vertical layers nTracer, & ! Number of tracers nTracerAer ! Number of aerosol tracers + integer,intent(in),dimension(:) :: & + idxday ! Daylit point indices real(kind_phys), dimension(:), intent(in) :: & lon, & ! Longitude lat, & ! Latitude @@ -63,6 +66,8 @@ subroutine rrtmgp_lw_aerosol_optics_run(doLWrad, nCol, nLev, nTracer, nTracerAer aerodp ! Vertical integrated optical depth for various aerosol species type(ty_optical_props_1scl),intent(inout) :: & lw_optical_props_aerosol ! RRTMGP DDT: Longwave aerosol optical properties (tau) + type(ty_optical_props_2str),intent(out) :: & + sw_optical_props_aerosol ! RRTMGP DDT: Shortwave aerosol optical properties (tau,ssa,omega) integer, intent(out) :: & errflg ! CCPP error flag character(len=*), intent(out) :: & @@ -82,23 +87,40 @@ subroutine rrtmgp_lw_aerosol_optics_run(doLWrad, nCol, nLev, nTracer, nTracerAer if (.not. doLWrad) return ! Call module_radiation_aerosols::setaer(),to setup aerosols property profile - call setaer(p_lev/100., p_lay/100., p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, ncol, nLev, & - nLev+1, .true., .true., aerosolssw, aerosolslw, aerodp) + call setaer(p_lev/100., p_lay/100., p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, & + lat, ncol, nLev, nLev+1, .true., .true., aerosolssw, aerosolslw, aerodp) - ! Copy aerosol optical information to RRTMGP DDT + ! Copy aerosol optical information to RRTMGP DDTs + ! + ! LW + ! lw_optical_props_aerosol%tau = aerosolslw(:,:,:,1) * (1. - aerosolslw(:,:,:,2)) - lw_optical_props_aerosol%band_lims_wvn = lw_gas_props%get_band_lims_wavenumber() do iBand=1,lw_gas_props%get_nband() lw_optical_props_aerosol%band2gpt(1:2,iBand) = iBand lw_optical_props_aerosol%gpt2band(iBand) = iBand end do + ! + ! SW + ! + if (nDay .gt. 0) then + ! Allocate RRTMGP DDT + call check_error_msg('rrtmgp_sw_aerosol_optics_run',sw_optical_props_aerosol%alloc_2str( & + nDay, nlev, sw_gas_props%get_band_lims_wavenumber())) + ! Copy + sw_optical_props_aerosol%tau(:,:,1) = aerosolssw(idxday(1:nDay),:,sw_gas_props%get_nband(), 1) + sw_optical_props_aerosol%tau(:,:,2:sw_gas_props%get_nband()-1) = aerosolssw(idxday(1:nDay),:,1:sw_gas_props%get_nband()-1,1) + sw_optical_props_aerosol%ssa(:,:,1) = aerosolssw(idxday(1:nDay),:,sw_gas_props%get_nband(), 2) + sw_optical_props_aerosol%ssa(:,:,2:sw_gas_props%get_nband()-1) = aerosolssw(idxday(1:nDay),:,1:sw_gas_props%get_nband()-1,2) + sw_optical_props_aerosol%g(:,:,1) = aerosolssw(idxday(1:nDay),:,sw_gas_props%get_nband(), 3) + sw_optical_props_aerosol%g(:,:,2:sw_gas_props%get_nband()-1) = aerosolssw(idxday(1:nDay),:,1:sw_gas_props%get_nband()-1,3) + endif - end subroutine rrtmgp_lw_aerosol_optics_run + end subroutine GFS_rrtmgp_aerosol_optics_run ! ######################################################################################### - ! SUBROUTINE rrtmgp_lw_aerosol_optics_finalize() + ! SUBROUTINE GFS_rrtmgp_aerosol_optics_finalize() ! ######################################################################################### - subroutine rrtmgp_lw_aerosol_optics_finalize() - end subroutine rrtmgp_lw_aerosol_optics_finalize -end module rrtmgp_lw_aerosol_optics + subroutine GFS_rrtmgp_aerosol_optics_finalize() + end subroutine GFS_rrtmgp_aerosol_optics_finalize +end module GFS_rrtmgp_aerosol_optics diff --git a/physics/rrtmgp_sw_aerosol_optics.meta b/physics/GFS_rrtmgp_aerosol_optics.meta similarity index 91% rename from physics/rrtmgp_sw_aerosol_optics.meta rename to physics/GFS_rrtmgp_aerosol_optics.meta index f4909c794..aa7f6d4a5 100644 --- a/physics/rrtmgp_sw_aerosol_optics.meta +++ b/physics/GFS_rrtmgp_aerosol_optics.meta @@ -1,15 +1,15 @@ [ccpp-table-properties] - name = rrtmgp_sw_aerosol_optics + name = GFS_rrtmgp_aerosol_optics type = scheme dependencies = iounitdef.f,machine.F,radiation_aerosols.f,radiation_tools.F90 ######################################################################## [ccpp-arg-table] - name = rrtmgp_sw_aerosol_optics_run + name = GFS_rrtmgp_aerosol_optics_run type = scheme -[doSWrad] - standard_name = flag_to_calc_sw - long_name = logical flags for sw radiation calls +[doLWrad] + standard_name = flag_to_calc_lw + long_name = logical flags for lw radiation calls units = flag dimensions = () type = logical @@ -23,6 +23,22 @@ type = integer intent = in optional = F +[nday] + standard_name = daytime_points_dimension + long_name = daytime points dimension + units = count + dimensions = () + type = integer + intent = in + optional = F +[idxday] + standard_name = daytime_points + long_name = daytime points + units = index + dimensions = (horizontal_loop_extent) + type = integer + intent = in + optional = F [nLev] standard_name = vertical_dimension long_name = number of vertical levels @@ -47,22 +63,6 @@ type = integer intent = in optional = F -[nday] - standard_name = daytime_points_dimension - long_name = daytime points dimension - units = count - dimensions = () - type = integer - intent = in - optional = F -[idxday] - standard_name = daytime_points - long_name = daytime points - units = index - dimensions = (horizontal_loop_extent) - type = integer - intent = in - optional = F [p_lev] standard_name = air_pressure_at_interface_for_RRTMGP_in_hPa long_name = air pressure at vertical interface for radiation calculation @@ -162,6 +162,14 @@ kind = kind_phys intent = inout optional = F +[lw_optical_props_aerosol] + standard_name = longwave_optical_properties_for_aerosols + long_name = Fortran DDT containing RRTMGP optical properties + units = DDT + dimensions = () + type = ty_optical_props_1scl + intent = inout + optional = F [sw_optical_props_aerosol] standard_name = shortwave_optical_properties_for_aerosols long_name = Fortran DDT containing RRTMGP optical properties diff --git a/physics/rrtmgp_lw_aerosol_optics.meta b/physics/rrtmgp_lw_aerosol_optics.meta deleted file mode 100644 index ad68fd546..000000000 --- a/physics/rrtmgp_lw_aerosol_optics.meta +++ /dev/null @@ -1,173 +0,0 @@ -[ccpp-table-properties] - name = rrtmgp_lw_aerosol_optics - type = scheme - dependencies = iounitdef.f,machine.F,radiation_aerosols.f,radiation_tools.F90 - -######################################################################## -[ccpp-arg-table] - name = rrtmgp_lw_aerosol_optics_run - type = scheme -[doLWrad] - standard_name = flag_to_calc_lw - long_name = logical flags for lw radiation calls - units = flag - dimensions = () - type = logical - intent = in - optional = F -[ncol] - standard_name = horizontal_loop_extent - long_name = horizontal dimension - units = count - dimensions = () - type = integer - intent = in - optional = F -[nLev] - standard_name = vertical_dimension - long_name = number of vertical levels - units = count - dimensions = () - type = integer - intent = in - optional = F -[nTracer] - standard_name = number_of_tracers - long_name = number of tracers - units = count - dimensions = () - type = integer - intent = in - optional = F -[nTracerAer] - standard_name = number_of_aerosol_tracers_MG - long_name = number of aerosol tracers for Morrison Gettelman MP - units = count - dimensions = () - type = integer - intent = in - optional = F -[p_lev] - standard_name = air_pressure_at_interface_for_RRTMGP_in_hPa - long_name = air pressure at vertical interface for radiation calculation - units = hPa - dimensions = (horizontal_loop_extent,vertical_dimension_plus_one) - type = real - kind = kind_phys - intent = in - optional = F -[p_lay] - standard_name = air_pressure_at_layer_for_RRTMGP_in_hPa - long_name = air pressure at vertical layer for radiation calculation - units = hPa - dimensions = (horizontal_loop_extent,vertical_dimension) - type = real - kind = kind_phys - intent = in - optional = F -[p_lk] - standard_name = dimensionless_exner_function_at_model_layers - long_name = dimensionless Exner function at model layer centers - units = none - dimensions = (horizontal_loop_extent,vertical_dimension) - type = real - kind = kind_phys - intent = in - optional = F -[tv_lay] - standard_name = virtual_temperature - long_name = layer virtual temperature - units = K - dimensions = (horizontal_loop_extent,vertical_dimension) - type = real - kind = kind_phys - intent = in - optional = F -[relhum] - standard_name = relative_humidity - long_name = layer relative humidity - units = frac - dimensions = (horizontal_loop_extent,vertical_dimension) - type = real - kind = kind_phys - intent = in - optional = F -[lsmask] - standard_name = sea_land_ice_mask_real - long_name = landmask: sea/land/ice=0/1/2 - units = flag - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = in - optional = F -[tracer] - standard_name = chemical_tracers - long_name = chemical tracers - units = g g-1 - dimensions = (horizontal_loop_extent,vertical_dimension,number_of_tracers) - type = real - kind = kind_phys - intent = in - optional = F -[aerfld] - standard_name = aerosol_number_concentration_from_gocart_aerosol_climatology - long_name = GOCART aerosol climatology number concentration - units = kg-1 - dimensions = (horizontal_loop_extent,vertical_dimension,number_of_aerosol_tracers_MG) - type = real - kind = kind_phys - intent = in - optional = F -[lon] - standard_name = longitude - long_name = longitude - units = radian - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = in - optional = F -[lat] - standard_name = latitude - long_name = latitude - units = radian - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = in - optional = F -[aerodp] - standard_name = atmosphere_optical_thickness_due_to_ambient_aerosol_particles - long_name = vertical integrated optical depth for various aerosol species - units = none - dimensions = (horizontal_loop_extent,number_of_species_for_aerosol_optical_depth) - type = real - kind = kind_phys - intent = inout - optional = F -[lw_optical_props_aerosol] - standard_name = longwave_optical_properties_for_aerosols - long_name = Fortran DDT containing RRTMGP optical properties - units = DDT - dimensions = () - type = ty_optical_props_1scl - intent = inout - optional = F -[errmsg] - standard_name = ccpp_error_message - long_name = error message for error handling in CCPP - units = none - dimensions = () - type = character - kind = len=* - intent = out - optional = F -[errflg] - standard_name = ccpp_error_flag - long_name = error flag for error handling in CCPP - units = flag - dimensions = () - type = integer - intent = out - optional = F diff --git a/physics/rrtmgp_sw_aerosol_optics.F90 b/physics/rrtmgp_sw_aerosol_optics.F90 deleted file mode 100644 index 3a74771b7..000000000 --- a/physics/rrtmgp_sw_aerosol_optics.F90 +++ /dev/null @@ -1,119 +0,0 @@ -module rrtmgp_sw_aerosol_optics - use machine, only: kind_phys - use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp - use mo_optical_props, only: ty_optical_props_2str - use radiation_tools, only: check_error_msg - use rrtmgp_sw_gas_optics, only: sw_gas_props - use rrtmgp_lw_gas_optics, only: lw_gas_props - use module_radiation_aerosols, only: & - NF_AESW, & ! Number of optical-fields in SW output (3=tau+g+omega) - NF_AELW, & ! Number of optical-fields in LW output (3=tau+g+omega) - setaer, & ! Routine to compute aerosol radiative properties (tau,g,omega) - NSPC1 ! Number of species for vertically integrated aerosol optical-depth - use netcdf - - implicit none - - public rrtmgp_sw_aerosol_optics_init, rrtmgp_sw_aerosol_optics_run, rrtmgp_sw_aerosol_optics_finalize - -contains - - ! ######################################################################################### - ! SUBROUTINE rrtmgp_sw_aerosol_optics_init() - ! ######################################################################################### - subroutine rrtmgp_sw_aerosol_optics_init() - end subroutine rrtmgp_sw_aerosol_optics_init - - ! ######################################################################################### - ! SUBROUTINE rrtmgp_sw_aerosol_optics_run() - ! ######################################################################################### -!! \section arg_table_rrtmgp_sw_aerosol_optics_run -!! \htmlinclude rrtmgp_sw_aerosol_optics_run.html -!! - subroutine rrtmgp_sw_aerosol_optics_run(doSWrad, nCol, nLev, nTracer, nTracerAer, nDay, & - idxday, p_lev, p_lay, p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, & - aerodp, sw_optical_props_aerosol, errmsg, errflg ) - - ! Inputs - logical, intent(in) :: & - doSWrad ! Logical flag for shortwave radiation call - integer, intent(in) :: & - nCol, & ! Number of horizontal grid points - nDay, & ! Number of daylit points - nLev, & ! Number of vertical layers - nTracer, & ! Number of tracers - nTracerAer ! Number of aerosol tracers - integer,intent(in),dimension(:) :: & - idxday ! Indices for daylit points. - real(kind_phys), dimension(:), intent(in) :: & - lon, & ! Longitude - lat, & ! Latitude - lsmask ! Land/sea/sea-ice mask - real(kind_phys), dimension(:,:),intent(in) :: & - p_lay, & ! Pressure @ layer-centers (Pa) - tv_lay, & ! Virtual-temperature @ layer-centers (K) - relhum, & ! Relative-humidity @ layer-centers - p_lk ! Exner function @ layer-centers (1) - real(kind_phys), dimension(:, :,:),intent(in) :: & - tracer ! trace gas concentrations - real(kind_phys), dimension(:, :,:),intent(in) :: & - aerfld ! aerosol input concentrations - real(kind_phys), dimension(:,:),intent(in) :: & - p_lev ! Pressure @ layer-interfaces (Pa) - - ! Outputs - real(kind_phys), dimension(:,:), intent(inout) :: & - aerodp ! Vertical integrated optical depth for various aerosol species - type(ty_optical_props_2str),intent(out) :: & - sw_optical_props_aerosol ! RRTMGP DDT: Longwave aerosol optical properties (tau) - integer, intent(out) :: & - errflg ! CCPP error flag - character(len=*), intent(out) :: & - errmsg ! CCPP error message - - ! Local variables - real(kind_phys), dimension(nCol, nLev, lw_gas_props%get_nband(), NF_AELW) :: & - aerosolslw ! - real(kind_phys), dimension(nCol, nLev, sw_gas_props%get_nband(), NF_AESW) :: & - aerosolssw, aerosolssw2 - - ! Initialize CCPP error handling variables - errmsg = '' - errflg = 0 - - if (.not. doSWrad) return - if (nDay .gt. 0) then - - ! Call module_radiation_aerosols::setaer(),to setup aerosols property profile - call setaer(p_lev/100., p_lay/100., p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, nCol, nLev, & - nLev+1, .true., .true., aerosolssw2, aerosolslw, aerodp) - - ! Store aerosol optical properties - ! SW. - ! For RRTMGP SW the bands are now ordered from [IR(band) -> nIR -> UV], in RRTMG the - ! band ordering was [nIR -> UV -> IR(band)] - aerosolssw(1:nCol,:,1,1) = aerosolssw2(1:nCol,:,sw_gas_props%get_nband(),1) - aerosolssw(1:nCol,:,1,2) = aerosolssw2(1:nCol,:,sw_gas_props%get_nband(),2) - aerosolssw(1:nCol,:,1,3) = aerosolssw2(1:nCol,:,sw_gas_props%get_nband(),3) - aerosolssw(1:nCol,:,2:sw_gas_props%get_nband(),1) = aerosolssw2(1:nCol,:,1:sw_gas_props%get_nband()-1,1) - aerosolssw(1:nCol,:,2:sw_gas_props%get_nband(),2) = aerosolssw2(1:nCol,:,1:sw_gas_props%get_nband()-1,2) - aerosolssw(1:nCol,:,2:sw_gas_props%get_nband(),3) = aerosolssw2(1:nCol,:,1:sw_gas_props%get_nband()-1,3) - - ! Allocate RRTMGP DDT: Aerosol optics [nCol,nlev,nBands] - call check_error_msg('rrtmgp_sw_aerosol_optics_run',sw_optical_props_aerosol%alloc_2str( & - nDay, nlev, sw_gas_props%get_band_lims_wavenumber())) - - ! Copy aerosol optical information to RRTMGP DDT - sw_optical_props_aerosol%tau = aerosolssw(idxday(1:nday),:,:,1) - sw_optical_props_aerosol%ssa = aerosolssw(idxday(1:nday),:,:,2) - sw_optical_props_aerosol%g = aerosolssw(idxday(1:nday),:,:,3) - endif - - end subroutine rrtmgp_sw_aerosol_optics_run - - ! ######################################################################################### - ! SUBROUTINE rrtmgp_sw_aerosol_optics_finalize() - ! ######################################################################################### - subroutine rrtmgp_sw_aerosol_optics_finalize() - end subroutine rrtmgp_sw_aerosol_optics_finalize -end module rrtmgp_sw_aerosol_optics From 7bc877dd5e69cf0f171b1267cf0776fa33d6fb48 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Tue, 25 May 2021 16:47:28 +0000 Subject: [PATCH 07/14] Revert "Reorganized RRTMGP aerosol optics." This reverts commit 2f6e70814fdaae4ba33b00c7d5c7ec421a8e69e8. --- ...ptics.F90 => rrtmgp_lw_aerosol_optics.F90} | 66 +++---- physics/rrtmgp_lw_aerosol_optics.meta | 173 ++++++++++++++++++ physics/rrtmgp_sw_aerosol_optics.F90 | 119 ++++++++++++ ...ics.meta => rrtmgp_sw_aerosol_optics.meta} | 50 +++-- 4 files changed, 335 insertions(+), 73 deletions(-) rename physics/{GFS_rrtmgp_aerosol_optics.F90 => rrtmgp_lw_aerosol_optics.F90} (59%) create mode 100644 physics/rrtmgp_lw_aerosol_optics.meta create mode 100644 physics/rrtmgp_sw_aerosol_optics.F90 rename physics/{GFS_rrtmgp_aerosol_optics.meta => rrtmgp_sw_aerosol_optics.meta} (91%) diff --git a/physics/GFS_rrtmgp_aerosol_optics.F90 b/physics/rrtmgp_lw_aerosol_optics.F90 similarity index 59% rename from physics/GFS_rrtmgp_aerosol_optics.F90 rename to physics/rrtmgp_lw_aerosol_optics.F90 index 194c33ba1..df0e77163 100644 --- a/physics/GFS_rrtmgp_aerosol_optics.F90 +++ b/physics/rrtmgp_lw_aerosol_optics.F90 @@ -1,8 +1,8 @@ -module GFS_rrtmgp_aerosol_optics +module rrtmgp_lw_aerosol_optics use machine, only: kind_phys use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp - use mo_optical_props, only: ty_optical_props_1scl,ty_optical_props_2str - use radiation_tools, only: check_error_msg + use mo_optical_props, only: ty_optical_props_1scl + use radiation_tools, only: check_error_msg use rrtmgp_sw_gas_optics, only: sw_gas_props use rrtmgp_lw_gas_optics, only: lw_gas_props use module_radiation_aerosols, only: & @@ -14,37 +14,34 @@ module GFS_rrtmgp_aerosol_optics implicit none - public GFS_rrtmgp_aerosol_optics_init, GFS_rrtmgp_aerosol_optics_run, GFS_rrtmgp_aerosol_optics_finalize + public rrtmgp_lw_aerosol_optics_init, rrtmgp_lw_aerosol_optics_run, rrtmgp_lw_aerosol_optics_finalize contains ! ######################################################################################### - ! SUBROUTINE GFS_rrtmgp_aerosol_optics_init() + ! SUBROUTINE rrtmgp_lw_aerosol_optics_init() ! ######################################################################################### - subroutine GFS_rrtmgp_aerosol_optics_init() - end subroutine GFS_rrtmgp_aerosol_optics_init + subroutine rrtmgp_lw_aerosol_optics_init() + end subroutine rrtmgp_lw_aerosol_optics_init ! ######################################################################################### - ! SUBROUTINE GFS_rrtmgp_aerosol_optics_run() + ! SUBROUTINE rrtmgp_lw_aerosol_optics_run() ! ######################################################################################### -!! \section arg_table_GFS_rrtmgp_aerosol_optics_run -!! \htmlinclude GFS_rrtmgp_aerosol_optics.html +!! \section arg_table_rrtmgp_lw_aerosol_optics_run +!! \htmlinclude rrtmgp_lw_aerosol_optics.html !! - subroutine GFS_rrtmgp_aerosol_optics_run(doLWrad, nCol, nLev, nDay, nTracer, nTracerAer, & - idxday, p_lev, p_lay, p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, & - aerodp, lw_optical_props_aerosol, sw_optical_props_aerosol, errmsg, errflg) + subroutine rrtmgp_lw_aerosol_optics_run(doLWrad, nCol, nLev, nTracer, nTracerAer,& + p_lev, p_lay, p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, & + aerodp, lw_optical_props_aerosol, errmsg, errflg) ! Inputs logical, intent(in) :: & doLWrad ! Logical flag for longwave radiation call integer, intent(in) :: & nCol, & ! Number of horizontal grid points - nDay, & ! Number of daylit points nLev, & ! Number of vertical layers nTracer, & ! Number of tracers nTracerAer ! Number of aerosol tracers - integer,intent(in),dimension(:) :: & - idxday ! Daylit point indices real(kind_phys), dimension(:), intent(in) :: & lon, & ! Longitude lat, & ! Latitude @@ -66,8 +63,6 @@ subroutine GFS_rrtmgp_aerosol_optics_run(doLWrad, nCol, nLev, nDay, nTracer, nTr aerodp ! Vertical integrated optical depth for various aerosol species type(ty_optical_props_1scl),intent(inout) :: & lw_optical_props_aerosol ! RRTMGP DDT: Longwave aerosol optical properties (tau) - type(ty_optical_props_2str),intent(out) :: & - sw_optical_props_aerosol ! RRTMGP DDT: Shortwave aerosol optical properties (tau,ssa,omega) integer, intent(out) :: & errflg ! CCPP error flag character(len=*), intent(out) :: & @@ -87,40 +82,23 @@ subroutine GFS_rrtmgp_aerosol_optics_run(doLWrad, nCol, nLev, nDay, nTracer, nTr if (.not. doLWrad) return ! Call module_radiation_aerosols::setaer(),to setup aerosols property profile - call setaer(p_lev/100., p_lay/100., p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, & - lat, ncol, nLev, nLev+1, .true., .true., aerosolssw, aerosolslw, aerodp) + call setaer(p_lev/100., p_lay/100., p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, ncol, nLev, & + nLev+1, .true., .true., aerosolssw, aerosolslw, aerodp) - ! Copy aerosol optical information to RRTMGP DDTs - ! - ! LW - ! + ! Copy aerosol optical information to RRTMGP DDT lw_optical_props_aerosol%tau = aerosolslw(:,:,:,1) * (1. - aerosolslw(:,:,:,2)) + lw_optical_props_aerosol%band_lims_wvn = lw_gas_props%get_band_lims_wavenumber() do iBand=1,lw_gas_props%get_nband() lw_optical_props_aerosol%band2gpt(1:2,iBand) = iBand lw_optical_props_aerosol%gpt2band(iBand) = iBand end do - ! - ! SW - ! - if (nDay .gt. 0) then - ! Allocate RRTMGP DDT - call check_error_msg('rrtmgp_sw_aerosol_optics_run',sw_optical_props_aerosol%alloc_2str( & - nDay, nlev, sw_gas_props%get_band_lims_wavenumber())) - ! Copy - sw_optical_props_aerosol%tau(:,:,1) = aerosolssw(idxday(1:nDay),:,sw_gas_props%get_nband(), 1) - sw_optical_props_aerosol%tau(:,:,2:sw_gas_props%get_nband()-1) = aerosolssw(idxday(1:nDay),:,1:sw_gas_props%get_nband()-1,1) - sw_optical_props_aerosol%ssa(:,:,1) = aerosolssw(idxday(1:nDay),:,sw_gas_props%get_nband(), 2) - sw_optical_props_aerosol%ssa(:,:,2:sw_gas_props%get_nband()-1) = aerosolssw(idxday(1:nDay),:,1:sw_gas_props%get_nband()-1,2) - sw_optical_props_aerosol%g(:,:,1) = aerosolssw(idxday(1:nDay),:,sw_gas_props%get_nband(), 3) - sw_optical_props_aerosol%g(:,:,2:sw_gas_props%get_nband()-1) = aerosolssw(idxday(1:nDay),:,1:sw_gas_props%get_nband()-1,3) - endif - end subroutine GFS_rrtmgp_aerosol_optics_run + end subroutine rrtmgp_lw_aerosol_optics_run ! ######################################################################################### - ! SUBROUTINE GFS_rrtmgp_aerosol_optics_finalize() + ! SUBROUTINE rrtmgp_lw_aerosol_optics_finalize() ! ######################################################################################### - subroutine GFS_rrtmgp_aerosol_optics_finalize() - end subroutine GFS_rrtmgp_aerosol_optics_finalize -end module GFS_rrtmgp_aerosol_optics + subroutine rrtmgp_lw_aerosol_optics_finalize() + end subroutine rrtmgp_lw_aerosol_optics_finalize +end module rrtmgp_lw_aerosol_optics diff --git a/physics/rrtmgp_lw_aerosol_optics.meta b/physics/rrtmgp_lw_aerosol_optics.meta new file mode 100644 index 000000000..ad68fd546 --- /dev/null +++ b/physics/rrtmgp_lw_aerosol_optics.meta @@ -0,0 +1,173 @@ +[ccpp-table-properties] + name = rrtmgp_lw_aerosol_optics + type = scheme + dependencies = iounitdef.f,machine.F,radiation_aerosols.f,radiation_tools.F90 + +######################################################################## +[ccpp-arg-table] + name = rrtmgp_lw_aerosol_optics_run + type = scheme +[doLWrad] + standard_name = flag_to_calc_lw + long_name = logical flags for lw radiation calls + units = flag + dimensions = () + type = logical + intent = in + optional = F +[ncol] + standard_name = horizontal_loop_extent + long_name = horizontal dimension + units = count + dimensions = () + type = integer + intent = in + optional = F +[nLev] + standard_name = vertical_dimension + long_name = number of vertical levels + units = count + dimensions = () + type = integer + intent = in + optional = F +[nTracer] + standard_name = number_of_tracers + long_name = number of tracers + units = count + dimensions = () + type = integer + intent = in + optional = F +[nTracerAer] + standard_name = number_of_aerosol_tracers_MG + long_name = number of aerosol tracers for Morrison Gettelman MP + units = count + dimensions = () + type = integer + intent = in + optional = F +[p_lev] + standard_name = air_pressure_at_interface_for_RRTMGP_in_hPa + long_name = air pressure at vertical interface for radiation calculation + units = hPa + dimensions = (horizontal_loop_extent,vertical_dimension_plus_one) + type = real + kind = kind_phys + intent = in + optional = F +[p_lay] + standard_name = air_pressure_at_layer_for_RRTMGP_in_hPa + long_name = air pressure at vertical layer for radiation calculation + units = hPa + dimensions = (horizontal_loop_extent,vertical_dimension) + type = real + kind = kind_phys + intent = in + optional = F +[p_lk] + standard_name = dimensionless_exner_function_at_model_layers + long_name = dimensionless Exner function at model layer centers + units = none + dimensions = (horizontal_loop_extent,vertical_dimension) + type = real + kind = kind_phys + intent = in + optional = F +[tv_lay] + standard_name = virtual_temperature + long_name = layer virtual temperature + units = K + dimensions = (horizontal_loop_extent,vertical_dimension) + type = real + kind = kind_phys + intent = in + optional = F +[relhum] + standard_name = relative_humidity + long_name = layer relative humidity + units = frac + dimensions = (horizontal_loop_extent,vertical_dimension) + type = real + kind = kind_phys + intent = in + optional = F +[lsmask] + standard_name = sea_land_ice_mask_real + long_name = landmask: sea/land/ice=0/1/2 + units = flag + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[tracer] + standard_name = chemical_tracers + long_name = chemical tracers + units = g g-1 + dimensions = (horizontal_loop_extent,vertical_dimension,number_of_tracers) + type = real + kind = kind_phys + intent = in + optional = F +[aerfld] + standard_name = aerosol_number_concentration_from_gocart_aerosol_climatology + long_name = GOCART aerosol climatology number concentration + units = kg-1 + dimensions = (horizontal_loop_extent,vertical_dimension,number_of_aerosol_tracers_MG) + type = real + kind = kind_phys + intent = in + optional = F +[lon] + standard_name = longitude + long_name = longitude + units = radian + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[lat] + standard_name = latitude + long_name = latitude + units = radian + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[aerodp] + standard_name = atmosphere_optical_thickness_due_to_ambient_aerosol_particles + long_name = vertical integrated optical depth for various aerosol species + units = none + dimensions = (horizontal_loop_extent,number_of_species_for_aerosol_optical_depth) + type = real + kind = kind_phys + intent = inout + optional = F +[lw_optical_props_aerosol] + standard_name = longwave_optical_properties_for_aerosols + long_name = Fortran DDT containing RRTMGP optical properties + units = DDT + dimensions = () + type = ty_optical_props_1scl + intent = inout + optional = F +[errmsg] + standard_name = ccpp_error_message + long_name = error message for error handling in CCPP + units = none + dimensions = () + type = character + kind = len=* + intent = out + optional = F +[errflg] + standard_name = ccpp_error_flag + long_name = error flag for error handling in CCPP + units = flag + dimensions = () + type = integer + intent = out + optional = F diff --git a/physics/rrtmgp_sw_aerosol_optics.F90 b/physics/rrtmgp_sw_aerosol_optics.F90 new file mode 100644 index 000000000..3a74771b7 --- /dev/null +++ b/physics/rrtmgp_sw_aerosol_optics.F90 @@ -0,0 +1,119 @@ +module rrtmgp_sw_aerosol_optics + use machine, only: kind_phys + use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp + use mo_optical_props, only: ty_optical_props_2str + use radiation_tools, only: check_error_msg + use rrtmgp_sw_gas_optics, only: sw_gas_props + use rrtmgp_lw_gas_optics, only: lw_gas_props + use module_radiation_aerosols, only: & + NF_AESW, & ! Number of optical-fields in SW output (3=tau+g+omega) + NF_AELW, & ! Number of optical-fields in LW output (3=tau+g+omega) + setaer, & ! Routine to compute aerosol radiative properties (tau,g,omega) + NSPC1 ! Number of species for vertically integrated aerosol optical-depth + use netcdf + + implicit none + + public rrtmgp_sw_aerosol_optics_init, rrtmgp_sw_aerosol_optics_run, rrtmgp_sw_aerosol_optics_finalize + +contains + + ! ######################################################################################### + ! SUBROUTINE rrtmgp_sw_aerosol_optics_init() + ! ######################################################################################### + subroutine rrtmgp_sw_aerosol_optics_init() + end subroutine rrtmgp_sw_aerosol_optics_init + + ! ######################################################################################### + ! SUBROUTINE rrtmgp_sw_aerosol_optics_run() + ! ######################################################################################### +!! \section arg_table_rrtmgp_sw_aerosol_optics_run +!! \htmlinclude rrtmgp_sw_aerosol_optics_run.html +!! + subroutine rrtmgp_sw_aerosol_optics_run(doSWrad, nCol, nLev, nTracer, nTracerAer, nDay, & + idxday, p_lev, p_lay, p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, & + aerodp, sw_optical_props_aerosol, errmsg, errflg ) + + ! Inputs + logical, intent(in) :: & + doSWrad ! Logical flag for shortwave radiation call + integer, intent(in) :: & + nCol, & ! Number of horizontal grid points + nDay, & ! Number of daylit points + nLev, & ! Number of vertical layers + nTracer, & ! Number of tracers + nTracerAer ! Number of aerosol tracers + integer,intent(in),dimension(:) :: & + idxday ! Indices for daylit points. + real(kind_phys), dimension(:), intent(in) :: & + lon, & ! Longitude + lat, & ! Latitude + lsmask ! Land/sea/sea-ice mask + real(kind_phys), dimension(:,:),intent(in) :: & + p_lay, & ! Pressure @ layer-centers (Pa) + tv_lay, & ! Virtual-temperature @ layer-centers (K) + relhum, & ! Relative-humidity @ layer-centers + p_lk ! Exner function @ layer-centers (1) + real(kind_phys), dimension(:, :,:),intent(in) :: & + tracer ! trace gas concentrations + real(kind_phys), dimension(:, :,:),intent(in) :: & + aerfld ! aerosol input concentrations + real(kind_phys), dimension(:,:),intent(in) :: & + p_lev ! Pressure @ layer-interfaces (Pa) + + ! Outputs + real(kind_phys), dimension(:,:), intent(inout) :: & + aerodp ! Vertical integrated optical depth for various aerosol species + type(ty_optical_props_2str),intent(out) :: & + sw_optical_props_aerosol ! RRTMGP DDT: Longwave aerosol optical properties (tau) + integer, intent(out) :: & + errflg ! CCPP error flag + character(len=*), intent(out) :: & + errmsg ! CCPP error message + + ! Local variables + real(kind_phys), dimension(nCol, nLev, lw_gas_props%get_nband(), NF_AELW) :: & + aerosolslw ! + real(kind_phys), dimension(nCol, nLev, sw_gas_props%get_nband(), NF_AESW) :: & + aerosolssw, aerosolssw2 + + ! Initialize CCPP error handling variables + errmsg = '' + errflg = 0 + + if (.not. doSWrad) return + if (nDay .gt. 0) then + + ! Call module_radiation_aerosols::setaer(),to setup aerosols property profile + call setaer(p_lev/100., p_lay/100., p_lk, tv_lay, relhum, lsmask, tracer, aerfld, lon, lat, nCol, nLev, & + nLev+1, .true., .true., aerosolssw2, aerosolslw, aerodp) + + ! Store aerosol optical properties + ! SW. + ! For RRTMGP SW the bands are now ordered from [IR(band) -> nIR -> UV], in RRTMG the + ! band ordering was [nIR -> UV -> IR(band)] + aerosolssw(1:nCol,:,1,1) = aerosolssw2(1:nCol,:,sw_gas_props%get_nband(),1) + aerosolssw(1:nCol,:,1,2) = aerosolssw2(1:nCol,:,sw_gas_props%get_nband(),2) + aerosolssw(1:nCol,:,1,3) = aerosolssw2(1:nCol,:,sw_gas_props%get_nband(),3) + aerosolssw(1:nCol,:,2:sw_gas_props%get_nband(),1) = aerosolssw2(1:nCol,:,1:sw_gas_props%get_nband()-1,1) + aerosolssw(1:nCol,:,2:sw_gas_props%get_nband(),2) = aerosolssw2(1:nCol,:,1:sw_gas_props%get_nband()-1,2) + aerosolssw(1:nCol,:,2:sw_gas_props%get_nband(),3) = aerosolssw2(1:nCol,:,1:sw_gas_props%get_nband()-1,3) + + ! Allocate RRTMGP DDT: Aerosol optics [nCol,nlev,nBands] + call check_error_msg('rrtmgp_sw_aerosol_optics_run',sw_optical_props_aerosol%alloc_2str( & + nDay, nlev, sw_gas_props%get_band_lims_wavenumber())) + + ! Copy aerosol optical information to RRTMGP DDT + sw_optical_props_aerosol%tau = aerosolssw(idxday(1:nday),:,:,1) + sw_optical_props_aerosol%ssa = aerosolssw(idxday(1:nday),:,:,2) + sw_optical_props_aerosol%g = aerosolssw(idxday(1:nday),:,:,3) + endif + + end subroutine rrtmgp_sw_aerosol_optics_run + + ! ######################################################################################### + ! SUBROUTINE rrtmgp_sw_aerosol_optics_finalize() + ! ######################################################################################### + subroutine rrtmgp_sw_aerosol_optics_finalize() + end subroutine rrtmgp_sw_aerosol_optics_finalize +end module rrtmgp_sw_aerosol_optics diff --git a/physics/GFS_rrtmgp_aerosol_optics.meta b/physics/rrtmgp_sw_aerosol_optics.meta similarity index 91% rename from physics/GFS_rrtmgp_aerosol_optics.meta rename to physics/rrtmgp_sw_aerosol_optics.meta index aa7f6d4a5..f4909c794 100644 --- a/physics/GFS_rrtmgp_aerosol_optics.meta +++ b/physics/rrtmgp_sw_aerosol_optics.meta @@ -1,15 +1,15 @@ [ccpp-table-properties] - name = GFS_rrtmgp_aerosol_optics + name = rrtmgp_sw_aerosol_optics type = scheme dependencies = iounitdef.f,machine.F,radiation_aerosols.f,radiation_tools.F90 ######################################################################## [ccpp-arg-table] - name = GFS_rrtmgp_aerosol_optics_run + name = rrtmgp_sw_aerosol_optics_run type = scheme -[doLWrad] - standard_name = flag_to_calc_lw - long_name = logical flags for lw radiation calls +[doSWrad] + standard_name = flag_to_calc_sw + long_name = logical flags for sw radiation calls units = flag dimensions = () type = logical @@ -23,22 +23,6 @@ type = integer intent = in optional = F -[nday] - standard_name = daytime_points_dimension - long_name = daytime points dimension - units = count - dimensions = () - type = integer - intent = in - optional = F -[idxday] - standard_name = daytime_points - long_name = daytime points - units = index - dimensions = (horizontal_loop_extent) - type = integer - intent = in - optional = F [nLev] standard_name = vertical_dimension long_name = number of vertical levels @@ -63,6 +47,22 @@ type = integer intent = in optional = F +[nday] + standard_name = daytime_points_dimension + long_name = daytime points dimension + units = count + dimensions = () + type = integer + intent = in + optional = F +[idxday] + standard_name = daytime_points + long_name = daytime points + units = index + dimensions = (horizontal_loop_extent) + type = integer + intent = in + optional = F [p_lev] standard_name = air_pressure_at_interface_for_RRTMGP_in_hPa long_name = air pressure at vertical interface for radiation calculation @@ -162,14 +162,6 @@ kind = kind_phys intent = inout optional = F -[lw_optical_props_aerosol] - standard_name = longwave_optical_properties_for_aerosols - long_name = Fortran DDT containing RRTMGP optical properties - units = DDT - dimensions = () - type = ty_optical_props_1scl - intent = inout - optional = F [sw_optical_props_aerosol] standard_name = shortwave_optical_properties_for_aerosols long_name = Fortran DDT containing RRTMGP optical properties From eb7837d7f32e83fb6a151dbae50172bda5fd83db Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Tue, 25 May 2021 20:05:16 +0000 Subject: [PATCH 08/14] Bug fix. Add transition to HR adjustment. --- physics/dcyc2.f | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index 09c80a97e..8e2f86e5a 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -424,7 +424,8 @@ subroutine dcyc2t3_run & else lfnc = 1. endif - dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + htrlw(i,k)*lfnc + dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + & + & htrlw(i,k)*lfnc + (1.-lfnc)*hlw(i,k) enddo enddo else From f2d55708f4f8d1162284145707e0c53a2677e8f8 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Wed, 26 May 2021 10:33:24 +0000 Subject: [PATCH 09/14] Removed exp(1) from scaling --- physics/dcyc2.f | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index 8e2f86e5a..796c36f12 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -420,7 +420,7 @@ subroutine dcyc2t3_run & ! using a logistic function if (damp_LW_fluxadj) then dp = p_lev(i,k) - lfnc_p0 - lfnc = L / (1+exp(-lfnc_k_grad*exp(1.)*dp/lfnc_p0)) + lfnc = L / (1+exp(-lfnc_k_grad*dp/lfnc_p0)) else lfnc = 1. endif From ea0113929c9ec254e15275ab55f094360d0ffb55 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Wed, 26 May 2021 11:18:19 +0000 Subject: [PATCH 10/14] Housekeeping. Added comments. --- physics/dcyc2.f | 34 ++++++++++++++++++---------------- 1 file changed, 18 insertions(+), 16 deletions(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index 796c36f12..cfe7c75a8 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -260,14 +260,14 @@ subroutine dcyc2t3_run & & rstl, solang, dT real(kind=kind_phys), dimension(im,levs+1) :: flxlwup_adj, & & flxlwdn_adj, t_lev2 - real(kind=kind_phys) :: fluxlwnet_adj,fluxlwnet,c1,c2,c3,c4,c5, & - & dP,lfnc + real(kind=kind_phys) :: fluxlwnet_adj,fluxlwnet,dT_sfc, & + &fluxlwDOWN_jac,dP,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 :: & - & c0 = 0.2 + & gamma = 0.2 ! !===> ... begin here ! @@ -400,19 +400,21 @@ subroutine dcyc2t3_run & ! L = 1, fix scale between 0-1. ! k (steepness) and p0 (midpoint) are controlled via namelist do i = 1, im - c1 = fluxlwUP_jac(i,iSFC) - c2 = fluxlwUP_jac(i,iTOA) / c1 - c3 = t_lev2(i,iSFC) - t_lev(i,iSFC) + c1 = fluxlwUP_jac(i,iTOA) / fluxlwUP_jac(i,iSFC) + dT_sfc = t_lev2(i,iSFC) - t_lev(i,iSFC) do k = 1, levs - ! Only apply the Jacobian adjustment below plim_fluxAdj_upper - c4 = fluxlwUP_jac(i,k)/c1 - fluxlwnet = (flux2D_lwUP(i,k+1) - flux2D_lwUP(i,k) - & - & flux2D_lwDOWN(i,k+1) + flux2D_lwDOWN(i,k)) - ! (Eq. 9) - c5 = c0 * (c4 - c2) / (1 - c2) - ! (Eq. 10) - fluxlwnet_adj = fluxlwnet + c3*(c4-c5) - ! Compute adjusted heating rate + ! 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))) @@ -420,7 +422,7 @@ subroutine dcyc2t3_run & ! using a logistic function if (damp_LW_fluxadj) then dp = p_lev(i,k) - lfnc_p0 - lfnc = L / (1+exp(-lfnc_k_grad*dp/lfnc_p0)) + lfnc = L / (1+exp(-(lfnc_k_grad/lfnc_p0)*dp)) else lfnc = 1. endif From d93adbeb519a58e0aedbbe92a840186ee312cb54 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Wed, 26 May 2021 22:00:59 +0000 Subject: [PATCH 11/14] Further cleanup of dcyc2 --- physics/dcyc2.f | 14 +++++++------- physics/dcyc2.meta | 12 ++++++------ 2 files changed, 13 insertions(+), 13 deletions(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index cfe7c75a8..8700e2bfb 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -179,7 +179,7 @@ subroutine dcyc2t3_run & & sfcnirbmd,sfcnirdfd,sfcvisbmd,sfcvisdfd, & & im, levs, deltim, fhswr, & & dry, icy, wet, minGPtemp, maxGPtemp, damp_LW_fluxadj, & - & lfnc_k_grad, lfnc_p0, minGPpres, use_LW_jacobian, sfculw, & + & 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, & ! & dry, icy, wet, lprnt, ipr, & @@ -217,7 +217,7 @@ subroutine dcyc2t3_run & & pert_radtend logical, intent(in) :: do_sppt,ca_global real(kind=kind_phys), intent(in) :: solhr, slag, cdec, sdec, & - & deltim, fhswr, minGPpres, minGPtemp, maxGPtemp, lfnc_k_grad, & + & deltim, fhswr, minGPpres, minGPtemp, maxGPtemp, lfnc_k, & & lfnc_p0 real(kind=kind_phys), dimension(:), intent(in) :: & @@ -261,7 +261,7 @@ subroutine dcyc2t3_run & 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,dP,lfnc,c1 + &fluxlwDOWN_jac,lfnc,c1 ! Length scale for flux-adjustment scaling real(kind=kind_phys), parameter :: & & L = 1. @@ -397,8 +397,9 @@ subroutine dcyc2t3_run & ! ! 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. - ! k (steepness) and p0 (midpoint) are controlled via namelist + ! 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) @@ -421,8 +422,7 @@ subroutine dcyc2t3_run & ! Add radiative heating rates to physics heating rate. Optionally, scaled w/ height ! using a logistic function if (damp_LW_fluxadj) then - dp = p_lev(i,k) - lfnc_p0 - lfnc = L / (1+exp(-(lfnc_k_grad/lfnc_p0)*dp)) + lfnc = L / (1+exp(-lfnc_k*(p_lev(i,k) - lfnc_p0))) else lfnc = 1. endif diff --git a/physics/dcyc2.meta b/physics/dcyc2.meta index dceb9ce77..5ba718c2e 100644 --- a/physics/dcyc2.meta +++ b/physics/dcyc2.meta @@ -396,18 +396,18 @@ type = logical intent = in optional = F -[lfnc_k_grad] - standard_name = steepness_of_flux_damping - long_name = steepness of logistic function for damping the LW flux adjustment - units = none +[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 = midpoint_used_for_flux_damping - long_name = midpoint for damping the LW flux adjustment + standard_name = transition_pressure_for_flux_damping + long_name = transition pressure for LW flux adjustment damping units = Pa dimensions = () type = real From fe44d62031c77b4232f331ce27661aef7847fe96 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Wed, 26 May 2021 22:04:26 +0000 Subject: [PATCH 12/14] Revert "Added more safeguards against out-of-bounds temperature to GP inputs." This reverts commit 6961b10546035b55021132e62ce139333eed9cc0. --- physics/GFS_rrtmgp_pre.F90 | 12 ++++-------- physics/GFS_rrtmgp_pre.meta | 9 --------- physics/dcyc2.f | 15 +++++++-------- physics/dcyc2.meta | 18 ------------------ physics/radiation_tools.F90 | 20 ++++++-------------- physics/rrtmgp_lw_gas_optics.F90 | 4 +--- physics/rrtmgp_lw_gas_optics.meta | 9 --------- 7 files changed, 18 insertions(+), 69 deletions(-) diff --git a/physics/GFS_rrtmgp_pre.F90 b/physics/GFS_rrtmgp_pre.F90 index af7e5f1a0..88e534595 100644 --- a/physics/GFS_rrtmgp_pre.F90 +++ b/physics/GFS_rrtmgp_pre.F90 @@ -98,8 +98,8 @@ end subroutine GFS_rrtmgp_pre_init !! subroutine GFS_rrtmgp_pre_run(nCol, nLev, nTracers, i_o3, lsswr, lslwr, fhswr, fhlwr, & xlat, xlon, prsl, tgrs, prslk, prsi, qgrs, tsfc, con_eps, con_epsm1, con_fvirt, & - con_epsqs, minGPpres, minGPtemp, maxGPtemp, raddt, p_lay, t_lay, p_lev, t_lev, tsfg, & - tsfa, qs_lay, q_lay, tv_lay, relhum, tracer, gas_concentrations, errmsg, errflg) + con_epsqs, minGPpres, minGPtemp, raddt, p_lay, t_lay, p_lev, t_lev, tsfg, tsfa, & + qs_lay, q_lay, tv_lay, relhum, tracer, gas_concentrations, errmsg, errflg) ! Inputs integer, intent(in) :: & @@ -112,7 +112,6 @@ subroutine GFS_rrtmgp_pre_run(nCol, nLev, nTracers, i_o3, lsswr, lslwr, fhswr, f lslwr ! Call LW radiation real(kind_phys), intent(in) :: & minGPtemp, & ! Minimum temperature allowed in RRTMGP. - maxGPtemp, & ! Maximum temperature allowed in RRTMGP. minGPpres, & ! Minimum pressure allowed in RRTMGP. fhswr, & ! Frequency of SW radiation call. fhlwr ! Frequency of LW radiation call. @@ -209,14 +208,11 @@ subroutine GFS_rrtmgp_pre_run(nCol, nLev, nTracers, i_o3, lsswr, lslwr, fhswr, f if (t_lay(iCol,iLay) .le. minGPtemp) then t_lay(iCol,iLay) = minGPtemp + epsilon(minGPtemp) endif - if (t_lay(iCol,iLay) .ge. maxGPtemp) then - t_lay(iCol,iLay) = maxGPtemp - epsilon(maxGPtemp) - endif enddo enddo ! Temperature at layer-interfaces - call cmp_tlev(nCol,nLev,minGPpres,minGPtemp,maxGPtemp,p_lay,t_lay,p_lev,tsfc,t_lev) + call cmp_tlev(nCol,nLev,minGPpres,p_lay,t_lay,p_lev,tsfc,t_lev) ! Compute a bunch of thermodynamic fields needed by the cloud microphysics schemes. ! Relative humidity, saturation mixing-ratio, vapor mixing-ratio, virtual temperature, @@ -277,7 +273,7 @@ subroutine GFS_rrtmgp_pre_run(nCol, nLev, nTracers, i_o3, lsswr, lslwr, fhswr, f ! ####################################################################################### ! Setup surface ground temperature and ground/air skin temperature if required. ! ####################################################################################### - tsfg(1:NCOL) = t_lev(1:NCOL,iSFC) + tsfg(1:NCOL) = tsfc(1:NCOL) tsfa(1:NCOL) = t_lay(1:NCOL,iSFC) end subroutine GFS_rrtmgp_pre_run diff --git a/physics/GFS_rrtmgp_pre.meta b/physics/GFS_rrtmgp_pre.meta index 895bbc630..8096aef2a 100644 --- a/physics/GFS_rrtmgp_pre.meta +++ b/physics/GFS_rrtmgp_pre.meta @@ -239,15 +239,6 @@ kind = kind_phys intent = in optional = F -[maxGPtemp] - standard_name = maximum_temperature_in_RRTMGP - long_name = maximum temperature allowed in RRTMGP - units = K - dimensions = () - type = real - kind = kind_phys - intent = in - optional = F [raddt] standard_name = time_step_for_radiation long_name = radiation time step diff --git a/physics/dcyc2.f b/physics/dcyc2.f index 8700e2bfb..6247f360f 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -178,10 +178,10 @@ subroutine dcyc2t3_run & & sfcnirbmu,sfcnirdfu,sfcvisbmu,sfcvisdfu, & & sfcnirbmd,sfcnirdfd,sfcvisbmd,sfcvisdfd, & & im, levs, deltim, fhswr, & - & dry, icy, wet, minGPtemp, maxGPtemp, 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, & + & 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, & ! & dry, icy, wet, lprnt, ipr, & ! --- input/output: & dtdt,dtdtnp,htrlw, & @@ -217,8 +217,7 @@ subroutine dcyc2t3_run & & pert_radtend logical, intent(in) :: do_sppt,ca_global real(kind=kind_phys), intent(in) :: solhr, slag, cdec, sdec, & - & deltim, fhswr, minGPpres, minGPtemp, maxGPtemp, lfnc_k, & - & lfnc_p0 + & deltim, fhswr, minGPpres, minGPtemp, lfnc_k, lfnc_p0 real(kind=kind_phys), dimension(:), intent(in) :: & & sinlat, coslat, xlon, coszen, tf, tsflw, sfcdlw, & @@ -384,8 +383,8 @@ subroutine dcyc2t3_run & ! ! Compute temperatute at level interfaces. ! - call cmp_tlev(im, levs, minGPpres, minGPtemp, maxGPtemp, p_lay,& - & t_lay, p_lev, tsfc, t_lev2) + 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 diff --git a/physics/dcyc2.meta b/physics/dcyc2.meta index 5ba718c2e..91e01a2d2 100644 --- a/physics/dcyc2.meta +++ b/physics/dcyc2.meta @@ -353,24 +353,6 @@ type = logical intent = in optional = F -[minGPtemp] - standard_name = minimum_temperature_in_RRTMGP - long_name = minimum temperature allowed in RRTMGP - units = K - dimensions = () - type = real - kind = kind_phys - intent = in - optional = F -[maxGPtemp] - standard_name = maximum_temperature_in_RRTMGP - long_name = maximum temperature allowed in RRTMGP - units = K - dimensions = () - type = real - kind = kind_phys - intent = in - optional = F [minGPpres] standard_name = minimum_pressure_in_RRTMGP long_name = minimum pressure allowed in RRTMGP diff --git a/physics/radiation_tools.F90 b/physics/radiation_tools.F90 index a8d3f5457..c6524aab6 100644 --- a/physics/radiation_tools.F90 +++ b/physics/radiation_tools.F90 @@ -2,16 +2,20 @@ module radiation_tools use machine, only: & kind_phys ! Working type implicit none + + real(kind_phys) :: & + rrtmgp_minP, & ! Minimum pressure allowed in RRTMGP + rrtmgp_minT ! Minimum temperature allowed in RRTMGP contains ! ######################################################################################### ! ######################################################################################### - subroutine cmp_tlev(nCol,nLev,minP,minT,maxT,p_lay,t_lay,p_lev,tsfc,t_lev) + subroutine cmp_tlev(nCol,nLev,minP,p_lay,t_lay,p_lev,tsfc,t_lev) ! Inputs integer, intent(in) :: & nCol,nLev real(kind_phys),intent(in) :: & - minP,minT,maxT + minP real(kind_phys),dimension(nCol),intent(in) :: & tsfc real(kind_phys),dimension(nCol,nLev),intent(in) :: & @@ -74,18 +78,6 @@ subroutine cmp_tlev(nCol,nLev,minP,minT,maxT,p_lay,t_lay,p_lev,tsfc,t_lev) t_lev(1:NCOL,iTOA+1) = t_lay(1:NCOL,iTOA) endif - ! Bound temperature at layer interfaces - do iCol=1,NCOL - do iLay=1,nLev+1 - if (t_lev(iCol,iLay) .le. minT) then - t_lev(iCol,iLay) = minT + epsilon(minT) - endif - if (t_lev(iCol,iLay) .ge. maxT) then - t_lev(iCol,iLay) = maxT - epsilon(maxT) - endif - enddo - enddo - end subroutine cmp_tlev ! ######################################################################################### diff --git a/physics/rrtmgp_lw_gas_optics.F90 b/physics/rrtmgp_lw_gas_optics.F90 index d7201e026..a116ad772 100644 --- a/physics/rrtmgp_lw_gas_optics.F90 +++ b/physics/rrtmgp_lw_gas_optics.F90 @@ -76,7 +76,7 @@ module rrtmgp_lw_gas_optics !! \htmlinclude rrtmgp_lw_gas_optics_init.html !! subroutine rrtmgp_lw_gas_optics_init(rrtmgp_root_dir, rrtmgp_lw_file_gas, mpicomm, & - mpirank, mpiroot, gas_concentrations, minGPpres, minGPtemp, maxGPtemp, errmsg, errflg) + mpirank, mpiroot, gas_concentrations, minGPpres, minGPtemp, errmsg, errflg) ! Inputs type(ty_gas_concs), intent(inout) :: & @@ -96,7 +96,6 @@ subroutine rrtmgp_lw_gas_optics_init(rrtmgp_root_dir, rrtmgp_lw_file_gas, mpicom errflg ! CCPP error code real(kind_phys), intent(out) :: & minGPtemp, & ! Minimum temperature allowed by RRTMGP. - maxGPtemp, & ! Maximum temperature allowed by RRTMG. minGPpres ! Minimum pressure allowed by RRTMGP. ! Local variables @@ -451,7 +450,6 @@ subroutine rrtmgp_lw_gas_optics_init(rrtmgp_root_dir, rrtmgp_lw_file_gas, mpicom ! temperature (GFS_rrtmgp_pre.F90) minGPpres = lw_gas_props%get_press_min() minGPtemp = lw_gas_props%get_temp_min() - maxGPtemp = lw_gas_props%get_temp_max() end subroutine rrtmgp_lw_gas_optics_init diff --git a/physics/rrtmgp_lw_gas_optics.meta b/physics/rrtmgp_lw_gas_optics.meta index 823501cfa..c92567e14 100644 --- a/physics/rrtmgp_lw_gas_optics.meta +++ b/physics/rrtmgp_lw_gas_optics.meta @@ -92,15 +92,6 @@ kind = kind_phys intent = out optional = F -[maxGPtemp] - standard_name = maximum_temperature_in_RRTMGP - long_name = maximum temperature allowed in RRTMGP - units = K - dimensions = () - type = real - kind = kind_phys - intent = out - optional = F ######################################################################## [ccpp-arg-table] From 6a0e904eadd35bd201e9848f68d025a3b8b7db51 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Thu, 27 May 2021 15:13:38 +0000 Subject: [PATCH 13/14] Omission from previous revert. --- physics/dcyc2.f | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index 6247f360f..0e3a4db42 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -217,7 +217,7 @@ subroutine dcyc2t3_run & & pert_radtend logical, intent(in) :: do_sppt,ca_global real(kind=kind_phys), intent(in) :: solhr, slag, cdec, sdec, & - & deltim, fhswr, minGPpres, minGPtemp, lfnc_k, lfnc_p0 + & deltim, fhswr, minGPpres, lfnc_k, lfnc_p0 real(kind=kind_phys), dimension(:), intent(in) :: & & sinlat, coslat, xlon, coszen, tf, tsflw, sfcdlw, & From 41782f1170b332a8c5e0c9324c83df48fb6df1b2 Mon Sep 17 00:00:00 2001 From: Dustin Swales Date: Thu, 27 May 2021 17:03:13 +0000 Subject: [PATCH 14/14] Change from PR review. --- physics/dcyc2.f | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/physics/dcyc2.f b/physics/dcyc2.f index 0e3a4db42..dfa9f02ed 100644 --- a/physics/dcyc2.f +++ b/physics/dcyc2.f @@ -421,7 +421,7 @@ subroutine dcyc2t3_run & ! 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(-lfnc_k*(p_lev(i,k) - lfnc_p0))) + lfnc = L / (1+exp(-(p_lev(i,k) - lfnc_p0)/lfnc_k)) else lfnc = 1. endif