CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_stochastics.F90
1
3
5
6 contains
7
16 subroutine gfs_stochastics_init (si,vfact_ca,km,do_ca,ca_global, errmsg, errflg)
17
18 use machine, only: kind_phys
19
20 implicit none
21 real(kind_phys), dimension(:), intent(in) :: si
22 real(kind_phys), dimension(:), intent(inout) :: vfact_ca
23 integer, intent(in) :: km
24 logical, intent(in) :: do_ca
25 logical, intent(in) :: ca_global
26 character(len=*), intent(out) :: errmsg
27 integer, intent(out) :: errflg
28 integer :: k,nz
29
30 errmsg = ''
31 errflg = 0
32 if (do_ca .and. ca_global) then
33 nz=min(km,size(vfact_ca))
34 vfact_ca(:)=0.0
35 do k=1,nz
36 if (si(k) .lt. 0.1 .and. si(k) .gt. 0.025) then
37 vfact_ca(k) = (si(k)-0.025)/(0.1-0.025)
38 else if (si(k) .lt. 0.025) then
39 vfact_ca(k) = 0.0
40 else
41 vfact_ca(k) = 1.0
42 endif
43 enddo
44 vfact_ca(2)=vfact_ca(3)*0.5
45 vfact_ca(1)=0.0
46 endif
47 end subroutine gfs_stochastics_init
48
60 subroutine gfs_stochastics_run (im, km, kdt, delt, do_sppt, pert_mp, use_zmtnblck, &
61 do_shum ,do_skeb, do_ca,ca_global,ca1,vfact_ca, &
62 zmtnblck, sppt_wts, skebu_wts, skebv_wts, shum_wts,&
63 diss_est, ugrs, vgrs, tgrs, qgrs_wv, &
64 qgrs_cw, qgrs_rw, qgrs_sw, qgrs_iw, qgrs_gl, &
65 gu0, gv0, gt0, gq0_wv, dtdtnp, &
66 gq0_cw, gq0_rw, gq0_sw, gq0_iw, gq0_gl, &
67 rain, rainc, tprcp, totprcp, cnvprcp, &
68 totprcpb, cnvprcpb, cplflx, cpllnd, &
69 rain_cpl, snow_cpl, drain_cpl, dsnow_cpl, &
70 ntcw,ntrw,ntsw,ntiw,ntgl, &
71 errmsg, errflg)
72
73 use machine, only: kind_phys
74
75 implicit none
76
77 integer, intent(in) :: im
78 integer, intent(in) :: km
79 integer, intent(in) :: kdt
80 real(kind_phys), intent(in) :: delt
81 logical, intent(in) :: do_sppt
82 logical, intent(in) :: pert_mp
83 logical, intent(in) :: do_ca
84 logical, intent(in) :: ca_global
85 logical, intent(in) :: use_zmtnblck
86 logical, intent(in) :: do_shum
87 logical, intent(in) :: do_skeb
88 real(kind_phys), dimension(:), intent(in) :: zmtnblck
89 ! sppt_wts only allocated if do_sppt == .true.
90 real(kind_phys), dimension(:,:), intent(inout), optional :: sppt_wts
91 ! skebu_wts, skebv_wts only allocated if do_skeb == .true.
92 real(kind_phys), dimension(:,:), intent(in), optional :: skebu_wts
93 real(kind_phys), dimension(:,:), intent(in), optional :: skebv_wts
94 ! shum_wts only allocated if do_shum == .true.
95 real(kind_phys), dimension(:,:), intent(in), optional :: shum_wts
96 real(kind_phys), dimension(:,:), intent(in) :: diss_est
97 real(kind_phys), dimension(:,:), intent(in) :: ugrs
98 real(kind_phys), dimension(:,:), intent(in) :: vgrs
99 real(kind_phys), dimension(:,:), intent(in) :: tgrs
100 real(kind_phys), dimension(:,:), intent(in) :: qgrs_wv
101 real(kind_phys), dimension(:,:), intent(in) :: qgrs_cw
102 real(kind_phys), dimension(:,:), intent(in) :: qgrs_rw
103 real(kind_phys), dimension(:,:), intent(in) :: qgrs_sw
104 real(kind_phys), dimension(:,:), intent(in) :: qgrs_iw
105 real(kind_phys), dimension(:,:), intent(in) :: qgrs_gl
106 real(kind_phys), dimension(:,:), intent(inout) :: gu0
107 real(kind_phys), dimension(:,:), intent(inout) :: gv0
108 real(kind_phys), dimension(:,:), intent(inout) :: gt0
109 real(kind_phys), dimension(:,:), intent(inout) :: gq0_wv
110 real(kind_phys), dimension(:,:), intent(inout) :: gq0_cw
111 real(kind_phys), dimension(:,:), intent(inout) :: gq0_rw
112 real(kind_phys), dimension(:,:), intent(inout) :: gq0_sw
113 real(kind_phys), dimension(:,:), intent(inout) :: gq0_iw
114 real(kind_phys), dimension(:,:), intent(inout) :: gq0_gl
115 integer, intent(in) :: ntcw
116 integer, intent(in) :: ntrw
117 integer, intent(in) :: ntsw
118 integer, intent(in) :: ntiw
119 integer, intent(in) :: ntgl
120 real(kind_phys), dimension(:,:), intent(inout), optional :: dtdtnp
121 real(kind_phys), dimension(:), intent(in) :: rain
122 real(kind_phys), dimension(:), intent(in) :: rainc
123 real(kind_phys), dimension(:), intent(inout) :: tprcp
124 real(kind_phys), dimension(:), intent(inout) :: totprcp
125 real(kind_phys), dimension(:), intent(inout) :: cnvprcp
126 real(kind_phys), dimension(:), intent(inout) :: totprcpb
127 real(kind_phys), dimension(:), intent(inout) :: cnvprcpb
128 logical, intent(in) :: cplflx
129 logical, intent(in) :: cpllnd
130 ! rain_cpl only allocated if cplflx == .true. or cplchm == .true. or cpllnd == .true.
131 real(kind_phys), dimension(:), intent(inout), optional :: rain_cpl
132 ! snow_cpl only allocated if cplflx == .true. or cplchm == .true.
133 real(kind_phys), dimension(:), intent(inout), optional :: snow_cpl
134 ! drain_cpl, dsnow_cpl only allocated if cplflx == .true. or cplchm == .true.
135 real(kind_phys), dimension(:), intent(in), optional :: drain_cpl
136 real(kind_phys), dimension(:), intent(in), optional :: dsnow_cpl
137 real(kind_phys), dimension(:), intent(in) :: vfact_ca
138 real(kind_phys), dimension(:), intent(in), optional :: ca1
139 character(len=*), intent(out) :: errmsg
140 integer, intent(out) :: errflg
141
142 !--- local variables
143 integer :: k, i
144 real(kind=kind_phys) :: upert, vpert, tpert, qpert, qnew, sppt_vwt
145 real(kind=kind_phys), dimension(1:im,1:km) :: ca
146
147 ! Initialize CCPP error handling variables
148 errmsg = ''
149 errflg = 0
150
151 if (do_sppt) then
152 do k=1,km
153 do i=1,im
154 sppt_vwt=1.0
155 if (zmtnblck(i).EQ.0.0) then
156 sppt_vwt=1.0
157 else
158 if (k.GT.zmtnblck(i)+2) then
159 sppt_vwt=1.0
160 endif
161 if (k.LE.zmtnblck(i)) then
162 sppt_vwt=0.0
163 endif
164 if (k.EQ.zmtnblck(i)+1) then
165 sppt_vwt=0.333333
166 endif
167 if (k.EQ.zmtnblck(i)+2) then
168 sppt_vwt=0.666667
169 endif
170 endif
171 if (use_zmtnblck)then
172 sppt_wts(i,k)=(sppt_wts(i,k)-1)*sppt_vwt+1.0
173 endif
174
175 upert = (gu0(i,k) - ugrs(i,k)) * sppt_wts(i,k)
176 vpert = (gv0(i,k) - vgrs(i,k)) * sppt_wts(i,k)
177 tpert = (gt0(i,k) - tgrs(i,k) - (delt*dtdtnp(i,k))) * sppt_wts(i,k)
178 qpert = (gq0_wv(i,k) - qgrs_wv(i,k)) * sppt_wts(i,k)
179
180 gu0(i,k) = ugrs(i,k)+upert
181 gv0(i,k) = vgrs(i,k)+vpert
182
183 !negative humidity check
184 qnew = qgrs_wv(i,k)+qpert
185 if (qnew >= 1.0e-10) then
186 gq0_wv(i,k) = qnew
187 gt0(i,k) = tgrs(i,k) + tpert + (delt*dtdtnp(i,k))
188 endif
189 if (pert_mp) then
190 if (ntcw>0) then
191 qpert = (gq0_cw(i,k) - qgrs_cw(i,k)) * sppt_wts(i,k)
192 qnew = qgrs_cw(i,k)+qpert
193 gq0_cw(i,k) = qnew
194 if (qnew < 0.0) then
195 gq0_cw(i,k) = 0.0
196 endif
197 endif
198 if (ntrw>0) then
199 qpert = (gq0_rw(i,k) - qgrs_rw(i,k)) * sppt_wts(i,k)
200 qnew = qgrs_rw(i,k)+qpert
201 gq0_rw(i,k) = qnew
202 if (qnew < 0.0) then
203 gq0_rw(i,k) = 0.0
204 endif
205 endif
206 if (ntsw>0) then
207 qpert = (gq0_sw(i,k) - qgrs_sw(i,k)) * sppt_wts(i,k)
208 qnew = qgrs_sw(i,k)+qpert
209 gq0_sw(i,k) = qnew
210 if (qnew < 0.0) then
211 gq0_sw(i,k) = 0.0
212 endif
213 endif
214 if (ntiw>0) then
215 qpert = (gq0_iw(i,k) - qgrs_iw(i,k)) * sppt_wts(i,k)
216 qnew = qgrs_iw(i,k)+qpert
217 gq0_iw(i,k) = qnew
218 if (qnew < 0.0) then
219 gq0_iw(i,k) = 0.0
220 endif
221 endif
222 if (ntgl>0) then
223 qpert = (gq0_gl(i,k) - qgrs_gl(i,k)) * sppt_wts(i,k)
224 qnew = qgrs_gl(i,k)+qpert
225 gq0_gl(i,k) = qnew
226 if (qnew < 0.0) then
227 gq0_gl(i,k) = 0.0
228 endif
229 endif
230 endif
231 enddo
232 enddo
233
234 ! instantaneous precip rate going into land model at the next time step
235 tprcp(:) = sppt_wts(:,15)*tprcp(:)
236 totprcp(:) = totprcp(:) + (sppt_wts(:,15) - 1 )*rain(:)
237 ! acccumulated total and convective preciptiation
238 cnvprcp(:) = cnvprcp(:) + (sppt_wts(:,15) - 1 )*rainc(:)
239 ! bucket precipitation adjustment due to sppt
240 totprcpb(:) = totprcpb(:) + (sppt_wts(:,15) - 1 )*rain(:)
241 cnvprcpb(:) = cnvprcpb(:) + (sppt_wts(:,15) - 1 )*rainc(:)
242
243 if (cplflx .or. cpllnd) then
244 rain_cpl(:) = rain_cpl(:) + (sppt_wts(:,15) - 1.0)*drain_cpl(:)
245 endif
246 if (cplflx) then
247 snow_cpl(:) = snow_cpl(:) + (sppt_wts(:,15) - 1.0)*dsnow_cpl(:)
248 endif
249 !zero out radiative heating tendency for next physics step
250 dtdtnp(:,:)=0.0
251
252 endif
253
254 if (do_ca .and. ca_global) then
255
256 !if(kdt == 1)then
257 !endif
258
259 do k = 1,km
260 do i = 1,im
261 sppt_vwt=1.0
262 if (zmtnblck(i).EQ.0.0) then
263 sppt_vwt=1.0
264 else
265 if (k.GT.zmtnblck(i)+2) then
266 sppt_vwt=1.0
267 endif
268 if (k.LE.zmtnblck(i)) then
269 sppt_vwt=0.0
270 endif
271 if (k.EQ.zmtnblck(i)+1) then
272 sppt_vwt=0.333333
273 endif
274 if (k.EQ.zmtnblck(i)+2) then
275 sppt_vwt=0.666667
276 endif
277 endif
278
279 ca(i,k)=((ca1(i)-1.)*sppt_vwt*vfact_ca(k))+1.0
280
281 upert = (gu0(i,k) - ugrs(i,k)) * ca(i,k)
282 vpert = (gv0(i,k) - vgrs(i,k)) * ca(i,k)
283 tpert = (gt0(i,k) - tgrs(i,k) - (delt*dtdtnp(i,k))) * ca(i,k)
284 qpert = (gq0_wv(i,k) - qgrs_wv(i,k)) * ca(i,k)
285 gu0(i,k) = ugrs(i,k)+upert
286 gv0(i,k) = vgrs(i,k)+vpert
287 !negative humidity check
288 qnew = qgrs_wv(i,k)+qpert
289 if (qnew >= 1.0e-10) then
290 gq0_wv(i,k) = qnew
291 gt0(i,k) = tgrs(i,k) + tpert + (delt*dtdtnp(i,k))
292 endif
293 if (pert_mp) then
294 if (ntcw>0) then
295 qpert = (gq0_cw(i,k) - qgrs_cw(i,k)) * ca(i,k)
296 qnew = qgrs_cw(i,k)+qpert
297 gq0_cw(i,k) = qnew
298 if (qnew < 0.0) then
299 gq0_cw(i,k) = 0.0
300 endif
301 endif
302 if (ntrw>0) then
303 qpert = (gq0_rw(i,k) - qgrs_rw(i,k)) * ca(i,k)
304 qnew = qgrs_rw(i,k)+qpert
305 gq0_rw(i,k) = qnew
306 if (qnew < 0.0) then
307 gq0_rw(i,k) = 0.0
308 endif
309 endif
310 if (ntsw>0) then
311 qpert = (gq0_sw(i,k) - qgrs_sw(i,k)) * ca(i,k)
312 qnew = qgrs_sw(i,k)+qpert
313 gq0_sw(i,k) = qnew
314 if (qnew < 0.0) then
315 gq0_sw(i,k) = 0.0
316 endif
317 endif
318 if (ntiw>0) then
319 qpert = (gq0_iw(i,k) - qgrs_iw(i,k)) * ca(i,k)
320 qnew = qgrs_iw(i,k)+qpert
321 gq0_iw(i,k) = qnew
322 if (qnew < 0.0) then
323 gq0_iw(i,k) = 0.0
324 endif
325 endif
326 if (ntgl>0) then
327 qpert = (gq0_gl(i,k) - qgrs_gl(i,k)) * ca(i,k)
328 qnew = qgrs_gl(i,k)+qpert
329 gq0_gl(i,k) = qnew
330 if (qnew < 0.0) then
331 gq0_gl(i,k) = 0.0
332 endif
333 endif
334 endif
335 enddo
336 enddo
337
338 ! instantaneous precip rate going into land model at the next time step
339 tprcp(:) = ca(:,15)*tprcp(:)
340 totprcp(:) = totprcp(:) + (ca(:,15) - 1 )*rain(:)
341 ! acccumulated total and convective preciptiation
342 cnvprcp(:) = cnvprcp(:) + (ca(:,15) - 1 )*rainc(:)
343 ! bucket precipitation adjustment due to sppt
344 totprcpb(:) = totprcpb(:) + (ca(:,15) - 1 )*rain(:)
345 cnvprcpb(:) = cnvprcpb(:) + (ca(:,15) - 1 )*rainc(:)
346
347 if (cplflx .or. cpllnd) then
348 rain_cpl(:) = rain_cpl(:) + (ca(:,15) - 1.0)*drain_cpl(:)
349 endif
350 if (cplflx) then
351 snow_cpl(:) = snow_cpl(:) + (ca(:,15) - 1.0)*dsnow_cpl(:)
352 endif
353 !zero out radiative heating tendency for next physics step
354 dtdtnp(:,:)=0.0
355
356
357 endif
358
359 if (do_shum) then
360 do k=1,km
361 gq0_wv(:,k) = gq0_wv(:,k)*(1.0 + shum_wts(:,k))
362 end do
363 endif
364
365 if (do_skeb) then
366 do k=1,km
367 gu0(:,k) = gu0(:,k)+skebu_wts(:,k)*(diss_est(:,k))
368 gv0(:,k) = gv0(:,k)+skebv_wts(:,k)*(diss_est(:,k))
369 enddo
370 endif
371
372 end subroutine gfs_stochastics_run
373 end module gfs_stochastics