C $Header: /u/gcmpack/MITgcm/pkg/bbl/bbl_calc_rhs.F,v 1.8 2013/08/27 21:34:21 dimitri Exp $ C $Name: $ #include "BBL_OPTIONS.h" CBOP C !ROUTINE: BBL_CALC_RHS C !INTERFACE: SUBROUTINE BBL_CALC_RHS( I myTime, myIter, myThid ) C !DESCRIPTION: C Calculate tendency of tracers due to bottom boundary layer. C !USES: IMPLICIT NONE #include "SIZE.h" #include "EEPARAMS.h" #include "PARAMS.h" #include "GRID.h" #include "DYNVARS.h" #include "BBL.h" C !INPUT PARAMETERS: C myTime :: Current time in simulation C myIter :: Current time-step number C myThid :: my Thread Id number _RL myTime INTEGER myIter, myThid C !OUTPUT PARAMETERS: C !LOCAL VARIABLES: C bi,bj :: Tile indices C i,j :: Loop indices C d,r :: Donnor/Receiver indices C kBot :: k index of bottommost wet grid box C kLowC1 :: k index of bottommost (i,j) cell C kLowC2 :: k index of bottommost (i+1,j) or (i,j+1) cell C kl :: k index at which to compare 2 cells C thk_d :: thickness of donnor bottommost wet grid cell C thk_r :: thickness of receiver bottommost wet grid cell C bblEta_d :: bbl_eta of donnor cell C bblEta_r :: bbl_eta of receiver cell C resThk_r :: thk_r - bblEta_r C Theta_r :: Theta of receiver cell C bblTheta_d:: Theta of donnor bbl C bblTheta_r:: Theta of receiver bbl C resTheta_r:: Theta of resThk_r C Salt_r :: Salt of receiver cell C bblSalt_d :: Salt of donnor bbl C bblSalt_r :: Salt of receiver bbl C resSalt_r :: Salt of resThk_r C deltaRho :: density change C deltaDpt :: depth change C dVol :: horizontal volume transport C bbl_tend :: temporary variable for tendency terms C sloc :: salinity of bottommost wet grid box C tloc :: temperature of bottommost wet grid box C rholoc :: in situ density of bottommost wet grid box C rhoBBL :: in situ density of bottom boundary layer C bbl_rho1 :: local (i,j) density C bbl_rho2 :: local (i+1, j) or (i,j+1) density INTEGER bi, bj INTEGER i, j, d, r, kBot, kLowC1, kLowC2, kl _RL thk_d, thk_r, bblEta_d, bblEta_r, resThk_r _RL Theta_r, bblTheta_d, bblTheta_r, resTheta_r _RL Salt_r, bblSalt_d, bblSalt_r, resSalt_r _RL deltaRho, deltaDpt, dVol, bbl_tend _RL sloc ( 0:sNx+1, 0:sNy+1 ) _RL tloc ( 0:sNx+1, 0:sNy+1 ) _RL rholoc ( 0:sNx+1, 0:sNy+1 ) _RL rhoBBL ( 0:sNx+1, 0:sNy+1 ) _RL bbl_rho1, bbl_rho2 CEOP C-- Loops on tile indices bi,bj DO bj=myByLo(myThid),myByHi(myThid) DO bi=myBxLo(myThid),myBxHi(myThid) C Initialize tendency terms, make local copy of C bottomost temperature, salinity, in-situ density C and in-situ BBL density. DO j=0,sNy+1 DO i=0,sNx+1 bbl_TendTheta(i,j,bi,bj) = 0. _d 0 bbl_TendSalt (i,j,bi,bj) = 0. _d 0 kBot = max(1,kLowC(i,j,bi,bj)) tLoc(i,j) = theta(i,j,kBot,bi,bj) sLoc(i,j) = salt (i,j,kBot,bi,bj) rholoc(i,j) = rhoInSitu(i,j,kBot,bi,bj) IF ( kBot .EQ. Nr ) THEN rhoBBL(i,j) = bbl_rho_nr(i,j,bi,bj) ELSE rhoBBL(i,j) = rhoInSitu(i,j,kBot+1,bi,bj) ENDIF ENDDO ENDDO C==== Compute and apply vertical exchange between BBL and C residual volume of botommost wet grid box. C This operation does not change total tracer quantity C in botommost wet grid box. DO j=0,sNy+1 DO i=0,sNx+1 c DO j=-oly,sNy+oly c DO i=-olx,sNx+olx kBot = kLowC(i,j,bi,bj) IF ( kBot .GT. 0 ) THEN C If model density is lower than BBL, slowly diffuse upward. IF ( rhoLoc(i,j) .LT. rhoBBL(i,j) ) & bbl_eta(i,j,bi,bj) = max ( 0. _d 0 , & bbl_eta(i,j,bi,bj) - bbl_wvel * dTtracerLev(kBot) ) C If model density is higher than BBL then mix instantly. IF ( rhoLoc(i,j) .GE. rhoBBL(i,j) .OR. & bbl_eta(i,j,bi,bj) .EQ. 0. _d 0 ) THEN bbl_theta(i,j,bi,bj) = tLoc(i,j) bbl_salt (i,j,bi,bj) = sLoc(i,j) bbl_eta (i,j,bi,bj) = 0. _d 0 ENDIF ENDIF ENDDO ENDDO C==== Compute meridional bbl exchange at northern edge. j=sNy DO i=0,sNx+1 kLowC1 = kLowC(i,j ,bi,bj) kLowC2 = kLowC(i,j+1,bi,bj) IF ((kLowC1.GT.0).AND.(kLowC2.GT.0)) THEN C Compare the bbl densities at the higher pressure C (highest possible density of given t,s) C bbl in situ density is stored in k > kLowC indices kl = MAX(kLowC1,kLowC2) + 1 deltaDpt = R_low(i,j,bi,bj) + bbl_eta(i,j,bi,bj) - & R_low(i,j+1,bi,bj) - bbl_eta(i,j+1,bi,bj) IF ( deltaDpt .GT. 0. ) THEN IF ( kl .GT. Nr ) THEN bbl_rho1 = bbl_rho_nr(i,j,bi,bj) ELSE bbl_rho1 = rhoInSitu(i,j,kl,bi,bj) ENDIF bbl_rho2 = rhoInSitu(i,j+1,kLowC2,bi,bj) ELSE bbl_rho1 = rhoInSitu(i,j,kLowC1,bi,bj) IF ( kl .GT. Nr ) THEN bbl_rho2 = bbl_rho_nr(i,j+1,bi,bj) ELSE bbl_rho2 = rhoInSitu(i,j+1,kl,bi,bj) ENDIF ENDIF deltaRho = bbl_rho2 - bbl_rho1 IF ( (deltaRho*deltaDpt) .LT. 0. ) THEN C If heavy BBL water is higher than light BBL water, C exchange properties laterally. C Determine donnor and receiver cells. IF ( bbl_rho1 .GT. bbl_rho2 ) THEN d = j r = j + 1 ELSE d = j + 1 r = j ENDIF C Replenish thickness of donor cell, if needed. thk_d = drF(kLowC(i,d,bi,bj)) * & hFacC(i,d,kLowC(i,d,bi,bj),bi,bj) IF ( ( bbl_theta(i,d,bi,bj) .EQ. tloc(i,d) ) .AND. & ( bbl_salt (i,d,bi,bj) .EQ. sloc(i,d) ) .AND. & ( bbl_eta (i,d,bi,bj) .LT. bbl_initEta ) ) & bbl_eta(i,d,bi,bj) = min ( bbl_initEta, thk_d ) C Compute some donnor and receiver cell properties. thk_r = drF(kLowC(i,r,bi,bj)) * & hFacC(i,r,kLowC(i,r,bi,bj),bi,bj) Theta_r = tLoc(i,r) Salt_r = sLoc(i,r) bblTheta_d = bbl_theta(i,d,bi,bj) bblTheta_r = bbl_theta(i,r,bi,bj) bblSalt_d = bbl_salt (i,d,bi,bj) bblSalt_r = bbl_salt (i,r,bi,bj) bblEta_d = bbl_eta (i,d,bi,bj) bblEta_r = bbl_eta (i,r,bi,bj) resThk_r = thk_r - bblEta_r resTheta_r = (Theta_r*thk_r-bblTheta_r*bblEta_r)/resThk_r resSalt_r = (Salt_r *thk_r-bblSalt_r *bblEta_r)/resThk_r C Compute volume transport from donnor to receiver. dVol = min ( bblEta_d * rA(i,d,bi,bj) / 2. _d 0, & resThk_r * rA(i,r,bi,bj) / 2. _d 0, & dxG(i,j+1,bi,bj) * bblEta_d * bbl_hvel * deltaT ) C Compute temperature tracer tendencies for donor and receiver cell. bbl_tend = dVol * (bblTheta_d - resTheta_r) / deltaT bbl_TendTheta(i,d,bi,bj) = bbl_TendTheta(i,d,bi,bj) - & bbl_tend / rA(i,d,bi,bj) / thk_d bbl_TendTheta(i,r,bi,bj) = bbl_TendTheta(i,r,bi,bj) + & bbl_tend / rA(i,r,bi,bj) / thk_r C Compute salinity tracer tendencies for donor and receiver cell. bbl_tend = dVol * (bblSalt_d - resSalt_r) / deltaT bbl_TendSalt(i,d,bi,bj) = bbl_TendSalt(i,d,bi,bj) - & bbl_tend / rA(i,d,bi,bj) / thk_d bbl_TendSalt(i,r,bi,bj) = bbl_TendSalt(i,r,bi,bj) + & bbl_tend / rA(i,r,bi,bj) / thk_r C Adjust pbl thickness and tracer properties. bbl_eta(i,d,bi,bj) = bblEta_d - dVol / rA(i,d,bi,bj) IF ( bbl_eta(i,d,bi,bj) .LT. 0.0001 ) THEN bbl_eta(i,d,bi,bj) = 0. _d 0 bbl_theta(i,d,bi,bj) = tLoc(i,d) bbl_salt (i,d,bi,bj) = sLoc(i,d) ENDIF bbl_eta(i,r,bi,bj) = bblEta_r + dVol / rA(i,r,bi,bj) bbl_theta(i,r,bi,bj) = ( dVol * bblTheta_d + & bblEta_r * rA(i,r,bi,bj) * bblTheta_r ) / & ( bbl_eta(i,r,bi,bj) * rA(i,r,bi,bj) ) bbl_salt(i,r,bi,bj) = ( dVol * bblSalt_d + & bblEta_r * rA(i,r,bi,bj) * bblSalt_r ) / & ( bbl_eta(i,r,bi,bj) * rA(i,r,bi,bj) ) ENDIF ENDIF ENDDO C==== Compute meridional bbl exchange inside tile. DO j=0,sNy-1 DO i=0,sNx+1 kLowC1 = kLowC(i,j ,bi,bj) kLowC2 = kLowC(i,j+1,bi,bj) IF ((kLowC1.GT.0).AND.(kLowC2.GT.0)) THEN C Compare the bbl densities at the higher pressure C (highest possible density of given t,s) C bbl in situ density is stored in k > kLowC indices kl = MAX(kLowC1,kLowC2) + 1 deltaDpt = R_low(i,j,bi,bj) + bbl_eta(i,j,bi,bj) - & R_low(i,j+1,bi,bj) - bbl_eta(i,j+1,bi,bj) IF ( deltaDpt .GT. 0. ) THEN IF ( kl .GT. Nr ) THEN bbl_rho1 = bbl_rho_nr(i,j,bi,bj) ELSE bbl_rho1 = rhoInSitu(i,j,kl,bi,bj) ENDIF bbl_rho2 = rhoInSitu(i,j+1,kLowC2,bi,bj) ELSE bbl_rho1 = rhoInSitu(i,j,kLowC1,bi,bj) IF ( kl .GT. Nr ) THEN bbl_rho2 = bbl_rho_nr(i,j+1,bi,bj) ELSE bbl_rho2 = rhoInSitu(i,j+1,kl,bi,bj) ENDIF ENDIF deltaRho = bbl_rho2 - bbl_rho1 IF ( (deltaRho*deltaDpt) .LT. 0. ) THEN C If heavy BBL water is higher than light BBL water, C exchange properties laterally. C Determine donnor and receiver cells. IF ( bbl_rho1 .GT. bbl_rho2 ) THEN d = j r = j + 1 ELSE d = j + 1 r = j ENDIF C Replenish thickness of donor cell, if needed. thk_d = drF(kLowC(i,d,bi,bj)) * & hFacC(i,d,kLowC(i,d,bi,bj),bi,bj) IF ( ( bbl_theta(i,d,bi,bj) .EQ. tloc(i,d) ) .AND. & ( bbl_salt (i,d,bi,bj) .EQ. sloc(i,d) ) .AND. & ( bbl_eta (i,d,bi,bj) .LT. bbl_initEta ) ) & bbl_eta(i,d,bi,bj) = min ( bbl_initEta, thk_d ) C Compute some donnor and receiver cell properties. thk_r = drF(kLowC(i,r,bi,bj)) * & hFacC(i,r,kLowC(i,r,bi,bj),bi,bj) Theta_r = tLoc(i,r) Salt_r = sLoc(i,r) bblTheta_d = bbl_theta(i,d,bi,bj) bblTheta_r = bbl_theta(i,r,bi,bj) bblSalt_d = bbl_salt (i,d,bi,bj) bblSalt_r = bbl_salt (i,r,bi,bj) bblEta_d = bbl_eta (i,d,bi,bj) bblEta_r = bbl_eta (i,r,bi,bj) resThk_r = thk_r - bblEta_r resTheta_r = (Theta_r*thk_r-bblTheta_r*bblEta_r)/resThk_r resSalt_r = (Salt_r *thk_r-bblSalt_r *bblEta_r)/resThk_r C Compute volume transport from donnor to receiver. dVol = min ( bblEta_d * rA(i,d,bi,bj) / 2. _d 0, & resThk_r * rA(i,r,bi,bj) / 2. _d 0, & dxG(i,j+1,bi,bj) * bblEta_d * bbl_hvel * deltaT ) C Compute temperature tracer tendencies for donor and receiver cell. bbl_tend = dVol * (bblTheta_d - resTheta_r) / deltaT bbl_TendTheta(i,d,bi,bj) = bbl_TendTheta(i,d,bi,bj) - & bbl_tend / rA(i,d,bi,bj) / thk_d bbl_TendTheta(i,r,bi,bj) = bbl_TendTheta(i,r,bi,bj) + & bbl_tend / rA(i,r,bi,bj) / thk_r C Compute salinity tracer tendencies for donor and receiver cell. bbl_tend = dVol * (bblSalt_d - resSalt_r) / deltaT bbl_TendSalt(i,d,bi,bj) = bbl_TendSalt(i,d,bi,bj) - & bbl_tend / rA(i,d,bi,bj) / thk_d bbl_TendSalt(i,r,bi,bj) = bbl_TendSalt(i,r,bi,bj) + & bbl_tend / rA(i,r,bi,bj) / thk_r C Adjust pbl thickness and tracer properties. bbl_eta(i,d,bi,bj) = bblEta_d - dVol / rA(i,d,bi,bj) IF ( bbl_eta(i,d,bi,bj) .LT. 0.0001 ) THEN bbl_eta(i,d,bi,bj) = 0. _d 0 bbl_theta(i,d,bi,bj) = tLoc(i,d) bbl_salt (i,d,bi,bj) = sLoc(i,d) ENDIF bbl_eta(i,r,bi,bj) = bblEta_r + dVol / rA(i,r,bi,bj) bbl_theta(i,r,bi,bj) = ( dVol * bblTheta_d + & bblEta_r * rA(i,r,bi,bj) * bblTheta_r ) / & ( bbl_eta(i,r,bi,bj) * rA(i,r,bi,bj) ) bbl_salt(i,r,bi,bj) = ( dVol * bblSalt_d + & bblEta_r * rA(i,r,bi,bj) * bblSalt_r ) / & ( bbl_eta(i,r,bi,bj) * rA(i,r,bi,bj) ) ENDIF ENDIF ENDDO ENDDO C==== Compute zonal bbl exchange at Eastern edge. i=sNx DO j=1,sNy kLowC1 = kLowC(i ,j,bi,bj) kLowC2 = kLowC(i+1,j,bi,bj) IF ((kLowC1.GT.0).AND.(kLowC2.GT.0)) THEN C Compare the bbl densities at the higher pressure C (highest possible density of given t,s) C bbl in situ density is stored in k > kLowC indices kl = MAX(kLowC1,kLowC2) + 1 deltaDpt = R_low(i,j,bi,bj) + bbl_eta(i,j,bi,bj) - & R_low(i+1,j,bi,bj) - bbl_eta(i+1,j,bi,bj) IF ( deltaDpt .GT. 0. ) THEN IF ( kl .GT. Nr ) THEN bbl_rho1 = bbl_rho_nr(i,j,bi,bj) ELSE bbl_rho1 = rhoInSitu(i,j,kl,bi,bj) ENDIF bbl_rho2 = rhoInSitu(i+1,j,kLowC2,bi,bj) ELSE bbl_rho1 = rhoInSitu(i,j,kLowC1,bi,bj) IF ( kl .GT. Nr ) THEN bbl_rho2 = bbl_rho_nr(i+1,j,bi,bj) ELSE bbl_rho2 = rhoInSitu(i+1,j,kl,bi,bj) ENDIF ENDIF deltaRho = bbl_rho2 - bbl_rho1 IF ( (deltaRho*deltaDpt) .LT. 0. ) THEN C If heavy BBL water is higher than light BBL water, C exchange properties laterally. C Determine donnor and receiver cells. IF ( bbl_rho1 .GT. bbl_rho2 ) THEN d = i r = i + 1 ELSE d = i + 1 r = i ENDIF C Replenish thickness of donor cell, if needed. thk_d = drF(kLowC(d,j,bi,bj)) * & hFacC(d,j,kLowC(d,j,bi,bj),bi,bj) IF ( ( bbl_theta(d,j,bi,bj) .EQ. tloc(d,j) ) .AND. & ( bbl_salt (d,j,bi,bj) .EQ. sloc(d,j) ) .AND. & ( bbl_eta (d,j,bi,bj) .LT. bbl_initEta ) ) & bbl_eta(d,j,bi,bj) = min ( bbl_initEta, thk_d ) C Compute some donnor and receiver cell properties. thk_r = drF(kLowC(r,j,bi,bj)) * & hFacC(r,j,kLowC(r,j,bi,bj),bi,bj) Theta_r = tLoc(r,j) Salt_r = sLoc(r,j) bblTheta_d = bbl_theta(d,j,bi,bj) bblTheta_r = bbl_theta(r,j,bi,bj) bblSalt_d = bbl_salt (d,j,bi,bj) bblSalt_r = bbl_salt (r,j,bi,bj) bblEta_d = bbl_eta (d,j,bi,bj) bblEta_r = bbl_eta (r,j,bi,bj) resThk_r = thk_r - bblEta_r resTheta_r = (Theta_r*thk_r-bblTheta_r*bblEta_r)/resThk_r resSalt_r = (Salt_r *thk_r-bblSalt_r *bblEta_r)/resThk_r C Compute volume transport from donnor to receiver. dVol = min ( bblEta_d * rA(d,j,bi,bj) / 2. _d 0, & resThk_r * rA(r,j,bi,bj) / 2. _d 0, & dxG(i+1,j,bi,bj) * bblEta_d * bbl_hvel * deltaT ) C Compute temperature tracer tendencies for donor and receiver cell. bbl_tend = dVol * (bblTheta_d - resTheta_r) / deltaT bbl_TendTheta(d,j,bi,bj) = bbl_TendTheta(d,j,bi,bj) - & bbl_tend / rA(d,j,bi,bj) / thk_d bbl_TendTheta(r,j,bi,bj) = bbl_TendTheta(r,j,bi,bj) + & bbl_tend / rA(r,j,bi,bj) / thk_r C Compute salinity tracer tendencies for donor and receiver cell. bbl_tend = dVol * (bblSalt_d - resSalt_r) / deltaT bbl_TendSalt(d,j,bi,bj) = bbl_TendSalt(d,j,bi,bj) - & bbl_tend / rA(d,j,bi,bj) / thk_d bbl_TendSalt(r,j,bi,bj) = bbl_TendSalt(r,j,bi,bj) + & bbl_tend / rA(r,j,bi,bj) / thk_r C Adjust pbl thickness and tracer properties. bbl_eta(d,j,bi,bj) = bblEta_d - dVol / rA(d,j,bi,bj) IF ( bbl_eta(d,j,bi,bj) .LT. 0.0001 ) THEN bbl_eta(d,j,bi,bj) = 0. _d 0 bbl_theta(d,j,bi,bj) = tLoc(d,j) bbl_salt (d,j,bi,bj) = sLoc(d,j) ENDIF bbl_eta(r,j,bi,bj) = bblEta_r + dVol / rA(r,j,bi,bj) bbl_theta(r,j,bi,bj) = ( dVol * bblTheta_d + & bblEta_r * rA(r,j,bi,bj) * bblTheta_r ) / & ( bbl_eta(r,j,bi,bj) * rA(r,j,bi,bj) ) bbl_salt(r,j,bi,bj) = ( dVol * bblSalt_d + & bblEta_r * rA(r,j,bi,bj) * bblSalt_r ) / & ( bbl_eta(r,j,bi,bj) * rA(r,j,bi,bj) ) ENDIF ENDIF ENDDO C==== Compute zonal bbl exchange inside tile. DO j=1,sNy DO i=0,sNx-1 kLowC1 = kLowC(i ,j,bi,bj) kLowC2 = kLowC(i+1,j,bi,bj) IF ((kLowC1.GT.0).AND.(kLowC2.GT.0)) THEN C Compare the bbl densities at the higher pressure C (highest possible density of given t,s) C bbl in situ density is stored in k > kLowC indices kl = MAX(kLowC1,kLowC2) + 1 deltaDpt = R_low(i,j,bi,bj) + bbl_eta(i,j,bi,bj) - & R_low(i+1,j,bi,bj) - bbl_eta(i+1,j,bi,bj) IF ( deltaDpt .GT. 0. ) THEN IF ( kl .GT. Nr ) THEN bbl_rho1 = bbl_rho_nr(i,j,bi,bj) ELSE bbl_rho1 = rhoInSitu(i,j,kl,bi,bj) ENDIF bbl_rho2 = rhoInSitu(i+1,j,kLowC2,bi,bj) ELSE bbl_rho1 = rhoInSitu(i,j,kLowC1,bi,bj) IF ( kl .GT. Nr ) THEN bbl_rho2 = bbl_rho_nr(i+1,j,bi,bj) ELSE bbl_rho2 = rhoInSitu(i+1,j,kl,bi,bj) ENDIF ENDIF deltaRho = bbl_rho2 - bbl_rho1 IF ( (deltaRho*deltaDpt) .LT. 0. ) THEN C If heavy BBL water is higher than light BBL water, C exchange properties laterally. C Determine donnor and receiver cells. IF ( bbl_rho1 .GT. bbl_rho2 ) THEN d = i r = i + 1 ELSE d = i + 1 r = i ENDIF C Replenish thickness of donor cell, if needed. thk_d = drF(kLowC(d,j,bi,bj)) * & hFacC(d,j,kLowC(d,j,bi,bj),bi,bj) IF ( ( bbl_theta(d,j,bi,bj) .EQ. tloc(d,j) ) .AND. & ( bbl_salt (d,j,bi,bj) .EQ. sloc(d,j) ) .AND. & ( bbl_eta (d,j,bi,bj) .LT. bbl_initEta ) ) & bbl_eta(d,j,bi,bj) = min ( bbl_initEta, thk_d ) C Compute some donnor and receiver cell properties. thk_r = drF(kLowC(r,j,bi,bj)) * & hFacC(r,j,kLowC(r,j,bi,bj),bi,bj) Theta_r = tLoc(r,j) Salt_r = sLoc(r,j) bblTheta_d = bbl_theta(d,j,bi,bj) bblTheta_r = bbl_theta(r,j,bi,bj) bblSalt_d = bbl_salt (d,j,bi,bj) bblSalt_r = bbl_salt (r,j,bi,bj) bblEta_d = bbl_eta (d,j,bi,bj) bblEta_r = bbl_eta (r,j,bi,bj) resThk_r = thk_r - bblEta_r resTheta_r = (Theta_r*thk_r-bblTheta_r*bblEta_r)/resThk_r resSalt_r = (Salt_r *thk_r-bblSalt_r *bblEta_r)/resThk_r C Compute volume transport from donnor to receiver. dVol = min ( bblEta_d * rA(d,j,bi,bj) / 2. _d 0, & resThk_r * rA(r,j,bi,bj) / 2. _d 0, & dxG(i+1,j,bi,bj) * bblEta_d * bbl_hvel * deltaT ) C Compute temperature tracer tendencies for donor and receiver cell. bbl_tend = dVol * (bblTheta_d - resTheta_r) / deltaT bbl_TendTheta(d,j,bi,bj) = bbl_TendTheta(d,j,bi,bj) - & bbl_tend / rA(d,j,bi,bj) / thk_d bbl_TendTheta(r,j,bi,bj) = bbl_TendTheta(r,j,bi,bj) + & bbl_tend / rA(r,j,bi,bj) / thk_r C Compute salinity tracer tendencies for donor and receiver cell. bbl_tend = dVol * (bblSalt_d - resSalt_r) / deltaT bbl_TendSalt(d,j,bi,bj) = bbl_TendSalt(d,j,bi,bj) - & bbl_tend / rA(d,j,bi,bj) / thk_d bbl_TendSalt(r,j,bi,bj) = bbl_TendSalt(r,j,bi,bj) + & bbl_tend / rA(r,j,bi,bj) / thk_r C Adjust pbl thickness and tracer properties. bbl_eta(d,j,bi,bj) = bblEta_d - dVol / rA(d,j,bi,bj) IF ( bbl_eta(d,j,bi,bj) .LT. 0.0001 ) THEN bbl_eta(d,j,bi,bj) = 0. _d 0 bbl_theta(d,j,bi,bj) = tLoc(d,j) bbl_salt (d,j,bi,bj) = sLoc(d,j) ENDIF bbl_eta(r,j,bi,bj) = bblEta_r + dVol / rA(r,j,bi,bj) bbl_theta(r,j,bi,bj) = ( dVol * bblTheta_d + & bblEta_r * rA(r,j,bi,bj) * bblTheta_r ) / & ( bbl_eta(r,j,bi,bj) * rA(r,j,bi,bj) ) bbl_salt(r,j,bi,bj) = ( dVol * bblSalt_d + & bblEta_r * rA(r,j,bi,bj) * bblSalt_r ) / & ( bbl_eta(r,j,bi,bj) * rA(r,j,bi,bj) ) ENDIF ENDIF ENDDO ENDDO C-- end bi,bj loops. ENDDO ENDDO CALL EXCH_XY_RL( bbl_eta , myThid ) CALL EXCH_XY_RL( bbl_theta , myThid ) CALL EXCH_XY_RL( bbl_salt , myThid ) CALL EXCH_XY_RL( bbl_TendTheta, myThid ) CALL EXCH_XY_RL( bbl_TendSalt , myThid ) RETURN END