CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
cu_c3_driver.F90
1
3
4
6
7 ! DH* TODO: replace constants with arguments to cu_c3_driver_run
8 !use physcons , g => con_g, cp => con_cp, xlv => con_hvap, r_v => con_rv
9 use machine , only: kind_phys
11 use cu_c3_sh , only: cu_c3_sh_run
12 use progsigma , only: progsigma_calc
13
14 implicit none
15
16 private
17
18 public :: cu_c3_driver_init, cu_c3_driver_run, progsigma_calc
19
20contains
21
32 subroutine cu_c3_driver_init(imfshalcnv, imfshalcnv_c3, imfdeepcnv, &
33 imfdeepcnv_c3,progsigma, cnx, mpirank, mpiroot, &
34 errmsg, errflg)
35
36 implicit none
37
38 integer, intent(in) :: imfshalcnv, imfshalcnv_c3
39 integer, intent(in) :: imfdeepcnv, imfdeepcnv_c3
40 integer, intent(in) :: mpirank
41 integer, intent(in) :: mpiroot
42 integer, intent(in) :: cnx
43 logical, intent(inout) :: progsigma
44 character(len=*), intent( out) :: errmsg
45 integer, intent( out) :: errflg
46
47 ! initialize ccpp error handling variables
48 errmsg = ''
49 errflg = 0
50
51 if(progsigma)then
52 if(cnx < 384)then
53 progsigma=.false.
54 write(*,*)'Forcing prognostic closure to .false. due to coarse resolution'
55 endif
56 endif
57
58 end subroutine cu_c3_driver_init
59
60!
61! t2di is temp after advection, but before physics
62! t = current temp (t2di + physics up to now)
63!===================
64
70 subroutine cu_c3_driver_run(ntracer,garea,im,km,dt,flag_init,flag_restart,&
71 do_ca,progsigma,cactiv,cactiv_m,g,cp,fv,r_d,xlv,r_v,forcet, &
72 forceqv_spechum,phil,delp,raincv,tmf,qmicro,sigmain, &
73 betascu,betamcu,betadcu,qv_spechum,t,cld1d,us,vs,t2di,w, &
74 qv2di_spechum,p2di,psuri, &
75 hbot,htop,kcnv,xland,hfx2,qfx2,aod_gf,cliw,clcw,ca_deep,rainevap,&
76 pbl,ud_mf,dd_mf,dt_mf,cnvw_moist,cnvc,imfshalcnv, &
77 flag_for_scnv_generic_tend,flag_for_dcnv_generic_tend, &
78 dtend,dtidx,ntqv,ntiw,ntcw,index_of_temperature,index_of_x_wind, &
79 index_of_y_wind,index_of_process_scnv,index_of_process_dcnv, &
80 fhour,fh_dfi_radar,ix_dfi_radar,num_dfi_radar,cap_suppress, &
81 dfi_radar_max_intervals,ldiag3d,qci_conv,do_cap_suppress, &
82 sigmaout,maxupmf,maxMF,do_mynnedmf,ichoice_in,ichoicem_in, &
83 ichoice_s_in,errmsg,errflg)
84!-------------------------------------------------------------
85 implicit none
86 integer, parameter :: maxiens=1
87 integer, parameter :: maxens=1
88 integer, parameter :: maxens2=1
89 integer, parameter :: maxens3=16
90 integer, parameter :: ensdim=16
91 integer :: imid_gf=1 ! gf congest conv.
92 integer, parameter :: ideep=1
93 integer :: ichoice=0 ! 0 2 5 13 8
94 integer :: ichoicem=13 ! 0 2 5 13
95 integer :: ichoice_s=3 ! 0 1 2 3
96
97 logical, intent(in) :: do_cap_suppress
98 real(kind=kind_phys), parameter :: aodc0=0.14
99 real(kind=kind_phys), parameter :: aodreturn=30.
100 real(kind=kind_phys) :: dts,fpi,fp
101 integer, parameter :: dicycle=0 ! diurnal cycle flag
102 integer, parameter :: dicycle_m=0 !- diurnal cycle flag
103 integer :: ishallow_g3 ! depend on imfshalcnv
104!-------------------------------------------------------------
105 integer :: its,ite, jts,jte, kts,kte
106 integer, intent(in ) :: im,km,ntracer
107 integer, intent(in ) :: ichoice_in,ichoicem_in,ichoice_s_in
108 logical, intent(in ) :: flag_init, flag_restart, do_mynnedmf
109 logical, intent(in ) :: flag_for_scnv_generic_tend,flag_for_dcnv_generic_tend, &
110 do_ca
111 real (kind=kind_phys), intent(in) :: g,cp,fv,r_d,xlv,r_v,betascu,betamcu,betadcu
112 logical, intent(in ) :: ldiag3d
113 logical, intent(in ) :: progsigma
114 real(kind=kind_phys), intent(inout), optional :: dtend(:,:,:)
115!$acc declare copy(dtend)
116 integer, intent(in) :: dtidx(:,:), &
117 index_of_x_wind, index_of_y_wind, index_of_temperature, &
118 index_of_process_scnv, index_of_process_dcnv, ntqv, ntcw, ntiw
119!$acc declare copyin(dtidx)
120 real(kind=kind_phys), dimension( : , : ), intent(in ), optional :: forcet,forceqv_spechum
121 real(kind=kind_phys), dimension( : , : ), intent(in ) :: w,phil,delp
122 real(kind=kind_phys), dimension ( : , : ), intent(in ), optional :: sigmain,qmicro
123 real(kind=kind_phys), dimension( : , : ), intent(inout ) :: t,us,vs
124 real(kind=kind_phys), dimension( : , : ), intent(inout ), optional :: qci_conv
125 real(kind=kind_phys), dimension( : , : ), intent(out ) :: cnvw_moist,cnvc
126 real(kind=kind_phys), dimension ( : , : ), intent(out ), optional :: sigmaout
127 real(kind=kind_phys), dimension( : , : ), intent(inout ) :: cliw, clcw
128 real(kind=kind_phys), dimension ( : , : , :), intent(in ) :: tmf
129!$acc declare copyin(forcet,forceqv_spechum,w,phil)
130!$acc declare copy(t,us,vs,qci_conv,cliw, clcw)
131!$acc declare copyout(cnvw_moist,cnvc)
132
133 real(kind=kind_phys), allocatable :: clcw_save(:,:), cliw_save(:,:)
134
135 integer, intent(in) :: dfi_radar_max_intervals
136 real(kind=kind_phys), intent(in) :: fhour, fh_dfi_radar(:)
137 integer, intent(in) :: num_dfi_radar, ix_dfi_radar(:)
138 real(kind=kind_phys), intent(in), optional :: cap_suppress(:,:)
139!$acc declare copyin(fh_dfi_radar,ix_dfi_radar,cap_suppress)
140
141 integer, dimension (:), intent(out) :: hbot,htop,kcnv
142 integer, dimension (:), intent(in) :: xland
143 real(kind=kind_phys), dimension (:), intent(in) :: pbl
144 real(kind=kind_phys), dimension (:), intent(in), optional :: maxmf
145!$acc declare copyout(hbot,htop,kcnv)
146!$acc declare copyin(xland,pbl)
147 integer, dimension (im) :: tropics
148!$acc declare create(tropics)
149! ruc variable
150 real(kind=kind_phys), dimension (:), intent(in) :: hfx2,qfx2,psuri
151 real(kind=kind_phys), dimension (:), intent(in), optional :: ca_deep
152 real(kind=kind_phys), dimension (:,:), intent(out), optional :: ud_mf
153 real(kind=kind_phys), dimension (:,:), intent(out) :: dd_mf,dt_mf
154 real(kind=kind_phys), dimension (:), intent(out) :: raincv,cld1d,rainevap
155 real(kind=kind_phys), dimension (:), intent(out), optional :: maxupmf
156 real(kind=kind_phys), dimension (:,:), intent(in) :: t2di,p2di
157!$acc declare copyin(hfx2,qfx2,psuri,t2di,p2di)
158!$acc declare copyout(ud_mf,dd_mf,dt_mf,raincv,cld1d)
159 ! Specific humidity from FV3
160 real(kind=kind_phys), dimension (:,:), intent(in) :: qv2di_spechum
161 real(kind=kind_phys), dimension (:,:), intent(inout) :: qv_spechum
162 real(kind=kind_phys), dimension (:), intent(inout), optional :: aod_gf
163!$acc declare copyin(qv2di_spechum) copy(qv_spechum,aod_gf)
164 ! Local water vapor mixing ratios and cloud water mixing ratios
165 real(kind=kind_phys), dimension (im,km) :: qv2di, qv, forceqv, cnvw
166!$acc declare create(qv2di, qv, forceqv, cnvw)
167 !
168 real(kind=kind_phys), dimension(:),intent(in) :: garea
169!$acc declare copyin(garea)
170 real(kind=kind_phys), intent(in ) :: dt
171
172 integer, intent(in ) :: imfshalcnv
173 integer, dimension(:), intent(inout), optional :: cactiv,cactiv_m
174!$acc declare copy(cactiv,cactiv_m)
175
176 character(len=*), intent(out) :: errmsg
177 integer, intent(out) :: errflg
178
179! local variables
180 integer, dimension(im) :: k22_shallow,kbcon_shallow,ktop_shallow
181 real(kind=kind_phys), dimension (im) :: rand_mom,rand_vmas
182 real(kind=kind_phys), dimension (im,4) :: rand_clos
183 real(kind=kind_phys), dimension (im,km,11) :: gdc,gdc2
184 real(kind=kind_phys), dimension (im) :: ht
185 real(kind=kind_phys), dimension (im) :: ccn_gf,ccn_m
186 real(kind=kind_phys) :: ccnclean
187 real(kind=kind_phys), dimension (im) :: dx
188 real(kind=kind_phys), dimension (im) :: frhm,frhd
189 real(kind=kind_phys), dimension (im,km) :: outt,outq,outqc,phh,subm,cupclw,cupclws
190 real(kind=kind_phys), dimension (im,km) :: dhdt,zu,zus,zd,phf,zum,zdm,outum,outvm
191 real(kind=kind_phys), dimension (im,km) :: outts,outqs,outqcs,outu,outv,outus,outvs
192 real(kind=kind_phys), dimension (im,km) :: outtm,outqm,outqcm,submm,cupclwm
193 real(kind=kind_phys), dimension (im,km) :: cnvwt,cnvwts,cnvwtm
194 real(kind=kind_phys), dimension (im,km) :: hco,hcdo,zdo,zdd,hcom,hcdom,zdom,tmfq
195 real(kind=kind_phys), dimension (km) :: zh
196 real(kind=kind_phys), dimension (im) :: tau_ecmwf,edt,edtm,edtd,ter11,aa0,xlandi
197 real(kind=kind_phys), dimension (im) :: pret,prets,pretm,hexec
198 real(kind=kind_phys), dimension (im,10) :: forcing,forcing2
199
200 integer, dimension (im) :: kbcon, ktop,ierr,ierrs,ierrm,kpbli
201 integer, dimension (im) :: k22s,kbcons,ktops,k22,jmin,jminm
202 integer, dimension (im) :: kbconm,ktopm,k22m
203!$acc declare create(k22_shallow,kbcon_shallow,ktop_shallow,rand_mom,rand_vmas, &
204!$acc rand_clos,gdc,gdc2,ht,ccn_gf,ccn_m,dx,frhm,frhd, &
205!$acc outt,outq,outqc,phh,subm,cupclw,cupclws, &
206!$acc dhdt,zu,zus,zd,phf,zum,zdm,outum,outvm, &
207!$acc outts,outqs,outqcs,outu,outv,outus,outvs, &
208!$acc outtm,outqm,outqcm,submm,cupclwm, &
209!$acc cnvwt,cnvwts,cnvwtm,hco,hcdo,zdo,zdd,hcom,hcdom,zdom, &
210!$acc tau_ecmwf,edt,edtm,edtd,ter11,aa0,xlandi, &
211!$acc pret,prets,pretm,hexec,forcing,forcing2, &
212!$acc kbcon, ktop,ierr,ierrs,ierrm,kpbli, &
213!$acc k22s,kbcons,ktops,k22,jmin,jminm,kbconm,ktopm,k22m)
214
215 integer :: iens,ibeg,iend,jbeg,jend,n
216 integer :: ibegh,iendh,jbegh,jendh
217 integer :: ibegc,iendc,jbegc,jendc,kstop
218 real(kind=kind_phys), dimension(im,km) :: rho_dryar
219!$acc declare create(rho_dryar)
220 real(kind=kind_phys) :: pten,pqen,paph,zrho,pahfs,pqhfl,zkhvfl,pgeoh
221 integer, parameter :: ipn = 0
222
223!
224! basic environmental input includes moisture convergence (mconv)
225! omega (omeg), windspeed (us,vs), and a flag (ierr) to turn off
226! convection for this call only and at that particular gridpoint
227!
228 real(kind=kind_phys), dimension (im,km) :: qcheck,zo,t2d,q2d,po,p2d,rhoi,clw_ten
229 real(kind=kind_phys), dimension (im,km) :: tn,qo,tshall,qshall,dz8w,omeg
230 real(kind=kind_phys), dimension (im) :: z1,psur,cuten,cutens,cutenm
231 real(kind=kind_phys), dimension (im) :: umean,vmean,pmean
232 real(kind=kind_phys), dimension (im) :: xmbs,xmbs2,xmb,xmbm,xmb_dumm,mconv
233!$acc declare create(qcheck,zo,t2d,q2d,po,p2d,rhoi,clw_ten,tn,qo,tshall,qshall,dz8w,omeg, &
234!$acc z1,psur,cuten,cutens,cutenm,umean,vmean,pmean, &
235!$acc xmbs,xmbs2,xmb,xmbm,xmb_dumm,mconv)
236
237 integer :: i,j,k,icldck,ipr,jpr,jpr_deep,ipr_deep,uidx,vidx,tidx,qidx
238 integer :: itf,jtf,ktf,iss,jss,nbegin,nend,cliw_idx,clcw_idx
239 integer :: high_resolution
240 real(kind=kind_phys) :: clwtot,clwtot1,excess,tcrit,tscl_kf,dp,dq,sub_spread,subcenter
241 real(kind=kind_phys) :: dsubclw,dsubclws,dsubclwm,dtime_max,ztm,ztq,hfm,qfm,rkbcon,rktop
242 real(kind=kind_phys), dimension(km) :: massflx,trcflx_in1,clw_in1,po_cup
243! real(kind=kind_phys), dimension(km) :: trcflx_in2,clw_in2,clw_ten2
244 real(kind=kind_phys), dimension (im) :: flux_tun,tun_rad_mid,tun_rad_shall,tun_rad_deep
245!$acc declare create(flux_tun,tun_rad_mid,tun_rad_shall,tun_rad_deep)
246 character*50 :: ierrc(im),ierrcm(im)
247 character*50 :: ierrcs(im)
248! ruc variable
249! hfx2 -- sensible heat flux (k m/s), positive upward from sfc
250! qfx2 -- latent heat flux (kg/kg m/s), positive upward from sfc
251! gf needs them in w/m2. define hfx and qfx after simple unit conversion
252 real(kind=kind_phys), dimension (im) :: hfx,qfx
253!$acc declare create(hfx,qfx)
254 real(kind=kind_phys) tem,tem1,tf,tcr,tcrf,psum
255 real(kind=kind_phys) :: cliw_shal,clcw_shal,tem_shal, cliw_both, weight_sum
256 real(kind=kind_phys) :: cliw_deep,clcw_deep,tem_deep, clcw_both
257 integer :: cliw_deep_idx, clcw_deep_idx, cliw_shal_idx, clcw_shal_idx
258
259 real(kind=kind_phys) :: cap_suppress_j(im)
260!$acc declare create(cap_suppress_j)
261 integer :: itime, do_cap_suppress_here
262 logical :: exit_func
263
264 !parameter (tf=243.16, tcr=270.16, tcrf=1.0/(tcr-tf)) ! FV3 original
265 !parameter (tf=263.16, tcr=273.16, tcrf=1.0/(tcr-tf))
266 !parameter (tf=233.16, tcr=263.16, tcrf=1.0/(tcr-tf))
267 parameter(tf=258.16, tcr=273.16, tcrf=1.0/(tcr-tf)) ! as fim, HCB tuning
268 ! initialize ccpp error handling variables
269 errmsg = ''
270 errflg = 0
271
272 ichoice = ichoice_in
273 ichoicem = ichoicem_in
274 ichoice_s = ichoice_s_in
275
276 if(do_cap_suppress) then
277!$acc serial
278 do itime=1,num_dfi_radar
279 if(ix_dfi_radar(itime)<1) cycle
280 if(fhour<fh_dfi_radar(itime)) cycle
281 if(fhour>=fh_dfi_radar(itime+1)) cycle
282 exit
283 enddo
284!$acc end serial
285 endif
286 if(do_cap_suppress .and. itime<=num_dfi_radar) then
287 do_cap_suppress_here = 1
288!$acc kernels
289 cap_suppress_j(:) = cap_suppress(:,itime)
290!$acc end kernels
291 else
292 do_cap_suppress_here = 0
293!$acc kernels
294 cap_suppress_j(:) = 0
295!$acc end kernels
296 endif
297
298 if(ldiag3d) then
299 if(flag_for_dcnv_generic_tend) then
300 cliw_deep_idx=0
301 clcw_deep_idx=0
302 else
303 cliw_deep_idx=dtidx(100+ntiw,index_of_process_dcnv)
304 clcw_deep_idx=dtidx(100+ntcw,index_of_process_dcnv)
305 endif
306 if(flag_for_scnv_generic_tend) then
307 cliw_shal_idx=0
308 clcw_shal_idx=0
309 else
310 cliw_shal_idx=dtidx(100+ntiw,index_of_process_scnv)
311 clcw_shal_idx=dtidx(100+ntcw,index_of_process_scnv)
312 endif
313 if(cliw_deep_idx>=1 .or. clcw_deep_idx>=1 .or. &
314 cliw_shal_idx>=1 .or. clcw_shal_idx>=1) then
315 allocate(clcw_save(im,km), cliw_save(im,km))
316!$acc enter data create(clcw_save,cliw_save)
317!$acc kernels
318 clcw_save(:,:)=clcw(:,:)
319 cliw_save(:,:)=cliw(:,:)
320!$acc end kernels
321 endif
322 endif
323
324!
325! Scale specific humidity to dry mixing ratio
326!
327!$acc kernels
328 ! state in before physics
329 qv2di = qv2di_spechum/(1.0_kind_phys-qv2di_spechum)
330 ! forcing by dynamics, based on state in
331 forceqv = forceqv_spechum/(1.0_kind_phys-qv2di_spechum)
332 ! current state (updated by preceeding physics)
333 qv = qv_spechum/(1.0_kind_phys-qv_spechum)
334!
335!
336! these should be coming in from outside
337!
338! cactiv(:) = 0
339 rand_mom(:) = 0.
340 rand_vmas(:) = 0.
341 rand_clos(:,:) = 0.
342!$acc end kernels
343!
344 its=1
345 ite=im
346 itf=ite
347 jts=1
348 jte=1
349 jtf=jte
350 kts=1
351 kte=km
352 ktf=kte-1
353!$acc kernels
354!
355 tropics(:)=0
356!
358!
359 tun_rad_shall(:)=.012
360 tun_rad_mid(:)=.15 !.02
361 tun_rad_deep(:)=.3 !.065
362 edt(:)=0.
363 edtm(:)=0.
364 edtd(:)=0.
365 zdd(:,:)=0.
366 flux_tun(:)=5.
367! dx for scale awareness
368!$acc end kernels
369
370 if (imfshalcnv == 5) then
371 ishallow_g3 = 1
372 else
373 ishallow_g3 = 0
374 end if
375 high_resolution=0
376 subcenter=0.
377 iens=1
378!
379! these can be set for debugging
380!
381 ipr=0
382 jpr=0
383 ipr_deep=0
384 jpr_deep= 0 !53322 ! 528196 !0 ! 1136 !0 !421755 !3536
385!
386!
387 ibeg=its
388 iend=ite
389 tcrit=258.
390
391 ztm=0.
392 ztq=0.
393 hfm=0.
394 qfm=0.
395!$acc kernels
396 ud_mf(:,:) =0.
397 dd_mf(:,:) =0.
398 dt_mf(:,:) =0.
399 tau_ecmwf(:)=0.
400!$acc end kernels
401!
402 j=1
403!$acc kernels
404 ht(:)=phil(:,1)/g
405!$acc loop private(zh)
406 do i=its,ite
407 cld1d(i)=0.
408 zo(i,:)=phil(i,:)/g
409 dz8w(i,1)=zo(i,2)-zo(i,1)
410 zh(1)=0.
411 kpbli(i)=2
412 do k=kts+1,ktf
413 dz8w(i,k)=zo(i,k+1)-zo(i,k)
414 enddo
415!$acc loop seq
416 do k=kts+1,ktf
417 zh(k)=zh(k-1)+dz8w(i,k-1)
418 if(zh(k).gt.pbl(i))then
419 kpbli(i)=max(2,k)
420 exit
421 endif
422 enddo
423 enddo
424!$acc end kernels
425
426!$acc kernels
427 do i= its,itf
428 forcing(i,:)=0.
429 forcing2(i,:)=0.
430 ccn_gf(i) = 0.
431 ccn_m(i) = 0.
432
433 ! set aod and ccn
434 if (flag_init .and. .not.flag_restart) then
435 aod_gf(i)=aodc0
436 else
437 if((cactiv(i).eq.0) .and. (cactiv_m(i).eq.0))then
438 if(aodc0>aod_gf(i)) aod_gf(i)=aod_gf(i)+((aodc0-aod_gf(i))*(dt/(aodreturn*60)))
439 if(aod_gf(i)>aodc0) aod_gf(i)=aodc0
440 endif
441 endif
442
443 ccn_gf(i)=max(5., (aod_gf(i)/0.0027)**(1/0.640))
444 ccn_m(i)=ccn_gf(i)
445
446 ccnclean=max(5., (aodc0/0.0027)**(1/0.640))
447
448 hbot(i) =kte
449 htop(i) =kts
450 raincv(i)=0.
451 xlandi(i)=real(xland(i))
452! if(abs(xlandi(i)-1.).le.1.e-3) tun_rad_shall(i)=.15
453! if(abs(xlandi(i)-1.).le.1.e-3) flux_tun(i)=1.5
454 enddo
455 do i= its,itf
456 mconv(i)=0.
457 enddo
458 do k=kts,kte
459 do i= its,itf
460 omeg(i,k)=0.
461 zu(i,k)=0.
462 zum(i,k)=0.
463 zus(i,k)=0.
464 zd(i,k)=0.
465 zdm(i,k)=0.
466 enddo
467 enddo
468
469 psur(:)=0.01*psuri(:)
470 do i=its,itf
471 ter11(i)=max(0.,ht(i))
472 enddo
473 do k=kts,kte
474 do i=its,ite
475 cnvw(i,k)=0.
476 cnvc(i,k)=0.
477 gdc(i,k,1)=0.
478 gdc(i,k,2)=0.
479 gdc(i,k,3)=0.
480 gdc(i,k,4)=0.
481 gdc(i,k,7)=0.
482 gdc(i,k,8)=0.
483 gdc(i,k,9)=0.
484 gdc(i,k,10)=0.
485 gdc2(i,k,1)=0.
486 enddo
487 enddo
488
489 do k=kts,kte
490 do i=its,ite
491 tmfq(i,k)=tmf(i,k,1)
492 enddo
493 enddo
494
495 ierr(:)=0
496 ierrm(:)=0
497 ierrs(:)=0
498 cuten(:)=0.
499 cutenm(:)=0.
500 cutens(:)=0.
501!$acc end kernels
502 ierrc(:)=" "
503!$acc kernels
504
505
506 kbcon(:)=0
507 kbcons(:)=0
508 kbconm(:)=0
509
510 ktop(:)=0
511 ktops(:)=0
512 ktopm(:)=0
513
514 xmb(:)=0.
515 xmb_dumm(:)=0.
516 xmbm(:)=0.
517 xmbs(:)=0.
518 xmbs2(:)=0.
519
520 k22s(:)=0
521 k22m(:)=0
522 k22(:)=0
523
524 jmin(:)=0
525 jminm(:)=0
526
527 pret(:)=0.
528 prets(:)=0.
529 pretm(:)=0.
530
531 umean(:)=0.
532 vmean(:)=0.
533 pmean(:)=0.
534
535 cupclw(:,:)=0.
536 cupclwm(:,:)=0.
537 cupclws(:,:)=0.
538
539 cnvwt(:,:)=0.
540 cnvwts(:,:)=0.
541 cnvwtm(:,:)=0.
542
543 hco(:,:)=0.
544 hcom(:,:)=0.
545 hcdo(:,:)=0.
546 hcdom(:,:)=0.
547
548 outt(:,:)=0.
549 outts(:,:)=0.
550 outtm(:,:)=0.
551
552 outu(:,:)=0.
553 outus(:,:)=0.
554 outum(:,:)=0.
555
556 outv(:,:)=0.
557 outvs(:,:)=0.
558 outvm(:,:)=0.
559
560 outq(:,:)=0.
561 outqs(:,:)=0.
562 outqm(:,:)=0.
563
564 outqc(:,:)=0.
565 outqcs(:,:)=0.
566 outqcm(:,:)=0.
567
568 subm(:,:)=0.
569 dhdt(:,:)=0.
570
571 frhm(:)=0.
572 frhd(:)=0.
573
574 do k=kts,ktf
575 do i=its,itf
576 p2d(i,k)=0.01*p2di(i,k)
577 po(i,k)=p2d(i,k) !*.01
578 rhoi(i,k) = 100.*p2d(i,k)/(287.04*(t2di(i,k)*(1.+0.608*qv2di(i,k))))
579 qcheck(i,k)=qv(i,k)
580 tn(i,k)=t(i,k)!+forcet(i,k)*dt
581 qo(i,k)=max(1.e-16,qv(i,k))!+forceqv(i,k)*dt
582 t2d(i,k)=t2di(i,k)-forcet(i,k)*dt
583 q2d(i,k)=max(1.e-16,qv2di(i,k)-forceqv(i,k)*dt)
584 if(qo(i,k).lt.1.e-16)qo(i,k)=1.e-16
585 tshall(i,k)=t2d(i,k)
586 qshall(i,k)=q2d(i,k)
587 enddo
588 enddo
589!$acc end kernels
590123 format(1x,i2,1x,2(1x,f8.0),1x,2(1x,f8.3),3(1x,e13.5))
591!$acc kernels
592 do i=its,itf
593 do k=kts,kpbli(i)
594 tshall(i,k)=t(i,k)
595 qshall(i,k)=max(1.e-16,qv(i,k))
596 enddo
597 enddo
598!
599! converting hfx2 and qfx2 to w/m2
600! hfx=cp*rho*hfx2
601! qfx=xlv*qfx2
602 do i=its,itf
603 hfx(i)=hfx2(i)*cp*rhoi(i,1)
604 qfx(i)=qfx2(i)*xlv*rhoi(i,1)
605 dx(i) = sqrt(garea(i))
606 enddo
607
608 do i=its,itf
609 do k=kts,kpbli(i)
610 tn(i,k)=t(i,k)
611 qo(i,k)=max(1.e-16,qv(i,k))
612 enddo
613 enddo
614 nbegin=0
615 nend=0
616 do i=its,itf
617 do k=kts,kpbli(i)
618 dhdt(i,k)=cp*(forcet(i,k)+(t(i,k)-t2di(i,k))/dt) + &
619 xlv*(forceqv(i,k)+(qv(i,k)-qv2di(i,k))/dt)
620! tshall(i,k)=t(i,k)
621! qshall(i,k)=qv(i,k)
622 enddo
623 enddo
624!$acc loop collapse(2) independent private(dp)
625 do k= kts+1,ktf-1
626 do i = its,itf
627 if((p2d(i,1)-p2d(i,k)).gt.150.and.p2d(i,k).gt.300)then
628 dp=-.5*(p2d(i,k+1)-p2d(i,k-1))
629!$acc atomic
630 umean(i)=umean(i)+us(i,k)*dp
631!$acc atomic
632 vmean(i)=vmean(i)+vs(i,k)*dp
633!$acc atomic
634 pmean(i)=pmean(i)+dp
635 endif
636 enddo
637 enddo
638 do i = its,itf
639 psum=0.
640 do k=kts,ktf-3
641 if (clcw(i,k) .gt. -999.0 .and. clcw(i,k+1) .gt. -999.0 )then
642 dp=(p2d(i,k)-p2d(i,k+1))
643 psum=psum+dp
644 clwtot = cliw(i,k) + clcw(i,k)
645 if(clwtot.lt.1.e-32)clwtot=0.
646 forcing(i,7)=forcing(i,7)+clwtot*dp
647 endif
648 enddo
649 if(psum.gt.0)forcing(i,7)=forcing(i,7)/psum
650 forcing2(i,7)=forcing(i,7)
651 enddo
652 do k=kts,ktf-1
653 do i = its,itf
654 omeg(i,k)= w(i,k) !-g*rhoi(i,k)*w(i,k)
655 enddo
656 enddo
657 do i = its,itf
658 if(mconv(i).lt.0.)mconv(i)=0.
659 if((dx(i)<6500.).and.do_mynnedmf.and.(maxmf(i).gt.0.))ierr(i)=555
660 enddo
661!$acc end kernels
662 if (dx(its)<6500.) then
663 imid_gf=0
664 endif
665!
666!---- call cumulus parameterization
667!
668 if(ishallow_g3.eq.1)then
669
670!$acc kernels
671 do i=its,ite
672 ierrs(i)=0
673 ierrm(i)=0
674 enddo
675!$acc end kernels
676!
678!
679 call cu_c3_sh_run (us,vs, &
680! input variables, must be supplied
681 zo,t2d,q2d,ter11,tshall,qshall,p2d,psur,dhdt,kpbli, &
682 rhoi,hfx,qfx,xlandi,ichoice_s,tcrit,dt, &
683! input variables. ierr should be initialized to zero or larger than zero for
684! turning off shallow convection for grid points
685 zus,xmbs,kbcons,ktops,k22s,ierrs,ierrcs, &
686! Prog closure
687 flag_init, flag_restart,fv,r_d,delp,tmfq,qmicro, &
688 forceqv_spechum,betascu,betamcu,betadcu,sigmain, &
689 sigmaout,progsigma,dx, &
690! output tendencies
691 outts,outqs,outqcs,outus,outvs,cnvwt,prets,cupclws, &
692! dimesnional variables
693 itf,ktf,its,ite, kts,kte,ipr,tropics)
694
695!$acc kernels
696 do i=its,itf
697 if(xmbs(i).gt.0.)then
698 cutens(i)=1.
699 endif
700 enddo
701!$acc end kernels
703 call neg_check('shallow',ipn,dt,qcheck,outqs,outts,outus,outvs, &
704 outqcs,prets,its,ite,kts,kte,itf,ktf,ktops)
705 endif
706
707 ipr=0
708 jpr_deep=0 !340765
710 if(imid_gf == 1)then
711 call cu_c3_deep_run( &
712 itf,ktf,its,ite, kts,kte &
713 ,flag_init &
714 ,flag_restart &
715 ,fv,r_d &
716 ,dicycle_m &
717 ,ichoicem &
718 ,ipr &
719 ,ccn_m &
720 ,ccnclean &
721 ,dt &
722 ,imid_gf &
723 ,kpbli &
724 ,dhdt &
725 ,xlandi &
726 ,delp &
727 ,zo &
728 ,forcing &
729 ,t2d &
730 ,q2d &
731 ,tmfq &
732 ,qmicro &
733 ,forceqv_spechum &
734 ,betascu &
735 ,betamcu &
736 ,betadcu &
737 ,sigmain &
738 ,sigmaout &
739 ,ter11 &
740 ,tshall &
741 ,qshall &
742 ,p2d &
743 ,psur &
744 ,us &
745 ,vs &
746 ,rhoi &
747 ,hfx &
748 ,qfx &
749 ,dx &
750 ,do_ca &
751 ,progsigma &
752 ,ca_deep &
753 ,mconv &
754 ,omeg &
755 ,cactiv_m &
756 ,cnvwtm &
757 ,zum &
758 ,zdm & ! hli
759 ,zdd &
760 ,edtm &
761 ,edtd & ! hli
762 ,xmbm &
763 ,xmb_dumm &
764 ,xmbs &
765 ,pretm &
766 ,outum &
767 ,outvm &
768 ,outtm &
769 ,outqm &
770 ,outqcm &
771 ,kbconm &
772 ,ktopm &
773 ,cupclwm &
774 ,frhm &
775 ,rainevap &
776 ,ierrm &
777 ,ierrcm &
778! the following should be set to zero if not available
779 ,rand_mom & ! for stochastics mom, if temporal and spatial patterns exist
780 ,rand_vmas & ! for stochastics vertmass, if temporal and spatial patterns exist
781 ,rand_clos & ! for stochastics closures, if temporal and spatial patterns exist
782 ,0 & ! flag to what you want perturbed
783 ! 1 = momentum transport
784 ! 2 = normalized vertical mass flux profile
785 ! 3 = closures
786 ! more is possible, talk to developer or
787 ! implement yourself. pattern is expected to be
788 ! betwee -1 and +1
789 ,do_cap_suppress_here,cap_suppress_j &
790 ,k22m &
791 ,jminm,tropics)
792!$acc kernels
793 do i=its,itf
794 do k=kts,ktf
795 qcheck(i,k)=qv(i,k) +outqs(i,k)*dt
796 enddo
797 enddo
798!$acc end kernels
800 call neg_check('mid',ipn,dt,qcheck,outqm,outtm,outum,outvm, &
801 outqcm,pretm,its,ite,kts,kte,itf,ktf,ktopm)
802 endif
804 if(ideep.eq.1)then
805 call cu_c3_deep_run( &
806 itf,ktf,its,ite, kts,kte &
807 ,flag_init &
808 ,flag_restart &
809 ,fv,r_d &
810 ,dicycle &
811 ,ichoice &
812 ,ipr &
813 ,ccn_gf &
814 ,ccnclean &
815 ,dt &
816 ,0 &
817 ,kpbli &
818 ,dhdt &
819 ,xlandi &
820 ,delp &
821 ,zo &
822 ,forcing2 &
823 ,t2d &
824 ,q2d &
825 ,tmfq &
826 ,qmicro &
827 ,forceqv_spechum &
828 ,betascu &
829 ,betamcu &
830 ,betadcu &
831 ,sigmain &
832 ,sigmaout &
833 ,ter11 &
834 ,tn &
835 ,qo &
836 ,p2d &
837 ,psur &
838 ,us &
839 ,vs &
840 ,rhoi &
841 ,hfx &
842 ,qfx &
843 ,dx &
844 ,do_ca &
845 ,progsigma &
846 ,ca_deep &
847 ,mconv &
848 ,omeg &
849 ,cactiv &
850 ,cnvwt &
851 ,zu &
852 ,zd &
853 ,zdm & ! hli
854 ,edt &
855 ,edtm & ! hli
856 ,xmb &
857 ,xmbm &
858 ,xmbs &
859 ,pret &
860 ,outu &
861 ,outv &
862 ,outt &
863 ,outq &
864 ,outqc &
865 ,kbcon &
866 ,ktop &
867 ,cupclw &
868 ,frhd &
869 ,rainevap &
870 ,ierr &
871 ,ierrc &
872! the following should be set to zero if not available
873 ,rand_mom & ! for stochastics mom, if temporal and spatial patterns exist
874 ,rand_vmas & ! for stochastics vertmass, if temporal and spatial patterns exist
875 ,rand_clos & ! for stochastics closures, if temporal and spatial patterns exist
876 ,0 & ! flag to what you want perturbed
877 ! 1 = momentum transport
878 ! 2 = normalized vertical mass flux profile
879 ! 3 = closures
880 ! more is possible, talk to developer or
881 ! implement yourself. pattern is expected to be
882 ! betwee -1 and +1
883 ,do_cap_suppress_here,cap_suppress_j &
884 ,k22 &
885 ,jmin,tropics)
886 jpr=0
887 ipr=0
888!$acc kernels
889 do i=its,itf
890 do k=kts,ktf
891 qcheck(i,k)=qv(i,k) +(outqs(i,k)+outqm(i,k))*dt
892 enddo
893 enddo
894!$acc end kernels
896 call neg_check('deep',ipn,dt,qcheck,outq,outt,outu,outv, &
897 outqc,pret,its,ite,kts,kte,itf,ktf,ktop)
898!
899 endif
900!$acc kernels
901 do i=its,itf
902 kcnv(i)=0
903 if(pretm(i).gt.0.)then
904 kcnv(i)= 1 !jmin(i)
905 cutenm(i)=1.
906 else
907 kbconm(i)=0
908 ktopm(i)=0
909 cutenm(i)=0.
910 endif ! pret > 0
911
912 if(pret(i).gt.0.)then
913 cuten(i)=1.
914 cutenm(i)=0.
915 pretm(i)=0.
916 kcnv(i)= 1 !jmin(i)
917 ktopm(i)=0
918 kbconm(i)=0
919 else
920 kbcon(i)=0
921 ktop(i)=0
922 cuten(i)=0.
923 endif ! pret > 0
924 enddo
925!$acc end kernels
926!
927!$acc parallel loop private(kstop,dtime_max,massflx,trcflx_in1,clw_in1,po_cup)
928 do i=its,itf
929 massflx(:)=0.
930 trcflx_in1(:)=0.
931 clw_in1(:)=0.
932 do k=kts,ktf
933 clw_ten(i, k)=0.
934 enddo
935 po_cup(:)=0.
936 kstop=kts
937 if(ktopm(i).gt.kts .or. ktop(i).gt.kts)kstop=max(ktopm(i),ktop(i))
938 if(ktops(i).gt.kts)kstop=max(kstop,ktops(i))
939 if(kstop.gt.2)then
940 htop(i)=kstop
941 if(kbcon(i).gt.2 .or. kbconm(i).gt.2)then
942 hbot(i)=max(kbconm(i),kbcon(i)) !jmin(i)
943 endif
944
945 dtime_max=dt
946 forcing2(i,3)=0.
947 do k=kts,kstop
948 cnvc(i,k) = 0.04 * log(1. + 675. * zu(i,k) * xmb(i)) + &
949 0.04 * log(1. + 675. * zum(i,k) * xmbm(i)) + &
950 0.04 * log(1. + 675. * zus(i,k) * xmbs(i))
951 cnvc(i,k) = min(cnvc(i,k), 0.6)
952 cnvc(i,k) = max(cnvc(i,k), 0.0)
953 cnvw(i,k)=cnvwt(i,k)*xmb(i)*dt+cnvwts(i,k)*xmbs(i)*dt+cnvwtm(i,k)*xmbm(i)*dt
954 ud_mf(i,k)=cuten(i)*zu(i,k)*xmb(i)*dt
955 dd_mf(i,k)=cuten(i)*zd(i,k)*edt(i)*xmb(i)*dt
956 t(i,k)=t(i,k)+dt*(cutens(i)*outts(i,k)+cutenm(i)*outtm(i,k)+outt(i,k)*cuten(i))
957 qv(i,k)=max(1.e-16,qv(i,k)+dt*(cutens(i)*outqs(i,k)+cutenm(i)*outqm(i,k)+outq(i,k)*cuten(i)))
958 gdc(i,k,7)=sqrt(us(i,k)**2 +vs(i,k)**2)
959 us(i,k)=us(i,k)+outu(i,k)*cuten(i)*dt +outum(i,k)*cutenm(i)*dt +outus(i,k)*cutens(i)*dt
960 vs(i,k)=vs(i,k)+outv(i,k)*cuten(i)*dt +outvm(i,k)*cutenm(i)*dt +outvs(i,k)*cutens(i)*dt
961
962 gdc(i,k,1)= max(0.,tun_rad_shall(i)*cupclws(i,k)*cutens(i)) ! my mod
963 !gdc2(i,k,1)=max(0.,tun_rad_deep(i)*(cupclwm(i,k)*cutenm(i)+cupclw(i,k)*cuten(i)))
964 gdc2(i,k,1)=max(0.,tun_rad_mid(i)*cupclwm(i,k)*cutenm(i)+tun_rad_deep(i)*cupclw(i,k)*cuten(i)+tun_rad_shall(i)*cupclws(i,k)*cutens(i))
965 !gdc2(i,k,1) = min(0.1, max(0.01, tun_rad_mid(i)*frhm(i)))*cupclwm(i,k)*cutenm(i) + min(0.1, max(0.01, tun_rad_deep(i)*(frhd(i))))*cupclw(i,k)*cuten(i) + tun_rad_shall(i)*cupclws(i,k)*cutens(i)
966 qci_conv(i,k)=gdc2(i,k,1)
967 gdc(i,k,2)=(outt(i,k))*86400.
968 gdc(i,k,3)=(outtm(i,k))*86400.
969 gdc(i,k,4)=(outts(i,k))*86400.
970 gdc(i,k,7)=-(gdc(i,k,7)-sqrt(us(i,k)**2 +vs(i,k)**2))/dt
971 !gdc(i,k,8)=(outq(i,k))*86400.*xlv/cp
972 gdc(i,k,8)=(outqm(i,k)+outqs(i,k)+outq(i,k))*86400.*xlv/cp
973 gdc(i,k,9)=gdc(i,k,2)+gdc(i,k,3)+gdc(i,k,4)
974
976 dp=100.*(p2d(i,k)-p2d(i,k+1))
977 dtime_max=min(dtime_max,.5*dp)
978 po_cup(k)=.5*(p2d(i,k)+p2d(i,k+1))
979 if (clcw(i,k) .gt. -999.0 .and. clcw(i,k+1) .gt. -999.0 )then
980 clwtot = cliw(i,k) + clcw(i,k)
981 if(clwtot.lt.1.e-32)clwtot=0.
982 clwtot1= cliw(i,k+1) + clcw(i,k+1)
983 if(clwtot1.lt.1.e-32)clwtot1=0.
984 clw_in1(k)=clwtot
985 massflx(k)=-(xmb(i) *( zu(i,k)- edt(i)* zd(i,k))) &
986 -(xmbm(i)*(zdm(i,k)-edtm(i)*zdm(i,k))) &
987 -(xmbs(i)*zus(i,k))
988 trcflx_in1(k)=massflx(k)*.5*(clwtot+clwtot1)
989 forcing2(i,3)=forcing2(i,3)+clwtot
990 endif
991 enddo
992
993 massflx(1)=0.
994 trcflx_in1(1)=0.
995 call fct1d3 (kstop,kte,dtime_max,po_cup, &
996 clw_in1,massflx,trcflx_in1,clw_ten(i,:),g)
997
998 do k=1,kstop
999 tem = dt*(outqcs(i,k)*cutens(i)+outqc(i,k)*cuten(i) &
1000 +outqcm(i,k)*cutenm(i) &
1001 +clw_ten(i,k) &
1002 )
1003 tem1 = max(0.0, min(1.0, (tcr-t(i,k))*tcrf))
1004 if (clcw(i,k) .gt. -999.0) then
1005 cliw(i,k) = max(0.,cliw(i,k) + tem * tem1) ! ice
1006 clcw(i,k) = max(0.,clcw(i,k) + tem *(1.0-tem1)) ! water
1007 else
1008 cliw(i,k) = max(0.,cliw(i,k) + tem)
1009 endif
1010
1011 enddo
1012
1013 gdc(i,1,10)=forcing(i,1)
1014 gdc(i,2,10)=forcing(i,2)
1015 gdc(i,3,10)=forcing(i,3)
1016 gdc(i,4,10)=forcing(i,4)
1017 gdc(i,5,10)=forcing(i,5)
1018 gdc(i,6,10)=forcing(i,6)
1019 gdc(i,7,10)=forcing(i,7)
1020 gdc(i,8,10)=forcing(i,8)
1021 gdc(i,10,10)=xmb(i)
1022 gdc(i,11,10)=xmbm(i)
1023 gdc(i,12,10)=xmbs(i)
1024 gdc(i,13,10)=hfx(i)
1025 gdc(i,15,10)=qfx(i)
1026 gdc(i,16,10)=pret(i)*3600.
1027
1028 maxupmf(i)=0.
1029 if(forcing2(i,6).gt.0.)then
1030 maxupmf(i)=maxval(xmb(i)*zu(i,kts:ktf)/forcing2(i,6))
1031 endif
1032
1033 if(ktop(i).gt.2 .and.pret(i).gt.0.)dt_mf(i,ktop(i)-1)=ud_mf(i,ktop(i))
1034 endif
1035 enddo
1036!$acc end parallel
1037!$acc kernels
1038 do i=its,itf
1039 if(pret(i).gt.0.)then
1040 cactiv(i)=1
1041 raincv(i)=.001*(cutenm(i)*pretm(i)+cutens(i)*prets(i)+cuten(i)*pret(i))*dt
1042 else
1043 cactiv(i)=0
1044 if(pretm(i).gt.0)raincv(i)=.001*cutenm(i)*pretm(i)*dt
1045 endif ! pret > 0
1046
1047 if(pretm(i).gt.0)then
1048 cactiv_m(i)=1
1049 else
1050 cactiv_m(i)=0
1051 endif
1052
1053 ! Unify ccn
1054 if(ccn_m(i).lt.ccn_gf(i))then
1055 ccn_gf(i)=ccn_m(i)
1056 endif
1057
1058 if(ccn_gf(i)<0) ccn_gf(i)=0
1059
1060 ! Convert ccn back to aod
1061 aod_gf(i)=0.0027*(ccn_gf(i)**0.64)
1062 if(aod_gf(i)<0.007)then
1063 aod_gf(i)=0.007
1064 ccn_gf(i)=(aod_gf(i)/0.0027)**(1/0.640)
1065 elseif(aod_gf(i)>aodc0)then
1066 aod_gf(i)=aodc0
1067 ccn_gf(i)=(aod_gf(i)/0.0027)**(1/0.640)
1068 endif
1069 enddo
1070!$acc end kernels
1071 100 continue
1072!
1073! Scale dry mixing ratios for water wapor and cloud water to specific humidy / moist mixing ratios
1074!
1075!$acc kernels
1076 qv_spechum = qv/(1.0_kind_phys+qv)
1077 cnvw_moist = cnvw/(1.0_kind_phys+qv)
1078!$acc end kernels
1079!
1080! Diagnostic tendency updates
1081!
1082 if(ldiag3d) then
1083 if(ishallow_g3.eq.1 .and. .not.flag_for_scnv_generic_tend) then
1084 uidx=dtidx(index_of_x_wind,index_of_process_scnv)
1085 vidx=dtidx(index_of_y_wind,index_of_process_scnv)
1086 tidx=dtidx(index_of_temperature,index_of_process_scnv)
1087 qidx=dtidx(100+ntqv,index_of_process_scnv)
1088 if(uidx>=1) then
1089!$acc kernels
1090 do k=kts,ktf
1091 dtend(:,k,uidx) = dtend(:,k,uidx) + cutens(:)*outus(:,k) * dt
1092 enddo
1093!$acc end kernels
1094 endif
1095 if(vidx>=1) then
1096!$acc kernels
1097 do k=kts,ktf
1098 dtend(:,k,vidx) = dtend(:,k,vidx) + cutens(:)*outvs(:,k) * dt
1099 enddo
1100!$acc end kernels
1101 endif
1102 if(tidx>=1) then
1103!$acc kernels
1104 do k=kts,ktf
1105 dtend(:,k,tidx) = dtend(:,k,tidx) + cutens(:)*outts(:,k) * dt
1106 enddo
1107!$acc end kernels
1108 endif
1109 if(qidx>=1) then
1110!$acc kernels
1111 do k=kts,ktf
1112 do i=its,itf
1113 tem = cutens(i)*outqs(i,k)* dt
1114 tem = tem/(1.0_kind_phys+tem)
1115 dtend(i,k,qidx) = dtend(i,k,qidx) + tem
1116 enddo
1117 enddo
1118!$acc end kernels
1119 endif
1120 endif
1121 if((ideep.eq.1. .or. imid_gf.eq.1) .and. .not.flag_for_dcnv_generic_tend) then
1122 uidx=dtidx(index_of_x_wind,index_of_process_dcnv)
1123 vidx=dtidx(index_of_y_wind,index_of_process_dcnv)
1124 tidx=dtidx(index_of_temperature,index_of_process_dcnv)
1125 if(uidx>=1) then
1126!$acc kernels
1127 do k=kts,ktf
1128 dtend(:,k,uidx) = dtend(:,k,uidx) + (cuten*outu(:,k)+cutenm*outum(:,k)) * dt
1129 enddo
1130!$acc end kernels
1131 endif
1132 if(vidx>=1) then
1133!$acc kernels
1134 do k=kts,ktf
1135 dtend(:,k,vidx) = dtend(:,k,vidx) + (cuten*outv(:,k)+cutenm*outvm(:,k)) * dt
1136 enddo
1137!$acc end kernels
1138 endif
1139 if(tidx>=1) then
1140!$acc kernels
1141 do k=kts,ktf
1142 dtend(:,k,tidx) = dtend(:,k,tidx) + (cuten*outt(:,k)+cutenm*outtm(:,k)) * dt
1143 enddo
1144!$acc end kernels
1145 endif
1146
1147 qidx=dtidx(100+ntqv,index_of_process_dcnv)
1148 if(qidx>=1) then
1149!$acc kernels
1150 do k=kts,ktf
1151 do i=its,itf
1152 tem = (cuten(i)*outq(i,k) + cutenm(i)*outqm(i,k))* dt
1153 tem = tem/(1.0_kind_phys+tem)
1154 dtend(i,k,qidx) = dtend(i,k,qidx) + tem
1155 enddo
1156 enddo
1157!$acc end kernels
1158 endif
1159 endif
1160 if(allocated(clcw_save)) then
1161!$acc parallel loop collapse(2) private(tem_shal,tem_deep,tem,tem1,weight_sum,cliw_both,clcw_both)
1162 do k=kts,ktf
1163 do i=its,itf
1164 tem_shal = dt*(outqcs(i,k)*cutens(i)+outqcm(i,k)*cutenm(i))
1165 tem_deep = dt*(outqc(i,k)*cuten(i)+clw_ten(i,k))
1166 tem = tem_shal+tem_deep
1167 tem1 = max(0.0, min(1.0, (tcr-t(i,k))*tcrf))
1168 weight_sum = abs(tem_shal)+abs(tem_deep)
1169 if(weight_sum<1e-12) then
1170 cycle
1171 endif
1172
1173 if (clcw_save(i,k) .gt. -999.0) then
1174 cliw_both = max(0.,cliw_save(i,k) + tem * tem1) - cliw_save(i,k)
1175 clcw_both = max(0.,clcw_save(i,k) + tem) - clcw_save(i,k)
1176 else if(cliw_idx>=1) then
1177 cliw_both = max(0.,cliw_save(i,k) + tem) - cliw_save(i,k)
1178 clcw_both = 0
1179 endif
1180 if(cliw_deep_idx>=1) then
1181 dtend(i,k,cliw_deep_idx) = dtend(i,k,cliw_deep_idx) + abs(tem_deep)/weight_sum*cliw_both
1182 endif
1183 if(clcw_deep_idx>=1) then
1184 dtend(i,k,clcw_deep_idx) = dtend(i,k,clcw_deep_idx) + abs(tem_deep)/weight_sum*clcw_both
1185 endif
1186 if(cliw_shal_idx>=1) then
1187 dtend(i,k,cliw_shal_idx) = dtend(i,k,cliw_shal_idx) + abs(tem_shal)/weight_sum*cliw_both
1188 endif
1189 if(clcw_shal_idx>=1) then
1190 dtend(i,k,clcw_shal_idx) = dtend(i,k,clcw_shal_idx) + abs(tem_shal)/weight_sum*clcw_both
1191 endif
1192 enddo
1193 enddo
1194!$acc end parallel
1195 endif
1196 endif
1197 end subroutine cu_c3_driver_run
1199end module cu_c3_driver
subroutine fct1d3(ktop, n, dt, z, tracr, massflx, trflx_in, dellac, g)
Calculates tracer fluxes due to subsidence, only up-stream differencing is currently used but flux co...
subroutine neg_check(name, j, dt, q, outq, outt, outu, outv, outqc, pret, its, ite, kts, kte, itf, ktf, ktop)
Checks for negative or excessive tendencies and corrects in a mass conversing way by adjusting the cl...
subroutine cu_c3_deep_run(itf, ktf, its, ite, kts, kte, flag_init, flag_restart, fv, r_d, dicycle, ichoice, ipr, ccn, ccnclean, dtime, imid, kpbl, dhdt, xland, delp, zo, forcing, t, q, tmf, qmicro, forceqv_spechum, betascu, betamcu, betadcu, sigmain, sigmaout, z1, tn, qo, po, psur, us, vs, rho, hfx, qfx, dx, do_ca, progsigma, ca_deep, mconv, omeg, csum, cnvwt, zuo, zdo, zdm, edto, edtm, xmb_out, xmbm_in, xmbs_in, pre, outu, outv, outt, outq, outqc, kbcon, ktop, cupclw, frh_out, rainevap, ierr, ierrc, rand_mom, rand_vmas, rand_clos, nranflag, do_capsuppress, cap_suppress_j, k22, jmin, tropics)
Driver for the deep or congestus routine.
subroutine, public cu_c3_driver_init(imfshalcnv, imfshalcnv_c3, imfdeepcnv, imfdeepcnv_c3, progsigma, cnx, mpirank, mpiroot, errmsg, errflg)
subroutine, public cu_c3_driver_run(ntracer, garea, im, km, dt, flag_init, flag_restart, do_ca, progsigma, cactiv, cactiv_m, g, cp, fv, r_d, xlv, r_v, forcet, forceqv_spechum, phil, delp, raincv, tmf, qmicro, sigmain, betascu, betamcu, betadcu, qv_spechum, t, cld1d, us, vs, t2di, w, qv2di_spechum, p2di, psuri, hbot, htop, kcnv, xland, hfx2, qfx2, aod_gf, cliw, clcw, ca_deep, rainevap, pbl, ud_mf, dd_mf, dt_mf, cnvw_moist, cnvc, imfshalcnv, flag_for_scnv_generic_tend, flag_for_dcnv_generic_tend, dtend, dtidx, ntqv, ntiw, ntcw, index_of_temperature, index_of_x_wind, index_of_y_wind, index_of_process_scnv, index_of_process_dcnv, fhour, fh_dfi_radar, ix_dfi_radar, num_dfi_radar, cap_suppress, dfi_radar_max_intervals, ldiag3d, qci_conv, do_cap_suppress, sigmaout, maxupmf, maxmf, do_mynnedmf, ichoice_in, ichoicem_in, ichoice_s_in, errmsg, errflg)
This is the Grell-Freitas convection scheme driver module.
subroutine cu_c3_sh_run(us, vs, zo, t, q, z1, tn, qo, po, psur, dhdt, kpbl, rho, hfx, qfx, xland, ichoice, tcrit, dtime, zuo, xmb_out, kbcon, ktop, k22, ierr, ierrc, flag_init, flag_restart, fv, r_d, delp, tmf, qmicro, forceqv_spechum, betascu, betamcu, betadcu, sigmain, sigmaout, progsigma, dx, outt, outq, outqc, outu, outv, cnvwt, pre, cupclw, itf, ktf, its, ite, kts, kte, ipr, tropics)
Definition cu_c3_sh.F90:75