CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
mynnedmf_wrapper.F90
1
4
9
10 contains
11
15 subroutine mynnedmf_wrapper_init ( &
16 & con_cp, con_grav, con_rd, con_rv, &
17 & con_cpv, con_cliq, con_cice, con_rcp, &
18 & con_XLV, con_XLF, con_p608, con_ep2, &
19 & con_karman, con_t0c, &
20 & do_mynnedmf, &
21 & errmsg, errflg )
22
23 use machine, only : kind_phys
25
26 implicit none
27
28 logical, intent(in) :: do_mynnedmf
29 character(len=*),intent(out):: errmsg
30 integer, intent(out) :: errflg
31
32 real(kind_phys),intent(in) :: con_xlv
33 real(kind_phys),intent(in) :: con_xlf
34 real(kind_phys),intent(in) :: con_rv
35 real(kind_phys),intent(in) :: con_rd
36 real(kind_phys),intent(in) :: con_ep2
37 real(kind_phys),intent(in) :: con_grav
38 real(kind_phys),intent(in) :: con_cp
39 real(kind_phys),intent(in) :: con_cpv
40 real(kind_phys),intent(in) :: con_rcp
41 real(kind_phys),intent(in) :: con_p608
42 real(kind_phys),intent(in) :: con_cliq
43 real(kind_phys),intent(in) :: con_cice
44 real(kind_phys),intent(in) :: con_karman
45 real(kind_phys),intent(in) :: con_t0c
46
47 ! Initialize CCPP error handling variables
48 errmsg = ''
49 errflg = 0
50
51 xlv = con_xlv
52 xlf = con_xlf
53 r_v = con_rv
54 r_d = con_rd
55 ep_2 = con_ep2
56 grav = con_grav
57 cp = con_cp
58 cpv = con_cpv
59 rcp = con_rcp
60 p608 = con_p608
61 cliq = con_cliq
62 cice = con_cice
63 karman = con_karman
64 t0c = con_t0c
65
66 xls = xlv+xlf != 2.85E6 (J/kg) sublimation
67 rvovrd = r_v/r_d != 1.608
68 ep_3 = 1.-ep_2 != 0.378
69 gtr = grav/tref
70 rk = cp/r_d
71 tv0 = p608*tref
72 tv1 = (1.+p608)*tref
73 xlscp = (xlv+xlf)/cp
74 xlvcp = xlv/cp
75 g_inv = 1./grav
76
77 ! Consistency checks
78 if (.not. do_mynnedmf) then
79 errmsg = 'Logic error: do_mynnedmf = .false.'
80 errflg = 1
81 return
82 end if
83
84 end subroutine mynnedmf_wrapper_init
85
91SUBROUTINE mynnedmf_wrapper_run( &
92 & im,levs, &
93 & flag_init,flag_restart, &
94 & lssav, ldiag3d, qdiag3d, &
95 & lsidea, cplflx, &
96 & delt,dtf,dx,zorl, &
97 & phii,u,v,omega,t3d, &
98 & qgrs_water_vapor, &
99 & qgrs_liquid_cloud, &
100 & qgrs_ice, &
101 & qgrs_snow, &
102 & qgrs_cloud_droplet_num_conc, &
103 & qgrs_cloud_ice_num_conc, &
104 & qgrs_ozone, &
105 & qgrs_water_aer_num_conc, &
106 & qgrs_ice_aer_num_conc, &
107 & qgrs_cccn, &
108 & prsl,prsi,exner, &
109 & slmsk,tsurf,qsfc,ps, &
110 & ust,ch,hflx,qflx,wspd,rb, &
111 & dtsfc1,dqsfc1, &
112 & dusfc1,dvsfc1, &
113 & dusfci_diag,dvsfci_diag, &
114 & dtsfci_diag,dqsfci_diag, &
115 & dusfc_diag,dvsfc_diag, &
116 & dtsfc_diag,dqsfc_diag, &
117 & dusfc_cice,dvsfc_cice, &
118 & dtsfc_cice,dqsfc_cice, &
119 & hflx_wat,qflx_wat,stress_wat, &
120 & oceanfrac,fice,wet,icy,dry, &
121 & dusfci_cpl,dvsfci_cpl, &
122 & dtsfci_cpl,dqsfci_cpl, &
123 & dusfc_cpl,dvsfc_cpl, &
124 & dtsfc_cpl,dqsfc_cpl, &
125 & recmol, &
126 & qke,qke_adv,Tsq,Qsq,Cov, &
127 & el_pbl,sh3d,sm3d,exch_h,exch_m, &
128 & dqke,qwt,qshear,qbuoy,qdiss, &
129 & Pblh,kpbl, &
130 & qc_bl,qi_bl,cldfra_bl, &
131 & edmf_a,edmf_w,edmf_qt, &
132 & edmf_thl,edmf_ent,edmf_qc, &
133 & sub_thl,sub_sqv,det_thl,det_sqv,&
134 & maxwidth,maxMF,ztop_plume, &
135 & ktop_plume, &
136 & dudt, dvdt, dtdt, &
137 & dqdt_water_vapor, dqdt_liquid_cloud, & ! <=== ntqv, ntcw
138 & dqdt_ice, dqdt_snow, & ! <=== ntiw, ntsw
139 & dqdt_ozone, & ! <=== ntoz
140 & dqdt_cloud_droplet_num_conc, dqdt_ice_num_conc, & ! <=== ntlnc, ntinc
141 & dqdt_water_aer_num_conc, dqdt_ice_aer_num_conc,& ! <=== ntwa, ntia
142 & dqdt_cccn, & ! <=== ntccn
143 & flag_for_pbl_generic_tend, &
144 & dtend, dtidx, index_of_temperature, &
145 & index_of_x_wind, index_of_y_wind, ntke, &
146 & ntqv, ntcw, ntiw, ntsw, &
147 & ntoz, ntlnc, ntinc, ntwa, ntia, &
148 & index_of_process_pbl, htrsw, htrlw, xmu, &
149 & tke_budget, bl_mynn_tkeadvect, &
150 & bl_mynn_cloudpdf, bl_mynn_mixlength, &
151 & bl_mynn_edmf, &
152 & bl_mynn_edmf_mom, bl_mynn_edmf_tke, &
153 & bl_mynn_cloudmix, bl_mynn_mixqt, &
154 & bl_mynn_output, bl_mynn_closure, &
155 & icloud_bl, do_mynnsfclay, &
156 & imp_physics, imp_physics_gfdl, &
157 & imp_physics_thompson, imp_physics_wsm6, &
158 & imp_physics_fa, &
159 & chem3d, frp, mix_chem, rrfs_sd, enh_mix, &
160 & nchem, ndvel, vdep, smoke_dbg, &
161 & imp_physics_nssl, nssl_ccn_on, &
162 & ltaerosol, mraerosol, spp_wts_pbl, spp_pbl, &
163 & lprnt, huge, errmsg, errflg )
164
165! should be moved to inside the mynn:
166 use machine, only: kind_phys
167 use bl_mynn_common, only: cp, r_d, grav, g_inv, zero, &
168 xlv, xlvcp, xlscp, p608
170
171!-------------------------------------------------------------------
172 implicit none
173!-------------------------------------------------------------------
174
175 real(kind_phys), intent(in) :: huge
176 character(len=*), intent(out) :: errmsg
177 integer, intent(out) :: errflg
178
179 logical, intent(in) :: lssav, ldiag3d, lsidea, qdiag3d
180 logical, intent(in) :: cplflx
181
182 !smoke/chem
183 integer, intent(in) :: nchem, ndvel
184 integer, parameter :: kdvel=1
185 logical, intent(in) :: smoke_dbg
186
187! NAMELIST OPTIONS (INPUT):
188 logical, intent(in) :: &
189 & bl_mynn_tkeadvect, &
190 & ltaerosol, &
191 & mraerosol, &
192 & lprnt, &
193 & do_mynnsfclay, &
194 & flag_for_pbl_generic_tend, &
195 & nssl_ccn_on
196 integer, intent(in) :: &
197 & bl_mynn_cloudpdf, &
198 & bl_mynn_mixlength, &
199 & icloud_bl, &
200 & bl_mynn_edmf, &
201 & bl_mynn_edmf_mom, &
202 & bl_mynn_edmf_tke, &
203 & bl_mynn_cloudmix, &
204 & bl_mynn_mixqt, &
205 & bl_mynn_output, &
206 & imp_physics, imp_physics_wsm6, &
207 & imp_physics_thompson, imp_physics_gfdl, &
208 & imp_physics_nssl, imp_physics_fa, &
209 & spp_pbl, &
210 & tke_budget
211 real(kind_phys), intent(in) :: &
212 & bl_mynn_closure
213
214!TENDENCY DIAGNOSTICS
215 real(kind_phys), intent(inout), optional :: dtend(:,:,:)
216 integer, intent(in) :: dtidx(:,:)
217 integer, intent(in) :: index_of_temperature, index_of_x_wind
218 integer, intent(in) :: index_of_y_wind, index_of_process_pbl
219 integer, intent(in) :: ntoz, ntqv, ntcw, ntiw, ntsw, ntlnc
220 integer, intent(in) :: ntinc, ntwa, ntia, ntke
221
222!MISC CONFIGURATION OPTIONS
223 INTEGER, PARAMETER :: &
224 & bl_mynn_mixscalars=1
225 LOGICAL :: &
226 & FLAG_QI, FLAG_QNI, FLAG_QC, FLAG_QS, FLAG_QNC, &
227 & FLAG_QNWFA, FLAG_QNIFA, FLAG_QNBCA, FLAG_OZONE
228 ! Define locally until needed from CCPP
229 LOGICAL, PARAMETER :: cycling = .false.
230
231!MYNN-1D
232 REAL(kind_phys), intent(in) :: delt, dtf
233 INTEGER, intent(in) :: im, levs
234 LOGICAL, intent(in) :: flag_init, flag_restart
235 INTEGER :: initflag, k, i
236 INTEGER :: IDS,IDE,JDS,JDE,KDS,KDE, &
237 & IMS,IME,JMS,JME,KMS,KME, &
238 & ITS,ITE,JTS,JTE,KTS,KTE
239
240 REAL(kind_phys) :: tem
241
242!MYNN-3D
243 real(kind_phys), dimension(:,:), intent(in) :: phii
244 real(kind_phys), dimension(:,:), intent(inout) :: &
245 & dtdt, dudt, dvdt, &
246 & dqdt_water_vapor, dqdt_liquid_cloud, dqdt_ice, &
247 & dqdt_snow, dqdt_ice_num_conc, dqdt_ozone
248 real(kind_phys), dimension(:,:), intent(inout), optional :: &
249 & dqdt_cloud_droplet_num_conc, dqdt_water_aer_num_conc, &
250 & dqdt_ice_aer_num_conc
251 real(kind_phys), dimension(:,:), intent(inout), optional :: qke, &
252 & EL_PBL, Sh3D, Sm3D, qc_bl, qi_bl, cldfra_bl, dqdt_cccn
253 real(kind_phys), dimension(:,:), intent(inout) :: &
254 & qke_adv
255 !These 10 arrays are only allocated when bl_mynn_output > 0
256 real(kind_phys), dimension(:,:), intent(inout), optional :: &
257 & edmf_a,edmf_w,edmf_qt, &
258 & edmf_thl,edmf_ent,edmf_qc, &
259 & sub_thl,sub_sqv,det_thl,det_sqv
260 real(kind_phys), dimension(:,:), intent(inout) :: &
261 & t3d,qgrs_water_vapor,qgrs_liquid_cloud,qgrs_ice, &
262 & qgrs_snow
263 real(kind_phys), dimension(:,:), intent(in) :: &
264 & qgrs_cloud_ice_num_conc, &
265 & u,v,omega, &
266 & exner,prsl,prsi, &
267 & qgrs_ozone
268 real(kind_phys), dimension(:,:), intent(in), optional :: &
269 & qgrs_water_aer_num_conc, &
270 & qgrs_cloud_droplet_num_conc, &
271 & qgrs_ice_aer_num_conc
272 real(kind_phys), dimension(:,:), intent(in), optional :: qgrs_cccn
273 real(kind_phys), dimension(:,:), intent(out), optional :: &
274 & Tsq, Qsq, Cov, exch_h, exch_m, dqke, qWT, qSHEAR, qBUOY, &
275 & qDISS
276 real(kind_phys), dimension(:), intent(in) :: xmu
277 real(kind_phys), dimension(:,:), intent(in) :: htrsw, htrlw
278 ! spp_wts_pbl only allocated if spp_pbl == 1
279 real(kind_phys), dimension(:,:), intent(in), optional :: spp_wts_pbl
280
281 !LOCAL
282 real(kind_phys), dimension(im,levs) :: &
283 & sqv,sqc,sqi,sqs,qnc,qni,ozone,qnwfa,qnifa,qnbca, &
284 & dz, w, p, rho, th, qv, delp, &
285 & RUBLTEN, RVBLTEN, RTHBLTEN, RQVBLTEN, &
286 & RQCBLTEN, RQNCBLTEN, RQIBLTEN, RQNIBLTEN, RQSBLTEN, &
287 & RQNWFABLTEN, RQNIFABLTEN, RQNBCABLTEN
288 real(kind_phys), allocatable :: old_ozone(:,:)
289
290!smoke/chem arrays
291 real(kind_phys), dimension(:), intent(inout), optional :: frp
292 logical, intent(in) :: mix_chem, enh_mix, rrfs_sd
293 real(kind_phys), dimension(:,:,:), intent(inout), optional :: chem3d
294 real(kind_phys), dimension(:,: ), intent(in), optional :: vdep
295 real(kind_phys), dimension(im) :: emis_ant_no
296
297!MYNN-2D
298 real(kind_phys), dimension(:), intent(in) :: &
299 & dx,zorl,slmsk,tsurf,qsfc,ps, &
300 & hflx,qflx,ust,wspd,rb,recmol
301
302 real(kind_phys), dimension(:), intent(in), optional :: &
303 & dusfc_cice,dvsfc_cice,dtsfc_cice,dqsfc_cice
304 real(kind_phys), dimension(:), intent(in) :: &
305 & stress_wat,hflx_wat,qflx_wat, &
306 & oceanfrac,fice
307
308 logical, dimension(:), intent(in) :: &
309 & wet, dry, icy
310
311 real(kind_phys), dimension(:), intent(inout) :: &
312 & pblh,dusfc_diag,dvsfc_diag,dtsfc_diag,dqsfc_diag
313 real(kind_phys), dimension(:), intent(out) :: &
314 & ch,dtsfc1,dqsfc1,dusfc1,dvsfc1, &
315 & dtsfci_diag,dqsfci_diag,dusfci_diag,dvsfci_diag
316 real(kind_phys), dimension(:), intent(out), optional :: &
317 & maxMF,maxwidth,ztop_plume
318 integer, dimension(:), intent(inout) :: &
319 & kpbl
320 integer, dimension(:), intent(inout), optional :: &
321 & ktop_plume
322
323 real(kind_phys), dimension(:), intent(inout), optional :: &
324 & dusfc_cpl,dvsfc_cpl,dtsfc_cpl,dqsfc_cpl
325 real(kind_phys), dimension(:), intent(inout), optional :: &
326 & dusfci_cpl,dvsfci_cpl,dtsfci_cpl,dqsfci_cpl
327
328 !LOCAL
329 real(kind_phys), dimension(im) :: &
330 & hfx,qfx,rmol,xland,uoce,voce,znt,ts
331 integer :: idtend
332 real(kind_phys), dimension(im) :: dusfci1,dvsfci1,dtsfci1,dqsfci1
333 real(kind_phys), allocatable :: save_qke_adv(:,:)
334 real(kind_phys), dimension(levs) :: kzero
335
336 ! Initialize CCPP error handling variables
337 errmsg = ''
338 errflg = 0
339
340 if (lprnt) then
341 write(0,*)"=============================================="
342 write(0,*)"in mynn wrapper..."
343 write(0,*)"flag_init=",flag_init
344 write(0,*)"flag_restart=",flag_restart
345 endif
346
347 if (.not. flag_for_pbl_generic_tend .and. ldiag3d) then
348 idtend = dtidx(ntke+100,index_of_process_pbl)
349 if (idtend>=1) then
350 allocate(save_qke_adv(im,levs))
351 save_qke_adv=qke_adv
352 endif
353 endif
354
355 ! DH* TODO: Use flag_restart to distinguish which fields need
356 ! to be initialized and which are read from restart files
357 if (flag_init) then
358 initflag=1
359 !print*,"in MYNN, initflag=",initflag
360 else
361 initflag=0
362 !print*,"in MYNN, initflag=",initflag
363 endif
364
365 kzero = zero !generic zero array
366 !initialize arrays for test
367 emis_ant_no = 0.
368
369 flag_ozone = ntoz>0
370
371 ! Assign variables for each microphysics scheme
372 if (imp_physics == imp_physics_wsm6 .or. imp_physics == imp_physics_fa) then
373 ! WSM6 or Ferrier-Aligo
374 flag_qi = .true.
375 flag_qni= .false.
376 flag_qc = .true.
377 flag_qnc= .false.
378 flag_qs = .false.
379 flag_qnwfa= .false.
380 flag_qnifa= .false.
381 flag_qnbca= .false.
382 do k=1,levs
383 do i=1,im
384 sqv(i,k) = qgrs_water_vapor(i,k)
385 sqc(i,k) = qgrs_liquid_cloud(i,k)
386 sqi(i,k) = qgrs_ice(i,k)
387 sqs(i,k) = 0.
388 ozone(i,k) = qgrs_ozone(i,k)
389 qnc(i,k) = 0.
390 qni(i,k) = 0.
391 qnwfa(i,k) = 0.
392 qnifa(i,k) = 0.
393 qnbca(i,k) = 0.
394 enddo
395 enddo
396 elseif (imp_physics == imp_physics_nssl ) then
397 ! NSSL
398 flag_qi = .true.
399 flag_qni= .true.
400 flag_qc = .true.
401 flag_qnc= .true.
402 flag_qs = .true.
403 flag_qnwfa= nssl_ccn_on ! ERM: Perhaps could use this field for CCN field?
404 flag_qnifa= .false.
405 flag_qnbca= .false.
406 do k=1,levs
407 do i=1,im
408 sqv(i,k) = qgrs_water_vapor(i,k)
409 sqc(i,k) = qgrs_liquid_cloud(i,k)
410 sqi(i,k) = qgrs_ice(i,k)
411 sqs(i,k) = qgrs_snow(i,k)
412 ozone(i,k) = qgrs_ozone(i,k)
413 qnc(i,k) = qgrs_cloud_droplet_num_conc(i,k)
414 qni(i,k) = qgrs_cloud_ice_num_conc(i,k)
415 qnwfa(i,k) = 0.
416 IF ( nssl_ccn_on ) THEN
417 qnwfa(i,k) = qgrs_cccn(i,k)
418 ENDIF
419 qnifa(i,k) = 0.
420 qnbca(i,k) = 0.
421 enddo
422 enddo
423 elseif (imp_physics == imp_physics_thompson) then
424 ! Thompson
425 if(ltaerosol) then
426 flag_qi = .true.
427 flag_qni= .true.
428 flag_qc = .true.
429 flag_qs = .true. !pipe it in, but do not mix
430 flag_qnc= .true.
431 flag_qnwfa= .true.
432 flag_qnifa= .true.
433 flag_qnbca= .false.
434 do k=1,levs
435 do i=1,im
436 sqv(i,k) = qgrs_water_vapor(i,k)
437 sqc(i,k) = qgrs_liquid_cloud(i,k)
438 sqi(i,k) = qgrs_ice(i,k)
439 sqs(i,k) = qgrs_snow(i,k)
440 qnc(i,k) = qgrs_cloud_droplet_num_conc(i,k)
441 qni(i,k) = qgrs_cloud_ice_num_conc(i,k)
442 ozone(i,k) = qgrs_ozone(i,k)
443 qnwfa(i,k) = qgrs_water_aer_num_conc(i,k)
444 qnifa(i,k) = qgrs_ice_aer_num_conc(i,k)
445 qnbca(i,k) = 0.
446 enddo
447 enddo
448 else if(mraerosol) then
449 flag_qi = .true.
450 flag_qni= .true.
451 flag_qc = .true.
452 flag_qs = .true.
453 flag_qnc= .true.
454 flag_qnwfa= .false.
455 flag_qnifa= .false.
456 flag_qnbca= .false.
457 do k=1,levs
458 do i=1,im
459 sqv(i,k) = qgrs_water_vapor(i,k)
460 sqc(i,k) = qgrs_liquid_cloud(i,k)
461 sqi(i,k) = qgrs_ice(i,k)
462 sqs(i,k) = qgrs_snow(i,k)
463 qnc(i,k) = qgrs_cloud_droplet_num_conc(i,k)
464 qni(i,k) = qgrs_cloud_ice_num_conc(i,k)
465 ozone(i,k) = qgrs_ozone(i,k)
466 qnwfa(i,k) = 0.
467 qnifa(i,k) = 0.
468 qnbca(i,k) = 0.
469 enddo
470 enddo
471 else
472 flag_qi = .true.
473 flag_qni= .true.
474 flag_qc = .true.
475 flag_qs = .true.
476 flag_qnc= .false.
477 flag_qnwfa= .false.
478 flag_qnifa= .false.
479 flag_qnbca= .false.
480 do k=1,levs
481 do i=1,im
482 sqv(i,k) = qgrs_water_vapor(i,k)
483 sqc(i,k) = qgrs_liquid_cloud(i,k)
484 sqi(i,k) = qgrs_ice(i,k)
485 sqs(i,k) = qgrs_snow(i,k)
486 qnc(i,k) = 0.
487 qni(i,k) = qgrs_cloud_ice_num_conc(i,k)
488 ozone(i,k) = qgrs_ozone(i,k)
489 qnwfa(i,k) = 0.
490 qnifa(i,k) = 0.
491 qnbca(i,k) = 0.
492 enddo
493 enddo
494 endif
495 elseif (imp_physics == imp_physics_gfdl) then
496 ! GFDL MP
497 flag_qi = .true.
498 flag_qni= .false.
499 flag_qc = .true.
500 flag_qnc= .false.
501 flag_qs = .false.
502 flag_qnwfa= .false.
503 flag_qnifa= .false.
504 flag_qnbca= .false.
505 do k=1,levs
506 do i=1,im
507 sqv(i,k) = qgrs_water_vapor(i,k)
508 sqc(i,k) = qgrs_liquid_cloud(i,k)
509 sqi(i,k) = qgrs_ice(i,k)
510 qnc(i,k) = 0.
511 qni(i,k) = 0.
512 sqs(i,k) = 0.
513 qnwfa(i,k) = 0.
514 qnifa(i,k) = 0.
515 qnbca(i,k) = 0.
516 ozone(i,k) = qgrs_ozone(i,k)
517 enddo
518 enddo
519 else
520 print*,"In MYNN wrapper. Unknown microphysics scheme, imp_physics=",imp_physics
521 print*,"Defaulting to qc and qv species only..."
522 flag_qi = .false.
523 flag_qni= .false.
524 flag_qc = .true.
525 flag_qnc= .false.
526 flag_qs = .false.
527 flag_qnwfa= .false.
528 flag_qnifa= .false.
529 flag_qnbca= .false.
530 do k=1,levs
531 do i=1,im
532 sqv(i,k) = qgrs_water_vapor(i,k)
533 sqc(i,k) = qgrs_liquid_cloud(i,k)
534 sqi(i,k) = 0.
535 sqs(i,k) = 0.
536 qnc(i,k) = 0.
537 qni(i,k) = 0.
538 qnwfa(i,k) = 0.
539 qnifa(i,k) = 0.
540 qnbca(i,k) = 0.
541 ozone(i,k) = qgrs_ozone(i,k)
542 enddo
543 enddo
544 endif
545 if(ldiag3d .and. dtidx(100+ntoz,index_of_process_pbl)>1) then
546 allocate(old_ozone(im,levs))
547 old_ozone = ozone
548 endif
549
550 do k=1,levs
551 do i=1,im
552 th(i,k)=t3d(i,k)/exner(i,k)
553 rho(i,k)=prsl(i,k)/(r_d*t3d(i,k)*(1.+p608*max(sqv(i,k),1e-8)))
554 w(i,k) = -omega(i,k)/(rho(i,k)*grav)
555 enddo
556 enddo
557
558 ! Check incoming moist species to ensure non-negative values
559 ! First, create height difference (dz)
560 do k=1,levs
561 do i=1,im
562 dz(i,k)=(phii(i,k+1) - phii(i,k))*g_inv
563 enddo
564 enddo
565
566 do i=1,im
567 do k=1,levs
568 delp(i,k) = prsi(i,k) - prsi(i,k+1)
569 enddo
570 enddo
571
572 do i=1,im
573 call moisture_check2(levs, delt, &
574 delp(i,:), exner(i,:), &
575 sqv(i,:), sqc(i,:), &
576 sqi(i,:), kzero(:), &
577 t3d(i,:) )
578 enddo
579
580 !intialize more variables
581 do i=1,im
582 if (slmsk(i)==1. .or. slmsk(i)==2.) then !sea/land/ice mask (=0/1/2) in FV3
583 xland(i)=1.0 !but land/water = (1/2) in SFCLAY_mynn
584 else
585 xland(i)=2.0
586 endif
587 uoce(i)=0.0
588 voce(i)=0.0
589 !ust(i) = sqrt(stress(i))
590 ch(i)=0.0
591 hfx(i)=hflx(i)*rho(i,1)*cp
592 qfx(i)=qflx(i)*rho(i,1)
593 !filter bad incoming fluxes
594 if (hfx(i) > 1200.)hfx(i) = 1200.
595 if (hfx(i) < -500.)hfx(i) = -500.
596 if (qfx(i) > .0005)qfx(i) = 0.0005
597 if (qfx(i) < -.0002)qfx(i) = -0.0002
598
599 dtsfc1(i) = hfx(i)
600 dqsfc1(i) = qfx(i)*xlv
601 dusfc1(i) = -1.*rho(i,1)*ust(i)*ust(i)*u(i,1)/wspd(i)
602 dvsfc1(i) = -1.*rho(i,1)*ust(i)*ust(i)*v(i,1)/wspd(i)
603
604 !BWG: diagnostic surface fluxes for scalars & momentum
605 dtsfci_diag(i) = dtsfc1(i)
606 dqsfci_diag(i) = dqsfc1(i)
607 dtsfc_diag(i) = dtsfc_diag(i) + dtsfc1(i)*delt
608 dqsfc_diag(i) = dqsfc_diag(i) + dqsfc1(i)*delt
609 dusfci_diag(i) = dusfc1(i)
610 dvsfci_diag(i) = dvsfc1(i)
611 dusfc_diag(i) = dusfc_diag(i) + dusfci_diag(i)*delt
612 dvsfc_diag(i) = dvsfc_diag(i) + dvsfci_diag(i)*delt
613
614 znt(i)=zorl(i)*0.01 !cm -> m?
615 if (do_mynnsfclay) then
616 rmol(i)=recmol(i)
617 else
618 if (hfx(i) .ge. 0.)then
619 rmol(i)=-hfx(i)/(200.*dz(i,1)*0.5)
620 else
621 rmol(i)=abs(rb(i))*1./(dz(i,1)*0.5)
622 endif
623 endif
624 ts(i)=tsurf(i)/exner(i,1) !theta
625! qsfc(i)=qss(i)
626! ps(i)=pgr(i)
627! wspd(i)=wind(i)
628 enddo
629
630 ! BWG: Coupling insertion
631 if (cplflx) then
632 do i=1,im
633 if (oceanfrac(i) > zero) then ! Ocean only, NO LAKES
634 if ( .not. wet(i)) then ! no open water, use results from CICE
635 dusfci_cpl(i) = dusfc_cice(i)
636 dvsfci_cpl(i) = dvsfc_cice(i)
637 dtsfci_cpl(i) = dtsfc_cice(i)
638 dqsfci_cpl(i) = dqsfc_cice(i)
639 elseif (icy(i) .or. dry(i)) then ! use stress_ocean for opw component at mixed point
640 if (wspd(i) > zero) then
641 dusfci_cpl(i) = -1.*rho(i,1)*stress_wat(i)*u(i,1)/wspd(i) ! U-momentum flux
642 dvsfci_cpl(i) = -1.*rho(i,1)*stress_wat(i)*v(i,1)/wspd(i) ! V-momentum flux
643 else
644 dusfci_cpl(i) = zero
645 dvsfci_cpl(i) = zero
646 endif
647 dtsfci_cpl(i) = cp*rho(i,1)*hflx_wat(i) ! sensible heat flux over open ocean
648 dqsfci_cpl(i) = xlv*rho(i,1)*qflx_wat(i) ! latent heat flux over open ocean
649 else ! use results from this scheme for 100% open ocean
650 dusfci_cpl(i) = dusfci_diag(i)
651 dvsfci_cpl(i) = dvsfci_diag(i)
652 dtsfci_cpl(i) = dtsfci_diag(i)
653 dqsfci_cpl(i) = dqsfci_diag(i)
654 endif
655!
656 dusfc_cpl(i) = dusfc_cpl(i) + dusfci_cpl(i) * delt
657 dvsfc_cpl(i) = dvsfc_cpl(i) + dvsfci_cpl(i) * delt
658 dtsfc_cpl(i) = dtsfc_cpl(i) + dtsfci_cpl(i) * delt
659 dqsfc_cpl(i) = dqsfc_cpl(i) + dqsfci_cpl(i) * delt
660 else ! If no ocean
661 dusfc_cpl(i) = huge
662 dvsfc_cpl(i) = huge
663 dtsfc_cpl(i) = huge
664 dqsfc_cpl(i) = huge
665 endif ! Ocean only, NO LAKES
666 enddo
667 endif ! End coupling insertion
668
669 if (lprnt) then
670 print*
671 write(0,*)"===CALLING mynn_bl_driver; input:"
672 print*,"tke_budget=",tke_budget," bl_mynn_tkeadvect=",bl_mynn_tkeadvect
673 print*,"bl_mynn_cloudpdf=",bl_mynn_cloudpdf," bl_mynn_mixlength=",bl_mynn_mixlength
674 print*,"bl_mynn_edmf=",bl_mynn_edmf," bl_mynn_edmf_mom=",bl_mynn_edmf_mom
675 print*,"bl_mynn_edmf_tke=",bl_mynn_edmf_tke
676 print*,"bl_mynn_cloudmix=",bl_mynn_cloudmix," bl_mynn_mixqt=",bl_mynn_mixqt
677 print*,"icloud_bl=",icloud_bl
678 print*,"T:",t3d(1,1),t3d(1,2),t3d(1,levs)
679 print*,"TH:",th(1,1),th(1,2),th(1,levs)
680 print*,"rho:",rho(1,1),rho(1,2),rho(1,levs)
681 print*,"exner:",exner(1,1),exner(1,2),exner(1,levs)
682 print*,"prsl:",prsl(1,1),prsl(1,2),prsl(1,levs)
683 print*,"dz:",dz(1,1),dz(1,2),dz(1,levs)
684 print*,"u:",u(1,1),u(1,2),u(1,levs)
685 print*,"v:",v(1,1),v(1,2),v(1,levs)
686 print*,"sqv:",sqv(1,1),sqv(1,2),sqv(1,levs)
687 print*,"sqc:",sqc(1,1),sqc(1,2),sqc(1,levs)
688 print*,"sqi:",sqi(1,1),sqi(1,2),sqi(1,levs)
689 print*,"rmol:",rmol(1)," ust:",ust(1)
690 print*," dx=",dx(1),"initflag=",initflag
691 print*,"Tsurf:",tsurf(1)," Thetasurf:",ts(1)
692 print*,"HFX:",hfx(1)," qfx",qfx(1)
693 print*,"qsfc:",qsfc(1)," ps:",ps(1)
694 print*,"wspd:",wspd(1)," rb=",rb(1)
695 print*,"znt:",znt(1)," delt=",delt
696 print*,"im=",im," levs=",levs
697 print*,"PBLH=",pblh(1)," KPBL=",kpbl(1)," xland=",xland(1)
698 print*,"ch=",ch(1)
699 !print*,"TKE:",TKE_PBL(1,1),TKE_PBL(1,2),TKE_PBL(1,levs)
700 print*,"qke:",qke(1,1),qke(1,2),qke(1,levs)
701 print*,"el_pbl:",el_pbl(1,1),el_pbl(1,2),el_pbl(1,levs)
702 print*,"Sh3d:",sh3d(1,1),sh3d(1,2),sh3d(1,levs)
703 !print*,"exch_h:",exch_h(1,1),exch_h(1,2),exch_h(1,levs) ! - intent(out)
704 !print*,"exch_m:",exch_m(1,1),exch_m(1,2),exch_m(1,levs) ! - intent(out)
705 print*,"max cf_bl:",maxval(cldfra_bl(1,:))
706 endif
707
708
709 CALL mynn_bl_driver( &
710 & initflag=initflag,restart=flag_restart, &
711 & cycling=cycling, &
712 & delt=delt,dz=dz,dx=dx,znt=znt, &
713 & u=u,v=v,w=w,th=th,sqv3d=sqv,sqc3d=sqc, &
714 & sqi3d=sqi,sqs3d=sqs,qnc=qnc,qni=qni, &
715 & qnwfa=qnwfa,qnifa=qnifa,qnbca=qnbca,ozone=ozone, &
716 & p=prsl,exner=exner,rho=rho,t3d=t3d, &
717 & xland=xland,ts=ts,qsfc=qsfc,ps=ps, &
718 & ust=ust,ch=ch,hfx=hfx,qfx=qfx,rmol=rmol, &
719 & wspd=wspd,uoce=uoce,voce=voce, & !input
720 & qke=qke,qke_adv=qke_adv, & !output
721 & sh3d=sh3d,sm3d=sm3d, &
722!chem/smoke
723 & nchem=nchem,kdvel=kdvel,ndvel=ndvel, &
724 & chem3d=chem3d,vdep=vdep,smoke_dbg=smoke_dbg, &
725 & frp=frp,emis_ant_no=emis_ant_no, &
726 & mix_chem=mix_chem,enh_mix=enh_mix, &
727 & rrfs_sd=rrfs_sd, &
728!-----
729 & tsq=tsq,qsq=qsq,cov=cov, & !output
730 & rublten=rublten,rvblten=rvblten,rthblten=rthblten, & !output
731 & rqvblten=rqvblten,rqcblten=rqcblten, &
732 & rqiblten=rqiblten,rqncblten=rqncblten, & !output
733 & rqsblten=rqsblten, & !output
734 & rqniblten=rqniblten,rqnwfablten=rqnwfablten, & !output
735 & rqnifablten=rqnifablten,rqnbcablten=rqnbcablten, & !output
736 & dozone=dqdt_ozone, & !output
737 & exch_h=exch_h,exch_m=exch_m, & !output
738 & pblh=pblh,kpbl=kpbl, & !output
739 & el_pbl=el_pbl, & !output
740 & dqke=dqke, & !output
741 & qwt=qwt,qshear=qshear,qbuoy=qbuoy,qdiss=qdiss, & !output
742 & bl_mynn_tkeadvect=bl_mynn_tkeadvect, &
743 & tke_budget=tke_budget, & !input parameter
744 & bl_mynn_cloudpdf=bl_mynn_cloudpdf, & !input parameter
745 & bl_mynn_mixlength=bl_mynn_mixlength, & !input parameter
746 & icloud_bl=icloud_bl, & !input parameter
747 & qc_bl=qc_bl,qi_bl=qi_bl,cldfra_bl=cldfra_bl, & !output
748 & closure=bl_mynn_closure,bl_mynn_edmf=bl_mynn_edmf, & !input parameter
749 & bl_mynn_edmf_mom=bl_mynn_edmf_mom, & !input parameter
750 & bl_mynn_edmf_tke=bl_mynn_edmf_tke, & !input parameter
751 & bl_mynn_mixscalars=bl_mynn_mixscalars, & !input parameter
752 & bl_mynn_output=bl_mynn_output, & !input parameter
753 & bl_mynn_cloudmix=bl_mynn_cloudmix, & !input parameter
754 & bl_mynn_mixqt=bl_mynn_mixqt, & !input parameter
755 & edmf_a=edmf_a,edmf_w=edmf_w,edmf_qt=edmf_qt, & !output
756 & edmf_thl=edmf_thl,edmf_ent=edmf_ent,edmf_qc=edmf_qc,& !output
757 & sub_thl3d=sub_thl,sub_sqv3d=sub_sqv, &
758 & det_thl3d=det_thl,det_sqv3d=det_sqv, &
759 & maxwidth=maxwidth,maxmf=maxmf,ztop_plume=ztop_plume,& !output
760 & ktop_plume=ktop_plume, & !output
761 & spp_pbl=spp_pbl,pattern_spp_pbl=spp_wts_pbl, & !input
762 & rthraten=htrlw, & !input
763 & flag_qi=flag_qi,flag_qni=flag_qni, & !input
764 & flag_qc=flag_qc,flag_qnc=flag_qnc,flag_qs=flag_qs, & !input
765 & flag_qnwfa=flag_qnwfa,flag_qnifa=flag_qnifa, & !input
766 & flag_qnbca=flag_qnbca,flag_ozone=flag_ozone, & !input
767 & ids=1,ide=im,jds=1,jde=1,kds=1,kde=levs, & !input
768 & ims=1,ime=im,jms=1,jme=1,kms=1,kme=levs, & !input
769 & its=1,ite=im,jts=1,jte=1,kts=1,kte=levs ) !input
770
771
772 ! POST MYNN (INTERSTITIAL) WORK:
773 !update/save MYNN-only variables
774 !do k=1,levs
775 ! do i=1,im
776 ! gq0(i,k,4)=qke(i,k,1) !tke*2
777 ! enddo
778 !enddo
779 !For MYNN, convert TH-tend to T-tend
780 do k = 1, levs
781 do i = 1, im
782 dtdt(i,k) = dtdt(i,k) + rthblten(i,k)*exner(i,k)
783 dudt(i,k) = dudt(i,k) + rublten(i,k)
784 dvdt(i,k) = dvdt(i,k) + rvblten(i,k)
785 enddo
786 enddo
787 accum_duvt3dt: if(ldiag3d .or. lsidea) then
788 call dtend_helper(index_of_x_wind,rublten)
789 call dtend_helper(index_of_y_wind,rvblten)
790 call dtend_helper(index_of_temperature,rthblten,exner)
791 if(ldiag3d) then
792 call dtend_helper(100+ntoz,dqdt_ozone)
793 ! idtend = dtidx(100+ntoz,index_of_process_pbl)
794 ! if(idtend>=1) then
795 ! dtend(:,:,idtend) = dtend(:,:,idtend) + (ozone-old_ozone)
796 ! deallocate(old_ozone)
797 ! endif
798 endif
799 endif accum_duvt3dt
800 !Update T, U and V:
801 !do k = 1, levs
802 ! do i = 1, im
803 ! T3D(i,k) = T3D(i,k) + RTHBLTEN(i,k)*exner(i,k)*delt
804 ! u(i,k) = u(i,k) + RUBLTEN(i,k)*delt
805 ! v(i,k) = v(i,k) + RVBLTEN(i,k)*delt
806 ! enddo
807 !enddo
808
809 !DO moist/scalar/tracer tendencies:
810 if (imp_physics == imp_physics_wsm6 .or. imp_physics == imp_physics_fa) then
811 ! WSM6
812 do k=1,levs
813 do i=1,im
814 dqdt_water_vapor(i,k) = rqvblten(i,k) !/(1.0 + qv(i,k))
815 dqdt_liquid_cloud(i,k) = rqcblten(i,k) !/(1.0 + qv(i,k))
816 dqdt_ice(i,k) = rqiblten(i,k) !/(1.0 + qv(i,k))
817 dqdt_snow(i,k) = rqsblten(i,k) !/(1.0 + qv(i,k))
818 !dqdt_ozone(i,k) = 0.0
819 enddo
820 enddo
821 if(ldiag3d .and. .not. flag_for_pbl_generic_tend) then
822 call dtend_helper(100+ntqv,rqvblten)
823 call dtend_helper(100+ntcw,rqcblten)
824 call dtend_helper(100+ntiw,rqiblten)
825 endif
826 !Update moist species:
827 !do k=1,levs
828 ! do i=1,im
829 ! qgrs_water_vapor(i,k) = qgrs_water_vapor(i,k) + (RQVBLTEN(i,k)/(1.0+RQVBLTEN(i,k)))*delt
830 ! qgrs_liquid_cloud(i,k) = qgrs_liquid_cloud(i,k) + RQCBLTEN(i,k)*delt
831 ! qgrs_ice(i,k) = qgrs_ice(i,k) + RQIBLTEN(i,k)*delt
832 ! !dqdt_ozone(i,k) = 0.0
833 ! enddo
834 !enddo
835 elseif (imp_physics == imp_physics_thompson) then
836 ! Thompson-Aerosol
837 if(ltaerosol) then
838 do k=1,levs
839 do i=1,im
840 dqdt_water_vapor(i,k) = rqvblten(i,k) !/(1.0 + qv(i,k))
841 dqdt_liquid_cloud(i,k) = rqcblten(i,k) !/(1.0 + qv(i,k))
842 dqdt_cloud_droplet_num_conc(i,k) = rqncblten(i,k)
843 dqdt_ice(i,k) = rqiblten(i,k) !/(1.0 + qv(i,k))
844 dqdt_ice_num_conc(i,k) = rqniblten(i,k)
845 dqdt_snow(i,k) = rqsblten(i,k) !/(1.0 + qv(i,k))
846 !dqdt_ozone(i,k) = 0.0
847 dqdt_water_aer_num_conc(i,k) = rqnwfablten(i,k)
848 dqdt_ice_aer_num_conc(i,k) = rqnifablten(i,k)
849 enddo
850 enddo
851 if(ldiag3d .and. .not. flag_for_pbl_generic_tend) then
852 call dtend_helper(100+ntqv,rqvblten)
853 call dtend_helper(100+ntcw,rqcblten)
854 call dtend_helper(100+ntlnc,rqncblten)
855 call dtend_helper(100+ntiw,rqiblten)
856 call dtend_helper(100+ntinc,rqniblten)
857 call dtend_helper(100+ntwa,rqnwfablten)
858 call dtend_helper(100+ntia,rqnifablten)
859 endif
860 !do k=1,levs
861 ! do i=1,im
862 ! qgrs_water_vapor(i,k) = qgrs_water_vapor(i,k) + (RQVBLTEN(i,k)/(1.0+RQVBLTEN(i,k)))*delt
863 ! qgrs_liquid_cloud(i,k) = qgrs_liquid_cloud(i,k) + RQCBLTEN(i,k)*delt
864 ! qgrs_ice(i,k) = qgrs_ice(i,k) + RQIBLTEN(i,k)*delt
865 ! qgrs_cloud_droplet_num_conc(i,k) = qgrs_cloud_droplet_num_conc(i,k) + RQNCBLTEN(i,k)*delt
866 ! qgrs_cloud_ice_num_conc(i,k) = qgrs_cloud_ice_num_conc(i,k) + RQNIBLTEN(i,k)*delt
867 ! !dqdt_ozone(i,k) = 0.0
868 ! !qgrs_water_aer_num_conc(i,k) = qgrs_water_aer_num_conc(i,k) + RQNWFABLTEN(i,k)*delt
869 ! !qgrs_ice_aer_num_conc(i,k) = qgrs_ice_aer_num_conc(i,k) + RQNIFABLTEN(i,k)*delt
870 ! enddo
871 !enddo
872 else if(mraerosol) then
873 do k=1,levs
874 do i=1,im
875 dqdt_water_vapor(i,k) = rqvblten(i,k) !/(1.0 + qv(i,k))
876 dqdt_liquid_cloud(i,k) = rqcblten(i,k) !/(1.0 + qv(i,k))
877 dqdt_cloud_droplet_num_conc(i,k) = rqncblten(i,k)
878 dqdt_ice(i,k) = rqiblten(i,k) !/(1.0 + qv(i,k))
879 dqdt_ice_num_conc(i,k) = rqniblten(i,k)
880 dqdt_snow(i,k) = rqsblten(i,k) !/(1.0 + qv(i,k))
881 enddo
882 enddo
883 if(ldiag3d .and. .not. flag_for_pbl_generic_tend) then
884 call dtend_helper(100+ntqv,rqvblten)
885 call dtend_helper(100+ntcw,rqcblten)
886 call dtend_helper(100+ntlnc,rqncblten)
887 call dtend_helper(100+ntiw,rqiblten)
888 call dtend_helper(100+ntinc,rqniblten)
889 endif
890 else
891 !Thompson (2008)
892 do k=1,levs
893 do i=1,im
894 dqdt_water_vapor(i,k) = rqvblten(i,k) !/(1.0 + qv(i,k))
895 dqdt_liquid_cloud(i,k) = rqcblten(i,k) !/(1.0 + qv(i,k))
896 dqdt_ice(i,k) = rqiblten(i,k) !/(1.0 + qv(i,k))
897 dqdt_ice_num_conc(i,k) = rqniblten(i,k)
898 dqdt_snow(i,k) = rqsblten(i,k) !/(1.0 + qv(i,k))
899 !dqdt_ozone(i,k) = 0.0
900 enddo
901 enddo
902 if(ldiag3d .and. .not. flag_for_pbl_generic_tend) then
903 call dtend_helper(100+ntqv,rqvblten)
904 call dtend_helper(100+ntcw,rqcblten)
905 call dtend_helper(100+ntiw,rqiblten)
906 call dtend_helper(100+ntinc,rqniblten)
907 !call dtend_helper(100+ntsw,RQSBLTEN)
908 endif
909 !do k=1,levs
910 ! do i=1,im
911 ! qgrs_water_vapor(i,k) = qgrs_water_vapor(i,k) + (RQVBLTEN(i,k)/(1.0+RQVBLTEN(i,k)))*delt
912 ! qgrs_liquid_cloud(i,k) = qgrs_liquid_cloud(i,k) + RQCBLTEN(i,k)*delt
913 ! qgrs_ice(i,k) = qgrs_ice(i,k) + RQIBLTEN(i,k)*delt
914 ! qgrs_cloud_ice_num_conc(i,k) = qgrs_cloud_ice_num_conc(i,k) + RQNIBLTEN(i,k)*delt
915 ! !dqdt_ozone(i,k) = 0.0
916 ! enddo
917 !enddo
918 endif !end thompson choice
919 elseif (imp_physics == imp_physics_nssl) then
920 ! NSSL
921 do k=1,levs
922 do i=1,im
923 dqdt_water_vapor(i,k) = rqvblten(i,k) !/(1.0 + qv(i,k))
924 dqdt_liquid_cloud(i,k) = rqcblten(i,k) !/(1.0 + qv(i,k))
925 dqdt_cloud_droplet_num_conc(i,k) = rqncblten(i,k)
926 dqdt_ice(i,k) = rqiblten(i,k) !/(1.0 + qv(i,k))
927 dqdt_ice_num_conc(i,k) = rqniblten(i,k)
928 dqdt_snow(i,k) = rqsblten(i,k) !/(1.0 + qv(i,k))
929 IF ( nssl_ccn_on ) THEN !
930 dqdt_cccn(i,k) = rqnwfablten(i,k)
931 ENDIF
932 enddo
933 enddo
934
935 elseif (imp_physics == imp_physics_gfdl) then
936 ! GFDL MP
937 do k=1,levs
938 do i=1,im
939 dqdt_water_vapor(i,k) = rqvblten(i,k) !/(1.0 + qv(i,k))
940 dqdt_liquid_cloud(i,k) = rqcblten(i,k) !/(1.0 + qv(i,k))
941 dqdt_ice(i,k) = rqiblten(i,k) !/(1.0 + qv(i,k))
942 !dqdt_rain(i,k) = 0.0
943 !dqdt_snow(i,k) = 0.0
944 !dqdt_graupel(i,k) = 0.0
945 !dqdt_ozone(i,k) = 0.0
946 enddo
947 enddo
948 if(ldiag3d .and. .not. flag_for_pbl_generic_tend) then
949 call dtend_helper(100+ntqv,rqvblten)
950 call dtend_helper(100+ntcw,rqcblten)
951 call dtend_helper(100+ntiw,rqiblten)
952 endif
953 !do k=1,levs
954 ! do i=1,im
955 ! qgrs_water_vapor(i,k) = qgrs_water_vapor(i,k) + (RQVBLTEN(i,k)/(1.0+RQVBLTEN(i,k)))*delt
956 ! qgrs_liquid_cloud(i,k) = qgrs_liquid_cloud(i,k) + RQCBLTEN(i,k)*delt
957 ! qgrs_ice(i,k) = qgrs_ice(i,k) + RQIBLTEN(i,k)*delt
958 ! !dqdt_ozone(i,k) = 0.0
959 ! enddo
960 !enddo
961 else
962! print*,"In MYNN wrapper. Unknown microphysics scheme, imp_physics=",imp_physics
963 do k=1,levs
964 do i=1,im
965 dqdt_water_vapor(i,k) = rqvblten(i,k) !/(1.0 + qv(i,k))
966 dqdt_liquid_cloud(i,k) = rqcblten(i,k) !/(1.0 + qv(i,k))
967 dqdt_ice(i,k) = 0.0
968 !dqdt_rain(i,k) = 0.0
969 !dqdt_snow(i,k) = 0.0
970 !dqdt_graupel(i,k) = 0.0
971 !dqdt_ozone(i,k) = 0.0
972 enddo
973 enddo
974 if(ldiag3d .and. .not. flag_for_pbl_generic_tend) then
975 call dtend_helper(100+ntqv,rqvblten)
976 call dtend_helper(100+ntcw,rqcblten)
977 call dtend_helper(100+ntiw,rqiblten)
978 endif
979 endif
980
981 if (lprnt) then
982 print*
983 print*,"===Finished with mynn_bl_driver; output:"
984 print*,"T:",t3d(1,1),t3d(1,2),t3d(1,levs)
985 print*,"TH:",th(1,1),th(1,2),th(1,levs)
986 print*,"rho:",rho(1,1),rho(1,2),rho(1,levs)
987 print*,"exner:",exner(1,1),exner(1,2),exner(1,levs)
988 print*,"prsl:",prsl(1,1),prsl(1,2),prsl(1,levs)
989 print*,"dz:",dz(1,1),dz(1,2),dz(1,levs)
990 print*,"u:",u(1,1),u(1,2),u(1,levs)
991 print*,"v:",v(1,1),v(1,2),v(1,levs)
992 print*,"sqv:",sqv(1,1),sqv(1,2),sqv(1,levs)
993 print*,"sqc:",sqc(1,1),sqc(1,2),sqc(1,levs)
994 print*,"sqi:",sqi(1,1),sqi(1,2),sqi(1,levs)
995 print*,"rmol:",rmol(1)," ust:",ust(1)
996 print*,"dx(1)=",dx(1),"initflag=",initflag
997 print*,"Tsurf:",tsurf(1)," Thetasurf:",ts(1)
998 print*,"HFX:",hfx(1)," qfx",qfx(1)
999 print*,"qsfc:",qsfc(1)," ps:",ps(1)
1000 print*,"wspd:",wspd(1)," rb=",rb(1)
1001 print*,"znt:",znt(1)," delt=",delt
1002 print*,"im=",im," levs=",levs
1003 print*,"PBLH=",pblh(1)," KPBL=",kpbl(1)," xland=",xland(1)
1004 print*,"ch=",ch(1)
1005 print*,"qke:",qke(1,1),qke(1,2),qke(1,levs)
1006 print*,"el_pbl:",el_pbl(1,1),el_pbl(1,2),el_pbl(1,levs)
1007 print*,"Sh3d:",sh3d(1,1),sh3d(1,2),sh3d(1,levs)
1008 print*,"exch_h:",exch_h(1,1),exch_h(1,2),exch_h(1,levs)
1009 print*,"exch_m:",exch_m(1,1),exch_m(1,2),exch_m(1,levs)
1010 print*,"max cf_bl:",maxval(cldfra_bl(1,:))
1011 print*,"max qc_bl:",maxval(qc_bl(1,:))
1012 print*,"dtdt:",dtdt(1,1),dtdt(1,2),dtdt(1,levs)
1013 print*,"dudt:",dudt(1,1),dudt(1,2),dudt(1,levs)
1014 print*,"dvdt:",dvdt(1,1),dvdt(1,2),dvdt(1,levs)
1015 print*,"dqdt:",dqdt_water_vapor(1,1),dqdt_water_vapor(1,2),dqdt_water_vapor(1,levs)
1016 print*,"ztop_plume:",ztop_plume(1)," maxmf:",maxmf(1)
1017 print*,"maxwidth:",maxwidth(1)
1018 print*
1019 endif
1020
1021 if(allocated(save_qke_adv)) then
1022 if(ldiag3d .and. .not. flag_for_pbl_generic_tend) then
1023 idtend = dtidx(100+ntke,index_of_process_pbl)
1024 if(idtend>=1) then
1025 dtend(:,:,idtend) = dtend(:,:,idtend) + qke_adv-save_qke_adv
1026 endif
1027 endif
1028 deallocate(save_qke_adv)
1029 endif
1030
1031 CONTAINS
1032
1033 SUBROUTINE dtend_helper(itracer,field,mult)
1034 real(kind_phys), intent(in) :: field(im,levs)
1035 real(kind_phys), intent(in), optional :: mult(im,levs)
1036 integer, intent(in) :: itracer
1037 integer :: idtend
1038
1039 idtend=dtidx(itracer,index_of_process_pbl)
1040 if(idtend>=1) then
1041 if(present(mult)) then
1042 dtend(:,:,idtend) = dtend(:,:,idtend) + field*dtf*mult
1043 else
1044 dtend(:,:,idtend) = dtend(:,:,idtend) + field*dtf
1045 endif
1046 endif
1047 END SUBROUTINE dtend_helper
1048
1049! ==================================================================
1050 SUBROUTINE moisture_check2(kte, delt, dp, exner, &
1051 qv, qc, qi, qs, th )
1052 !
1053 ! If qc < qcmin, qi < qimin, or qv < qvmin happens in any layer,
1054 ! force them to be larger than minimum value by (1) condensating
1055 ! water vapor into liquid or ice, and (2) by transporting water vapor
1056 ! from the very lower layer.
1057 !
1058 ! We then update the final state variables and tendencies associated
1059 ! with this correction. If any condensation happens, update theta/temperature too.
1060 ! Note that (qv,qc,qi,th) are the final state variables after
1061 ! applying corresponding input tendencies and corrective tendencies.
1062
1063 implicit none
1064 integer, intent(in) :: kte
1065 real(kind_phys), intent(in) :: delt
1066 real(kind_phys), dimension(kte), intent(in) :: dp, exner
1067 real(kind_phys), dimension(kte), intent(inout) :: qv, qc, qi, qs, th
1068 integer k
1069 real :: dqc2, dqi2, dqs2, dqv2, sum, aa, dum
1070 real, parameter :: qvmin1= 1e-8, & !min at k=1
1071 qvmin = 1e-20, & !min above k=1
1072 qcmin = 0.0, &
1073 qimin = 0.0
1074
1075 do k = kte, 1, -1 ! From the top to the surface
1076 dqc2 = max(0.0, qcmin-qc(k)) !qc deficit (>=0)
1077 dqi2 = max(0.0, qimin-qi(k)) !qi deficit (>=0)
1078 dqs2 = max(0.0, qimin-qs(k)) !qs deficit (>=0)
1079
1080 !update species
1081 qc(k) = qc(k) + dqc2
1082 qi(k) = qi(k) + dqi2
1083 qs(k) = qs(k) + dqs2
1084 qv(k) = qv(k) - dqc2 - dqi2 - dqs2
1085 !for theta
1086 !th(k) = th(k) + xlvcp/exner(k)*dqc2 + &
1087 ! xlscp/exner(k)*dqi2
1088 !for temperature
1089 th(k) = th(k) + xlvcp*dqc2 + &
1090 xlscp*(dqi2+dqs2)
1091
1092 !then fix qv if lending qv made it negative
1093 if (k .eq. 1) then
1094 dqv2 = max(0.0, qvmin1-qv(k)) !qv deficit (>=0)
1095 qv(k) = qv(k) + dqv2
1096 qv(k) = max(qv(k),qvmin1)
1097 dqv2 = 0.0
1098 else
1099 dqv2 = max(0.0, qvmin-qv(k)) !qv deficit (>=0)
1100 qv(k) = qv(k) + dqv2
1101 qv(k-1)= qv(k-1) - dqv2*dp(k)/dp(k-1)
1102 qv(k) = max(qv(k),qvmin)
1103 endif
1104 qc(k) = max(qc(k),qcmin)
1105 qi(k) = max(qi(k),qimin)
1106 qs(k) = max(qs(k),qimin)
1107 end do
1108
1109 ! Extra moisture used to satisfy 'qv(1)>=qvmin' is proportionally
1110 ! extracted from all the layers that has 'qv > 2*qvmin'. This fully
1111 ! preserves column moisture.
1112 if( dqv2 .gt. 1.e-20 ) then
1113 sum = 0.0
1114 do k = 1, kte
1115 if( qv(k) .gt. 2.0*qvmin ) sum = sum + qv(k)*dp(k)
1116 enddo
1117 aa = dqv2*dp(1)/max(1.e-20,sum)
1118 if( aa .lt. 0.5 ) then
1119 do k = 1, kte
1120 if( qv(k) .gt. 2.0*qvmin ) then
1121 dum = aa*qv(k)
1122 qv(k) = qv(k) - dum
1123 endif
1124 enddo
1125 else
1126 ! For testing purposes only (not yet found in any output):
1127 ! write(*,*) 'Full moisture conservation is impossible'
1128 endif
1129 endif
1130
1131 return
1132
1133 END SUBROUTINE moisture_check2
1134
1135 END SUBROUTINE mynnedmf_wrapper_run
1136
1137!###=================================================================
1138
1139END MODULE mynnedmf_wrapper
subroutine mynn_bl_driver(initflag, restart, cycling, delt, dz, dx, znt, u, v, w, th, sqv3d, sqc3d, sqi3d, sqs3d, qnc, qni, qnwfa, qnifa, qnbca, ozone, p, exner, rho, t3d, xland, ts, qsfc, ps, ust, ch, hfx, qfx, rmol, wspd, uoce, voce, qke, qke_adv, sh3d, sm3d, nchem, kdvel, ndvel, chem3d, vdep, smoke_dbg, frp, emis_ant_no, mix_chem, enh_mix, rrfs_sd, tsq, qsq, cov, rublten, rvblten, rthblten, rqvblten, rqcblten, rqiblten, rqncblten, rqniblten, rqsblten, rqnwfablten, rqnifablten, rqnbcablten, dozone, exch_h, exch_m, pblh, kpbl, el_pbl, dqke, qwt, qshear, qbuoy, qdiss, qc_bl, qi_bl, cldfra_bl, bl_mynn_tkeadvect, tke_budget, bl_mynn_cloudpdf, bl_mynn_mixlength, icloud_bl, closure, bl_mynn_edmf, bl_mynn_edmf_mom, bl_mynn_edmf_tke, bl_mynn_mixscalars, bl_mynn_output, bl_mynn_cloudmix, bl_mynn_mixqt, edmf_a, edmf_w, edmf_qt, edmf_thl, edmf_ent, edmf_qc, sub_thl3d, sub_sqv3d, det_thl3d, det_sqv3d, maxwidth, maxmf, ztop_plume, ktop_plume, spp_pbl, pattern_spp_pbl, rthraten, flag_qc, flag_qi, flag_qnc, flag_qni, flag_qs, flag_qnwfa, flag_qnifa, flag_qnbca, flag_ozone, ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
This subroutine is the MYNN-EDNF PBL driver routine,which encompassed the majority of the subroutines...
This module defines model-specific constants/parameters.
This module contains the entity of MYNN-EDMF PBL scheme.
The following references best describe the code within Olson et al. (2019, NOAA Technical Memorandum)...