← Back to Documentation Home
CCPP SciDoc v8.0.0  v8.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
canopy_driver.F
1
2
4 module canopy_driver
5 use mfpbltq_mod
6 use tridi_mod
7 use mfscuq_mod
8 use canopy_utils_mod
9 use satmedmfvdifq_can_mod
10
11 use canopy_mask_mod !Feb13: , only : canopy_mask_init, canopy_mask_run
12 use canopy_levs_mod !Feb13: , only : canopy_levs_init, canopy_levs_run
13 use canopy_transfer_mod !Feb13: , only : canopy_transfer_init,
14!Feb13: & canopy_transfer_run
15
16 contains
17
21 subroutine canopy_driver_init (satmedmf,
22 & isatmedmf,isatmedmf_vdifq,
23 & errmsg,errflg)
24
25 logical, intent(in ) :: satmedmf
26 integer, intent(in) :: isatmedmf,isatmedmf_vdifq
27
28 character(len=*), intent(out) :: errmsg
29 integer, intent(out) :: errflg
30
31! Initialize CCPP error handling variables
32 errmsg = ''
33 errflg = 0
34
35! Consistency checks
36 if (.not. satmedmf) then
37 write(errmsg,fmt='(*(a))') 'Logic error: satmedmf = .false.'
38 errflg = 1
39 return
40 end if
41
42 if (.not. isatmedmf==isatmedmf_vdifq) then
43 write(errmsg,fmt='(*(a))') 'Logic error: satmedmfvdif is ',
44 & 'called, but isatmedmf/=isatmedmf_vdifq.'
45 errflg = 1
46 return
47 end if
48
49 end subroutine canopy_driver_init
50
54 subroutine canopy_driver_run(im,km,ntrac,ntcw,ntrw,
55 & ntiw,ntke,ntchm,ntchs,ntche, ntoz,nto3,
56 & ndtend, ! in
57 & con_rocp,
58 & grav,pi,rd,cp,rv,hvap,hfus,fv,eps,epsm1,
59!The following three variables are for SA-3D-TKE
60 & def_1,def_2,def_3,sa3dtke,dku3d_h,dku3d_e,
61 & dv,du,tdt,rtg, ! inout: dv,du,tdt,rtg
62 & u1,v1,t1,q1,usfco,vsfco,use_oceanuv, ! in
63 & swh,hlw,xmu,garea,zvfun,sigmaf, ! in
64 & psk,rbsoil,zorl,u10m,v10m,t2m, q2m,fm,fh,
65 & tsea,heat,evap,stress,spd1,
66 & kpbl, ! in: kpbl
67 & pgr,
68 & prsi,del,prsl,prslk,phii,phil,delt,tte_edmf,
69 & dspheat, ! in
70 & dusfc,dvsfc,dtsfc,dqsfc,hpbl, ! in: dusfc,dvsfc,dtsfc,dqsfc,hpbl
71 & dkt,dku,tkeh, ! inout: dkt,dku, tkeh
72 & dkt_can,dku_can, ! out
73 & kinver,xkzm_m,xkzm_h,xkzm_s,dspfac,bl_upfr,bl_dnfr, ! in
74 & rlmx,elmx,sfc_rlm,tc_pbl,use_lpt, ! in
75 & do_canopy, cplaqm, claie, cfch, cfrt, cclu, cpopu, ! in
76 & ntqv,dtend,dtidx,index_of_temperature,index_of_x_wind, ! inout: dtend (.ldiag3d.)
77 & index_of_y_wind,index_of_process_pbl,gen_tend,ldiag3d,
78 & errmsg,errflg)
79
80!
81 use machine , only : kind_phys
82 use funcphys , only : fpvs
83!
84 implicit none
85!
86!----------------------------------------------------------------------
87 integer, intent(in) :: im, km,
88 & ntrac, ntcw, ntrw, ntiw, ntke, ntqv,
89 & ntchm,ntchs,ntche, ntoz,nto3, ndtend
90 integer, intent(in) :: sfc_rlm
91 integer, intent(in) :: tc_pbl
92 integer, intent(in) :: use_lpt
93 integer, intent(in) :: kinver(:)
94 integer, intent(in) :: kpbl(:)
95 logical, intent(in) :: gen_tend,ldiag3d
96!
97 real(kind=kind_phys), intent(in) :: grav,pi,rd,cp,rv,hvap,hfus,fv,
98 & eps,epsm1,
99 & con_rocp
100 real(kind=kind_phys), intent(in) :: delt, xkzm_m, xkzm_h, xkzm_s
101 real(kind=kind_phys), intent(in) :: dspfac, bl_upfr, bl_dnfr
102 real(kind=kind_phys), intent(in) :: rlmx, elmx
103!PCC CANOPY------------------------------------
104 logical, intent(in) :: do_canopy, cplaqm
105
106 real(kind=kind_phys), optional, intent(in) ::
107! 2D
108 & claie(:) , cfch(:),
109 & cfrt(:) , cclu(:), cpopu(:)
110
111 !----------------------------------------------
112 real(kind=kind_phys), intent(inout) ::
113 & dv(:,:), du(:,:),
114 & tdt(:,:), tkeh(:,:),
115 & rtg(:,:,:)
116 real(kind=kind_phys), intent(in) ::
117 & u1(:,:), v1(:,:),
118 & usfco(:), vsfco(:),
119 & t1(:,:), q1(:,:,:),
120!The following two variables are for SA-3D-TKE
121 & def_1(:,:), def_2(:,:), def_3(:,:),
122 & swh(:,:), hlw(:,:),
123 & xmu(:), garea(:),
124 & zvfun(:), sigmaf(:),
125 & psk(:), rbsoil(:),
126 & zorl(:), tsea(:),
127 & u10m(:), v10m(:),
128 & t2m(:), q2m(:),
129 & fm(:), fh(:),
130 & evap(:), heat(:),
131 & stress(:), spd1(:),
132 & pgr(:),
133 & prsi(:,:), del(:,:),
134 & prsl(:,:), prslk(:,:),
135 & phii(:,:), phil(:,:)
136
137 real(kind=kind_phys), intent(inout), dimension(:,:,:), optional ::
138 & dtend
139
140 integer, intent(in) :: dtidx(:,:), index_of_temperature,
141 & index_of_x_wind, index_of_y_wind, index_of_process_pbl
142 logical, intent(in) :: use_oceanuv
143 real(kind=kind_phys), intent(in) ::
144 & dusfc(:), dvsfc(:),
145 & dtsfc(:), dqsfc(:)
146 real(kind=kind_phys), intent(in) ::
147 & hpbl(:) ! use resolved hpbl in non-canopy columns
148 real(kind=kind_phys), intent(inout) ::
149 & dkt(:,:), dku(:,:)
150 real(kind=kind_phys), optional, intent(out) ::
151 & dkt_can(:,:), dku_can(:,:)
152
153 logical, intent(in) :: sa3dtke !flag for SA-3D-TKE scheme
154!
155!
156! flag for tke dissipative heating
157 logical, intent(in) :: dspheat
158! flag for TTE-EDMF scheme
159 logical, intent(in) :: tte_edmf
160
161 character(len=*), intent(out) :: errmsg
162 integer, intent(out) :: errflg
163
164!For passing dku to the dyn_core (SA-3D-TKE scheme)
165 real(kind=kind_phys), intent(in) ::
166 & dku3d_h(:,:),dku3d_e(:,:)
167
168!
169!----------------------------------------------------------------------
170!***
171!*** local variables
172
173 integer i,is,k,n,ndt,km1,kmpbl,kmscu,ntrac1,
174 & idtend
175
176 real(kind=kind_phys) dtdz1(im), gdx(im),
177 & phih(im), phim(im),
178 & phims(im), prn(im,km-1),
179 & rbdn(im), rbup(im), thermal(im),
180 & ustar(im), wstar(im),
181 & ust3(im), wst3(im),
182 & z0(im), crb(im), tkemean(im),
183 & hgamt(im), hgamq(im),
184 & wscale(im),vpert(im),
185 & zol(im), sflux(im),
186 & sumx(im), tx1(im), tx2(im)
187!
188 real(kind=kind_phys) zi(im,km+1), zl(im,km), zm(im,km),
189 & xkzo(im,km-1),xkzmo(im,km-1),
190 & xkzm_hx(im), xkzm_mx(im), tkmnz(im,km-1),
191 & rdzt(im,km-1),rlmnz(im,km),
192 & al(im,km-1), ad(im,km), au(im,km-1),
193 & f1(im,km), f2(im,km*(ntrac-1))
194!
195 real(kind=kind_phys) elm(im,km), ele(im,km),
196 & ckz(im,km), chz(im,km),
197 & diss(im,km-1),prod(im,km-1),
198 & bf(im,km-1), shr2(im,km-1),
199 & xlamue(im,km-1), xlamde(im,km-1),
200 & gotvx(im,km), rlam(im,km-1)
201!
202 real(kind=kind_phys) aphi16, aphi5,
203 & wfac, cfac,
204 & gamcrt, gamcrq, sfcfrac,
205! & conq, cont, conw,
206 & dsdz2, dsdzt, dkmax,
207 & dsig, dt2, dtodsd,
208 & dtodsu, g, factor, dz,
209 & gocp, gravi, zol1, zolcru,
210 & concmin,
211 & buop, shrp, dtn,
212 & prnum, prmax, prmin, prtke,
213 & prscu, pr0, ri,
214 & dw2, dw2min, zk,
215 & elmfac, elefac, dspmax,
216 & alp, clwt, cql,
217 & f0, robn, crbmin, crbmax,
218 & es, qs, value, onemrh,
219 & cfh, gamma, elocp, el2orc,
220 & epsi, beta, chx, cqx,
221 & rdt, rdz, qmin, qlmin,
222 & rimin, rbcr, rbint, tdzmin,
223 & rlmn, rlmn0, rlmn1, rlmn2,
224 & ttend, utend, vtend, qtend,
225 & zfac, zfmin, vk, spdk2,
226 & tkmin, tkbmx, xkgdx,
227 & xkinv1, xkinv2,
228 & zlup, zldn, cs0, csmf,
229 & tem, tem1, tem2, tem3,
230 & ptem, ptem0, ptem1, ptem2
231!
232 real(kind=kind_phys) slfac
233!
234 real(kind=kind_phys) vegflo, vegfup, z0lo, z0up, vc0, zc0
235!
236 real(kind=kind_phys) ck0, ck1, ch0, ch1, ce0, rchck
237!
238 real(kind=kind_phys) qlcr, zstblmax, hcrinv
239!
240 real(kind=kind_phys) h1
241
242 real(kind=kind_phys) bfac, mffac
243
244 real(kind=kind_phys) q1_new(im,km,ntrac-1)
245
246 integer kount
247
248!PCC_CANOPY------------------------------------
249 integer COUNTCAN,KCAN
250 ! Declare local maximum canopy layers
251 integer, parameter :: MAXCAN = 1000
252
253 real(kind=kind_phys) :: fch,
254 & zi05, zfl, dzfl, botcan, fzi05, ttcorr,
255 & uucan, vvcan, ttcan, tkecan, tkehcan,
256 & uu1, vv1, tt1, tke1, tkeh1,
257 & uu_int,vv_int,tt_int,tke_int, tkeh_int,
258 & tt_sum, zz_int
259
260 real(kind=kind_phys) :: zcanx(maxcan), zooox(maxcan),
261 & uux(maxcan), vvx(maxcan), ttx(maxcan),
262 & tkex(maxcan), tkehx(maxcan)
263!PCC_CANOPY------------------------------------
264
265 real(kind=kind_phys) ::
266 & dv_can(im,km), du_can(im,km),
267 & duv_can(im,km),
268 & tdt_can(im,km), rtg_can(im,km,ntrac)
269
270 real(kind=kind_phys) ::
271 & dum3d_h(im, km) , dum3d_e(im, km),
272 & tkeh_mod(im, km)
273
274 real(kind=kind_phys) :: dtend_can(im, km , ndtend)
275
276! Out: list sat_canopy call
277 real(kind=kind_phys) :: dusfc_can(im), dvsfc_can(im),
278 & dtsfc_can(im), dqsfc_can(im),
279 & hpbl_can(im)
280
281 integer :: kpbl_can(im)
282
283 real(kind=kind_phys) ::
284 & rho1(im,km),
285 & t2(im,km),
286 & u2(im,km), v2(im,km), ws2(im,km),
287 & rho2(im,km), duv(im,km)
288
289! Number of canopy layers
290 integer, parameter :: nkc = 3
291 integer :: nkt ! # of resolved model layers plus canopy layers
292
293 integer
294 & kcan1, kc_can,
295 & kc, nkt1 ,
296 & kmod(im, km) ,
297 & kcan3(im, nkc) ,
298 & nfrct(km+nkc, im) , !nkt
299 & ifrct(km+nkc, 2, im) !nkt
300
301 real(kind=kind_phys) ::
302! 2D arrays
303 & frt_mask(im) ,
304! met2d arrays
305 & u2_mod(im, km) , v2_mod(im, km),
306 & ws2_mod(im, km) , t2_mod(im, km),
307 & du_mod(im, km) , dv_mod(im, km),
308 & duv_mod(im, km) , tdt_mod(im, km),
309 & u10m_can(im) ,
310 & v10m_can(im) ,
311! all gas-phase species array
312! NB. mfpbltq_mod: q1(ix,km,ntrac1) with ntrac1 = ntrac - 1
313 & rtg1_mod(im, km, ntrac), ! before diffusion
314 & q1_mod(im, km, ntrac), ! before diffusion
315 & q2_mod(im, km, ntrac), ! after diffusion
316 & q2(im, km, ntrac), ! after diffusion
317 & rtg_mod(im, km, ntrac), ! after diffusion
318! sat_can inputs
319 & q1_can(im, km , ntrac), ! size (km) before diffusion
320 & q2_can(im, km , ntrac), ! size (km) after diffusion
321 & u1_can(im, km ) , ! size (km)
322 & v1_can(im, km ) , ! size (km)
323 & ws1_can(im, km ) , ! size (km)
324 & t1_can(im, km ) , ! size (km)
325 & tkeh_can(im, km ) ,
326!
327 & swh_can(im, km) ,
328 & hlw_can(im, km) ,
329! sat_can inputs
330 & phii_can(im, km+1) ,
331 & prsi_can(im, km+1) ,
332 & prsl_can(im, km) ,
333 & del_can(im, km) ,
334 & prslk_can(im, km) ,
335 & phil_can(im, km) ,
336! Canopy layers
337 & dkt_can3(im, km+nkc) ,
338 & dku_can3(im, km+nkc) ,
339 & qv_can3(im, km+nkc) , ! nkt before diffusion
340 & q1_can3(im, km+nkc, ntrac), ! nkt before diffusion
341 & q2_can3(im, km+nkc, ntrac), ! nkt after diffsion
342! before diffusion
343 & ws1_can3(im, km+nkc) , ! using km for now only
344 & u1_can3(im, km+nkc) , ! using km for now only
345 & v1_can3(im, km+nkc) , ! using km for now only
346 & t1_can3(im, km+nkc) , ! using km for now only
347 & rho1_can3(im, km+nkc) ,
348 & tkeh_can3(im, km+nkc) ,
349! after diffusion
350 & ws2_can3(im, km+nkc) ,
351 & u2_can3(im, km+nkc) ,
352 & v2_can3(im, km+nkc) ,
353 & t2_can3(im, km+nkc) ,
354 & rho2_can3(im, km+nkc) ,
355! 2D
356 & q1_2m(im, ntrac), ! before diffusion
357 & q2_2m(im, ntrac), ! after diffusion
358! met3d arrays
359 & phii_can3(im, km+nkc+1) , !nkt
360 & phil_can3(im, km+nkc) ,
361 & prsi_can3(im, km+nkc+1) ,
362 & prsl_can3(im, km+nkc) ,
363 & prslk_can3(im, km+nkc) ,
364 & del_can3(im, km+nkc) ,
365!
366 & zi_can3(im, km+nkc+1) , !nkt
367 & zl_can3(im, km+nkc) , ! zl_can is ZH_CAN
368 & zm_can3(im, km+nkc) , ! zm_can is ZF_CAN
369 & dz_can3(im, km+nkc) ,
370!
371! model layers
372 & wind_dir_to_rad(im, km) ,
373 & ws1(im, km) ,
374 & wdir(im, km) ,
375
376! layer height arrays !layers are in reverse order!
377! 1 is top resolved layer
378! km is bottom model hybrid layer
379! km+nkc=nkt is bottom canopy layer
380 & zmom_can3(im, km+nkc+1) , ! dim(im, nkt+1)
381 & zmid_can3(im, km+nkc) ,
382 & sigmom_can3(im, km+nkc+1) , ! ~zm (nkt) or ~zi (nkt+1)
383 & sigmid_can3(im, km+nkc) , ! ~zl
384 & massair_can3(im, km+nkc) ,
385 & massair(im, km) ,
386 & mmr_o3_can3(im, km+nkc) ,
387 & frctr2c(km+nkc, 2, im) ,
388 & frctc2r(km+nkc, 2, im)
389
390!!
391 parameter(bfac=100.)
392 parameter(wfac=7.0,cfac=4.5)
393 parameter(gamcrt=3.,gamcrq=0.,sfcfrac=0.1)
394 parameter(vk=0.4,rimin=-100.,slfac=0.1)
395 parameter(rbcr=0.25,zolcru=-0.02,tdzmin=1.e-3)
396 parameter(rlmn=30.,rlmn0=5.,rlmn1=5.,rlmn2=10.)
397 parameter(prmin=0.25,prmax=4.0)
398 parameter(pr0=1.0,prtke=1.0,prscu=0.67)
399 parameter(f0=1.e-4,crbmin=0.15,crbmax=0.35)
400 parameter(tkmin=1.e-9,tkbmx=0.2,dspmax=10.0)
401 parameter(qmin=1.e-8,qlmin=1.e-12,zfmin=1.e-8)
402 parameter(concmin = 1.0e-30) ! Minimum conc "rbdriver"
403 parameter(aphi5=5.,aphi16=16.)
404 parameter(elmfac=1.0,elefac=1.0,cql=100.)
405 parameter(dw2min=1.e-4,dkmax=1000.,xkgdx=1000.)
406 parameter(qlcr=3.5e-5,zstblmax=2500.)
407 parameter(xkinv1=0.15,xkinv2=0.3)
408 parameter(h1=0.33333333,hcrinv=250.)
409 parameter(vegflo=0.1,vegfup=1.0,z0lo=0.1,z0up=1.0)
410 parameter(vc0=1.0,zc0=1.0)
411 parameter(ck1=0.15,ch1=0.15)
412 parameter(cs0=0.4,csmf=0.5)
413 parameter(rchck=1.5,ndt=20)
414 real(kind=kind_phys), parameter :: epsilon = 1.e-10
415
416! Consistency checks
417 if (.not. (do_canopy .and. cplaqm ) ) then
418 write(errmsg,fmt='(*(a))') 'Logic error: do_canopy = .false.'
419 return
420 end if
421
422! Number of combined canopy plus resolved model layers
423 nkt = km + nkc ! # of resolved model layers plus canopy layers
424 nkt1 = nkt - 1
425
426 if (tc_pbl == 0) then
427 ck0 = 0.4
428 ch0 = 0.4
429 ce0 = 0.4
430 else if (tc_pbl == 1) then
431 ck0 = 0.55
432 ch0 = 0.55
433 ce0 = 0.12
434 endif
435 gravi = 1.0 / grav
436 g = grav
437 gocp = g / cp
438! cont=cp/g
439! conq=hvap/g
440! conw=1.0/g ! for del in pa
441!! parameter(cont=1000.*cp/g,conq=1000.*hvap/g,conw=1000./g) !kpa
442 elocp = hvap / cp
443 el2orc = hvap * hvap / (rv * cp)
444
445!
446!************************************************************************
447! Initialize CCPP error handling variables
448 errmsg = ''
449 errflg = 0
450
452 dt2 = delt
453 rdt = 1. / dt2
454!
455! the code is written assuming ntke=ntrac
456! if ntrac > ntke, the code needs to be modified
457!
458 ntrac1 = ntrac - 1
459 km1 = km - 1
460 kmpbl = km / 2
461 kmscu = km / 2
464 do k=1,km
465 do i=1,im
466 zi(i,k) = phii(i,k) * gravi
467 zl(i,k) = phil(i,k) * gravi
468 ckz(i,k) = ck1
469 chz(i,k) = ch1
470 rlmnz(i,k) = rlmn0
471 enddo
472 enddo
473 do i=1,im
474 zi(i,km+1) = phii(i,km+1) * gravi
475 enddo
476 do k=1,km
477 do i=1,im
478 zm(i,k) = zi(i,k+1)
479 enddo
480 enddo
482 do i=1,im
483 gdx(i) = sqrt(garea(i))
484 enddo
485
487! do i = 1,im
488! z0(i) = 0.01 * zorl(i)
489! rho_a(i) = prsl(i,1)/(rd*t1(i,1)*(1.+fv*max(q1(i,1,1),qmin)))
490
491! pblflg(i)= .true.
492! sfcflg(i)= .true.
493! if(rbsoil(i) > 0.) sfcflg(i) = .false.
494! pcnvflg(i)= .false.
495! scuflg(i)= .true.
496! if(scuflg(i)) then
497! radmin(i)= 0.
498! mrad(i) = km1
499! krad(i) = 1
500! lcld(i) = km1
501! kcld(i) = km1
502! endif
503! enddo
504
505! Initialize canopy layers concentrations with values before diffusion
506 if (do_canopy .and. cplaqm) then
507
508 do k = 1,km
509 do i = 1,im
510 rho1(i,k) = prsl(i,k)/
511 & (rd*t1(i,k)*
512 & (1.+fv*max(q1(i,k, 1),qmin)))
513
514 enddo
515 enddo
516
517 CALL canopy_mask_init( im, km, nkc, nkt,
518 & claie, cfch, cfrt, cclu, cpopu, ! in:
519 & frt_mask, ! out
520 & errmsg, errflg)
521
522 if (errflg /= 0) return
523
524 CALL canopy_levs_init( im, km, nkc, nkt,
525 & ntrac, ntqv, ntke, ! ndtend, ! in
526 & zi, zl, zm, ! in 3D
527 & prsl, prsi, ! in 3D
528 & dv, du, tdt, rtg, ! in 3D
529 & u1, v1, t1, q1, ! in 3D / 4D
530 & rho1, dkt, dku, ! in 3D
531! & dtend,
532 & zmom_can3, zmid_can3, !out 3D
533 & sigmom_can3, sigmid_can3, !out 3D
534 & zl_can3, zm_can3, !out 3D ZL=half- and ZM=full-layer height
535 & prsl_can3, prsi_can3, !out 3D set to zero
536 & dv_can, du_can, tdt_can, rtg_can, !out 3D size (km)
537 & t1_can3, qv_can3,rho1_can3, !out 3D set to zero
538 & ws1_can3, dkt_can3, dku_can3, !out 3D set to WS DKT
539 & q1_can3, q1_2m, !out 4D set to Q1
540! & DTEND_CAN,
541 & errmsg, errflg)
542
543! Zero in-canopy tendencies
544 dtend_can(:, :, : ) = 0.0
545
546 if (errflg /= 0) return
547
548! ================
549! In; Q1 (:,km ,NTRAC)
550! Out: Q1_CAN3(:,1:3,NTRAC) <= Q1(:,1,NTRAC) ! ALL CANOPY & NON-CANOPY COLUMNS !
551! ===============
552
553 CALL canopy_transfer_init(im, km, nkc, nkt, !in
554 & massair_can3, massair, !out
555 & mmr_o3_can3, !inout
556 & nfrct, ifrct, !out
557 & frctr2c, frctc2r, !out
558 & errmsg, errflg )
559
560 if (errflg /= 0) return
561
562
563 endif ! (do_canopy .and. cplaqm)
564
565
566 !PCC_CANOPY------------------------------------
567 kount=0
568 if (do_canopy .and. cplaqm) then
569
570! NB. Call canopy routines after eddy diffusivities are calculated!!!
571
572 CALL canopy_mask_run( im, km, nkc, nkt, !in
573 & claie, cfch, cfrt, cclu, cpopu, !in
574 & frt_mask, !out
575 & errmsg, errflg)
576
577! if (errflg /= 0) return
578
579! Wind direction, degrees
580! ATAN2(Y, X) computes the principal value of the argument function of the complex number X + i Y.
581! This function can be used to transform from Cartesian into polar coordinates and allows to determine the angle in the correct quadrant.
582
583 do k=1,km
584 do i=1,im
585
586 ws1(i,k) = sqrt(u1(i,k)**2+v1(i,k)**2)
587
588 wind_dir_to_rad(i,k) =
589 & atan2(u1(i,k),v1(i,k)) ! to radians
590! & atan2(u1(i,k)/ws1(i,k),v1(i,k)/ws1(i,k)) ! to radians
591 enddo
592 enddo
593
594! ===============
595! In canopy_levs_run, the vertical layers are going down (reversed to phot)
596!
597! k = 1, NLAYT from top to bottom of augmented canopy+resolved model layers
598! 1 = top model layer
599! NLAYS = bottom resolved model layer
600! NLAYC = 3 canopy layers
601! NLAYT = NLAYS + NLAYC, augmented canopy plus resolved model layers
602! NLAYS+1= top canopy layer at 1.0*CH
603! NLAYS+2= middle canopy layer at 0.5*CH
604! NLAYS+3= bottom canopy layer at 0.2*CH
605! =================
606
607 CALL canopy_levs_run(im, km, nkc, nkt, ! in
608 & ntrac, ntqv, ntke, ! in
609 & rd, pi, ! in gry gas constant
610 & zi, zl, zm, ! in
611 & prsl, prsi, pgr, ! in (Pa)
612 & cfch, ! in: canopy data input
613 & garea, u10m,v10m, fm, fh, ! in: 2D
614 & rbsoil, ! in: 2D
615 & t2m, q2m, ! in 2D
616 & stress, spd1, ! in: 2D
617 & dv, du, tdt, rtg, ! in: 3D
618 & u1, v1, t1, q1, ! in: 3D " 4D
619 & rho1, dkt, dku, ! in 3D
620 & frt_mask, ! in 2D canopy_mask
621 & kmod, kcan3, ! out
622 & zmom_can3, zmid_can3, ! out
623 & sigmom_can3, sigmid_can3, ! out
624 & zl_can3, zm_can3, ! out: zl_can=half- and zm_can=full-layer height
625 & prsl_can3, prsi_can3, ! out: mean layer pressure; air pressure at model layer interfaces
626 & dv_can, du_can, tdt_can, rtg_can, ! out: 3D size (km)
627 & t1_can3, qv_can3, rho1_can3, ! out 3D: 2-m interpolated T1 QV rho1
628 & ws1_can3, dkt_can3, dku_can3, ! out 3D: 10-m interpolated WS1
629 & q1_can3, q1_2m, ! inout kg kg-1
630 & errmsg, errflg) ! out
631
632 if (errflg /= 0) return
633
634! ================
635! Out:
636! T1_CAN3 (:,:,NLAYT)
637! QV_CAN3 (:,:,NLAYT) Q2m interpolated
638! PRSI_CAN3 ( NLAYT+1)
639! PRSL_CAN3 (:, NLAYT)
640! rho1_CAN3 (:, NLAYT)
641! Q1_CAN3(1,2,3) <= Q1(1) ! ALL CANOPY & NON-CANOPY COLUMNS !
642
643! !Layers in reverse order!
644! 1 is top resolved layer
645! km is bottom model hybrid layer
646! nkt is bottom canopy layer
647! zmid_can3 (:,:, NLAYT) layers are in reverse order!
648! zmom_can3 (:,:, NLAYT+1) layers are in reverse order!
649! massair_can3(:,:, NLAYT) : mass of air in canopy layers (kg)
650! massair (:,:, NLAYS) : mass of air in model layers (kg)
651! (gathered canopy + resolved scale columns)
652! nfrct (NLAYT, :,:) : Number of original model levels contributing to canopy level k
653! ifrct (NLAYT,2,:,:) : Index of the original model level contributing to canopy level k
654! frctr2c(NLAYT,2,:,:) : Fractional contribution of the original model level to canopy level k
655! frctc2r(NLAYT,2,:,:) : Fractional contribution of the canopy level to the original model level
656! ================
657
658 do i = 1,im
659! rho1_can_sfc(i) = prsl_can3(i,1)/
660! & (rd*t1_can3(i,1)*(1.+fv*max(q1_can3(i,1, ntqv),qmin))) ! ntqv=1
661 enddo
662
663 do k = 1,km ! ntk to do
664 do i = 1,im
665
666 q1_can3(i,k, ntqv) = qv_can3(i,k) ! 2-m interpolated humidity
667
668 rho1_can3(i,k) = prsl_can3(i,k)/
669 & (rd*t1_can3(i,k)*
670 & (1.+fv*max(q1_can3(i,k, ntqv),qmin))) ! ntqv=1 2-m interpolated
671
672 enddo
673 enddo
674
675! Humidity on canopy layers, interpolated from q2m qv_can3
676! Humidity overwritten in "resolved_to_canopy" mass transfer
677! q1_can3(:,1:km, ntqv ) = qv_can3(:,1:km)
678
679! Above canopy layers
680 do k=1,km
681 do i=1,im
682 ! kc = 4,5,6.. 67
683 kc = nkc + k
684 del_can3(i,kc) = prsi_can3(i,kc) - prsi_can3(i,kc+1)
685
686! Exner function canopy layers
687! !< exner function = (p/p0)**rocp
688 prslk_can3(i,kc) = (prsl_can3(i,kc) /pgr(i)) ** con_rocp
689
690 enddo
691 enddo
692
693! Canopy layers: kc = 1, 2, 3
694 do kc = 1, nkc
695 do i = 1,im
696
697! Canopy columns
698 IF (frt_mask(i) > 0.) THEN
699 del_can3(i,kc) = prsi_can3(i, kc) - prsi_can3(i, kc+1)
700
701! Non-canopy columns set to del(1)
702 ELSE IF (frt_mask(i) <= 0.) THEN
703 del_can3(i,kc) = del(i,1)
704
705 ENDIF
706
707! Exner function canopy layers
708! !< exner function = (p/p0)**rocp
709 prslk_can3(i,kc) = (prsl_can3(i,kc) /pgr(i)) ** con_rocp
710
711 end do
712 end do
713
714!--------
715! zmid(1) = zl(km) is top height model layer centers
716! zmid(km) = zl(1) is bottom height model layer centers
717! NB.
718! ZH_CAN (i, nkt + 1 - k) = zmid_can(i, k) (k = 1, nkt) combined layer centers => rename zl_can
719! zi height model layer interfaces
720! dz_can = zl_can3(i,k+1) - zl_can3(i,k)
721! so zm(i,k) = zi(i,k+1) = zl_can3(i,k) + dz_can3(i,k)/2
722! zm (:,k) = zi(:,k+1), so zm_can (i, k) = zi_can3(i,k+1) (k=1,km)
723
724! In-Canopy layers: kc = 1, 2, 3
725 zi_can3(:, 1) = 0.
726 do kc = 1, nkc
727 do i = 1,im
728 ! kc+1 = 2, 3, 4
729 zi_can3(i,kc+1) = zm_can3(i, kc)
730 dz_can3(i,kc) = zi_can3(i, kc+1) - zi_can3(i, kc)
731 end do
732 end do
733
734! Above canopy layers
735 do k = 1,km
736 do i = 1,im
737 ! kc = 4,5,6.. 67
738 kc = nkc + k
739
740 ! kc+1 = 5,6,7 ...68
741 zi_can3(i,kc+1) = zm_can3(i, kc) ! upper interface
742
743! print*,'canopy_driver: zm_can3= ',i, kc, ! 5
744! & zm_can3(i, kc),
745! & zi_can3(i, kc)
746
747 dz_can3(i,kc) = zi_can3(i, kc+1) - zi_can3(i, kc)
748 end do
749 end do
750
751
752
755 do k=1,nkt
756 do i=1,im
757 phil_can3(i,k) = zl_can3(i,k) * grav
758 phii_can3(i,k) = zi_can3(i,k) * grav
759 enddo
760 enddo
761 do i=1,im
762 phii_can3(i,nkt+1) = zi_can3(i,nkt+1) * grav
763 enddo
764
765 do k = 1,km
766
767 ! nkt is top
768 ! nkc+1 is bottom
769 kc=nkc+k
770
771 do i=1,im
772
773 u1_can3(i,kc) = ws1_can3(i,kc) * sin(wind_dir_to_rad(i,k)) ! m/s
774 v1_can3(i,kc) = ws1_can3(i,kc) * cos(wind_dir_to_rad(i,k)) ! m/s
775
776 end do
777 end do
778
779! Canopy Layers: use 1hy resolved model layer wind direction
780 do kc = 1,nkc
781 do i=1,im
782
783 u1_can3(i,kc) = ws1_can3(i,kc) * sin(wind_dir_to_rad(i,1)) ! m/s
784 v1_can3(i,kc) = ws1_can3(i,kc) * cos(wind_dir_to_rad(i,1)) ! m/s
785
786 end do
787 end do
788
789! All columns & all layers
790 q1_mod(:,1:km, :) = q1(:,1:km, :) ! before "resolved_to_canopy"
791! Q1_2M ... move out of "canopy_levs"
792
793!===============================================================================
794! Distribute tracer concentration from model resolved layers into canopy layers
795! flag = 0 "resolved_to_canopy"
796!===============================================================================
797
798 CALL canopy_transfer_run(im, km, nkc, nkt, !in
799 & ntrac, ntoz, !in
800 & garea, !in
801 & zi, zl, zm, !in
802 & q1, rho1, !in kg kg-1
803 & 0, !in 0 = "resolved_to_canopy"
804 & frt_mask, !in
805 & kmod, kcan3, !in
806 & zmom_can3, zmid_can3, !in
807 & prsl_can3, rho1_can3, !in: before diffusion
808 & q1_mod, q1_can3, q1_2m, !inout: kg kg-1 before diffusion
809 & massair_can3, massair, !inout
810 & mmr_o3_can3, !inout
811 & nfrct, ifrct, !inout
812 & frctr2c, frctc2r, !inout
813 & errmsg, errflg )
814
815 if (errflg /= 0) return
816
817! ==============
818! Input:
819! Q1 (:,:, NLAYS, ntrac) : Chemical tracers conc. ppmv on model levels
820!
821! Output:
822! Q1_CAN3(:,:, NLAYT, ntrac) : Chemical tracers conc. ppmv on combined canopy+resolved layers
823! ! CANOPY COLUMNS ONLY !
824! Q1_2M (:,: , ntrac) : 2M Chemical tracers conc. ppmv Diagnostics
825!
826! ================================
827
828! Comment out to use 2-m interpolated value from "canopy_levs" qv_can3(:,:)
829! q1_can3 (:,1:km, ntqv) = qv_can3(:,1:km) ! ntqv=1
830
831! Subset the canopy tracers/arrays for input to "sat_can", since routine is coded on dim(km)
832 prsi_can(:,1:km+1) = prsi_can3(:,1:km+1)
833 prsl_can(:,1:km) = prsl_can3(:,1:km)
834 prslk_can(:,1:km) = prslk_can3(:,1:km)
835
836 do i = 1,im
837 IF (frt_mask(i) > 0.) THEN
838 del_can(i,1:km) = del_can3(i,1:km)
839 phii_can(i,1:km+1) = phii_can3(i,1:km+1)
840 phil_can(i,1:km) = phil_can3(i,1:km)
841! Set Non-Canopy columns to resolved layer thickness
842 ELSE IF (frt_mask(i) <= 0.) THEN
843 del_can(i,1:km) = del(i,1:km)
844 phii_can(i,1:km+1) = phii(i,1:km+1)
845 phil_can(i,1:km) = phil(i,1:km)
846 ENDIF
847 enddo
848
849! NB. Using 10-m interpolated values creates shear and gives very high TKE tendencies
850 u1_can(:,1:km) = u1_can3(:,1:km)
851 v1_can(:,1:km) = v1_can3(:,1:km)
852
853 ws1_can(:,1:km) = sqrt(u1_can(:,1:km)**2+v1_can(:,1:km)**2)
854
855 t1_can(:,1:km) = t1_can3(:,1:km)
856
857 dku_can(:,1:km) = dku_can3(:,1:km) ! "canopy_levs"
858 dkt_can(:,1:km) = dkt_can3(:,1:km) ! "canopy_levs"
859
860! Prepare arrays for "Resolved_to_Canopy" transfer
861! Mass tracers (ntrac1) except TKE
862 q1_can(:,1:km, 1:ntrac1) = q1_can3(:,1:km, 1:ntrac1) ! ntrac1 "resolved_to_canopy"
863! TKE tracer
864 q1_can(:,1:km, ntke ) = q1_can3(:,1:km, ntke ) ! ntke "resolved_to_canopy"
865
866 endif !do_canopy .and. cplaqm
867
868 if (do_canopy .and. cplaqm) then
869
870! 3D array on combined canopy plus resolved model layers
871
872! Test with km combined canopy plus resolved layers, so skip the top combined 3 layers
873! This should be nkt layers... dv_can(:,nkc+1:nkt) = dv(:,1:km) ! nkt combined canopy plus resolved layers
874
875! Sub-Canopy
876 swh_can(:, nkc+1:km ) = swh(:,1:km-nkc)
877 swh_can(:, 3 ) = swh(:,1 )
878 swh_can(:, 2 ) = swh(:,1 )
879 swh_can(:, 1 ) = swh(:,1 )
880
881 hlw_can(:, nkc+1:km ) = hlw(:,1:km-nkc)
882 hlw_can(:, 3 ) = hlw(:,1 )
883 hlw_can(:, 2 ) = hlw(:,1 )
884 hlw_can(:, 1 ) = hlw(:,1 )
885
886! Subset combined layers (minus top nkc layers)
887 do k = km-nkc, 1, -1 ! top to 1hy model layer
888 ! km is top combined subset
889 ! nkc+1 is bot combined
890 kc= nkc+k ! 4th from top (nkt) to nkc+1 combined canopy plus resolved model layer
891
892! Above-canopy TKE tracer set to hybrid model layers ("resolved_to_canopy" only does mass tranfer of mass tracers)
893 q1_can(:,kc, ntke ) = q1(:,k, ntke ) ! ntke always on
894
895! Above-canopy wind components set to hybrid model layers
896! u1_can (:,kc) = u1 (:,k)
897! v1_can (:,kc) = v1 (:,k)
898
899! ws1_can(:,kc) = sqrt(u1_can(:,kc)**2+v1_can(:,kc)**2)
900
901 end do
902
903 do kc = 1, nkc ! 3-nkc canopy layers
904
905! Sub-canopy values of TKE set to 1hy model layer
906! ("canopy_transfer" only does mass transfer to mass conc. tracers)
907 q1_can(:,kc, ntke ) = q1(:,1, ntke ) ! ntke always on
908
909!Sub-canopy values of wind components set to 1hy model layer
910! u1_can (:,kc) = u1 (:,1)
911! v1_can (:,kc) = v1 (:,1)
912
913! ws1_can(:,kc) = sqrt(u1_can(:,kc)**2+v1_can(:,kc)**2)
914 end do
915
916!
919 CALL satmedmfvdifq_can(im,km,ntrac,ntcw,ntrw,
920 & ntiw,ntke,grav,pi,rd,cp,rv,hvap,hfus,fv,eps,epsm1,
921!The following three variables are for SA-3D-TKE
922 & def_1,def_2,def_3,sa3dtke,dum3d_h,dum3d_e,
923 & dv_can,du_can,tdt_can,rtg_can, ! InOut
924 & u1_can,v1_can,t1_can, q1_can, ! In: canopy inputs
925 & usfco,vsfco,use_oceanuv,
926 & swh_can,hlw_can, ! In: canopy inputs
927 & xmu,garea,zvfun,sigmaf,
928 & psk,rbsoil,zorl,u10m,v10m,fm,fh,
929 & tsea,heat,evap,stress,spd1,
930 & kpbl_can, ! Out
931 & prsi_can,del_can,prsl_can,prslk_can,phii_can,phil_can, ! In: canopy inputs
932 & delt,tte_edmf,
933 & dspheat,
934 & dusfc_can,dvsfc_can,dtsfc_can,dqsfc_can,hpbl_can, ! Out
935 & dkt, dku, tkeh_can, ! Out/Out:tkeh_can
936 & dkt_can,dku_can, ! In: canopy inputs
937 & kinver,xkzm_m,xkzm_h,xkzm_s,dspfac,bl_upfr,bl_dnfr,
938 & rlmx,elmx,sfc_rlm,tc_pbl,use_lpt,
939 & do_canopy, cplaqm,
940 & ntqv, dtend_can, dtidx,index_of_temperature,index_of_x_wind, !inout: dtend (.ldiag3d.)
941 & index_of_y_wind,index_of_process_pbl,gen_tend,ldiag3d,
942 & errmsg,errflg)
943
944! Set non-canopy columns to resolved values
945! NB. Only vars not ALREADY defined in non-canopy columns
946 do i = 1,im
947 IF (frt_mask(i) <= 0.) THEN
948
949! Non-canopy columns
950 kpbl_can(i) = kpbl(i) ! kpbl zero before original sat call
951 hpbl_can(i) = hpbl(i)
952
953! Non-canopy columns
954 dusfc_can(i) = dusfc(i) ! dusfc is zero before the main sat call
955 dvsfc_can(i) = dvsfc(i) ! dvsfc is zero before the main sat call
956 dtsfc_can(i) = dtsfc(i) ! dtsfc is zero before the main sat call
957 dqsfc_can(i) = dqsfc(i) ! dqsfc is zero before the main sat call
958
959 END IF !(FRT_MASK)
960 end do
961
962c
964! NB. before doing "canopy_to_resolved" mass transfer
965c
966
967!U-wind/V-wind on original model layers after diffusion
968 u2(:,1:km) = u1(:,1:km) +
969 & du(:,1:km) * dt2 ! before "canopy_to_resolved"
970 v2(:,1:km) = v1(:,1:km) +
971 & dv(:,1:km) * dt2 ! before "canopy_to_resolved"
972
973 ws2(:,1:km) = sqrt(u2(:,1:km)**2+v2(:,1:km)**2)
974 duv(:,1:km) = (ws2(:,1:km) - ws1(:,1:km) )*rdt ! before "canopy_to_resolved"
975
976! Air temperature on original model layers after diffusion
977 t2(:,1:km) = t1(:,1:km) +
978 & tdt(:,1:km) * dt2 ! before "canopy_to_resolved"
979
980! All mass tracers (excepr TKE) on original model layers after diffusion for use in "canopy_to_resolved"
981 q2(:,1:km, 1:ntrac1) = q1(:,1:km, 1:ntrac1) +
982 & rtg(:,1:km, 1:ntrac1) * dt2 ! before "canopy_to_resolved"
983
984! TKE tracers on original model layers after diffusion for use in "canopy_to_resolved"
985 q2(:,1:km, ntke ) = q1(:,1:km, ntke ) +
986 & rtg(:,1:km, ntke ) * dt2 ! before "canopy_to_resolved"
987 q2(:,:, ntke) = max(q2(:,:, ntke), tkmin) ! before "canopy_to_resolved"
988
989! Temperature & wind components
990 t2_mod(:,1:km) = t2(:,1:km)
991 u2_mod(:,1:km) = u2(:,1:km)
992 v2_mod(:,1:km) = v2(:,1:km)
993
994! Wind speed
995 ws2_mod(:,1:km) = sqrt(u2_mod(:,1:km)**2+v2_mod(:,1:km)**2) ! after diffusion
996
997! Tracers
998 q2_mod(:,1:km,:) = q2(:,1:km, :) ! before "canopy_to_resolved"
999 tkeh_mod(:,1:km) = tkeh(:,1:km) ! before "canopy_to_resolved"
1000
1001! Tendencies
1002 rtg_mod(:,1:km, :) = rtg(:,1:km, :) ! before "canopy_to_resolved" (km)
1003
1004 tdt_mod(:,1:km) = tdt(:,1:km) ! before "canopy_to_resolved" (km)
1005
1006 du_mod(:,1:km) = du(:,1:km) ! before "canopy_to_resolved" (km)
1007 dv_mod(:,1:km) = dv(:,1:km) ! before "canopy_to_resolved" (km)
1008 duv_mod(:,1:km) = duv(:,1:km) ! before "canopy_to_resolved" (km)
1009
1010! Air Density after diffusion model layers
1011 rho2(:,1:km) = prsl(:,1:km)/
1012 & (rd*t2(:,1:km)*
1013 & (1.+fv*max(q2(:,1:km, ntqv),qmin))) ! ntqv=1 before "canopy_to_resolved"
1014
1015! Set non-canopy columns to resolved values
1016! NB. Only vars not ALREADY defined in non-canopy columns
1017
1018! Above-canopy layers (canopy layers below) non-canopy columns
1019 do k = 1, km-nkc
1020 ! km is top combined subset
1021 ! nkc+1 is bot combined
1022 kc= nkc+k ! 4th from top (nkt) to nkc+1 combined canopy plus resolved model layer
1023
1024 do i = 1,im
1025 IF (frt_mask(i) <= 0.) THEN
1026
1027! Above-canopy U-Wind/V-Wind/Temp non-canopy columns
1028 u2_can3(i,kc) = u2(i,k)
1029 v2_can3(i,kc) = v2(i,k)
1030 t2_can3(i,kc) = t2(i,k)
1031
1032! Above-canopy TKE tracer half layers non-canopy columns
1033 tkeh_can(i,kc) = tkeh(i,k) ! J s-1 after diffusion (km)
1034
1035! Tendencies
1036 tdt_can(i,kc) = tdt(i,k) ! K s-1 after diffusion (km)
1037 du_can(i,kc) = du(i,k) ! m s-2 after diffusion (km)
1038 dv_can(i,kc) = dv(i,k) ! m s-2 after diffusion (km)
1039 duv_can(i,kc) = duv(i,k) ! m s-2 after diffusion (km)
1040
1041! All tendencies except TKE non-canopy columns
1042 rtg_can(i,kc, 1:ntrac1) = rtg(i,k, 1:ntrac1) ! kg kg-1 s-1
1043! TKE Tendency non-canopy columns
1044 rtg_can(i,kc, ntke) = rtg(i,k, ntke ) ! s-1
1045
1046 END IF ! (FRT_MASK)
1047 end do ! i=1,im
1048 end do ! k = 1, km-nkc
1049
1050! Canopy layers non-canopy columns
1051 do kc = 1, nkc ! 3-nkc canopy layers
1052 do i = 1,im
1053 IF (frt_mask(i) <= 0.) THEN
1054
1055! Canopy layers U-wind/V-wind/Temp non-canopy columns
1056 u2_can3(i,kc) = u2(i,1)
1057 v2_can3(i,kc) = v2(i,1)
1058 t2_can3(i,kc) = t2(i,1)
1059
1060! Canopy layers TKE half layers non-canopy columns
1061 tkeh_can(i,kc) = tkeh(i,1) ! J S-1 (km)
1062
1063! Canopy layers Tendencies momentum and heat non-canopy columns
1064 tdt_can(i,kc) = tdt(i,1) ! K s-1 (km)
1065 du_can(i,kc) = du(i,1) ! m s-2 (km)
1066 dv_can(i,kc) = dv(i,1) ! m s-2 (km)
1067 duv_can(i,kc) = duv(i,1) ! m s-2 (km)
1068
1069! Canopy layers Tendencies mass tracers non-canopy columns
1070 rtg_can(i,kc, 1:ntrac1) = rtg(i,1, 1:ntrac1) ! kg kg-1 s-1 (km)
1071! Canopy layers Tendency TKE tracer non-canopy columns
1072 rtg_can(i,kc, ntke ) = rtg(i,1, ntke) ! J s-1 s-1 (km)
1073
1074 ENDIF ! (FRT_MASK)
1075 end do ! do i=1,im
1076 end do ! kc=1,nkc
1077
1078! Combined layers (1,km) Humidity after diffusion
1079! Apply minimum value on humidity qmin before "canopy_to_resolved"
1080 q2_can(:,:, ntqv) = q1_can(:,:, ntqv) +
1081 & rtg_can(:,:, ntqv) * dt2 ! after diffusion (km)
1082 q2_can(:,:, ntqv) = max(q2_can(:,:, ntqv),qmin) ! after diffusion (km)
1083
1084! Winds & temperature on combined layers after diffusion
1085
1086! Top 3 combined layers NOT calculated in "sat_can": fill with top 3 model layers
1087 u2_can3(:,km+nkc) = u2(:,km ) ! after diffusion (nkt)
1088 u2_can3(:,km+2 ) = u2(:,km-1) ! after diffusion (nkt)
1089 u2_can3(:,km+1 ) = u2(:,km-2) ! after diffusion (nkt)
1090
1091! Top 3 combined layers NOT calculated in "sat_can": fill with top 3 model layers
1092 v2_can3(:,km+nkc) = v2(:,km ) ! after diffusion (nkt)
1093 v2_can3(:,km+2 ) = v2(:,km-1) ! after diffusion (nkt)
1094 v2_can3(:,km+1 ) = v2(:,km-2) ! after diffusion (nkt)
1095
1096! Top 3 combined layers NOT calculated in "sat_can": fill with top 3 model layers
1097 t2_can3(:,km+nkc) = t2(:,km ) ! after diffusion (nkt)
1098 t2_can3(:,km+2 ) = t2(:,km-1) ! after diffusion (nkt)
1099 t2_can3(:,km+1 ) = t2(:,km-2) ! after diffusion (nkt)
1100
1101! All tracers on combined layers after diffusion, for use in "canopy_to_resolved"
1102
1103! Top 3 combined layers NOT calculated in "sat_can": fill with top 3 model layers
1104 q2_can3(:,km+nkc , ntke ) = q2(:,km , ntke ) ! after diffusion (nkt)
1105 q2_can3(:,km+2 , ntke ) = q2(:,km-1 , ntke ) ! after diffusion (nkt)
1106 q2_can3(:,km+1 , ntke ) = q2(:,km-2 , ntke ) ! after diffusion (nkt)
1107
1108! Top 3 combined layers NOT calculated in "sat_can": fill with top 3 model layers
1109 q2_can3(:,km+nkc , 1:ntrac1) = q2(:,km , 1:ntrac1) ! nkt after diffusion (nkt)
1110 q2_can3(:,km+2 , 1:ntrac1) = q2(:,km-1 , 1:ntrac1) ! nkt after diffusion (nkt)
1111 q2_can3(:,km+1 , 1:ntrac1) = q2(:,km-2 , 1:ntrac1) ! nkt after diffusion (nkt)
1112
1113! Wind and temperature after diffusion
1114
1115! Combined layers (1,km) U-Wind/V-Wind/Temp after diffusion
1116 u2_can3(:,1:km) = u1_can(:,1:km) + du_can(:,1:km) * dt2 ! after diffusion (nkt)
1117 v2_can3(:,1:km) = v1_can(:,1:km) + dv_can(:,1:km) * dt2 ! after diffusion (nkt)
1118
1119! Combined layers (1,km) Temperature after diffusion
1120 t2_can3(:,1:km) = t1_can(:,1:km) + tdt_can(:,1:km) * dt2 ! after diffusion (nkt)
1121! t2_can (:,1:km) = t1_can (:,1:km) + tdt_can(:,1:km) * dt2 ! after diffusion (km)
1122
1123! Wind Speed after diffusion on canopy layers
1124 ws2_can3(:,1:km) = sqrt(u2_can3(:,1:km)**2+v2_can3(:,1:km)**2)
1125
1126! 10m-interpolated wind ws1_can3
1127! duv_can (:,1:km) = (ws2_can3(:,1:km) - ws1_can3(:,1:km)) * rdt ! 10m-interpolated ws1_can3
1128! 1hy model layer wind ws1_can (1-4cy)
1129 duv_can(:,1:km) = (ws2_can3(:,1:km) - ws1_can(:,1:km)) * rdt ! ws1_can is using 1hy model layer u1&v1(:,1)
1130
1131! Tracers after diffusion
1132! Subset (1:km) combined layers out of total ntk layers (NB. dim(:,nkt,:) <= dim(:,km,:)
1133 q2_can3(:,1:km, 1:ntrac1) = q1_can(:,1:km, 1:ntrac1) +
1134 & rtg_can(:,1:km, 1:ntrac1) * dt2 ! after diffusion (nkt)
1135! TKE after diffusion
1136 q2_can3(:,1:km, ntke ) = q1_can(:,1:km, ntke ) +
1137 & rtg_can(:,1:km, ntke ) * dt2 ! after diffusion (nkt)
1138 q2_can(:,1:km, ntke ) = q1_can(:,1:km, ntke ) +
1139 & rtg_can(:,1:km, ntke ) * dt2 ! after diffusion (km)
1140
1141! Apply minimum value on TKE tracer before "canopy_to_resolved"
1142 q2_can3(:,:, ntke) = max(q2_can3(:,:, ntke), tkmin) ! after diffusion (nkt)
1143 q2_can(:,:, ntke) = max(q2_can(:,:, ntke), tkmin) ! after diffusion (km)
1144
1145! Apply minimum value on humidity qmin before "canopy_to_resolved" and the tendency update
1146 q2_can3(:,:, ntqv) = max(q2_can3(:,:, ntqv), qmin ) ! ntqv=1
1147
1148! Apply minimum value on chemical conc before "canopy_to_resolved"
1149 q2_can3(:,:, ntoz:ntrac1 ) =
1150 & max(q2_can3(:,:, ntoz:ntrac1), concmin)
1151
1152! Top 3 combined layers set to resolved
1153! NB. Q2_can3 tracers array & t2_can3 after diffusion only updated 1:km
1154 rho2_can3(:,km+nkc) = prsl(:,km )/ ! after diffusion
1155 & (rd*t2(:,km )*
1156 & (1.+fv*max(q2(:,km , ntqv),qmin))) ! ntqv=1
1157 rho2_can3(:,km+2 ) = prsl(:,km-1)/ ! after diffusion
1158 & (rd*t2(:,km-1)*
1159 & (1.+fv*max(q2(:,km-1, ntqv),qmin))) ! ntqv=1
1160
1161 rho2_can3(:,km+1 ) = prsl(:,km-2)/ ! after diffusion
1162 & (rd*t2(:,km-2)*
1163 & (1.+fv*max(q2(:,km-2, ntqv),qmin))) ! ntqv=1
1164
1165! Air density after diffusion on canopy layers
1166 rho2_can3(:,1:km ) = prsl_can3(:,1:km)/ ! after diffusion
1167 & (rd*t2_can3(:,1:km)*
1168 & (1.+fv*max(q2_can3(:,1:km, ntqv),qmin))) ! ntqv=1
1169
1170 endif !do_canopy .and. cplaqm
1171
1172 if (do_canopy) then
1173
1174!===============================================================================
1175! Gather tracer concentration from canopy layers into model resolved layers (flag = 1)
1176!===============================================================================
1177
1178 CALL canopy_transfer_run(im, km, nkc, nkt, !in
1179 & ntrac1, ntoz, !in
1180 & garea, !in
1181 & zi, zl, zm, !in
1182 & q2, rho2, !in kg kg-1
1183 & 1, !in 1 = "canopy_to_resolved"
1184 & frt_mask, !in
1185 & kmod, kcan3, !in
1186 & zmom_can3, zmid_can3, !in
1187 & prsl_can3, rho2_can3, !in
1188 & q2_mod, q2_can3, q2_2m, !inout kg kg-1 after diffusion
1189 & massair_can3, massair, !inout
1190 & mmr_o3_can3, !inout
1191 & nfrct, ifrct, !inout
1192 & frctr2c, frctc2r, !inout
1193 & errmsg, errflg )
1194
1195 if (errflg /= 0) return
1196
1197! ==============
1198! Input:
1199! Q2_CAN3(:,:, NLAYT, NSPCSD) : Chemical tracers mass conc. kg kg-1 on combined canopy+resolved layers after diffusion
1200!
1201! InOutput
1202! Q2_MOD(:,:, NLAYS, NSPCSD) : Chemical tracers mass conc. kg kg-1 on model levels after diffusion
1203! Q2 Canopy columns only!!!
1204!
1205! ================================
1206
1207
1208!
1209! Canopy columns calculated in "canopy_to_resolved"
1210!
1211! Non-Canopy columns filled with resolved values
1212!
1213 do k = 1, km
1214 do i = 1,im
1215
1216! Non-canopy columns set to resolved values
1217 IF (frt_mask(i) <= 0.) THEN
1218
1219! Non-canopy columns after "canopy_to_resolved" Line 1843: forrtl: error (78): process killed (SIGTERM)
1220 q2_mod(i,k, 1:ntrac1) = q2(i,k, 1:ntrac1) ! after diffusion & after "canopy_to_resolved"
1221 q2_mod(i,k, ntke ) = q2(i,k, ntke) ! after diffusion
1222
1223 t2_mod(i,k) = t2(i,k) ! after diffusion
1224 u2_mod(i,k) = u2(i,k) ! after diffusion
1225 v2_mod(i,k) = v2(i,k) ! after diffusion
1226
1227 ws2_mod(i,k) = sqrt(u2_mod(i,k)**2+v2_mod(i,k)**2) ! after diffusion
1228
1229 END IF ! (FRT_MASK)
1230 end do ! i=1,im
1231 end do ! k = 1, km-nkc
1232
1233! Apply minimum value on chemical conc after "canopy_to_resolved"
1234 q2_mod(:,:, ntchs-1:ntrac1 ) =
1235 & max(q2_mod(:,:, ntchs-1:ntrac1), concmin) ! after "canopy_to_resolved"
1236
1237! Apply minimum value on humidity qmin before doing canopy_transfer & update tendency
1238 q2_mod(:,:, ntqv) = max(q2_mod(:,:, ntqv),qmin)
1239
1240!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1241 do k = 1, km1-1 ! from bottom to top resolved model levels
1242
1243 kount=0
1244 zooox(:) = 1.
1245 do i = 1, im
1246
1247! Canopy columns/grid cells
1248! There is a contiguous forest canopy, calculated integrated canopy correction over model layers
1249 IF ( frt_mask(i) > 0. ) THEN
1250
1251 fch = cfch(i)
1252
1253! Determine if canopy inside the model layer (kcan=1) or not (kcan=0)
1254 IF (k .EQ. 1) THEN !use model layer interfaces
1255 kcan = 1
1256 ELSE
1257 IF ( cfch(i) > zi(i,k) .AND. cfch(i) <= zi(i,k+1) ) THEN
1258 kcan = 1
1259 ELSE
1260 kcan= 0
1261 END IF
1262 END IF
1263
1264 IF (kcan == 1 ) THEN !canopy inside model layer k
1265
1266 zi05 = zi(i,k+1) ! Initialize each model layer top that contains canopy (m)
1267 ! Integrate across total model interface
1268 zfl = zi(i,k+1) ! Set ZFL = ZI05
1269 countcan = 0 ! Initialize canopy layers
1270
1271 IF (k .EQ. 1) THEN !Find bottom in each model layer
1272 botcan = 0.0
1273 ELSE
1274 botcan = zi(i,k)
1275 END IF
1276
1277 dzfl = zi(i,k+1) - botcan
1278
1279! Canopy inside model layer
1280 DO WHILE (zi05 .GE. botcan)
1281
1282!!!!!!!!!!!!!!!!!!!!!!!!
1283! Steping down in-canopy zi_can3(:,1)= 0.
1284
1285 do kc = nkt-1, 1, -1 ! from top to bottom combined canopy plus resolved model layers
1286
1287!!!!!!!!!!!!!!!!!!!!!!!!
1288
1289! Between two canopy layers
1290 IF ( zi05 > zi_can3(i,kc ) .and.
1291 & zi05 <= zi_can3(i,kc+1) ) THEN
1292
1293 kc_can = kc
1294 fzi05 = 1./
1295! max(zi_can3 (i,kc+1) -
1296! & zi_can3 (i,kc), epsilon) *
1297 & (zi_can3(i,kc+1) - zi_can3(i,kc)) *
1298 & (zi05 - zi_can3(i,kc))
1299
1300 ttcorr = (t2_can3(i,kc+1) - t2_can3(i,kc))/
1301! & max(zi_can3 (i,kc+1) -
1302! & zi_can3 (i,kc), epsilon) *
1303 & (zi_can3(i,kc+1) - zi_can3(i,kc)) *
1304 & (zi05 - zi_can3(i,kc))
1305
1306! U-Wind/V-Wind after diffusion on canopy layers
1307 uucan = u2_can3(i,kc ) +
1308 & (u2_can3(i,kc+1) - u2_can3(i,kc))/
1309! & max(zi_can3 (i,kc+1) -
1310! & zi_can3 (i,kc), epsilon) *
1311 & (zi_can3(i,kc+1) - zi_can3(i,kc)) *
1312 & (zi05 - zi_can3(i,kc))
1313
1314 vvcan = v2_can3(i,kc ) +
1315 & (v2_can3(i,kc+1) - v2_can3(i,kc))/
1316! & max(zi_can3 (i,kc+1) -
1317! & zi_can3 (i,kc), epsilon) *
1318 & (zi_can3(i,kc+1) - zi_can3(i,kc)) *
1319 & (zi05 - zi_can3(i,kc))
1320! Temperature after diffusion on canopy layers
1321 ttcan = t2_can3(i,kc ) +
1322 & (t2_can3(i,kc+1) - t2_can3(i,kc))/
1323! & max(zi_can3 (i,kc+1) -
1324! & zi_can3 (i,kc), epsilon) *
1325 & (zi_can3(i,kc+1) - zi_can3(i,kc)) *
1326 & (zi05 - zi_can3(i,kc))
1327
1328
1329! TKE after diffusion on canopy layers
1330 tkecan = q2_can3(i,kc, ntke) +
1331 & (q2_can3(i,kc+1, ntke) -
1332 & q2_can3(i,kc , ntke))/
1333! & max(zi_can3 (i,kc+1) -
1334! & zi_can3 (i,kc), epsilon) *
1335 & (zi_can3(i,kc+1) - zi_can3(i,kc)) *
1336 & (zi05 - zi_can3(i,kc))
1337
1338! TKE half layers
1339 tkehcan = tkeh_can(i,kc ) +
1340 & (tkeh_can(i,kc+1) - tkeh_can(i,kc))/
1341! & max(zi_can3(i,kc+1) -
1342! & zi_can3(i,kc), epsilon) *
1343 & ( zi_can3(i,kc+1) - zi_can3(i,kc)) *
1344 & (zi05 - zi_can3(i,kc))
1345
1346 END IF ! "zi_can3(kc) < ZI05 <= zi_can3(kc+1)"
1347
1348!!!!!!!!!!!!!!
1349! End steping down in-canopy
1350!
1351 end do ! kc = nkt-1, 1, -1 ! from top to bottom combined canopy plus resolved model layers
1352!
1353!!!!!!!!!!!!!!
1354
1355! IF ( ZI05 .LE. FCH ) THEN ! in-canopy layers (before Dec 2025)
1356 IF ( zi05 .LE. zfl ) THEN ! Model layers
1357
1358 countcan = countcan + 1
1359 zcanx(countcan) = zi05
1360
1361! U-Wind/V-Wind after diffusion on model layer
1362 uux(countcan) = uucan
1363 vvx(countcan) = vvcan
1364! Temperature after diffusion on model layers
1365 ttx(countcan) = ttcan
1366! TKE on model layers
1367 tkex(countcan) = tkecan
1368! TKE half layers on model layers
1369 tkehx(countcan) = tkehcan
1370
1371 END IF ! ( ZI05 .LE. ZFL )
1372
1373 zi05 = zi05-0.5 !step down in-canopy resolution of 0.5m
1374
1375 END DO ! DO WHILE (ZI05.GE.BOTCAN)
1376
1377! IntegrateTrapezoid = sum((y(1+1:n-0) + y(1+0:n-1))*
1378! & (x(1+1:n-0) - x(1+0:n-1)))/2
1379 zz_int= integratetrapezoid(
1380 & zcanx(countcan:1:-1), zooox(countcan:1:-1))
1381
1382! U-wind
1383 uu_int = integratetrapezoid(
1384 & zcanx(countcan:1:-1) ,
1385 & uux(countcan:1:-1) ) /
1386! & max(ZZ_INT, epsilon)
1387 & zz_int
1388
1389! V-wind
1390 vv_int = integratetrapezoid(
1391 & zcanx(countcan:1:-1) ,
1392 & vvx(countcan:1:-1) ) /
1393! & max(ZZ_INT, epsilon)
1394 & zz_int
1395
1396! Temp
1397 tt_int = integratetrapezoid(
1398 & zcanx(countcan:1:-1) ,
1399 & ttx(countcan:1:-1) ) /
1400! & max(ZZ_INT, epsilon)
1401 & zz_int
1402
1403! TKE
1404 tke_int = integratetrapezoid(
1405 & zcanx(countcan:1:-1) ,
1406 & tkex(countcan:1:-1) ) /
1407! & max(ZZ_INT, epsilon)
1408 & zz_int
1409
1410! TKEH
1411 tkeh_int= integratetrapezoid(
1412 & zcanx(countcan:1:-1) ,
1413 & tkehx(countcan:1:-1) ) /
1414! & max(ZZ_INT, epsilon)
1415 & zz_int
1416
1417! Sum
1418 tt_sum = sum( ttx(countcan:1:-1))/countcan
1419
1420! U-wind/V-wind Canopy Columns
1421 u2_mod(i,k) = uu_int ! after "canopy-to-resolved"
1422 v2_mod(i,k) = vv_int ! after "canopy-to-resolved"
1423
1424 ws2_mod(i,k) = sqrt(u2_mod(i,k)**2+v2_mod(i,k)**2) ! after "canopy-to-resolved"
1425
1426! Temperature Canopy Columns
1427 t2_mod(i,k) = tt_int ! after "canopy-to-resolved"
1428
1429! TKE Canopy Columns
1430 q2_mod(i,k, ntke) = tke_int ! after "canopy-to-resolved"
1431 q2_mod(i,k, ntke) = max(q2_mod(i,k, ntke), tkmin) ! after "canopy-to-resolved"
1432
1433! TKEH Canopy Columns
1434 tkeh_mod(i,k) = tkeh_int ! after "canopy-to-resolved"
1435 tkeh_mod(i,k) = max(tkeh_mod(i,k), tkmin)
1436! Apply minimum value on TKE tracer before "canopy_to_resolved"
1437
1438 END IF ! (KCAN .EQ. 1) model layer(s) containing canopy
1439
1440 END IF ! contiguous canopy conditions
1441 enddo ! i = 1, im
1442
1443 kount = kount + 1
1444
1445 enddo ! k = 1, km1-1 ! from bottom to top resolved model levels
1446
1447!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1448! Update ALL tracers with in-canopy tendencies (average sub-canopy values )
1449! Here just wind components, temperature TKE, and interstitial tracers
1450! The chemical tracers are updated below with values from "canopy_to_resolved" transfer
1451!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1452
1453! Update ALL model layer
1454 do k = 1,km
1455
1456 do i = 1,im
1457
1458 ! Update Canopy columns only
1459 IF (frt_mask(i) > 0.) THEN
1460
1461! Canopy Columns
1462! U-Wind/V-Wind after sub-canopy diffusion
1463 du_mod(i,k) = (u2_mod(i,k) - u1(i,k) )*rdt ! after "canopy_to_resolved" dim(km)
1464 dv_mod(i,k) = (v2_mod(i,k) - v1(i,k) )*rdt ! after "canopy_to_resolved" dim(km)
1465
1466 duv_mod(i,k) = (ws2_mod(i,k) - ws1(i,k) )*rdt ! after "canopy_to_resolved" dim(km)
1467
1468! Temperature after sub-canopy diffusion
1469 tdt_mod(i,k) = (t2_mod(i,k) - t1(i,k) )*rdt ! after "canopy_to_resolved" dim(km)
1470
1471! Tendency mass tracers after sub-canopy diffusion
1472 rtg_mod(i,k, 1:ntrac1) = (q2_mod(i,k, 1:ntrac1) -
1473 & q1(i,k, 1:ntrac1))*rdt ! after "canopy_to_resolved"
1474! Tendency TKE tracer after sub-canopy diffusion
1475 rtg_mod(i,k, ntke) = (q2_mod(i,k, ntke) -
1476 & q1(i,k, ntke) )*rdt ! after "canopy_to_resolved"
1477
1478!!!!!!!!!!!!!!!!!!!!!
1479
1480 ENDIF ! (FRT_MASK)
1481 enddo ! do i=1,im
1482 enddo ! k=1,km
1483
1484!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1485! Update Met & TKE & MP (microphysics) cloud fields
1486!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1487! Update all model layers
1488 do k = 1,km
1489 do i = 1,im
1490
1491 IF (frt_mask(i) > 0.) THEN
1492! Tendency Wind components
1493! du (i,k) = du_mod (i,k) ! <<<<<<<========== UPDATE UU TEND =========>>>>>>>
1494! dv (i,k) = dv_mod (i,k) ! <<<<<<<========== UPDATE VV TEND =========>>>>>>>
1495! Tendency Temperature
1496! tdt(i,k) = tdt_mod(i,k) ! <<<<<<<========== UPDATE TT TEND =========>>>>>>>
1497! Tendency TKE (ntke=198)
1498! rtg(i,k, ntke) = rtg_mod(i,k, ntke) ! <<<<<<<========== UPDATE TKE TEND =========>>>>>>>
1499
1500! TKE half layers
1501! tkeh(i,k) = tkeh_mod(i,k) ! <<<<<<<========== UPDATE TKEH =========>>>>>>>
1502
1503 ENDIF ! Contiguous canopy
1504 enddo ! i
1505 enddo !k
1506
1507! cloud/rain
1508! ------------------------
1509! n=1 (ntqv)
1510! n=1 (ntcw)
1511! n=3 ...
1512! n=7 "o3mr"
1513! -------------
1514!NTOZ do n = 1, ntoz
1515! do k = 1,km
1516! do i = 1,im
1517! IF (FRT_MASK(i) > 0.) THEN
1518! Humidity & Clouds
1519! rtg(i,k, n) = rtg_mod(i,k, n) ! <<<<<<<========== UPDATE MET TEND =========>>>>>>>
1520! rtg(i,k, ntqv) = rtg_mod(i,k, ntqv) ! <<<<<<<========== UPDATE VAP TEND =========>>>>>>>
1521! ENDIF ! Contiguous canopy
1522! enddo ! i
1523! enddo !k
1524!NTOZ enddo !n
1525
1526
1527
1528!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1529! Update ONLY chemical tracers (n=8, ntrac1)
1530!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1531
1532! All chemical tracers (n= ntchs-1 , ntche-1 (ntrac1)
1533! ntqv = 1-7 cloud/rain and "sgs_tke"
1534! sgs_tke = 8 N/A
1535!
1536! ntchs-1 = 8 "no2" GFS_v16 GFDL MP
1537! ntchs = 9 "no" GFS_v16 GFDL MP
1538! n = 10 "o3" GFS_v16 GFDL MP
1539! ...
1540! ntche-1= 196 NTRAC1 GFS_v16 GFDL MP
1541! ntche = 197 "cld_amt" GFS_v16 GFDL MP
1542! ntke = 198 "tke" GFS_v16 GFDL MP
1543! ---------------------------------------
1544! ntqv = 1-9 cloud/rain
1545! sgs_tke = 10 "sgs_tke" N/A
1546!
1547! ntchs-1 = 10 "no2" GFS_v17p8
1548! ntchs = 11 "no" GFS_v17p8
1549! n = 12 "o3" GFS_v17p8
1550! ...
1551! ntche-1 = 198 NTRAC1 GFS_v17p8
1552! ntke = 199 "tke" GFS_v17p8
1553! ---------------------------------------
1554 DO n = ntchs-1, ntche-1 ! 10, ntche-1 (same as NTRAC1)
1555
1556! Update all model layers
1557 do k = 1,km
1558
1559 do i = 1,im
1560 IF (frt_mask(i) > 0.) THEN
1561
1562!!!!!!!!!!!!!!!!!!!!!!!!!
1563! Select all chemical tracers (selected tracers above)
1564!!!!!!!!!!!!!!!!!!!!!!!!!
1565 rtg(i,k, n) = rtg_mod(i,k, n) ! <<<<<<<========== UPDATE TEND =========>>>>>>>
1566 END IF
1567 end do
1568
1569 end do ! k = 1,km
1570 end do ! n = 1, NTRAC1
1571
1572 endif !if(do_canopy)
1573
1574!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1575
1576 return
1577 end subroutine canopy_driver_run
1579 end module canopy_driver
satmedmfvdifq_can_mod::satmedmfvdifq_can
subroutine satmedmfvdifq_can(im, km, ntrac, ntcw, ntrw,
Definition satmedmfvdifq_can.F:51
canopy_driver
This file contains ...
Definition canopy_driver.F:4
mfpbltq_mod
This module contains the subroutine that calculates mass flux and updraft parcel properties for therm...
Definition mfpbltq.f:9
mfscuq_mod
This module contains the mass flux and downdraft parcel properties parameterization for stratocumulus...
Definition mfscuq.f:6
satmedmfvdifq_can_mod
This file contains the CCPP-compliant SATMEDMF scheme (updated version) which computes subgrid vertic...
Definition satmedmfvdifq_can.F:10
tridi_mod
This module contains routine to compute tridiagonal matrix elements for TKE, heat,...
Definition tridi.f:6