C $Header: /home/ubuntu/mnt/e9_copy/MITgcm_contrib/atnguyen/code_21Dec2012_saltplume/salt_plume_apply.F,v 1.7 2014/05/04 09:25:16 atn Exp $ C $Name: $ #include "SALT_PLUME_OPTIONS.h" CBOP C !ROUTINE: SALT_PLUME_APPLY C !INTERFACE: SUBROUTINE SALT_PLUME_APPLY( I trIdentity, bi, bj, I recip_hFac_arg, I tracer,trApplyFlag, U trStar, I myTime, myIter, myThid ) C !DESCRIPTION: \bv C *==========================================================* C | SUBROUTINE SALT_PLUME_APPLY C | o Apply the salt_pume-transport to tracer field C *==========================================================* C \ev C !USES: IMPLICIT NONE C === Global variables === #include "SIZE.h" #include "GRID.h" #include "EEPARAMS.h" #include "PARAMS.h" #include "DYNVARS.h" #include "SALT_PLUME.h" #ifdef ALLOW_GENERIC_ADVDIFF # include "GAD.h" #endif C !INPUT/OUTPUT PARAMETERS: C === Routine arguments === C trIdentity :: tracer identification number C bi,bj :: Tile indices C recip_drF :: Reciprol of cell thickness C recip_hFac_arg :: Reciprol of cell open-depth factor C tracer :: tracer field at current time (input) C trApplyFlag:: [0]=update saltplume forcing T/S terms C :: [1]=update gTr tendency C trStar :: future tracer field (modified) C myTime :: Current time in simulation C myIter :: Current time-step number C myThid :: my Thread Id. number INTEGER trIdentity, trApplyFlag INTEGER bi, bj _RS recip_hFac_arg(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy) _RL tracer (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy) _RL trStar (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy) _RL myTime INTEGER myIter, myThid #ifdef ALLOW_SALT_PLUME #ifdef SALT_PLUME_VOLUME C !LOCAL VARIABLES: C === Local variables === C msgBuf :: Informational/error message buffer C plumetend :: forcing terms [W/m2 or kg/m2/s*psu] C work :: working array C CHARACTER*(MAX_LEN_MBUF) msgBuf INTEGER i, j, k INTEGER upward LOGICAL onOffFlag _RL gTr_Surf2kLev, gTr_Below2kLev, gTr_kLev2Above, & dSPvolBelow2kLev, gTr_totSurf2Below, & SurfVal, SEAICE_Tfrz, ConvertFac, recip_ConvertFac integer kp1, Nrp1 _RL plumetend(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr) _RL work(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1) #ifdef ALLOW_DIAGNOSTICS CHARACTER*8 diagName CHARACTER*4 diagSufx LOGICAL doDiagSPtend C- Functions: LOGICAL DIAGNOSTICS_IS_ON EXTERNAL DIAGNOSTICS_IS_ON #ifdef ALLOW_GENERIC_ADVDIFF CHARACTER*5 GAD_DIAG_SUFX EXTERNAL GAD_DIAG_SUFX #endif /* ALLOW_GENERIC_ADVDIFF */ #endif /* ALLOW_DIAGNOSTICS */ CEOP IF ( trApplyFlag.LT.0 .OR. trApplyFlag.GT.1) THEN STOP 'S/R SALT_PLUME_APPLY: incorrect setting of trApplyFlag!' ELSE SEAICE_Tfrz = -1.96 _d 0 onOffFlag = .FALSE. #ifdef ALLOW_GENERIC_ADVDIFF IF ( trIdentity.EQ.GAD_TEMPERATURE ) onOffFlag = .TRUE. IF ( trIdentity.EQ.GAD_SALINITY ) onOffFlag = .TRUE. #endif IF ( onOffFlag ) THEN C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| c upward = rkSign*NINT(-gravitySign) upward = 1 IF (usingZCoords) upward = -1 IF ( trIdentity.EQ.1 ) THEN SurfVal = SEAICE_Tfrz ConvertFac = HeatCapacity_Cp*rhoConst recip_ConvertFac = recip_Cp*mass2rUnit #ifdef ALLOW_DIAGNOSTICS IF ( useDiagnostics ) diagSufx = '_TH ' #endif /* ALLOW_DIAGNOSTICS */ ENDIF IF ( trIdentity.EQ.2 ) THEN SurfVal = SPbrineSconst ConvertFac = rhoConst recip_ConvertFac = mass2rUnit #ifdef ALLOW_DIAGNOSTICS IF ( useDiagnostics ) diagSufx = '_SLT' #endif /* ALLOW_DIAGNOSTICS */ ENDIF #ifdef ALLOW_DIAGNOSTICS doDiagSPtend = .FALSE. diagName = 'SPtd' IF ( useDiagnostics ) THEN C-- Set diagnostic suffix for the current tracer #ifdef ALLOW_GENERIC_ADVDIFF diagSufx = GAD_DIAG_SUFX( trIdentity, myThid ) #endif /* ALLOW_GENERIC_ADVDIFF */ diagName = 'SPtd'//diagSufx doDiagSPtend = DIAGNOSTICS_IS_ON(diagName,myThid) C WRITE(msgBuf,'(3A)') 'diagSufx,diagName: ', C & diagSufx,diagName C CALL PRINT_MESSAGE(msgBuf,standardMessageUnit, C & SQUEEZE_RIGHT, myThid ) ENDIF #endif /* ALLOW_DIAGNOSTICS */ C-- initializing: Nrp1=Nr+1 DO k=1,Nr DO j=1,OLy DO i=1,OLx plumetend(i,j,k) = 0. _d 0 work(i,j,k) = tracer(i,j,k,bi,bj) ENDDO ENDDO ENDDO DO j=1,OLy DO i=1,OLx work(i,j,Nrp1) = 0. _d 0 ENDDO ENDDO C----------------- Catn: After discussion with JM, it's cleaner to remove the negative c buoyancy associated with saltplumeflux from the surface lev c here instead of inside salt_plume_forcing_surf.F and c kpp_forcing_surf.F: DO k=Nr,1,-1 kp1=k+1 DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx C IF(trIdentity.EQ.GAD_SALINITY) SurfVal=SPbrineSalt(i,j,bi,bj) Catn: m/s*[degC,psu] gTr_totSurf2Below = SPbrineVolFlux(i,j,bi,bj)*SurfVal C dSPvolBelow2kLev = -dSPvolkLev2Above(i,j,kp1,bi,bj) gTr_Surf2kLev = dSPvolSurf2kLev(i,j,k,bi,bj) * SurfVal gTr_Below2kLev= dSPvolBelow2kLev * work(i,j,kp1) Catn: gTr_kLev2Above works even for kLev=1 because this is how much C volume of original [salinity,heat] associated with [SSS,SST] C was replaced by same volume of brine [salt,heat(from SEAICE_Tfrz)]. C Note: by design, dSPvolkLev2Above already is negative gTr_kLev2Above= dSPvolkLev2Above(i,j,k,bi,bj) * work(i,j,k) C salt: [m/s * psu * kg/m3] = [kg/s/m2 psu] = unit of saltPlumeFlux C theta:[m/s * kg/m3 * J/kg/degC * degC] = [W/m2] plumetend(i,j,k) = ConvertFac * & ( gTr_Surf2kLev + gTr_Below2kLev + gTr_kLev2Above ) IF(k.EQ.1) THEN plumetend(i,j,k) = plumetend(i,j,k) & - ConvertFac*gTr_totSurf2Below ENDIF IF(trApplyFlag.EQ.0) THEN Catn: report T/S SPforcing[T,S] related to saltplumeflux for kpp C and return zero to do_oceanic_phys.F; unit: [g/m2/s or W/m2] trStar(i,j,k,bi,bj) = 0. _d 0 IF (trIdentity.EQ.GAD_SALINITY) THEN SPforcingS(i,j,k,bi,bj)=plumetend(i,j,k) ENDIF IF (trIdentity.EQ.GAD_TEMPERATURE) THEN SPforcingT(i,j,k,bi,bj)=plumetend(i,j,k) ENDIF ELSE Catn: updating tendency gTr (gS,gT); unit: [psu/s or degC/s] trStar(i,j,k,bi,bj)=trStar(i,j,k,bi,bj)+plumetend(i,j,k) & *recip_drF(k)*recip_hFac_arg(i,j,k,bi,bj) & *recip_ConvertFac ENDIF ENDDO ENDDO ENDDO #ifdef ALLOW_DIAGNOSTICS IF ( doDiagSPtend ) & CALL DIAGNOSTICS_FILL(plumetend, diagName, 0,Nr,2,bi,bj,myThid) #endif /* ALLOW_DIAGNOSTICS */ C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| C-- end if on-off-flag ENDIF C-- end trApplyFlag ENDIF #endif /* SALT_PLUME_VOLUME */ #endif /* ALLOW_SALT_PLUME */ RETURN END