@@ -35,7 +35,7 @@ module satmedmfvdifq
3535!!
3636!! Incorporate the TTE- EDMF; if (tte_edmf= .true. ),
3737!! TKE- EDMF scheme becomes TTE- EDMF scheme and the variable ' te'
38- !! is read as TTE; if (tte_edmf= .false. ), the variable ' te' is
38+ !! is read as TTE; if (tte_edmf= .false. ), the variable ' te' is
3939!! read as TKE, 5 / 22 / 2025
4040!!
4141!!
@@ -85,25 +85,25 @@ end subroutine satmedmfvdifq_init
8585!! (mfpbltq.f), is represented using a mass flux approach (Siebesma et al.(2007 ) \cite Siebesma_2007 ).
8686!! - # A mass- flux approach is also used to represent the stratocumulus- top- induced turbulence
8787!! (mfscuq.f).
88- !! \section detail_satmedmfvidfq GFS satmedmfvdifq Detailed Algorithm
88+ !! \section detail_satmedmfvidfq GFS satmedmfvdifq Detailed Algorithm
8989 subroutine satmedmfvdifq_run (im ,km ,ntrac ,ntcw ,ntrw , &
9090 & ntiw ,ntke ,grav ,pi ,rd ,cp ,rv ,hvap ,hfus ,fv ,eps ,epsm1 , &
9191!The following three variables are for SA-3D-TKE
9292 & def_1 ,def_2 ,def_3 ,sa3dtke ,dku3d_h ,dku3d_e , &
93- & dv ,du ,tdt ,rtg ,u1 ,v1 ,t1 ,q1 ,usfco ,vsfco ,icplocn2atm , &
93+ & dv ,du ,tdt ,rtg ,u1 ,v1 ,t1 ,q1 ,usfco ,vsfco ,use_oceanuv , &
9494 & swh ,hlw ,xmu ,garea ,zvfun ,sigmaf , &
9595 & psk ,rbsoil ,zorl ,u10m ,v10m ,fm ,fh , &
9696 & tsea ,heat ,evap ,stress ,spd1 ,kpbl , &
9797 & prsi ,del ,prsl ,prslk ,phii ,phil ,delt ,tte_edmf , &
9898 & dspheat ,dusfc ,dvsfc ,dtsfc ,dqsfc ,hpbl ,dkt ,dku ,tkeh , &
9999 & kinver ,xkzm_m ,xkzm_h ,xkzm_s ,dspfac ,bl_upfr ,bl_dnfr , &
100100 & rlmx ,elmx ,sfc_rlm ,tc_pbl ,use_lpt , &
101- !IVAI: canopy inputs from AQM
101+ !IVAI: canopy inputs from AQM
102102 & do_canopy , cplaqm , claie , cfch , cfrt , cclu , cpopu , &
103103!IVAI
104104 & ntqv ,dtend ,dtidx ,index_of_temperature ,index_of_x_wind , &
105105 & index_of_y_wind ,index_of_process_pbl ,gen_tend ,ldiag3d , &
106- & errmsg ,errflg )
106+ & errmsg ,errflg )
107107!IVAI: aux arrays
108108! & naux2d,naux3d,aux2d,aux3d)
109109
@@ -158,7 +158,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
158158 & dtend
159159 integer , intent (in ) :: dtidx(:,:), index_of_temperature, &
160160 & index_of_x_wind, index_of_y_wind, index_of_process_pbl
161- integer , intent (in ) :: icplocn2atm
161+ logical , intent (in ) :: use_oceanuv
162162 real (kind= kind_phys), intent (out ) :: &
163163 & dusfc(:), dvsfc(:), &
164164 & dtsfc(:), dqsfc(:), &
@@ -192,7 +192,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
192192 integer lcld(im),kcld(im),krad(im),mrad(im)
193193 integer kx1(im), kb1(im), kpblx(im)
194194!
195- real (kind= kind_phys) te(im,km), tei(im,km-1 ), tke(im,km),
195+ real (kind= kind_phys) te(im,km), tei(im,km-1 ), tke(im,km),
196196 & tteh(im,km), tesq(im,km-1 ),e2(im,0 :km)
197197!
198198 real (kind= kind_phys) theta(im,km),thvx(im,km), thlvx(im,km),
@@ -227,7 +227,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
227227!
228228 real (kind= kind_phys) elm(im,km), ele(im,km),
229229 & ckz(im,km), chz(im,km),
230- & diss(im,km-1 ),prod(im,km-1 ),
230+ & diss(im,km-1 ),prod(im,km-1 ),
231231 & bf(im,km-1 ), shr2(im,km-1 ), wush(im,km),
232232 & xlamue(im,km-1 ), xlamde(im,km-1 ),
233233 & gotvx(im,km), rlam(im,km-1 )
@@ -317,7 +317,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
317317 real (kind= kind_phys) h1
318318
319319 real (kind= kind_phys) bfac, mffac
320-
320+
321321 real (kind= kind_phys) qice(im,km),qliq(im,km)
322322
323323!PCC_CANOPY------------------------------------
@@ -438,7 +438,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
438438 km1 = km - 1
439439 kmpbl = km / 2
440440 kmscu = km / 2
441- !> - Compute physical height of the layer centers and interfaces from
441+ !> - Compute physical height of the layer centers and interfaces from
442442!! the geopotential height (\p zi and \p zl)
443443 do k= 1 ,km
444444 do i= 1 ,im
@@ -466,7 +466,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
466466 do i= 1 ,im
467467 gdx(i) = sqrt (garea(i))
468468 enddo
469- !> - Initialize tke value at vertical layer centers and interfaces
469+ !> - Initialize tke value at vertical layer centers and interfaces
470470!! from tracer (\p tke and \p tkeh)
471471 do k= 1 ,km
472472 do i= 1 ,im
@@ -583,11 +583,11 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
583583 kcld(i) = km1
584584 endif
585585 enddo
586- !
586+ !
587587!> - Compute a function for green vegetation fraction and surface roughness.
588588!! Entrainment rate in updraft is a function of vegetation fraction and surface
589589!! roughness length
590- !
590+ !
591591 do i = 1 ,im
592592 tem = (sigmaf(i) - vegflo) / (vegfup - vegflo)
593593 tem = min (max (tem, 0 .), 1 .)
@@ -598,7 +598,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
598598 enddo
599599!
600600!> - Compute \f$\theta\f$(theta), and \f$q_l\f$(qlx), \f$\theta_e\f$(thetae),
601- !! \f$\theta_v\f$(thvx),\f$\theta_{l,v}\f$ (thlvx) including ice water
601+ !! \f$\theta_v\f$(thvx),\f$\theta_{l,v}\f$ (thlvx) including ice water
602602 do k= 1 ,km
603603 do i= 1 ,im
604604 pix(i,k) = psk(i) / prslk(i,k)
@@ -637,7 +637,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
637637 enddo
638638 enddo
639639!
640- !> - Compute an empirical cloud fraction based on
640+ !> - Compute an empirical cloud fraction based on
641641!! Xu and Randall (1996 ) \cite xu_and_randall_1996
642642 do k = 1 , km
643643 do i = 1 , im
@@ -688,7 +688,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
688688!
689689!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
690690!
691- !> - Initialize diffusion coefficients to 0 and calculate the total
691+ !> - Initialize diffusion coefficients to 0 and calculate the total
692692!! radiative heating rate (dku, dkt, radx)
693693 do k= 1 ,km
694694 do i= 1 ,im
@@ -705,7 +705,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
705705 enddo
706706!> - Compute stable/ unstable PBL flag (pblflg) based on the total
707707!! surface energy flux (\e false if the total surface energy flux
708- !! is into the surface)
708+ !! is into the surface)
709709 do i = 1 ,im
710710 sflux(i) = heat(i) + evap(i)* fv* theta(i,1 )
711711 if (.not. sfcflg(i) .or. sflux(i) <= 0 .) pblflg(i)= .false.
@@ -716,12 +716,12 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
716716!! - Compute critical bulk Richardson number (\f$Rb_{cr}\f$) (crb)
717717!! - For the unstable PBL, crb is a constant (0.25 )
718718!! - For the stable boundary layer (SBL), \f$Rb_{cr}\f$ varies
719- !! with the surface Rossby number, \f$R_{0 }\f$, as given by
719+ !! with the surface Rossby number, \f$R_{0 }\f$, as given by
720720!! Vickers and Mahrt (2004 ) \cite Vickers_2004
721721!! \f[
722722!! Rb_{cr}= 0.16 (10 ^{- 7 }R_{0 })^{- 0.18 }
723723!! \f]
724- !! \f[
724+ !! \f[
725725!! R_{0 }= \frac{U_{10 }}{f_{0 }z_{0 }}
726726!! \f]
727727!! where \f$U_{10 }\f$ is the wind speed at 10m above the ground surface,
@@ -753,7 +753,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
753753 ustar(i) = sqrt (stress(i))
754754 enddo
755755!
756- !> - Compute buoyancy \f$\frac{\partial \theta_v}{\partial z}\f$ (bf)
756+ !> - Compute buoyancy \f$\frac{\partial \theta_v}{\partial z}\f$ (bf)
757757!! and the wind shear squared (shr2)
758758!
759759 do k = 1 , km1
@@ -980,7 +980,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
980980!
981981!> - The \f$z/L\f$ (zol) is used as the stability criterion for the PBL.Currently,
982982!! strong unstable (convective) PBL for \f$z/L < -0.02\f$ and weakly and moderately
983- !! unstable PBL for \f$0>z/L>-0.02\f$
983+ !! unstable PBL for \f$0>z/L>-0.02\f$
984984!> - Compute the velocity scale \f$w_s\f$ (wscale) (eqn 22 of Han et al. 2019). It
985985!! is represented by the value scaled at the top of the surface layer:
986986!! \f[
@@ -1172,11 +1172,11 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
11721172 do i = 1, im
11731173 if(scuflg(i) .and. kcld(i)==km1) scuflg(i)=.false.
11741174 enddo
1175- !> - Starting at the PBL top and going downward, if the level is less
1175+ !> - Starting at the PBL top and going downward, if the level is less
11761176!! than the cloud top, find the level of the minimum radiative heating
11771177!! rate wihin the cloud. If the level of the minimum is the lowest model
11781178!! level or the minimum radiative heating rate is positive, then set
1179- !! scuflg to F.
1179+ !! scuflg to F.
11801180 do i = 1, im
11811181 flg(i)=scuflg(i)
11821182 enddo
@@ -1202,7 +1202,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
12021202!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
12031203!> ## Compute components for mass flux mixing by large thermals
12041204!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1205- !> - If the PBL is convective, the updraft properties are initialized
1205+ !> - If the PBL is convective, the updraft properties are initialized
12061206!! to be the same as the state variables.
12071207 do k = 1, km
12081208 do i = 1, im
@@ -1236,7 +1236,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
12361236 & pcnvflg,zl,zm,q1,t1,u1,v1,plyr,pix,thlx,thvx,
12371237 & gdx,hpbl,kpbl,vpert,buou,wush,tkemean,vez0fun,xmf,
12381238 & tcko,qcko,ucko,vcko,xlamue,bl_upfr)
1239- !> - Call mfscuq(), which is a new mass-flux parameterization for
1239+ !> - Call mfscuq(), which is a new mass-flux parameterization for
12401240!! stratocumulus-top-induced turbulence mixing. For details of the mfscuq subroutine, step into its documentation ::mfscuq
12411241 call mfscuq(im,im,km,kmscu,ntcw,ntrac1,dt2,
12421242 & scuflg,zl,zm,q1,t1,u1,v1,plyr,pix,
@@ -1389,7 +1389,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
13891389!! l_2=min(l_{up},l_{down})
13901390!!\f]
13911391!! and dissipation length scale \f$l_d\f$ is given by:
1392- !!\f[
1392+ !!\f[
13931393!! l_d=(l_{up}l_{down})^{1/2}
13941394!!\f]
13951395!! where \f$l_{up}\f$ and \f$l_{down}\f$ are the distances that a parcel
@@ -1413,7 +1413,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
14131413!
14141414 enddo
14151415 enddo
1416- !> - Compute the surface layer length scale (\f$l_1\f$) following
1416+ !> - Compute the surface layer length scale (\f$l_1\f$) following
14171417!! Nakanishi (2001) \cite Nakanish_2001 (eqn 9 of Han et al.(2019) \cite Han_2019)
14181418 do k = 1, km1
14191419 do i = 1, im
@@ -1581,7 +1581,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
15811581 enddo
15821582 enddo
15831583!
1584- ! eddy diffusivities at model interface (zm level) in LES scale
1584+ ! eddy diffusivities at model interface (zm level) in LES scale
15851585!
15861586 do k = 1, km1
15871587 do i = 1, im
@@ -1614,7 +1614,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
16141614 do k = 2, kmpbl
16151615 do i = 1, im
16161616 if(tke(i,k) > tkemax(i)) then
1617- tkemax(i) = tke(i,k)
1617+ tkemax(i) = tke(i,k)
16181618 ktkemax(i) = k
16191619 endif
16201620 enddo
@@ -1673,7 +1673,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
16731673 enddo
16741674 enddo
16751675!
1676- ! eddy diffusivities at model layer (zl level) in LES scale
1676+ ! eddy diffusivities at model layer (zl level) in LES scale
16771677!
16781678 do k = 1, km1
16791679 do i = 1, im
@@ -1717,7 +1717,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
17171717
17181718!PCC_CANOPY------------------------------------
17191719 kount=0 !IVAI
1720- if (do_canopy .and. cplaqm) then
1720+ if (do_canopy .and. cplaqm) then
17211721
17221722!IVAI
17231723! print*, ' SATMEDMFVDIFQ_RUN: CLAIE = ' , claie(:)
@@ -1774,7 +1774,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
17741774! & .OR. MAX(0.0, 1.0 - XCANOPYFRT) .GT. 0.5
17751775! & .OR. POPU .GT. 10000.0
17761776! & .OR. EXP(-0.5*XCANOPYLAI*XCANOPYCLU).GT. 0.45
1777- ! & .AND. FCH .LT. 18.0 ) THEN
1777+ ! & .AND. FCH .LT. 18.0 ) THEN
17781778
17791779 dkt(i,k)= dkt(i,k)
17801780 dkq(i,k)= dkq(i,k)
@@ -1811,9 +1811,9 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
18111811 IF ( ZCAN/FCH .GT. 1.25 ) THEN !SIGMACAN = Eulerian vertical velocity variance
18121812 SIGMACAN = 1.25*ustar(i)
18131813 END IF
1814- IF ( ZCAN/FCH .GE. 0.175
1814+ IF ( ZCAN/FCH .GE. 0.175
18151815 & .AND. ZCAN/FCH .LE. 1.25 ) THEN
1816- SIGMACAN = ustar(i) * ( 0.75 +
1816+ SIGMACAN = ustar(i) * ( 0.75 +
18171817 & (0.5 * COS((PI/1.06818) *
18181818 & (1.25 - (ZCAN/FCH)))) )
18191819 END IF
@@ -1825,9 +1825,9 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
18251825 IF ( ZCAN/FCH .GT. 1.25 ) THEN
18261826 SIGMACAN = 1.0*ustar(i)
18271827 END IF
1828- IF ( ZCAN/FCH .GE. 0.175
1828+ IF ( ZCAN/FCH .GE. 0.175
18291829 & .AND. ZCAN/FCH .LE. 1.25 ) THEN
1830- SIGMACAN = ustar(i) * ( 0.625 +
1830+ SIGMACAN = ustar(i) * ( 0.625 +
18311831 & (0.375* COS((PI/1.06818) *
18321832 & (1.25 - (ZCAN/FCH)))) )
18331833 END IF
@@ -1839,12 +1839,12 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
18391839 IF ( ZCAN/FCH .GT. 1.25 ) THEN
18401840 SIGMACAN = 0.25*(4.375 - (3.75*HOL))*ustar(i)
18411841 END IF
1842- IF ( ZCAN/FCH .GE. 0.175
1842+ IF ( ZCAN/FCH .GE. 0.175
18431843 & .AND. ZCAN/FCH .LE. 1.25 ) THEN
18441844 RRCAN=4.375-(3.75*HOL)
18451845 AACAN=(0.125*RRCAN) + 0.125
18461846 BBCAN=(0.125*RRCAN) - 0.125
1847- SIGMACAN = ustar(i) * ( AACAN +
1847+ SIGMACAN = ustar(i) * ( AACAN +
18481848 & (BBCAN * COS((PI/1.06818) *
18491849 & (1.25 - (ZCAN/FCH)))) )
18501850 END IF
@@ -2051,7 +2051,7 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
20512051 enddo
20522052!
20532053!----------------------------------------------------------------------
2054- !> - First predict te due to te production & dissipation(diss)
2054+ !> - First predict te due to te production & dissipation(diss)
20552055!
20562056 if(sa3dtke) then
20572057!The following is for SA-3D-TKE
@@ -3052,10 +3052,10 @@ subroutine satmedmfvdifq_run(im,km,ntrac,ntcw,ntrw, &
30523052 enddo
30533053 enddo
30543054 do i = 1,im
3055- if(icplocn2atm == 0 ) then
3055+ if(.not. use_oceanuv ) then
30563056 dusfc(i) = -1.*rho_a(i)*stress(i)*u1(i,1)/spd1(i)
30573057 dvsfc(i) = -1.*rho_a(i)*stress(i)*v1(i,1)/spd1(i)
3058- else if (icplocn2atm ==1 ) then
3058+ else if (use_oceanuv ) then
30593059 spd1_m=sqrt( (u1(i,1)-usfco(i))**2+(v1(i,1)-vsfco(i))**2 )
30603060 dusfc(i) = -1.*rho_a(i)*stress(i)*(u1(i,1)-usfco(i))/spd1_m
30613061 dvsfc(i) = -1.*rho_a(i)*stress(i)*(v1(i,1)-vsfco(i))/spd1_m
0 commit comments