@@ -8,7 +8,8 @@ module samfdeepcnv
88
99 use samfcnv_aerosols, only : samfdeepcnv_aerosols
1010 use progsigma, only : progsigma_calc
11-
11+ use progomega, only : progomega_calc
12+
1213 contains
1314
1415 subroutine samfdeepcnv_init (imfdeepcnv ,imfdeepcnv_samf , &
@@ -77,14 +78,15 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
7778 & eps ,epsm1 ,fv ,grav ,hvap ,rd ,rv , &
7879 & t0c ,delt ,ntk ,ntr ,delp , &
7980 & prslp ,psp ,phil ,tkeh ,qtr ,prevsq ,q ,q1 ,t1 ,u1 ,v1 ,fscav , &
80- & hwrf_samfdeep ,progsigma ,cldwrk ,rn ,kbot ,ktop ,kcnv , &
81+ & hwrf_samfdeep ,progsigma ,progomega , cldwrk ,rn ,kbot ,ktop ,kcnv , &
8182 & islimsk ,garea ,dot ,ncloud ,hpbl ,ud_mf ,dd_mf ,dt_mf ,cnvw ,cnvc , &
8283 & QLCN , QICN , w_upi , cf_upi , CNV_MFD , &
8384 & CNV_DQLDT ,CLCN ,CNV_FICE ,CNV_NDROP ,CNV_NICE ,mp_phys ,mp_phys_mg ,&
8485 & clam ,c0s ,c1 ,betal ,betas ,evef ,pgcon ,asolfac , &
8586 & do_ca , ca_closure , ca_entr , ca_trigger , nthresh ,ca_deep , &
86- & rainevap ,sigmain ,sigmaout ,betadcu ,betamcu ,betascu , &
87- & maxMF , do_mynnedmf ,sigmab_coldstart ,errmsg ,errflg )
87+ & rainevap ,sigmain ,sigmaout ,omegain ,omegaout ,betadcu ,betamcu , &
88+ & betascu ,maxMF ,do_mynnedmf ,sigmab_coldstart ,errmsg ,errflg )
89+
8890!
8991 use machine , only : kind_phys
9092 use funcphys , only : fpvs
@@ -100,16 +102,17 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
100102 & prslp(:,:), garea(:), hpbl(:), dot(:,:), phil(:,:)
101103 real (kind= kind_phys), dimension (:), intent (in ) :: fscav
102104 logical , intent (in ) :: first_time_step,restart,hwrf_samfdeep, &
103- & progsigma,do_mynnedmf,sigmab_coldstart
105+ & progsigma,progomega, do_mynnedmf,sigmab_coldstart
104106 real (kind= kind_phys), intent (in ) :: nthresh,betadcu,betamcu, &
105107 & betascu
106108 real (kind= kind_phys), intent (in ), optional :: ca_deep(:)
107109 real (kind= kind_phys), intent (in ), optional :: sigmain(:,:), &
108- & qmicro(:,:), prevsq(:,:)
110+ & qmicro(:,:), prevsq(:,:), omegain(:,:)
109111 real (kind= kind_phys), intent (in ) :: tmf(:,:,:),q(:,:)
110112 real (kind= kind_phys), dimension (:), intent (in ), optional :: maxMF
111113 real (kind= kind_phys), intent (out ) :: rainevap(:)
112- real (kind= kind_phys), intent (out ), optional :: sigmaout(:,:)
114+ real (kind= kind_phys), intent (inout ), optional :: sigmaout(:,:), &
115+ & omegaout(:,:)
113116 logical , intent (in ) :: do_ca,ca_closure,ca_entr,ca_trigger
114117 integer , intent (inout ) :: kcnv(:)
115118 ! DH* TODO - check dimensions of qtr, ntr+2 correct? * DH
@@ -216,7 +219,7 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
216219! parameters for prognostic sigma closure
217220 real (kind= kind_phys) omega_u(im,km),zdqca(im,km),tmfq(im,km),
218221 & omegac(im),zeta(im,km),dbyo1(im,km),sigmab(im),qadv(im,km),
219- & sigmaoutx(im)
222+ & sigmaoutx(im),tentr(im,km)
220223 real (kind= kind_phys) gravinv,invdelt,sigmind,sigminm,sigmins
221224 parameter (sigmind= 0.01 ,sigmins= 0.03 ,sigminm= 0.01 )
222225 logical flag_shallow, flag_mid
@@ -333,6 +336,7 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
333336c -----------------------------------------------------------------------
334337!> ## Compute preliminary quantities needed for static, dynamic, and feedback control portions of the algorithm.
335338!> - Convert input pressure terms to centibar units.
339+
336340!************************************************************************
337341! convert input Pa terms to Cb terms -- Moorthi
338342 ps = psp * 0.001
@@ -1130,7 +1134,8 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
11301134 do k = 2 , km1
11311135 do i= 1 ,im
11321136 if (cnvflg(i) .and.
1133- & (k > kbcon(i) .and. k < kmax(i))) then
1137+ & (k > kbcon(i) .and. k < kmax(i))) then
1138+ tentr(i,k)= xlamue(i,k)* fent1(i,k)
11341139 tem = cxlamet(i) * frh(i,k) * fent2(i,k)
11351140 xlamue(i,k) = xlamue(i,k)* fent1(i,k) + tem
11361141 tem1 = cxlamdt(i) * frh(i,k)
@@ -1760,43 +1765,64 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
17601765 enddo
17611766!
17621767! compute updraft velocity square(wu2)
1763- !> - Calculate updraft velocity square(wu2) according to Han et al.' s (2017 ) \cite han_et_al_2017 equation 7 .
1764- !
1768+ !> - Calculate diagnostic updraft velocity square(wu2) according to Han et al.' s (2017 ) \cite han_et_al_2017 equation 7 .
1769+ !> - if progomega = true, calculate prognostic updraft velocity (Pa/ s) according to progomega routine.
1770+
17651771 if (hwrf_samfdeep) then
1766- do i = 1 , im
1767- if (cnvflg(i)) then
1768- k = kbcon1(i)
1769- tem = po(i,k) / (rd * to (i,k))
1770- wucb = - 0.01 * dot(i,k) / (tem * grav)
1771- if (wucb.gt. 0 .) then
1772- wu2(i,k) = wucb * wucb
1773- else
1774- wu2(i,k) = 0 .
1775- endif
1776- endif
1777- enddo
1778- endif
1779- !
1780- do k = 2 , km1
1781- do i = 1 , im
1782- if (cnvflg(i)) then
1783- if (k > kbcon1(i) .and. k < ktcon(i)) then
1784- dz = zi(i,k) - zi(i,k-1 )
1785- tem = 0.25 * bb1 * (drag(i,k-1 )+ drag(i,k)) * dz
1786- tem1 = 0.5 * bb2 * (buo(i,k-1 )+ buo(i,k))
1787- tem2 = wush(i,k) * sqrt (wu2(i,k-1 ))
1788- tem2 = (tem1 - tem2) * dz
1789- ptem = (1 . - tem) * wu2(i,k-1 )
1790- ptem1 = 1 . + tem
1791- wu2(i,k) = (ptem + tem2) / ptem1
1792- wu2(i,k) = max (wu2(i,k), 0 .)
1772+ do i = 1 , im
1773+ if (cnvflg(i)) then
1774+ k = kbcon1(i)
1775+ tem = po(i,k) / (rd * to (i,k))
1776+ wucb = - 0.01 * dot(i,k) / (tem * grav)
1777+ if (wucb.gt. 0 .) then
1778+ wu2(i,k) = wucb * wucb
1779+ else
1780+ wu2(i,k) = 0 .
1781+ endif
17931782 endif
1794- endif
1795- enddo
1796- enddo
1797-
1798- if (progsigma)then
1799- do k = 2 , km1
1783+ enddo
1784+ endif
1785+ !
1786+ if (progomega) then
1787+ call progomega_calc(first_time_step,restart,im,km,
1788+ & kbcon1,ktcon,omegain,delt,del,zi,cnvflg,omegaout,
1789+ & grav,buo,drag,wush,tentr,bb1,bb2)
1790+ do k = 1 , km
1791+ do i = 1 , im
1792+ if (cnvflg(i)) then
1793+ if (k > kbcon1(i) .and. k < ktcon(i)) then
1794+ omega_u(i,k)= omegaout(i,k)
1795+ omega_u(i,k)= MAX (omega_u(i,k),- 80 .)
1796+ ! Convert to m/ s for use in convective time- scale:
1797+ rho = po(i,k)* 100 . / (rd * to (i,k))
1798+ tem = (- omega_u(i,k)) / ((rho * grav))
1799+ wu2(i,k) = tem** 2
1800+ wu2(i,k) = max (wu2(i,k), 0 .)
1801+ endif
1802+ endif
1803+ enddo
1804+ enddo
1805+ else
1806+ ! diagnostic method:
1807+ do k = 2 , km1
1808+ do i = 1 , im
1809+ if (cnvflg(i)) then
1810+ if (k > kbcon1(i) .and. k < ktcon(i)) then
1811+ dz = zi(i,k) - zi(i,k-1 )
1812+ tem = 0.25 * bb1 * (drag(i,k-1 )+ drag(i,k)) * dz
1813+ tem1 = 0.5 * bb2 * (buo(i,k-1 )+ buo(i,k))
1814+ tem2 = wush(i,k) * sqrt (wu2(i,k-1 ))
1815+ tem2 = (tem1 - tem2) * dz
1816+ ptem = (1 . - tem) * wu2(i,k-1 )
1817+ ptem1 = 1 . + tem
1818+ wu2(i,k) = (ptem + tem2) / ptem1
1819+ wu2(i,k) = max (wu2(i,k), 0 .)
1820+ endif
1821+ endif
1822+ enddo
1823+ enddo
1824+ ! convert to Pa/ s for use in closure
1825+ do k = 1 , km
18001826 do i = 1 , im
18011827 if (cnvflg(i)) then
18021828 if (k > kbcon1(i) .and. k < ktcon(i)) then
@@ -1807,10 +1833,11 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
18071833 endif
18081834 enddo
18091835 enddo
1810- endif
1836+
1837+ endif !progomega
1838+
18111839!
18121840! compute updraft velocity average over the whole cumulus
1813- !
18141841!> - Calculate the mean updraft velocity within the cloud (wc).
18151842 do i = 1 , im
18161843 wc(i) = 0 .
@@ -1838,11 +1865,10 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
18381865 val = 1.e-4
18391866 if (wc(i) < val) cnvflg(i)= .false.
18401867 endif
1841- enddo
1868+ enddo
18421869c
1843-
18441870!> - For progsigma = T, calculate the mean updraft velocity within the cloud (omegac),cast in pressure coordinates.
1845- if (progsigma)then
1871+ if (progsigma)then
18461872 do i = 1 , im
18471873 omegac(i) = 0 .
18481874 sumx(i) = 0 .
@@ -2929,10 +2955,9 @@ subroutine samfdeepcnv_run (im,km,first_time_step,restart, &
29292955 advfac(i) = min (advfac(i), 1 .)
29302956 endif
29312957 enddo
2932-
2958+
29332959!> - From Bengtsson et al. (2022 ) \cite Bengtsson_2022 prognostic closure scheme, equation 8 , call progsigma_calc() to compute updraft area fraction based on a moisture budget
29342960 if (progsigma)then
2935-
29362961!Initial computations, dynamic q- tendency
29372962 if (first_time_step .and. (.not. restart
29382963 & .or. sigmab_coldstart))then
0 commit comments