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

Diff of /MITgcm/pkg/monitor/mon_ke.F

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

revision 1.23 by jmc, Wed Dec 26 23:57:46 2012 UTC revision 1.24 by jmc, Sun Feb 17 04:07:30 2013 UTC
# Line 34  C     !LOCAL VARIABLES: Line 34  C     !LOCAL VARIABLES:
34        INTEGER i,j,k        INTEGER i,j,k
35        INTEGER ks, kp1        INTEGER ks, kp1
36        _RL numPnts,theVol,tmpVal, mskp1, msk_1        _RL numPnts,theVol,tmpVal, mskp1, msk_1
37        _RL weight0, weight1        _RL abFac1, abFac2, R_drK, cosLat
38        _RL theMax,theMean,theVolMean,potEnMean        _RL theMax,theMean,theVolMean,potEnMean
39        _RL uBarC, vBarC, totAMu, totAMs        _RL totAMu, totAMs
40        _RL tileMean(nSx,nSy)        _RL tileMean(nSx,nSy)
41        _RL tileVlAv(nSx,nSy)        _RL tileVlAv(nSx,nSy)
42        _RL tilePEav(nSx,nSy)        _RL tilePEav(nSx,nSy)
43        _RL tileVol (nSx,nSy)        _RL tileVol (nSx,nSy)
44        _RL tileAMu (nSx,nSy)        _RL tileAMu (nSx,nSy)
45        _RL tileAMs (nSx,nSy)        _RL tileAMs (nSx,nSy)
46        _RL radDist(1:sNx,1:sNy)        _RL tmpFld(1:sNx,1:sNy)
47          _RS cos2LatG(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
48  #ifdef ALLOW_NONHYDROSTATIC  #ifdef ALLOW_NONHYDROSTATIC
49        _RL tmpWke        _RL tmpWke
50  #endif  #endif
51    #ifdef ALLOW_ADAMSBASHFORTH_3
52          INTEGER m1, m2
53    #endif
54    
55  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
56    
# Line 165  C--   Compute total angular momentum Line 169  C--   Compute total angular momentum
169        IF ( mon_output_AM ) THEN        IF ( mon_output_AM ) THEN
170         DO bj=myByLo(myThid),myByHi(myThid)         DO bj=myByLo(myThid),myByHi(myThid)
171          DO bi=myBxLo(myThid),myBxHi(myThid)          DO bi=myBxLo(myThid),myBxHi(myThid)
172  C-    calculate radial distance  C-    Calculate contribution from zonal velocity
173             abFac1 = 0. _d 0
174             abFac2 = 0. _d 0
175    #ifdef ALLOW_ADAMSBASHFORTH_3
176              m1 = 1 + mod(myIter+1,2)
177              m2 = 1 + mod( myIter ,2)
178              IF ( myIter.GE.2 ) abFac2 = beta_AB
179              IF ( myIter.GE.1 ) abFac1 = -( alph_AB + abFac2 )
180    #else
181              IF ( myIter.GE.1 ) abFac1 = -( 0.5 _d 0 + abEps )
182    #endif
183    C-    contribution from uVel component: 1rst integrate vertically
184           DO j=1,sNy           DO j=1,sNy
185            DO i=1,sNx            DO i=1,sNx
186              radDist(i,j) = rSphere*COS( deg2rad*yC(i,j,bi,bj) )              tmpFld(i,j) = 0. _d 0
      &                            *maskInC(i,j,bi,bj)  
187            ENDDO            ENDDO
188           ENDDO           ENDDO
 C-    calculate contribution from zonal velocity  
          tileAMu(bi,bj) = 0. _d 0  
          tileAMs(bi,bj) = 0. _d 0  
189           DO k=1,Nr           DO k=1,Nr
190              R_drK = rSphere*deepFacC(k)*deepFac2C(k)
191         &                   *rhoFacC(k)*drF(k)
192            DO j=1,sNy            DO j=1,sNy
193             DO i=1,sNx             DO i=1,sNx
194              uBarC = (uVel(i,j,k,bi,bj)+uVel(i+1,j,k,bi,bj))*0.5 _d 0  #ifdef ALLOW_ADAMSBASHFORTH_3
195              vBarC = (vVel(i,j,k,bi,bj)+vVel(i,j+1,k,bi,bj))*0.5 _d 0              tmpVal = abFac1*guNm(i,j,k,bi,bj,m1)
196              tmpVal = ( angleCosC(i,j,bi,bj)*uBarC       &             + abFac2*guNm(i,j,k,bi,bj,m2)
197       &                -angleSinC(i,j,bi,bj)*vBarC  #else
198       &               )*radDist(i,j)*deepFacC(k)              tmpVal = abFac1*guNm1(i,j,k,bi,bj)
199              tileAMu(bi,bj) = tileAMu(bi,bj)  #endif
200       &             + tmpVal*rA(i,j,bi,bj)*deepFac2C(k)              tmpVal = tmpVal*deltaTMom + uVel(i,j,k,bi,bj)
201       &                     *rhoFacC(k)*drF(k)*_hFacC(i,j,k,bi,bj)              tmpFld(i,j) = tmpFld(i,j)
202         &             + R_drK*tmpVal*_hFacW(i,j,k,bi,bj)
203             ENDDO             ENDDO
204            ENDDO            ENDDO
205           ENDDO           ENDDO
206  C-    and contribution from mass distribution anomaly (i.e., free-surface)  C-    and then integrate horizontally over this tile
207  c        IF ( .FALSE. ) THEN           DO j=1,sNy
208           IF ( exactConserv ) THEN            DO i=1,sNx
209  C-    To improve balance between zonal-wind (AMu) and mass (AMs) AM, need              cosLat = COS( deg2rad*
210  C     a consistent time-stepping of meridional mass transport in Coriolis       &             ( yG(i,j,bi,bj) + yG(i,j+1,bi,bj) )*halfRL )
211  C     (zonal momentum, go through AB) and mass transport divergence;              tmpFld(i,j) = tmpFld(i,j)*u2zonDir(i,j,bi,bj)
212  C     try to account for AB-2 in AMs the same way it applies to f*v:       &                   *cosLat*rAw(i,j,bi,bj)
213         &                   *maskInW(i,j,bi,bj)
214              ENDDO
215             ENDDO
216             tileAMu(bi,bj) = 0. _d 0
217             DO j=1,sNy
218              DO i=1,sNx
219                tileAMu(bi,bj) = tileAMu(bi,bj) + tmpFld(i,j)
220              ENDDO
221             ENDDO
222    C-    contribution from vVel component: 1rst integrate vertically
223             DO j=1,sNy
224              DO i=1,sNx
225                tmpFld(i,j) = 0. _d 0
226              ENDDO
227             ENDDO
228             DO k=1,Nr
229              R_drK = rSphere*deepFacC(k)*deepFac2C(k)
230         &                   *rhoFacC(k)*drF(k)
231              DO j=1,sNy
232               DO i=1,sNx
233  #ifdef ALLOW_ADAMSBASHFORTH_3  #ifdef ALLOW_ADAMSBASHFORTH_3
234            weight1 = alph_AB              tmpVal = abFac1*gvNm(i,j,k,bi,bj,m1)
235         &             + abFac2*gvNm(i,j,k,bi,bj,m2)
236  #else  #else
237            weight1 = 0.5 _d 0 + abEps              tmpVal = abFac1*gvNm1(i,j,k,bi,bj)
238  #endif  #endif
239            weight0 = 1.0 _d 0 - weight1              tmpVal = tmpVal*deltaTMom + vVel(i,j,k,bi,bj)
240                tmpFld(i,j) = tmpFld(i,j)
241         &             + R_drK*tmpVal*_hFacS(i,j,k,bi,bj)
242               ENDDO
243              ENDDO
244             ENDDO
245    C-    and then integrate horizontally over this tile
246             DO j=1,sNy
247              DO i=1,sNx
248                cosLat = COS( deg2rad*
249         &             ( yG(i,j,bi,bj) + yG(i+1,j,bi,bj) )*halfRL )
250                tmpFld(i,j) = tmpFld(i,j)*v2zonDir(i,j,bi,bj)
251         &                   *cosLat*rAs(i,j,bi,bj)
252         &                   *maskInS(i,j,bi,bj)
253              ENDDO
254             ENDDO
255             DO j=1,sNy
256              DO i=1,sNx
257                tileAMu(bi,bj) = tileAMu(bi,bj) + tmpFld(i,j)
258              ENDDO
259             ENDDO
260    C-    Calculate contribution from mass distribution anomaly (i.e., free-surface)
261             IF ( exactConserv ) THEN
262            DO j=1,sNy            DO j=1,sNy
263             DO i=1,sNx             DO i=1,sNx
             ks = kSurfC(i,j,bi,bj)  
             tmpVal = weight1*etaH(i,j,bi,bj)  
264  #ifdef EXACT_CONSERV  #ifdef EXACT_CONSERV
265       &             + weight0*etaHnm1(i,j,bi,bj)              tmpFld(i,j) = etaHnm1(i,j,bi,bj)
266    #else
267                tmpFld(i,j) = 0.
268  #endif  #endif
             tmpVal = omega*tmpVal  
      &             * radDist(i,j)*radDist(i,j)*deepFac2F(ks)  
             tileAMs(bi,bj) = tileAMs(bi,bj)  
      &             + tmpVal*rA(i,j,bi,bj)*deepFac2F(ks)*rhoFacF(ks)  
269             ENDDO             ENDDO
270            ENDDO            ENDDO
271           ELSE           ELSE
272            DO j=1,sNy            DO j=1,sNy
273             DO i=1,sNx             DO i=1,sNx
274              ks = kSurfC(i,j,bi,bj)              tmpFld(i,j) = etaN(i,j,bi,bj)
             tmpVal = omega*etaN(i,j,bi,bj)  
      &             * radDist(i,j)*radDist(i,j)*deepFac2F(ks)  
             tileAMs(bi,bj) = tileAMs(bi,bj)  
      &             + tmpVal*rA(i,j,bi,bj)*deepFac2F(ks)*rhoFacF(ks)  
275             ENDDO             ENDDO
276            ENDDO            ENDDO
277           ENDIF           ENDIF
278    C-    calculate angular momentum from mass-distribution anomaly
279    C     using square of radial distance (averaged @ center point)
280             DO j=1-OLy,sNy+OLy
281              DO i=1-OLx,sNx+OLx
282                cosLat = COS( deg2rad*yG(i,j,bi,bj) )
283                cos2LatG(i,j) = cosLat*cosLat
284              ENDDO
285             ENDDO
286             DO j=1,sNy
287              DO i=1,sNx
288                tmpFld(i,j) = tmpFld(i,j)
289         &        *omega*rSphere*rSphere
290         &        *( ( cos2LatG(i,j) + cos2LatG(i+1,j+1) )
291         &         + ( cos2LatG(i+1,j) + cos2LatG(i,j+1) )
292         &         )*0.25 _d 0
293              ENDDO
294             ENDDO
295             DO j=1,sNy
296              DO i=1,sNx
297                ks = kSurfC(i,j,bi,bj)
298                tmpFld(i,j) = tmpFld(i,j)
299         &             *maskInC(i,j,bi,bj)*deepFac2F(ks)
300         &             *rA(i,j,bi,bj)*deepFac2F(ks)*rhoFacF(ks)
301              ENDDO
302             ENDDO
303             tileAMs(bi,bj) = 0. _d 0
304             DO j=1,sNy
305              DO i=1,sNx
306                tileAMs(bi,bj) = tileAMs(bi,bj) + tmpFld(i,j)
307              ENDDO
308             ENDDO
309  C- end bi,bj loops  C- end bi,bj loops
310          ENDDO          ENDDO
311         ENDDO         ENDDO
# Line 242  C--   Print stats for total Angular Mome Line 322  C--   Print stats for total Angular Mome
322       &         '_eta_mean', myThid )       &         '_eta_mean', myThid )
323         CALL MON_OUT_RL( mon_string_none, totAMu,         CALL MON_OUT_RL( mon_string_none, totAMu,
324       &         '_uZo_mean', myThid )       &         '_uZo_mean', myThid )
325         totAMu = totAMu + totAMs         totAMu = totAMu + freeSurfFac*totAMs
326         CALL MON_OUT_RL( mon_string_none, totAMu,         CALL MON_OUT_RL( mon_string_none, totAMu,
327       &         '_tot_mean', myThid )       &         '_tot_mean', myThid )
328    
Legend:
Removed from v.1.23  
changed lines
  Added in v.1.24
  ViewVC Help
Powered by ViewVC 1.1.22