/[MITgcm]/MITgcm/pkg/pp81/pp81_ri_number.F
ViewVC logotype

Diff of /MITgcm/pkg/pp81/pp81_ri_number.F

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

revision 1.2 by jmc, Thu Aug 23 19:13:10 2007 UTC revision 1.3 by jmc, Mon Aug 11 22:28:06 2008 UTC
# Line 8  C !ROUTINE: PP81_RI_NUMBER Line 8  C !ROUTINE: PP81_RI_NUMBER
8    
9  C !INTERFACE: ==========================================================  C !INTERFACE: ==========================================================
10        subroutine PP81_RI_NUMBER(        subroutine PP81_RI_NUMBER(
11       I     bi, bj, K, iMin, iMax, jMin, jMax,       I     bi, bj, K, iMin, iMax, jMin, jMax,
12       O     RiNumber,       O     RiNumber,
13       I     myTime, myThid )       I     myTime, myThid )
14    
# Line 41  C !USES: =============================== Line 41  C !USES: ===============================
41  C !INPUT PARAMETERS: ===================================================  C !INPUT PARAMETERS: ===================================================
42  C Routine arguments  C Routine arguments
43  C     bi, bj - array indices on which to apply calculations  C     bi, bj - array indices on which to apply calculations
44  C     iMin, iMax, jMin, jMax  C     iMin, iMax, jMin, jMax
45  C            - array boundaries  C            - array boundaries
46  C     k      - depth level  C     k      - depth level
47  C     myTime - Current time in simulation  C     myTime - Current time in simulation
# Line 59  C     p0-125       - averaging coefficie Line 59  C     p0-125       - averaging coefficie
59  C     tempu, tempv - temporary variables  C     tempu, tempv - temporary variables
60  C     rhoK, rhoKm1 - Density below and above current interface  C     rhoK, rhoKm1 - Density below and above current interface
61  C     epsilon      - small number  C     epsilon      - small number
62  C     RiFlux       - denominator of Richardson number  C     RiFlux       - denominator of Richardson number
63  C     BuoyFreq     - buoyancy frequency  C     BuoyFreq     - buoyancy frequency
64        INTEGER I,J,Km1        INTEGER I,J,Km1
65        _RL        p5    , p125              _RL        p5    , p125
66        PARAMETER( p5=0.5, p125=0.125 )        PARAMETER( p5=0.5, p125=0.125 )
67        _RL tempu, tempv        _RL tempu, tempv
68        _RL epsilon        _RL epsilon
69        PARAMETER    (  epsilon = 1.D-10 )        PARAMETER    (  epsilon = 1.D-10 )
70        _RL RiFlux          _RL RiFlux
71        _RL buoyFreq        _RL buoyFreq
72        _RL rhoKm1   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhoKm1   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
73        _RL rhoK     (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rhoK     (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
74  #ifdef PP81_SMOOTH_RI  #ifdef PP81_SMOOTH_RI
75        _RL RiTmp   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL RiTmp   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
76  #endif /* PP81_SMOOTH_RI */  #endif /* PP81_SMOOTH_RI */
77  CEOP      CEOP
78    
79        Km1   = MAX(1,K-1)        Km1   = MAX(1,K-1)
80  C     Not sure what is correct for pressure coordinates:  C     Not sure what is correct for pressure coordinates:
# Line 92  CML      ELSE Line 92  CML      ELSE
92  CML       STOP 'PP81_RI_NUMBER: We should never reach this point'  CML       STOP 'PP81_RI_NUMBER: We should never reach this point'
93  CML      ENDIF  CML      ENDIF
94  C     preparation: find density a kth and k-1st level  C     preparation: find density a kth and k-1st level
95        CALL FIND_RHO(        CALL FIND_RHO_2D(
96       I     bi, bj, iMin, iMax, jMin, jMax, Km1, K,       I     iMin, iMax, jMin, jMax, K,
97       I     theta, salt,       I     theta(1-OLx,1-OLy,Km1,bi,bj), salt(1-OLx,1-OLy,Km1,bi,bj),
98       O     rhoKm1,       O     rhoKm1,
99       I     myThid )       I     Km1, bi, bj, myThid )
100        CALL FIND_RHO(        CALL FIND_RHO_2D(
101       I     bi, bj, iMin, iMax, jMin, jMax, K, K,       I     iMin, iMax, jMin, jMax, K,
102       I     theta, salt,       I     theta(1-OLx,1-OLy,K,bi,bj), salt(1-OLx,1-OLy,K,bi,bj),
103       O     rhoK,       O     rhoK,
104       I     myThid )       I     K, bi, bj, myThid )
105    
106  C     first step:  calculate flux Richardson number.  C     first step:  calculate flux Richardson number.
107  C     calculate (du/dz)^2 + (dv/dz)^2 on W (between T) points.  C     calculate (du/dz)^2 + (dv/dz)^2 on W (between T) points.
# Line 120  C     calculate - g*(drho/dz)/rho0= N^2 Line 120  C     calculate - g*(drho/dz)/rho0= N^2
120  C  C
121          buoyFreq = - gravity*mass2rUnit *          buoyFreq = - gravity*mass2rUnit *
122       &       (rhoKm1(I,J) - rhoK(I,J))*recip_drC(K)       &       (rhoKm1(I,J) - rhoK(I,J))*recip_drC(K)
123  C      C
124  C     calculates gradient Richardson number, bounded by  C     calculates gradient Richardson number, bounded by
125  C     a very large negative value.  C     a very large negative value.
126  C      C
127          RiNumber(I,J) = buoyFreq/max(RiFlux,epsilon)          RiNumber(I,J) = buoyFreq/max(RiFlux,epsilon)
128         ENDDO         ENDDO
129        ENDDO        ENDDO
130          
131  #ifdef PP81_SMOOTH_RI  #ifdef PP81_SMOOTH_RI
132  C     average Richardson number horizontally  C     average Richardson number horizontally
133        DO J=jMin,jMax        DO J=jMin,jMax
134         DO I=iMin,iMax         DO I=iMin,iMax
135          RiTmp(I,J) = RiNumber(I,J)          RiTmp(I,J) = RiNumber(I,J)
136         ENDDO         ENDDO
137        ENDDO            ENDDO
138        DO J=1-Oly+1,sNy+Oly-1        DO J=1-Oly+1,sNy+Oly-1
139         DO I=1-Olx+1,sNx+Olx-1         DO I=1-Olx+1,sNx+Olx-1
140          RiNumber(I,J) = p5*RiTmp(I,J)          RiNumber(I,J) = p5*RiTmp(I,J)
141       &       + p125*RiTmp(I-1,J) + p125*RiTmp(I+1,J)       &       + p125*RiTmp(I-1,J) + p125*RiTmp(I+1,J)
142       &       + p125*RiTmp(I,J-1) + p125*RiTmp(I,J+1)       &       + p125*RiTmp(I,J-1) + p125*RiTmp(I,J+1)
143         ENDDO         ENDDO
144        ENDDO            ENDDO
145  #endif /* PP81_SMOOTH_RI */  #endif /* PP81_SMOOTH_RI */
146    
147  #endif /* ALLOW_PP81 */  #endif /* ALLOW_PP81 */
Legend:
Removed from v.1.2  
changed lines
  Added in v.1.3
  ViewVC Help
Powered by ViewVC 1.1.22