CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
cires_ugwp_initialize.F90
1
3! initialization of ugwp_common_v0
4! init gw-solvers (1,2) .. no UFS-funds for (3,4) tests
5! init gw-source specifications
6! init gw-background dissipation
7!===============================
8
11!
12 use machine, only: kind_phys
13 use physcons, only : pi => con_pi, grav => con_g, rd => con_rd, &
14 rv => con_rv, cpd => con_cp, fv => con_fvirt,&
15 arad => con_rerth
16 implicit none
17
18 real(kind=kind_phys), parameter :: grcp = grav/cpd, rgrav = 1.0d0/grav, &
19 rdi = 1.0d0/rd, &
20 gor = grav/rd, gr2 = grav*gor, gocp = grav/cpd, &
21 rcpd = 1./cpd, rcpd2 = 0.5*rcpd, &
22 pi2 = pi + pi, omega1 = pi2/86400.0, &
23 omega2 = omega1+omega1, &
24 rad_to_deg=180.0/pi, deg_to_rad=pi/180.0, &
25 dw2min=1.0, bnv2min=1.e-6, velmin=sqrt(dw2min)
26
27
28 end module ugwp_common_v0
29!
30!
31!===================================================
32!
33!Part-1 init => wave dissipation + RFriction
34!
35!===================================================
37 subroutine init_global_gwdis_v0(levs, zkm, pmb, kvg, ktg, krad, kion)
38 implicit none
39
40 integer :: levs
41 real, intent(in) :: zkm(levs), pmb(levs)
42 real, intent(out), dimension(levs+1) :: kvg, ktg, krad, kion
43!
44!locals + data
45!
46 integer :: k
47 real, parameter :: vusurf = 2.e-5
48 real, parameter :: musurf = vusurf/1.95
49 real, parameter :: hpmol = 8.5
50!
51 real, parameter :: kzmin = 0.1
52 real, parameter :: kturbo = 100.
53 real, parameter :: zturbo = 130.
54 real, parameter :: zturw = 30.
55 real, parameter :: inv_pra = 3. !kt/kv =inv_pr
56!
57 real, parameter :: alpha = 1./86400./15.
58!
59 real, parameter :: kdrag = 1./86400./10.
60 real, parameter :: zdrag = 100.
61 real, parameter :: zgrow = 50.
62!
63 real :: vumol, mumol, keddy, ion_drag
64!
65 do k=1, levs
66 vumol = vusurf*exp(-zkm(k)/hpmol)
67 mumol = musurf*exp(-zkm(k)/hpmol)
68
69 keddy = kturbo*exp(-((zkm(k)-zturbo) /zturw)**2)
70
71 kvg(k) = vumol + keddy
72 ktg(k) = mumol + keddy*inv_pra
73
74 krad(k) = alpha
75!
76 ion_drag = kdrag
77!
78 kion(k) = ion_drag
79 enddo
80
81 k= levs+1
82 kion(k) = kion(k-1)
83 krad(k) = krad(k-1)
84 kvg(k) = kvg(k-1)
85 ktg(k) = ktg(k-1)
86!
87 end subroutine init_global_gwdis_v0
88
89
90! ========================================================================
91! Part 2 - sources
92! wave sources
93! ========================================================================
94!
95! ugwpv0_oro_init
96!
97!=========================================================================
99
100 use ugwp_common_v0, only : bnv2min, grav, grcp, fv, grav, cpd, grcp, pi
101
102 implicit none
103!
104! constants and "crirtical" values to run oro-mtb_gw physics
105!
106! choice of oro-scheme: strver = 'vay_2018' , 'gfs_2018', 'kdn_2005', 'smc_2000'
107!
108 character(len=8) :: strver = 'gfs_2018'
109 character(len=8) :: strbase = 'gfs_2018'
110 real, parameter :: rimin=-10., ric=0.25
111
112!
113 real, parameter :: efmin=0.5, efmax=10.0
114 real, parameter :: hpmax=2400.0, hpmin=25.0
115 real, parameter :: sigma_std=1./100., gamm_std=1.0
116
117 real, parameter :: frmax=10., frc =1.0, frmin =0.01
118!
119
120 real, parameter :: ce=0.8, ceofrc=ce/frc, cg=0.5
121 real, parameter :: gmax=1.0, veleps=1.0, factop=0.5
122!
123 real, parameter :: rlolev=50000.0
124!
125 real, parameter :: hncrit=9000. ! max value in meters for elvmax
126
127! hncrit set to 8000m and sigfac added to enhance elvmax mtn hgt
128
129 real, parameter :: sigfac=4.0 ! mb3a expt test for elvmax factor
130 real, parameter :: hminmt=50. ! min mtn height (*j*)
131 real, parameter :: minwnd=1.0 ! min wind component (*j*)
132 real, parameter :: dpmin=5000.0 ! minimum thickness of the reference layer in pa
133
134 real, parameter :: kxoro=6.28e-3/200. !
135 real, parameter :: coro = 0.0
136 integer, parameter :: nridge=2
137
138 real :: cdmb ! scale factors for mtb
139 real :: cleff ! scale factors for orogw
140 integer :: nworo ! number of waves
141 integer :: nazoro ! number of azimuths
142 integer :: nstoro ! flag for stochastic launch above SG-peak
143
144 integer, parameter :: mdir = 8
145 real, parameter :: fdir=.5*mdir/pi
146
147 integer nwdir(mdir)
148 data nwdir/6,7,5,8,2,3,1,4/
149 save nwdir
150
151 real, parameter :: odmin = 0.1, odmax = 10.0
152!------------------------------------------------------------------------------
153! small-scale orography parameters for TOFD of Beljaars et al., 2004, QJRMS
154!------------------------------------------------------------------------------
155
156 integer, parameter :: n_tofd = 2 ! depth of SSO for TOFD compared with Zpbl
157 real, parameter :: const_tofd = 0.0759 ! alpha*beta*Cmd*Ccorr*2.109 = 12.*1.*0.005*0.6*2.109 = 0.0759
158 real, parameter :: ze_tofd = 1500.0 ! BJ's z-decay in meters
159 real, parameter :: a12_tofd = 0.0002662*0.005363 ! BJ's k-spect const for sigf2 * a1*a2*exp(-[z/zdec]**1.5]
160 real, parameter :: ztop_tofd = 10.*ze_tofd ! no TOFD > this height too higher 15 km
161!------------------------------------------------------------------------------
162!
163 real, parameter :: fcrit_sm = 0.7, fcrit_sm2 = fcrit_sm * fcrit_sm
164 real, parameter :: fcrit_gfs = 0.7
165 real, parameter :: fcrit_mtb = 0.7
166
167 real, parameter :: lzmax = 18.e3 ! 18 km
168 real, parameter :: mkzmin = 6.28/lzmax
169 real, parameter :: mkz2min = mkzmin*mkzmin
170 real, parameter :: zbr_pi = (3.0/2.0)*pi
171 real, parameter :: zbr_ifs = 0.5*pi
172
173 contains
174!
175 subroutine init_oro_gws_v0(nwaves, nazdir, nstoch, effac, &
176 lonr, kxw, cdmbgwd )
177!
178!
179 integer :: nwaves, nazdir, nstoch
180 integer :: lonr
181 real :: cdmbgwd(2) ! scaling factors for MTb (1) & (2) for cleff = cleff * cdmbgwd(2)
182 ! high res-n "larger" MTB and "less-active" cleff in GFS-2018
183 real :: cdmbX
184 real :: kxw
185 real :: effac ! it is analog of cdmbgwd(2) for GWs, off for now
186!-----------------------------! GFS-setup for cdmb & cleff
187! cdmb = 4.0 * (192.0/IMX)
188! cleff = 0.5E-5 / SQRT(IMX/192.0) = 0.5E-5*SQRT(192./IMX)
189!
190 real, parameter :: lonr_refmb = 4.0 * 192.0
191 real, parameter :: lonr_refgw = 192.0
192
193! copy to "ugwp_oro_init" => nwaves, nazdir, nstoch
194
195 nworo = nwaves
196 nazoro = nazdir
197 nstoro = nstoch
198
199 cdmbx = lonr_refmb/float(lonr)
200 cdmb = cdmbx
201 if (cdmbgwd(1) >= 0.0) cdmb = cdmb * cdmbgwd(1)
202
203 cleff = 0.5e-5 * sqrt(lonr_refgw/float(lonr)) !* effac
204
205!!! cleff = kxw * sqrt(lonr_refgw/float(lonr)) !* effac
206
207 if (cdmbgwd(2) >= 0.0) cleff = cleff * cdmbgwd(2)
208!
209!....................................................................
210! higher res => smaller h' ..&.. higher kx
211! flux_gwd ~ 'u'^2*kx/kz ~kxu/n ~1/dx *u/n tau ~ h'*h'*kx*kx = const (h'-less kx-grow)
212!....................................................................
213!
214! print *, ' init_oro_gws 2-1cdmb', cdmbgwd(2), cdmbgwd(1)
215 end subroutine init_oro_gws_v0
216!
217
218 end module ugwpv0_oro_init
219!=============================== end of GW sources
220!
221! init specific gw-solvers (1,2,3,4)
222!
223
224!===============================
225! Part -3 init wave solvers
226!===============================
227
229 implicit none
230
231 integer :: nwav, nazd
232 integer :: nst
233 real :: eff
234 integer, parameter :: incdim = 4, iazdim = 4
235!
236 contains
237
238 subroutine initsolv_lsatdis_v0(me, master, nwaves, nazdir, nstoch, effac, do_physb, kxw)
239
240 implicit none
241!
242 integer :: me, master
243 integer :: nwaves, nazdir
244 integer :: nstoch
245 real :: effac
246 logical :: do_physb
247 real :: kxw
248!
249!locals: define azimuths and Ch(nwaves) - domain when physics-based soureces
250! are not actibve
251!
252 integer :: inc, jk, jl, iazi, i, j, k
253
254 if( nwaves == 0 .or. nstoch == 1 ) then
255! redefine from the default
256 nwav = incdim
257 nazd = iazdim
258 nst = 0
259 eff = 1.0
260 else
261! from input_nml multi-wave spectra
262 nwav = nwaves
263 nazd = nazdir
264 nst = nstoch
265 eff = effac
266 endif
267!
268 end subroutine initsolv_lsatdis_v0
269!
270 end module ugwpv0_lsatdis_init
271!
272!
274
275 use ugwp_common_v0, only : pi, pi2
276 implicit none
277
278 real, parameter :: maxdudt = 250.e-5
279
280 real, parameter :: hpscale= 7000., rhp2 = 0.5/hpscale
281 real, parameter :: omega2 = 2.*6.28/86400
282 real, parameter :: gptwo=2.0
283
284 real, parameter :: dked_min =0.01
285 real, parameter :: gssec = (6.28/30.)**2 ! max-value for bn2
286 real, parameter :: bv2min = (6.28/60./120.)**2 ! min-value for bn2 7.6(-7) 2 hrs
287 real, parameter :: minvel = 0.5
288
289!
290! make parameter list that will be passed to SOLVER
291!
292
293 real, parameter :: v_kxw = 6.28e-3/200.
294 real, parameter :: v_kxw2 = v_kxw*v_kxw
295 real, parameter :: tamp_mpa = 30.e-3
296 real, parameter :: zfluxglob= 3.75e-3
297
298 real , parameter :: nslope=1 ! the GW sprctral slope at small-m
299
300 integer , parameter :: iazidim=4 ! number of azimuths
301 integer , parameter :: incdim=25 ! number of discrete cx - spectral elements in launch spectrum
302 real , parameter :: ucrit2=0.5
303
304 real , parameter :: zcimin = ucrit2
305 real , parameter :: zcimax = 125.0
306 real , parameter :: zgam = 0.25
307 real , parameter :: zms_l = 2000.0, zms = pi2 / zms_l, zmsi = 1.0 / zms
308
309 integer :: ilaunch
310 real :: gw_eff
311
312!===========================================================================
313 integer :: nwav, nazd, nst
314 real :: eff
315
316 real :: zaz_fct
317 real, allocatable :: zci(:), zci4(:), zci3(:),zci2(:), zdci(:)
318 real, allocatable :: zcosang(:), zsinang(:)
319 contains
320!============================================================================
321 subroutine initsolv_wmsdis_v0(me, master, nwaves, nazdir, nstoch, effac, do_physb, kxw)
322
323 implicit none
324!
325!input -control for solvers:
326! nwaves, nazdir, nstoch, effac, do_physb, kxw
327!
328!
329 integer :: me, master, nwaves, nazdir, nstoch
330 real :: effac, kxw
331 logical :: do_physb
332!
333!locals
334!
335 integer :: inc, jk, jl, iazi
336!
337 real :: zang, zang1, znorm
338 real :: zx1, zx2, ztx, zdx, zxran, zxmin, zxmax, zx, zpexp
339
340 if( nwaves == 0) then
341!
342! redefine from the deafault
343!
344 nwav = incdim
345 nazd = iazidim
346 nst = 0
347 eff = 1.0
348 gw_eff = eff
349 else
350!
351! from input.nml
352!
353 nwav = nwaves
354 nazd = nazdir
355 nst = nstoch
356 gw_eff = effac
357 endif
358
359 allocate ( zci(nwav), zci4(nwav), zci3(nwav),zci2(nwav), zdci(nwav) )
360 allocate ( zcosang(nazd), zsinang(nazd) )
361
362 if (me == master) then
363 print *, 'ugwp_v0: init_gw_wmsdis_control '
364! print *, 'ugwp_v0: WMSDIS launch layer ', klaunch
365 print *, 'ugwp_v0: WMSDIS launch layer ', ilaunch
366 print *, 'ugwp_v0: WMSDID tot_mflux in mpa', tamp_mpa*1000.
367 endif
368
369 zpexp = gptwo * 0.5 ! gptwo=2 , zpexp = 1.
370
371!
372! set up azimuth directions and some trig factors
373!
374!
375 zang = pi2 / float(nazd)
376
377! get normalization factor to ensure that the same amount of momentum
378! flux is directed (n,s,e,w) no mater how many azimuths are selected.
379!
380 znorm = 0.0
381 do iazi=1, nazd
382 zang1 = (iazi-1)*zang
383 zcosang(iazi) = cos(zang1)
384 zsinang(iazi) = sin(zang1)
385 znorm = znorm + abs(zcosang(iazi))
386 enddo
387! zaz_fct = 1.0
388 zaz_fct = 2.0 / znorm ! correction factor for azimuthal sums
389
390! define coordinate transform for "Ch" ....x = 1/c stretching transform
391! -----------------------------------------------
392! note that this is expresed in terms of the intrinsic phase speed
393! at launch ci=c-u_o so that the transformation is identical
394! see eq. 28-30 of scinocca 2003. x = 1/c stretching transform
395!
396 zxmax = 1.0 / zcimin
397 zxmin = 1.0 / zcimax
398 zxran = zxmax - zxmin
399 zdx = zxran / real(nwav-1) ! dkz
400!
401 zx1 = zxran/(exp(zxran/zgam)-1.0 ) ! zgam =1./4.
402 zx2 = zxmin - zx1
403
404!
405! computations for zci =1/zx
406! if(lgacalc) zgam=(zxmax-zxmin)/log(zxmax/zxmin)
407! zx1=zxran/(exp(zxran/zgam)-1.0_jprb)
408! zx2=zxmin-zx1
409! zms = pi2 / zms_l
410 do inc=1, nwav
411 ztx = real(inc-1)*zdx+zxmin
412 zx = zx1*exp((ztx-zxmin)/zgam)+zx2 !eq. 29 of scinocca 2003
413 zci(inc) = 1.0 /zx !eq. 28 of scinocca 2003
414 zdci(inc) = zci(inc)**2*(zx1/zgam)*exp((ztx-zxmin)/zgam)*zdx !eq. 30 of scinocca 2003
415 zci4(inc) = (zms*zci(inc))**4
416 zci2(inc) = (zms*zci(inc))**2
417 zci3(inc) = (zms*zci(inc))**3
418 enddo
419!
420!
421! all done and print-out
422!
423!
424 if (me == master) then
425 print *
426 print *, 'ugwp_v0: zcimin=' , zcimin
427 print *, 'ugwp_v0: zcimax=' , zcimax
428 print *, 'ugwp_v0: cd_crit=', zgam ! m/s precision for crit-level
429 print *, 'ugwp_v0: launch_level', ilaunch
430 print *, ' ugwp_v0 zms_l=', zms_l
431 print *, ' ugwp_vgw nslope=', nslope
432
433 print *
434 endif
435
436 end subroutine initsolv_wmsdis_v0
437!
438 end module ugwpv0_wmsdis_init
Define constants.