File indexing completed on 2024-10-18 05:11:24 UTC
view on githubraw file Latest commit 5bb179dd on 2024-10-17 18:00:27 UTC
6e2f4e58fa Mart*0001 #include "SEAICE_OPTIONS.h"
fec75090d3 Jean*0002 #ifdef ALLOW_OBCS
0003 # include "OBCS_OPTIONS.h"
0004 #else
0005 # define OBCS_UVICE_OLD
0006 #endif
772b2ed80e Gael*0007 #ifdef ALLOW_AUTODIFF
0008 # include "AUTODIFF_OPTIONS.h"
0009 #endif
6e2f4e58fa Mart*0010
c9a6744012 Jean*0011
0012
0013
6e2f4e58fa Mart*0014 SUBROUTINE SEAICE_EVP( myTime, myIter, myThid )
c9a6744012 Jean*0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
6e2f4e58fa Mart*0029 IMPLICIT NONE
0030
0031
0032 #include "SIZE.h"
0033 #include "EEPARAMS.h"
0034 #include "PARAMS.h"
bcf4255045 Dimi*0035 #include "DYNVARS.h"
6e2f4e58fa Mart*0036 #include "GRID.h"
ccaa3c61f4 Patr*0037 #include "SEAICE_SIZE.h"
6e2f4e58fa Mart*0038 #include "SEAICE_PARAMS.h"
ccaa3c61f4 Patr*0039 #include "SEAICE.h"
6e2f4e58fa Mart*0040
0041 #ifdef ALLOW_AUTODIFF_TAMC
0042 # include "tamc.h"
0043 #endif
0044
c9a6744012 Jean*0045
6e2f4e58fa Mart*0046
c9a6744012 Jean*0047
0048
0049
6e2f4e58fa Mart*0050 _RL myTime
0051 INTEGER myIter
0052 INTEGER myThid
0053
45315406aa Mart*0054 #if ( defined SEAICE_CGRID && defined SEAICE_ALLOW_EVP )
6e2f4e58fa Mart*0055
0056
c9a6744012 Jean*0057
0058
0059
0060
0061
a9588c2504 Mart*0062
0063
0064
0065
0066
0067
c9a6744012 Jean*0068
0069
0070
0071
0072
0073
0074
0075
e4120ef5ff Jean*0076
85f5225996 Mart*0077
a6389affa6 Mart*0078
0079
0080
0081
85f5225996 Mart*0082
0083
e4120ef5ff Jean*0084
85f5225996 Mart*0085
e4120ef5ff Jean*0086
85f5225996 Mart*0087
e4120ef5ff Jean*0088
0089
85f5225996 Mart*0090
e4120ef5ff Jean*0091
0092
85f5225996 Mart*0093
0094
b4949dd6db Jean*0095
6e2f4e58fa Mart*0096 INTEGER i, j, bi, bj
282fbc7b4c Mart*0097 INTEGER kSrf
6e2f4e58fa Mart*0098 INTEGER nEVPstep, iEVPstep
8377b8ee87 Mart*0099 #ifndef ALLOW_AUTODIFF
0100 INTEGER nEVPstepMax
0101 #endif
452bde753c Mart*0102 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0103
0104
0105 INTEGER evpkey, tkey
cf9fa44a59 Patr*0106 #endif
6e2f4e58fa Mart*0107
35fda33b05 Jean*0108 _RL COSWAT
0109 _RS SINWAT
a9588c2504 Mart*0110 _RL ecc2, recip_ecc2, recip_evpAlpha, recip_deltaT
0912318009 Mart*0111 _RL betaFacP1, betaFac, evpStarFac, evpRevFac, recip_evpRevFac
6e2f4e58fa Mart*0112
b7030e3a29 Mart*0113 _RL seaice_div (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0114 _RL seaice_tension(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0115 _RL seaice_shear (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
68ff576152 Mart*0116 _RL sig11 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0117 _RL sig22 (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
70e078b38a Mart*0118
4d65bb090e Mart*0119 _RL areaW (1:sNx,1:sNy,nSx,nSy)
0120 _RL areaS (1:sNx,1:sNy,nSx,nSy)
037351a1a6 Mart*0121
da4a30d88f Mart*0122 _RL ep (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0123 _RL em (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
d2082aecba Mart*0124 _RL e12Csq (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
1dd53d2d6a Mart*0125 _RL pressC (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0126 _RL zetaC (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0127 _RL deltaZ (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0fd6b2de1a Mart*0128
8ce59fd29e Mart*0129 #ifdef SEAICE_ALLOW_MOM_ADVECTION
0130
0131 _RL gUmom (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0132 _RL gVmom (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0133 #endif /* SEAICE_ALLOW_MOM_ADVECTION */
e4120ef5ff Jean*0134 _RL deltaCreg, deltaSq, deltaMinSq, tmp
452bde753c Mart*0135 #ifdef SEAICE_ALLOW_TEM
b8e39f4242 Mart*0136 _RL etaDenC, zetaMaxC, etaDenZ, zetaMaxZ
452bde753c Mart*0137 #endif /* SEAICE_ALLOW_TEM */
fc27138d73 Mart*0138 #ifdef SEAICE_ALLOW_CLIPZETA
0139 _RL zMaxZ, zMinZ, fac
0140 #endif /* SEAICE_ALLOW_CLIPZETA */
e4120ef5ff Jean*0141 _RL denom1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
a6389affa6 Mart*0142 _RL denom2 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
967b0b9f87 Mart*0143 _RL sumNorm, denomU, denomV
fec75090d3 Jean*0144 _RL locMaskU, locMaskV
a6389affa6 Mart*0145 _RL EVPcFac
85f5225996 Mart*0146 _RL evpAlphaC(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
a6389affa6 Mart*0147 _RL evpAlphaZ(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
0148 _RL evpBetaU (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
0149 _RL evpBetaV (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
85f5225996 Mart*0150 _RL betaFacP1U, betaFacP1V
0151 _RL betaFacU, betaFacV
a6389affa6 Mart*0152 LOGICAL useAdaptiveEVP
85f5225996 Mart*0153
e4120ef5ff Jean*0154 #ifdef ALLOW_SEAICE_EVP_RESIDUAL
85f5225996 Mart*0155 _RL resTile(nSx,nSy)
0156 _RL resLoc
4669aaa4e4 Mart*0157 _RL uIcePm1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
0158 _RL vIcePm1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
0159 _RL sig11pm1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
0160 _RL sig22pm1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
0161 _RL sig12pm1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
85f5225996 Mart*0162
0163 LOGICAL printResidual
0164
0165
0166 LOGICAL DIFFERENT_MULTIPLE
0167 EXTERNAL DIFFERENT_MULTIPLE
0168 #endif /* ALLOW_SEAICE_EVP_RESIDUAL */
0169
6e2f4e58fa Mart*0170
a6389affa6 Mart*0171
0172 useAdaptiveEVP = .FALSE.
0173 IF ( SEAICEaEvpCoeff .NE. UNSET_RL ) useAdaptiveEVP = .TRUE.
0174 EVPcFac = 0. _d 0
0175 IF ( useAdaptiveEVP )
0176 & EVPcFac = SEAICE_deltaTdyn*SEAICEaEVPcStar
0177 & * (SEAICEaEvpCoeff * PI)**2
85f5225996 Mart*0178
e4120ef5ff Jean*0179 #ifdef ALLOW_SEAICE_EVP_RESIDUAL
85f5225996 Mart*0180 printResidual = debugLevel.GE.debLevA
0181 & .AND. DIFFERENT_MULTIPLE( SEAICE_monFreq, myTime, deltaTClock )
0182 #endif /* ALLOW_SEAICE_EVP_RESIDUAL */
c9a6744012 Jean*0183
0320e25227 Mart*0184 IF ( usingPCoords ) THEN
0185 kSrf = Nr
0186 ELSE
0187 kSrf = 1
0188 ENDIF
6e2f4e58fa Mart*0189
0190 SINWAT=SIN(SEAICE_waterTurnAngle*deg2rad)
0191 COSWAT=COS(SEAICE_waterTurnAngle*deg2rad)
0192
0193
0194 ecc2=SEAICE_eccen**2
037351a1a6 Mart*0195 recip_ecc2 = 0. _d 0
0196 IF ( ecc2 .NE. 0. _d 0 ) recip_ecc2=ONE/ecc2
3daf25222c Mart*0197 deltaMinSq = SEAICE_deltaMin**2
dc26f158aa Mart*0198
0199 nEVPstep = SEAICEnEVPstarSteps
fb2ddc4f23 Mart*0200
a6389affa6 Mart*0201 DO bj=myByLo(myThid),myByHi(myThid)
0202 DO bi=myBxLo(myThid),myBxHi(myThid)
8377b8ee87 Mart*0203 DO j=1-OLy,sNy+OLy
0204 DO i=1-OLx,sNx+OLx
0205 denom1(i,j,bi,bj) = 1. _d 0 / ( SEAICE_evpAlpha + 1. _d 0 )
0206 denom2(i,j,bi,bj) = 1. _d 0 / ( SEAICE_evpAlpha + ecc2 )
a6389affa6 Mart*0207 ENDDO
0208 ENDDO
0209 ENDDO
0210 ENDDO
2bf4ed5c39 Mart*0211 recip_deltaT = 1. _d 0 / SEAICE_deltaTdyn
0212 recip_evpAlpha = 0. _d 0
e4120ef5ff Jean*0213 IF ( SEAICE_evpAlpha .GT. 0. _d 0 )
2bf4ed5c39 Mart*0214 & recip_evpAlpha = 1. _d 0 / SEAICE_evpAlpha
0215 evpStarFac = 0. _d 0
0912318009 Mart*0216 evpRevFac = 0. _d 0
0217 recip_evpRevFac = 1. _d 0
2bf4ed5c39 Mart*0218 IF ( SEAICEuseEVPstar ) evpStarFac = 1. _d 0
0912318009 Mart*0219 IF ( SEAICEuseEVPrev ) THEN
e4120ef5ff Jean*0220
d2082aecba Mart*0221
0912318009 Mart*0222 evpRevFac = 1. _d 0
0223 evpStarFac = 1. _d 0
0224 recip_evpRevFac = recip_ecc2
a6389affa6 Mart*0225 DO bj=myByLo(myThid),myByHi(myThid)
0226 DO bi=myBxLo(myThid),myBxHi(myThid)
8377b8ee87 Mart*0227 DO j=1-OLy,sNy+OLy
0228 DO i=1-OLx,sNx+OLx
0229 denom1(i,j,bi,bj) = 1. _d 0 / SEAICE_evpAlpha
0230 denom2(i,j,bi,bj) = denom1(i,j,bi,bj)
a6389affa6 Mart*0231 ENDDO
0232 ENDDO
0233 ENDDO
0234 ENDDO
0912318009 Mart*0235 ENDIF
85f5225996 Mart*0236
0237 betaFac = SEAICE_evpBeta*recip_deltaT
0238 betaFacU = betaFac
0239 betaFacV = betaFac
0240
0241 betaFacP1 = betaFac + evpStarFac*recip_deltaT
0242 betaFacP1U = betaFacP1
0243 betaFacP1V = betaFacP1
8377b8ee87 Mart*0244 #ifndef ALLOW_AUTODIFF
cf9fa44a59 Patr*0245 nEVPstepMax = nEVPstep
0246 #endif
6e2f4e58fa Mart*0247
0248 DO bj=myByLo(myThid),myByHi(myThid)
0249 DO bi=myBxLo(myThid),myBxHi(myThid)
9e706e0219 Jean*0250 DO j=1-OLy,sNy+OLy
0251 DO i=1-OLx,sNx+OLx
a3f6dcaab3 Mart*0252
8377b8ee87 Mart*0253 uIceNm1(i,j,bi,bj) = uIce(i,j,bi,bj)
0254 vIceNm1(i,j,bi,bj) = vIce(i,j,bi,bj)
c9a6744012 Jean*0255
8377b8ee87 Mart*0256 e11 (i,j,bi,bj) = 0. _d 0
0257 e22 (i,j,bi,bj) = 0. _d 0
0258 e12 (i,j,bi,bj) = 0. _d 0
a6389affa6 Mart*0259
8377b8ee87 Mart*0260 evpAlphaC(i,j,bi,bj) = SEAICE_evpAlpha
0261 evpAlphaZ(i,j,bi,bj) = SEAICE_evpAlpha
0262 evpBetaU (i,j,bi,bj) = SEAICE_evpBeta
0263 evpBetaV (i,j,bi,bj) = SEAICE_evpBeta
6e2f4e58fa Mart*0264 ENDDO
0265 ENDDO
4d65bb090e Mart*0266
0267 IF ( SEAICEscaleSurfStress ) THEN
8377b8ee87 Mart*0268 DO j=1,sNy
0269 DO i=1,sNx
0270 areaW(i,j,bi,bj) =
0271 & 0.5 _d 0*(AREA(i,j,bi,bj)+AREA(i-1,j,bi,bj))
0272 areaS(i,j,bi,bj) =
0273 & 0.5 _d 0*(AREA(i,j,bi,bj)+AREA(i,j-1,bi,bj))
4d65bb090e Mart*0274 ENDDO
0275 ENDDO
0276 ELSE
8377b8ee87 Mart*0277 DO j=1,sNy
0278 DO i=1,sNx
0279 areaW(i,j,bi,bj) = 1. _d 0
0280 areaS(i,j,bi,bj) = 1. _d 0
4d65bb090e Mart*0281 ENDDO
0282 ENDDO
0283 ENDIF
cf9fa44a59 Patr*0284 ENDDO
0285 ENDDO
fc27138d73 Mart*0286 #ifdef SEAICE_ALLOW_CLIPZETA
8a89485793 Mart*0287
0288 IF ( SEAICE_evpDampC .GT. 0. _d 0 ) THEN
2bf4ed5c39 Mart*0289
0290
0291 fac = 0.25 _d 0 * SEAICE_evpDampC * SEAICE_evpAlpha
0292 & /SEAICE_deltaTevp
8a89485793 Mart*0293 DO bj=myByLo(myThid),myByHi(myThid)
0294 DO bi=myBxLo(myThid),myBxHi(myThid)
9e706e0219 Jean*0295 DO j=1-OLy,sNy+OLy
0296 DO i=1-OLx,sNx+OLx
8e32c48b8f Mart*0297 SEAICE_zMax(i,j,bi,bj) = _rA(i,j,bi,bj) * fac
8a89485793 Mart*0298 ENDDO
0299 ENDDO
0300 ENDDO
0301 ENDDO
0302 ENDIF
fc27138d73 Mart*0303 #endif /* SEAICE_ALLOW_CLIPZETA */
b4949dd6db Jean*0304
6e2f4e58fa Mart*0305
cf9fa44a59 Patr*0306 DO iEVPstep = 1, nEVPstepMax
0307 IF (iEVPstep.LE.nEVPstep) THEN
f6a7d7e353 Patr*0308
6e2f4e58fa Mart*0309 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0310 evpkey = iEVPstep + (ikey_dynamics-1)*nEVPstepMax
0311
0312
0313
0314
0315
6e2f4e58fa Mart*0316 #endif /* ALLOW_AUTODIFF_TAMC */
0317
0318
b4949dd6db Jean*0319
0320 CALL SEAICE_CALC_STRAINRATES(
a3f6dcaab3 Mart*0321 I uIce, vIce,
037351a1a6 Mart*0322 O e11, e22, e12,
2e75568507 Mart*0323 I iEVPstep, myTime, myIter, myThid )
037351a1a6 Mart*0324
88601fb2cb Gael*0325 #ifdef ALLOW_AUTODIFF_TAMC
5bb179ddc2 Mart*0326
edb6656069 Mart*0327
5bb179ddc2 Mart*0328
88601fb2cb Gael*0329 #endif /* ALLOW_AUTODIFF_TAMC */
0330
037351a1a6 Mart*0331 DO bj=myByLo(myThid),myByHi(myThid)
0332 DO bi=myBxLo(myThid),myBxHi(myThid)
77fc19a7e9 Mart*0333 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0334 tkey = bi + (bj-1)*nSx + (evpkey-1)*nSx*nSy
77fc19a7e9 Mart*0335 #endif /* ALLOW_AUTODIFF_TAMC */
85f5225996 Mart*0336 #ifdef ALLOW_SEAICE_EVP_RESIDUAL
0337
0338 IF ( printResidual ) THEN
0339 DO j=1,sNy
0340 DO i=1,sNx
8377b8ee87 Mart*0341 sig11Pm1(i,j,bi,bj) = seaice_sigma1(i,j,bi,bj)
0342 sig22Pm1(i,j,bi,bj) = seaice_sigma2(i,j,bi,bj)
0343 sig12Pm1(i,j,bi,bj) = seaice_sigma12(i,j,bi,bj)
0344 uIcePm1 (i,j,bi,bj) = uIce(i,j,bi,bj)
0345 vIcePm1 (i,j,bi,bj) = vIce(i,j,bi,bj)
85f5225996 Mart*0346 ENDDO
0347 ENDDO
0348 ENDIF
0349 #endif /* ALLOW_SEAICE_EVP_RESIDUAL */
9e706e0219 Jean*0350 DO j=1-OLy,sNy+OLy
0351 DO i=1-OLx,sNx+OLx
8377b8ee87 Mart*0352 seaice_div (i,j) = 0. _d 0
0353 seaice_tension(i,j) = 0. _d 0
0354 seaice_shear (i,j) = 0. _d 0
0355 pressC (i,j) = 0. _d 0
0356 e12Csq (i,j) = 0. _d 0
0357 zetaC (i,j) = 0. _d 0
0358 deltaZ (i,j) = 0. _d 0
0359 zetaZ (i,j,bi,bj) = 0. _d 0
0360 deltaC (i,j,bi,bj) = 0. _d 0
b7030e3a29 Mart*0361 ENDDO
0362 ENDDO
9e706e0219 Jean*0363 DO j=1-OLy,sNy+OLy
0364 DO i=1-OLx,sNx+OLx
da4a30d88f Mart*0365 ep(i,j) = e11(i,j,bi,bj) + e22(i,j,bi,bj)
0366 em(i,j) = e11(i,j,bi,bj) - e22(i,j,bi,bj)
0367 ENDDO
0368 ENDDO
d2082aecba Mart*0369
0370
0371
0372 IF ( SEAICEetaZmethod .EQ. 0 ) THEN
0373 DO j=1-OLy+1,sNy+OLy-1
0374 DO i=1-OLx+1,sNx+OLx-1
0375 tmp = 0.25 *
8377b8ee87 Mart*0376 & ( e12(i,j ,bi,bj) + e12(i+1,j ,bi,bj)
0377 & + e12(i,j+1,bi,bj) + e12(i+1,j+1,bi,bj) )
d2082aecba Mart*0378 e12Csq(i,j) = tmp*tmp
0379 ENDDO
0380 ENDDO
0381 ELSEIF ( SEAICEetaZmethod .EQ. 3 ) THEN
0382 DO j=1-OLy+1,sNy+OLy-1
0383 DO i=1-OLx+1,sNx+OLx-1
0384
85f5225996 Mart*0385
8377b8ee87 Mart*0386 e12Csq(i,j) = 0.25 _d 0 * recip_rA(i,j,bi,bj) *
0387 & ( rAz(i ,j ,bi,bj)*e12(i ,j ,bi,bj)**2
0388 & + rAz(i+1,j ,bi,bj)*e12(i+1,j ,bi,bj)**2
0389 & + rAz(i ,j+1,bi,bj)*e12(i ,j+1,bi,bj)**2
0390 & + rAz(i+1,j+1,bi,bj)*e12(i+1,j+1,bi,bj)**2 )
d2082aecba Mart*0391 ENDDO
0392 ENDDO
0393 ENDIF
68ff576152 Mart*0394 DO j=0,sNy+1
0395 DO i=0,sNx+1
8377b8ee87 Mart*0396 deltaSq = ep(i,j)**2 + recip_ecc2 * em(i,j)**2
0397 & + recip_ecc2 * 4. _d 0 * e12Csq(i,j)
0398 #ifdef ALLOW_AUTODIFF
d2082aecba Mart*0399
8377b8ee87 Mart*0400 deltaC(i,j,bi,bj) = 0. _d 0
d2082aecba Mart*0401 IF ( deltaSq .GT. 0. _d 0 )
8377b8ee87 Mart*0402 & deltaC(i,j,bi,bj) = SQRT(deltaSq)
d2082aecba Mart*0403 #else
8377b8ee87 Mart*0404 deltaC(i,j,bi,bj) = SQRT(deltaSq)
0405 #endif /* ALLOW_AUTODIFF */
d2082aecba Mart*0406 #ifdef SEAICE_DELTA_SMOOTHREG
0407
0408
85f5225996 Mart*0409
8377b8ee87 Mart*0410 deltaCreg = deltaC(i,j,bi,bj) + SEAICE_deltaMin
d2082aecba Mart*0411 #else
8377b8ee87 Mart*0412 deltaCreg = MAX(deltaC(i,j,bi,bj),SEAICE_deltaMin)
d2082aecba Mart*0413 #endif /* SEAICE_DELTA_SMOOTHREG */
8377b8ee87 Mart*0414 zetaC(i,j) = HALF*( press0(i,j,bi,bj)
0415 & * ( 1. _d 0 + tensileStrFac(i,j,bi,bj) )
45c6bb8414 Mart*0416 & )/deltaCreg
d2082aecba Mart*0417 ENDDO
0418 ENDDO
a6389affa6 Mart*0419 IF ( useAdaptiveEVP ) THEN
0420 DO j=0,sNy+1
0421 DO i=0,sNx+1
e4120ef5ff Jean*0422
85f5225996 Mart*0423
af61e5eb16 Mart*0424 #ifdef ALLOW_AUTODIFF
0425 evpAlphaC(i,j,bi,bj) = SEAICE_evpAlpha
0426 IF ( zetaC(i,j) .GT. 0. _d 0 ) THEN
0427 #endif
8377b8ee87 Mart*0428 evpAlphaC(i,j,bi,bj) = SQRT(zetaC(i,j)
0429 & * EVPcFac / MAX(seaiceMassC(i,j,bi,bj), 1.D-04)
0430 & * recip_rA(i,j,bi,bj) ) * HEFFM(i,j,bi,bj)
af61e5eb16 Mart*0431 #ifdef ALLOW_AUTODIFF
0432 ENDIF
0433 #endif
8377b8ee87 Mart*0434 evpAlphaC(i,j,bi,bj) =
0435 & MAX(evpAlphaC(i,j,bi,bj),SEAICEaEVPalphaMin)
a6389affa6 Mart*0436 ENDDO
85f5225996 Mart*0437 ENDDO
a6389affa6 Mart*0438 ENDIF
e4120ef5ff Jean*0439
d2082aecba Mart*0440
8377b8ee87 Mart*0441 DO j=1,sNy+1
0442 DO i=1,sNx+1
0443 sumNorm = HEFFM(i,j, bi,bj)+HEFFM(i-1,j, bi,bj)
0444 & + HEFFM(i,j-1,bi,bj)+HEFFM(i-1,j-1,bi,bj)
967b0b9f87 Mart*0445 IF ( sumNorm.GT.0. _d 0 ) sumNorm = 1. _d 0 / sumNorm
8377b8ee87 Mart*0446 zetaZ(i,j,bi,bj) = sumNorm *
0447 & ( zetaC(i, j) + zetaC(i-1,j-1)
0448 & + zetaC(i-1,j) + zetaC(i, j-1) )
0449 deltaZ(i,j) = sumNorm *
0450 & ( deltaC(i, j,bi,bj) + deltaC(i-1,j-1,bi,bj)
0451 & + deltaC(i-1,j,bi,bj) + deltaC(i, j-1,bi,bj) )
d2082aecba Mart*0452 ENDDO
0453 ENDDO
fc27138d73 Mart*0454 #ifdef SEAICE_ALLOW_CLIPZETA
4619635089 Mart*0455 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0456
0457
4619635089 Mart*0458 #endif /* ALLOW_AUTODIFF_TAMC */
d8ae85c68c Mart*0459
0460 DO j=0,sNy+1
0461 DO i=0,sNx+1
8e32c48b8f Mart*0462 zetaC(i,j) = MAX( SEAICE_zMin(i,j,bi,bj),
0463 & MIN( SEAICE_zMax(i,j,bi,bj), zetaC(i,j) ) )
ec0d7df165 Mart*0464
d8ae85c68c Mart*0465
282fbc7b4c Mart*0466
0467
1cd6ce9d27 Mart*0468 zMaxZ = MAX(
8e32c48b8f Mart*0469 & MAX(SEAICE_zMax(i, j,bi,bj),SEAICE_zMax(i, j-1,bi,bj)),
0470 & MAX(SEAICE_zMax(i-1,j,bi,bj),SEAICE_zMax(i-1,j-1,bi,bj)) )
1cd6ce9d27 Mart*0471 zMinZ = MAX(
8e32c48b8f Mart*0472 & MAX(SEAICE_zMin(i, j,bi,bj),SEAICE_zMin(i, j-1,bi,bj)),
0473 & MAX(SEAICE_zMin(i-1,j,bi,bj),SEAICE_zMin(i-1,j-1,bi,bj)) )
8377b8ee87 Mart*0474 zetaZ(i,j,bi,bj) = MAX(zMinZ,MIN(zMaxZ,zetaZ(i,j,bi,bj)))
1cd6ce9d27 Mart*0475 ENDDO
0476 ENDDO
fc27138d73 Mart*0477 #endif /* SEAICE_ALLOW_CLIPZETA */
d8ae85c68c Mart*0478
0479 DO j=0,sNy+1
0480 DO i=0,sNx+1
8377b8ee87 Mart*0481 pressC(i,j) =
0482 & ( press0(i,j,bi,bj) * ( 1. _d 0 - SEAICEpressReplFac )
0483 & + TWO*zetaC(i,j)*deltaC(i,j,bi,bj)*SEAICEpressReplFac
0484 & /(1. _d 0 + tensileStrFac(i,j,bi,bj))
0485 & ) * (1. _d 0 - tensileStrFac(i,j,bi,bj))
d8ae85c68c Mart*0486 ENDDO
0487 ENDDO
0beea9c759 Mart*0488 #ifdef ALLOW_DIAGNOSTICS
1cd6ce9d27 Mart*0489 IF ( useDiagnostics ) THEN
a3f6dcaab3 Mart*0490
1cd6ce9d27 Mart*0491 DO j=1,sNy
0492 DO i=1,sNx
8377b8ee87 Mart*0493 press(i,j,bi,bj) = pressC(i,j)
0494 zeta (i,j,bi,bj) = zetaC(i,j)
0495 eta (i,j,bi,bj) = zetaC(i,j)*recip_ecc2
1cd6ce9d27 Mart*0496 ENDDO
0497 ENDDO
0498 ENDIF
0beea9c759 Mart*0499 #endif /* ALLOW_DIAGNOSTICS */
b7030e3a29 Mart*0500
0501
0502
0503
0504
dadd13178c Mart*0505 #ifdef SEAICE_ALLOW_TEM
1cd6ce9d27 Mart*0506 IF ( SEAICEuseTEM ) THEN
d2082aecba Mart*0507 DO j=0,sNy
0508 DO i=0,sNx
8377b8ee87 Mart*0509 etaDenC = em(i,j)**2 + 4. _d 0 * e12Csq(i,j)
d2082aecba Mart*0510 etaDenC = SQRT(MAX(deltaMinSq,etaDenC))
8377b8ee87 Mart*0511 zetaMaxC = ecc2*zetaC(i,j)
0512 & *(deltaC(i,j,bi,bj)-ep(i,j))/etaDenC
d2082aecba Mart*0513 #ifdef ALLOW_DIAGNOSTICS
0514
8377b8ee87 Mart*0515 eta(i,j,bi,bj) = MIN(zetaC(i,j),zetaMaxC)*recip_ecc2
d2082aecba Mart*0516 #endif /* ALLOW_DIAGNOSTICS */
8377b8ee87 Mart*0517 seaice_div (i,j) =
0518 & ( 2. _d 0 *zetaC(i,j)*ep(i,j) - pressC(i,j)
0519 & ) * HEFFM(i,j,bi,bj)
0520 seaice_tension(i,j) = 2. _d 0*MIN(zetaC(i,j),zetaMaxC)
0521 & * em(i,j) * HEFFM(i,j,bi,bj)
d2082aecba Mart*0522 ENDDO
0523 ENDDO
0524 DO j=1,sNy+1
0525 DO i=1,sNx+1
967b0b9f87 Mart*0526 sumNorm = 0.25 _d 0
0527
ec0d7df165 Mart*0528
0529
0530
0531
da4a30d88f Mart*0532
0533
c9a6744012 Jean*0534 etaDenZ =
8377b8ee87 Mart*0535 & sumNorm * recip_rAz(i,j,bi,bj) *
0536 & ( _rA(i ,j ,bi,bj) * em(i, j )**2
0537 & + _rA(i-1,j-1,bi,bj) * em(i-1,j-1)**2
0538 & + _rA(i-1,j ,bi,bj) * em(i-1,j )**2
0539 & + _rA(i ,j-1,bi,bj) * em(i, j-1)**2 )
0540 & + 4. _d 0*e12(i,j,bi,bj)**2
3daf25222c Mart*0541 etaDenZ = SQRT(MAX(deltaMinSq,etaDenZ))
8377b8ee87 Mart*0542 zetaMaxZ = ecc2*zetaZ(i,j,bi,bj) * ( deltaZ(i,j)
0543 & - sumNorm * ( ep(i,j ) + ep(i-1,j )
0544 & + ep(i,j-1) + ep(i-1,j-1) )
da4a30d88f Mart*0545 & )/etaDenZ
8377b8ee87 Mart*0546 seaice_shear (i,j) =
0547 & 2. _d 0*MIN(zetaZ(i,j,bi,bj),zetaMaxZ)
0548 & * 2. _d 0*e12(i,j,bi,bj)
1cd6ce9d27 Mart*0549 ENDDO
0550 ENDDO
0551 ELSE
0552 #else
0553 IF (.TRUE. ) THEN
dadd13178c Mart*0554 #endif /* SEAICE_ALLOW_TEM */
8377b8ee87 Mart*0555 DO j=0,sNy
0556 DO i=0,sNx
0557 seaice_div (i,j) =
0558 & ( 2. _d 0 *zetaC(i,j)*ep(i,j) - pressC(i,j)
0559 & ) * HEFFM(i,j,bi,bj)
0560 seaice_tension(i,j) = 2. _d 0*zetaC(i,j)
0561 & * em(i,j) * HEFFM(i,j,bi,bj)
d2082aecba Mart*0562 ENDDO
0563 ENDDO
8377b8ee87 Mart*0564 DO j=1,sNy+1
0565 DO i=1,sNx+1
0566 seaice_shear (i,j) =
0567 & 2. _d 0*zetaZ(i,j,bi,bj)*e12(i,j,bi,bj)
1cd6ce9d27 Mart*0568 ENDDO
037351a1a6 Mart*0569 ENDDO
1cd6ce9d27 Mart*0570 ENDIF
af61e5eb16 Mart*0571 #ifdef ALLOW_AUTODIFF_TAMC
0572
0573
0574
0575
0576 #endif /* ALLOW_AUTODIFF_TAMC */
1fb8ac1204 Mart*0577
b7030e3a29 Mart*0578
1fb8ac1204 Mart*0579
a6389affa6 Mart*0580 IF ( useAdaptiveEVP ) THEN
0581 DO j=0,sNy
0582 DO i=0,sNx
8377b8ee87 Mart*0583 denom1(i,j,bi,bj) = 1. _d 0 / evpAlphaC(i,j,bi,bj)
0584 denom2(i,j,bi,bj) = denom1(i,j,bi,bj)
a6389affa6 Mart*0585 ENDDO
0586 ENDDO
0587 ENDIF
af61e5eb16 Mart*0588 #ifdef ALLOW_AUTODIFF_TAMC
0589
0590
0591 #endif /* ALLOW_AUTODIFF_TAMC */
1fb8ac1204 Mart*0592 DO j=0,sNy
0593 DO i=0,sNx
b7030e3a29 Mart*0594
8377b8ee87 Mart*0595 seaice_sigma1 (i,j,bi,bj) = ( seaice_sigma1 (i,j,bi,bj)
0596 & * ( evpAlphaC(i,j,bi,bj) - evpRevFac )
0597 & + seaice_div(i,j)
0598 & ) * denom1(i,j,bi,bj)
0599 & *HEFFM(i,j,bi,bj)
0600 seaice_sigma2 (i,j,bi,bj) = ( seaice_sigma2 (i,j,bi,bj)
0601 & * ( evpAlphaC(i,j,bi,bj) - evpRevFac )
0602 & + seaice_tension(i,j)*recip_evpRevFac
0603 & ) * denom2(i,j,bi,bj)
0604 & *HEFFM(i,j,bi,bj)
dc7b968f4a Mart*0605 #ifdef SEAICE_EVP_ELIMINATE_UNDERFLOWS
0606
9e706e0219 Jean*0607
dc7b968f4a Mart*0608
8377b8ee87 Mart*0609 seaice_sigma1(i,j,bi,bj) = SIGN(MAX(
0610 & ABS( seaice_sigma1(i,j,bi,bj) ), SEAICE_EPS ),
0611 & seaice_sigma1(i,j,bi,bj) )
0612 seaice_sigma2(i,j,bi,bj) = SIGN(MAX(
0613 & ABS( seaice_sigma2(i,j,bi,bj) ), SEAICE_EPS ),
0614 & seaice_sigma2(i,j,bi,bj) )
dc7b968f4a Mart*0615 #endif /* SEAICE_EVP_ELIMINATE_UNDERFLOWS */
b7030e3a29 Mart*0616
8377b8ee87 Mart*0617 sig11(i,j) = 0.5 _d 0 *
0618 & ( seaice_sigma1(i,j,bi,bj)+seaice_sigma2(i,j,bi,bj) )
0619 sig22(i,j) = 0.5 _d 0 *
0620 & ( seaice_sigma1(i,j,bi,bj)-seaice_sigma2(i,j,bi,bj) )
b7030e3a29 Mart*0621 ENDDO
1fb8ac1204 Mart*0622 ENDDO
0623
0624
a6389affa6 Mart*0625 IF ( useAdaptiveEVP ) THEN
0626 DO j=1,sNy+1
0627 DO i=1,sNx+1
8377b8ee87 Mart*0628 evpAlphaZ(i,j,bi,bj) = 0.25 _d 0 *
0629 & ( evpAlphaC(i, j,bi,bj)+evpAlphaC(i-1,j-1,bi,bj)
0630 & + evpAlphaC(i-1,j,bi,bj)+evpAlphaC(i, j-1,bi,bj) )
0631 denom2(i,j,bi,bj) = 1. _d 0 / evpAlphaZ(i,j,bi,bj)
a6389affa6 Mart*0632 ENDDO
0633 ENDDO
0634 ENDIF
af61e5eb16 Mart*0635 #ifdef ALLOW_AUTODIFF_TAMC
0636
0637
0638 #endif /* ALLOW_AUTODIFF_TAMC */
1fb8ac1204 Mart*0639 DO j=1,sNy+1
0640 DO i=1,sNx+1
8377b8ee87 Mart*0641 seaice_sigma12(i,j,bi,bj) = ( seaice_sigma12(i,j,bi,bj)
0642 & * ( evpAlphaZ(i,j,bi,bj) - evpRevFac )
0643 & + seaice_shear(i,j)*recip_evpRevFac
0644 & ) * denom2(i,j,bi,bj)
dc7b968f4a Mart*0645 #ifdef SEAICE_EVP_ELIMINATE_UNDERFLOWS
8377b8ee87 Mart*0646 seaice_sigma12(i,j,bi,bj) = SIGN(MAX(
0647 & ABS( seaice_sigma12(i,j,bi,bj) ), SEAICE_EPS ),
0648 & seaice_sigma12(i,j,bi,bj) )
dc7b968f4a Mart*0649 #endif /* SEAICE_EVP_ELIMINATE_UNDERFLOWS */
1fb8ac1204 Mart*0650 ENDDO
b4949dd6db Jean*0651 ENDDO
b7030e3a29 Mart*0652
0653
0654
8377b8ee87 Mart*0655 DO j=1,sNy
0656 DO i=1,sNx
0657 stressDivergenceX(i,j,bi,bj) =
0658 & ( sig11(i ,j ) * _dyF(i ,j,bi,bj)
0659 & - sig11(i-1,j ) * _dyF(i-1,j,bi,bj)
0660 & + seaice_sigma12(i,j+1,bi,bj) * _dxV(i,j+1,bi,bj)
0661 & - seaice_sigma12(i,j ,bi,bj) * _dxV(i,j ,bi,bj)
0662 & ) * recip_rAw(i,j,bi,bj)
0663 stressDivergenceY(i,j,bi,bj) =
0664 & ( sig22(i,j ) * _dxF(i,j ,bi,bj)
0665 & - sig22(i,j-1) * _dxF(i,j-1,bi,bj)
0666 & + seaice_sigma12(i+1,j,bi,bj) * _dyU(i+1,j,bi,bj)
0667 & - seaice_sigma12(i ,j,bi,bj) * _dyU(i ,j,bi,bj)
0668 & ) * recip_rAs(i,j,bi,bj)
b7030e3a29 Mart*0669 ENDDO
b4949dd6db Jean*0670 ENDDO
4669aaa4e4 Mart*0671 ENDDO
0672 ENDDO
0673
e4120ef5ff Jean*0674 #ifdef ALLOW_SEAICE_EVP_RESIDUAL
4669aaa4e4 Mart*0675 IF ( printResidual ) THEN
0676 DO bj=myByLo(myThid),myByHi(myThid)
0677 DO bi=myBxLo(myThid),myBxHi(myThid)
85f5225996 Mart*0678 resTile(bi,bj) = 0. _d 0
0679 DO j=1,sNy
0680 DO i=1,sNx
e4120ef5ff Jean*0681 sig11Pm1(i,j,bi,bj) =
4669aaa4e4 Mart*0682 & seaice_sigma1(i,j,bi,bj)-sig11pm1(i,j,bi,bj)
e4120ef5ff Jean*0683 sig22Pm1(i,j,bi,bj) =
4669aaa4e4 Mart*0684 & seaice_sigma2(i,j,bi,bj)-sig22pm1(i,j,bi,bj)
e4120ef5ff Jean*0685 sig12Pm1(i,j,bi,bj) =
4669aaa4e4 Mart*0686 & seaice_sigma12(i,j,bi,bj)-sig12Pm1(i,j,bi,bj)
e4120ef5ff Jean*0687 sig11Pm1(i,j,bi,bj) =
8377b8ee87 Mart*0688 & evpAlphaC(i,j,bi,bj) * sig11Pm1(i,j,bi,bj)
e4120ef5ff Jean*0689 sig22Pm1(i,j,bi,bj) =
8377b8ee87 Mart*0690 & evpAlphaC(i,j,bi,bj) * sig22Pm1(i,j,bi,bj)
e4120ef5ff Jean*0691 sig12Pm1(i,j,bi,bj) =
8377b8ee87 Mart*0692 & evpAlphaZ(i,j,bi,bj) * sig12Pm1(i,j,bi,bj)
85f5225996 Mart*0693 ENDDO
0694 ENDDO
859b587af3 Mart*0695 IF ( .NOT. SEAICEscaleSurfStress ) THEN
e4120ef5ff Jean*0696
859b587af3 Mart*0697 DO j=1,sNy
0698 DO i=1,sNx
0699 resTile(bi,bj) = resTile(bi,bj) + AREA(i,j,bi,bj) *
0700 & ( sig11Pm1(i,j,bi,bj)*sig11Pm1(i,j,bi,bj)
0701 & + sig22Pm1(i,j,bi,bj)*sig22Pm1(i,j,bi,bj)
0702 & + sig12Pm1(i,j,bi,bj)*sig12Pm1(i,j,bi,bj) )
0703 ENDDO
0704 ENDDO
0705 ELSE
0706
0707 DO j=1,sNy
0708 DO i=1,sNx
0709 resTile(bi,bj) = resTile(bi,bj)
0710 & + sig11Pm1(i,j,bi,bj)*sig11Pm1(i,j,bi,bj)
0711 & + sig22Pm1(i,j,bi,bj)*sig22Pm1(i,j,bi,bj)
0712 & + sig12Pm1(i,j,bi,bj)*sig12Pm1(i,j,bi,bj)
0713 ENDDO
0714 ENDDO
0715 ENDIF
4669aaa4e4 Mart*0716 ENDDO
88601fb2cb Gael*0717 ENDDO
85f5225996 Mart*0718 resloc = 0. _d 0
0719 CALL GLOBAL_SUM_TILE_RL( resTile, resloc, myThid )
0720 resloc = SQRT(resloc)
0721 WRITE(standardMessageUnit,'(A,1X,I6,1PE16.8)')
a6389affa6 Mart*0722 & ' SEAICE_EVP: iEVPstep, residual sigma = ', iEVPstep, resLoc
85f5225996 Mart*0723 ENDIF
0724 #endif /* ALLOW_SEAICE_EVP_RESIDUAL */
88601fb2cb Gael*0725 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0726
0727
a6389affa6 Mart*0728 # ifdef SEAICE_DYN_STABLE_ADJOINT
c9a6744012 Jean*0729
4f95e6bec9 Gael*0730 CALL ZERO_ADJ( 1, stressDivergenceX, myThid)
0731 CALL ZERO_ADJ( 1, stressDivergenceY, myThid)
a6389affa6 Mart*0732 # endif /* SEAICE_DYN_STABLE_ADJOINT */
4f95e6bec9 Gael*0733 #endif /* ALLOW_AUTODIFF_TAMC */
0734
6e2f4e58fa Mart*0735
037351a1a6 Mart*0736
6e2f4e58fa Mart*0737
3f31c7d5de Mart*0738 CALL SEAICE_OCEANDRAG_COEFFS(
ec0d7df165 Mart*0739 I uIce, vIce, HEFFM,
3f31c7d5de Mart*0740 O DWATN,
0741 I iEVPstep, myTime, myIter, myThid )
df1dac8b7b Mart*0742 #ifdef SEAICE_ALLOW_BOTTOMDRAG
d5254b4e3d Mart*0743 CALL SEAICE_BOTTOMDRAG_COEFFS(
ec0d7df165 Mart*0744 I uIce, vIce, HEFFM,
df1dac8b7b Mart*0745 #ifdef SEAICE_ITD
0746 I HEFFITD, AREAITD, AREA,
0747 #else
0748 I HEFF, AREA,
e4120ef5ff Jean*0749 #endif
df1dac8b7b Mart*0750 O CbotC,
0751 I iEVPstep, myTime, myIter, myThid )
0752 #endif /* SEAICE_ALLOW_BOTTOMDRAG */
5bb179ddc2 Mart*0753 #ifdef SEAICE_ALLOW_SIDEDRAG
0754 IF ( SEAICEsideDrag .NE. 0. _d 0 ) CALL SEAICE_SIDEDRAG_STRESS(
0755 I uIce, vIce, coastRoughU, coastRoughV, AREA,
0756 O sideDragU, sideDragV,
0757 I iEVPstep, myTime, myIter, myThid )
0758 #endif /* SEAICE_ALLOW_SIDEDRAG */
0759
af61e5eb16 Mart*0760 #ifdef ALLOW_AUTODIFF_TAMC
0761
0762 # ifdef SEAICE_ALLOW_BOTTOMDRAG
0763
0764 # endif
5bb179ddc2 Mart*0765 # ifdef SEAICE_ALLOW_SIDEDRAG
0766
0767
0768 # endif
af61e5eb16 Mart*0769 #endif
88601fb2cb Gael*0770 DO bj=myByLo(myThid),myByHi(myThid)
0771 DO bi=myBxLo(myThid),myBxHi(myThid)
af61e5eb16 Mart*0772 #ifdef ALLOW_AUTODIFF_TAMC
0773 tkey = bi + (bj-1)*nSx + (evpkey-1)*nSx*nSy
0774
0775
0776 #endif /* ALLOW_AUTODIFF_TAMC */
8377b8ee87 Mart*0777 DO j=1,sNy
0778 DO i=1,sNx
e8b26f9f39 Mart*0779
0780
e4120ef5ff Jean*0781
0782
e8b26f9f39 Mart*0783
0784
b8e39f4242 Mart*0785
e8b26f9f39 Mart*0786
8377b8ee87 Mart*0787 locMaskU = seaiceMassU(i,j,bi,bj)
0788 locMaskV = seaiceMassV(i,j,bi,bj)
e8b26f9f39 Mart*0789 IF ( locMaskU .NE. 0. _d 0 ) locMaskU = 1. _d 0
0790 IF ( locMaskV .NE. 0. _d 0 ) locMaskV = 1. _d 0
0791
0792
0793
037351a1a6 Mart*0794
6e2f4e58fa Mart*0795
8377b8ee87 Mart*0796 FORCEX(i,j,bi,bj)=FORCEX0(i,j,bi,bj)+
0797 & ( 0.5 _d 0 * ( DWATN(i,j,bi,bj)+DWATN(i-1,j,bi,bj) ) *
0798 & COSWAT * uVel(i,j,kSrf,bi,bj)
0799 & - SIGN(SINWAT, _fCori(i,j,bi,bj))* 0.5 _d 0 *
0800 & ( DWATN(i ,j,bi,bj) * 0.5 _d 0 *
0801 & (vVel(i ,j ,kSrf,bi,bj)-vIce(i ,j ,bi,bj)
0802 & +vVel(i ,j+1,kSrf,bi,bj)-vIce(i ,j+1,bi,bj))
0803 & + DWATN(i-1,j,bi,bj) * 0.5 _d 0 *
0804 & (vVel(i-1,j ,kSrf,bi,bj)-vIce(i-1,j ,bi,bj)
0805 & +vVel(i-1,j+1,kSrf,bi,bj)-vIce(i-1,j+1,bi,bj))
0806 & )*locMaskU ) * areaW(i,j,bi,bj)
0807 FORCEY(i,j,bi,bj)=FORCEY0(i,j,bi,bj)+
0808 & ( 0.5 _d 0 * ( DWATN(i,j,bi,bj)+DWATN(i,j-1,bi,bj) ) *
0809 & COSWAT * vVel(i,j,kSrf,bi,bj)
0810 & + SIGN(SINWAT, _fCori(i,j,bi,bj)) * 0.5 _d 0 *
0811 & ( DWATN(i,j ,bi,bj) * 0.5 _d 0 *
0812 & (uVel(i ,j ,kSrf,bi,bj)-uIce(i ,j ,bi,bj)
0813 & +uVel(i+1,j ,kSrf,bi,bj)-uIce(i+1,j ,bi,bj))
0814 & + DWATN(i,j-1,bi,bj) * 0.5 _d 0 *
0815 & (uVel(i ,j-1,kSrf,bi,bj)-uIce(i ,j-1,bi,bj)
0816 & +uVel(i+1,j-1,kSrf,bi,bj)-uIce(i+1,j-1,bi,bj))
0817 & )*locMaskV ) * areaS(i,j,bi,bj)
037351a1a6 Mart*0818
8377b8ee87 Mart*0819 FORCEX(i,j,bi,bj)=FORCEX(i,j,bi,bj) + 0.5 _d 0*(
0820 & seaiceMassC(i ,j,bi,bj) * _fCori(i ,j,bi,bj)
0821 & * 0.5 _d 0*( vIce(i ,j,bi,bj)+vIce(i ,j+1,bi,bj) )
0822 & + seaiceMassC(i-1,j,bi,bj) * _fCori(i-1,j,bi,bj)
0823 & * 0.5 _d 0*( vIce(i-1,j,bi,bj)+vIce(i-1,j+1,bi,bj) )
b7030e3a29 Mart*0824 & )
8377b8ee87 Mart*0825 FORCEY(i,j,bi,bj)=FORCEY(i,j,bi,bj) - 0.5 _d 0*(
0826 & seaiceMassC(i,j ,bi,bj) * _fCori(i,j ,bi,bj)
0827 & * 0.5 _d 0*( uIce(i,j ,bi,bj)+uIce(i+1, j,bi,bj) )
0828 & + seaiceMassC(i,j-1,bi,bj) * _fCori(i,j-1,bi,bj)
0829 & * 0.5 _d 0*( uIce(i,j-1,bi,bj)+uIce(i+1,j-1,bi,bj) )
b7030e3a29 Mart*0830 & )
452bde753c Mart*0831 ENDDO
0832 ENDDO
af61e5eb16 Mart*0833 #ifdef ALLOW_AUTODIFF_TAMC
0834
0835
0836 #endif /* ALLOW_AUTODIFF_TAMC */
8ce59fd29e Mart*0837 #ifdef SEAICE_ALLOW_MOM_ADVECTION
ec0d7df165 Mart*0838 IF ( SEAICEmomAdvection ) THEN
8377b8ee87 Mart*0839 DO j=1-OLy,sNy+OLy
0840 DO i=1-OLx,sNx+OLx
0841 gUmom(i,j) = 0. _d 0
0842 gVmom(i,j) = 0. _d 0
8ce59fd29e Mart*0843 ENDDO
0844 ENDDO
0845 CALL SEAICE_MOM_ADVECTION(
0846 I bi,bj,1,sNx,1,sNy,
0847 I uIce, vIce,
0848 O gUmom, gVmom,
0849 I myTime, myIter, myThid )
8377b8ee87 Mart*0850 DO j=1,sNy
0851 DO i=1,sNx
0852 FORCEX(i,j,bi,bj) = FORCEX(i,j,bi,bj) + gUmom(i,j)
0853 FORCEY(i,j,bi,bj) = FORCEY(i,j,bi,bj) + gVmom(i,j)
8ce59fd29e Mart*0854 ENDDO
0855 ENDDO
0856 ENDIF
0857 #endif /* SEAICE_ALLOW_MOM_ADVECTION */
6e2f4e58fa Mart*0858
0859
0860
a6389affa6 Mart*0861 IF ( useAdaptiveEVP ) THEN
8377b8ee87 Mart*0862 DO j=1,sNy
0863 DO i=1,sNx
85f5225996 Mart*0864
8377b8ee87 Mart*0865 evpBetaU(i,j,bi,bj) = 0.5 _d 0*(evpAlphaC(i-1,j,bi,bj)
0866 & +evpAlphaC(i, j,bi,bj))
0867 evpBetaV(i,j,bi,bj) = 0.5 _d 0*(evpAlphaC(i,j-1,bi,bj)
0868 & +evpAlphaC(i,j, bi,bj))
85f5225996 Mart*0869
a6389affa6 Mart*0870 ENDDO
0871 ENDDO
0872 ENDIF
af61e5eb16 Mart*0873 #ifdef ALLOW_AUTODIFF_TAMC
0874
0875
0876
0877
0878 #endif /* ALLOW_AUTODIFF_TAMC */
8377b8ee87 Mart*0879 DO j=1,sNy
0880 DO i=1,sNx
0881 betaFacU = evpBetaU(i,j,bi,bj)*recip_deltaT
0882 betaFacV = evpBetaV(i,j,bi,bj)*recip_deltaT
85f5225996 Mart*0883 tmp = evpStarFac*recip_deltaT
0884 betaFacP1V = betaFacV + tmp
0885 betaFacP1U = betaFacU + tmp
8377b8ee87 Mart*0886 denomU = seaiceMassU(i,j,bi,bj)*betaFacP1U
0887 & + 0.5 _d 0*( DWATN(i,j,bi,bj) + DWATN(i-1,j,bi,bj) )
0888 & * COSWAT * areaW(i,j,bi,bj)
0889 denomV = seaiceMassV(i,j,bi,bj)*betaFacP1V
0890 & + 0.5 _d 0*( DWATN(i,j,bi,bj) + DWATN(i,j-1,bi,bj) )
0891 & * COSWAT * areaS(i,j,bi,bj)
df1dac8b7b Mart*0892 #ifdef SEAICE_ALLOW_BOTTOMDRAG
8377b8ee87 Mart*0893 denomU = denomU + areaW(i,j,bi,bj)
0894 & * 0.5 _d 0*( CbotC(i,j,bi,bj) + CbotC(i-1,j,bi,bj) )
0895 denomV = denomV + areaS(i,j,bi,bj)
0896 & * 0.5 _d 0*( CbotC(i,j,bi,bj) + CbotC(i,j-1,bi,bj) )
df1dac8b7b Mart*0897 #endif /* SEAICE_ALLOW_BOTTOMDRAG */
5bb179ddc2 Mart*0898 #ifdef SEAICE_ALLOW_SIDEDRAG
0899 denomU = denomU + sideDragU(i,j,bi,bj)
0900 denomV = denomV + sideDragV(i,j,bi,bj)
0901 #endif
70e078b38a Mart*0902 IF ( denomU .EQ. 0. _d 0 ) denomU = 1. _d 0
0903 IF ( denomV .EQ. 0. _d 0 ) denomV = 1. _d 0
8377b8ee87 Mart*0904 uIce(i,j,bi,bj) = seaiceMaskU(i,j,bi,bj) *
0905 & ( seaiceMassU(i,j,bi,bj)*betaFacU
0906 & * uIce(i,j,bi,bj)
0907 & + seaiceMassU(i,j,bi,bj)*recip_deltaT*evpStarFac
0908 & * uIceNm1(i,j,bi,bj)
0909 & + FORCEX(i,j,bi,bj)
0910 & + stressDivergenceX(i,j,bi,bj) ) / denomU
0911 vIce(i,j,bi,bj) = seaiceMaskV(i,j,bi,bj) *
0912 & ( seaiceMassV(i,j,bi,bj)*betaFacV
0913 & * vIce(i,j,bi,bj)
0914 & + seaiceMassV(i,j,bi,bj)*recip_deltaT*evpStarFac
0915 & * vIceNm1(i,j,bi,bj)
0916 & + FORCEY(i,j,bi,bj)
0917 & + stressDivergenceY(i,j,bi,bj) ) / denomV
c9a6744012 Jean*0918
a3f6dcaab3 Mart*0919
0920
0921
0922
037351a1a6 Mart*0923 ENDDO
0924 ENDDO
fec75090d3 Jean*0925 #ifndef OBCS_UVICE_OLD
af61e5eb16 Mart*0926 #ifdef ALLOW_AUTODIFF_TAMC
0927
0928
0929 #endif /* ALLOW_AUTODIFF_TAMC */
fec75090d3 Jean*0930 DO j=1,sNy
0931 DO i=1,sNx
0932 locMaskU = maskInC(i,j,bi,bj)*maskInC(i-1,j,bi,bj)
0933 locMaskV = maskInC(i,j,bi,bj)*maskInC(i,j-1,bi,bj)
a3f6dcaab3 Mart*0934 uIce(i,j,bi,bj) = uIce(i,j,bi,bj)*locMaskU
0935 & + uIceNm1(i,j,bi,bj)*(ONE-locMaskU)
0936 vIce(i,j,bi,bj) = vIce(i,j,bi,bj)*locMaskV
0937 & + vIceNm1(i,j,bi,bj)*(ONE-locMaskV)
fec75090d3 Jean*0938 ENDDO
0939 ENDDO
0940 #endif /* OBCS_UVICE_OLD */
6e2f4e58fa Mart*0941 ENDDO
0942 ENDDO
b4949dd6db Jean*0943
a3f6dcaab3 Mart*0944 CALL EXCH_UV_XY_RL(uIce,vIce,.TRUE.,myThid)
b4949dd6db Jean*0945
e4120ef5ff Jean*0946 #ifdef ALLOW_SEAICE_EVP_RESIDUAL
85f5225996 Mart*0947 IF ( printResidual ) THEN
0948 DO bj=myByLo(myThid),myByHi(myThid)
0949 DO bi=myBxLo(myThid),myBxHi(myThid)
0950 resTile(bi,bj) = 0. _d 0
8377b8ee87 Mart*0951 DO j=1,sNy
0952 DO i=1,sNx
0953 uIcePm1(i,j,bi,bj) = seaiceMaskU(i,j,bi,bj) *
0954 & ( uIce(i,j,bi,bj)-uIcePm1(i,j,bi,bj) )
0955 & * evpBetaU(i,j,bi,bj)
0956 vIcePm1(i,j,bi,bj) = seaiceMaskV(i,j,bi,bj) *
0957 & ( vIce(i,j,bi,bj)-vIcePm1(i,j,bi,bj) )
0958 & * evpBetaV(i,j,bi,bj)
85f5225996 Mart*0959 ENDDO
0960 ENDDO
859b587af3 Mart*0961 IF ( .NOT. SEAICEscaleSurfStress ) THEN
e4120ef5ff Jean*0962
859b587af3 Mart*0963 DO j=1,sNy
0964 DO i=1,sNx
0965 resTile(bi,bj) = resTile(bi,bj) + AREA(i,j,bi,bj) *
8377b8ee87 Mart*0966 & ( uIcePm1(i,j,bi,bj)*uIcePm1(i,j,bi,bj)
0967 & + vIcePm1(i,j,bi,bj)*vIcePm1(i,j,bi,bj) )
859b587af3 Mart*0968 ENDDO
0969 ENDDO
0970 ELSE
0971 DO j=1,sNy
0972 DO i=1,sNx
0973 resTile(bi,bj) = resTile(bi,bj)
8377b8ee87 Mart*0974 & + uIcePm1(i,j,bi,bj)*uIcePm1(i,j,bi,bj)
0975 & + vIcePm1(i,j,bi,bj)*vIcePm1(i,j,bi,bj)
859b587af3 Mart*0976 ENDDO
0977 ENDDO
0978 ENDIF
85f5225996 Mart*0979 ENDDO
0980 ENDDO
0981 resloc = 0. _d 0
0982 CALL GLOBAL_SUM_TILE_RL( resTile, resloc, myThid )
0983 resloc = SQRT(resloc)
0984 WRITE(standardMessageUnit,'(A,1X,I6,1PE16.8)')
a6389affa6 Mart*0985 & ' SEAICE_EVP: iEVPstep, residual U = ', iEVPstep, resLoc
85f5225996 Mart*0986 ENDIF
0987
0988
0989
0990
0991
0992
0993
0994
0995
0996
0997
0998
0999
1000
1001 #endif /* ALLOW_SEAICE_EVP_RESIDUAL */
1002
cf9fa44a59 Patr*1003 ENDIF
6e2f4e58fa Mart*1004
b4949dd6db Jean*1005 ENDDO
6e2f4e58fa Mart*1006
45315406aa Mart*1007 #endif /* SEAICE_CGRID and SEAICE_ALLOW_EVP */
6e2f4e58fa Mart*1008
1009 RETURN
1010 END