pkg/generic_advdiff/gad_dst2u1_adv_x.F

76b8386f6e Jean*0001 #include "GAD_OPTIONS.h"
                0002
                0003 CBOP
                0004 C !ROUTINE: GAD_DST2U1_ADV_X
                0005
                0006 C !INTERFACE: ==========================================================
0ef5a3ea3b Jean*0007       SUBROUTINE GAD_DST2U1_ADV_X(
692dd30681 Jean*0008      I           bi,bj,k, advectionScheme, calcCFL,
0af3073e4e Jean*0009      I           deltaTloc, uTrans, uFld,
76b8386f6e Jean*0010      I           tracer,
                0011      O           uT,
                0012      I           myThid )
                0013
                0014 C !DESCRIPTION:
                0015 C  Calculates the area integrated zonal flux due to advection
0ef5a3ea3b Jean*0016 C  of a tracer using second-order Direct Space and Time (DST-2)
76b8386f6e Jean*0017 C  interpolation (=Lax-Wendroff) or simple 1rst order upwind scheme.
                0018
                0019 C !USES: ===============================================================
                0020       IMPLICIT NONE
                0021 #include "SIZE.h"
                0022 #include "GRID.h"
                0023 #include "GAD.h"
                0024
                0025 C !INPUT PARAMETERS: ===================================================
                0026 C  bi,bj             :: tile indices
                0027 C  k                 :: vertical level
0ef5a3ea3b Jean*0028 C  advectionScheme   :: advection scheme to use: either 2nd Order DST
76b8386f6e Jean*0029 C                                                or 1rst Order Upwind
692dd30681 Jean*0030 C  calcCFL           :: =T: calculate CFL number ; =F: take uFld as CFL
                0031 C  deltaTloc         :: local time-step (s)
76b8386f6e Jean*0032 C  uTrans            :: zonal volume transport
692dd30681 Jean*0033 C  uFld              :: zonal flow / CFL number
76b8386f6e Jean*0034 C  tracer            :: tracer field
                0035 C  myThid            :: thread number
                0036       INTEGER bi,bj,k
                0037       INTEGER advectionScheme
692dd30681 Jean*0038       LOGICAL calcCFL
76b8386f6e Jean*0039       _RL deltaTloc
                0040       _RL uTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
0af3073e4e Jean*0041       _RL uFld  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
76b8386f6e Jean*0042       _RL tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0043       INTEGER myThid
                0044
                0045 C !OUTPUT PARAMETERS: ==================================================
                0046 C  uT                :: zonal advective flux
                0047       _RL uT    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0048
                0049 C !LOCAL VARIABLES: ====================================================
                0050 C  i,j               :: loop indices
                0051 C  rLimit            :: centered (vs upwind) fraction
                0052 C  uCFL              :: Courant-Friedrich-Levy number
                0053       INTEGER i,j
692dd30681 Jean*0054       _RL uCFL, xLimit, uAbs
76b8386f6e Jean*0055 CEOP
                0056
                0057       xLimit = 0. _d 0
                0058       IF ( advectionScheme.EQ.ENUM_DST2 ) xLimit = 1. _d 0
                0059
                0060       DO j=1-Oly,sNy+Oly
                0061        uT(1-Olx,j)=0.
360ad14abb Mart*0062       ENDDO
                0063       DO j=1-Oly,sNy+Oly
76b8386f6e Jean*0064        DO i=1-Olx+1,sNx+Olx
                0065
692dd30681 Jean*0066         uCFL = uFld(i,j)
                0067         IF ( calcCFL ) uCFL = ABS( uFld(i,j)*deltaTloc
62743db5c5 Jean*0068      &                  *recip_dxC(i,j,bi,bj)*recip_deepFacC(k) )
76b8386f6e Jean*0069
0ef5a3ea3b Jean*0070 c       uT(i,j) =
                0071 c    &     uTrans(i,j)*(tracer(i-1,j)+tracer(i,j))*0.5 _d 0
                0072 c    &   + ( 1. _d 0 - xLimit*(1. _d 0 - uCFL) )*ABS(uTrans(i,j))
                0073 c    &                *(tracer(i-1,j)-tracer(i,j))*0.5 _d 0
                0074 C-- above formulation produces large truncation error when:
                0075 C    1rst.O upWind and   u > 0 & |tracer(i-1,j)| << |tracer(i,j)|
                0076 C                   or   u < 0 & |tracer(i-1,j)| >> |tracer(i,j)|
                0077 C-- change to a more robust expression:
                0078         uAbs = ABS(uTrans(i,j))
                0079      &       *( 1. _d 0 - xLimit*(1. _d 0 - uCFL) )
                0080         uT(i,j) = ( uTrans(i,j)+uAbs )* 0.5 _d 0 * tracer(i-1,j)
                0081      &          + ( uTrans(i,j)-uAbs )* 0.5 _d 0 * tracer(i,j)
76b8386f6e Jean*0082        ENDDO
                0083       ENDDO
                0084
                0085       RETURN
                0086       END