/[MITgcm]/MITgcm/pkg/generic_advdiff/gad_dst3_adv_r.F
ViewVC logotype

Diff of /MITgcm/pkg/generic_advdiff/gad_dst3_adv_r.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.3 by jmc, Wed Mar 6 01:29:36 2002 UTC revision 1.4 by heimbach, Sun Apr 3 16:05:34 2005 UTC
# Line 38  C     wFld     :: velocity, vertical com Line 38  C     wFld     :: velocity, vertical com
38        _RL Rjm,Rj,Rjp,cfl,d0,d1        _RL Rjm,Rj,Rjp,cfl,d0,d1
39        _RL psiP,psiM,thetaP,thetaM        _RL psiP,psiM,thetaP,thetaM
40        _RL wFld        _RL wFld
41          _RL smallNo
42          _RL Rjjm,Rjjp
43    
44        IF (.NOT. multiDimAdvection) THEN        IF (.NOT. multiDimAdvection) THEN
45  C      If using the standard time-stepping/advection schemes (ie. AB-II)  C      If using the standard time-stepping/advection schemes (ie. AB-II)
# Line 52  C      for maskC(...) and wVel(...) Line 54  C      for maskC(...) and wVel(...)
54         bj=1         bj=1
55        ENDIF        ENDIF
56    
57          IF (inAdMode) THEN
58           smallNo = 1.0D-20
59          ELSE
60           smallNo = 1.0D-20
61          ENDIF
62    
63        km2=MAX(1,k-2)        km2=MAX(1,k-2)
64        km1=MAX(1,k-1)        km1=MAX(1,k-1)
65        kp1=MIN(Nr,k+1)        kp1=MIN(Nr,k+1)
# Line 70  c       wFld = wVel(i,j,k,bi_arg,bj_arg) Line 78  c       wFld = wVel(i,j,k,bi_arg,bj_arg)
78          cfl=abs(wFld*dTarg*recip_drC(k))          cfl=abs(wFld*dTarg*recip_drC(k))
79          d0=(2.-cfl)*(1.-cfl)*oneSixth          d0=(2.-cfl)*(1.-cfl)*oneSixth
80          d1=(1.-cfl*cfl)*oneSixth          d1=(1.-cfl*cfl)*oneSixth
81  c       thetaP=0.          IF ( ABS(Rj).LT.smallNo .OR.
82  c       IF (Rj.NE.0.) thetaP=Rjm/Rj       &       ABS(Rjm).LT.smallNo ) THEN
83          thetaP=Rjm/(1.D-20+Rj)           thetaP=0.
84          psiP=d0+d1*thetaP           psiP=0.
85  c       psiP=max(0.,min(min(1.,psiP),(1.-cfl)/(1.D-20+cfl)*thetaP))          ELSE
86          thetaM=Rjp/(1.D-20+Rj)           thetaP=(Rjm+smallNo)/(smallNo+Rj)
87  c       thetaM=0.           psiP=d0+d1*thetaP
88  c       IF (Rj.NE.0.) thetaM=Rjp/Rj          ENDIF
89          psiM=d0+d1*thetaM          IF ( ABS(Rj).LT.smallNo .OR.
90  c       psiM=max(0.,min(min(1.,psiM),(1.-cfl)/(1.D-20+cfl)*thetaM))       &       ABS(Rjp).LT.smallNo ) THEN
91          wT(i,j)=           thetaM=0.
92       &   0.5*(rTrans(i,j)+abs(rTrans(i,j)))           psiM=0.
93       &      *( Tracer(i,j, k ,bi,bj) + psiM*Rj )          ELSE
94       &  +0.5*(rTrans(i,j)-abs(rTrans(i,j)))           thetaM=(Rjp+smallNo)/(smallNo+Rj)
95       &      *( Tracer(i,j,km1,bi,bj) - psiP*Rj )           psiM=d0+d1*thetaM
96            ENDIF
97             wT(i,j)=
98         &    0.5*(rTrans(i,j)+abs(rTrans(i,j)))
99         &       *( Tracer(i,j, k ,bi,bj) + psiM*Rj )
100         &   +0.5*(rTrans(i,j)-abs(rTrans(i,j)))
101         &       *( Tracer(i,j,km1,bi,bj) - psiP*Rj )
102         ENDDO         ENDDO
103        ENDDO        ENDDO
104    
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.4
  ViewVC Help
Powered by ViewVC 1.1.22