75 |
C rMixingLength:: inverse of mixing length |
C rMixingLength:: inverse of mixing length |
76 |
C totalDepth :: thickness of water column (inverse of recip_Rcol) |
C totalDepth :: thickness of water column (inverse of recip_Rcol) |
77 |
C TKEPrandtlNumber :: here, an empirical function of the Richardson number |
C TKEPrandtlNumber :: here, an empirical function of the Richardson number |
|
C rhoK, rhoKm1 :: density at layer k and km1 (relative to k) |
|
78 |
INTEGER iMin ,iMax ,jMin ,jMax |
INTEGER iMin ,iMax ,jMin ,jMax |
79 |
INTEGER i, j, k, kp1, km1, kSurf, kBottom |
INTEGER i, j, k, kp1, km1, kSurf, kBottom |
80 |
_RL explDissFac, implDissFac |
_RL explDissFac, implDissFac |
95 |
_RL rMixingLength (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
_RL rMixingLength (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
96 |
_RL mxLength_Dn (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
_RL mxLength_Dn (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
97 |
_RL KappaE (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
_RL KappaE (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
|
_RL rhoK (1-OLx:sNx+OLx,1-OLy:sNy+OLy) |
|
|
_RL rhoKm1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy) |
|
98 |
_RL totalDepth (1-OLx:sNx+OLx,1-OLy:sNy+OLy) |
_RL totalDepth (1-OLx:sNx+OLx,1-OLy:sNy+OLy) |
99 |
_RL GGL90visctmp (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
_RL GGL90visctmp (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) |
100 |
C- tri-diagonal matrix |
C- tri-diagonal matrix |
149 |
ENDDO |
ENDDO |
150 |
DO j=1-OLy,sNy+OLy |
DO j=1-OLy,sNy+OLy |
151 |
DO i=1-OLx,sNx+OLx |
DO i=1-OLx,sNx+OLx |
|
rhoK(i,j) = 0. _d 0 |
|
|
rhoKm1(i,j) = 0. _d 0 |
|
152 |
totalDepth(i,j) = Ro_surf(i,j,bi,bj) - R_low(i,j,bi,bj) |
totalDepth(i,j) = Ro_surf(i,j,bi,bj) - R_low(i,j,bi,bj) |
153 |
rMixingLength(i,j,1) = 0. _d 0 |
rMixingLength(i,j,1) = 0. _d 0 |
154 |
mxLength_Dn(i,j,1) = GGL90mixingLengthMin |
mxLength_Dn(i,j,1) = GGL90mixingLengthMin |
158 |
|
|
159 |
C start k-loop |
C start k-loop |
160 |
DO k = 2, Nr |
DO k = 2, Nr |
161 |
km1 = k-1 |
c km1 = k-1 |
162 |
c kp1 = MIN(Nr,k+1) |
c kp1 = MIN(Nr,k+1) |
|
CALL FIND_RHO_2D( |
|
|
I iMin, iMax, jMin, jMax, k, |
|
|
I theta(1-OLx,1-OLy,km1,bi,bj), salt(1-OLx,1-OLy,km1,bi,bj), |
|
|
O rhoKm1, |
|
|
I km1, bi, bj, myThid ) |
|
|
|
|
|
CALL FIND_RHO_2D( |
|
|
I iMin, iMax, jMin, jMax, k, |
|
|
I theta(1-OLx,1-OLy,k,bi,bj), salt(1-OLx,1-OLy,k,bi,bj), |
|
|
O rhoK, |
|
|
I k, bi, bj, myThid ) |
|
163 |
DO j=jMin,jMax |
DO j=jMin,jMax |
164 |
DO i=iMin,iMax |
DO i=iMin,iMax |
165 |
SQRTTKE(i,j,k)=SQRT( GGL90TKE(i,j,k,bi,bj) ) |
SQRTTKE(i,j,k)=SQRT( GGL90TKE(i,j,k,bi,bj) ) |
167 |
C buoyancy frequency |
C buoyancy frequency |
168 |
Nsquare(i,j,k) = gravity*gravitySign*recip_rhoConst |
Nsquare(i,j,k) = gravity*gravitySign*recip_rhoConst |
169 |
& * sigmaR(i,j,k) |
& * sigmaR(i,j,k) |
|
Nsquare(i,j,k) = - gravity*recip_rhoConst*recip_drC(k) |
|
|
& * ( rhoKm1(i,j) - rhoK(i,j) )*maskC(i,j,k,bi,bj) |
|
170 |
cC vertical shear term (dU/dz)^2+(dV/dz)^2 |
cC vertical shear term (dU/dz)^2+(dV/dz)^2 |
171 |
c tempU= .5 _d 0*( uVel(i,j,km1,bi,bj)+uVel(i+1,j,km1,bi,bj) |
c tempU= .5 _d 0*( uVel(i,j,km1,bi,bj)+uVel(i+1,j,km1,bi,bj) |
172 |
c & -( uVel(i,j,k ,bi,bj)+uVel(i+1,j,k ,bi,bj)) ) |
c & -( uVel(i,j,k ,bi,bj)+uVel(i+1,j,k ,bi,bj)) ) |