/[MITgcm]/MITgcm/pkg/fizhi/fizhi_step_diag.F
ViewVC logotype

Annotation of /MITgcm/pkg/fizhi/fizhi_step_diag.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (hide annotations) (download)
Thu Jun 24 15:06:51 2004 UTC (20 years ago) by molod
Branch: MAIN
Code to fill some fizhi diagnostics from do_fizhi

1 molod 1.1 C $Header: $
2     C $Name: $
3    
4     #include "CPP_OPTIONS.h"
5     subroutine fizhi_step_diag(myThid,p,uphy,vphy,thphy,sphy,qq,
6     . radswt,radswg,swgclr,osr,osrclr,st4,dst4,tgz,tg0,radlwg,lwgclr,
7     . turbu,turbv,turbt,turbq,moistu,moistv,moistt,moistq,
8     . lwdt,swdt,lwdtclr,swdtclr,im1,im2,jm1,jm2,Nrphys,Nsx,Nsy,bi,bj)
9     C***********************************************************************
10     implicit none
11    
12     #include "diagnostics.h"
13    
14     integer myThid,im1,im2,jm1,jm2,Nrphys,Nsx,Nsy,bi,bj
15     real radswt(im2,jm2)
16    
17     integer i,j,L,m
18     real getcon
19     real cp, pstd, tstd, akap, pkstd, thstd, grav, delp
20    
21     cp = getcon('CP')
22     pstd = getcon('PSTD')
23     tstd = getcon('TSTD')
24     akap = getcon('KAPPA')
25     pkstd = pstd**akap
26     thstd = tstd/pkstd
27     grav = getcon('GRAVITY')
28    
29     C **********************************************************************
30     C **** Compute 2-D Diagnostics ****
31     C **********************************************************************
32    
33     do j=jm1,jm2
34     do i=im1,im2
35     pinv(i,j) = 1.0 / p(i,j)
36     enddo
37     enddo
38    
39     c Analysis Increment of Surface Pressure (mb/day)
40     c -----------------------------------------------
41     if( ipiau.ne.0 ) then
42     do j=jm1,jm2
43     do i=im1,im2
44     qdiag(i,j,ipiau) = qdiag(i,j,ipiau) + tend%iau%dp(i,j)*86400
45     enddo
46     enddo
47     endif
48    
49     c Incident Solar Radiation (W/m**2)
50     c ---------------------------------
51     if (iradswt.ne.0) then
52     do j=jm1,jm2
53     do i=im1,im2
54     qdiag(i,j,iradswt) = qdiag(i,j,iradswt) + radswt(i,j)
55     enddo
56     enddo
57     endif
58    
59     c Net Solar Radiation at the Ground (W/m**2)
60     c ------------------------------------------
61     if (iradswg.ne.0) then
62     do j=jm1,jm2
63     do i=im1,im2
64     qdiag(i,j,iradswg) = qdiag(i,j,iradswg) + coup%sw%radswg(i,j)*radswt(i,j)
65     enddo
66     enddo
67     endif
68    
69     c Net Clear Sky Solar Radiation at the Ground (W/m**2)
70     c ----------------------------------------------------
71     if (iswgclr.ne.0) then
72     do j=jm1,jm2
73     do i=im1,im2
74     qdiag(i,j,iswgclr) = qdiag(i,j,iswgclr) + coup%sw%swgclr(i,j)*radswt(i,j)
75     enddo
76     enddo
77     endif
78    
79     c Outgoing Solar Radiation at Ptop (W/m**2)
80     c -----------------------------------------
81     if (iosr.ne.0) then
82     do j=jm1,jm2
83     do i=im1,im2
84     qdiag(i,j,iosr) = qdiag(i,j,iosr) + (1.0-coup%sw%osr(i,j))*radswt(i,j)
85     enddo
86     enddo
87     endif
88    
89     c Outgoing Clear Sky Solar Radiation at Ptop (W/m**2)
90     c ---------------------------------------------------
91     if (iosrclr.ne.0) then
92     do j=jm1,jm2
93     do i=im1,im2
94     qdiag(i,j,iosrclr) = qdiag(i,j,iosrclr) + (1.0-coup%sw%osrclr(i,j))*radswt(i,j)
95     enddo
96     enddo
97     endif
98    
99     c Upward Longwave Flux at the Ground (W/m**2)
100     c -------------------------------------------
101     if (ilwgup.ne.0) then
102     do j=jm1,jm2
103     do i=im1,im2
104     qdiag(i,j,ilwgup) = qdiag(i,j,ilwgup) + coup%lw%st4(i,j)
105     . + coup%lw%dst4(i,j)*(coup%land%tgz(i,j)-coup%lw%tg0(i,j))
106     enddo
107     enddo
108     endif
109    
110     c Net Longwave Flux at the Ground (W/m**2)
111     c ----------------------------------------
112     if (iradlwg.ne.0) then
113     do j=jm1,jm2
114     do i=im1,im2
115     qdiag(i,j,iradlwg) = qdiag(i,j,iradlwg) + coup%lw%radlwg(i,j)
116     . + coup%lw%dst4(i,j)*(coup%land%tgz(i,j)-coup%lw%tg0(i,j))
117     enddo
118     enddo
119     endif
120    
121     c Net Longwave Flux at the Ground Clear Sky (W/m**2)
122     c --------------------------------------------------
123     if (ilwgclr.ne.0) then
124     do j=jm1,jm2
125     do i=im1,im2
126     qdiag(i,j,ilwgclr) = qdiag(i,j,ilwgclr) + coup%lw%lwgclr(i,j)
127     . + coup%lw%dst4(i,j)*(coup%land%tgz(i,j)-coup%lw%tg0(i,j))
128     enddo
129     enddo
130     endif
131    
132     c Total Surface Pressure Tendency (mb/day)
133     c ----------------------------------------
134     if( idpdt.ne.0 ) then
135     do j=jm1,jm2
136     do i=im1,im2
137     qdiag(i,j,idpdt) = qdiag(i,j,idpdt) + (updt%p(i,j)-curr%p(i,j))*fact
138     enddo
139     enddo
140     endif
141    
142     c Averaged P-Field (mb)
143     c ---------------------
144     if( ips.ne.0 ) then
145     do j=jm1,jm2
146     do i=im1,im2
147     qdiag(i,j,ips) = qdiag(i,j,ips) + curr%p(i,j)+ptop
148     enddo
149     enddo
150     endif
151    
152     c Averaged SLP-Field (mb)
153     c ----------------------
154     if( islp.ne.0 ) then
155     do L=1,Nrphys
156     do j=jm1,jm2
157     do i=im1,im2
158     tmp2(i,j,L) = curr%t(i,j,L) * (1.+0.609*curr%q(i,j,L,1))
159     enddo
160     enddo
161     enddo
162     call slprs ( tmp1(1,1,2),curr%p,ptop,coup%earth%phis_cmp,tmp2,dsig,coup%earth%lw_cmp,im,jm,lm )
163     do j=jm1,jm2
164     do i=im1,im2
165     qdiag(i,j,islp) = qdiag(i,j,islp) + tmp1(i,j,2)
166     enddo
167     enddo
168     endif
169    
170     npiau = npiau + 1
171     nradswt = nradswt + 1
172     nradswg = nradswg + 1
173     nswgclr = nswgclr + 1
174     nosr = nosr + 1
175     nosrclr = nosrclr + 1
176     nradlwg = nradlwg + 1
177     nlwgclr = nlwgclr + 1
178     nlwgup = nlwgup + 1
179     ndpdt = ndpdt + 1
180     nps = nps + 1
181     nslp = nslp + 1
182    
183     C **********************************************************************
184     C **** Compute 3-D Diagnostics ****
185     C **********************************************************************
186    
187     do L=1,Nrphys
188    
189     if( itiau.ne.0 .or. ivavetiau.ne.0 .or. idiabt.ne.0 ) then
190     do j=jm1,jm2
191     do i=im1,im2
192     dtiau(i,j,L) = ( tend%iau%dt(i,j,L) + curr%t(i,j,L)*tend%iau%dp(i,j)
193     . * ( sig(L)*akap*curr%p(i,j)/(sig(L)*curr%p(i,j)+ptop) - 1.0) )
194     enddo
195     enddo
196     endif
197    
198     if( iqiau.ne.0 .or. ivaveqiau.ne.0 .or. idiabq.ne.0 ) then
199     do j=jm1,jm2
200     do i=im1,im2
201     dqiau(i,j,L) = tend%iau%dq(i,j,L) - curr%q(i,j,L,1)*tend%iau%dp(i,j)
202     enddo
203     enddo
204     endif
205    
206     c Total Diabatic U-Tendency (m/sec/day)
207     c -------------------------------------
208     if( idiabu.ne.0 ) then
209     do j=jm1,jm2
210     do i=im1,im2
211     qdiag(i,j,idiabu+L-1) = qdiag (i,j,idiabu+L-1)
212     . + ( tend%iau%du (i,j,L)
213     . + tend%turb%du(i,j,L) )*86400
214     enddo
215     enddo
216     endif
217    
218     c Total Diabatic V-Tendency (m/sec/day)
219     c -------------------------------------
220     if( idiabv.ne.0 ) then
221     do j=jm1,jm2
222     do i=im1,im2
223     qdiag(i,j,idiabv+L-1) = qdiag (i,j,idiabv+L-1)
224     . + ( tend%iau%dv (i,j,L)
225     . + tend%turb%dv(i,j,L) )*86400
226     enddo
227     enddo
228     endif
229    
230     c Total Diabatic T-Tendency (deg/day)
231     c -----------------------------------
232     if( idiabt.ne.0 ) then
233     do j=jm1,jm2
234     do i=im1,im2
235     qdiag(i,j,idiabt+L-1) = qdiag(i,j,idiabt+L-1)
236     . + ( dtiau(i,j,L) + tend%turb%dt(i,j,L) + tend%lw%dt(i,j,L)
237     . + coup%lw%dlwdtg(i,j,L)*(coup%land%tgc(i,j)-coup%lw%tg0c(i,j))
238     . + tend%sw%dt(i,j,L)*radswc(i,j)
239     . + tend%moist%dt(i,j,L) )
240     . * pkn(i,j,L)*pinv(i,j)*86400
241     enddo
242     enddo
243     endif
244    
245     c Total Diabatic Q-Tendency (g/kg/day)
246     c ------------------------------------
247     if( idiabq.ne.0 ) then
248     do j=jm1,jm2
249     do i=im1,im2
250     qdiag(i,j,idiabq+L-1) = qdiag(i,j,idiabq+L-1)
251     . + ( dqiau(i,j,L) + tend%turb%dq(i,j,L,1) + tend%moist%dq(i,j,L,1) )
252     . * pinv(i,j)*86400*1000
253     enddo
254     enddo
255     endif
256    
257     c Analysis U-Wind Increment (m/sec/day)
258     c -------------------------------------
259     if( iuiau.ne.0 ) then
260     do j=jm1,jm2
261     do i=im1,im2
262     qdiag(i,j,iuiau+L-1) = qdiag (i,j,iuiau+L-1) + tend%iau%du(i,j,L)*86400
263     enddo
264     enddo
265     endif
266    
267     c Analysis V-Wind Increment (m/sec/day)
268     c -------------------------------------
269     if( iviau.ne.0 ) then
270     do j=jm1,jm2
271     do i=im1,im2
272     qdiag(i,j,iviau+L-1) = qdiag (i,j,iviau+L-1) + tend%iau%dv(i,j,L)*86400
273     enddo
274     enddo
275     endif
276    
277     c Analysis Temperature Tendency (deg/day)
278     c ---------------------------------------
279     if( itiau.ne.0 ) then
280     do j=jm1,jm2
281     do i=im1,im2
282     qdiag(i,j,itiau+L-1) = qdiag(i,j,itiau+L-1)
283     . + dtiau(i,j,L) * pkn(i,j,L)*pinv(i,j)*86400
284     enddo
285     enddo
286     endif
287    
288     c Analysis Moisture Tendency (g/kg/day)
289     c -------------------------------------
290     if( iqiau.ne.0 ) then
291     do j=jm1,jm2
292     do i=im1,im2
293     qdiag(i,j,iqiau+L-1) = qdiag(i,j,iqiau+L-1)
294     . + dqiau(i,j,L) * pinv(i,j)*86400*1000
295     enddo
296     enddo
297     endif
298    
299     c Longwave Heating (deg/day)
300     c --------------------------
301     if (iradlw.ne.0) then
302     do j=jm1,jm2
303     do i=im1,im2
304     qdiag(i,j,iradlw+l-1) = qdiag(i,j,iradlw+l-1)
305     . + ( tend%lw%dt(i,j,l)
306     . + coup%lw%dlwdtg (i,j,L)*(coup%land%tgc(i,j)-coup%lw%tg0c(i,j)))
307     . * pkn(i,j,l)*pinv(i,j)*86400
308     enddo
309     enddo
310     endif
311    
312     c Longwave Heating Clear-Sky (deg/day)
313     c ------------------------------------
314     if (ilwclr.ne.0) then
315     do j=jm1,jm2
316     do i=im1,im2
317     qdiag(i,j,ilwclr+l-1) = qdiag(i,j,ilwclr+l-1)
318     . + ( tend%lw%dtclr(i,j,l)
319     . + coup%lw%dlwdtg(i,j,L)*(coup%land%tgc(i,j)-coup%lw%tg0c(i,j)))
320     . * pkn(i,j,l)*pinv(i,j)*86400
321     enddo
322     enddo
323     endif
324    
325     c Solar Radiative Heating (deg/day)
326     c ---------------------------------
327     if (iradsw.ne.0) then
328     do j=jm1,jm2
329     do i=im1,im2
330     qdiag(i,j,iradsw+l-1) = qdiag(i,j,iradsw+l-1)
331     . + tend%sw%dt(i,j,l)*radswc(i,j)
332     . * pkn(i,j,l)*pinv(i,j)*86400
333     enddo
334     enddo
335     endif
336    
337     c Clear Sky Solar Radiative Heating (deg/day)
338     c -------------------------------------------
339     if (iswclr.ne.0) then
340     do j=jm1,jm2
341     do i=im1,im2
342     qdiag(i,j,iswclr+l-1) = qdiag(i,j,iswclr+l-1)
343     . + tend%sw%dtclr(i,j,l)*radswc(i,j)
344     . * pkn(i,j,l)*pinv(i,j)*86400
345     enddo
346     enddo
347     endif
348    
349     c Total U-Tendency (m/sec/day)
350     c ----------------------------
351     if( idudt.ne.0 ) then
352     do j=jm1,jm2
353     do i=im1,im2
354     qdiag(i,j,idudt+L-1) = qdiag(i,j,idudt+L-1)
355     . + ( updt%u(i,j,L)-curr%u(i,j,L) )*fact
356     enddo
357     enddo
358     endif
359    
360     c Total V-Tendency (m/sec/day)
361     c ----------------------------
362     if( idvdt.ne.0 ) then
363     do j=jm1,jm2
364     do i=im1,im2
365     qdiag(i,j,idvdt+L-1) = qdiag(i,j,idvdt+L-1)
366     . + ( updt%v(i,j,L)-curr%v(i,j,L) )*fact
367     enddo
368     enddo
369     endif
370    
371     c Total T-Tendency (deg/day)
372     c --------------------------
373     if( idtdt.ne.0 ) then
374     do j=jm1,jm2
375     do i=im1,im2
376     qdiag(i,j,idtdt+L-1) = qdiag(i,j,idtdt+L-1)
377     . + ( updt%t(i,j,L)*pknp1(i,j,L) - curr%t(i,j,L)*pkn(i,j,L) )*fact
378     enddo
379     enddo
380     endif
381    
382     c Total Q-Tendency (g/kg/day)
383     c ---------------------------
384     if( idqdt.ne.0 ) then
385     do j=jm1,jm2
386     do i=im1,im2
387     qdiag(i,j,idqdt+L-1) = qdiag(i,j,idqdt+L-1)
388     . + ( updt%q(i,j,L,1)-curr%q(i,j,L,1) )*fact*1000
389     enddo
390     enddo
391     endif
392    
393     c Averaged U-Field (m/sec)
394     c ------------------------
395     if( iuwnd.ne.0 ) then
396     do j=jm1,jm2
397     do i=im1,im2
398     qdiag(i,j,iuwnd+L-1) = qdiag(i,j,iuwnd+L-1) + curr%u(i,j,L)
399     enddo
400     enddo
401     endif
402    
403     c Averaged V-Field (m/sec)
404     c ------------------------
405     if( ivwnd.ne.0 ) then
406     do j=jm1,jm2
407     do i=im1,im2
408     qdiag(i,j,ivwnd+L-1) = qdiag(i,j,ivwnd+L-1) + curr%v(i,j,L)
409     enddo
410     enddo
411     endif
412    
413     c Averaged T-Field (deg)
414     c ----------------------
415     if( itmpu.ne.0 ) then
416     do j=jm1,jm2
417     do i=im1,im2
418     qdiag(i,j,itmpu+L-1) = qdiag(i,j,itmpu+L-1) + curr%t(i,j,L)*pkn(i,j,L)
419     enddo
420     enddo
421     endif
422    
423     c Averaged QQ-Field (m/sec)**2
424     c ----------------------------
425     if( itke.ne.0 ) then
426     do j=jm1,jm2
427     do i=im1,im2
428     qdiag(i,j,itke+L-1) = qdiag(i,j,itke+L-1) + coup%turb%qq(i,j,L)
429     enddo
430     enddo
431     endif
432    
433     c Averaged Q-Field (g/kg)
434     c -----------------------
435     if( isphu.ne.0 ) then
436     do j=jm1,jm2
437     do i=im1,im2
438     qdiag(i,j,isphu+L-1) = qdiag(i,j,isphu+L-1) + curr%q(i,j,L,1)*1000
439     enddo
440     enddo
441     endif
442    
443     enddo ! End Level Loop
444    
445     ndiabu = ndiabu + 1
446     ndiabv = ndiabv + 1
447     ndiabt = ndiabt + 1
448     nuiau = nuiau + 1
449     nviau = nviau + 1
450     ntiau = ntiau + 1
451     nradlw = nradlw + 1
452     nlwclr = nlwclr + 1
453     nradsw = nradsw + 1
454     nswclr = nswclr + 1
455     ndudt = ndudt + 1
456     ndvdt = ndvdt + 1
457     ndtdt = ndtdt + 1
458     ndiabq = ndiabq + 1
459     nqiau = nqiau + 1
460     ndqdt = ndqdt + 1
461     nuwnd = nuwnd + 1
462     nvwnd = nvwnd + 1
463     ntmpu = ntmpu + 1
464     ntke = ntke + 1
465     nsphu = nsphu + 1
466    
467     C **********************************************************************
468     C **** Compute Vertically Integrated Diagnostics ****
469     C **********************************************************************
470    
471     c Compute Moisture and Temperature Convergence Diagnostic
472     c -------------------------------------------------------
473     if( ivaveut.ne.0 .or. ivavevt.ne.0 .or.
474     . ivaveuq.ne.0 .or. ivavevq.ne.0 ) then
475     do j=jm1,jm2
476     do i=im1,im2
477     vintut(i,j) = 0.0
478     vintvt(i,j) = 0.0
479     vintuq(i,j) = 0.0
480     vintvq(i,j) = 0.0
481     enddo
482     enddo
483    
484     call ctoa ( curr%t,curr%t,dlam,dphi,im,jm,lm,0,grid%lattice )
485     call ctoa ( curr%q,curr%q,dlam,dphi,im,jm,lm,0,grid%lattice )
486     call ctoa_winds ( curr%u,curr%v,tmp1,tmp2,dlam,dphi,im,jm,lm,grid%lattice )
487     do L=1,Nrphys
488     do j=jm1,jm2
489     do i=im1,im2
490     vintut(i,j) = vintut(i,j) + tmp1(i,j,L)*curr%t(i,j,L) *dsig(L)* pkn(i,j,L)
491     vintvt(i,j) = vintvt(i,j) + tmp2(i,j,L)*curr%t(i,j,L) *dsig(L)* pkn(i,j,L)
492     vintuq(i,j) = vintuq(i,j) + tmp1(i,j,L)*curr%q(i,j,L,1)*dsig(L)
493     vintvq(i,j) = vintvq(i,j) + tmp2(i,j,L)*curr%q(i,j,L,1)*dsig(L)
494     enddo
495     enddo
496     enddo
497    
498     if( ivaveut.ne.0 ) then
499     do j=jm1,jm2
500     do i=im1,im2
501     qdiag(i,j,ivaveut) = qdiag(i,j,ivaveut) + vintut(i,j)
502     enddo
503     enddo
504     endif
505     if( ivavevt.ne.0 ) then
506     do j=jm1,jm2
507     do i=im1,im2
508     qdiag(i,j,ivavevt) = qdiag(i,j,ivavevt) + vintvt(i,j)
509     enddo
510     enddo
511     endif
512     if( ivaveuq.ne.0 ) then
513     do j=jm1,jm2
514     do i=im1,im2
515     qdiag(i,j,ivaveuq) = qdiag(i,j,ivaveuq) + vintuq(i,j)*1000
516     enddo
517     enddo
518     endif
519     if( ivavevq.ne.0 ) then
520     do j=jm1,jm2
521     do i=im1,im2
522     qdiag(i,j,ivavevq) = qdiag(i,j,ivavevq) + vintvq(i,j)*1000
523     enddo
524     enddo
525     endif
526    
527     endif ! End Convergence Diagnostic
528    
529     c Total Precipitable Water (gm/cm**2)
530     c -----------------------------------
531     if( itpw.ne.0 ) then
532     do j=jm1,jm2
533     do i=im1,im2
534     qbar(i,j) = 0.0
535     enddo
536     enddo
537     do L=1,Nrphys
538     do j=jm1,jm2
539     do i=im1,im2
540     qbar(i,j) = qbar(i,j) + curr%q(i,j,L,1)*dsig(L)
541     enddo
542     enddo
543     enddo
544     do j=jm1,jm2
545     do i=im1,im2
546     qdiag(i,j,itpw) = qdiag(i,j,itpw) + qbar(i,j)*10*curr%p(i,j)/grav
547     enddo
548     enddo
549     endif
550    
551     c Total Precipitable Analysis Increment (mm/day)
552     c ----------------------------------------------
553     if( ivaveqiau.ne.0 ) then
554     do j=jm1,jm2
555     do i=im1,im2
556     qbar(i,j) = 0.0
557     enddo
558     enddo
559     do L=1,Nrphys
560     do j=jm1,jm2
561     do i=im1,im2
562     qbar(i,j) = qbar(i,j) + dqiau(i,j,L)*dsig(L)
563     enddo
564     enddo
565     enddo
566     do j=jm1,jm2
567     do i=im1,im2
568     qdiag(i,j,ivaveqiau) = qdiag(i,j,ivaveqiau) + qbar(i,j)*(100*86400/grav)
569     enddo
570     enddo
571     endif
572    
573     c Vertically Averaged Analysis Temperature Increment (K/day)
574     c ----------------------------------------------------------
575     if( ivavetiau.ne.0 ) then
576     do j=jm1,jm2
577     do i=im1,im2
578     qbar(i,j) = 0.0
579     enddo
580     enddo
581     do L=1,Nrphys
582     do j=jm1,jm2
583     do i=im1,im2
584     qbar(i,j) = qbar(i,j) + dtiau(i,j,L)*pkn(i,j,l)*dsig(L)
585     enddo
586     enddo
587     enddo
588     do j=jm1,jm2
589     do i=im1,im2
590     qdiag(i,j,ivavetiau) = qdiag(i,j,ivavetiau) + qbar(i,j)*pinv(i,j)*86400
591     enddo
592     enddo
593     endif
594    
595     c Vertically Averaged Moist-T Increment (K/day)
596     c ---------------------------------------------
597     if( ivdtmoist.ne.0 ) then
598     do j=jm1,jm2
599     do i=im1,im2
600     qbar(i,j) = 0.0
601     enddo
602     enddo
603     do L=1,Nrphys
604     do j=jm1,jm2
605     do i=im1,im2
606     qbar(i,j) = qbar(i,j) + tend%moist%dt(i,j,L)*pkn(i,j,l)*dsig(L)
607     enddo
608     enddo
609     enddo
610     do j=jm1,jm2
611     do i=im1,im2
612     qdiag(i,j,ivdtmoist) = qdiag(i,j,ivdtmoist) + qbar(i,j)*pinv(i,j)*86400
613     enddo
614     enddo
615     endif
616    
617     c Vertically Averaged Turb-T Increment (K/day)
618     c --------------------------------------------
619     if( ivdtturb.ne.0 ) then
620     do j=jm1,jm2
621     do i=im1,im2
622     qbar(i,j) = 0.0
623     enddo
624     enddo
625     do L=1,Nrphys
626     do j=jm1,jm2
627     do i=im1,im2
628     qbar(i,j) = qbar(i,j) + tend%turb%dt(i,j,L)*pkn(i,j,l)*dsig(L)
629     enddo
630     enddo
631     enddo
632     do j=jm1,jm2
633     do i=im1,im2
634     qdiag(i,j,ivdtturb) = qdiag(i,j,ivdtturb) + qbar(i,j)*pinv(i,j)*86400
635     enddo
636     enddo
637     endif
638    
639     c Vertically Averaged RADLW Temperature Increment (K/day)
640     c -------------------------------------------------------
641     if( ivdtradlw.ne.0 ) then
642     do j=jm1,jm2
643     do i=im1,im2
644     qbar(i,j) = 0.0
645     enddo
646     enddo
647     do L=1,Nrphys
648     do j=jm1,jm2
649     do i=im1,im2
650     qbar(i,j) = qbar(i,j) + ( tend%lw%dt(i,j,L)
651     . + coup%lw%dlwdtg(i,j,L)*(coup%land%tgc(i,j)-coup%lw%tg0c(i,j)) )*pkn(i,j,l)*dsig(L)
652     enddo
653     enddo
654     enddo
655     do j=jm1,jm2
656     do i=im1,im2
657     qdiag(i,j,ivdtradlw) = qdiag(i,j,ivdtradlw) + qbar(i,j)*pinv(i,j)*86400
658     enddo
659     enddo
660     endif
661    
662     c Vertically Averaged RADSW Temperature Increment (K/day)
663     c -------------------------------------------------------
664     if( ivdtradsw.ne.0 ) then
665     do j=jm1,jm2
666     do i=im1,im2
667     qbar(i,j) = 0.0
668     enddo
669     enddo
670     do L=1,Nrphys
671     do j=jm1,jm2
672     do i=im1,im2
673     qbar(i,j) = qbar(i,j) + tend%sw%dt(i,j,L)*pkn(i,j,l)*dsig(L)
674     enddo
675     enddo
676     enddo
677     do j=jm1,jm2
678     do i=im1,im2
679     qdiag(i,j,ivdtradsw) = qdiag(i,j,ivdtradsw) + qbar(i,j)*radswc(i,j)*pinv(i,j)*86400
680     enddo
681     enddo
682     endif
683    
684     nvaveut = nvaveut + 1
685     nvavevt = nvavevt + 1
686     nvaveuq = nvaveuq + 1
687     nvavevq = nvavevq + 1
688     ntpw = ntpw + 1
689     nvaveqiau = nvaveqiau + 1
690     nvavetiau = nvavetiau + 1
691     nvdtmoist = nvdtmoist + 1
692     nvdtturb = nvdtturb + 1
693     nvdtradlw = nvdtradlw + 1
694     nvdtradsw = nvdtradsw + 1
695    
696     C *****************************************************************
697     C **** Release Workspace ****
698     C *****************************************************************
699    
700     deallocate ( dlam )
701     deallocate ( dphi )
702     deallocate ( dsig )
703     deallocate ( sig )
704     deallocate ( sige )
705    
706     deallocate ( pinv )
707     deallocate ( tmp1 )
708     deallocate ( tmp2 )
709     deallocate ( qbar )
710     deallocate ( vintuq )
711     deallocate ( vintvq )
712     deallocate ( vintut )
713     deallocate ( vintvt )
714     deallocate ( pkn )
715     deallocate ( pknp1 )
716     deallocate ( dtiau )
717     deallocate ( dqiau )
718    
719     call timeend (' step_diag')
720     return
721     end

  ViewVC Help
Powered by ViewVC 1.1.22