18 use machine,
only : kind_phys
22 real(kind=kind_phys) :: pi, pi2, pih, rad_to_deg, deg_to_rad
23 real(kind=kind_phys) :: arad, p0s
24 real(kind=kind_phys) :: grav, grav2, rgrav, rgrav2
25 real(kind=kind_phys) :: cpd, rd, rv, fv
26 real(kind=kind_phys) :: rdi, rcpd, rcpd2
28 real(kind=kind_phys) :: gor, gr2, grcp, gocp, rcpdl, grav2cpd
29 real(kind=kind_phys) :: bnv2min, bnv2max
30 real(kind=kind_phys) :: dw2min, velmin, minvel
31 real(kind=kind_phys) :: omega1, omega2, omega3
32 real(kind=kind_phys) :: hpscale, rhp, rhp2, rh4, rhp4, khp, hpskm
33 real(kind=kind_phys) :: mkzmin, mkz2min, mkzmax, mkz2max, cdmin
34 real(kind=kind_phys) :: rcpdt
65 subroutine init_nazdir(naz, xaz, yaz)
67 use machine,
only : kind_phys
73 real(kind=kind_phys),
dimension(naz) :: xaz, yaz
75 real(kind=kind_phys) :: phic, drad
80 phic = drad*(float(idir)-1.0)
95 end subroutine init_nazdir
103 subroutine init_global_gwdis(levs, zkm, pmb, kvg, ktg, krad, kion, me, master)
105 use machine ,
only : kind_phys
109 integer ,
intent(in) :: me, master
110 integer ,
intent(in) :: levs
111 real(kind=kind_phys),
intent(in) :: zkm(levs), pmb(levs)
112 real(kind=kind_phys),
intent(out),
dimension(levs+1) :: kvg, ktg, krad, kion
117 real(kind=kind_phys),
parameter :: vusurf = 2.e-5
118 real(kind=kind_phys),
parameter :: musurf = vusurf/1.95
119 real(kind=kind_phys),
parameter :: hpmol = 7.0
121 real(kind=kind_phys),
parameter :: kzmin = 0.1
122 real(kind=kind_phys),
parameter :: kturbo = 100.
123 real(kind=kind_phys),
parameter :: zturbo = 130.
124 real(kind=kind_phys),
parameter :: zturw = 30.
125 real(kind=kind_phys),
parameter :: inv_pra = 3.
127 real(kind=kind_phys),
parameter :: alpha = 1./86400./15.
128 real(kind=kind_phys) :: pa_alp = 750.
129 real(kind=kind_phys) :: tau_alp = 10.
131 real(kind=kind_phys),
parameter :: kdrag = 1./86400./30.
132 real(kind=kind_phys),
parameter :: zdrag = 100.
133 real(kind=kind_phys),
parameter :: zgrow = 50.
135 real(kind=kind_phys) :: vumol, mumol, keddy, ion_drag
136 real(kind=kind_phys) :: rf_fv3, rtau_fv3, ptop, pih_dlog
138 real(kind=kind_phys) :: ae1 ,ae2
142 rtau_fv3 = 1./86400./tau_alp
143 pih_dlog = pih/log(pa_alp/ptop)
147 vumol = vusurf*exp(ae1)
148 mumol = musurf*exp(ae1)
149 ae2 = -((zkm(k)-zturbo) /zturw)**2
150 keddy = kturbo*exp(ae2)
152 kvg(k) = vumol + keddy
153 ktg(k) = mumol + keddy*inv_pra
162 if (pmb(k) .le. pa_alp)
then
163 rf_fv3=rtau_fv3*sin(pih_dlog*log(pa_alp/pmb(k)))**2
164 krad(k) = krad(k) + rf_fv3
165 kion(k) = kion(k) + rf_fv3
187 end subroutine init_global_gwdis
198 use machine ,
only : kind_phys
199 use ugwp_common,
only : bnv2min, grav, grcp, fv, grav, cpd, grcp, pi
206 real(kind=kind_phys),
parameter :: hncrit=9000.
207 real(kind=kind_phys),
parameter :: hminmt=50.
208 real(kind=kind_phys),
parameter :: sigfac=4.0
209 real(kind=kind_phys),
parameter :: hpmax=2500.0
210 real(kind=kind_phys),
parameter :: hpmin=25.0
213 real(kind=kind_phys),
parameter :: minwnd=1.0
214 real(kind=kind_phys),
parameter :: dpmin=5000.0
217 character(len=8) :: strver =
'gfs_2018'
218 character(len=8) :: strbase =
'gfs_2018'
220 real(kind=kind_phys),
parameter :: rimin=-10., ric=0.25
222 real(kind=kind_phys),
parameter :: frmax=10., frc =1.0, frmin =0.01
223 real(kind=kind_phys),
parameter :: ce=0.8, ceofrc=ce/frc, cg=0.5
224 real(kind=kind_phys),
parameter :: gmax=1.0, veleps=1.0, factop=0.5
225 real(kind=kind_phys),
parameter :: efmin=0.5, efmax=10.0
227 real(kind=kind_phys),
parameter :: rlolev=50000.0
228 integer,
parameter :: mdir = 8
229 real(kind=kind_phys),
parameter :: fdir=mdir/(8.*atan(1.0))
230 real(kind=kind_phys),
parameter :: zpgeo=2.*atan(1.0)
233 data nwdir/6,7,5,8,2,3,1,4/
236 real(kind=kind_phys),
parameter :: odmin = 0.1, odmax = 10.0
237 real(kind=kind_phys),
parameter :: fcrit_sm = 0.7, fcrit_sm2 = fcrit_sm * fcrit_sm
238 real(kind=kind_phys),
parameter :: fcrit_gfs = 0.7, fcrit_v1 = 0.7
239 real(kind=kind_phys),
parameter :: fcrit_mtb = 0.7
241 real(kind=kind_phys),
parameter :: zbr_pi = zpgeo
242 real(kind=kind_phys),
parameter :: zbr_ifs = zpgeo
246 real(kind=kind_phys),
parameter :: kxoro=6.28e-3/200.
247 real(kind=kind_phys),
parameter :: coro = 0.0
248 integer,
parameter :: nridge=2
249 real(kind=kind_phys),
parameter :: sigma_std=1./100., gamm_std=1.0
251 real(kind=kind_phys) :: cdmb
252 real(kind=kind_phys) :: cleff
266 real(kind=kind_phys),
parameter :: kmax = 6.28/(10.*25.)
267 real(kind=kind_phys),
parameter :: k1tr = 6.28/(2100)
268 real(kind=kind_phys),
parameter :: kflt = 6.28/(18.e3)
269 real(kind=kind_phys),
parameter :: k0tr = 6.28/(10.e3)
270 real(kind=kind_phys),
parameter :: nk1tr = 2.8
271 real(kind=kind_phys),
parameter :: nk0tr = 1.9
272 real(kind=kind_phys),
parameter :: a1_tofd = kflt ** nk1tr *1.e3
273 real(kind=kind_phys),
parameter :: a2_tofd = k1tr ** (nk0tr-nk1tr)
274 real(kind=kind_phys),
parameter :: fix_tofd = 2.* 0.005 * 12 *0.6
282 integer,
parameter :: n_tofd = 2
283 real(kind=kind_phys),
parameter :: const_tofd = 0.0759
284 real(kind=kind_phys),
parameter :: ze_tofd = 1500.0
285 real(kind=kind_phys),
parameter :: a12_tofd = 0.0002662*0.005363
286 real(kind=kind_phys),
parameter :: ztop_tofd = 3.*ze_tofd
292 subroutine init_oro_gws(nwaves, nazdir, nstoch, effac, &
296 integer :: nwaves, nazdir, nstoch
298 real(kind=kind_phys) :: cdmbx
299 real(kind=kind_phys) :: kxw
300 real(kind=kind_phys) :: effac
305 real(kind=kind_phys),
parameter :: lonr_refmb = 4.0 * 192.0
306 real(kind=kind_phys),
parameter :: lonr_refgw = 192.0
307 real(kind=kind_phys),
parameter :: cleff_ref = 0.5e-5
315 cdmbx = lonr_refmb/float(lonr)
318 cleff = cleff_ref * sqrt(lonr_refgw/float(lonr))
320 end subroutine init_oro_gws
331 use machine ,
only : kind_phys
335 real(kind=kind_phys) :: eff_con
339 real(kind=kind_phys) :: con_dlength
340 real(kind=kind_phys) :: con_cldf
342 real(kind=kind_phys),
parameter :: cmin = 5
343 real(kind=kind_phys),
parameter :: cmax = 95.
344 real(kind=kind_phys),
parameter :: cmid = 22.5
345 real(kind=kind_phys),
parameter :: cwid = cmid
346 real(kind=kind_phys),
parameter :: bns = 2.e-2, bns2 = bns*bns, bns4=bns2*bns2
347 real(kind=kind_phys),
parameter :: mstar = 6.28e-3/2.
348 real(kind=kind_phys) :: dc
350 real(kind=kind_phys),
allocatable :: ch_conv(:), spf_conv(:)
351 real(kind=kind_phys),
allocatable :: xaz_conv(:), yaz_conv(:)
354 subroutine init_conv_gws(nwaves, nazdir, nstoch, effac, lonr, kxw)
361 integer :: nwaves, nazdir, nstoch
367 real(kind=kind_phys) :: kxw, effac
368 real(kind=kind_phys) :: work1 = 0.5
369 real(kind=kind_phys) :: chk, tn4, snorm
377 con_dlength = pi2*arad/float(lonr)
381 if (.not.
allocated(ch_conv))
allocate (ch_conv(nwaves))
382 if (.not.
allocated(spf_conv))
allocate (spf_conv(nwaves))
383 if (.not.
allocated(xaz_conv))
allocate (xaz_conv(nazdir))
384 if (.not.
allocated(yaz_conv))
allocate (yaz_conv(nazdir))
386 dc = (cmax-cmin)/float(nwaves-1)
393 chk = cmin + (k-1)*dc
396 spf_conv(k) = bns4*chk/(bns4+tn4)
399 snorm = sum(spf_conv)
400 spf_conv = spf_conv/snorm*1.5
402 call init_nazdir(nazdir, xaz_conv, yaz_conv)
403 end subroutine init_conv_gws
414 use machine ,
only : kind_phys
419 real(kind=kind_phys) :: eff_fj
424 real(kind=kind_phys),
parameter :: fjet_trig=0.
427 real(kind=kind_phys),
parameter :: cmin = 2.5
428 real(kind=kind_phys),
parameter :: cmax = 67.5
429 real(kind=kind_phys) :: dc
430 real(kind=kind_phys),
allocatable :: ch_fjet(:) , spf_fjet(:)
431 real(kind=kind_phys),
allocatable :: xaz_fjet(:), yaz_fjet(:)
434 subroutine init_fjet_gws(nwaves, nazdir, nstoch, effac,lonr, kxw)
440 integer :: nwaves, nazdir, nstoch
442 real(kind=kind_phys) :: kxw, effac , chk
451 if (.not.
allocated(ch_fjet))
allocate (ch_fjet(nwaves))
452 if (.not.
allocated(spf_fjet))
allocate (spf_fjet(nwaves))
453 if (.not.
allocated(xaz_fjet))
allocate (xaz_fjet(nazdir))
454 if (.not.
allocated(yaz_fjet))
allocate (yaz_fjet(nazdir))
456 dc = (cmax-cmin)/float(nwaves-1)
458 chk = cmin + (k-1)*dc
462 call init_nazdir(nazdir, xaz_fjet, yaz_fjet)
464 end subroutine init_fjet_gws
473 use machine ,
only : kind_phys
477 real(kind=kind_phys) :: eff_okw
482 real(kind=kind_phys),
parameter :: okw_trig=0.
484 real(kind=kind_phys),
parameter :: cmin = 2.5
485 real(kind=kind_phys),
parameter :: cmax = 67.5
486 real(kind=kind_phys) :: dc
487 real(kind=kind_phys),
allocatable :: ch_okwp(:), spf_okwp(:)
488 real(kind=kind_phys),
allocatable :: xaz_okwp(:), yaz_okwp(:)
494 subroutine init_okw_gws(nwaves, nazdir, nstoch, effac, lonr, kxw)
499 integer :: nwaves, nazdir, nstoch
501 real(kind=kind_phys) :: kxw, effac , chk
510 if (.not.
allocated(ch_okwp))
allocate (ch_okwp(nwaves))
511 if (.not.
allocated(spf_okwp))
allocate (spf_okwp(nwaves))
512 if (.not.
allocated(xaz_okwp))
allocate (xaz_okwp(nazdir))
513 if (.not.
allocated(yaz_okwp))
allocate (yaz_okwp(nazdir))
514 dc = (cmax-cmin)/float(nwaves-1)
516 chk = cmin + (k-1)*dc
521 call init_nazdir(nazdir, xaz_okwp, yaz_okwp)
523 end subroutine init_okw_gws
536 use machine ,
only : kind_phys
539 integer :: nwav, nazd
541 real(kind=kind_phys) :: eff
542 integer,
parameter :: incdim = 4, iazdim = 4
546 subroutine initsolv_lsatdis(me, master, nwaves, nazdir, nstoch, effac, do_physb, kxw)
550 integer :: me, master
551 integer :: nwaves, nazdir
553 real(kind=kind_phys) :: effac
555 real(kind=kind_phys) :: kxw
560 integer :: inc, jk, jl, iazi, i, j, k
562 if( nwaves == 0 .or. nstoch == 1 )
then
576 end subroutine initsolv_lsatdis
583 use machine ,
only : kind_phys
584 use ugwp_common,
only : arad, pi, pi2, hpscale, rhp, rhp2, rh4, omega2
586 use ugwp_common,
only : mkzmin, mkz2min, mkzmax, mkz2max, ucrit => cdmin
590 real(kind=kind_phys),
parameter :: maxdudt = 250.e-5, maxdtdt=15.e-2
591 real(kind=kind_phys),
parameter :: dked_min =0.01, dked_max=250.0
593 real(kind=kind_phys),
parameter :: gptwo=2.0
595 real(kind=kind_phys) ,
parameter :: bnfix = 6.28/300., bnfix2= bnfix * bnfix
596 real(kind=kind_phys) ,
parameter :: bnfix4 = bnfix2 * bnfix2
597 real(kind=kind_phys) ,
parameter :: bnfix3 = bnfix2 * bnfix
601 integer ,
parameter :: iazidim=4
602 integer ,
parameter :: incdim=25
604 real(kind=kind_phys) ,
parameter :: zcimin = 2.5
605 real(kind=kind_phys) ,
parameter :: zcimax = 125.0
606 real(kind=kind_phys) ,
parameter :: zgam = 0.25
610 real(kind=kind_phys) ,
parameter :: pind_wd = 5./3.
611 real(kind=kind_phys) ,
parameter :: sind_kz = 1.
612 real(kind=kind_phys) ,
parameter :: tind_kz = 3.
613 real(kind=kind_phys) ,
parameter :: stind_kz = sind_kz + tind_kz
617 real(kind=kind_phys) :: nslope
618 real(kind=kind_phys) :: lzstar
619 real(kind=kind_phys) :: lzmin
620 real(kind=kind_phys) :: lzmax
621 real(kind=kind_phys) :: lhmet
622 real(kind=kind_phys) :: tamp_mpa
623 real(kind=kind_phys) :: tau_min
625 real(kind=kind_phys) :: gw_eff
627 real(kind=kind_phys) :: v_kxw, rv_kxw, v_kxw2
632 integer :: nwav, nazd, nst
633 real(kind=kind_phys) :: eff
635 real(kind=kind_phys) :: zaz_fct, zms
636 real(kind=kind_phys),
allocatable :: zci(:), zci4(:), zci3(:),zci2(:), zdci(:)
637 real(kind=kind_phys),
allocatable :: zcosang(:), zsinang(:)
638 real(kind=kind_phys),
allocatable :: lzmet(:), czmet(:), mkzmet(:), dczmet(:), dmkz(:)
644 real(kind=kind_phys),
parameter :: ipr_pt = 0.5
645 real(kind=kind_phys),
parameter :: lturb = 30., sc2 = lturb*lturb
646 real(kind=kind_phys),
parameter :: ulturb=150., sc2u = ulturb* ulturb
647 real(kind=kind_phys),
parameter :: ric =0.25
648 real(kind=kind_phys),
parameter :: rimin = -10., prmin = 0.25
649 real(kind=kind_phys),
parameter :: prmax = 4.0
653 subroutine initsolv_wmsdis(me, master, nwaves, nazdir, nstoch, effac, do_physb, kxw, version)
664 integer,
intent(in) :: me, master, nwaves, nazdir, nstoch
665 integer,
intent(in) :: version
667 real(kind=kind_phys),
intent(in) :: effac, kxw
668 logical,
intent(in) :: do_physb
673 real(kind=kind_phys) :: dlzmet
674 real(kind=kind_phys) :: cstar,rcstar, nslope3, fnorm, zcin
676 integer :: inc, jk, jl, iazi
678 real(kind=kind_phys) :: zang, zang1, znorm
679 real(kind=kind_phys) :: zx1, zx2, ztx, zdx, zxran, zxmin, zxmax, zx, zpexp
680 real(kind=kind_phys) :: fpc, fpc_dc
681 real(kind=kind_phys) :: ae1,ae2
682 if( nwaves == 0)
then
702 v_kxw = kxw ; v_kxw2 = v_kxw*v_kxw
705 allocate ( zci(nwav), zci4(nwav), zci3(nwav),zci2(nwav), zdci(nwav) )
706 allocate ( zcosang(nazd), zsinang(nazd) )
707 allocate (lzmet(nwav), czmet(nwav), mkzmet(nwav), dczmet(nwav), dmkz(nwav) )
723 zang = pi2 / float(nazd)
730 zang1 = (iazi-1)*zang
731 zcosang(iazi) = cos(zang1)
732 zsinang(iazi) = sin(zang1)
733 znorm = znorm + abs(zcosang(iazi))
736 zaz_fct = 2.0 / znorm
746 zxran = zxmax - zxmin
747 zdx = zxran / real(nwav-1)
750 zx1 = zxran/(exp(ae1)-1.0 )
760 ztx = real(inc-1)*zdx+zxmin
761 ae1 = (ztx-zxmin)/zgam
762 zx = zx1*exp(ae1)+zx2
764 zdci(inc) = zci(inc)**2*(zx1/zgam)*exp(ae1)*zdx
765 zci4(inc) = (zms*zci(inc))**4
766 zci2(inc) = (zms*zci(inc))**2
767 zci3(inc) = (zms*zci(inc))**3
774 if (version == 1)
then
776 dlzmet = (lzmax-lzmin)/ real(nwav-1)
778 lzmet(inc) = lzmin + (inc-1)*dlzmet
779 mkzmet(inc) = pi2/lzmet(inc)
780 zci(inc) =lzmet(inc)/(pi2/bnfix)
781 zci4(inc) = (zms*zci(inc))**4
782 zci2(inc) = (zms*zci(inc))**2
783 zci3(inc) = (zms*zci(inc))**3
787 zdx = (zci(nwav)-zci(1))/ real(nwav-1)
798 if (me == master)
then
800 print *,
'ugwp_v0: zcimin=' , zcimin
801 print *,
'ugwp_v0: zcimax=' , zcimax
802 print *,
'ugwp_v0: zgam= ', zgam
807 print *,
'ugwp_v1: zcimin/zci=' , maxval(zci)
808 print *,
'ugwp_v1: zcimax/zci=' , minval(zci)
809 print *,
'ugwp_v1: cd_crit=', ucrit
810 print *,
'ugwp_v1: launch_level', ilaunch
811 print *,
' ugwp_v1 lzstar=', lzstar
812 print *,
' ugwp_v1 nslope=', nslope
817 zcin =zci(inc)*rcstar
818 fpc = rcstar*(zcin*zcin)/(1.+ zcin**nslope3)
819 fpc_dc = fpc * zdci(inc)
820 write(6,111) inc, zci(inc), zdci(inc),ucrit, fpc, fpc_dc, 6.28e-3/bnfix*zci(inc)
826 111
format(
'wms-zci', i4, 7 (3x, f8.3))
828 end subroutine initsolv_wmsdis