CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_suite_interstitial_3.F90
1
3
5
6 contains
7
11 subroutine gfs_suite_interstitial_3_run (otsptflag, &
12 im, levs, nn, cscnv,imfshalcnv, imfdeepcnv, &
13 imfshalcnv_samf, imfdeepcnv_samf, imfdeepcnv_c3, &
14 imfshalcnv_c3,progsigma, &
15 first_time_step, restart, &
16 satmedmf, trans_trac, do_shoc, ltaerosol, ntrac, ntcw, &
17 ntiw, ntclamt, ntrw, ntsw, ntrnc, ntsnc, ntgl, ntgnc, &
18 xlon, xlat, gt0, gq0, sigmain,sigmaout,qmicro, &
19 imp_physics, imp_physics_mg, &
20 imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, &
21 imp_physics_gfdl, imp_physics_thompson, dtidx, ntlnc, &
22 imp_physics_wsm6, imp_physics_fer_hires, prsi, ntinc, &
23 imp_physics_nssl, &
24 prsl, prslk, rhcbot,rhcpbl, rhctop, rhcmax, islmsk, &
25 work1, work2, kpbl, kinver, ras, me, save_lnc, save_inc, &
26 ldiag3d, qdiag3d, index_of_process_conv_trans, &
27 clw, rhc, save_qc, save_qi, save_tcp, errmsg, errflg)
28
29 use machine, only: kind_phys
30
31 implicit none
32
33 ! interface variables
34 logical, intent(in) :: otsptflag(:)! on/off switch for tracer transport (size ntrac)
35 integer, intent(in ) :: im, levs, nn, ntrac, ntcw, ntiw, ntclamt, ntrw, ntsw,&
36 ntrnc, ntsnc, ntgl, ntgnc, imp_physics, imp_physics_mg, imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, &
37 imp_physics_gfdl, imp_physics_thompson, imp_physics_wsm6,imp_physics_fer_hires, &
38 imp_physics_nssl, me, index_of_process_conv_trans
39 integer, intent(in ), dimension(:) :: islmsk, kpbl, kinver
40 logical, intent(in ) :: cscnv, satmedmf, trans_trac, do_shoc, ltaerosol, ras, progsigma
41 logical, intent(in ) :: first_time_step, restart
42 integer, intent(in ) :: imfshalcnv, imfdeepcnv, imfshalcnv_samf,imfdeepcnv_samf
43 integer, intent(in ) :: imfshalcnv_c3,imfdeepcnv_c3
44 integer, intent(in) :: ntinc, ntlnc
45 logical, intent(in) :: ldiag3d, qdiag3d
46 integer, dimension(:,:), intent(in) :: dtidx
47 real, dimension(:,:), intent(out) :: save_lnc, save_inc
48
49 real(kind=kind_phys), intent(in ) :: rhcbot, rhcmax, rhcpbl, rhctop
50 real(kind=kind_phys), intent(in ), dimension(:) :: work1, work2
51 real(kind=kind_phys), intent(in ), dimension(:,:) :: prsl, prslk
52 real(kind=kind_phys), intent(in ), dimension(:,:) :: prsi
53 real(kind=kind_phys), intent(in ), dimension(:) :: xlon, xlat
54 real(kind=kind_phys), intent(in ), dimension(:,:) :: gt0
55 real(kind=kind_phys), intent(in ), dimension(:,:,:) :: gq0
56
57 real(kind=kind_phys), intent(inout ), dimension(:,:), optional :: sigmain
58 real(kind=kind_phys), intent(inout ), dimension(:,:), optional :: sigmaout, qmicro
59 real(kind=kind_phys), intent(inout), dimension(:,:) :: rhc, save_qc
60 ! save_qi is not allocated for Zhao-Carr MP
61 real(kind=kind_phys), intent(inout), dimension(:,:) :: save_qi
62 real(kind=kind_phys), intent(inout), dimension(:,:) :: save_tcp
63 real(kind=kind_phys), intent(inout), dimension(:,:,:) :: clw
64
65 character(len=*), intent( out) :: errmsg
66 integer, intent( out) :: errflg
67
68 ! local variables
69 integer :: i,k,n,tracers,kk
70 real(kind=kind_phys) :: tem, tem1, tem2
71 real(kind=kind_phys), dimension(im) :: tx1, tx2, tx3, tx4
72
73 !real(kind=kind_phys),parameter :: slope_mg = 0.02, slope_upmg = 0.04, &
74 ! turnrhcrit = 0.900, turnrhcrit_upper = 0.150
75 ! in the following inverse of slope_mg and slope_upmg are specified
76 real(kind=kind_phys), parameter :: zero = 0.0_kind_phys, one = 1.0_kind_phys
77 real(kind=kind_phys), parameter :: slope_mg = 50.0_kind_phys, &
78 slope_upmg = 25.0_kind_phys
79
80 ! Initialize CCPP error handling variables
81 errmsg = ''
82 errflg = 0
83
84 ! In case of using prognostic updraf area fraction, initialize area fraction here
85 ! since progsigma_calc is called from both deep and shallow schemes.
86 if(((imfshalcnv == imfshalcnv_samf) .or. (imfdeepcnv == imfdeepcnv_samf) &
87 .or. (imfshalcnv == imfshalcnv_c3) .or. (imfdeepcnv == imfdeepcnv_c3)) &
88 .and. progsigma)then
89 if(first_time_step .and. .not. restart)then
90 do k=1,levs
91 do i=1,im
92 sigmain(i,k)=0.0
93 sigmaout(i,k)=0.0
94 qmicro(i,k)=0.0
95 enddo
96 enddo
97 endif
98 do k=1,levs
99 do i=1,im
100 sigmaout(i,k)=0.0
101 enddo
102 enddo
103 endif
104
105
106 if (cscnv .or. satmedmf .or. trans_trac .or. ras) then
107 tracers = 2
108 do n=2,ntrac
109! if ( n /= ntcw .and. n /= ntiw .and. n /= ntclamt .and. &
110! n /= ntrw .and. n /= ntsw .and. n /= ntrnc .and. &
111! n /= ntsnc .and. n /= ntgl .and. n /= ntgnc) then
112 IF ( otsptflag(n) ) THEN
113 tracers = tracers + 1
114 do k=1,levs
115 do i=1,im
116 clw(i,k,tracers) = gq0(i,k,n)
117 enddo
118 enddo
119 endif
120 enddo
121 endif ! end if_ras or cfscnv or samf
122
123 if (ntcw > 0) then
124 if (imp_physics == imp_physics_mg .and. rhcpbl < 0.5_kind_phys) then ! compute rhc for GMAO macro physics cloud pdf
125 do i=1,im
126 tx1(i) = one / prsi(i,1)
127 tx2(i) = one - rhcmax*work1(i)-rhcbot*work2(i)
128
129 kk = min(kinver(i), max(2,kpbl(i)))
130 tx3(i) = prsi(i,kk)*tx1(i)
131 tx4(i) = rhcpbl - rhctop*abs(cos(xlat(i)))
132 enddo
133 do k = 1, levs
134 do i = 1, im
135 tem = prsl(i,k) * tx1(i)
136 tem1 = min(max((tem-tx3(i))*slope_mg, -20.0_kind_phys), 20.0_kind_phys)
137 ! Using rhcpbl and rhctop from the namelist instead of 0.3 and 0.2
138 ! and rhcbot represents pbl top critical relative humidity
139 tem2 = min(max((tx4(i)-tem)*slope_upmg, -20.0_kind_phys), 20.0_kind_phys) ! Anning
140 if (islmsk(i) > 0) then
141 tem1 = one / (one+exp(tem1+tem1))
142 else
143 tem1 = 2.0_kind_phys / (one+exp(tem1+tem1))
144 endif
145 tem2 = one / (one+exp(tem2))
146
147 rhc(i,k) = min(rhcmax, max(0.7_kind_phys, one-tx2(i)*tem1*tem2))
148 enddo
149 enddo
150 else
151 do k=1,levs
152 do i=1,im
153 kk = max(10,kpbl(i))
154#ifdef SINGLE_PREC
155 if (k < kk) then
156 tem = rhcbot - (rhcbot-rhcpbl) * (one-prslk(i,k)) / max(one-prslk(i,kk),1e-7)
157 else
158 tem = rhcpbl - (rhcpbl-rhctop) * (prslk(i,kk)-prslk(i,k)) / max(prslk(i,kk),1e-7)
159 endif
160#else
161 if (k < kk) then
162 tem = rhcbot - (rhcbot-rhcpbl) * (one-prslk(i,k)) / (one-prslk(i,kk))
163 else
164 tem = rhcpbl - (rhcpbl-rhctop) * (prslk(i,kk)-prslk(i,k)) / prslk(i,kk)
165 endif
166#endif
167 tem = rhcmax * work1(i) + tem * work2(i)
168 rhc(i,k) = max(zero, min(one,tem))
169 enddo
170 enddo
171 endif
172 else
173 rhc(:,:) = 1.0
174 endif
175
176 if (imp_physics == imp_physics_zhao_carr .or. imp_physics == imp_physics_zhao_carr_pdf) then ! zhao-carr microphysics
177 !GF* move to GFS_MP_generic_pre (from gscond/precpd)
178 ! do i=1,im
179 ! psautco_l(i) = Model%psautco(1)*work1(i) + Model%psautco(2)*work2(i)
180 ! prautco_l(i) = Model%prautco(1)*work1(i) + Model%prautco(2)*work2(i)
181 ! enddo
182 !*GF
183 do k=1,levs
184 do i=1,im
185 clw(i,k,1) = gq0(i,k,ntcw)
186 enddo
187 enddo
188 elseif (imp_physics == imp_physics_gfdl) then
189 clw(1:im,:,1) = gq0(1:im,:,ntcw)
190 elseif (imp_physics == imp_physics_thompson) then
191 do k=1,levs
192 do i=1,im
193 clw(i,k,1) = gq0(i,k,ntiw) ! ice
194 clw(i,k,2) = gq0(i,k,ntcw) ! water
195 save_tcp(i,k) = gt0(i,k)
196 enddo
197 enddo
198 if(ltaerosol) then
199 save_qi(:,:) = clw(:,:,1)
200 save_qc(:,:) = clw(:,:,2)
201 else
202 save_qi(:,:) = clw(:,:,1)
203 endif
204 else if (imp_physics == imp_physics_nssl ) then
205 do k=1,levs
206 do i=1,im
207 clw(i,k,1) = gq0(i,k,ntiw) ! cloud ice
208 clw(i,k,2) = gq0(i,k,ntcw) ! cloud droplets
209 enddo
210 enddo
211 save_qi(:,:) = clw(:,:,1)
212 save_qc(:,:) = clw(:,:,2)
213 elseif (imp_physics == imp_physics_wsm6 .or. imp_physics == imp_physics_mg .or. imp_physics == imp_physics_fer_hires) then
214 do k=1,levs
215 do i=1,im
216 clw(i,k,1) = gq0(i,k,ntiw) ! ice
217 clw(i,k,2) = gq0(i,k,ntcw) ! water
218 enddo
219 enddo
220 endif
221
222 if(imp_physics == imp_physics_thompson .and. ldiag3d .and. qdiag3d) then
223 if(dtidx(100+ntlnc,index_of_process_conv_trans)>0) then
224 save_lnc = gq0(:,:,ntlnc)
225 endif
226 if(dtidx(100+ntinc,index_of_process_conv_trans)>0) then
227 save_inc = gq0(:,:,ntinc)
228 endif
229 endif
230
231 end subroutine gfs_suite_interstitial_3_run
232
233 end module gfs_suite_interstitial_3