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

Contents of /MITgcm/pkg/generic_advdiff/gad_dst3_adv_y.F

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


Revision 1.10 - (show annotations) (download)
Sun Oct 22 01:08:04 2006 UTC (17 years, 6 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint58r_post
Changes since 1.9: +40 -24 lines
use Samar's version as the default (instead of only if useMatrix=T)

1 C $Header: /u/gcmpack/MITgcm/pkg/generic_advdiff/gad_dst3_adv_y.F,v 1.9 2006/06/19 14:40:43 jmc Exp $
2 C $Name: $
3
4 #include "GAD_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: GAD_DST3_ADV_Y
8
9 C !INTERFACE: ==========================================================
10 SUBROUTINE GAD_DST3_ADV_Y(
11 I bi,bj,k,deltaTloc,
12 I vTrans, vFld,
13 I maskLocS, tracer,
14 O vT,
15 I myThid )
16 C !DESCRIPTION:
17 C Calculates the area integrated Meridional flux due to advection of a
18 C tracer using 3rd-order Direct Space and Time (DST-3) Advection Scheme
19
20 C !USES: ===============================================================
21 IMPLICIT NONE
22
23 C == GLobal variables ==
24 #include "SIZE.h"
25 #include "GRID.h"
26 #include "EEPARAMS.h"
27 #include "PARAMS.h"
28 #include "GAD.h"
29
30 C == Routine arguments ==
31 C !INPUT PARAMETERS: ===================================================
32 C bi,bj :: tile indices
33 C k :: vertical level
34 C deltaTloc :: local time-step (s)
35 C vTrans :: meridional volume transport
36 C vFld :: meridional flow
37 C tracer :: tracer field
38 C myThid :: thread number
39 INTEGER bi,bj,k
40 _RL deltaTloc
41 _RL vTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
42 _RL vFld (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
43 _RS maskLocS(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
44 _RL tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
45 INTEGER myThid
46
47 C !OUTPUT PARAMETERS: ==================================================
48 C vT :: meridional advective flux
49 _RL vT (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
50
51 C == Local variables ==
52 C !LOCAL VARIABLES: ====================================================
53 C i,j :: loop indices
54 C vLoc :: velocity [m/s], meridional component
55 C cfl :: Courant-Friedrich-Levy number
56 INTEGER i,j
57 _RL Rjm,Rj,Rjp,cfl,d0,d1
58 _RL vLoc
59 #ifdef OLD_DST3_FORMULATION
60 _RL psiP,psiM,thetaP,thetaM
61 _RL smallNo
62 c _RL Rjjm,Rjjp
63
64 IF (inAdMode) THEN
65 smallNo = 1.0D-20
66 ELSE
67 smallNo = 1.0D-20
68 ENDIF
69 #endif
70
71 DO i=1-Olx,sNx+Olx
72 vT(i,1-Oly)=0.
73 vT(i,2-Oly)=0.
74 vT(i,sNy+Oly)=0.
75 ENDDO
76 DO j=1-Oly+2,sNy+Oly-1
77 DO i=1-Olx,sNx+Olx
78 Rjp=(tracer(i,j+1)-tracer(i, j ))*maskLocS(i,j+1)
79 Rj =(tracer(i, j )-tracer(i,j-1))*maskLocS(i, j )
80 Rjm=(tracer(i,j-1)-tracer(i,j-2))*maskLocS(i,j-1)
81
82 vLoc = vFld(i,j)
83 c vLoc = vTrans(i,j)*recip_dxG(i,j,bi,bj)
84 c & *recip_drF(k)*_recip_hFacS(i,j,k,bi,bj)
85 cfl=ABS(vLoc*deltaTloc*recip_dyC(i,j,bi,bj))
86 d0=(2.-cfl)*(1.-cfl)*oneSixth
87 d1=(1.-cfl*cfl)*oneSixth
88 #ifdef OLD_DST3_FORMULATION
89 IF ( ABS(Rj).LT.smallNo .OR.
90 & ABS(Rjm).LT.smallNo ) THEN
91 thetaP=0.
92 psiP=0.
93 ELSE
94 thetaP=(Rjm+smallNo)/(smallNo+Rj)
95 psiP=d0+d1*thetaP
96 ENDIF
97 IF ( ABS(Rj).LT.smallNo .OR.
98 & ABS(Rjp).LT.smallNo ) THEN
99 thetaM=0.
100 psiM=0.
101 ELSE
102 thetaM=(Rjp+smallNo)/(smallNo+Rj)
103 psiM=d0+d1*thetaM
104 ENDIF
105 vT(i,j)=
106 & 0.5*(vTrans(i,j)+ABS(vTrans(i,j)))
107 & *( Tracer(i,j-1) + psiP*Rj )
108 & +0.5*(vTrans(i,j)-ABS(vTrans(i,j)))
109 & *( Tracer(i, j ) - psiM*Rj )
110 #else /* OLD_DST3_FORMULATION */
111 vT(i,j)=
112 & 0.5*(vTrans(i,j)+ABS(vTrans(i,j)))
113 & *( Tracer(i,j-1) + (d0*Rj+d1*Rjm) )
114 & +0.5*(vTrans(i,j)-ABS(vTrans(i,j)))
115 & *( Tracer(i, j ) - (d0*Rj+d1*Rjp) )
116 #endif /* OLD_DST3_FORMULATION */
117
118 ENDDO
119 ENDDO
120
121 RETURN
122 END

  ViewVC Help
Powered by ViewVC 1.1.22