C $Header: /u/gcmpack/MITgcm/pkg/generic_advdiff/salt_fill.F,v 1.4 2006/06/29 16:14:36 dimitri Exp $
C $Name:  $

#include "GAD_OPTIONS.h"

C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
CBOP
C !ROUTINE: SALT_FILL
C !INTERFACE: ==========================================================
      SUBROUTINE SALT_FILL(
     .     uVel, vVel,
     .     salt,
     .     flag,myTime,myIter,myThid)
C
C !DESCRIPTION:
C Fills in negatives for the salt (specific humidity) field
C
C The algorithm is as follows:
C 
C Simplest scheme (flag = 1) -> borrow from below and create
C                               salt if needed at bottom level
C 'Get it back'   (flag = 2) -> Fill negative of salt by getting it 
C                               back from where it went
C     If no immediate surrounding value is large enough to fill negative, 
C     the sum of immediate surrounding positive values is tried.  
C     If sum is not large enough, salt is simply set to zero.
C  NOTE AS OF 6/2/06 -- DO NOT USE FLAG=2 OPTION - NOT WORKING
C
C !INPUT PARAMETERS: ===================================================
C  uVel              :: zonal velocity component
C  vVel              :: meridional velocity component
C  salt              :: salt field
C  flag              :: integer flag telling scheme how to fill
C  myTime            :: current time
C  myIter            :: iteration number
C  myThid            :: thread number
C !OUTPUT PARAMETERS: ==================================================
C  salt             :: salt array is replaced 
CEOP
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
      IMPLICIT NONE
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
      _RL uVel  (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
      _RL vVel  (1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
      _RL salt(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
      INTEGER flag
      _RL myTime
      INTEGER myIter
      INTEGER myThid

#ifdef ALLOW_DIAGNOSTICS
      logical  diagnostics_is_on
      external 
#endif

      INTEGER bi,bj,i,j,L,LM1
      _RL dpratio, scale

      DO bj=myByLo(myThid),myByHi(myThid)
       DO bi=myBxLo(myThid),myBxHi(myThid)
C Fill diagnostic for filling with negative of salt
#ifdef ALLOW_DIAGNOSTICS
        IF (useDiagnostics)THEN
        if(diagnostics_is_on('SALTFILL',myThid) ) then
        scale = -1. _d 0
        CALL DIAGNOSTICS_SCALE_FILL(salt,scale,1,'SALTFILL',
     .                                   0,Nr,-1,bi,bj,myThid)
        endif
        ENDIF
#endif

c Flag = 1: 
c ---------------------------------
        if(flag.eq.1) then

        do L=Nr,2,-1
         LM1 = L-1
         dpratio= rC(L)/rC(LM1)
         do j=1,sNy
         do i=1,sNx
          if( salt(i,j,L,bi,bj).lt.0.0  _d 0) then
           salt(i,j,LM1,bi,bj) = salt(i,j,LM1,bi,bj) +
     .          salt(i,j,L,bi,bj)*dpratio
           salt(i,j,L,bi,bj) = 0.0 _d 0
          endif
         enddo
         enddo
        enddo

        do j=1,sNy
        do i=1,sNx
         if(salt(i,j,1,bi,bj).lt.0.0 _d 0)
     .                  salt(i,j,1,bi,bj) = 0.0 _d 0
        enddo
        enddo

        else
         print *,'Invalid Flag in salt_fill - nothing done '
        endif

C Fill diagnostic for filling with salt - get tendency
#ifdef ALLOW_DIAGNOSTICS
        IF (useDiagnostics)THEN
        if(diagnostics_is_on('SALTFILL',myThid) ) then
        CALL DIAGNOSTICS_FILL(salt,'SALTFILL',0,Nr,1,bi,bj,myThid)
        endif
        ENDIF
#endif
       ENDDO
      ENDDO

      RETURN
      END