/[MITgcm]/MITgcm/pkg/monitor/mon_vort3.F
ViewVC logotype

Annotation of /MITgcm/pkg/monitor/mon_vort3.F

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


Revision 1.15 - (hide annotations) (download)
Sun Jun 28 01:05:41 2009 UTC (14 years, 11 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint64, checkpoint65, checkpoint62, checkpoint63, checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, checkpoint62c, checkpoint62b, checkpoint62a, checkpoint62g, checkpoint62f, checkpoint62e, checkpoint62d, checkpoint62k, checkpoint62j, checkpoint62i, checkpoint62h, checkpoint62o, checkpoint62n, checkpoint62m, checkpoint62l, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62p, checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62z, checkpoint62y, checkpoint62x, checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61s, checkpoint61z, checkpoint61x, checkpoint61y, HEAD
Changes since 1.14: +2 -2 lines
add bj in exch2 arrays and S/R

1 jmc 1.15 C $Header: /u/gcmpack/MITgcm/pkg/monitor/mon_vort3.F,v 1.14 2009/05/12 19:56:36 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4 adcroft 1.2 #include "MONITOR_OPTIONS.h"
5 jmc 1.1
6 edhill 1.8 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7     CBOP
8     C !ROUTINE: MON_VORT3
9    
10     C !INTERFACE:
11 jmc 1.1 SUBROUTINE MON_VORT3(
12 edhill 1.8 I myIter, myThid )
13    
14     C !DESCRIPTION:
15     C Calculates stats for Vorticity (z-component).
16    
17     C !USES:
18 jmc 1.1 IMPLICIT NONE
19     #include "SIZE.h"
20     #include "EEPARAMS.h"
21     #include "PARAMS.h"
22     #include "DYNVARS.h"
23     #include "MONITOR.h"
24     #include "GRID.h"
25 adcroft 1.6 #ifdef ALLOW_EXCH2
26 jmc 1.14 #include "W2_EXCH2_SIZE.h"
27 afe 1.3 #include "W2_EXCH2_TOPOLOGY.h"
28 adcroft 1.6 #endif /* ALLOW_EXCH2 */
29 jmc 1.1
30 edhill 1.8 C !INPUT PARAMETERS:
31 jmc 1.1 INTEGER myIter, myThid
32 edhill 1.8 CEOP
33 jmc 1.1
34 edhill 1.8 C !LOCAL VARIABLES:
35 jmc 1.1 INTEGER bi,bj,i,j,k
36     INTEGER iMax,jMax
37 jmc 1.9 _RL theVol, theArea, tmpVal, tmpAre, tmpVol
38 jmc 1.1 _RL theMin, theMax, theMean, theVar, volMean, volVar, theSD
39 jmc 1.12 c _RL areaTile, volTile, sumTile, sqsTile, vSumTile, vSqsTile
40     _RL tileArea(nSx,nSy)
41     _RL tileVol (nSx,nSy)
42     _RL tileSum (nSx,nSy)
43     _RL tileVar (nSx,nSy)
44     _RL tileVSum(nSx,nSy)
45     _RL tileVSq (nSx,nSy)
46 jmc 1.1 _RL vort3(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
47     _RS hFacZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
48     _RL AZcorner
49     #ifdef MONITOR_TEST_HFACZ
50     _RL tmpFac
51     _RL etaFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
52     #endif
53 jmc 1.12 LOGICAL northWestCorner, northEastCorner
54     LOGICAL southWestCorner, southEastCorner
55 jmc 1.11 INTEGER iG
56     #ifdef ALLOW_EXCH2
57     INTEGER myTile
58     #endif
59 jmc 1.1
60     theMin = 1. _d 20
61     theMax =-1. _d 20
62     theArea= 0. _d 0
63     theMean= 0. _d 0
64     theVar = 0. _d 0
65     theVol = 0. _d 0
66     volMean= 0. _d 0
67     volVar = 0. _d 0
68     theSD = 0. _d 0
69     AZcorner = 1. _d 0
70    
71     DO bj=myByLo(myThid),myByHi(myThid)
72     DO bi=myBxLo(myThid),myBxHi(myThid)
73 jmc 1.12 tileArea(bi,bj)= 0. _d 0
74     tileVol(bi,bj) = 0. _d 0
75     tileSum(bi,bj) = 0. _d 0
76     tileVar(bi,bj) = 0. _d 0
77     tileVSum(bi,bj)= 0. _d 0
78     tileVSq(bi,bj) = 0. _d 0
79 jmc 1.1 #ifdef MONITOR_TEST_HFACZ
80     tmpFac = 0.
81 jmc 1.12 IF( implicDiv2Dflow.GT.0 .AND. abEps.GT.-0.5 )
82 jmc 1.1 & tmpFac = (0.5+abEps) / implicDiv2Dflow
83     DO j=1-Oly,sNy+Oly
84     DO i=1-Olx,sNx+Olx
85     etaFld(i,j) = etaH(i,j,bi,bj)
86     & + tmpFac*(etaN(i,j,bi,bj)-etaH(i,j,bi,bj))
87     ENDDO
88     ENDDO
89     #endif
90     DO k=1,Nr
91    
92     iMax = sNx
93     jMax = sNy
94     DO j=1,sNy
95     DO i=1,sNx
96     #ifdef MONITOR_TEST_HFACZ
97     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
98     C- Test various definitions of hFacZ (for 1 layer, flat bottom ocean):
99     c hFacZ(i,j) = 1. +
100 jmc 1.12 c & 0.25 _d 0*( etaFld(i-1,j-1)
101 jmc 1.1 c & + etaFld( i ,j-1)
102 jmc 1.12 c & + etaFld(i-1, j )
103     c & + etaFld( i , j )
104 jmc 1.1 c & )*recip_drF(k)
105     c hFacZ(i,j) = 1. +
106     c & 0.25 _d 0*( etaFld(i-1,j-1)*rA(i-1,j-1,bi,bj)
107 jmc 1.12 c & + etaFld( i ,j-1)*rA( i ,j-1,bi,bj)
108 jmc 1.1 c & + etaFld(i-1, j )*rA(i-1, j ,bi,bj)
109 jmc 1.12 c & + etaFld( i , j )*rA( i , j ,bi,bj)
110 jmc 1.1 c & )*recip_drF(k)*recip_rAz(i,j,bi,bj)
111     hFacZ(i,j) = 1. + 0.125 _d 0*
112     & ( ( etaFld(i-1,j-1)*rA(i-1,j-1,bi,bj)
113     & +etaFld( i ,j-1)*rA( i ,j-1,bi,bj)
114 jmc 1.12 & )*recip_rAw(i,j-1,bi,bj)
115 jmc 1.1 & + ( etaFld(i-1, j )*rA(i-1, j ,bi,bj)
116 jmc 1.12 & +etaFld( i , j )*rA( i , j ,bi,bj)
117 jmc 1.1 & )*recip_rAw(i, j ,bi,bj)
118     & + ( etaFld(i-1,j-1)*rA(i-1,j-1,bi,bj)
119     & +etaFld(i-1, j )*rA(i-1, j ,bi,bj)
120 jmc 1.12 & )*recip_rAs(i-1,j,bi,bj)
121 jmc 1.1 & + ( etaFld( i ,j-1)*rA( i ,j-1,bi,bj)
122 jmc 1.12 & + etaFld( i , j )*rA( i , j ,bi,bj)
123 jmc 1.1 & )*recip_rAs( i ,j,bi,bj)
124     & )*recip_drF(k)
125     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
126     #else
127     C- Standard definition of hFac at vorticity point:
128     hFacZ(i,j) =
129 jmc 1.12 & 0.25 _d 0*( _hFacW(i,j-1,k,bi,bj)
130     & + _hFacW(i, j ,k,bi,bj)
131     & + _hFacS(i-1,j,k,bi,bj)
132 jmc 1.1 & + _hFacS( i ,j,k,bi,bj)
133     & )
134 jmc 1.12 #endif /* MONITOR_TEST_HFACZ */
135 jmc 1.1 vort3(i,j) = recip_rAz(i,j,bi,bj)*(
136     & vVel( i ,j,k,bi,bj)*dyC( i ,j,bi,bj)
137     & -vVel(i-1,j,k,bi,bj)*dyC(i-1,j,bi,bj)
138     & -uVel(i, j ,k,bi,bj)*dxC(i, j ,bi,bj)
139     & +uVel(i,j-1,k,bi,bj)*dxC(i,j-1,bi,bj)
140 jmc 1.12 & )
141 jmc 1.1 ENDDO
142     ENDDO
143    
144     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
145     C Special stuff for Cubed Sphere:
146     IF (useCubedSphereExchange) THEN
147     c AZcorner = 0.75 _d 0
148     iMax = sNx+1
149     jMax = sNy+1
150     DO i=1,iMax
151     hFacZ(i,jMax)=0.
152     vort3(i,jMax)=0.
153     ENDDO
154     DO j=1,jMax
155     hFacZ(iMax,j)=0.
156     vort3(iMax,j)=0.
157     ENDDO
158 jmc 1.5
159     southWestCorner = .TRUE.
160     southEastCorner = .TRUE.
161     northWestCorner = .TRUE.
162     northEastCorner = .TRUE.
163     iG = bi+(myXGlobalLo-1)/sNx
164 adcroft 1.6 #ifdef ALLOW_EXCH2
165 jmc 1.15 myTile = W2_myTileList(bi,bj)
166 jmc 1.5 iG = exch2_myFace(myTile)
167 jmc 1.12 southWestCorner = exch2_isWedge(myTile).EQ.1
168 jmc 1.5 & .AND. exch2_isSedge(myTile).EQ.1
169     southEastCorner = exch2_isEedge(myTile).EQ.1
170     & .AND. exch2_isSedge(myTile).EQ.1
171     northEastCorner = exch2_isEedge(myTile).EQ.1
172     & .AND. exch2_isNedge(myTile).EQ.1
173     northWestCorner = exch2_isWedge(myTile).EQ.1
174     & .AND. exch2_isNedge(myTile).EQ.1
175 adcroft 1.6 #endif /* ALLOW_EXCH2 */
176 jmc 1.5
177     C-- avoid to count 3 times the same corner:
178     southEastCorner = southEastCorner .AND. iG.EQ.2
179     northWestCorner = northWestCorner .AND. iG.EQ.1
180     northEastCorner = .FALSE.
181    
182 afe 1.3 C-- S.W. corner:
183     IF ( southWestCorner ) THEN
184 jmc 1.1 i=1
185     j=1
186     vort3(i,j)=
187     & +recip_rAz(i,j,bi,bj)/AZcorner*(
188     & vVel(i,j,k,bi,bj)*dyC(i,j,bi,bj)
189     & -uVel(i,j,k,bi,bj)*dxC(i,j,bi,bj)
190     & +uVel(i,j-1,k,bi,bj)*dxC(i,j-1,bi,bj)
191     & )
192 jmc 1.12 hFacZ(i,j) = ( _hFacW(i,j-1,k,bi,bj)
193     & + _hFacW(i, j ,k,bi,bj)
194 jmc 1.1 & + _hFacS( i ,j,k,bi,bj)
195     & )/3. _d 0
196 afe 1.3 ENDIF
197     IF ( southEastCorner ) THEN
198 jmc 1.1 C-- S.E. corner:
199     i=iMax
200     j=1
201 jmc 1.12 vort3(i,j)=
202     & +recip_rAz(i,j,bi,bj)/AZcorner*(
203 jmc 1.1 & -vVel(i-1,j,k,bi,bj)*dyC(i-1,j,bi,bj)
204     & -uVel(i,j,k,bi,bj)*dxC(i,j,bi,bj)
205     & +uVel(i,j-1,k,bi,bj)*dxC(i,j-1,bi,bj)
206     & )
207 jmc 1.12 hFacZ(i,j) = ( _hFacW(i,j-1,k,bi,bj)
208     & + _hFacW(i, j ,k,bi,bj)
209     & + _hFacS(i-1,j,k,bi,bj)
210 jmc 1.1 & )/3. _d 0
211 afe 1.3 ENDIF
212     IF ( northWestCorner ) THEN
213 jmc 1.1 C-- N.W. corner:
214     i=1
215     j=jMax
216     vort3(i,j)=
217     & +recip_rAz(i,j,bi,bj)/AZcorner*(
218     & vVel(i,j,k,bi,bj)*dyC(i,j,bi,bj)
219     & -uVel(i,j,k,bi,bj)*dxC(i,j,bi,bj)
220     & +uVel(i,j-1,k,bi,bj)*dxC(i,j-1,bi,bj)
221     & )
222 jmc 1.12 hFacZ(i,j) = ( _hFacW(i,j-1,k,bi,bj)
223     & + _hFacW(i, j ,k,bi,bj)
224 jmc 1.1 & + _hFacS( i ,j,k,bi,bj)
225     & )/3. _d 0
226 afe 1.3 ENDIF
227     IF ( northEastCorner ) THEN
228 jmc 1.1 C-- N.E. corner:
229     i=iMax
230     j=jMax
231     vort3(i,j)=
232     & +recip_rAz(i,j,bi,bj)/AZcorner*(
233     & -vVel(i-1,j,k,bi,bj)*dyC(i-1,j,bi,bj)
234     & -uVel(i,j,k,bi,bj)*dxC(i,j,bi,bj)
235     & +uVel(i,j-1,k,bi,bj)*dxC(i,j-1,bi,bj)
236     & )
237 jmc 1.12 hFacZ(i,j) = ( _hFacW(i,j-1,k,bi,bj)
238     & + _hFacW(i, j ,k,bi,bj)
239     & + _hFacS(i-1,j,k,bi,bj)
240 jmc 1.1 & )/3. _d 0
241 afe 1.3 ENDIF
242 jmc 1.1 ENDIF
243    
244     C- Special stuff for North & South Poles, LatLon grid
245     IF ( usingSphericalPolarGrid ) THEN
246     IF (yG(1,sNy+1,bi,bj).EQ.90.) THEN
247     jMax = sNy+1
248     j = jMax
249     DO i=1,sNx
250     vort3(i,j) = 0.
251     vort3(1,j) = vort3(1,j)
252     & + uVel(i,j-1,k,bi,bj)*dxC(i,j-1,bi,bj)
253     hFacZ(i,j) = 0.
254     #ifndef MONITOR_TEST_HFACZ
255 jmc 1.12 hFacZ(1,j) = hFacZ(1,j) + _hFacW(i,j-1,k,bi,bj)
256 jmc 1.1 ENDDO
257     #else
258 jmc 1.12 hFacZ(1,j) = hFacZ(1,j) + etaFld(i,j-1)
259 jmc 1.1 ENDDO
260     hFacZ(1,j) = sNx + hFacZ(1,j)*recip_drF(k)
261     #endif
262     hFacZ(1,j) = hFacZ(1,j) / FLOAT(sNx)
263     vort3(1,j) = vort3(1,j)*recip_rAz(1,j,bi,bj)
264     ENDIF
265     IF (yG(1,1,bi,bj).EQ.-90.) THEN
266     j = 1
267     DO i=1,sNx
268     vort3(i,j) = 0.
269     vort3(1,j) = vort3(1,j)
270     & - uVel(i,j,k,bi,bj)*dxC(i,j,bi,bj)
271     hFacZ(i,j) = 0.
272     #ifndef MONITOR_TEST_HFACZ
273 jmc 1.12 hFacZ(1,j) = hFacZ(1,j) + _hFacW(i,j,k,bi,bj)
274 jmc 1.1 ENDDO
275     #else
276 jmc 1.12 hFacZ(1,j) = hFacZ(1,j) + etaFld(i,j)
277 jmc 1.1 ENDDO
278     hFacZ(1,j) = sNx + hFacZ(1,j)*recip_drF(k)
279     #endif
280     hFacZ(1,j) = hFacZ(1,j) / FLOAT(sNx)
281     vort3(1,j) = vort3(1,j)*recip_rAz(1,j,bi,bj)
282     ENDIF
283     ENDIF
284    
285     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
286    
287 jmc 1.12 DO j=1,jMax
288     DO i=1,iMax
289 jmc 1.1 IF (hFacZ(i,j).GT.0. _d 0) THEN
290     tmpVal = vort3(i,j)
291     tmpAre = rAz(i,j,bi,bj)*drF(k)
292     tmpVol = rAz(i,j,bi,bj)*drF(k)*hFacZ(i,j)
293 jmc 1.12 tileArea(bi,bj) = tileArea(bi,bj) + tmpAre
294 jmc 1.1 C- min,max of relative vorticity ("r")
295     theMin = MIN(theMin,tmpVal)
296     theMax = MAX(theMax,tmpVal)
297     C- average & std.dev of absolute vorticity ("a")
298     tmpVal = tmpVal + fCoriG(i,j,bi,bj)
299 jmc 1.12 tileSum(bi,bj) = tileSum(bi,bj) + tmpAre*tmpVal
300     tileVar(bi,bj) = tileVar(bi,bj) + tmpAre*tmpVal*tmpVal
301 jmc 1.1 C- average & std.dev of potential vorticity ("p")
302     tmpVal = tmpVal / hFacZ(i,j)
303 jmc 1.12 tileVol(bi,bj) = tileVol(bi,bj) + tmpVol
304     tileVSum(bi,bj)= tileVSum(bi,bj)+ tmpVol*tmpVal
305     tileVSq(bi,bj) = tileVSq(bi,bj) + tmpVol*tmpVal*tmpVal
306 jmc 1.1 ENDIF
307     ENDDO
308     ENDDO
309    
310     ENDDO
311 jmc 1.12 c theArea= theArea + tileArea(bi,bj)
312     c theVol = theVol + tileVol (bi,bj)
313     c theMean= theMean + tileSum(bi,bj)
314     c theVar = theVar + tileVar(bi,bj)
315     c volMean= volMean + tileVSum(bi,bj)
316     c volVar = volVar + tileVSq(bi,bj)
317 jmc 1.1 ENDDO
318     ENDDO
319    
320     theMin = -theMin
321 jmc 1.13 _GLOBAL_MAX_RL(theMin, myThid)
322     _GLOBAL_MAX_RL(theMax, myThid)
323 jmc 1.1 theMin = -theMin
324 jmc 1.13 c _GLOBAL_SUM_RL(theArea,myThid)
325     c _GLOBAL_SUM_RL(theVol, myThid)
326     c _GLOBAL_SUM_RL(theMean,myThid)
327     c _GLOBAL_SUM_RL(theVar, myThid)
328     c _GLOBAL_SUM_RL(volMean,myThid)
329     c _GLOBAL_SUM_RL(volVar ,myThid)
330 jmc 1.12 CALL GLOBAL_SUM_TILE_RL( tileArea, theArea, myThid )
331     CALL GLOBAL_SUM_TILE_RL( tileVol, theVol, myThid )
332     CALL GLOBAL_SUM_TILE_RL( tileSum, theMean, myThid )
333     CALL GLOBAL_SUM_TILE_RL( tileVar, theVar, myThid )
334     CALL GLOBAL_SUM_TILE_RL( tileVSum, volMean, myThid )
335     CALL GLOBAL_SUM_TILE_RL( tileVSq, volVar, myThid )
336     IF (theArea.GT.0.) THEN
337 jmc 1.1 theMean= theMean/theArea
338     theVar = theVar /theArea
339     theVar = theVar - theMean*theMean
340     c IF (theVar.GT.0.) theSD = SQRT(theVar)
341     IF (theVar.GT.0.) theVar = SQRT(theVar)
342     ENDIF
343     IF (theVol.GT.0.) THEN
344     volMean= volMean/theVol
345     volVar = volVar /theVol
346     volVar = volVar - volMean*volMean
347     IF (volVar.GT.0.) theSD = SQRT(volVar)
348     ENDIF
349    
350     C- Print stats for (relative/absolute) Vorticity AND Pot.Vort.
351     CALL MON_SET_PREF('vort',myThid)
352     CALL MON_OUT_RL(mon_string_none,theMin, '_r_min', myThid)
353     CALL MON_OUT_RL(mon_string_none,theMax, '_r_max', myThid)
354     CALL MON_OUT_RL(mon_string_none,theMean,'_a_mean', myThid)
355     CALL MON_OUT_RL(mon_string_none,theVar, '_a_sd', myThid)
356     CALL MON_OUT_RL(mon_string_none,volMean,'_p_mean', myThid)
357     CALL MON_OUT_RL(mon_string_none,theSD, '_p_sd', myThid)
358     c CALL MON_OUT_RL(mon_string_none,theVol,mon_foot_vol,myThid)
359    
360     RETURN
361     END

  ViewVC Help
Powered by ViewVC 1.1.22