CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_PBL_generic_post.F90
1
3
5
6 contains
7
11 subroutine gfs_pbl_generic_post_run (im, levs, nvdiff, ntrac, &
12 ntqv, ntcw, ntiw, ntrw, ntsw, ntlnc, ntinc, ntrnc, ntsnc, ntgnc, ntwa, ntia, ntgl, ntoz, ntke, ntkev,nqrimef, &
13 trans_aero, ntchs, ntchm, ntccn, nthl, nthnc, ntgv, nthv, ntrz, ntgz, nthz, &
14 imp_physics, imp_physics_gfdl, imp_physics_thompson, imp_physics_wsm6, imp_physics_zhao_carr, imp_physics_mg, &
15 imp_physics_fer_hires, imp_physics_nssl, nssl_ccn_on, ltaerosol, mraerosol, nssl_hail_on, nssl_3moment, &
16 cplflx, cplaqm, cplchm, lssav, flag_for_pbl_generic_tend, ldiag3d, lsidea, hybedmf, do_shoc, satmedmf, &
17 shinhong, do_ysu, dvdftra, dusfc1, dvsfc1, dtsfc1, dqsfc1, dtf, dudt, dvdt, dtdt, htrsw, htrlw, xmu, &
18 dqdt, dusfc_cpl, dvsfc_cpl, dtsfc_cpl, dtend, dtidx, index_of_temperature, index_of_x_wind, index_of_y_wind, &
19 index_of_process_pbl, dqsfc_cpl, dusfci_cpl, dvsfci_cpl, dtsfci_cpl, dqsfci_cpl, dusfc_diag, dvsfc_diag, dtsfc_diag, &
20 dqsfc_diag, dusfci_diag, dvsfci_diag, dtsfci_diag, dqsfci_diag, &
21 rd, cp, fvirt, hvap, t1, q1, prsl, hflx, ushfsfci, oceanfrac, kdt, dusfc_cice, dvsfc_cice, &
22 dtsfc_cice, dqsfc_cice, use_med_flux, dtsfc_med, dqsfc_med, dusfc_med, dvsfc_med, wet, dry, icy, wind, stress_wat, &
23 hflx_wat, evap_wat, ugrs1, vgrs1, hffac, ugrs, vgrs, tgrs, qgrs, save_u, save_v, save_t, save_q, huge, errmsg, errflg)
24
25 use machine, only : kind_phys
26 use gfs_pbl_generic_common, only : set_aerosol_tracer_index
27
28 implicit none
29
30 integer, parameter :: kp = kind_phys
31 integer, intent(in) :: im, levs, nvdiff, ntrac, ntchs, ntchm, kdt
32 integer, intent(in) :: ntqv, ntcw, ntiw, ntrw, ntsw, ntlnc, ntinc, ntrnc, ntsnc, ntgnc, ntwa, ntia, ntgl, ntoz, ntke, ntkev, nqrimef
33 integer, intent(in) :: ntccn, nthl, nthnc, ntgv, nthv, ntrz, ntgz, nthz
34 logical, intent(in) :: trans_aero
35 integer, intent(in) :: imp_physics, imp_physics_gfdl, imp_physics_thompson, imp_physics_wsm6
36 integer, intent(in) :: imp_physics_zhao_carr, imp_physics_mg, imp_physics_fer_hires
37 integer, intent(in) :: imp_physics_nssl
38 logical, intent(in) :: nssl_ccn_on, nssl_hail_on, nssl_3moment
39 logical, intent(in) :: ltaerosol, cplflx, cplaqm, cplchm, lssav, ldiag3d, lsidea, use_med_flux, mraerosol
40 logical, intent(in) :: hybedmf, do_shoc, satmedmf, shinhong, do_ysu
41
42 logical, intent(in) :: flag_for_pbl_generic_tend
43 real(kind=kind_phys), dimension(:,:), intent(in) :: save_u, save_v, save_t
44 real(kind=kind_phys), dimension(:,:, :), intent(in) :: save_q
45
46 real(kind=kind_phys), intent(in) :: dtf
47 real(kind=kind_phys), intent(in) :: rd, cp, fvirt, hvap, huge
48 real(kind=kind_phys), dimension(:), intent(in) :: t1, q1, hflx, oceanfrac
49 real(kind=kind_phys), dimension(:,:), intent(in) :: prsl
50 real(kind=kind_phys), dimension(:), intent(in), optional :: dusfc_cice, dvsfc_cice, dtsfc_cice, dqsfc_cice, &
51 dtsfc_med, dqsfc_med, dusfc_med, dvsfc_med
52 real(kind=kind_phys), dimension(:), intent(in) :: wind, stress_wat, hflx_wat, evap_wat, ugrs1, vgrs1
53
54 real(kind=kind_phys), dimension(:,:, :), intent(in) :: qgrs
55 real(kind=kind_phys), dimension(:,:), intent(in) :: ugrs, vgrs, tgrs
56
57 real(kind=kind_phys), dimension(:,:, :), intent(in) :: dvdftra
58 real(kind=kind_phys), dimension(:), intent(in) :: dusfc1, dvsfc1, dtsfc1, dqsfc1, xmu
59 real(kind=kind_phys), dimension(:,:), intent(in) :: dudt, dvdt, dtdt, htrsw, htrlw
60
61 real(kind=kind_phys), dimension(:,:, :), intent(inout) :: dqdt
62
63 ! The following arrays may not be allocated, depending on certain flags (cplflx, ...).
64 ! Since Intel 15 crashes when passing unallocated arrays to arrays defined with explicit shape,
65 ! use assumed-shape arrays. Note that Intel 18 and GNU 6.2.0-8.1.0 tolerate explicit-shape arrays
66 ! as long as these do not get used when not allocated.
67 real(kind=kind_phys), dimension(:), intent(inout), optional :: dusfc_cpl, dvsfc_cpl, &
68 dtsfc_cpl, dqsfc_cpl, dusfci_cpl, dvsfci_cpl, dtsfci_cpl, dqsfci_cpl
69 real(kind=kind_phys), dimension(:), intent(inout) :: dusfc_diag, dvsfc_diag, &
70 dtsfc_diag, dqsfc_diag, dusfci_diag, dvsfci_diag, dtsfci_diag, dqsfci_diag
71 real(kind=kind_phys), intent(inout), optional :: dtend(:,:,:)
72 integer, intent(in) :: dtidx(:,:)
73 integer, intent(in) :: index_of_temperature, index_of_x_wind, index_of_y_wind, index_of_process_pbl
74
75 logical, dimension(:),intent(in) :: wet, dry, icy
76 real(kind=kind_phys), dimension(:), intent(out), optional :: ushfsfci
77
78 ! From canopy heat storage - reduction factors in latent/sensible heat flux due to surface roughness
79 real(kind=kind_phys), dimension(:), intent(in) :: hffac
80
81 character(len=*), intent(out) :: errmsg
82 integer, intent(out) :: errflg
83
84 real(kind=kind_phys), parameter :: zero = 0.0_kp, one = 1.0_kp
85! real(kind=kind_phys), parameter :: huge = 9.9692099683868690E36 ! NetCDF float FillValue, same as in GFS_typedefs.F90
86 real(kind=kind_phys), parameter :: qmin = 1.0e-8_kp
87 integer :: i, k, kk, k1, n
88 real(kind=kind_phys) :: tem, rho
89 integer :: idtend
90
91 ! Initialize CCPP error handling variables
92 errmsg = ''
93 errflg = 0
94!GJF: dvdftra is only used if nvdiff != ntrac or (nvdiff == ntrac .and. )
95 if (nvdiff == ntrac .and. (hybedmf .or. do_shoc .or. satmedmf)) then
96 dqdt = dvdftra
97 elseif (nvdiff /= ntrac .and. .not. shinhong .and. .not. do_ysu) then
98!
99 if (ntke>0) then
100 do k=1,levs
101 do i=1,im
102 dqdt(i,k,ntke) = dvdftra(i,k,ntkev)
103 enddo
104 enddo
105 endif
106!
107 if (trans_aero) then
108 ! Set kk if chemistry-aerosol tracers are diffused
109 call set_aerosol_tracer_index(imp_physics, imp_physics_wsm6, &
110 imp_physics_thompson, ltaerosol,mraerosol, &
111 imp_physics_mg, ntgl, imp_physics_gfdl, &
112 imp_physics_zhao_carr, imp_physics_nssl,&
113 nssl_hail_on, nssl_ccn_on, kk, &
114 errmsg, errflg)
115 if (errflg /= 0) return
116 !
117 k1 = kk
118 do n=ntchs,ntchm+ntchs-1
119 k1 = k1 + 1
120 do k=1,levs
121 do i=1,im
122 dqdt(i,k,n) = dvdftra(i,k,k1)
123 enddo
124 enddo
125 enddo
126 endif
127!
128 if (imp_physics == imp_physics_wsm6) then
129 ! WSM6
130 do k=1,levs
131 do i=1,im
132 dqdt(i,k,ntqv) = dvdftra(i,k,1)
133 dqdt(i,k,ntcw) = dvdftra(i,k,2)
134 dqdt(i,k,ntiw) = dvdftra(i,k,3)
135 dqdt(i,k,ntoz) = dvdftra(i,k,4)
136 enddo
137 enddo
138
139 elseif (imp_physics == imp_physics_fer_hires) then
140 ! Ferrier-Aligo
141 do k=1,levs
142 do i=1,im
143 dqdt(i,k,ntqv) = dvdftra(i,k,1)
144 dqdt(i,k,ntcw) = dvdftra(i,k,2)
145 dqdt(i,k,ntiw) = dvdftra(i,k,3)
146 dqdt(i,k,ntrw) = dvdftra(i,k,4)
147 dqdt(i,k,nqrimef) = dvdftra(i,k,5)
148 dqdt(i,k,ntoz) = dvdftra(i,k,6)
149 enddo
150 enddo
151
152 elseif (imp_physics == imp_physics_thompson) then
153 ! Thompson
154 if(ltaerosol) then
155 do k=1,levs
156 do i=1,im
157 dqdt(i,k,ntqv) = dvdftra(i,k,1)
158 dqdt(i,k,ntcw) = dvdftra(i,k,2)
159 dqdt(i,k,ntiw) = dvdftra(i,k,3)
160 dqdt(i,k,ntrw) = dvdftra(i,k,4)
161 dqdt(i,k,ntsw) = dvdftra(i,k,5)
162 dqdt(i,k,ntgl) = dvdftra(i,k,6)
163 dqdt(i,k,ntlnc) = dvdftra(i,k,7)
164 dqdt(i,k,ntinc) = dvdftra(i,k,8)
165 dqdt(i,k,ntrnc) = dvdftra(i,k,9)
166 dqdt(i,k,ntoz) = dvdftra(i,k,10)
167 dqdt(i,k,ntwa) = dvdftra(i,k,11)
168 dqdt(i,k,ntia) = dvdftra(i,k,12)
169 enddo
170 enddo
171 else if(mraerosol) then
172 do k=1,levs
173 do i=1,im
174 dqdt(i,k,ntqv) = dvdftra(i,k,1)
175 dqdt(i,k,ntcw) = dvdftra(i,k,2)
176 dqdt(i,k,ntiw) = dvdftra(i,k,3)
177 dqdt(i,k,ntrw) = dvdftra(i,k,4)
178 dqdt(i,k,ntsw) = dvdftra(i,k,5)
179 dqdt(i,k,ntgl) = dvdftra(i,k,6)
180 dqdt(i,k,ntlnc) = dvdftra(i,k,7)
181 dqdt(i,k,ntinc) = dvdftra(i,k,8)
182 dqdt(i,k,ntrnc) = dvdftra(i,k,9)
183 dqdt(i,k,ntoz) = dvdftra(i,k,10)
184 enddo
185 enddo
186 else
187 do k=1,levs
188 do i=1,im
189 dqdt(i,k,ntqv) = dvdftra(i,k,1)
190 dqdt(i,k,ntcw) = dvdftra(i,k,2)
191 dqdt(i,k,ntiw) = dvdftra(i,k,3)
192 dqdt(i,k,ntrw) = dvdftra(i,k,4)
193 dqdt(i,k,ntsw) = dvdftra(i,k,5)
194 dqdt(i,k,ntgl) = dvdftra(i,k,6)
195 dqdt(i,k,ntinc) = dvdftra(i,k,7)
196 dqdt(i,k,ntrnc) = dvdftra(i,k,8)
197 dqdt(i,k,ntoz) = dvdftra(i,k,9)
198 enddo
199 enddo
200 endif
201 elseif (imp_physics == imp_physics_mg) then ! MG3/2
202 if (ntgl > 0) then ! MG
203 do k=1,levs
204 do i=1,im
205 dqdt(i,k,1) = dvdftra(i,k,1)
206 dqdt(i,k,ntcw) = dvdftra(i,k,2)
207 dqdt(i,k,ntiw) = dvdftra(i,k,3)
208 dqdt(i,k,ntrw) = dvdftra(i,k,4)
209 dqdt(i,k,ntsw) = dvdftra(i,k,5)
210 dqdt(i,k,ntgl) = dvdftra(i,k,6)
211 dqdt(i,k,ntlnc) = dvdftra(i,k,7)
212 dqdt(i,k,ntinc) = dvdftra(i,k,8)
213 dqdt(i,k,ntrnc) = dvdftra(i,k,9)
214 dqdt(i,k,ntsnc) = dvdftra(i,k,10)
215 dqdt(i,k,ntgnc) = dvdftra(i,k,11)
216 dqdt(i,k,ntoz) = dvdftra(i,k,12)
217 enddo
218 enddo
219 else ! MG2
220 do k=1,levs
221 do i=1,im
222 dqdt(i,k,1) = dvdftra(i,k,1)
223 dqdt(i,k,ntcw) = dvdftra(i,k,2)
224 dqdt(i,k,ntiw) = dvdftra(i,k,3)
225 dqdt(i,k,ntrw) = dvdftra(i,k,4)
226 dqdt(i,k,ntsw) = dvdftra(i,k,5)
227 dqdt(i,k,ntlnc) = dvdftra(i,k,6)
228 dqdt(i,k,ntinc) = dvdftra(i,k,7)
229 dqdt(i,k,ntrnc) = dvdftra(i,k,8)
230 dqdt(i,k,ntsnc) = dvdftra(i,k,9)
231 dqdt(i,k,ntoz) = dvdftra(i,k,10)
232 enddo
233 enddo
234 endif
235 elseif (imp_physics == imp_physics_gfdl) then ! GFDL MP
236 do k=1,levs
237 do i=1,im
238 dqdt(i,k,ntqv) = dvdftra(i,k,1)
239 dqdt(i,k,ntcw) = dvdftra(i,k,2)
240 dqdt(i,k,ntiw) = dvdftra(i,k,3)
241 dqdt(i,k,ntrw) = dvdftra(i,k,4)
242 dqdt(i,k,ntsw) = dvdftra(i,k,5)
243 dqdt(i,k,ntgl) = dvdftra(i,k,6)
244 dqdt(i,k,ntoz) = dvdftra(i,k,7)
245 enddo
246 enddo
247 elseif (imp_physics == imp_physics_zhao_carr) then
248 do k=1,levs
249 do i=1,im
250 dqdt(i,k,1) = dvdftra(i,k,1)
251 dqdt(i,k,ntcw) = dvdftra(i,k,2)
252 dqdt(i,k,ntoz) = dvdftra(i,k,3)
253 enddo
254 enddo
255 elseif (imp_physics == imp_physics_nssl ) then
256 ! nssl
257 IF ( nssl_hail_on ) THEN
258 do k=1,levs
259 do i=1,im
260 dqdt(i,k,ntqv) = dvdftra(i,k,1)
261 dqdt(i,k,ntcw) = dvdftra(i,k,2)
262 dqdt(i,k,ntiw) = dvdftra(i,k,3)
263 dqdt(i,k,ntrw) = dvdftra(i,k,4)
264 dqdt(i,k,ntsw) = dvdftra(i,k,5)
265 dqdt(i,k,ntgl) = dvdftra(i,k,6)
266 dqdt(i,k,nthl) = dvdftra(i,k,7)
267 dqdt(i,k,ntlnc) = dvdftra(i,k,8)
268 dqdt(i,k,ntinc) = dvdftra(i,k,9)
269 dqdt(i,k,ntrnc) = dvdftra(i,k,10)
270 dqdt(i,k,ntsnc) = dvdftra(i,k,11)
271 dqdt(i,k,ntgnc) = dvdftra(i,k,12)
272 dqdt(i,k,nthnc) = dvdftra(i,k,13)
273 dqdt(i,k,ntgv) = dvdftra(i,k,14)
274 dqdt(i,k,nthv) = dvdftra(i,k,15)
275 dqdt(i,k,ntoz) = dvdftra(i,k,16)
276 n = 16
277 IF ( nssl_ccn_on ) THEN
278 dqdt(i,k,ntccn) = dvdftra(i,k,n+1)
279 n = n+1
280 ENDIF
281 IF ( nssl_3moment ) THEN
282 dqdt(i,k,ntrz) = dvdftra(i,k,n+1)
283 dqdt(i,k,ntgz) = dvdftra(i,k,n+2)
284 dqdt(i,k,nthz) = dvdftra(i,k,n+3)
285 n = n+3
286 ENDIF
287 enddo
288 enddo
289
290 ELSE
291
292 do k=1,levs
293 do i=1,im
294 dqdt(i,k,ntqv) = dvdftra(i,k,1)
295 dqdt(i,k,ntcw) = dvdftra(i,k,2)
296 dqdt(i,k,ntiw) = dvdftra(i,k,3)
297 dqdt(i,k,ntrw) = dvdftra(i,k,4)
298 dqdt(i,k,ntsw) = dvdftra(i,k,5)
299 dqdt(i,k,ntgl) = dvdftra(i,k,6)
300 dqdt(i,k,ntlnc) = dvdftra(i,k,7)
301 dqdt(i,k,ntinc) = dvdftra(i,k,8)
302 dqdt(i,k,ntrnc) = dvdftra(i,k,9)
303 dqdt(i,k,ntsnc) = dvdftra(i,k,10)
304 dqdt(i,k,ntgnc) = dvdftra(i,k,11)
305 dqdt(i,k,ntgv) = dvdftra(i,k,12)
306 dqdt(i,k,ntoz) = dvdftra(i,k,13)
307 n = 13
308 IF ( nssl_ccn_on ) THEN
309 dqdt(i,k,ntccn) = dvdftra(i,k,n+1)
310 n = n+1
311 ENDIF
312 IF ( nssl_3moment ) THEN
313 dqdt(i,k,ntrz) = dvdftra(i,k,n+1)
314 dqdt(i,k,ntgz) = dvdftra(i,k,n+2)
315 n = n+2
316 ENDIF
317 enddo
318 enddo
319
320 ENDIF
321 endif
322
323 endif ! nvdiff == ntrac
324
325! --- ... coupling insertion
326
327 if (cplflx) then
328 do i=1,im
329 if (oceanfrac(i) > zero) then ! Ocean only, NO LAKES
330 if ( .not. wet(i)) then ! no open water
331 if (kdt > 1) then !use results from CICE
332 dusfci_cpl(i) = dusfc_cice(i)
333 dvsfci_cpl(i) = dvsfc_cice(i)
334 dtsfci_cpl(i) = dtsfc_cice(i)
335 dqsfci_cpl(i) = dqsfc_cice(i)
336 else !use PBL fluxes when CICE fluxes is unavailable
337 dusfci_cpl(i) = dusfc1(i)
338 dvsfci_cpl(i) = dvsfc1(i)
339 dtsfci_cpl(i) = dtsfc1(i)*hffac(i)
340 dqsfci_cpl(i) = dqsfc1(i)
341 end if
342 elseif (icy(i) .or. dry(i)) then ! use stress_ocean from sfc_diff for opw component at mixed point
343 rho = prsl(i,1) / (rd*t1(i)*(one+fvirt*max(q1(i), qmin)))
344 if (wind(i) > zero) then
345 tem = - rho * stress_wat(i) / wind(i)
346 dusfci_cpl(i) = tem * ugrs1(i) ! U-momentum flux
347 dvsfci_cpl(i) = tem * vgrs1(i) ! V-momentum flux
348 else
349 dusfci_cpl(i) = zero
350 dvsfci_cpl(i) = zero
351 endif
352 dtsfci_cpl(i) = cp * rho * hflx_wat(i) ! sensible heat flux over open ocean
353 dqsfci_cpl(i) = hvap * rho * evap_wat(i) ! latent heat flux over open ocean
354 else ! 100% open ocean
355 if (use_med_flux .and. kdt > 1) then ! use results from CMEPS mediator
356 dusfci_cpl(i) = dusfc_med(i)
357 dvsfci_cpl(i) = dvsfc_med(i)
358 dtsfci_cpl(i) = dtsfc_med(i)
359 dqsfci_cpl(i) = dqsfc_med(i)
360 else ! use results from PBL scheme
361 dusfci_cpl(i) = dusfc1(i)
362 dvsfci_cpl(i) = dvsfc1(i)
363 dtsfci_cpl(i) = dtsfc1(i)*hffac(i)
364 dqsfci_cpl(i) = dqsfc1(i)
365 end if
366 endif
367!
368 dusfc_cpl(i) = dusfc_cpl(i) + dusfci_cpl(i) * dtf
369 dvsfc_cpl(i) = dvsfc_cpl(i) + dvsfci_cpl(i) * dtf
370 dtsfc_cpl(i) = dtsfc_cpl(i) + dtsfci_cpl(i) * dtf
371 dqsfc_cpl(i) = dqsfc_cpl(i) + dqsfci_cpl(i) * dtf
372!
373 else
374 dusfc_cpl(i) = huge
375 dvsfc_cpl(i) = huge
376 dtsfc_cpl(i) = huge
377 dqsfc_cpl(i) = huge
378!!
379 endif ! Ocean only, NO LAKES
380 enddo
381 endif
382
383 if (cplchm) then
384 if (cplflx) then
385 do i = 1, im
386 if (oceanfrac(i) > zero) then
387 ushfsfci(i) = dtsfci_cpl(i)
388 else
389 rho = prsl(i,1) / (rd*t1(i)*(one+fvirt*max(q1(i), qmin)))
390 ushfsfci(i) = cp * rho * hflx(i)
391 end if
392 end do
393 else
394 do i = 1, im
395 rho = prsl(i,1) / (rd*t1(i)*(one+fvirt*max(q1(i), qmin)))
396 ushfsfci(i) = cp * rho * hflx(i)
397 end do
398 end if
399 end if
400
401 if (cplaqm) then
402 do i = 1, im
403 if (oceanfrac(i) > zero) then
404 if (.not.cplflx) then
405 dtsfci_cpl(i) = dtsfc1(i)*hffac(i)
406 dqsfci_cpl(i) = dqsfc1(i)
407 end if
408 else ! heat fluxes are required over land
409 dtsfci_cpl(i) = dtsfc1(i)*hffac(i)
410 dqsfci_cpl(i) = dqsfc1(i)
411 end if
412 end do
413 end if
414
415!-------------------------------------------------------lssav if loop ----------
416 if (lssav) then
417 do i=1,im
418 dusfc_diag(i) = dusfc_diag(i) + dusfc1(i) * dtf
419 dvsfc_diag(i) = dvsfc_diag(i) + dvsfc1(i) * dtf
420 dusfci_diag(i) = dusfc1(i)
421 dvsfci_diag(i) = dvsfc1(i)
422 dtsfci_diag(i) = dtsfc1(i)*hffac(i)
423 dqsfci_diag(i) = dqsfc1(i)
424 dtsfc_diag(i) = dtsfc_diag(i) + dtsfci_diag(i) * dtf
425 dqsfc_diag(i) = dqsfc_diag(i) + dqsfci_diag(i) * dtf
426 enddo
427
428 if (ldiag3d .and. flag_for_pbl_generic_tend) then
429 if (lsidea) then
430 idtend = dtidx(index_of_temperature, index_of_process_pbl)
431 if(idtend>=1) then
432 dtend(1:im,1:levs,idtend) = dtend(1:im,1:levs,idtend) + dtdt(1:im,1:levs)*dtf
433 endif
434 else
435 idtend = dtidx(index_of_temperature, index_of_process_pbl)
436 if(idtend>=1) then
437 dtend(1:im,1:levs,idtend) = dtend(1:im,1:levs,idtend) + (tgrs(1:im,1:levs) - save_t(1:im,1:levs))
438 endif
439 endif
440 idtend = dtidx(index_of_x_wind, index_of_process_pbl)
441 if(idtend>=1) then
442 dtend(1:im,1:levs,idtend) = dtend(1:im,1:levs,idtend) + (ugrs(1:im,1:levs) - save_u(1:im,1:levs))
443 endif
444 idtend = dtidx(index_of_y_wind, index_of_process_pbl)
445 if(idtend>=1) then
446 dtend(1:im,1:levs,idtend) = dtend(1:im,1:levs,idtend) + (vgrs(1:im,1:levs) - save_v(1:im,1:levs))
447 endif
448 idtend = dtidx(100+ntqv, index_of_process_pbl)
449 if(idtend>=1) then
450 dtend(1:im,1:levs,idtend) = dtend(1:im,1:levs,idtend) + qgrs(1:im,1:levs,ntqv) - save_q(1:im,1:levs,ntqv)
451 endif
452 idtend = dtidx(100+ntoz, index_of_process_pbl)
453 if(idtend>=1) then
454 dtend(1:im,1:levs,idtend) = dtend(1:im,1:levs,idtend) + qgrs(1:im,1:levs,ntoz) - save_q(1:im,1:levs,ntoz)
455 endif
456 idtend = dtidx(100+ntke, index_of_process_pbl)
457 if(idtend>=1) then
458 dtend(1:im,1:levs,idtend) = dtend(1:im,1:levs,idtend) + (qgrs(1:im,1:levs,ntke) - save_q(1:im,1:levs,ntke))
459 endif
460 endif
461
462 endif ! end if_lssav
463
464 end subroutine gfs_pbl_generic_post_run
465
466 end module gfs_pbl_generic_post