pkg/generic_advdiff/gad_fluxlimit_adv_x.F

31566b6684 Alis*0001 #include "GAD_OPTIONS.h"
                0002
527a84022c Alis*0003 CBOP
                0004 C !ROUTINE: GAD_FLUXLIMIT_ADV_X
                0005
                0006 C !INTERFACE: ==========================================================
0af3073e4e Jean*0007       SUBROUTINE GAD_FLUXLIMIT_ADV_X(
8f93275e63 Jean*0008      I           bi, bj, k, calcCFL, deltaTloc,
0af3073e4e Jean*0009      I           uTrans, uFld,
38a970fb5d Jean*0010      I           maskLocW, tracer,
31566b6684 Alis*0011      O           uT,
                0012      I           myThid )
                0013
527a84022c Alis*0014 C !DESCRIPTION:
                0015 C Calculates the area integrated zonal flux due to advection of a tracer
                0016 C using second-order interpolation with a flux limiter:
                0017 C \begin{equation*}
0af3073e4e Jean*0018 C F^x_{adv} = U \overline{ \theta }^i
527a84022c Alis*0019 C - \frac{1}{2} \left(
                0020 C     [ 1 - \psi(C_r) ] |U|
                0021 C    + U \frac{u \Delta t}{\Delta x_c} \psi(C_r)
                0022 C              \right) \delta_i \theta
                0023 C \end{equation*}
                0024 C where the $\psi(C_r)$ is the limiter function and $C_r$ is
                0025 C the slope ratio.
                0026
                0027 C !USES: ===============================================================
                0028       IMPLICIT NONE
31566b6684 Alis*0029 #include "SIZE.h"
8f93275e63 Jean*0030 #include "EEPARAMS.h"
31566b6684 Alis*0031 #include "GRID.h"
                0032
527a84022c Alis*0033 C !INPUT PARAMETERS: ===================================================
8f93275e63 Jean*0034 C  bi, bj            :: tile indices
692dd30681 Jean*0035 C  k                 :: vertical level
                0036 C  calcCFL           :: =T: calculate CFL number ; =F: take uFld as CFL
                0037 C  deltaTloc         :: local time-step (s)
                0038 C  uTrans            :: zonal volume transport
                0039 C  uFld              :: zonal flow / CFL number
                0040 C  tracer            :: tracer field
                0041 C  myThid            :: thread number
8f93275e63 Jean*0042       INTEGER bi, bj, k
692dd30681 Jean*0043       LOGICAL calcCFL
1816bb0a69 Patr*0044       _RL deltaTloc
31566b6684 Alis*0045       _RL uTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0af3073e4e Jean*0046       _RL uFld  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
38a970fb5d Jean*0047       _RS maskLocW(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
31566b6684 Alis*0048       _RL tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0049       INTEGER myThid
                0050
527a84022c Alis*0051 C !OUTPUT PARAMETERS: ==================================================
692dd30681 Jean*0052 C  uT                :: zonal advective flux
527a84022c Alis*0053       _RL uT    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0054
                0055 C !LOCAL VARIABLES: ====================================================
8f93275e63 Jean*0056 C  i, j              :: loop indices
692dd30681 Jean*0057 C  Cr                :: slope ratio
8f93275e63 Jean*0058 C  Rjm, Rj, Rjp      :: differences at i-1,i,i+1
31566b6684 Alis*0059       INTEGER i,j
8f93275e63 Jean*0060       _RL Cr, Rjm, Rj, Rjp
692dd30681 Jean*0061       _RL uCFL
8f93275e63 Jean*0062       _RL CrMax
                0063       PARAMETER( CrMax = 1.D+6 )
                0064
527a84022c Alis*0065 C Statement function provides Limiter(Cr)
31566b6684 Alis*0066 #include "GAD_FLUX_LIMITER.h"
527a84022c Alis*0067 CEOP
31566b6684 Alis*0068
8f93275e63 Jean*0069       DO j=1-OLy,sNy+OLy
                0070        uT(1-OLx,j) = 0. _d 0
                0071        uT(2-OLx,j) = 0. _d 0
                0072        uT(sNx+OLx,j) = 0. _d 0
360ad14abb Mart*0073       ENDDO
8f93275e63 Jean*0074       DO j=1-OLy,sNy+OLy
                0075        DO i=1-OLx+2,sNx+OLx-1
2b22d4bdc3 Jean*0076
692dd30681 Jean*0077         uCFL = uFld(i,j)
                0078         IF ( calcCFL ) uCFL = ABS( uFld(i,j)*deltaTloc
c6bc48bdb5 Jean*0079      &                  *recip_dxC(i,j,bi,bj)*recip_deepFacC(k) )
8f93275e63 Jean*0080         Rjp = (tracer(i+1,j)-tracer( i ,j))*maskLocW(i+1,j)
                0081         Rj  = (tracer( i ,j)-tracer(i-1,j))*maskLocW( i ,j)
                0082         Rjm = (tracer(i-1,j)-tracer(i-2,j))*maskLocW(i-1,j)
2b22d4bdc3 Jean*0083
8f93275e63 Jean*0084         IF ( uTrans(i,j).GT.zeroRL ) THEN
                0085           Cr = Rjm
                0086         ELSE
                0087           Cr = Rjp
                0088         ENDIF
                0089         IF ( ABS(Rj)*CrMax .LE. ABS(Cr) ) THEN
                0090           Cr = SIGN( CrMax, Cr )*SIGN( oneRL, Rj )
31566b6684 Alis*0091         ELSE
8f93275e63 Jean*0092           Cr = Cr/Rj
31566b6684 Alis*0093         ENDIF
8f93275e63 Jean*0094
                0095 C       calculate Limiter Function:
                0096         Cr = Limiter(Cr)
                0097
0af3073e4e Jean*0098         uT(i,j) =
8f93275e63 Jean*0099      &     uTrans(i,j)*(Tracer(i,j)+Tracer(i-1,j))*0.5 _d 0
                0100      &   - ABS(uTrans(i,j))*( (oneRL-Cr) + uCFL*Cr )
                0101      &                     *Rj*0.5 _d 0
31566b6684 Alis*0102        ENDDO
                0103       ENDDO
                0104
                0105       RETURN
                0106       END