CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_PBL_generic_pre.F90
1
3
5
6 contains
7
12 subroutine gfs_pbl_generic_pre_run (im, levs, nvdiff, ntrac, rtg_ozone_index, &
13 ntqv, ntcw, ntiw, ntrw, ntsw, ntlnc, ntinc, ntrnc, ntsnc, ntgnc, &
14 ntwa, ntia, ntgl, ntoz, ntke, ntkev, nqrimef, trans_aero, ntchs, ntchm, &
15 ntccn, nthl, nthnc, ntgv, nthv, ntrz, ntgz, nthz, &
16 imp_physics, imp_physics_gfdl, imp_physics_thompson, imp_physics_wsm6, &
17 imp_physics_zhao_carr, imp_physics_mg, imp_physics_fer_hires, imp_physics_nssl, &
18 ltaerosol, mraerosol, nssl_ccn_on, nssl_hail_on, nssl_3moment, &
19 hybedmf, do_shoc, satmedmf, qgrs, vdftra, save_u, save_v, save_t, save_q, &
20 flag_for_pbl_generic_tend, ldiag3d, qdiag3d, lssav, ugrs, vgrs, tgrs, errmsg, errflg)
21
22 use machine, only : kind_phys
23 use gfs_pbl_generic_common, only : set_aerosol_tracer_index
24
25 implicit none
26
27 integer, parameter :: kp = kind_phys
28 integer, intent(out) :: rtg_ozone_index
29 integer, intent(in) :: im, levs, nvdiff, ntrac
30 integer, intent(in) :: ntqv, ntcw, ntiw, ntrw, ntsw, ntlnc, ntinc, ntrnc, ntsnc, ntgnc
31 integer, intent(in) :: ntwa, ntia, ntgl, ntoz, ntke, ntkev, nqrimef,ntchs, ntchm
32 integer, intent(in) :: ntccn, nthl, nthnc, ntgv, nthv, ntrz, ntgz, nthz
33 logical, intent(in) :: trans_aero, ldiag3d, qdiag3d, lssav
34 integer, intent(in) :: imp_physics, imp_physics_gfdl, imp_physics_thompson, imp_physics_wsm6
35 integer, intent(in) :: imp_physics_zhao_carr, imp_physics_mg, imp_physics_fer_hires
36 logical, intent(in) :: ltaerosol, hybedmf, do_shoc, satmedmf, flag_for_pbl_generic_tend, mraerosol
37 integer, intent(in) :: imp_physics_nssl
38 logical, intent(in) :: nssl_hail_on, nssl_ccn_on, nssl_3moment
39
40 real(kind=kind_phys), dimension(:,:,:), intent(in) :: qgrs
41 real(kind=kind_phys), dimension(:,:), intent(in) :: ugrs, vgrs, tgrs
42 real(kind=kind_phys), dimension(:,:, :), intent(inout) :: vdftra
43 real(kind=kind_phys), dimension(:,:), intent(out) :: save_u, save_v, save_t
44 real(kind=kind_phys), dimension(:,:, :), intent(out) :: save_q
45
46 ! CCPP error handling variables
47 character(len=*), intent(out) :: errmsg
48 integer, intent(out) :: errflg
49
50 real (kind=kind_phys), parameter :: zero = 0.0_kp, one=1.0_kp
51
52 ! Local variables
53 integer :: i, k, kk, k1, n
54
55 ! Initialize CCPP error handling variables
56 errmsg = ''
57 errflg = 0
58
59 rtg_ozone_index=-1
60!DH: dvdftra is only used if nvdiff != ntrac or (nvdiff == ntrac .and. )
61 if (nvdiff == ntrac .and. (hybedmf .or. do_shoc .or. satmedmf)) then
62 vdftra = qgrs
63 rtg_ozone_index = ntoz
64 else
65 if (imp_physics == imp_physics_wsm6) then
66 ! WSM6
67 do k=1,levs
68 do i=1,im
69 vdftra(i,k,1) = qgrs(i,k,ntqv)
70 vdftra(i,k,2) = qgrs(i,k,ntcw)
71 vdftra(i,k,3) = qgrs(i,k,ntiw)
72 vdftra(i,k,4) = qgrs(i,k,ntoz)
73 enddo
74 enddo
75 rtg_ozone_index = 4
76
77 ! Ferrier-Aligo
78 elseif (imp_physics == imp_physics_fer_hires) then
79 do k=1,levs
80 do i=1,im
81 vdftra(i,k,1) = qgrs(i,k,ntqv)
82 vdftra(i,k,2) = qgrs(i,k,ntcw)
83 vdftra(i,k,3) = qgrs(i,k,ntiw)
84 vdftra(i,k,4) = qgrs(i,k,ntrw)
85 vdftra(i,k,5) = qgrs(i,k,nqrimef)
86 vdftra(i,k,6) = qgrs(i,k,ntoz)
87 enddo
88 enddo
89 rtg_ozone_index = 6
90
91 elseif (imp_physics == imp_physics_thompson) then
92 ! Thompson
93 if(ltaerosol) then
94 do k=1,levs
95 do i=1,im
96 vdftra(i,k,1) = qgrs(i,k,ntqv)
97 vdftra(i,k,2) = qgrs(i,k,ntcw)
98 vdftra(i,k,3) = qgrs(i,k,ntiw)
99 vdftra(i,k,4) = qgrs(i,k,ntrw)
100 vdftra(i,k,5) = qgrs(i,k,ntsw)
101 vdftra(i,k,6) = qgrs(i,k,ntgl)
102 vdftra(i,k,7) = qgrs(i,k,ntlnc)
103 vdftra(i,k,8) = qgrs(i,k,ntinc)
104 vdftra(i,k,9) = qgrs(i,k,ntrnc)
105 vdftra(i,k,10) = qgrs(i,k,ntoz)
106 vdftra(i,k,11) = qgrs(i,k,ntwa)
107 vdftra(i,k,12) = qgrs(i,k,ntia)
108 enddo
109 enddo
110 rtg_ozone_index = 10
111 elseif(mraerosol) then
112 do k=1,levs
113 do i=1,im
114 vdftra(i,k,1) = qgrs(i,k,ntqv)
115 vdftra(i,k,2) = qgrs(i,k,ntcw)
116 vdftra(i,k,3) = qgrs(i,k,ntiw)
117 vdftra(i,k,4) = qgrs(i,k,ntrw)
118 vdftra(i,k,5) = qgrs(i,k,ntsw)
119 vdftra(i,k,6) = qgrs(i,k,ntgl)
120 vdftra(i,k,7) = qgrs(i,k,ntlnc)
121 vdftra(i,k,8) = qgrs(i,k,ntinc)
122 vdftra(i,k,9) = qgrs(i,k,ntrnc)
123 vdftra(i,k,10) = qgrs(i,k,ntoz)
124 enddo
125 enddo
126 rtg_ozone_index = 10
127 else
128 do k=1,levs
129 do i=1,im
130 vdftra(i,k,1) = qgrs(i,k,ntqv)
131 vdftra(i,k,2) = qgrs(i,k,ntcw)
132 vdftra(i,k,3) = qgrs(i,k,ntiw)
133 vdftra(i,k,4) = qgrs(i,k,ntrw)
134 vdftra(i,k,5) = qgrs(i,k,ntsw)
135 vdftra(i,k,6) = qgrs(i,k,ntgl)
136 vdftra(i,k,7) = qgrs(i,k,ntinc)
137 vdftra(i,k,8) = qgrs(i,k,ntrnc)
138 vdftra(i,k,9) = qgrs(i,k,ntoz)
139 enddo
140 enddo
141 rtg_ozone_index = 9
142 endif
143 ! MG
144 elseif (imp_physics == imp_physics_mg) then ! MG3/2
145 if (ntgl > 0) then ! MG3
146 do k=1,levs
147 do i=1,im
148 vdftra(i,k,1) = qgrs(i,k,ntqv)
149 vdftra(i,k,2) = qgrs(i,k,ntcw)
150 vdftra(i,k,3) = qgrs(i,k,ntiw)
151 vdftra(i,k,4) = qgrs(i,k,ntrw)
152 vdftra(i,k,5) = qgrs(i,k,ntsw)
153 vdftra(i,k,6) = qgrs(i,k,ntgl)
154 vdftra(i,k,7) = qgrs(i,k,ntlnc)
155 vdftra(i,k,8) = qgrs(i,k,ntinc)
156 vdftra(i,k,9) = qgrs(i,k,ntrnc)
157 vdftra(i,k,10) = qgrs(i,k,ntsnc)
158 vdftra(i,k,11) = qgrs(i,k,ntgnc)
159 vdftra(i,k,12) = qgrs(i,k,ntoz)
160 enddo
161 enddo
162 rtg_ozone_index = 12
163 else ! MG2
164 do k=1,levs
165 do i=1,im
166 vdftra(i,k,1) = qgrs(i,k,ntqv)
167 vdftra(i,k,2) = qgrs(i,k,ntcw)
168 vdftra(i,k,3) = qgrs(i,k,ntiw)
169 vdftra(i,k,4) = qgrs(i,k,ntrw)
170 vdftra(i,k,5) = qgrs(i,k,ntsw)
171 vdftra(i,k,6) = qgrs(i,k,ntlnc)
172 vdftra(i,k,7) = qgrs(i,k,ntinc)
173 vdftra(i,k,8) = qgrs(i,k,ntrnc)
174 vdftra(i,k,9) = qgrs(i,k,ntsnc)
175 vdftra(i,k,10) = qgrs(i,k,ntoz)
176 enddo
177 enddo
178 rtg_ozone_index = 10
179 endif
180 elseif (imp_physics == imp_physics_gfdl) then
181 ! GFDL MP
182 do k=1,levs
183 do i=1,im
184 vdftra(i,k,1) = qgrs(i,k,ntqv)
185 vdftra(i,k,2) = qgrs(i,k,ntcw)
186 vdftra(i,k,3) = qgrs(i,k,ntiw)
187 vdftra(i,k,4) = qgrs(i,k,ntrw)
188 vdftra(i,k,5) = qgrs(i,k,ntsw)
189 vdftra(i,k,6) = qgrs(i,k,ntgl)
190 vdftra(i,k,7) = qgrs(i,k,ntoz)
191 enddo
192 enddo
193 rtg_ozone_index = 7
194 elseif (imp_physics == imp_physics_zhao_carr) then
195! Zhao/Carr/Sundqvist
196 do k=1,levs
197 do i=1,im
198 vdftra(i,k,1) = qgrs(i,k,ntqv)
199 vdftra(i,k,2) = qgrs(i,k,ntcw)
200 vdftra(i,k,3) = qgrs(i,k,ntoz)
201 enddo
202 enddo
203 rtg_ozone_index = 3
204 elseif (imp_physics == imp_physics_nssl ) then
205 ! nssl
206 IF ( nssl_hail_on ) THEN
207 do k=1,levs
208 do i=1,im
209 vdftra(i,k,1) = qgrs(i,k,ntqv)
210 vdftra(i,k,2) = qgrs(i,k,ntcw)
211 vdftra(i,k,3) = qgrs(i,k,ntiw)
212 vdftra(i,k,4) = qgrs(i,k,ntrw)
213 vdftra(i,k,5) = qgrs(i,k,ntsw)
214 vdftra(i,k,6) = qgrs(i,k,ntgl)
215 vdftra(i,k,7) = qgrs(i,k,nthl)
216 vdftra(i,k,8) = qgrs(i,k,ntlnc)
217 vdftra(i,k,9) = qgrs(i,k,ntinc)
218 vdftra(i,k,10) = qgrs(i,k,ntrnc)
219 vdftra(i,k,11) = qgrs(i,k,ntsnc)
220 vdftra(i,k,12) = qgrs(i,k,ntgnc)
221 vdftra(i,k,13) = qgrs(i,k,nthnc)
222 vdftra(i,k,14) = qgrs(i,k,ntgv)
223 vdftra(i,k,15) = qgrs(i,k,nthv)
224 vdftra(i,k,16) = qgrs(i,k,ntoz)
225 n = 16
226 IF ( nssl_ccn_on ) THEN
227 vdftra(i,k,n+1) = qgrs(i,k,ntccn)
228 n = n+1
229 ENDIF
230 IF ( nssl_3moment ) THEN
231 vdftra(i,k,n+1) = qgrs(i,k,ntrz)
232 vdftra(i,k,n+2) = qgrs(i,k,ntgz)
233 vdftra(i,k,n+3) = qgrs(i,k,nthz)
234 n = n+3
235 ENDIF
236 enddo
237 enddo
238
239 ELSE
240 ! no hail
241 do k=1,levs
242 do i=1,im
243 vdftra(i,k,1) = qgrs(i,k,ntqv)
244 vdftra(i,k,2) = qgrs(i,k,ntcw)
245 vdftra(i,k,3) = qgrs(i,k,ntiw)
246 vdftra(i,k,4) = qgrs(i,k,ntrw)
247 vdftra(i,k,5) = qgrs(i,k,ntsw)
248 vdftra(i,k,6) = qgrs(i,k,ntgl)
249 vdftra(i,k,7) = qgrs(i,k,ntlnc)
250 vdftra(i,k,8) = qgrs(i,k,ntinc)
251 vdftra(i,k,9) = qgrs(i,k,ntrnc)
252 vdftra(i,k,10) = qgrs(i,k,ntsnc)
253 vdftra(i,k,11) = qgrs(i,k,ntgnc)
254 vdftra(i,k,12) = qgrs(i,k,ntgv)
255 vdftra(i,k,13) = qgrs(i,k,ntoz)
256 n = 13
257 IF ( nssl_ccn_on ) THEN
258 vdftra(i,k,n+1) = qgrs(i,k,ntccn)
259 n = n+1
260 ENDIF
261 IF ( nssl_3moment ) THEN
262 vdftra(i,k,n+1) = qgrs(i,k,ntrz)
263 vdftra(i,k,n+2) = qgrs(i,k,ntgz)
264 n = n+2
265 ENDIF
266 enddo
267 enddo
268
269 ENDIF
270
271
272 endif
273!
274 if (trans_aero) then
275 call set_aerosol_tracer_index(imp_physics, imp_physics_wsm6, &
276 imp_physics_thompson, ltaerosol,mraerosol, &
277 imp_physics_mg, ntgl, imp_physics_gfdl, &
278 imp_physics_zhao_carr, imp_physics_nssl,&
279 nssl_hail_on, nssl_ccn_on, kk, &
280 errmsg, errflg)
281 if (errflg /= 0) return
282 !
283 k1 = kk
284 do n=ntchs,ntchm+ntchs-1
285 k1 = k1 + 1
286 do k=1,levs
287 do i=1,im
288 vdftra(i,k,k1) = qgrs(i,k,n)
289 enddo
290 enddo
291 enddo
292 endif
293!
294 if (ntke>0) then
295 do k=1,levs
296 do i=1,im
297 vdftra(i,k,ntkev) = qgrs(i,k,ntke)
298 enddo
299 enddo
300 endif
301!
302 endif
303
304 if(ldiag3d .and. lssav .and. flag_for_pbl_generic_tend) then
305 do k=1,levs
306 do i=1,im
307 save_t(i,k) = tgrs(i,k)
308 save_u(i,k) = ugrs(i,k)
309 save_v(i,k) = vgrs(i,k)
310 enddo
311 enddo
312 if(qdiag3d) then
313 do k=1,levs
314 do i=1,im
315 save_q(i,k,ntqv) = qgrs(i,k,ntqv)
316 save_q(i,k,ntoz) = qgrs(i,k,ntoz)
317 enddo
318 enddo
319 if(ntke>0) then
320 do k=1,levs
321 do i=1,im
322 save_q(i,k,ntke) = qgrs(i,k,ntke)
323 enddo
324 enddo
325 endif
326 endif
327 endif
328
329 end subroutine gfs_pbl_generic_pre_run
330
331 end module gfs_pbl_generic_pre