/[MITgcm]/MITgcm_contrib/submesoscale/code/calc_oce_mxlayer.F
ViewVC logotype

Annotation of /MITgcm_contrib/submesoscale/code/calc_oce_mxlayer.F

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


Revision 1.3 - (hide annotations) (download)
Fri Mar 12 18:31:00 2010 UTC (15 years, 4 months ago) by zhc
Branch: MAIN
CVS Tags: HEAD
Changes since 1.2: +47 -18 lines
updated with checkpoint62c

1 zhc 1.3 C $Header: /u/gcmpack/MITgcm/model/src/calc_oce_mxlayer.F,v 1.8 2009/10/16 17:27:33 dimitri Exp $
2 dimitri 1.1 C $Name: $
3    
4     #include "PACKAGES_CONFIG.h"
5     #include "CPP_OPTIONS.h"
6    
7     CBOP
8     C !ROUTINE: CALC_OCE_MXLAYER
9     C !INTERFACE:
10     SUBROUTINE CALC_OCE_MXLAYER(
11     I rhoSurf, sigmaR,
12     I bi, bj, myTime, myIter, myThid )
13     C !DESCRIPTION: \bv
14     C *==========================================================*
15     C | S/R CALC_OCE_MXLAYER
16     C | o Diagnose the Oceanic surface Mixed-Layer
17     C *==========================================================*
18     C \ev
19    
20     C !USES:
21     IMPLICIT NONE
22     C == Global variables ==
23     #include "SIZE.h"
24     #include "EEPARAMS.h"
25     #include "PARAMS.h"
26     #include "GRID.h"
27     #include "DYNVARS.h"
28     #ifdef ALLOW_GMREDI
29     # include "GMREDI.h"
30     #endif
31    
32     C !INPUT/OUTPUT PARAMETERS:
33     C == Routine Arguments ==
34     C rhoSurf :: Surface density anomaly
35     C sigmaR :: Vertical gradient of potential density
36     C bi,bj :: tile indices
37     C myTime :: Current time in simulation
38     C myIter :: Current iteration number in simulation
39     C myThid :: my Thread Id number
40     _RL rhoSurf(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
41     _RL sigmaR (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
42     INTEGER bi, bj
43     _RL myTime
44     INTEGER myIter
45     INTEGER myThid
46    
47     C === Functions ====
48     #ifdef ALLOW_DIAGNOSTICS
49     LOGICAL DIAGNOSTICS_IS_ON
50     EXTERNAL DIAGNOSTICS_IS_ON
51     #endif /* ALLOW_DIAGNOSTICS */
52    
53     C !LOCAL VARIABLES:
54     C == Local variables ==
55     C i,j :: Loop counters
56     INTEGER i,j,k
57     LOGICAL calcMixLayerDepth
58     INTEGER method
59 zhc 1.3 _RL rhoBigNb
60 dimitri 1.1 _RL rhoMxL(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
61     _RL rhoKm1(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
62     _RL rhoLoc(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
63     _RL tmpFac, sigmAv
64     CEOP
65    
66     calcMixLayerDepth = .FALSE.
67     #ifdef ALLOW_GMREDI
68 zhc 1.3 IF ( useGMRedi ) THEN
69     IF ( .NOT.useKPP ) calcMixLayerDepth = GM_taper_scheme.EQ.'fm07'
70     C SM(1)
71     IF ( GM_SM_Ce.GT.0.0 ) calcMixLayerDepth = .TRUE.
72     ENDIF
73 dimitri 1.1 #endif
74     #ifdef ALLOW_DIAGNOSTICS
75     IF ( useDiagnostics.AND. .NOT.calcMixLayerDepth ) THEN
76     calcMixLayerDepth = DIAGNOSTICS_IS_ON('MXLDEPTH',myThid)
77     ENDIF
78     #endif
79     IF ( calcMixLayerDepth ) THEN
80    
81     C-- Select which "method" to use:
82     method = 0
83     IF ( hMixCriteria.LT.0. ) method = 1
84     IF ( hMixCriteria.GT.1. ) method = 2
85    
86     IF ( method.EQ.1 ) THEN
87    
88     C-- First method :
89     C where the potential density (ref.lev=surface) is larger than
90     C surface density plus Delta_Rho = hMixCriteria * Alpha(surf)
91     C = density of water which is -hMixCriteria colder than surface water
92 zhc 1.3 C (see Kara, Rochford, and Hurlburt JGR 2000 for default criterion)
93 dimitri 1.1
94     c hMixCriteria = -0.8 _d 0
95 zhc 1.3 c dRhoSmall = 1. _d -6
96 dimitri 1.1 rhoBigNb = rhoConst*1. _d 10
97     CALL FIND_ALPHA(
98     I bi, bj, 1-OLx, sNx+OLx, 1-OLy, sNy+OLy, 1, 1,
99     O rhoMxL, myThid )
100    
101     DO j=1-Oly,sNy+Oly
102     DO i=1-Olx,sNx+Olx
103     rhoKm1(i,j) = rhoSurf(i,j)
104     rhoMxL(i,j) = rhoSurf(i,j)
105     & + MAX( rhoMxL(i,j)*hMixCriteria, dRhoSmall )
106 zhc 1.3 hMixLayer(i,j,bi,bj) = rF(1)-R_low(i,j,bi,bj)
107 dimitri 1.1 ENDDO
108     ENDDO
109     DO k = 2,Nr
110     C- potential density (reference level = surface level)
111 zhc 1.3 CALL FIND_RHO_2D(
112     I 1-OLx, sNx+OLx, 1-OLy, sNy+OLy, 1,
113     I theta(1-OLx,1-OLy,k,bi,bj), salt(1-OLx,1-OLy,k,bi,bj),
114     O rhoLoc,
115     I k, bi, bj, myThid )
116 dimitri 1.1
117     DO j=1-Oly,sNy+Oly
118     DO i=1-Olx,sNx+Olx
119     IF ( k.LE.klowC(i,j,bi,bj) .AND.
120     & rhoLoc(i,j).GE.rhoMxL(i,j) ) THEN
121     IF ( rhoLoc(i,j).GT.rhoKm1(i,j) ) THEN
122     tmpFac = ( rhoMxL(i,j) - rhoKm1(i,j) )
123     & / ( rhoLoc(i,j) - rhoKm1(i,j) )
124     ELSE
125     tmpFac = 0.
126     ENDIF
127     hMixLayer(i,j,bi,bj) = rF(1)-rC(k-1)+tmpFac*drC(k)
128     rhoMxL(i,j) = rhoBigNb
129     ELSE
130     rhoKm1(i,j) = rhoLoc(i,j)
131     ENDIF
132     ENDDO
133     ENDDO
134     ENDDO
135    
136     ELSEIF ( method.EQ.2 ) THEN
137    
138     C-- Second method :
139     C where the local stratification exceed the mean stratification above
140     C (from surface down to here) by factor hMixCriteria
141    
142 zhc 1.3 c hMixCriteria = 1.5 _d 0
143     c dRhoSmall = 1. _d -2
144 dimitri 1.1 DO j=1-Oly,sNy+Oly
145     DO i=1-Olx,sNx+Olx
146     IF ( klowC(i,j,bi,bj) .GT. 0 ) THEN
147     hMixLayer(i,j,bi,bj) = drF(1)
148     rhoMxL(i,j) = 1.
149     ELSE
150     hMixLayer(i,j,bi,bj) = rF(1)
151     rhoMxL(i,j) = -1.
152     ENDIF
153     ENDDO
154     ENDDO
155     DO k = 2,Nr-1
156     C- potential density (reference level = surface level)
157 zhc 1.3 CALL FIND_RHO_2D(
158     I 1-OLx, sNx+OLx, 1-OLy, sNy+OLy, 1,
159     I theta(1-OLx,1-OLy,k,bi,bj), salt(1-OLx,1-OLy,k,bi,bj),
160     O rhoLoc,
161     I k, bi, bj, myThid )
162 dimitri 1.1
163     DO j=1-Oly,sNy+Oly
164     DO i=1-Olx,sNx+Olx
165     IF ( k.LT.klowC(i,j,bi,bj)
166     & .AND. rhoMxL(i,j).GE.0. ) THEN
167     sigmAv = ( rhoLoc(i,j)-rhoSurf(i,j)+dRhoSmall )
168     & / ( rC(1)-rC(k) )
169     IF ( -sigmaR(i,j,k+1).GT.sigmAv*hMixCriteria ) THEN
170     tmpFac = 0. _d 0
171     IF ( sigmAv.GT.0. _d 0 ) THEN
172     tmpFac = hMixCriteria*sigmaR(i,j,k)/sigmaR(i,j,k+1)
173     IF ( tmpFac .GT. 1. _d 0 ) THEN
174     tmpFac = 1. _d 0
175     & + ( tmpFac - 1. _d 0 )/( hMixCriteria - 1. _d 0 )
176     ENDIF
177     tmpFac = MAX( 0. _d 0, MIN( tmpFac, 2. _d 0 ) )
178     ENDIF
179     hMixLayer(i,j,bi,bj) = rF(1)-rF(k+1)
180     & - drF(k)*tmpFac*0.5 _d 0
181     rhoMxL(i,j) = -1.
182     ENDIF
183     ENDIF
184     ENDDO
185     ENDDO
186     ENDDO
187    
188     ELSE
189     STOP 'S/R CALC_OCE_MXLAYER: invalid method'
190     ENDIF
191    
192 zhc 1.3 IF ( hMixSmooth .GT. 0. _d 0 ) THEN
193     tmpFac = (1. _d 0 - hMixSmooth ) / 4. _d 0
194     DO j=1-Oly+1,sNy+Oly-1
195     DO i=1-Olx+1,sNx+Olx-1
196     rhoLoc(i,j)=(hMixSmooth * hMixLayer(i,j,bi,bj) +
197     & tmpFac * ( hMixLayer(i-1,j,bi,bj) +
198     & hMixLayer(i+1,j,bi,bj) +
199     & hMixLayer(i,j-1,bi,bj) +
200     & hMixLayer(i,j+1,bi,bj) )
201     & )
202     & /(hMixSmooth +
203     & tmpFac * ( maskC(i-1,j,1,bi,bj) +
204     & maskC(i+1,j,1,bi,bj) +
205     & maskC(i,j-1,1,bi,bj) +
206     & maskC(i,j+1,1,bi,bj) )
207     & ) * maskC(i,j,1,bi,bj)
208     ENDDO
209     ENDDO
210     DO j=1-Oly+1,sNy+Oly-1
211     DO i=1-Olx+1,sNx+Olx-1
212     hMixLayer(i,j,bi,bj) = rhoLoc(i,j)
213     ENDDO
214     ENDDO
215     ENDIF
216    
217 dimitri 1.1 #ifdef ALLOW_DIAGNOSTICS
218     IF ( useDiagnostics ) THEN
219     CALL DIAGNOSTICS_FILL( hMixLayer, 'MXLDEPTH',
220     & 0, 1, 1, bi, bj, myThid )
221     ENDIF
222     #endif /* ALLOW_DIAGNOSTICS */
223    
224     C-- end if calcMixLayerDepth
225     ENDIF
226    
227     RETURN
228     END

  ViewVC Help
Powered by ViewVC 1.1.22