CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
cires_ugwpv1_triggers.F90
1
4
5 use machine, only: kind_phys
6
7contains
8
9!
10!
12 subroutine slat_geos5_tamp_v1(im, tau_amp, xlatdeg, tau_gw)
13
14 implicit none
15 integer :: im
16 real(kind=kind_phys) :: tau_amp, xlatdeg(im), tau_gw(im)
17 real(kind=kind_phys) :: latdeg, flat_gw, tem
18 integer :: i
19
20!
21! if-lat
22!
23 do i=1, im
24 latdeg = abs(xlatdeg(i))
25 if (latdeg < 15.3) then
26 tem = (latdeg-3.0) / 8.0
27 flat_gw = 0.75 * exp(-tem * tem)
28 if (flat_gw < 1.2 .and. latdeg <= 3.0) flat_gw = 0.75
29 elseif (latdeg < 31.0 .and. latdeg >= 15.3) then
30 flat_gw = 0.10
31 elseif (latdeg < 60.0 .and. latdeg >= 31.0) then
32 tem = (latdeg-60.0) / 23.0
33 flat_gw = 0.50 * exp(- tem * tem)
34 elseif (latdeg >= 60.0) then
35 tem = (latdeg-60.0) / 70.0
36 flat_gw = 0.50 * exp(- tem * tem)
37 endif
38 tau_gw(i) = tau_amp*flat_gw
39 enddo
40!
41 end subroutine slat_geos5_tamp_v1
42!
45 subroutine slat_geos5_2020(im, tau_amp, xlatdeg, tau_gw)
46
47 implicit none
48 integer :: im
49 real(kind=kind_phys) :: tau_amp, xlatdeg(im), tau_gw(im)
50 real(kind=kind_phys) :: latdeg, flat_gw, tem
51 real(kind=kind_phys), parameter :: fampqbo = 1.25 ! 1.5
52 real(kind=kind_phys), parameter :: famp60s = 1.0 ! 1.5
53 real(kind=kind_phys), parameter :: famp60n = 1.0 ! 1.0
54 real(kind=kind_phys), parameter :: famp30 = 0.25 ! 0.4
55
56 real(kind=kind_phys), parameter :: swid15 = 12.5
57 real(kind=kind_phys), parameter :: swid60s = 30.0 ! 40
58 real(kind=kind_phys), parameter :: swid60n = 25.0 ! 30
59 integer :: i
60!
61!
62!
63 do i=1, im
64
65 latdeg = abs(xlatdeg(i))
66 if (latdeg < 15.3) then
67 tem = (latdeg-3.0) / swid15
68 flat_gw = fampqbo * exp(-tem * tem)
69 if (latdeg <= 3.0) flat_gw = fampqbo
70 elseif (latdeg < 31.0 .and. latdeg >= 15.3) then
71 flat_gw = famp30
72 elseif (latdeg < 60.0 .and. latdeg >= 31.0) then
73 tem = (latdeg-60.0) / 23.0
74 flat_gw = famp60n* exp(- tem * tem)
75 elseif (latdeg >= 60.0) then
76 tem = (latdeg-60.0) /swid60n
77 flat_gw = famp60n * exp(- tem * tem)
78 endif
79
80 if (xlatdeg(i) <= -31.0) then
81!
82 if (latdeg < 60.0 .and. latdeg >= 31.0) then
83 tem = (latdeg-60.0) / 23.0
84 flat_gw = famp60s * exp(- tem * tem)
85 endif
86 if (latdeg >= 60.0) then
87 tem = (latdeg-60.0) /swid60s
88 flat_gw = famp60s * exp(- tem * tem)
89 endif
90
91 endif
92 tau_gw(i) = tau_amp*flat_gw
93 enddo
94!
95 end subroutine slat_geos5_2020
96
98 subroutine slat_geos5(im, xlatdeg, tau_gw)
99
100 implicit none
101 integer :: im
102 real(kind=kind_phys) :: xlatdeg(im)
103 real(kind=kind_phys) :: tau_gw(im)
104 real(kind=kind_phys) :: latdeg
105 real(kind=kind_phys), parameter :: tau_amp = 3.5e-3 ! 3.5 mPa
106 real(kind=kind_phys) :: trop_gw, flat_gw
107 integer :: i
108!
109! if-lat
110!
111 trop_gw = 0.75
112 do i=1, im
113 latdeg = xlatdeg(i)
114 if (-15.3 < latdeg .and. latdeg < 15.3) then
115 flat_gw = trop_gw*exp(-( (abs(latdeg)-3.)/8.0)**2)
116 if (flat_gw < 1.2 .and. abs(latdeg) <= 3.) flat_gw = trop_gw
117 else if (latdeg > -31. .and. latdeg <= -15.3) then
118 flat_gw = 0.10
119 else if (latdeg < 31. .and. latdeg >= 15.3) then
120 flat_gw = 0.10
121 else if (latdeg > -60. .and. latdeg <= -31.) then
122 flat_gw = 0.50*exp(-((abs(latdeg)-60.)/23.)**2)
123 else if (latdeg < 60. .and. latdeg >= 31.) then
124 flat_gw = 0.50*exp(-((abs(latdeg)-60.)/23.)**2)
125 else if (latdeg <= -60.) then
126 flat_gw = 0.50*exp(-((abs(latdeg)-60.)/70.)**2)
127 else if (latdeg >= 60.) then
128 flat_gw = 0.50*exp(-((abs(latdeg)-60.)/70.)**2)
129 end if
130 tau_gw(i) = tau_amp*flat_gw
131 enddo
132!
133 end subroutine slat_geos5
134
138 subroutine get_spectra_tau_convgw &
139 (nw, im, levs, dcheat, scheat, precip, icld, xlatd, sinlat, coslat,taub, klev, if_src, nf_src)
140!
141! temporarily can put GEOS-5/MERRA-2 GW-lat dependent function
142!
143 integer :: nw, im, levs
144 integer,dimension(im,3) :: icld
145 real(kind=kind_phys), dimension(im, levs) :: dcheat, scheat
146 real(kind=kind_phys), dimension(im) :: precip, xlatd, sinlat, coslat
147 real(kind=kind_phys), dimension(im) :: taub
148 integer, dimension(im) :: klev, if_src
149 integer :: nf_src
150!
151! locals
152 real(kind=kind_phys), parameter :: precip_max = 100. ! mm/day
153 real(kind=kind_phys), parameter :: tau_amp = 3.5e-3 ! 3.5 mPa
154
155 integer :: i, k, klow, ktop, kmid
156 real(kind=kind_phys) :: dtot, dmax, daver
157!
158 nf_src = 0
159 if_src(1:im) = 0
160 taub(1:im) = 0.0
161 do i=1, im
162 klow = icld(i,1)
163 ktop = icld(i,2)
164 kmid= icld(i,3)
165 if (klow == -99 .and. ktop == -99) then
166 cycle
167 else
168 klev(i) = ktop
169 k = klow
170 klev(i) = k
171 dmax = abs(dcheat(i,k) + scheat(i,k))
172 do k=klow+1, ktop
173 dtot =abs(dcheat(i,k) + scheat(i,k))
174 if ( dtot > dmax) then
175 klev(i) = k
176 dmax = dtot
177 endif
178 enddo
179!
180! klev as max( dcheat(i,k) + scheat)
181! vertical width of conv-heating
182!
183! counts/triiger=1 & taub(i)
184!
185 nf_src = nf_src +1
186 if_src(i) = 1
187 taub(i) = tau_amp* precip(i)/precip_max*coslat(i)
188 endif
189
190 enddo
191!
192!
193!
194 call slat_geos5(im, xlatd, taub)
195 nf_src =im
196 do i=1, im
197 if_src(i) = 1
198 klev(i) = 127-45
199 enddo
200
201! with info on precip/clouds/dc_heat create Bulk
202! taub(im), klev(im)
203!
204! print *, ' get_spectra_tau_convgw '
205 end subroutine get_spectra_tau_convgw
206
208 subroutine get_spectra_tau_nstgw(nw, im, levs, trig_fgf, xlatd, sinlat, coslat, taub, klev, if_src, nf_src)
209 integer :: nw, im, levs
210 real(kind=kind_phys), dimension(im, levs) :: trig_fgf
211! real(kind=kind_phys), dimension(im, levs+1) :: pint
212 real(kind=kind_phys), dimension(im) :: xlatd, sinlat, coslat
213 real(kind=kind_phys), dimension(im) :: taub
214 integer, dimension(im) :: klev, if_src
215 integer :: nf_src
216! locals
217 real(kind=kind_phys), parameter :: tlim_fgf = 100. ! trig_fgf > tlim_fgf, launch waves should scale-dependent
218 real(kind=kind_phys), parameter :: tau_amp = 3.5e-3 ! 3.5 mPa
219 real(kind=kind_phys), parameter :: pmax = 750.e2, pmin = 100.e2
220 integer, parameter :: klow =127-92, ktop=127-45
221 integer, parameter :: kwidth = ktop-klow+1
222 integer :: i, k, kex
223 real(kind=kind_phys) :: dtot, dmax, daver
224 real(kind=kind_phys) :: fnorm, tau_min
225 nf_src = 0
226 if_src(1:im) = 0
227 taub(1:im) = 0.0
228 fnorm = 1.0 / float(kwidth)
229 tau_min = tau_amp*fnorm
230 do i=1, im
231!
232! only trop-c fjets so find max(trig_fgf) => klev
233! use abs-values to scale tau_amp
234!
235
236 k = klow
237 klev(i) = k
238 dmax = abs(trig_fgf(i,k))
239 kex = 0
240 if (dmax >= tlim_fgf) kex = kex+1
241 do k=klow+1, ktop
242 dtot = abs(trig_fgf(i,k))
243 if (dtot >= tlim_fgf) kex = kex+1
244 if ( dtot > dmax) then
245 klev(i) = k
246 dmax = dtot
247 endif
248 enddo
249
250 if (dmax .ge. tlim_fgf) then
251 nf_src = nf_src +1
252 if_src(i) = 1
253 taub(i) = tau_min*float(kex) !* precip(i)/precip_max*coslat(i)
254 endif
255
256 enddo
257!
258! print *, ' get_spectra_tau_nstgw '
259 call slat_geos5(im, xlatd, taub)
260 nf_src =im
261 do i=1, im
262 if_src(i) = 1
263 klev(i) = 127-45 ! FV3-127L
264 enddo
265!
266 end subroutine get_spectra_tau_nstgw
267!
269 subroutine get_spectra_tau_okw(nw, im, levs, trig_okw, xlatd, sinlat, coslat, taub, klev, if_src, nf_src)
270 integer :: nw, im, levs
271 real(kind=kind_phys), dimension(im, levs) :: trig_okw
272! real(kind=kind_phys), dimension(im, levs+1) :: pint
273 real(kind=kind_phys), dimension(im) :: xlatd, sinlat, coslat
274 real(kind=kind_phys), dimension(im) :: taub
275 integer, dimension(im) :: klev, if_src
276 integer :: nf_src
277! locals
278 real(kind=kind_phys), parameter :: tlim_okw = 100. ! trig_fgf > tlim_fgf, launch GWs should scale-dependent
279 real(kind=kind_phys), parameter :: tau_amp = 35.e-3 ! 35 mPa
280 real(kind=kind_phys), parameter :: pmax = 750.e2, pmin = 100.e2
281 integer, parameter :: klow =127-92, ktop=127-45 ! for FV3-127L
282 integer, parameter :: kwidth = ktop-klow+1
283 integer :: i, k, kex
284 real(kind=kind_phys) :: dtot, dmax, daver
285 real(kind=kind_phys) :: fnorm, tau_min
286
287 nf_src = 0
288 if_src(1:im) = 0
289 taub(1:im) = 0.0
290 fnorm = 1./float(kwidth)
291 tau_min = tau_amp*fnorm
292 print *, ' get_spectra_tau_okwgw '
293 do i=1, im
294 k = klow
295 klev(i) = k
296 dmax = abs(trig_okw(i,k))
297 kex = 0
298 if (dmax >= tlim_okw) kex = kex+1
299 do k=klow+1, ktop
300 dtot = abs(trig_okw(i,k))
301 if (dtot >= tlim_fgf ) kex = kex+1
302 if ( dtot > dmax) then
303 klev(i) = k
304 dmax = dtot
305 endif
306 enddo
307!
308 if (dmax >= tlim_okw) then
309 nf_src = nf_src + 1
310 if_src(i) = 1
311 taub(i) = tau_min*float(kex) !* precip(i)/precip_max*coslat(i)
312 endif
313
314 enddo
315 print *, ' get_spectra_tau_okwgw '
316 end subroutine get_spectra_tau_okw
317
318end module cires_ugwpv1_triggers