File indexing completed on 2023-09-03 05:10:25 UTC
view on githubraw file Latest commit 74487008 on 2023-09-03 01:50:18 UTC
0c49347dc7 Alis*0001 #include "GMREDI_OPTIONS.h"
14e0496834 Jean*0002 #ifdef ALLOW_AUTODIFF
0003 # include "AUTODIFF_OPTIONS.h"
0004 #endif
0c49347dc7 Alis*0005
e2259a1942 Jean*0006
0007
0008
0c49347dc7 Alis*0009 SUBROUTINE GMREDI_SLOPE_LIMIT(
549d1a8d8c Jean*0010 O SlopeX, SlopeY,
e9fd580bc4 Jean*0011 O SlopeSqr, taperFct,
e2259a1942 Jean*0012 U hTransLay, baseSlope, recipLambda,
549d1a8d8c Jean*0013 U dSigmaDr,
0014 I dSigmaDx, dSigmaDy,
5b172de0d2 Jean*0015 I Lrho, hMixLay, rDepth, depthZ, kLow,
0016 I kPos, k, bi, bj, myTime, myIter, myThid )
e2259a1942 Jean*0017
0018
0019
0020
0021
0022
0023
5b172de0d2 Jean*0024
e2259a1942 Jean*0025
0026
0027
5b172de0d2 Jean*0028
0029
e2259a1942 Jean*0030
0031
0032
5b172de0d2 Jean*0033
e2259a1942 Jean*0034
0035
0036
0037
0038
0039
0040
0041
0042
0043
0c49347dc7 Alis*0044 IMPLICIT NONE
0045
0046
0047 #include "SIZE.h"
b6b11b9b2f Patr*0048 #include "GRID.h"
0c49347dc7 Alis*0049 #include "EEPARAMS.h"
0050 #include "GMREDI.h"
0051 #include "PARAMS.h"
e673dad091 Patr*0052 #ifdef ALLOW_AUTODIFF_TAMC
0053 #include "tamc.h"
0054 #endif /* ALLOW_AUTODIFF_TAMC */
0c49347dc7 Alis*0055
549d1a8d8c Jean*0056
e2259a1942 Jean*0057
0058
0059
0060
5b172de0d2 Jean*0061
e2259a1942 Jean*0062
0063
5b172de0d2 Jean*0064
e2259a1942 Jean*0065
0066
0067
0068
0069
549d1a8d8c Jean*0070
e2259a1942 Jean*0071
0072
0073
0074
0075
0076
0077
5b172de0d2 Jean*0078
14e0496834 Jean*0079 _RL SlopeX (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0080 _RL SlopeY (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0081 _RL SlopeSqr (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0082 _RL taperFct (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0083 _RL hTransLay (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0084 _RL baseSlope (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0085 _RL recipLambda(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0086 _RL dSigmaDr (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0087 _RL dSigmaDx (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0088 _RL dSigmaDy (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0089 _RL Lrho (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0090 _RL hMixLay (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
5b172de0d2 Jean*0091 _RL rDepth
e2259a1942 Jean*0092 _RS depthZ(*)
14e0496834 Jean*0093 INTEGER kLow (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
5b172de0d2 Jean*0094 INTEGER kPos, k, bi, bj
e2259a1942 Jean*0095 _RL myTime
0096 INTEGER myIter
0097 INTEGER myThid
0c49347dc7 Alis*0098
0099 #ifdef ALLOW_GMREDI
e2259a1942 Jean*0100
0c49347dc7 Alis*0101 _RL fpi
413f546150 Jean*0102 PARAMETER( fpi = PI )
9f5240b52a Jean*0103 INTEGER i, j
0104 _RL f1, Smod, f2, Rnondim
0105 _RL maxSlopeSqr
5b172de0d2 Jean*0106 _RL convSlopeUnit, loc_rMaxSlope
9f5240b52a Jean*0107 _RL dSigmMod (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0108 _RL dRdSigmaLtd(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0109 _RL tmpFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0110 #ifndef GM_EXCLUDE_FM07_TAP
25ef785e26 Jean*0111 _RL dTransLay, rLambMin, DoverLamb
ac0d9d59a8 Jean*0112 _RL taperFctLoc, taperFctHat
25ef785e26 Jean*0113 _RL minTransLay
9f5240b52a Jean*0114 _RL SlopeMod(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0115 _RL locVar (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0116 #endif
0117 #ifdef GMREDI_WITH_STABLE_ADJOINT
0118 _RL slopeSqTmp, slopeTmp, slopeMax
5a6997640b Mart*0119
0120
0121
7c50f07931 Mart*0122 #endif
df4facdcff Jean*0123
0124 _RL GM_bigSlope
0125 GM_bigSlope = 1. _d +02
5b172de0d2 Jean*0126
0127
0128 IF (kPos.EQ.3 ) THEN
0129 GM_bigSlope = GM_bigSlope*wUnit2rVel(k)
0130 maxSlopeSqr = GM_maxSlope*GM_maxSlope
0131 & *wUnit2rVel(k)*wUnit2rVel(k)
0132 convSlopeUnit = rVel2wUnit(k)
0133 ELSE
0134 GM_bigSlope = GM_bigSlope*z2rUnit(k)
0135 maxSlopeSqr = GM_maxSlope*GM_maxSlope
0136 & *z2rUnit(k)*z2rUnit(k)
0137 convSlopeUnit = rUnit2z(k)
0138 ENDIF
0139 loc_rMaxSlope = GM_rMaxSlope*convSlopeUnit
df4facdcff Jean*0140
e673dad091 Patr*0141 #ifdef ALLOW_AUTODIFF_TAMC
5b172de0d2 Jean*0142
9d44c9cacc Mart*0143
0144
5b172de0d2 Jean*0145
0146
113b21cd45 Mart*0147
0148
0149
0150
7448700841 Mart*0151
0152
0153
5a6997640b Mart*0154
0155
0156
0157
0158
0159
0160
0161
e673dad091 Patr*0162 #endif /* ALLOW_AUTODIFF_TAMC */
0163
14e0496834 Jean*0164 DO j=1-OLy+1,sNy+OLy-1
0165 DO i=1-OLx+1,sNx+OLx-1
e2259a1942 Jean*0166 dSigmMod(i,j) = 0. _d 0
ac0d9d59a8 Jean*0167 tmpFld(i,j) = 0. _d 0
b6b11b9b2f Patr*0168 ENDDO
0169 ENDDO
0170
e9fd580bc4 Jean*0171 IF (GM_taper_scheme.EQ.'orig' .OR.
0172 & GM_taper_scheme.EQ.'clipping') THEN
0c49347dc7 Alis*0173
bd3f833a36 Jean*0174 #ifdef GM_EXCLUDE_CLIPPING
0175
0176 STOP 'Need to compile without "#define GM_EXCLUDE_CLIPPING"'
0177
0178 #else /* GM_EXCLUDE_CLIPPING */
2092dbb101 Patr*0179
0c49347dc7 Alis*0180
0181
e9fd580bc4 Jean*0182
0183
14e0496834 Jean*0184 DO j=1-OLy+1,sNy+OLy-1
0185 DO i=1-OLx+1,sNx+OLx-1
ac0d9d59a8 Jean*0186 tmpFld(i,j) = dSigmaDx(i,j)*dSigmaDx(i,j)
e2259a1942 Jean*0187 & + dSigmaDy(i,j)*dSigmaDy(i,j)
5b172de0d2 Jean*0188 IF ( tmpFld(i,j) .EQ. zeroRL ) THEN
e2259a1942 Jean*0189 dSigmMod(i,j) = 0. _d 0
b6b11b9b2f Patr*0190 ELSE
ac0d9d59a8 Jean*0191 dSigmMod(i,j) = SQRT( tmpFld(i,j) )
b6b11b9b2f Patr*0192 ENDIF
e673dad091 Patr*0193 ENDDO
0194 ENDDO
0c49347dc7 Alis*0195
14e0496834 Jean*0196 DO j=1-OLy+1,sNy+OLy-1
0197 DO i=1-OLx+1,sNx+OLx-1
5b172de0d2 Jean*0198 IF ( dSigmMod(i,j) .NE. zeroRL ) THEN
0199 tmpFld(i,j) = dSigmMod(i,j)*loc_rMaxSlope
0200 IF ( dSigmaDr(i,j) .LE. tmpFld(i,j) )
ac0d9d59a8 Jean*0201 & dSigmaDr(i,j) = tmpFld(i,j)
b6b11b9b2f Patr*0202 ENDIF
e673dad091 Patr*0203 ENDDO
0204 ENDDO
0205
14e0496834 Jean*0206 DO j=1-OLy+1,sNy+OLy-1
0207 DO i=1-OLx+1,sNx+OLx-1
5b172de0d2 Jean*0208 IF ( dSigmMod(i,j) .EQ. zeroRL ) THEN
b6b11b9b2f Patr*0209 SlopeX(i,j) = 0. _d 0
0210 SlopeY(i,j) = 0. _d 0
0211 ELSE
e2259a1942 Jean*0212 dRdSigmaLtd(i,j) = 1. _d 0/( dSigmaDr(i,j) )
5b172de0d2 Jean*0213 SlopeX(i,j) = dSigmaDx(i,j)*dRdSigmaLtd(i,j)
0214 SlopeY(i,j) = dSigmaDy(i,j)*dRdSigmaLtd(i,j)
b6b11b9b2f Patr*0215 ENDIF
e673dad091 Patr*0216 ENDDO
0217 ENDDO
0218
14e0496834 Jean*0219 DO j=1-OLy+1,sNy+OLy-1
0220 DO i=1-OLx+1,sNx+OLx-1
8233d0ceb9 Jean*0221 SlopeSqr(i,j) = SlopeX(i,j)*SlopeX(i,j)
0222 & + SlopeY(i,j)*SlopeY(i,j)
0223 taperFct(i,j) = 1. _d 0
0c49347dc7 Alis*0224 ENDDO
0225 ENDDO
0226
bd3f833a36 Jean*0227 #endif /* GM_EXCLUDE_CLIPPING */
2092dbb101 Patr*0228
e2259a1942 Jean*0229 ELSEIF (GM_taper_scheme.EQ.'fm07' ) THEN
0230
0231
5755dbe390 Patr*0232 #ifdef GM_EXCLUDE_FM07_TAP
0233
0234 STOP 'Need to compile without "#define GM_EXCLUDE_FM07_TAP"'
0235
0236 #else /* GM_EXCLUDE_FM07_TAP */
0237
e2259a1942 Jean*0238
0239
14e0496834 Jean*0240 DO j=1-OLy+1,sNy+OLy-1
0241 DO i=1-OLx+1,sNx+OLx-1
e2259a1942 Jean*0242 IF ( k.GT.kLow(i,j) ) THEN
0243
0244 SlopeX (i,j) = 0. _d 0
0245 SlopeY (i,j) = 0. _d 0
0246 SlopeMod(i,j) = 0. _d 0
0247 taperFct(i,j) = 0. _d 0
0248 ELSE
0249
5b172de0d2 Jean*0250 IF ( dSigmaDr(i,j).LE. GM_Small_Number )
0251 & dSigmaDr(i,j) = GM_Small_Number
ac0d9d59a8 Jean*0252 tmpFld(i,j) = dSigmaDx(i,j)*dSigmaDx(i,j)
0253 & + dSigmaDy(i,j)*dSigmaDy(i,j)
5b172de0d2 Jean*0254 IF ( tmpFld(i,j).GT.zeroRL ) THEN
ac0d9d59a8 Jean*0255 locVar(i,j) = SQRT( tmpFld(i,j) )
5755dbe390 Patr*0256 SlopeX (i,j) = dSigmaDx(i,j)/locVar(i,j)
0257 SlopeY (i,j) = dSigmaDy(i,j)/locVar(i,j)
5b172de0d2 Jean*0258 SlopeMod(i,j) = locVar(i,j)/dSigmaDr(i,j)
e2259a1942 Jean*0259 taperFct(i,j) = 1. _d 0
0260 ELSE
0261 SlopeX (i,j) = 0. _d 0
0262 SlopeY (i,j) = 0. _d 0
0263 SlopeMod(i,j) = 0. _d 0
0264 taperFct(i,j) = 0. _d 0
0265 ENDIF
0266 ENDIF
0267 ENDDO
0268 ENDDO
0269
0270
ac0d9d59a8 Jean*0271 IF ( k.EQ.1 ) THEN
86e2f7bc1f Jean*0272 minTransLay = GM_facTrL2dz*( depthZ(k) - depthZ(k+1) )
ac0d9d59a8 Jean*0273 ELSE
86e2f7bc1f Jean*0274 minTransLay = GM_facTrL2dz*( depthZ(k-1) - depthZ(k) )
ac0d9d59a8 Jean*0275 ENDIF
14e0496834 Jean*0276 DO j=1-OLy+1,sNy+OLy-1
0277 DO i=1-OLx+1,sNx+OLx-1
e2259a1942 Jean*0278 IF ( hTransLay(i,j).LE.0. _d 0 ) THEN
0279
ac0d9d59a8 Jean*0280 tmpFld(i,j) = Lrho(i,j)*SlopeMod(i,j)
0281
0282
e2259a1942 Jean*0283 dTransLay =
ac0d9d59a8 Jean*0284 & MIN( MAX( tmpFld(i,j), minTransLay ),
86e2f7bc1f Jean*0285 & MAX( GM_facTrL2ML*hMixLay(i,j), GM_maxTransLay ) )
e2259a1942 Jean*0286 IF ( k.GE.kLow(i,j) ) THEN
0287
0288 recipLambda(i,j) = 0. _d 0
ac0d9d59a8 Jean*0289 baseSlope(i,j) = SlopeMod(i,j)
e2259a1942 Jean*0290
0291
0292
5b172de0d2 Jean*0293 ELSEIF ( dTransLay+hMixLay(i,j)+depthZ(k) .GE. zeroRL ) THEN
e2259a1942 Jean*0294
0295 hTransLay(i,j) = -depthZ(k+1)
0296
5b172de0d2 Jean*0297 IF ( baseSlope(i,j).GT.zeroRL ) THEN
ac0d9d59a8 Jean*0298 recipLambda(i,j) = recipLambda(i,j)
0299 & / MIN( baseSlope(i,j), GM_maxSlope )
e2259a1942 Jean*0300 ELSE
0301 recipLambda(i,j) = 0. _d 0
0302 ENDIF
25ef785e26 Jean*0303
5b172de0d2 Jean*0304 IF ( hMixLay(i,j)+depthZ(k+1).LT.zeroRL ) THEN
25ef785e26 Jean*0305 rLambMin = 1. _d 0 /( hMixLay(i,j)+depthZ(k+1) )
0306 recipLambda(i,j) = MAX( recipLambda(i,j), rLambMin )
0307 ENDIF
e2259a1942 Jean*0308 ELSE
0309
ac0d9d59a8 Jean*0310 recipLambda(i,j) = ( MIN( SlopeMod(i,j), GM_maxSlope )
0311 & - MIN( baseSlope(i,j), GM_maxSlope )
0312 & ) / ( depthZ(k) - depthZ(k+1) )
e2259a1942 Jean*0313 baseSlope(i,j) = SlopeMod(i,j)
0314 ENDIF
0315 ENDIF
0316 ENDDO
0317 ENDDO
0318
0319
0320
14e0496834 Jean*0321 DO j=1-OLy+1,sNy+OLy-1
0322 DO i=1-OLx+1,sNx+OLx-1
e2259a1942 Jean*0323 IF ( hTransLay(i,j).GT.0. _d 0 ) THEN
0324
0325
ac0d9d59a8 Jean*0326 DoverLamb = (hTransLay(i,j)-hMixLay(i,j))*recipLambda(i,j)
e2259a1942 Jean*0327 IF ( -depthZ(k).LE.hMixLay(i,j) ) THEN
0328
0329 taperFctLoc =
0330 & ( -depthZ(k)/(hTransLay(i,j)+hMixLay(i,j))
ac0d9d59a8 Jean*0331 & *( 2. _d 0 + DoverLamb )
e2259a1942 Jean*0332 & )
0333
0334 taperFctHat =
0335 & ( -depthZ(k)/(hTransLay(i,j)+hMixLay(i,j))
0336 & * 2. _d 0
ac0d9d59a8 Jean*0337 & *( 1. _d 0 + DoverLamb )
e2259a1942 Jean*0338 & )
0339 ELSE
0340
0341 taperFctLoc =
0342 & ( -depthZ(k)/(hTransLay(i,j)+hMixLay(i,j))
ac0d9d59a8 Jean*0343 & *( 2. _d 0 + DoverLamb )
e2259a1942 Jean*0344 & )
0345 & - ( (depthZ(k)+hMixLay(i,j))*(depthZ(k)+hMixLay(i,j))
0346 & /( hTransLay(i,j)*hTransLay(i,j)
0347 & - hMixLay(i,j)*hMixLay(i,j) )
0348 & *( 1. _d 0 + hTransLay(i,j)*recipLambda(i,j) )
0349 & )
0350
0351 taperFctHat =
0352 & ( -depthZ(k)/(hTransLay(i,j)+hMixLay(i,j))
0353 & * 2. _d 0
ac0d9d59a8 Jean*0354 & *( 1. _d 0 + DoverLamb )
e2259a1942 Jean*0355 & )
0356 & - ( (depthZ(k)+hMixLay(i,j))*(depthZ(k)+hMixLay(i,j))
0357 & /( hTransLay(i,j)*hTransLay(i,j)
0358 & - hMixLay(i,j)*hMixLay(i,j) )
0359 & *( 1. _d 0 + hTransLay(i,j)*recipLambda(i,j)*2. _d 0 )
0360 & )
0361 ENDIF
0362
0363
0364
0365 Smod = MIN( baseSlope(i,j), GM_maxSlope )
0366 SlopeX(i,j) = SlopeX(i,j)*Smod*taperFctLoc
0367 SlopeY(i,j) = SlopeY(i,j)*Smod*taperFctLoc
0368
df4facdcff Jean*0369
0370 SlopeSqr(i,j) = MIN( baseSlope(i,j), GM_bigSlope )
0371 & *Smod*taperFctHat
e2259a1942 Jean*0372
0373 ELSE
0374
0375
0376
0377 Smod = MIN( SlopeMod(i,j), GM_maxSlope )
0378 SlopeX(i,j) = SlopeX(i,j)*Smod
0379 SlopeY(i,j) = SlopeY(i,j)*Smod
0380
df4facdcff Jean*0381
0382 SlopeSqr(i,j) = MIN( SlopeMod(i,j), GM_bigSlope )
0383 & *Smod
e2259a1942 Jean*0384
0385
0386 ENDIF
0387
0388 ENDDO
0389 ENDDO
0390
5755dbe390 Patr*0391 #endif /* GM_EXCLUDE_FM07_TAP */
0392
8233d0ceb9 Jean*0393 ELSEIF (GM_taper_scheme.EQ.'ac02') THEN
2092dbb101 Patr*0394
bd3f833a36 Jean*0395 #ifdef GM_EXCLUDE_AC02_TAP
2092dbb101 Patr*0396
bd3f833a36 Jean*0397 STOP 'Need to compile without "#define GM_EXCLUDE_AC02_TAP"'
0398
0399 #else /* GM_EXCLUDE_AC02_TAP */
2092dbb101 Patr*0400
bd3f833a36 Jean*0401
e2259a1942 Jean*0402
bd3f833a36 Jean*0403
2092dbb101 Patr*0404
14e0496834 Jean*0405 DO j=1-OLy+1,sNy+OLy-1
0406 DO i=1-OLx+1,sNx+OLx-1
5b172de0d2 Jean*0407 dRdSigmaLtd(i,j) = ( dSigmaDx(i,j)*dSigmaDx(i,j)
0408 & + dSigmaDy(i,j)*dSigmaDy(i,j)
0409 & )*convSlopeUnit*convSlopeUnit
0410 & + dSigmaDr(i,j)*dSigmaDr(i,j)
20b8842d78 Patr*0411 taperFct(i,j) = 1. _d 0
e2259a1942 Jean*0412
5b172de0d2 Jean*0413 IF ( dRdSigmaLtd(i,j).NE.zeroRL ) THEN
8233d0ceb9 Jean*0414 dRdSigmaLtd(i,j) = 1. _d 0 / ( dRdSigmaLtd(i,j) )
0415 SlopeSqr(i,j) = ( dSigmaDx(i,j)*dSigmaDx(i,j)
0416 & + dSigmaDy(i,j)*dSigmaDy(i,j)
0417 & )*dRdSigmaLtd(i,j)
5b172de0d2 Jean*0418 SlopeX(i,j) = dSigmaDx(i,j)
0419 & *dRdSigmaLtd(i,j)*dSigmaDr(i,j)
0420 SlopeY(i,j) = dSigmaDy(i,j)
0421 & *dRdSigmaLtd(i,j)*dSigmaDr(i,j)
8233d0ceb9 Jean*0422 ELSE
0423 SlopeSqr(i,j) = 0. _d 0
0424 SlopeX(i,j) = 0. _d 0
0425 SlopeY(i,j) = 0. _d 0
2092dbb101 Patr*0426 ENDIF
0f5add5564 Jean*0427 #ifndef ALLOW_AUTODIFF_TAMC
04522666de Ed H*0428
9cb619cfcd Patr*0429 IF ( SlopeSqr(i,j) .GT. maxSlopeSqr .AND.
0430 & SlopeSqr(i,j) .LT. GM_slopeSqCutoff ) THEN
0431 taperFct(i,j) = maxSlopeSqr/SlopeSqr(i,j)
8233d0ceb9 Jean*0432 ELSEIF ( SlopeSqr(i,j) .GE. GM_slopeSqCutoff ) THEN
9cb619cfcd Patr*0433 taperFct(i,j) = 0. _d 0
0434 ENDIF
0435 #endif
2092dbb101 Patr*0436 ENDDO
0437 ENDDO
0438
bd3f833a36 Jean*0439 #endif /* GM_EXCLUDE_AC02_TAP */
2092dbb101 Patr*0440
e9fd580bc4 Jean*0441 ELSE
0c49347dc7 Alis*0442
bd3f833a36 Jean*0443 #ifdef GM_EXCLUDE_TAPERING
0444
0445 STOP 'Need to compile without "#define GM_EXCLUDE_TAPERING"'
0446
0447 #else /* GM_EXCLUDE_TAPERING */
2092dbb101 Patr*0448
b6b11b9b2f Patr*0449
0450
e673dad091 Patr*0451
5b172de0d2 Jean*0452
e673dad091 Patr*0453
0454 #ifdef ALLOW_AUTODIFF_TAMC
113b21cd45 Mart*0455
0456
0457
0458
e673dad091 Patr*0459 #endif
0460
14e0496834 Jean*0461 DO j=1-OLy+1,sNy+OLy-1
0462 DO i=1-OLx+1,sNx+OLx-1
5b172de0d2 Jean*0463 IF ( dSigmaDr(i,j) .NE. zeroRL ) THEN
0464 IF ( dSigmaDr(i,j) .LE. GM_Small_Number )
0465 & dSigmaDr(i,j) = GM_Small_Number
b6b11b9b2f Patr*0466 ENDIF
e673dad091 Patr*0467 ENDDO
0468 ENDDO
0469
14e0496834 Jean*0470 DO j=1-OLy+1,sNy+OLy-1
0471 DO i=1-OLx+1,sNx+OLx-1
5b172de0d2 Jean*0472 IF ( dSigmaDr(i,j) .EQ. zeroRL ) THEN
0473 IF ( dSigmaDx(i,j) .NE. zeroRL ) THEN
df4facdcff Jean*0474 SlopeX(i,j) = SIGN( GM_bigSlope, dSigmaDx(i,j) )
15dd2a0860 Patr*0475 ELSE
0476 SlopeX(i,j) = 0. _d 0
0477 ENDIF
5b172de0d2 Jean*0478 IF ( dSigmaDy(i,j) .NE. zeroRL ) THEN
df4facdcff Jean*0479 SlopeY(i,j) = SIGN( GM_bigSlope, dSigmaDy(i,j) )
15dd2a0860 Patr*0480 ELSE
0481 SlopeY(i,j) = 0. _d 0
0482 ENDIF
b6b11b9b2f Patr*0483 ELSE
549d1a8d8c Jean*0484 dRdSigmaLtd(i,j) = 1. _d 0 / dSigmaDr(i,j)
5b172de0d2 Jean*0485 SlopeX(i,j) = dSigmaDx(i,j)*dRdSigmaLtd(i,j)
0486 SlopeY(i,j) = dSigmaDy(i,j)*dRdSigmaLtd(i,j)
0487
0488
b6b11b9b2f Patr*0489 ENDIF
e673dad091 Patr*0490 ENDDO
0491 ENDDO
0c49347dc7 Alis*0492
e673dad091 Patr*0493 #ifdef ALLOW_AUTODIFF_TAMC
5a6997640b Mart*0494
0495
0496
0497
e673dad091 Patr*0498 #endif
0c49347dc7 Alis*0499
14e0496834 Jean*0500 DO j=1-OLy+1,sNy+OLy-1
0501 DO i=1-OLx+1,sNx+OLx-1
b6b11b9b2f Patr*0502 SlopeSqr(i,j) = SlopeX(i,j)*SlopeX(i,j)
0503 & +SlopeY(i,j)*SlopeY(i,j)
0504 taperFct(i,j) = 1. _d 0
5a6997640b Mart*0505 #ifdef ALLOW_AUTODIFF_TAMC
0506 ENDDO
0507 ENDDO
0508
0509
0510 DO j=1-OLy+1,sNy+OLy-1
0511 DO i=1-OLx+1,sNx+OLx-1
0512 #endif
8233d0ceb9 Jean*0513 IF ( SlopeSqr(i,j) .GE. GM_slopeSqCutoff ) THEN
0514 slopeSqr(i,j) = GM_slopeSqCutoff
0515 taperFct(i,j) = 0. _d 0
00e514b3a3 Jean*0516 ENDIF
0c49347dc7 Alis*0517 ENDDO
0518 ENDDO
0519
e9fd580bc4 Jean*0520
0c49347dc7 Alis*0521
e9fd580bc4 Jean*0522 IF (GM_taper_scheme.EQ.'linear') THEN
0523
0524
14e0496834 Jean*0525 DO j=1-OLy+1,sNy+OLy-1
0526 DO i=1-OLx+1,sNx+OLx-1
0c49347dc7 Alis*0527
5b172de0d2 Jean*0528 IF ( SlopeSqr(i,j) .EQ. zeroRL ) THEN
b6b11b9b2f Patr*0529 taperFct(i,j) = 1. _d 0
8233d0ceb9 Jean*0530 ELSEIF ( SlopeSqr(i,j) .GT. maxSlopeSqr .AND.
bd3f833a36 Jean*0531 & SlopeSqr(i,j) .LT. GM_slopeSqCutoff ) THEN
e2259a1942 Jean*0532 taperFct(i,j) = SQRT(maxSlopeSqr / SlopeSqr(i,j))
df4facdcff Jean*0533 slopeSqr(i,j) = MIN( slopeSqr(i,j),GM_bigSlope*GM_bigSlope )
b6b11b9b2f Patr*0534 ENDIF
0c49347dc7 Alis*0535
e9fd580bc4 Jean*0536 ENDDO
0537 ENDDO
0538
0539 ELSEIF (GM_taper_scheme.EQ.'gkw91') THEN
0540
0541
14e0496834 Jean*0542 DO j=1-OLy+1,sNy+OLy-1
0543 DO i=1-OLx+1,sNx+OLx-1
e9fd580bc4 Jean*0544
5b172de0d2 Jean*0545 IF ( SlopeSqr(i,j) .EQ. zeroRL ) THEN
b6b11b9b2f Patr*0546 taperFct(i,j) = 1. _d 0
8233d0ceb9 Jean*0547 ELSEIF ( SlopeSqr(i,j) .GT. maxSlopeSqr .AND.
bd3f833a36 Jean*0548 & SlopeSqr(i,j) .LT. GM_slopeSqCutoff ) THEN
b6b11b9b2f Patr*0549 taperFct(i,j) = maxSlopeSqr/SlopeSqr(i,j)
0550 ENDIF
0c49347dc7 Alis*0551
0552 ENDDO
0553 ENDDO
0554
e9fd580bc4 Jean*0555 ELSEIF (GM_taper_scheme.EQ.'dm95') THEN
0c49347dc7 Alis*0556
0557
14e0496834 Jean*0558 DO j=1-OLy+1,sNy+OLy-1
0559 DO i=1-OLx+1,sNx+OLx-1
0c49347dc7 Alis*0560
5b172de0d2 Jean*0561 IF ( SlopeSqr(i,j) .EQ. zeroRL ) THEN
b6b11b9b2f Patr*0562 taperFct(i,j) = 1. _d 0
8233d0ceb9 Jean*0563 ELSEIF ( SlopeSqr(i,j) .LT. GM_slopeSqCutoff ) THEN
5b172de0d2 Jean*0564 Smod = SQRT(SlopeSqr(i,j))*convSlopeUnit
8233d0ceb9 Jean*0565 taperFct(i,j) = op5*( oneRL + TANH( (GM_Scrit-Smod)/GM_Sd ) )
b6b11b9b2f Patr*0566 ENDIF
0c49347dc7 Alis*0567 ENDDO
0568 ENDDO
0569
e9fd580bc4 Jean*0570 ELSEIF (GM_taper_scheme.EQ.'ldd97') THEN
0c49347dc7 Alis*0571
0572
14e0496834 Jean*0573 DO j=1-OLy+1,sNy+OLy-1
0574 DO i=1-OLx+1,sNx+OLx-1
0c49347dc7 Alis*0575
5b172de0d2 Jean*0576 IF ( SlopeSqr(i,j) .EQ. zeroRL ) THEN
b6b11b9b2f Patr*0577 taperFct(i,j) = 1. _d 0
07c79baaea Jean*0578 ELSEIF ( SlopeSqr(i,j) .LT. GM_slopeSqCutoff ) THEN
8233d0ceb9 Jean*0579 Smod = SQRT(SlopeSqr(i,j))
5b172de0d2 Jean*0580 f1 = op5*( oneRL
0581 & + TANH( (GM_Scrit-Smod*convSlopeUnit)/GM_Sd ) )
0582 Rnondim = rDepth/(Lrho(i,j)*Smod)
07c79baaea Jean*0583 IF ( Rnondim.GE.1. _d 0 ) THEN
0584 f2 = 1. _d 0
0585 ELSE
8233d0ceb9 Jean*0586 f2 = op5*( 1. _d 0 + SIN( fpi*(Rnondim-op5) ) )
07c79baaea Jean*0587 ENDIF
8233d0ceb9 Jean*0588 taperFct(i,j) = f1*f2
b6b11b9b2f Patr*0589 ENDIF
0c49347dc7 Alis*0590
0591 ENDDO
0592 ENDDO
0593
d8ee9652bd Gael*0594 ELSEIF (GM_taper_scheme.EQ.'stableGmAdjTap') THEN
0595
0a637c269e Ou W*0596 #ifdef GMREDI_WITH_STABLE_ADJOINT
d8ee9652bd Gael*0597
0a637c269e Ou W*0598
0599
d8ee9652bd Gael*0600
5a6997640b Mart*0601
1b7098b402 Gael*0602 slopeMax= 2. _d -3
d8ee9652bd Gael*0603
20dee61641 Mart*0604 # ifdef ALLOW_AUTODIFF_TAMC
5a6997640b Mart*0605
0606
0607
0608
0609 #endif
8233d0ceb9 Jean*0610 DO j=1-OLy+1,sNy+OLy-1
0611 DO i=1-OLx+1,sNx+OLx-1
0612 slopeSqTmp = SlopeX(i,j)*SlopeX(i,j)
0613 & + SlopeY(i,j)*SlopeY(i,j)
0a637c269e Ou W*0614
1b7098b402 Gael*0615 IF ( slopeSqTmp .GT. slopeMax**2 ) then
8233d0ceb9 Jean*0616 slopeTmp = SQRT(slopeSqTmp)
0617 SlopeX(i,j) = SlopeX(i,j)*slopeMax/slopeTmp
0618 SlopeY(i,j) = SlopeY(i,j)*slopeMax/slopeTmp
f5fe6ba1f3 Gael*0619 ENDIF
0a637c269e Ou W*0620 ENDDO
0621 ENDDO
0622
5a6997640b Mart*0623
7c50f07931 Mart*0624
0a637c269e Ou W*0625
8233d0ceb9 Jean*0626 DO j=1-OLy+1,sNy+OLy-1
0627 DO i=1-OLx+1,sNx+OLx-1
d8ee9652bd Gael*0628 SlopeSqr(i,j) = SlopeX(i,j)*SlopeX(i,j)
8233d0ceb9 Jean*0629 & + SlopeY(i,j)*SlopeY(i,j)
d8ee9652bd Gael*0630 taperFct(i,j) = 1. _d 0
0631 ENDDO
0632 ENDDO
0633
0a637c269e Ou W*0634 #else /* GMREDI_WITH_STABLE_ADJOINT */
0635
0636 STOP 'Need to compile wth "#define GMREDI_WITH_STABLE_ADJOINT"'
0637
d8ee9652bd Gael*0638 #endif /* GMREDI_WITH_STABLE_ADJOINT */
0639
f42e64b3e7 Jean*0640 ELSEIF (GM_taper_scheme.NE.' ') THEN
0641 STOP 'GMREDI_SLOPE_LIMIT: Bad GM_taper_scheme'
e9fd580bc4 Jean*0642 ENDIF
0c49347dc7 Alis*0643
bd3f833a36 Jean*0644 #endif /* GM_EXCLUDE_TAPERING */
2092dbb101 Patr*0645
e9fd580bc4 Jean*0646 ENDIF
0c49347dc7 Alis*0647
0648 #endif /* ALLOW_GMREDI */
0649
0650 RETURN
0651 END