C $Header: /u/gcmpack/MITgcm/pkg/layers/layers_save.F,v 1.4 2015/06/10 20:51:17 jmc Exp $
C $Name:  $

#include "LAYERS_OPTIONS.h"
C--  File layers_save.F:
C--   Contents
C--   o LAYERS_FILL
C--   o LAYERS_FILL_FIELD

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CBOP

      SUBROUTINE LAYERS_FILL(
     I               df, trIdentity, fluxid,
     I               kLev, nLevs, bibjFlg, biArg, bjArg, myThid )
C     !DESCRIPTION: \bv
C     *==========================================================*
C     | SUBROUTINE LAYERS_FILL
C     | "Remember" the merid. advective flux for use later in layers_thermodynamics
C     *==========================================================*
       IMPLICIT NONE
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "LAYERS_SIZE.h"
#include "LAYERS.h"
C***********************************************************************
C   This is designed to look and work exactly like the a regular
C   diagnostics_fill call.
C***********************************************************************
C     surfflux  :: The surface temperature flux, the same as what is filled into
C                   the TFLUX and SFLUX diagnostics
C     trIdentity:: Index to let us know what tracer it is (1 for T, 2 for S)
C     kLev      :: Integer flag for vertical levels:
C                  > 0 (any integer): WHICH single level to increment in qdiag.
C                  0,-1 to increment "nLevs" levels in qdiag,
C                  0 : fill-in in the same order as the input array
C                  -1: fill-in in reverse order.
C     nLevs     :: indicates Number of levels of the input field array
C                  (whether to fill-in all the levels (kLev<1) or just one (kLev>0))
C     bibjFlg   :: Integer flag to indicate instructions for bi bj loop
C                  0 indicates that the bi-bj loop must be done here
C                  1 indicates that the bi-bj loop is done OUTSIDE
C                  2 indicates that the bi-bj loop is done OUTSIDE
C                     AND that we have been sent a local array (with overlap regions)
C                  3 indicates that the bi-bj loop is done OUTSIDE
C                     AND that we have been sent a local array
C                     AND that the array has no overlap region (interior only)
C                  NOTE - bibjFlg can be NEGATIVE to indicate not to increment counter
C     biArg     :: X-direction tile number - used for bibjFlg=1-3
C     bjArg     :: Y-direction tile number - used for bibjFlg=1-3

C       _RL df(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
       _RL df(*)
       INTEGER trIdentity, kLev, nLevs, bibjFlg, biArg, bjArg
       INTEGER myThid
       CHARACTER*(3) fluxid

#ifdef LAYERS_THERMODYNAMICS

C !LOCAL VARIABLES: ====================================================
C i,j              :: loop indices
C msgBuf           :: error message buffer
       CHARACTER*(MAX_LEN_MBUF) msgBuf

       IF ((trIdentity.EQ.1).OR.(trIdentity.EQ.2)) THEN

       IF (fluxid.EQ.'SUR') THEN
         CALL LAYERS_FILL_FIELD(df, trIdentity, 1, layers_surfflux,'M',
     &           klev, nLevs, bibjFlg, biArg, bjArg, myThid)
       ELSE IF (fluxid.EQ.'DFX') THEN
         CALL LAYERS_FILL_FIELD(df, trIdentity, Nr, layers_dfx,'U',
     &           kLev, nLevs, bibjFlg, biArg, bjArg, myThid)
       ELSE IF (fluxid.EQ.'DFY') THEN
         CALL LAYERS_FILL_FIELD(df, trIdentity, Nr, layers_dfy,'V',
     &           kLev, nLevs, bibjFlg, biArg, bjArg, myThid)
       ELSE IF (fluxid.EQ.'DFR') THEN
         CALL LAYERS_FILL_FIELD(df, trIdentity, Nr, layers_dfr,'M',
     &           kLev, nLevs, bibjFlg, biArg, bjArg, myThid)
       ELSE IF (fluxid.EQ.'AFX') THEN
         CALL LAYERS_FILL_FIELD(df, trIdentity, Nr, layers_afx,'U',
     &           kLev, nLevs, bibjFlg, biArg, bjArg, myThid)
       ELSE IF (fluxid.EQ.'AFY') THEN
         CALL LAYERS_FILL_FIELD(df, trIdentity, Nr, layers_afy,'V',
     &           kLev, nLevs, bibjFlg, biArg, bjArg, myThid)
       ELSE IF (fluxid.EQ.'AFR') THEN
         CALL LAYERS_FILL_FIELD(df, trIdentity, Nr, layers_afr,'M',
     &           kLev, nLevs, bibjFlg, biArg, bjArg, myThid)
       ELSE IF (fluxid.EQ.'TOT') THEN
         CALL LAYERS_FILL_FIELD(df, trIdentity, Nr, layers_tottend,'M',
     &           kLev, nLevs, bibjFlg, biArg, bjArg, myThid)
       ELSE
         WRITE(msgBuf,'(2A)')
     &          'S/R LAYERS_FILL: ',
     &          'invalid flux ID'
         CALL PRINT_ERROR( msgBuf, myThid )
         STOP 'ABNORMAL END: S/R LAYERS_FILL'
       ENDIF

       ELSE
C  ---- Cannot throw an error for different trIdentity
C       because subroutine also gets called for ptracers.
C       Just have to do nothing

C         WRITE(msgBuf,'(5A,I2)')
C     &          'S/R LAYERS_FILL: ',
C     &          'only works on THETA (1) or SALT (2)',
C     &          'fluxid=', fluxid, 'trId=',
C     &          trIdentity
C         CALL PRINT_ERROR( msgBuf, myThid )
C         STOP 'ABNORMAL END: S/R LAYERS_FILL'
       ENDIF

#endif /* LAYERS_THERMODYNAMICS */

      RETURN
      END


C end of S/R LAYERS_FILL SUBROUTINE LAYERS_FILL_FIELD( I df, trIdentity, myNr, U layers_saved_flux, I fldType, I kLev, nLevs, bibjFlg, biArg, bjArg, myThid ) IMPLICIT NONE #include "SIZE.h" #include "EEPARAMS.h" #include "PARAMS.h" #include "GRID.h" #include "LAYERS_SIZE.h" #include "LAYERS.h" INTEGER trIdentity, myNr, kLev, nLevs, bibjFlg, biArg, bjArg CHARACTER fldType _RL layers_saved_flux(1-OLx:sNx+OLx,1-OLy:sNy+OLy, & myNr,2,nSx,nSy) C _RL df(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) _RL df(*) INTEGER myThid #ifdef LAYERS_THERMODYNAMICS C !LOCAL VARIABLES: ==================================================== C i,j :: loop indices C msgBuf :: error message buffer INTEGER sizI1,sizI2,sizJ1,sizJ2 INTEGER sizTx,sizTy INTEGER iRun, jRun, k, bi, bj INTEGER kFirst, kLast INTEGER kd, kd0, ksgn C CHARACTER*(MAX_LEN_MBUF) msgBuf C- select range for 1rst & 2nd indices to accumulate C depending on variable location on C-grid, IF ( fldType.EQ.'U' ) THEN iRun = sNx+1 jRun = sNy ELSEIF ( fldType.EQ.'V' ) THEN iRun = sNx jRun = sNy+1 ELSE iRun = sNx jRun = sNy ENDIF C- Dimension of the input array: IF (abs(bibjFlg).EQ.3) THEN sizI1 = 1 sizI2 = sNx sizJ1 = 1 sizJ2 = sNy iRun = sNx jRun = sNy ELSE sizI1 = 1-OLx sizI2 = sNx+OLx sizJ1 = 1-OLy sizJ2 = sNy+OLy ENDIF IF (abs(bibjFlg).GE.2) THEN sizTx = 1 sizTy = 1 ELSE sizTx = nSx sizTy = nSy ENDIF C- Which part of inpFld to add : k = 3rd index, C and do the loop >> do k=kFirst,kLast << IF (kLev.LE.0) THEN kFirst = 1 kLast = nLevs ELSEIF ( nLevs.EQ.1 ) THEN kFirst = 1 kLast = 1 ELSEIF ( kLev.LE.nLevs ) THEN kFirst = kLev kLast = kLev ELSE STOP 'ABNORMAL END in LAYERS_SAVE: kLev > nLevs >0' ENDIF C- Which part of qdiag to update: kd = 3rd index, C and do the loop >> do k=kFirst,kLast ; kd = kd0 + k*ksgn << IF ( kLev.EQ.-1 ) THEN ksgn = -1 kd0 = 1 + nLevs ELSEIF ( kLev.EQ.0 ) THEN ksgn = 1 kd0 = 0 ELSE ksgn = 0 kd0 = kLev ENDIF IF ( bibjFlg.EQ.0 ) THEN DO bj=myByLo(myThid), myByHi(myThid) DO bi=myBxLo(myThid), myBxHi(myThid) DO k = kFirst,kLast kd = kd0 + ksgn*k CALL LAYERS_CUMULATE( U layers_saved_flux(1-OLx,1-OLy,kd,trIdentity,bi,bj), I df, I sizI1,sizI2,sizJ1,sizJ2,nLevs,sizTx,sizTy, I iRun,jRun,k,bi,bj, I myThid) ENDDO ENDDO ENDDO ELSE bi = MIN(biArg,sizTx) bj = MIN(bjArg,sizTy) DO k = kFirst,kLast kd = kd0 + ksgn*k CALL LAYERS_CUMULATE( U layers_saved_flux(1-OLx,1-OLy,kd,trIdentity,biArg,bjArg), I df, I sizI1,sizI2,sizJ1,sizJ2,nLevs,sizTx,sizTy, I iRun,jRun,k,bi,bj, I myThid) ENDDO ENDIF C IF (bibjFlg.EQ.2) THEN C C -- called INSIDE the bi-bj loop, with overlap present C DO k = kstart,kend C DO j = 1-OLy,sNy+OLy C DO i = 1-OLx,sNx+OLx C layers_saved_flux(i,j,k,trIdentity,biArg,bjArg) = C & layers_saved_flux(i,j,k,trIdentity,biArg,bjArg) + C & df(i,j,1,1) C ENDDO C ENDDO C ENDDO C ELSE IF (bibjFlg.EQ.0) THEN C C -- the bi-bj loop must be done here C DO bj=myByLo(myThid), myByHi(myThid) C DO bi=myBxLo(myThid), myBxHi(myThid) C DO k = kstart,kend C DO j = 1-OLy,sNy+OLy C DO i = 1-OLx,sNx+OLx C layers_saved_flux(i,j,k,trIdentity,bi,bj) = C & layers_saved_flux(i,j,k,trIdentity,bi,bj) + C & df(i,j,bi,bj) C ENDDO C ENDDO C ENDDO C ENDDO C ENDDO C ELSE C WRITE(msgBuf,'(2A)') C & 'S/R LAYERS_FILL_FIELD: ', C & 'got unexpected bibjFlg' C CALL PRINT_ERROR( msgBuf, myThid ) C STOP 'ABNORMAL END: S/R LAYERS_FILL_FIELD' C ENDIF #endif /* LAYERS_THERMODYNAMICS */ RETURN END


C end of S/R LAYERS_FILL_FIELD SUBROUTINE LAYERS_CUMULATE( U cumFld, I inpFld, I sizI1,sizI2,sizJ1,sizJ2,sizK,sizTx,sizTy, I iRun,jRun,k,bi,bj, I myThid ) C !DESCRIPTION: C Update array cumFld C by adding content of input field array inpFld C over the range [1:iRun],[1:jRun] C !USES: IMPLICIT NONE #include "EEPARAMS.h" #include "SIZE.h" C !INPUT/OUTPUT PARAMETERS: C == Routine Arguments == C cumFld :: cumulative array (updated) C inpFld :: input field array to add to cumFld C sizI1,sizI2 :: size of inpFld array: 1rst index range (min,max) C sizJ1,sizJ2 :: size of inpFld array: 2nd index range (min,max) C sizK :: size of inpFld array: 3rd dimension C sizTx,sizTy :: size of inpFld array: tile dimensions C iRun,jRun :: range of 1rst & 2nd index C k,bi,bj :: level and tile indices of inFld array to add to cumFld array C myThid :: my Thread Id number _RL cumFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy) INTEGER sizI1,sizI2,sizJ1,sizJ2 INTEGER sizK,sizTx,sizTy _RL inpFld(sizI1:sizI2,sizJ1:sizJ2,sizK,sizTx,sizTy) INTEGER iRun, jRun, k, bi, bj INTEGER myThid CEOP C !LOCAL VARIABLES: C i,j :: loop indices INTEGER i, j C _RL tmpFact C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| DO j = 1,jRun DO i = 1,iRun cumFld(i,j) = cumFld(i,j) + inpFld(i,j,k,bi,bj) ENDDO ENDDO RETURN END