CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
myjsfc_wrapper.F90
1
3
5
6 USE machine, only: kfpt => kind_phys, &
7 kind_phys
8
9 contains
10
11 subroutine myjsfc_wrapper_init (do_myjsfc, &
12 & errmsg,errflg)
13
14 logical, intent(in) :: do_myjsfc
15 character(len=*), intent(out) :: errmsg
16 integer, intent(out) :: errflg
17
18 ! Initialize CCPP error handling variables
19 errmsg = ''
20 errflg = 0
21
22 ! Consistency checks
23 if (.not. do_myjsfc) then
24 write(errmsg,fmt='(*(a))') 'Logic error: do_myjsfc = .false.'
25 errflg = 1
26 return
27 end if
28 end subroutine myjsfc_wrapper_init
29
30!!
35!###===================================================================
36 SUBROUTINE myjsfc_wrapper_run( &
37 restart, &
38 im,levs, &
39 kdt,ntrac,ntke, &
40 ntcw,ntiw,ntrw,ntsw,ntgl, &
41 iter,flag_iter, &
42 ugrs, vgrs, tgrs, qgrs, &
43 prsl, prsi, phii, &
44 prsik_1, prslk_1, tsfc, qsfc, &
45 phy_myj_qsfc, phy_myj_thz0, phy_myj_qz0, &
46 phy_myj_uz0, phy_myj_vz0, phy_myj_z0base, &
47 phy_myj_akhs, phy_myj_akms, &
48 phy_myj_chkqlm, phy_myj_elflx, &
49 phy_myj_a1u, phy_myj_a1t, phy_myj_a1q, &
50 pblh, slmsk, zorl, ustar, rib, &
51 cm,ch,stress,ffm,ffh,fm10,fh2, &
52 landfrac, oceanfrac,fice, &
53 z0rl_wat, z0rl_lnd, z0rl_ice, & ! intent(inout)
54 ustar_wat, ustar_lnd, ustar_ice, & ! intent(inout)
55 cm_wat, cm_lnd, cm_ice, & ! intent(inout)
56 ch_wat, ch_lnd, ch_ice, & ! intent(inout)
57 rb_wat, rb_lnd, rb_ice, & ! intent(inout)
58 stress_wat,stress_lnd,stress_ice, & ! intent(inout)
59 fm_wat, fm_lnd, fm_ice, & ! intent(inout)
60 fh_wat, fh_lnd, fh_ice, & ! intent(inout)
61 fm10_wat, fm10_lnd, fm10_ice, & ! intent(inout)
62 fh2_wat, fh2_lnd, fh2_ice, & ! intent(inout)
63 wind, con_cp, con_g, con_rd, &
64 me, lprnt, errmsg, errflg ) ! intent(inout)
65!
66
67 use module_sf_jsfc, only: jsfc_init,jsfc
68
69!-------------------------------------------------------------------
70 implicit none
71!-------------------------------------------------------------------
72
73! integer,parameter:: &
74! klog=4 & ! logical variables
75! ,kint=4 & ! integer variables
76! !,kfpt=4 & ! floating point variables
77! ,kfpt=8 & ! floating point variables
78! ,kdbl=8 ! double precision
79!
80! --- constant parameters:
81! real(kind=kind_phys), parameter :: karman = 0.4
82
83!-------------------------------------------------------------------
84!-------------------------------------------------------------------
85!For reference
86! real , parameter :: karman = 0.4
87! real , parameter :: g = 9.81
88! real , parameter :: r_d = 287.
89! real , parameter :: cp = 7.*r_d/2.
90! real , parameter :: r_v = 461.6
91! real , parameter :: cpv = 4.*r_v
92! real , parameter :: rcp = r_d/cp
93
94! real, parameter :: g_inv=1/g, cappa=r_d/cp
95
96 character(len=*), intent(out) :: errmsg
97 integer, intent(out) :: errflg
98
99!MYJ-1D
100 integer,intent(in) :: im, levs
101 integer,intent(in) :: kdt, iter, me
102 integer,intent(in) :: ntrac,ntke,ntcw,ntiw,ntrw,ntsw,ntgl
103 logical,intent(in) :: restart, lprnt
104 real(kind=kind_phys),intent(in) :: con_cp, con_g, con_rd
105
106!MYJ-2D
107 logical,dimension(:),intent(in) :: flag_iter
108 real(kind=kind_phys),dimension(:),intent(in) :: &
109 prsik_1, prslk_1, tsfc, qsfc, slmsk
110 real(kind=kind_phys),dimension(:),intent(inout),optional :: &
111 phy_myj_thz0, phy_myj_z0base, phy_myj_chkqlm, &
112 phy_myj_akhs, phy_myj_akms, phy_myj_qz0, &
113 phy_myj_qsfc, phy_myj_elflx, phy_myj_a1u, &
114 phy_myj_a1t, phy_myj_a1q, phy_myj_uz0, &
115 phy_myj_vz0
116 real(kind=kind_phys),dimension(:),intent(inout) :: &
117 pblh, zorl, ustar, rib
118 real(kind=kind_phys),dimension(:),intent(inout) :: &
119 cm, ch, stress, ffm, ffh, fm10, fh2
120 real(kind=kind_phys), dimension(:), intent(inout) :: &
121 landfrac, oceanfrac, fice
122 real(kind=kind_phys), dimension(:), intent(inout) :: &
123 z0rl_wat, z0rl_lnd, z0rl_ice, &
124 ustar_wat, ustar_lnd, ustar_ice, &
125 cm_wat, cm_lnd, cm_ice, &
126 ch_wat, ch_lnd, ch_ice, &
127 rb_wat, rb_lnd, rb_ice, &
128 stress_wat,stress_lnd,stress_ice, &
129 fm_wat, fm_lnd, fm_ice, &
130 fh_wat, fh_lnd, fh_ice, &
131 fm10_wat, fm10_lnd, fm10_ice, &
132 fh2_wat, fh2_lnd, fh2_ice, &
133 wind
134
135
136!MYJ-3D
137 real(kind=kind_phys),dimension(:,:),intent(in) :: &
138 phii, prsi
139 real(kind=kind_phys),dimension(:,:),intent(in) :: &
140 ugrs, vgrs, tgrs, prsl
141!MYJ-4D
142 real(kind=kind_phys),dimension(:,:,:),intent(in) :: &
143 qgrs
144
145!LOCAL
146 logical :: lprnt1, lprnt2
147 integer :: ntsd, k, k1, i, n, ide, jde, kde
148
149 real(kind=kind_phys) :: g, r_d, g_inv, cappa
150 real(kind=kfpt),dimension(levs) :: epsq2
151 real(kind=kfpt),dimension(im) :: &
152 sfcz,tsk,xland,mavail,rmol, &
153 ustar1,z0,rib1,sm,pblh_myj
154 real(kind=kfpt),dimension(im,13) :: &
155 phy_f2d_myj
156 real(kind=kfpt), dimension(im,levs) :: &
157 u_myj, v_myj, t_myj, q_myj, th_myj, &
158 cw, dz_myj, pmid, q2, exner
159 real(kind=kfpt), dimension(im,levs+1) :: pint
160 real(kind=kfpt),dimension(im) :: &
161 cm1,ch1,stress1,ffm1,ffh1,wind1,ffm10,ffh2
162! real(kind=kind_phys), dimension(im,levs,ntrac) :: &
163! & qgrs_myj
164
165 ! Initialize CCPP error handling variables
166 errmsg = ''
167 errflg = 0
168
169 ntsd = kdt-1
170
171 lprnt1 =.false.
172 lprnt2 =.false.
173
174 if (lprnt2) then
175 if(me.eq.0)then
176 print*,'in myj surface layer wrapper...'
177 print*,'ntsd,iter=',ntsd,iter
178 end if
179 endif
180
181 r_d = con_rd
182 g = con_g
183 g_inv = 1./con_g
184 cappa = con_rd/con_cp
185
186 if (ntsd==0.and.iter==1)then
187 do i=1,im
188 if(flag_iter(i))then
189 phy_myj_qsfc(i) = qgrs(i,1,1) ! qsfc(:)
190 phy_myj_thz0(i) = tsfc(i) ! thz0
191 phy_myj_qz0(i) = qgrs(i,1,1) ! qz0(:)
192 phy_myj_uz0(i) = 0. ! uz0(:)
193 phy_myj_vz0(i) = 0. ! vz0(:)
194 phy_myj_z0base(i) = zorl(i)*0.01 ! z0base
195 phy_myj_akhs(i) = 0.01 ! akhs(:)
196 phy_myj_akms(i) = 0.01 ! akms(:)
197 phy_myj_chkqlm(i) = 0. ! chkqlm(:)
198 phy_myj_elflx(i) = 0. ! elflx(:)
199 phy_myj_a1u(i) = 0. ! a1u
200 phy_myj_a1t(i) = 0. ! a1t
201 phy_myj_a1q(i) = 0. ! a1q
202 end if
203 end do
204 end if
205
206!prep MYJ-only variables
207 do i=1,im
208 sm(i)=1.; if(slmsk(i) > 0.5 ) sm(i)=0.
209 xland(i)=sm(i)+1.
210 sfcz(i)=phii(i,1)*g_inv
211 enddo
212
213 do k=1,levs
214 k1=levs+1-k
215 do i=1,im
216 u_myj(i,k)=ugrs(i,k1)
217 v_myj(i,k)=vgrs(i,k1)
218 t_myj(i,k)=tgrs(i,k1)
219 q_myj(i,k)=qgrs(i,k1,1)
220 cw(i,k) =qgrs(i,k1,ntcw)
221! if(ntrw.gt.0)cw(i,k) = cw(i,k) + qgrs(i,k1,ntrw)
222! if(ntiw.gt.0)cw(i,k) = cw(i,k) + qgrs(i,k1,ntiw)
223! if(ntsw.gt.0)cw(i,k) = cw(i,k) + qgrs(i,k1,ntsw)
224! if(ntgl.gt.0)cw(i,k) = cw(i,k) + qgrs(i,k1,ntgl)
225 if(ntke.gt.0)then
226 q2(i,k) =qgrs(i,k1,ntke)*2.
227 else
228 q2(i,k) =0.02
229 end if
230 pmid(i,k) =prsl(i,k1)
231 exner(i,k)=(prsl(i,k1)*1.e-5)**cappa
232 th_myj(i,k)=tgrs(i,k1)/exner(i,k)
233 end do
234 end do
235 do k=1,levs+1
236 k1=levs+2-k
237 do i=1,im
238 pint(i,k) =prsi(i,k1)
239 end do
240 end do
241
242 do k = 1, levs
243 k1 = levs-k+1
244 do i = 1, im
245 dz_myj(i,k) = (phii(i,k1+1)-phii(i,k1)) * g_inv
246 enddo
247 enddo
248
249 if (lprnt1) then
250 if(me==0.and.ntsd.lt.2)then
251 k=63
252 k1=levs+1-k
253 print*,'Qingfu starts MYJSFC'
254 print*,'ntsd,iter,me,1=',ntsd,iter,me
255 print*,'ntrac,ntcw,ntiw,ntrw,ntsw,ntgl,ntke=', &
256 ntrac,ntcw,ntiw,ntrw,ntsw,ntgl,ntke
257 print*,'im,levs,ntsd=',im,levs,ntsd
258 do i=10,40,40
259 print*,'Qingfu before MYJ surface kdt,i,k1=',kdt,i,k1
260 print*,'sfcz,dz_myj,th_myj,tsfc,qsfc=',sfcz(i),dz_myj(i,k), &
261 th_myj(i,k),tsfc(i),qsfc(i)
262 print*,'sm,z0,xland=', &
263 sm(i),z0(i),xland(i)
264! print*,'phy_f2d_myj(i,1:13)=', &
265! (phy_f2d_myj(i,n),n=1,13)
266 print*,'u_myj,v_myj=', &
267 u_myj(i,k),v_myj(i,k)
268 print*,'t_myj,q_myj,cw,q2=', &
269 t_myj(i,k),q_myj(i,k),cw(i,k),q2(i,k)
270 print*,'phii,pint,pmid', &
271 phii(i,k1),pint(i,k),pmid(i,k)
272 print*,'exner,th_myj=',exner(i,k),th_myj(i,k)
273 end do
274 end if
275 endif
276
277!-----------------------------------------------------------------------
278 ide=im+1
279 jde=2
280 kde=levs+1
281
282 do i = 1, im
283 epsq2(i)=0.02
284 mavail(i)=1.0
285 tsk(i)=tsfc(i)
286 phy_f2d_myj(i,1) = phy_myj_qsfc(i)
287 phy_f2d_myj(i,2) = phy_myj_thz0(i)
288 phy_f2d_myj(i,3) = phy_myj_qz0(i)
289 phy_f2d_myj(i,4) = phy_myj_uz0(i)
290 phy_f2d_myj(i,5) = phy_myj_vz0(i)
291 phy_f2d_myj(i,6) = phy_myj_z0base(i)
292 phy_f2d_myj(i,7) = phy_myj_akhs(i)
293 phy_f2d_myj(i,8) = phy_myj_akms(i)
294 phy_f2d_myj(i,9) = phy_myj_chkqlm(i)
295 phy_f2d_myj(i,10) = phy_myj_elflx(i)
296 phy_f2d_myj(i,11) = phy_myj_a1u(i)
297 phy_f2d_myj(i,12) = phy_myj_a1t(i)
298 phy_f2d_myj(i,13) = phy_myj_a1q(i)
299 z0(i)=zorl(i)*0.01
300 rmol(i)=0.
301 rib1(i)=rib(i)
302 pblh_myj(i)=pblh(i)
303 ustar1(i)=ustar(i)
304 cm1(i)=0.
305 ch1(i)=0.
306 stress1(i)=0.
307 ffm1(i)=0.
308 ffh1(i)=0.
309 wind1(i)=0.
310 ffm10(i)=0.
311 ffh2(i)=0.
312 end do
313
314 if((ntsd==0.and.iter.eq.1).or.restart)then
315 call jsfc_init(ustar1,restart &
316 ,1,ide,1,jde,1,kde &
317 ,1,im,1,1,1,levs &
318 ,1,im,1,1,1,levs)
319 end if
320
321 call jsfc(flag_iter,iter,me &
322 ,ntsd,epsq2,sfcz,dz_myj &
323 ,pmid,pint,th_myj,t_myj,q_myj,cw &
324 ,u_myj,v_myj,q2,tsk &
325 ,phy_f2d_myj(1:im,1),phy_f2d_myj(1:im,2) &
326 ,phy_f2d_myj(1:im,3),phy_f2d_myj(1:im,4) &
327 ,phy_f2d_myj(1:im,5),xland &
328 ,ustar1,z0,phy_f2d_myj(1:im,6) &
329 ,pblh_myj,mavail,rmol &
330 ,phy_f2d_myj(1:im,7),phy_f2d_myj(1:im,8) &
331 ,phy_f2d_myj(1:im,9),phy_f2d_myj(1:im,10) &
332 ,rib1,cm1,ch1,stress1,ffm1,ffh1,wind1,ffm10,ffh2 &
333 ,phy_f2d_myj(1:im,11),phy_f2d_myj(1:im,12) &
334 ,phy_f2d_myj(1:im,13) &
335 ,1,im,1,1,1,levs &
336 ,1,im,1,1,1,levs &
337 ,1,im,1,1,1,levs, errmsg, errflg)
338
339 do i = 1, im
340 if(flag_iter(i))then
341 zorl(i) = z0(i)*100.
342
343 phy_myj_qsfc(i) = phy_f2d_myj(i,1)
344 phy_myj_thz0(i) = phy_f2d_myj(i,2)
345 phy_myj_qz0(i) = phy_f2d_myj(i,3)
346 phy_myj_uz0(i) = phy_f2d_myj(i,4)
347 phy_myj_vz0(i) = phy_f2d_myj(i,5)
348 phy_myj_z0base(i) = phy_f2d_myj(i,6)
349 phy_myj_akhs(i) = phy_f2d_myj(i,7)
350 phy_myj_akms(i) = phy_f2d_myj(i,8)
351 phy_myj_chkqlm(i) = phy_f2d_myj(i,9)
352 phy_myj_elflx(i) = - phy_f2d_myj(i,10) ! change flux definition
353 phy_myj_a1u(i) = phy_f2d_myj(i,11)
354 phy_myj_a1t(i) = phy_f2d_myj(i,12)
355 phy_myj_a1q(i) = phy_f2d_myj(i,13)
356
357 rib(i)=rib1(i)
358 pblh(i)=pblh_myj(i)
359 cm(i)=cm1(i)
360 ch(i)=ch1(i)
361 stress(i)=stress1(i)
362 ffm(i)=ffm1(i)
363 ffh(i)=ffh1(i)
364 wind(i)=wind1(i)
365 fm10(i)=ffm10(i)
366 fh2(i)=ffh2(i)
367 ustar(i)=ustar1(i)
368 end if
369 end do
370
371 if (lprnt1) then
372
373 if(me==0.and.ntsd.lt.10)then
374 print*,'ntsd,iter,me,2=',ntsd,iter,me
375 do i=10,40,40
376 if(flag_iter(i))then
377 print*,'Qingfu after MYJ surface kdt,i,k1=',kdt,i,k1
378 print*,'xland,cm,ch=',xland(i),cm(i),ch(i)
379 print*,'ustar,z0,stress=',ustar(i),z0(i),stress(i)
380 print*,'ffm,ffh,wind,fm10,fh2=',ffm(i),ffh(i),wind(i),fm10(i),fh2(i)
381 print*,'phy_f2d_myj(9,1:13)=', &
382 (phy_f2d_myj(i,n),n=1,13)
383 print*,'u_myj,v_myj=', &
384 u_myj(i,k),v_myj(i,k)
385 print*,'t_myj,q_myj,cw,q2=', &
386 t_myj(i,k),q_myj(i,k),cw(i,k),q2(i,k)
387 print*,'phii,pint,pmid', &
388 phii(i,k1),pint(i,k),pmid(i,k)
389 print*,'exner,th_myj=',exner(i,k),th_myj(i,k)
390 print*,'Qingfu finish MYJSFC'
391 end if
392 end do
393 end if
394
395 do k=1,levs
396 k1=levs+1-k
397 do i=1,im
398 if(t_myj(i,k).gt.320..or.t_myj(i,k).lt.150.)then
399 print*,'xland,cm,ch=',xland(i),cm(i),ch(i)
400 print*,'ustar,z0,stress=',ustar(i),z0(i),stress(i)
401 print*,'ffm,ffh,wind,fm10,fh2=',ffm(i),ffh(i),wind(i),fm10(i),fh2(i)
402 print*,'phy_f2d_myj(9,1:13)=', &
403 (phy_f2d_myj(i,n),n=1,13)
404 print*,'u_myj,v_myj=', &
405 u_myj(i,k),v_myj(i,k)
406 print*,'t_myj,q_myj,cw,q2=', &
407 t_myj(i,k),q_myj(i,k),cw(i,k),q2(i,k)
408 print*,'phii,pint,pmid', &
409 phii(i,k1),pint(i,k),pmid(i,k)
410 print*,'exner,th_myj=',exner(i,k),th_myj(i,k)
411 print*,'Qingfu finish MYJSFC'
412 end if
413 end do
414 end do
415
416 end if
417
418 do i = 1, im
419 if(flag_iter(i))then
420 z0rl_wat(i) = zorl(i)
421 cm_wat(i) = cm(i)
422 ch_wat(i) = ch(i)
423 rb_wat(i) = rib(i)
424 stress_wat(i) = stress(i)
425 fm_wat(i) = ffm(i)
426 fh_wat(i) = ffh(i)
427 ustar_wat(i) = ustar(i)
428 fm10_wat(i) = fm10(i)
429 fh2_wat(i) = fh2(i)
430
431 z0rl_lnd(i) = zorl(i)
432 cm_lnd(i) = cm(i)
433 ch_lnd(i) = ch(i)
434 rb_lnd(i) = rib(i)
435 stress_lnd(i) = stress(i)
436 fm_lnd(i) = ffm(i)
437 fh_lnd(i) = ffh(i)
438 ustar_lnd(i) = ustar(i)
439 fm10_lnd(i) = fm10(i)
440 fh2_lnd(i) = fh2(i)
441
442 z0rl_ice(i) = zorl(i)
443 cm_ice(i) = cm(i)
444 ch_ice(i) = ch(i)
445 rb_ice(i) = rib(i)
446 stress_ice(i) = stress(i)
447 fm_ice(i) = ffm(i)
448 fh_ice(i) = ffh(i)
449 ustar_ice(i) = ustar(i)
450 fm10_ice(i) = fm10(i)
451 fh2_ice(i) = fh2(i)
452 end if
453 end do
454
455 if (lprnt2) then
456 if(me==0.and.ntsd.lt.10)then
457 print*,'ntsd,iter,me,3=',ntsd,iter,me
458 do i=10,40,40
459 if(flag_iter(i))then
460 print*,'Qingfu after MYJ surface kdt,i,k1,3=',kdt,i,k1
461 print*,'Qingfu test after MYJ surface kdt,i=',kdt,i,slmsk(i)
462 print*,'a1u,a1t,a1q=',(phy_f2d_myj(i,k),k=11,13)
463 print*,'zorl,cm,ch,rb,stress=',z0(i), &
464 cm(i),ch(i), &
465 rib(i),stress(i)
466 print*,'ffmm,ffhh,ustar,fm10,fh2,wind=', ffm(i), &
467 ffh(i),ustar(i),fm10(i),fh2(i),wind(i)
468 print*,'cm(i),ch(i)=', &
469 (0.4/ffm(i))**2,(0.4/ffm(i)*0.4/ffh(i))
470 end if
471 end do
472 endif
473 endif
474
475
476 END SUBROUTINE myjsfc_wrapper_run
477
478!###=================================================================
479
480END MODULE myjsfc_wrapper