Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Revisions to repair iovr=5 cloud overlap option #829

Closed
wants to merge 3 commits into from
Closed

Revisions to repair iovr=5 cloud overlap option #829

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
57d7533
Revisions to repair iovr=5 cloud overlap option
mjiacono Jan 18, 2022
f50f1f8
add fix for unbalanced parenthses
mjiacono Jan 19, 2022
5ec888f
Merge branch 'NCAR:main' into master
mjiacono Jan 20, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
185 changes: 154 additions & 31 deletions physics/radiation_clouds.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -876,6 +876,19 @@ subroutine progcld1 &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

!> - Call gethml() to compute low,mid,high,total, and boundary layer
!! cloud fractions and clouds top/bottom layer indices for low, mid,
!! and high clouds. The three cloud domain boundaries are defined by
Expand Down Expand Up @@ -1272,6 +1285,19 @@ subroutine progcld2 &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

!> - Call gethml() to compute low,mid,high,total, and boundary layer
!! cloud fractions and clouds top/bottom layer indices for low, mid,
!! and high clouds.
Expand Down Expand Up @@ -1699,6 +1725,19 @@ subroutine progcld3 &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

!> -# Call gethml() to compute low,mid,high,total, and boundary layer
!! cloud fractions and clouds top/bottom layer indices for low, mid,
!! and high clouds.
Expand Down Expand Up @@ -2062,6 +2101,19 @@ subroutine progcld4 &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

! --- compute low, mid, high, total, and boundary layer cloud fractions
! and clouds top/bottom layer indices for low, mid, and high clouds.
! The three cloud domain boundaries are defined by ptopc. The cloud
Expand Down Expand Up @@ -2416,6 +2468,19 @@ subroutine progcld4o &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

!> - Call gethml() to compute low, mid, high, total, and boundary layer cloud fractions
!! and clouds top/bottom layer indices for low, mid, and high clouds.
!! The three cloud domain boundaries are defined by ptopc. The cloud
Expand Down Expand Up @@ -2792,6 +2857,19 @@ subroutine progcld5 &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

!> - Call gethml() to compute low,mid,high,total, and boundary layer
!! cloud fractions and clouds top/bottom layer indices for low, mid,
!! and high clouds.
Expand Down Expand Up @@ -3128,8 +3206,7 @@ subroutine progcld6 &
enddo
enddo

! What portion of water and ice contents is associated with the
! partly cloudy boxes
! What portion of water and ice contents is associated with the partly cloudy boxes
do i = 1, IX
do k = 1, NLAY-1
if (cldtot(i,k).ge.climit .and. cldtot(i,k).lt.ovcst) then
Expand Down Expand Up @@ -3188,6 +3265,19 @@ subroutine progcld6 &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

!> - Call gethml() to compute low,mid,high,total, and boundary layer
!! cloud fractions and clouds top/bottom layer indices for low, mid,
!! and high clouds.
Expand Down Expand Up @@ -3274,7 +3364,7 @@ subroutine progcld_thompson &
! slmsk (IX) : sea/land mask array (sea:0,land:1,sea-ice:2) !
! dz (ix,nlay) : layer thickness (km) !
! delp (ix,nlay) : model layer pressure thickness in mb (100Pa) !
! gridkm (IX) : grid length in km !
! gridkm : grid length in km !
! IX : horizontal dimention !
! NLAY,NLP1 : vertical layer/level dimensions !
! uni_cld : logical - true for cloud fraction from shoc !
Expand Down Expand Up @@ -3333,8 +3423,8 @@ subroutine progcld_thompson &
real (kind=kind_phys), dimension(:), intent(in) :: xlat, xlon, &
& slmsk

real(kind=kind_phys), dimension(:), intent(in) :: latdeg, gridkm
real(kind=kind_phys), intent(in) :: julian
real(kind=kind_phys), dimension(:), intent(in) :: latdeg
real(kind=kind_phys), intent(in) :: julian, gridkm
integer, intent(in) :: yearlen

! --- outputs
Expand Down Expand Up @@ -3408,14 +3498,14 @@ subroutine progcld_thompson &
enddo

!> - Compute cloud liquid/ice condensate path in \f$ g/m^2 \f$ .
!> - Since using Thompson MP, assume 1 percent of snow is actually in
!> - Since using Thompson MP, assume 20 percent of snow is actually in
!! ice sizes.

do k = 1, NLAY-1
do i = 1, IX
cwp(i,k) = max(0.0, clw(i,k,ntcw) * dz(i,k)*1.E6)
crp(i,k) = 0.0
snow_mass_factor = 0.99
snow_mass_factor = 0.85
cip(i,k) = max(0.0, (clw(i,k,ntiw) &
& + (1.0-snow_mass_factor)*clw(i,k,ntsw))*dz(i,k)*1.E6)
if (re_snow(i,k) .gt. snow_max_radius)then
Expand Down Expand Up @@ -3481,7 +3571,7 @@ subroutine progcld_thompson &
endif

call cal_cldfra3(cldfra1d, qv1d, qc1d, qi1d, qs1d, dz1d, &
& p1d, t1d, xland, gridkm(i), &
& p1d, t1d, xland, gridkm, &
& .false., max_relh, 1, nlay, .false.)

do k = 1, NLAY
Expand Down Expand Up @@ -3555,6 +3645,19 @@ subroutine progcld_thompson &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

!> - Call gethml() to compute low,mid,high,total, and boundary layer
!! cloud fractions and clouds top/bottom layer indices for low, mid,
!! and high clouds.
Expand Down Expand Up @@ -3952,6 +4055,19 @@ subroutine progclduni &
alpha(:,:) = 0.
endif

! Revise alpha for exponential-random cloud overlap
! Decorrelate layers when a clear layer follows a cloudy layer to enforce
! random correlation between non-adjacent blocks of cloudy layers
if (iovr == 5) then
do k = 2, nLay
do i = 1, ix
if (clouds(i,k,1) == 0.0 .and. clouds(i,k-1,1) > 0.0) then
alpha(i,k) = 0.0
endif
enddo
enddo
endif

!> - Call gethml() to compute low,mid,high,total, and boundary layer
!! cloud fractions and clouds top/bottom layer indices for low, mid,
!! and high clouds.
Expand Down Expand Up @@ -4495,16 +4611,16 @@ SUBROUTINE cal_cldfra3(CLDFRA, qv, qc, qi, qs, dz, &
DO k = kts,kte

delz = MAX(100., dz(k))
RH_00L = 0.77+MIN(0.22,SQRT(1./(50.0+gridkm*gridkm*delz*0.01)))
RH_00O = 0.85+MIN(0.14,SQRT(1./(50.0+gridkm*gridkm*delz*0.01)))
RH_00L = 0.74+MIN(0.25,SQRT(1./(50.0+gridkm*gridkm*delz*0.01)))
RH_00O = 0.82+MIN(0.17,SQRT(1./(50.0+gridkm*gridkm*delz*0.01)))
RHUM = rh(k)

if (qc(k).ge.1.E-6 .or. qi(k).ge.1.E-7 &
& .or. (qs(k).gt.1.E-6 .and. t(k).lt.273.)) then
if (qc(k).ge.1.E-5 .or. qi(k).ge.1.E-5 &
& .or. (qs(k).gt.1.E-5 .and. t(k).lt.273.)) then
CLDFRA(K) = 1.0
elseif (((qc(k)+qi(k)).gt.1.E-10) .and. &
& ((qc(k)+qi(k)).lt.1.E-6)) then
CLDFRA(K) = MIN(0.99, 0.1*(11.0 + log10(qc(k)+qi(k))))
& ((qc(k)+qi(k)).lt.1.E-5)) then
CLDFRA(K) = MIN(0.99, 0.20*(10.0 + log10(qc(k)+qi(k))))
else

IF ((XLAND-1.5).GT.0.) THEN !--- Ocean
Expand All @@ -4513,7 +4629,7 @@ SUBROUTINE cal_cldfra3(CLDFRA, qv, qc, qi, qs, dz, &
RH_00 = RH_00L
ENDIF

tc = MAX(-80.0, t(k) - 273.15)
tc = t(k) - 273.15
if (tc .lt. -12.0) RH_00 = RH_00L

if (tc .gt. 20.0) then
Expand All @@ -4525,12 +4641,12 @@ SUBROUTINE cal_cldfra3(CLDFRA, qv, qc, qi, qs, dz, &
if (max_relh.gt.1.12 .or. (.NOT.(modify_qvapor)) ) then
!..For HRRR model, the following look OK.
RHUM = MIN(rh(k), 1.45)
RH_00 = RH_00 + (1.45-RH_00)*(-12.0-tc)/(-12.0+85.)
RH_00 = RH_00 + (1.45-RH_00)*(-12.0-tc)/(-12.0+112.)
CLDFRA(K) = MAX(0.,1.0-SQRT((1.46-RHUM)/(1.46-RH_00)))
else
!..but for the GFS model, RH is way lower.
RHUM = MIN(rh(k), 1.05)
RH_00 = RH_00 + (1.05-RH_00)*(-12.0-tc)/(-12.0+85.)
RH_00 = RH_00 + (1.05-RH_00)*(-12.0-tc)/(-12.0+112.)
CLDFRA(K) = MAX(0.,1.0-SQRT((1.06-RHUM)/(1.06-RH_00)))
endif
endif
Expand All @@ -4548,6 +4664,15 @@ SUBROUTINE cal_cldfra3(CLDFRA, qv, qc, qi, qs, dz, &

call adjust_cloudFinal(cldfra, qc, qi, rhoa, dz, kts,kte)

if (debug_flag .and. ndebug.lt.25) then
do k = kts,kte
write(6,'(a,i3,f9.2,f7.1,f7.2,f6.1,f6.3,f12.7,f12.7,f12.7)') &
& ' DEBUG-GT: ', k, p(k)*0.01, dz(k), t(k)-273.15, &
& rh(k)*100., cldfra(k), qc(k)*1.E3, qi(k)*1.E3, qs(k)*1.E3
enddo
ndebug = ndebug + 1
endif

!..Intended for cold start model runs, we use modify_qvapor to ensure that cloudy
!.. areas are actually saturated such that the inserted clouds do not evaporate a
!.. timestep later.
Expand Down Expand Up @@ -4689,9 +4814,9 @@ SUBROUTINE find_cloudLayers(qvs1d, cfr1d, T1d, P1d, Dz1d, entrmnt,&
k = k - 1
ENDDO

k_cldb = k_m12C + 3
k_cldb = k_m12C + 5
in_cloud = .false.
k = k_m12C + 2
k = k_m12C + 4
DO WHILE (.not. in_cloud .AND. k.gt.kbot)
k_cldt = 0
if (cfr1d(k).ge.0.01) then
Expand Down Expand Up @@ -4740,13 +4865,12 @@ SUBROUTINE adjust_cloudIce(cfr,qi,qs,qvs,T,dz,entr, k1,k2,kts,kte)
do k = k1, k2
tdz = tdz + dz(k)
enddo
! max_iwc = ABS(qvs(k2)-qvs(k1))
max_iwc = MAX(0.0, qvs(k1)-qvs(k2))
max_iwc = ABS(qvs(k2)-qvs(k1))

do k = k1, k2
max_iwc = MAX(1.E-6, max_iwc - (qi(k)+qs(k)))
max_iwc = MAX(1.E-5, max_iwc - (qi(k)+qs(k)))
enddo
max_iwc = MIN(1.E-4, max_iwc)
max_iwc = MIN(2.E-3, max_iwc)

this_dz = 0.0
do k = k1, k2
Expand All @@ -4756,7 +4880,7 @@ SUBROUTINE adjust_cloudIce(cfr,qi,qs,qvs,T,dz,entr, k1,k2,kts,kte)
this_dz = this_dz + dz(k)
endif
this_iwc = max_iwc*this_dz/tdz
iwc = MAX(1.E-6, this_iwc*(1.-entr))
iwc = MAX(5.E-6, this_iwc*(1.-entr))
if (cfr(k).gt.0.0.and.cfr(k).lt.1.0.and.T(k).ge.203.16) then
qi(k) = qi(k) + cfr(k)*cfr(k)*iwc
endif
Expand All @@ -4781,14 +4905,13 @@ SUBROUTINE adjust_cloudH2O(cfr, qc, qvs,T,dz,entr, k1,k2,kts,kte)
do k = k1, k2
tdz = tdz + dz(k)
enddo
! max_lwc = ABS(qvs(k2)-qvs(k1))
max_lwc = MAX(0.0, qvs(k1)-qvs(k2))
max_lwc = ABS(qvs(k2)-qvs(k1))
! print*, ' max_lwc = ', max_lwc, ' over DZ=',tdz

do k = k1, k2
max_lwc = MAX(1.E-6, max_lwc - qc(k))
max_lwc = MAX(1.E-5, max_lwc - qc(k))
enddo
max_lwc = MIN(1.E-4, max_lwc)
max_lwc = MIN(2.E-3, max_lwc)

this_dz = 0.0
do k = k1, k2
Expand All @@ -4798,8 +4921,8 @@ SUBROUTINE adjust_cloudH2O(cfr, qc, qvs,T,dz,entr, k1,k2,kts,kte)
this_dz = this_dz + dz(k)
endif
this_lwc = max_lwc*this_dz/tdz
lwc = MAX(1.E-6, this_lwc*(1.-entr))
if (cfr(k).gt.0.0.and.cfr(k).lt.1.0.and.T(k).ge.258.16) then
lwc = MAX(5.E-6, this_lwc*(1.-entr))
if (cfr(k).gt.0.0.and.cfr(k).lt.1.0.and.T(k).ge.253.16) then
qc(k) = qc(k) + cfr(k)*cfr(k)*lwc
endif
enddo
Expand Down Expand Up @@ -4851,6 +4974,6 @@ SUBROUTINE adjust_cloudFinal(cfr, qc, qi, Rho,dz, kts,kte)
END SUBROUTINE adjust_cloudFinal

!........................................!
end module module_radiation_clouds
end module module_radiation_clouds !
!! @}
!========================================!
Loading