CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_suite_interstitial_4.F90
1
3
5
6 contains
7
11 subroutine gfs_suite_interstitial_4_run (im, levs, ltaerosol, tracers_total, ntrac, ntcw, ntiw, ntclamt, &
12 ntrw, ntsw, ntrnc, ntsnc, ntgl, ntgnc, ntlnc, ntinc, ntccn, nn, imp_physics, imp_physics_gfdl, imp_physics_thompson, &
13 imp_physics_nssl, nssl_invertccn, nssl_ccn_on, &
14 imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, convert_dry_rho, dtf, save_qc, save_qi, con_pi, dtidx, dtend,&
15 index_of_process_conv_trans, gq0, clw, prsl, save_tcp, con_rd, con_eps, nssl_cccn, nwfa, spechum, ldiag3d, &
16 qdiag3d, save_lnc, save_inc, ntk, ntke, otsptflag, errmsg, errflg)
17
18 use machine, only: kind_phys
20
21 implicit none
22
23 ! interface variables
24
25 logical, intent(in) :: otsptflag(:)! on/off switch for tracer transport by updraft and
26 integer, intent(in ) :: im, levs, tracers_total, ntrac, ntcw, ntiw, ntclamt, ntrw, &
27 ntsw, ntrnc, ntsnc, ntgl, ntgnc, ntlnc, ntinc, ntccn, nn, imp_physics, imp_physics_gfdl, imp_physics_thompson, &
28 imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, imp_physics_nssl
29
30 logical, intent(in) :: ltaerosol, convert_dry_rho
31 logical, intent(in) :: nssl_ccn_on, nssl_invertccn
32
33 real(kind=kind_phys), intent(in ) :: con_pi, dtf
34 real(kind=kind_phys), intent(in ), dimension(:,:) :: save_qc
35 ! save_qi is not allocated for Zhao-Carr MP
36 real(kind=kind_phys), intent(in ), dimension(:,:) :: save_qi, save_lnc, save_inc
37
38 ! dtend and dtidx are only allocated if ldiag3d
39 logical, intent(in) :: ldiag3d, qdiag3d
40 real(kind=kind_phys), dimension(:,:,:), intent(inout), optional :: dtend
41 integer, dimension(:,:), intent(in) :: dtidx
42 integer, intent(in) :: index_of_process_conv_trans,ntk,ntke
43
44 real(kind=kind_phys), dimension(:,:,:), intent(inout) :: gq0
45 real(kind=kind_phys), dimension(:,:,:), intent(inout) :: clw
46 real(kind=kind_phys), dimension(:,:), intent(in) :: prsl
47 real(kind=kind_phys), intent(in) :: con_rd, con_eps, nssl_cccn
48 real(kind=kind_phys), dimension(:,:), intent(in), optional :: nwfa
49 real(kind=kind_phys), dimension(:,:), intent(in) :: save_tcp
50 real(kind=kind_phys), dimension(:,:), intent(in) :: spechum
51
52 character(len=*), intent( out) :: errmsg
53 integer, intent( out) :: errflg
54
55 ! local variables
56 real(kind=kind_phys), parameter :: zero = 0.0_kind_phys, one = 1.0_kind_phys
57 integer :: i,k,n,tracers,idtend
58 real(kind=kind_phys) :: liqm, icem, xccn, xcwmas, xccw, xcimas, qccn
59
60 real(kind=kind_phys) :: rho, orho
61 real(kind=kind_phys), dimension(im,levs) :: qv_mp
62 real(kind=kind_phys), dimension(im,levs) :: qc_mp
63 real(kind=kind_phys), dimension(im,levs) :: qi_mp
64 real(kind=kind_phys), dimension(im,levs) :: nc_mp
65 real(kind=kind_phys), dimension(im,levs) :: ni_mp
66
67 ! Initialize CCPP error handling variables
68 errmsg = ''
69 errflg = 0
70
71 ! This code was previously in GFS_SCNV_generic_post, but it really belongs
72 ! here, because it fixes the convective transportable_tracers mess for Zhao-Carr
73 ! and GFDL MP from GFS_suite_interstitial_3. This whole code around clw(:,:,2)
74 ! being set to -999 for Zhao-Carr MP (which doesn't have cloud ice) and GFDL-MP
75 ! (which does have cloud ice, but for some reason it was decided to code it up
76 ! in the same way as for Zhao-Carr, nowadays unnecessary and confusing) needs
77 ! to be cleaned up. The convection schemes doing something different internally
78 ! based on clw(i,k,2) being -999.0 or not is not a good idea.
79 do k=1,levs
80 do i=1,im
81 if (clw(i,k,2) <= -999.0) clw(i,k,2) = 0.0
82 enddo
83 enddo
84
85 if(ldiag3d) then
86 if(ntk>0 .and. ntk<=size(clw,3)) then
87 idtend=dtidx(100+ntke,index_of_process_conv_trans)
88 if(idtend>=1) then
89 dtend(:,:,idtend) = dtend(:,:,idtend) + clw(:,:,ntk)-gq0(:,:,ntk)
90 endif
91 endif
92 if(ntcw>0) then
93 if (imp_physics == imp_physics_zhao_carr .or. &
94 imp_physics == imp_physics_zhao_carr_pdf .or. &
95 imp_physics == imp_physics_gfdl) then
96 idtend=dtidx(100+ntcw,index_of_process_conv_trans)
97 if(idtend>=1) then
98 dtend(:,:,idtend) = dtend(:,:,idtend) + clw(:,:,1)+clw(:,:,2) - gq0(:,:,ntcw)
99 endif
100 else if(ntiw>0) then
101 idtend=dtidx(100+ntiw,index_of_process_conv_trans)
102 if(idtend>=1) then
103 dtend(:,:,idtend) = dtend(:,:,idtend) + clw(:,:,1)-gq0(:,:,ntiw)
104 endif
105 idtend=dtidx(100+ntcw,index_of_process_conv_trans)
106 if(idtend>=1) then
107 dtend(:,:,idtend) = dtend(:,:,idtend) + clw(:,:,2)-gq0(:,:,ntcw)
108 endif
109 else
110 idtend=dtidx(100+ntcw,index_of_process_conv_trans)
111 if(idtend>=1) then
112 dtend(:,:,idtend) = dtend(:,:,idtend) + clw(:,:,1)+clw(:,:,2) - gq0(:,:,ntcw)
113 endif
114 endif
115 endif
116 endif
117
118! --- update the tracers due to deep & shallow cumulus convective transport
119! (except for suspended water and ice)
120
121 if (tracers_total > 0) then
122 tracers = 2
123 do n=2,ntrac
124! if ( n /= ntcw .and. n /= ntiw .and. n /= ntclamt) then
125! if ( n /= ntcw .and. n /= ntiw .and. n /= ntclamt .and. &
126! n /= ntrw .and. n /= ntsw .and. n /= ntrnc .and. &
127! n /= ntsnc .and. n /= ntgl .and. n /= ntgnc &
128! .and. &
129! n /= nthl .and. n /= nthnc .and. n /= ntgv .and. &
130! n /= nthv .and. n /= ntccn &
131! ) then
132 IF ( otsptflag(n) ) THEN
133 tracers = tracers + 1
134 if(n/=ntk .and. n/=ntlnc .and. n/=ntinc .and. n /= ntcw .and. n /= ntiw) then
135 idtend=dtidx(100+n,index_of_process_conv_trans)
136 if(idtend>=1) then
137 dtend(:,:,idtend) = dtend(:,:,idtend) + clw(:,:,tracers)-gq0(:,:,n)
138 endif
139 endif
140 do k=1,levs
141 do i=1,im
142 gq0(i,k,n) = clw(i,k,tracers)
143 enddo
144 enddo
145 endif
146 enddo
147 endif
148
149 if (ntcw > 0) then
150
151! for microphysics
152 if (imp_physics == imp_physics_zhao_carr .or. &
153 imp_physics == imp_physics_zhao_carr_pdf .or. &
154 imp_physics == imp_physics_gfdl) then
155 gq0(1:im,:,ntcw) = clw(1:im,:,1) + clw(1:im,:,2)
156
157 elseif (ntiw > 0) then
158 do k=1,levs
159 do i=1,im
160 gq0(i,k,ntiw) = clw(i,k,1) ! ice
161 gq0(i,k,ntcw) = clw(i,k,2) ! water
162 enddo
163 enddo
164
165 if ( imp_physics == imp_physics_nssl ) then
166 liqm = con_pi/6.*1.e3*(18.e-6)**3 ! 4./3.*con_pi*1.e-12
167 icem = con_pi/6.*1.e3*(120.e-6)**3 ! 4./3.*con_pi*3.2768*1.e-14*890.
168 qccn = nssl_cccn/1.225 !1.225 is a reference air density and should match what is used in module_mp_nssl_2mom.F90 (rho00)
169 do k=1,levs
170 do i=1,im
171 ! check number of available ccn
172 IF ( nssl_ccn_on ) THEN
173 IF ( nssl_invertccn ) THEN
174 xccn = qccn - gq0(i,k,ntccn)
175 ELSE
176 xccn = gq0(i,k,ntccn)
177 ENDIF
178 ELSE
179 xccn = max(0.0, qccn - gq0(i,k,ntlnc))
180 ENDIF
181
182 IF ( gq0(i,k,ntlnc) > 0.0 .and. save_qc(i,k) > 0.0 ) THEN
183 xcwmas = max( liqm, clw(i,k,2)/gq0(i,k,ntlnc) )
184 ELSE
185 xcwmas = liqm
186 ENDIF
187
188 IF ( gq0(i,k,ntinc) > 0.0 .and. save_qi(i,k) > 0.0 ) THEN
189 xcimas = max( liqm, clw(i,k,1)/gq0(i,k,ntinc) )
190 ELSE
191 xcimas = icem
192 ENDIF
193
194 IF ( xccn > 0.0 ) THEN
195 xccw = min( xccn, max(0.0, (clw(i,k,2)-save_qc(i,k))) / xcwmas )
196 gq0(i,k,ntlnc) = gq0(i,k,ntlnc) + xccw
197 IF ( nssl_ccn_on ) THEN
198 IF ( nssl_invertccn ) THEN
199 ! ccn are activated CCN, so add
200 gq0(i,k,ntccn) = gq0(i,k,ntccn) + xccw
201 ELSE
202 ! ccn are unactivated CCN, so subtract
203 gq0(i,k,ntccn) = gq0(i,k,ntccn) - xccw
204 ENDIF
205 ENDIF
206 ENDIF
207
208 gq0(i,k,ntinc) = gq0(i,k,ntinc) &
209 + max(0.0, (clw(i,k,1)-save_qi(i,k))) / xcimas
210 enddo
211 enddo
212 endif
213
214 if (imp_physics == imp_physics_thompson .and. (ntlnc>0 .or. ntinc>0)) then
215 if_convert_dry_rho: if (convert_dry_rho) then
216 do k=1,levs
217 do i=1,im
219 qv_mp(i,k) = spechum(i,k) / (one-spechum(i,k))
221 rho = con_eps*prsl(i,k) / (con_rd*save_tcp(i,k)*(qv_mp(i,k)+con_eps))
222 orho = one/rho
223 if (ntlnc>0) then
225 qc_mp(i,k) = (clw(i,k,2)-save_qc(i,k)) / (one-spechum(i,k))
227 nc_mp(i,k) = gq0(i,k,ntlnc) / (one-spechum(i,k))
228 nc_mp(i,k) = max(zero, nc_mp(i,k) + make_dropletnumber(qc_mp(i,k) * rho, nwfa(i,k)*rho) * orho)
230 gq0(i,k,ntlnc) = nc_mp(i,k) / (one+qv_mp(i,k))
231 endif
232 if (ntinc>0) then
234 qi_mp(i,k) = (clw(i,k,1)-save_qi(i,k)) / (one-spechum(i,k))
236 ni_mp(i,k) = gq0(i,k,ntinc) / (one-spechum(i,k))
237 ni_mp(i,k) = max(zero, ni_mp(i,k) + make_icenumber(qi_mp(i,k) * rho, save_tcp(i,k)) * orho)
239 gq0(i,k,ntinc) = ni_mp(i,k) / (one+qv_mp(i,k))
240 endif
241 enddo
242 enddo
243 else
244 do k=1,levs
245 do i=1,im
247 rho = con_eps*prsl(i,k) / (con_rd*save_tcp(i,k)*(spechum(i,k)+con_eps))
248 orho = one/rho
249 if (ntlnc>0) then
251 qc_mp(i,k) = (clw(i,k,2)-save_qc(i,k))
253 gq0(i,k,ntlnc) = max(zero, gq0(i,k,ntlnc) + make_dropletnumber(qc_mp(i,k) * rho, nwfa(i,k)*rho) * orho)
254 endif
255 if (ntinc>0) then
257 qi_mp(i,k) = (clw(i,k,1)-save_qi(i,k))
259 gq0(i,k,ntinc) = max(zero, gq0(i,k,ntinc) + make_icenumber(qi_mp(i,k) * rho, save_tcp(i,k)) * orho)
260 endif
261 enddo
262 enddo
263 end if if_convert_dry_rho
264 if(ldiag3d .and. qdiag3d) then
265 idtend = dtidx(100+ntlnc,index_of_process_conv_trans)
266 if(idtend>0) then
267 dtend(:,:,idtend) = dtend(:,:,idtend) + gq0(:,:,ntlnc) - save_lnc
268 endif
269 idtend = dtidx(100+ntinc,index_of_process_conv_trans)
270 if(idtend>0) then
271 dtend(:,:,idtend) = dtend(:,:,idtend) + gq0(:,:,ntinc) - save_inc
272 endif
273 endif
274 endif
275
276 else
277 do k=1,levs
278 do i=1,im
279 gq0(i,k,ntcw) = clw(i,k,1) + clw(i,k,2)
280 enddo
281 enddo
282 endif ! end if_ntiw
283
284 else
285 do k=1,levs
286 do i=1,im
287 clw(i,k,1) = clw(i,k,1) + clw(i,k,2)
288 enddo
289 enddo
290 endif ! end if_ntcw
291
292 end subroutine gfs_suite_interstitial_4_run
293
294 end module gfs_suite_interstitial_4
elemental real function, public make_icenumber(q_ice, temp)
Table of lookup values of radiative effective radius of ice crystals as a function of Temperature fro...