File indexing completed on 2021-11-06 05:18:23 UTC
view on githubraw file Latest commit 016b84c4 on 2021-11-02 20:24:44 UTC
08be60903a Mart*0001 #include "MY82_OPTIONS.h"
0002
0003
0004
0005
0006
0007 subroutine MY82_CALC(
dc6107c029 Jean*0008 I bi, bj, sigmaR, myTime, myIter, myThid )
08be60903a Mart*0009
0010
5e48dccc42 Jean*0011
08be60903a Mart*0012
0013
5e48dccc42 Jean*0014
08be60903a Mart*0015
5e48dccc42 Jean*0016
08be60903a Mart*0017 IMPLICIT NONE
0018
0019
0020
0021
016b84c482 Mart*0022
0023
08be60903a Mart*0024
0025
0026
0027
0028 #include "SIZE.h"
0029 #include "EEPARAMS.h"
0030 #include "PARAMS.h"
0031 #include "GRID.h"
dbbea2be4e Jean*0032 #include "MY82.h"
0033 #ifdef ALLOW_3D_DIFFKR
0034 # include "DYNVARS.h"
0035 #endif
08be60903a Mart*0036
0037
dc6107c029 Jean*0038
0039
0040
0041
0042
0043
08be60903a Mart*0044 INTEGER bi, bj
dc6107c029 Jean*0045 _RL sigmaR(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
08be60903a Mart*0046 _RL myTime
dc6107c029 Jean*0047 INTEGER myIter
fb62f539dc Jean*0048 INTEGER myThid
08be60903a Mart*0049
0050
0051
0052
0053
0054
0055
0056
0057
0058 INTEGER I, J, K
0059 INTEGER iMin ,iMax ,jMin ,jMax
0060 _RL RiFlux, RiTmp
0061 _RL SHtmp, bTmp, tkesquare, tkel
0062 _RL RiNumber(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0063 _RL GH(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0064 _RL GM(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0065 _RL SH(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
0066 _RL SM(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
0067 _RL tke(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
0068
0069 iMin = 2-OLx
0070 iMax = sNx+OLx-1
0071 jMin = 2-OLy
0072 jMax = sNy+OLy-1
0073
0074
dc6107c029 Jean*0075 DO J=1-OLy,sNy+OLy
0076 DO I=1-OLx,sNx+OLx
fb62f539dc Jean*0077 GH(I,J) = 0. _d 0
0078 GM(I,J) = 0. _d 0
08be60903a Mart*0079 ENDDO
0080 ENDDO
0081 DO K = 1, Nr
dc6107c029 Jean*0082 DO J=1-OLy,sNy+OLy
0083 DO I=1-OLx,sNx+OLx
8fef4b1a57 Mart*0084 SH(I,J,K) = 0. _d 0
0085 SM(I,J,K) = 0. _d 0
0086 tke(I,J,K) = 0. _d 0
08be60903a Mart*0087 ENDDO
fb62f539dc Jean*0088 ENDDO
08be60903a Mart*0089 ENDDO
0090
0091
0092 DO K = 2, Nr
0093 CALL MY82_RI_NUMBER(
fb62f539dc Jean*0094 I bi, bj, K, iMin, iMax, jMin, jMax,
08be60903a Mart*0095 O RiNumber, GH, GM,
0096 I myTime, myThid )
0097 DO J=jMin,jMax
0098 DO I=iMin,iMax
0099 RiTmp = MIN(RiNumber(I,J),RiMax)
0100 btmp = beta1+beta4*RiTmp
0101
0102
fb62f539dc Jean*0103 RiFlux =( btmp - SQRT(btmp*btmp - 4. _d 0 *beta2*beta3*RiTmp) )
0104 & /(2. _d 0*beta2)
08be60903a Mart*0105
fb62f539dc Jean*0106 SHtmp = (alpha1-alpha2*RiFlux)/(1. _d 0-RiFlux)
08be60903a Mart*0107 SH(I,J,K) = SHtmp
0108 SM(I,J,K) = SHtmp*(beta1-beta2*RiFlux)/(beta3-beta4*RiFlux)
0109
0110 tkesquare = MAX( 0. _d 0,
0111 & b1*(SH(I,J,K)*GH(I,J) + SM(I,J,K)*GM(I,J)) )
0112 tke(I,J,K) = sqrt(tkesquare)
0113
0114
0115
0116
0117 ENDDO
0118 ENDDO
0119
fb62f539dc Jean*0120 ENDDO
08be60903a Mart*0121
fb62f539dc Jean*0122
08be60903a Mart*0123
0124 DO J=jMin,jMax
0125 DO I=iMin,iMax
8fef4b1a57 Mart*0126 GH(I,J) = 0. _d 0
0127 GM(I,J) = 0. _d 0
08be60903a Mart*0128 ENDDO
0129 ENDDO
0130
0131
0132
0133
0134 DO K = 2, Nr
0135 DO J=jMin,jMax
0136 DO I=iMin,iMax
0137 GM(I,J) = GM(I,J) + tke(I,J,K)*rF(K)
fb62f539dc Jean*0138 GH(I,J) = GH(I,J) + tke(I,J,K)
08be60903a Mart*0139 ENDDO
0140 ENDDO
0141
0142 END DO
0143
0144 DO J=jMin,jMax
0145 DO I=iMin,iMax
8fef4b1a57 Mart*0146 IF ( GH(I,J) .EQ. 0. _d 0 ) THEN
0147 MYhbl(I,J,bi,bj) = 0. _d 0
08be60903a Mart*0148 ELSE
0149 MYhbl(I,J,bi,bj) = -GM(I,J)/GH(I,J)*MYhblScale
0150 ENDIF
0151 ENDDO
fb62f539dc Jean*0152 ENDDO
08be60903a Mart*0153
0154 DO K = 1, Nr
0155
0156 DO J=jMin,jMax
0157 DO I=iMin,iMax
0158 tkel = MYhbl(I,J,bi,bj)*tke(I,J,K)
0159
0160 MYviscAr(I,J,K,bi,bj) = MYhbl(I,J,bi,bj)*tkel*SM(I,J,K)
0161 MYdiffKr(I,J,K,bi,bj) = MYhbl(I,J,bi,bj)*tkel*SH(I,J,K)
0162
5e48dccc42 Jean*0163 MYviscAr(I,J,K,bi,bj) = MAX(MYviscAr(I,J,K,bi,bj),
0164 & viscArnr(k) )
648edfaf30 Jean*0165 MYdiffKr(I,J,K,bi,bj) = MAX(MYdiffKr(I,J,K,bi,bj),
dbbea2be4e Jean*0166 #ifdef ALLOW_3D_DIFFKR
0167 & diffKr(i,j,k,bi,bj) )
0168 #else
78524d1402 Jean*0169 & diffKrNrS(k) )
dbbea2be4e Jean*0170 #endif
08be60903a Mart*0171
0172 MYviscAr(I,J,K,bi,bj) = MIN(MYviscAr(I,J,K,bi,bj),MYviscMax)
0173 & * maskC(I,J,K,bi,bj)
0174 MYdiffKr(I,J,K,bi,bj) = MIN(MYdiffKr(I,J,K,bi,bj),MYdiffMax)
0175 & * maskC(I,J,K,bi,bj)
0176 ENDDO
fb62f539dc Jean*0177 ENDDO
08be60903a Mart*0178
fb62f539dc Jean*0179 ENDDO
08be60903a Mart*0180
016b84c482 Mart*0181 #ifdef ALLOW_DIAGNOSTICS
0182 IF ( useDiagnostics ) THEN
0183 CALL DIAGNOSTICS_FILL(MYviscAr,'MYVISCAR',0,Nr,1,bi,bj,myThid)
0184 CALL DIAGNOSTICS_FILL(MYdiffKr,'MYDIFFKR',0,Nr,1,bi,bj,myThid)
0185 CALL DIAGNOSTICS_FILL(MYhbl ,'MYHBL ',0,1 ,1,bi,bj,myThid)
0186 ENDIF
0187 #endif /* ALLOW_DIAGNOSTICS */
08be60903a Mart*0188
0189 RETURN
0190 END