File indexing completed on 2024-03-30 05:10:15 UTC
view on githubraw file Latest commit 598aebfc on 2024-03-29 19:16:48 UTC
4794324d0b Jean*0001 #include "PACKAGES_CONFIG.h"
0002 #include "CPP_OPTIONS.h"
517dbdc414 Jean*0003 #ifdef ALLOW_AUTODIFF
0004 # include "AUTODIFF_OPTIONS.h"
0005 #endif
fb8688a588 Jean*0006 #ifdef ALLOW_GENERIC_ADVDIFF
0007 # include "GAD_OPTIONS.h"
0008 #endif
4794324d0b Jean*0009
0010
0011
0012
0013 SUBROUTINE SALT_INTEGRATE(
0b191c5f5a Jean*0014 I bi, bj, recip_hFac,
fb8688a588 Jean*0015 I uFld, vFld, wFld,
0016 U KappaRk,
4794324d0b Jean*0017 I myTime, myIter, myThid )
0018
0019
0020
bd27360393 Jean*0021
0022
0023
0024
4794324d0b Jean*0025
e874fa47e5 Jean*0026
4794324d0b Jean*0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042
fb8688a588 Jean*0043
0044
4794324d0b Jean*0045
0046
0047
0048
0049
0050 IMPLICIT NONE
0051
0052 #include "SIZE.h"
0053 #include "EEPARAMS.h"
0054 #include "PARAMS.h"
0b191c5f5a Jean*0055 #include "GRID.h"
0056 #include "DYNVARS.h"
4794324d0b Jean*0057 #include "RESTART.h"
0058 #ifdef ALLOW_GENERIC_ADVDIFF
fb8688a588 Jean*0059 # include "GAD.h"
0060 # include "GAD_SOM_VARS.h"
4794324d0b Jean*0061 #endif
7c50f07931 Mart*0062 #ifdef ALLOW_AUTODIFF_TAMC
4794324d0b Jean*0063 # include "tamc.h"
0064 #endif
0065
0066
0067
0b191c5f5a Jean*0068
0069
0070
0071
0072
0073
0074
0075
fb8688a588 Jean*0076 INTEGER bi, bj
0b191c5f5a Jean*0077 _RS recip_hFac(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
0078 _RL uFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
0079 _RL vFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
0080 _RL wFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
0081 _RL KappaRk (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
4794324d0b Jean*0082 _RL myTime
0083 INTEGER myIter
0084 INTEGER myThid
0085
0086
0087 #ifdef ALLOW_GENERIC_ADVDIFF
ecf3fdef88 Jean*0088 #ifdef ALLOW_DIAGNOSTICS
0089
0090 LOGICAL DIAGNOSTICS_IS_ON
0091 EXTERNAL DIAGNOSTICS_IS_ON
0092 #endif /* ALLOW_DIAGNOSTICS */
0093
4794324d0b Jean*0094
0b191c5f5a Jean*0095
0096
0097
0098
0099
0100
0101
0102
0103
0104
0105
0106
0107
0108
0109
0110
0111
1ece1facde Jean*0112
ecf3fdef88 Jean*0113
e874fa47e5 Jean*0114
fb8688a588 Jean*0115 INTEGER iMin, iMax, jMin, jMax
4794324d0b Jean*0116 INTEGER i, j, k
0117 INTEGER kUp, kDown, kM1
0118 _RS xA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0119 _RS yA (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0120 _RS maskUp (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0121 _RL uTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0122 _RL vTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0123 _RL rTrans (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0124 _RL rTransKp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0b191c5f5a Jean*0125 _RL fZon (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0126 _RL fMer (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0127 _RL fVer (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
1ece1facde Jean*0128 _RL gS_loc (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
ecf3fdef88 Jean*0129 _RL gsForc (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
4794324d0b Jean*0130 _RL gs_AB (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
ecf3fdef88 Jean*0131 #ifdef ALLOW_DIAGNOSTICS
0132 LOGICAL diagForcing, diagAB_tend
0133 #endif
4794324d0b Jean*0134 LOGICAL calcAdvection
0135 INTEGER iterNb
0136 #ifdef ALLOW_ADAMSBASHFORTH_3
972c0130ec Jean*0137 INTEGER m2
4794324d0b Jean*0138 #endif
7c50f07931 Mart*0139 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0140
0141
0142 INTEGER tkey, kkey
7c50f07931 Mart*0143 #endif
fb8688a588 Jean*0144
0145
0146 iterNb = myIter
0147 IF (staggerTimeStep) iterNb = myIter - 1
4794324d0b Jean*0148
d3548146a8 Jean*0149
0150
0151
0152
0153
0154
0155
0156
0157
0158 iMin = 0
0159 iMax = sNx+1
0160 jMin = 0
0161 jMax = sNy+1
0162
ecf3fdef88 Jean*0163 #ifdef ALLOW_DIAGNOSTICS
0164 diagForcing = .FALSE.
0165 diagAB_tend = .FALSE.
0166 IF ( useDiagnostics .AND. saltForcing )
0167 & diagForcing = DIAGNOSTICS_IS_ON( 'gS_Forc ', myThid )
0168 IF ( useDiagnostics .AND. AdamsBashforthGs )
0169 & diagAB_tend = DIAGNOSTICS_IS_ON( 'AB_gS ', myThid )
0170 #endif
0171
4794324d0b Jean*0172 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0173 tkey = bi + (bj-1)*nSx + (ikey_dynamics-1)*nSx*nSy
4794324d0b Jean*0174 #endif /* ALLOW_AUTODIFF_TAMC */
0175
ff1bd1eafd Jean*0176
0177 IF ( AdamsBashforth_S ) THEN
0178
0179 #ifdef ALLOW_ADAMSBASHFORTH_3
0180
0181
0182 CALL ADAMS_BASHFORTH3(
0183 I bi, bj, 0, Nr,
0184 I salt(1-OLx,1-OLy,1,bi,bj),
0185 U gsNm, gs_AB,
0186 I saltStartAB, iterNb, myThid )
0187 #else /* ALLOW_ADAMSBASHFORTH_3 */
0188 CALL ADAMS_BASHFORTH2(
0189 I bi, bj, 0, Nr,
0190 I salt(1-OLx,1-OLy,1,bi,bj),
0191 U gsNm1(1-OLx,1-OLy,1,bi,bj), gs_AB,
0192 I saltStartAB, iterNb, myThid )
0193 #endif /* ALLOW_ADAMSBASHFORTH_3 */
0194 ENDIF
0195
fb8688a588 Jean*0196
0197 DO k=1,Nr
0198 DO j=1-OLy,sNy+OLy
0199 DO i=1-OLx,sNx+OLx
1ece1facde Jean*0200 gS_loc(i,j,k) = 0. _d 0
fb8688a588 Jean*0201 ENDDO
0202 ENDDO
0203 ENDDO
4794324d0b Jean*0204 DO j=1-OLy,sNy+OLy
0205 DO i=1-OLx,sNx+OLx
0b191c5f5a Jean*0206 fVer(i,j,1) = 0. _d 0
0207 fVer(i,j,2) = 0. _d 0
4794324d0b Jean*0208 ENDDO
0209 ENDDO
fb8688a588 Jean*0210 #ifdef ALLOW_AUTODIFF
0211 DO k=1,Nr
0212 DO j=1-OLy,sNy+OLy
0213 DO i=1-OLx,sNx+OLx
0214 kappaRk(i,j,k) = 0. _d 0
0215 ENDDO
0216 ENDDO
0217 ENDDO
4e4ad91a39 Jean*0218 #endif /* ALLOW_AUTODIFF */
0219 #ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0220
0221
ff1bd1eafd Jean*0222 # ifdef ALLOW_ADAMSBASHFORTH_3
598aebfcee Mart*0223
0224
ff1bd1eafd Jean*0225 # else
598aebfcee Mart*0226
ff1bd1eafd Jean*0227 # endif
4e4ad91a39 Jean*0228 #endif /* ALLOW_AUTODIFF_TAMC */
4794324d0b Jean*0229
fb8688a588 Jean*0230 #ifdef INCLUDE_CALC_DIFFUSIVITY_CALL
0b191c5f5a Jean*0231 CALL CALC_3D_DIFFUSIVITY(
fb8688a588 Jean*0232 I bi, bj, iMin, iMax, jMin, jMax,
0233 I GAD_SALINITY, useGMredi, useKPP,
0234 O kappaRk,
0235 I myThid )
d9efc637ba Mart*0236 # ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0237
d9efc637ba Mart*0238 # endif
fb8688a588 Jean*0239 #endif /* INCLUDE_CALC_DIFFUSIVITY_CALL */
0240
0241 #ifndef DISABLE_MULTIDIM_ADVECTION
0242
0243
0244
0245
0246
0247
0248
0249
0250
0251
0252 #ifdef GAD_ALLOW_TS_SOM_ADV
0253 # ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0254
0255
fb8688a588 Jean*0256 # endif
0257 IF ( saltSOM_Advection ) THEN
0258 # ifdef ALLOW_DEBUG
0259 IF (debugMode) CALL DEBUG_CALL('GAD_SOM_ADVECT',myThid)
0260 # endif
0261 CALL GAD_SOM_ADVECT(
e874fa47e5 Jean*0262 I saltImplVertAdv,
0263 I saltAdvScheme, saltVertAdvScheme, GAD_SALINITY,
0264 I dTtracerLev, uFld, vFld, wFld, salt,
fb8688a588 Jean*0265 U som_S,
1ece1facde Jean*0266 O gS_loc,
fb8688a588 Jean*0267 I bi, bj, myTime, myIter, myThid )
0268 ELSEIF (saltMultiDimAdvec) THEN
0269 #else /* GAD_ALLOW_TS_SOM_ADV */
0270 IF (saltMultiDimAdvec) THEN
0271 #endif /* GAD_ALLOW_TS_SOM_ADV */
0272 # ifdef ALLOW_DEBUG
0273 IF (debugMode) CALL DEBUG_CALL('GAD_ADVECTION',myThid)
0274 # endif
0275 CALL GAD_ADVECTION(
e874fa47e5 Jean*0276 I saltImplVertAdv,
0277 I saltAdvScheme, saltVertAdvScheme, GAD_SALINITY,
0278 I dTtracerLev, uFld, vFld, wFld, salt,
1ece1facde Jean*0279 O gS_loc,
fb8688a588 Jean*0280 I bi, bj, myTime, myIter, myThid )
0281 ENDIF
0282 #endif /* DISABLE_MULTIDIM_ADVECTION */
0283
0284
0285
0286
0287 calcAdvection = saltAdvection .AND. .NOT.saltMultiDimAdvec
4794324d0b Jean*0288 DO k=Nr,1,-1
0289 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0290 kkey = k + (tkey-1)*Nr
4794324d0b Jean*0291 #endif
0292 kM1 = MAX(1,k-1)
0293 kUp = 1+MOD(k+1,2)
0294 kDown= 1+MOD(k,2)
0295
0296 #ifdef ALLOW_AUTODIFF_TAMC
d9efc637ba Mart*0297
0298
4794324d0b Jean*0299 # ifdef ALLOW_ADAMSBASHFORTH_3
d9efc637ba Mart*0300
0301
4794324d0b Jean*0302 # else
d9efc637ba Mart*0303
4794324d0b Jean*0304 # endif
0305 #endif /* ALLOW_AUTODIFF_TAMC */
0306 CALL CALC_ADV_FLOW(
0307 I uFld, vFld, wFld,
0308 U rTrans,
0309 O uTrans, vTrans, rTransKp,
0310 O maskUp, xA, yA,
0311 I k, bi, bj, myThid )
0312
ecf3fdef88 Jean*0313
0314 DO j=1-OLy,sNy+OLy
0315 DO i=1-OLx,sNx+OLx
0316 gsForc(i,j) = 0. _d 0
0317 ENDDO
0318 ENDDO
0319 IF ( saltForcing ) THEN
0320 CALL APPLY_FORCING_S(
0321 U gsForc,
0322 I iMin,iMax,jMin,jMax, k, bi,bj,
0323 I myTime, myIter, myThid )
0324 #ifdef ALLOW_DIAGNOSTICS
0325 IF ( diagForcing ) THEN
0326 CALL DIAGNOSTICS_FILL(gsForc,'gS_Forc ',k,1,2,bi,bj,myThid)
0327 ENDIF
0328 #endif /* ALLOW_DIAGNOSTICS */
0329 ENDIF
0330
4794324d0b Jean*0331 #ifdef ALLOW_ADAMSBASHFORTH_3
972c0130ec Jean*0332
4794324d0b Jean*0333 m2 = 1 + MOD( iterNb ,2)
0334 CALL GAD_CALC_RHS(
0335 I bi, bj, iMin,iMax,jMin,jMax, k, kM1, kUp, kDown,
0336 I xA, yA, maskUp, uFld(1-OLx,1-OLy,k),
0337 I vFld(1-OLx,1-OLy,k), wFld(1-OLx,1-OLy,k),
0338 I uTrans, vTrans, rTrans, rTransKp,
0339 I diffKhS, diffK4S, KappaRk(1-OLx,1-OLy,k), diffKr4S,
cb30c1718d Jean*0340 I salt(1-OLx,1-OLy,1,bi,bj),
0341 I gsNm(1-OLx,1-OLy,1,bi,bj,m2), dTtracerLev,
4794324d0b Jean*0342 I GAD_SALINITY, saltAdvScheme, saltVertAdvScheme,
0343 I calcAdvection, saltImplVertAdv, AdamsBashforth_S,
46918f1b26 Jean*0344 I saltVertDiff4, useGMRedi, useKPP, salt_stayPositive,
0b191c5f5a Jean*0345 O fZon, fMer,
1ece1facde Jean*0346 U fVer, gS_loc,
4794324d0b Jean*0347 I myTime, myIter, myThid )
0348 #else /* ALLOW_ADAMSBASHFORTH_3 */
0349 CALL GAD_CALC_RHS(
0350 I bi, bj, iMin,iMax,jMin,jMax, k, kM1, kUp, kDown,
0351 I xA, yA, maskUp, uFld(1-OLx,1-OLy,k),
0352 I vFld(1-OLx,1-OLy,k), wFld(1-OLx,1-OLy,k),
0353 I uTrans, vTrans, rTrans, rTransKp,
0354 I diffKhS, diffK4S, KappaRk(1-OLx,1-OLy,k), diffKr4S,
cb30c1718d Jean*0355 I salt(1-OLx,1-OLy,1,bi,bj),
0356 I gsNm1(1-OLx,1-OLy,1,bi,bj), dTtracerLev,
4794324d0b Jean*0357 I GAD_SALINITY, saltAdvScheme, saltVertAdvScheme,
0358 I calcAdvection, saltImplVertAdv, AdamsBashforth_S,
46918f1b26 Jean*0359 I saltVertDiff4, useGMRedi, useKPP, salt_stayPositive,
0b191c5f5a Jean*0360 O fZon, fMer,
1ece1facde Jean*0361 U fVer, gS_loc,
4794324d0b Jean*0362 I myTime, myIter, myThid )
0363 #endif /* ALLOW_ADAMSBASHFORTH_3 */
0364
0365
ecf3fdef88 Jean*0366 IF ( saltForcing .AND. tracForcingOutAB.NE.1 ) THEN
0367 DO j=1-OLy,sNy+OLy
0368 DO i=1-OLx,sNx+OLx
1ece1facde Jean*0369 gS_loc(i,j,k) = gS_loc(i,j,k) + gsForc(i,j)
ecf3fdef88 Jean*0370 ENDDO
0371 ENDDO
0372 ENDIF
4794324d0b Jean*0373
0374 IF ( AdamsBashforthGs ) THEN
0375 #ifdef ALLOW_ADAMSBASHFORTH_3
0376 CALL ADAMS_BASHFORTH3(
0377 I bi, bj, k, Nr,
1ece1facde Jean*0378 U gS_loc, gsNm,
0379 O gs_AB,
4794324d0b Jean*0380 I saltStartAB, iterNb, myThid )
0381 #else
0382 CALL ADAMS_BASHFORTH2(
0383 I bi, bj, k, Nr,
1ece1facde Jean*0384 U gS_loc, gsNm1(1-OLx,1-OLy,1,bi,bj),
0385 O gs_AB,
4794324d0b Jean*0386 I saltStartAB, iterNb, myThid )
0387 #endif
0388 #ifdef ALLOW_DIAGNOSTICS
ecf3fdef88 Jean*0389 IF ( diagAB_tend ) THEN
4794324d0b Jean*0390 CALL DIAGNOSTICS_FILL(gs_AB,'AB_gS ',k,1,2,bi,bj,myThid)
0391 ENDIF
0392 #endif /* ALLOW_DIAGNOSTICS */
0393 ENDIF
0394
0395
ecf3fdef88 Jean*0396 IF ( saltForcing .AND. tracForcingOutAB.EQ.1 ) THEN
0397 DO j=1-OLy,sNy+OLy
0398 DO i=1-OLx,sNx+OLx
1ece1facde Jean*0399 gS_loc(i,j,k) = gS_loc(i,j,k) + gsForc(i,j)
ecf3fdef88 Jean*0400 ENDDO
0401 ENDDO
0402 ENDIF
4794324d0b Jean*0403
0404 #ifdef NONLIN_FRSURF
d9efc637ba Mart*0405 # ifdef ALLOW_AUTODIFF_TAMC
0406
0407 # endif
4794324d0b Jean*0408 IF (nonlinFreeSurf.GT.0) THEN
0409 CALL FREESURF_RESCALE_G(
0410 I bi, bj, k,
1ece1facde Jean*0411 U gS_loc,
4794324d0b Jean*0412 I myThid )
0413 IF ( AdamsBashforthGs ) THEN
0414 #ifdef ALLOW_ADAMSBASHFORTH_3
0415 # ifdef ALLOW_AUTODIFF_TAMC
d9efc637ba Mart*0416
0417
4794324d0b Jean*0418 # endif
0419 CALL FREESURF_RESCALE_G(
0420 I bi, bj, k,
1ece1facde Jean*0421 U gsNm(1-OLx,1-OLy,1,bi,bj,1),
4794324d0b Jean*0422 I myThid )
0423 CALL FREESURF_RESCALE_G(
0424 I bi, bj, k,
1ece1facde Jean*0425 U gsNm(1-OLx,1-OLy,1,bi,bj,2),
4794324d0b Jean*0426 I myThid )
0427 #else
d9efc637ba Mart*0428 # ifdef ALLOW_AUTODIFF_TAMC
0429
0430 # endif
4794324d0b Jean*0431 CALL FREESURF_RESCALE_G(
0432 I bi, bj, k,
1ece1facde Jean*0433 U gsNm1(1-OLx,1-OLy,1,bi,bj),
4794324d0b Jean*0434 I myThid )
0435 #endif
0436 ENDIF
0437 ENDIF
0438 #endif /* NONLIN_FRSURF */
0439
0440
0441 ENDDO
0442
0b191c5f5a Jean*0443 #ifdef ALLOW_DOWN_SLOPE
0444 IF ( useDOWN_SLOPE ) THEN
d9efc637ba Mart*0445 # ifdef ALLOW_AUTODIFF_TAMC
0446
0447 # endif
0b191c5f5a Jean*0448 IF ( usingPCoords ) THEN
0449 CALL DWNSLP_APPLY(
0450 I GAD_SALINITY, bi, bj, kSurfC,
6fca64f237 Jean*0451 I salt(1-OLx,1-OLy,1,bi,bj),
1ece1facde Jean*0452 U gS_loc,
6fca64f237 Jean*0453 I recip_hFac, recip_rA, recip_drF,
0454 I dTtracerLev, myTime, myIter, myThid )
0b191c5f5a Jean*0455 ELSE
0456 CALL DWNSLP_APPLY(
0457 I GAD_SALINITY, bi, bj, kLowC,
6fca64f237 Jean*0458 I salt(1-OLx,1-OLy,1,bi,bj),
1ece1facde Jean*0459 U gS_loc,
6fca64f237 Jean*0460 I recip_hFac, recip_rA, recip_drF,
0461 I dTtracerLev, myTime, myIter, myThid )
0b191c5f5a Jean*0462 ENDIF
0463 ENDIF
0464 #endif /* ALLOW_DOWN_SLOPE */
0465
1ece1facde Jean*0466
6fca64f237 Jean*0467 CALL TIMESTEP_TRACER(
0468 I bi, bj, dTtracerLev,
0469 I salt(1-OLx,1-OLy,1,bi,bj),
1ece1facde Jean*0470 U gS_loc,
6fca64f237 Jean*0471 I myTime, myIter, myThid )
0472
0b191c5f5a Jean*0473
0474
d9efc637ba Mart*0475 #ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0476
d9efc637ba Mart*0477 #endif /* ALLOW_AUTODIFF_TAMC */
0b191c5f5a Jean*0478 #ifdef INCLUDE_IMPLVERTADV_CODE
5685a5c914 Jean*0479 IF ( saltImplVertAdv .OR. implicitDiffusion ) THEN
0480
0481
d9efc637ba Mart*0482 # ifdef ALLOW_AUTODIFF_TAMC
598aebfcee Mart*0483
0484
0485
d9efc637ba Mart*0486 # endif /* ALLOW_AUTODIFF_TAMC */
0b191c5f5a Jean*0487 CALL GAD_IMPLICIT_R(
0488 I saltImplVertAdv, saltVertAdvScheme, GAD_SALINITY,
1ece1facde Jean*0489 I dTtracerLev, kappaRk, recip_hFac, wFld,
0490 I salt(1-OLx,1-OLy,1,bi,bj),
0491 U gS_loc,
0b191c5f5a Jean*0492 I bi, bj, myTime, myIter, myThid )
0493 ELSEIF ( implicitDiffusion ) THEN
0494 #else /* INCLUDE_IMPLVERTADV_CODE */
0495 IF ( implicitDiffusion ) THEN
0496 #endif /* INCLUDE_IMPLVERTADV_CODE */
0497 CALL IMPLDIFF(
0498 I bi, bj, iMin, iMax, jMin, jMax,
0499 I GAD_SALINITY, kappaRk, recip_hFac,
1ece1facde Jean*0500 U gS_loc,
0b191c5f5a Jean*0501 I myThid )
0502 ENDIF
0503
bd27360393 Jean*0504 IF ( AdamsBashforth_S ) THEN
0505
ff1bd1eafd Jean*0506 #ifdef ALLOW_ADAMSBASHFORTH_3
bd27360393 Jean*0507 CALL CYCLE_AB_TRACER(
cb30c1718d Jean*0508 I bi, bj, gS_loc,
0509 U salt(1-OLx,1-OLy,1,bi,bj),
0510 O gsNm(1-OLx,1-OLy,1,bi,bj,m2),
bd27360393 Jean*0511 I myTime, myIter, myThid )
0512 #else /* ALLOW_ADAMSBASHFORTH_3 */
ff1bd1eafd Jean*0513 CALL CYCLE_AB_TRACER(
0514 I bi, bj, gS_loc,
0515 U salt(1-OLx,1-OLy,1,bi,bj),
0516 O gsNm1(1-OLx,1-OLy,1,bi,bj),
0517 I myTime, myIter, myThid )
bd27360393 Jean*0518 #endif /* ALLOW_ADAMSBASHFORTH_3 */
ff1bd1eafd Jean*0519 ELSE
bd27360393 Jean*0520
0521 CALL CYCLE_TRACER(
0522 I bi, bj,
cb30c1718d Jean*0523 O salt(1-OLx,1-OLy,1,bi,bj),
0524 I gS_loc, myTime, myIter, myThid )
bd27360393 Jean*0525 ENDIF
0526
4794324d0b Jean*0527 #endif /* ALLOW_GENERIC_ADVDIFF */
0528
0529 RETURN
0530 END