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

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

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

revision 1.31 by jmc, Sat Dec 4 00:20:27 2004 UTC revision 1.32 by jmc, Thu Dec 16 22:29:34 2004 UTC
# Line 107  C !LOCAL VARIABLES: ==================== Line 107  C !LOCAL VARIABLES: ====================
107  C i,j              :: loop indices  C i,j              :: loop indices
108  C df4              :: used for storing del^2 T for bi-harmonic term  C df4              :: used for storing del^2 T for bi-harmonic term
109  C fZon             :: zonal flux  C fZon             :: zonal flux
110  C fmer             :: meridional flux  C fMer             :: meridional flux
111  C af               :: advective flux  C af               :: advective flux
112  C df               :: diffusive flux  C df               :: diffusive flux
113  C localT           :: local copy of tracer field  C localT           :: local copy of tracer field
114    #ifdef ALLOW_DIAGNOSTICS
115          INTEGER     kk
116          CHARACTER*8 diagName
117          CHARACTER*4 GAD_DIAG_SUFX, diagSufx
118          EXTERNAL    GAD_DIAG_SUFX
119    #endif
120        INTEGER i,j        INTEGER i,j
121        _RL df4   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL df4   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
122        _RL fZon  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL fZon  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
# Line 127  C--   the kDown is still required Line 133  C--   the kDown is still required
133        fVerT(1,1,kDown) = fVerT(1,1,kDown)        fVerT(1,1,kDown) = fVerT(1,1,kDown)
134  #endif  #endif
135    
136    #ifdef ALLOW_DIAGNOSTICS
137    C--   Set diagnostic suffix for the current tracer
138          IF ( useDiagnostics ) THEN
139            diagSufx = GAD_DIAG_SUFX( tracerIdentity, myThid )
140          ENDIF
141    #endif
142    
143        advFac  = 0. _d 0        advFac  = 0. _d 0
144        IF (calcAdvection) advFac = 1. _d 0        IF (calcAdvection) advFac = 1. _d 0
145        rAdvFac = rkFac*advFac        rAdvFac = rkFac*advFac
# Line 176  C--   Initialize net flux in X direction Line 189  C--   Initialize net flux in X direction
189    
190  C-    Advective flux in X  C-    Advective flux in X
191        IF (calcAdvection) THEN        IF (calcAdvection) THEN
192        IF (advectionScheme.EQ.ENUM_CENTERED_2ND) THEN          IF (advectionScheme.EQ.ENUM_CENTERED_2ND) THEN
193         CALL GAD_C2_ADV_X(bi,bj,k,uTrans,localT,af,myThid)            CALL GAD_C2_ADV_X(bi,bj,k,uTrans,localT,af,myThid)
194        ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN          ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN
195         CALL GAD_FLUXLIMIT_ADV_X( bi,bj,k, dTtracerLev(k),            CALL GAD_FLUXLIMIT_ADV_X( bi,bj,k, dTtracerLev(k),
196       I          uTrans, uVel, maskW(1-Olx,1-Oly,k,bi,bj), localT,       I            uTrans, uVel, maskW(1-Olx,1-Oly,k,bi,bj), localT,
197       O          af, myThid )       O            af, myThid )
198        ELSEIF (advectionScheme.EQ.ENUM_UPWIND_3RD ) THEN          ELSEIF (advectionScheme.EQ.ENUM_UPWIND_3RD ) THEN
199         CALL GAD_U3_ADV_X(bi,bj,k,uTrans,localT,af,myThid)            CALL GAD_U3_ADV_X(bi,bj,k,uTrans,localT,af,myThid)
200        ELSEIF (advectionScheme.EQ.ENUM_CENTERED_4TH) THEN          ELSEIF (advectionScheme.EQ.ENUM_CENTERED_4TH) THEN
201         CALL GAD_C4_ADV_X(bi,bj,k,uTrans,localT,af,myThid)            CALL GAD_C4_ADV_X(bi,bj,k,uTrans,localT,af,myThid)
202        ELSEIF (advectionScheme.EQ.ENUM_DST3 ) THEN          ELSEIF (advectionScheme.EQ.ENUM_DST3 ) THEN
203         CALL GAD_DST3_ADV_X( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3_ADV_X( bi,bj,k, dTtracerLev(k),
204       I          uTrans, uVel, maskW(1-Olx,1-Oly,k,bi,bj), localT,       I            uTrans, uVel, maskW(1-Olx,1-Oly,k,bi,bj), localT,
205       O          af, myThid )       O            af, myThid )
206        ELSEIF (advectionScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN          ELSEIF (advectionScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN
207         CALL GAD_DST3FL_ADV_X( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3FL_ADV_X( bi,bj,k, dTtracerLev(k),
208       I          uTrans, uVel, maskW(1-Olx,1-Oly,k,bi,bj), localT,       I            uTrans, uVel, maskW(1-Olx,1-Oly,k,bi,bj), localT,
209       O          af, myThid )       O            af, myThid )
210        ELSE          ELSE
211         STOP 'GAD_CALC_RHS: Bad advectionScheme (X)'           STOP 'GAD_CALC_RHS: Bad advectionScheme (X)'
212        ENDIF          ENDIF
213        DO j=1-Oly,sNy+Oly          DO j=1-Oly,sNy+Oly
214         DO i=1-Olx,sNx+Olx           DO i=1-Olx,sNx+Olx
215          fZon(i,j) = fZon(i,j) + af(i,j)            fZon(i,j) = fZon(i,j) + af(i,j)
216         ENDDO           ENDDO
217        ENDDO          ENDDO
218    #ifdef ALLOW_DIAGNOSTICS
219            IF ( useDiagnostics ) THEN
220              diagName = 'ADVx'//diagSufx
221              kk = -k
222              CALL DIAGNOSTICS_FILL(af,diagName, kk,1, 2,bi,bj, myThid)
223            ENDIF
224    #endif
225        ENDIF        ENDIF
226    
227  C-    Diffusive flux in X  C-    Diffusive flux in X
# Line 215  C-    Diffusive flux in X Line 235  C-    Diffusive flux in X
235         ENDDO         ENDDO
236        ENDIF        ENDIF
237    
238    C-    Add bi-harmonic diffusive flux in X
239          IF (diffK4 .NE. 0.) THEN
240           CALL GAD_BIHARM_X(bi,bj,k,xA,df4,diffK4,df,myThid)
241          ENDIF
242    
243  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
244  C-    GM/Redi flux in X  C-    GM/Redi flux in X
245        IF (useGMRedi) THEN        IF (useGMRedi) THEN
# Line 232  C *note* should update GMREDI_XTRANSPORT Line 257  C *note* should update GMREDI_XTRANSPORT
257         ENDDO         ENDDO
258        ENDDO        ENDDO
259    
260  C-    Bi-harmonic duffusive flux in X  #ifdef ALLOW_DIAGNOSTICS
261        IF (diffK4 .NE. 0.) THEN  C-    Diagnostics of Tracer flux in X dir (mainly Diffusive term),
262         CALL GAD_BIHARM_X(bi,bj,k,xA,df4,diffK4,df,myThid)  C       excluding advective terms:
263         DO j=1-Oly,sNy+Oly        IF ( useDiagnostics .AND.
264          DO i=1-Olx,sNx+Olx       &    (diffKh.NE.0. .OR. diffK4 .NE.0. .OR. useGMRedi) ) THEN
265           fZon(i,j) = fZon(i,j) + df(i,j)            diagName = 'DIFx'//diagSufx
266          ENDDO            kk = -k
267         ENDDO            CALL DIAGNOSTICS_FILL(df,diagName, kk,1, 2,bi,bj, myThid)
268        ENDIF        ENDIF
269    #endif
270    
271  C--   Initialize net flux in Y direction  C--   Initialize net flux in Y direction
272        DO j=1-Oly,sNy+Oly        DO j=1-Oly,sNy+Oly
# Line 251  C--   Initialize net flux in Y direction Line 277  C--   Initialize net flux in Y direction
277    
278  C-    Advective flux in Y  C-    Advective flux in Y
279        IF (calcAdvection) THEN        IF (calcAdvection) THEN
280        IF (advectionScheme.EQ.ENUM_CENTERED_2ND) THEN          IF (advectionScheme.EQ.ENUM_CENTERED_2ND) THEN
281         CALL GAD_C2_ADV_Y(bi,bj,k,vTrans,localT,af,myThid)            CALL GAD_C2_ADV_Y(bi,bj,k,vTrans,localT,af,myThid)
282        ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN          ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN
283         CALL GAD_FLUXLIMIT_ADV_Y( bi,bj,k, dTtracerLev(k),            CALL GAD_FLUXLIMIT_ADV_Y( bi,bj,k, dTtracerLev(k),
284       I          vTrans, vVel, maskS(1-Olx,1-Oly,k,bi,bj), localT,       I            vTrans, vVel, maskS(1-Olx,1-Oly,k,bi,bj), localT,
285       O          af, myThid )       O            af, myThid )
286        ELSEIF (advectionScheme.EQ.ENUM_UPWIND_3RD ) THEN          ELSEIF (advectionScheme.EQ.ENUM_UPWIND_3RD ) THEN
287         CALL GAD_U3_ADV_Y(bi,bj,k,vTrans,localT,af,myThid)            CALL GAD_U3_ADV_Y(bi,bj,k,vTrans,localT,af,myThid)
288        ELSEIF (advectionScheme.EQ.ENUM_CENTERED_4TH) THEN          ELSEIF (advectionScheme.EQ.ENUM_CENTERED_4TH) THEN
289         CALL GAD_C4_ADV_Y(bi,bj,k,vTrans,localT,af,myThid)            CALL GAD_C4_ADV_Y(bi,bj,k,vTrans,localT,af,myThid)
290        ELSEIF (advectionScheme.EQ.ENUM_DST3 ) THEN          ELSEIF (advectionScheme.EQ.ENUM_DST3 ) THEN
291         CALL GAD_DST3_ADV_Y( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3_ADV_Y( bi,bj,k, dTtracerLev(k),
292       I          vTrans, vVel, maskS(1-Olx,1-Oly,k,bi,bj), localT,       I            vTrans, vVel, maskS(1-Olx,1-Oly,k,bi,bj), localT,
293       O          af, myThid )       O            af, myThid )
294        ELSEIF (advectionScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN          ELSEIF (advectionScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN
295         CALL GAD_DST3FL_ADV_Y( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3FL_ADV_Y( bi,bj,k, dTtracerLev(k),
296       I          vTrans, vVel, maskS(1-Olx,1-Oly,k,bi,bj), localT,       I            vTrans, vVel, maskS(1-Olx,1-Oly,k,bi,bj), localT,
297       O          af, myThid )       O            af, myThid )
298        ELSE          ELSE
299         STOP 'GAD_CALC_RHS: Bad advectionScheme (Y)'            STOP 'GAD_CALC_RHS: Bad advectionScheme (Y)'
300        ENDIF          ENDIF
301        DO j=1-Oly,sNy+Oly          DO j=1-Oly,sNy+Oly
302         DO i=1-Olx,sNx+Olx           DO i=1-Olx,sNx+Olx
303          fMer(i,j) = fMer(i,j) + af(i,j)            fMer(i,j) = fMer(i,j) + af(i,j)
304         ENDDO           ENDDO
305        ENDDO          ENDDO
306    #ifdef ALLOW_DIAGNOSTICS
307            IF ( useDiagnostics ) THEN
308              diagName = 'ADVy'//diagSufx
309              kk = -k
310              CALL DIAGNOSTICS_FILL(af,diagName, kk,1, 2,bi,bj, myThid)
311            ENDIF
312    #endif
313        ENDIF        ENDIF
314    
315  C-    Diffusive flux in Y  C-    Diffusive flux in Y
# Line 290  C-    Diffusive flux in Y Line 323  C-    Diffusive flux in Y
323         ENDDO         ENDDO
324        ENDIF        ENDIF
325    
326    C-    Add bi-harmonic flux in Y
327          IF (diffK4 .NE. 0.) THEN
328           CALL GAD_BIHARM_Y(bi,bj,k,yA,df4,diffK4,df,myThid)
329          ENDIF
330    
331  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
332  C-    GM/Redi flux in Y  C-    GM/Redi flux in Y
333        IF (useGMRedi) THEN        IF (useGMRedi) THEN
# Line 307  C *note* should update GMREDI_YTRANSPORT Line 345  C *note* should update GMREDI_YTRANSPORT
345         ENDDO         ENDDO
346        ENDDO        ENDDO
347    
348  C-    Bi-harmonic flux in Y  #ifdef ALLOW_DIAGNOSTICS
349        IF (diffK4 .NE. 0.) THEN  C-    Diagnostics of Tracer flux in Y dir (mainly Diffusive terms),
350         CALL GAD_BIHARM_Y(bi,bj,k,yA,df4,diffK4,df,myThid)  C       excluding advective terms:
351         DO j=1-Oly,sNy+Oly        IF ( useDiagnostics .AND.
352          DO i=1-Olx,sNx+Olx       &    (diffKh.NE.0. .OR. diffK4 .NE.0. .OR. useGMRedi) ) THEN
353           fMer(i,j) = fMer(i,j) + df(i,j)            diagName = 'DIFy'//diagSufx
354          ENDDO            kk = -k
355         ENDDO            CALL DIAGNOSTICS_FILL(df,diagName, kk,1, 2,bi,bj, myThid)
356        ENDIF        ENDIF
357    #endif
358    
359  C--   Compute vertical flux fVerT(kUp) at interface k (between k-1 & k):  C--   Compute vertical flux fVerT(kUp) at interface k (between k-1 & k):
360  C-    Advective flux in R  C-    Advective flux in R
# Line 328  C- a hack to prevent Water-Vapor vert.tr Line 367  C- a hack to prevent Water-Vapor vert.tr
367        IF (calcAdvection .AND. .NOT.implicitAdvection .AND. K.GE.2) THEN        IF (calcAdvection .AND. .NOT.implicitAdvection .AND. K.GE.2) THEN
368  #endif  #endif
369  C-    Compute vertical advective flux in the interior:  C-    Compute vertical advective flux in the interior:
370         IF (vertAdvecScheme.EQ.ENUM_CENTERED_2ND) THEN          IF (vertAdvecScheme.EQ.ENUM_CENTERED_2ND) THEN
371          CALL GAD_C2_ADV_R(bi,bj,k,rTrans,tracer,af,myThid)            CALL GAD_C2_ADV_R(bi,bj,k,rTrans,tracer,af,myThid)
372         ELSEIF (vertAdvecScheme.EQ.ENUM_FLUX_LIMIT) THEN          ELSEIF (vertAdvecScheme.EQ.ENUM_FLUX_LIMIT) THEN
373          CALL GAD_FLUXLIMIT_ADV_R(            CALL GAD_FLUXLIMIT_ADV_R(
374       &       bi,bj,k,dTtracerLev(k),rTrans,wVel,tracer,af,myThid)       &         bi,bj,k,dTtracerLev(k),rTrans,wVel,tracer,af,myThid)
375         ELSEIF (vertAdvecScheme.EQ.ENUM_UPWIND_3RD ) THEN          ELSEIF (vertAdvecScheme.EQ.ENUM_UPWIND_3RD ) THEN
376          CALL GAD_U3_ADV_R(bi,bj,k,rTrans,tracer,af,myThid)            CALL GAD_U3_ADV_R(bi,bj,k,rTrans,tracer,af,myThid)
377         ELSEIF (vertAdvecScheme.EQ.ENUM_CENTERED_4TH) THEN          ELSEIF (vertAdvecScheme.EQ.ENUM_CENTERED_4TH) THEN
378          CALL GAD_C4_ADV_R(bi,bj,k,rTrans,tracer,af,myThid)            CALL GAD_C4_ADV_R(bi,bj,k,rTrans,tracer,af,myThid)
379         ELSEIF (vertAdvecScheme.EQ.ENUM_DST3 ) THEN          ELSEIF (vertAdvecScheme.EQ.ENUM_DST3 ) THEN
380          CALL GAD_DST3_ADV_R(            CALL GAD_DST3_ADV_R(
381       &       bi,bj,k,dTtracerLev(k),rTrans,wVel,tracer,af,myThid)       &         bi,bj,k,dTtracerLev(k),rTrans,wVel,tracer,af,myThid)
382         ELSEIF (vertAdvecScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN          ELSEIF (vertAdvecScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN
383          CALL GAD_DST3FL_ADV_R(            CALL GAD_DST3FL_ADV_R(
384       &       bi,bj,k,dTtracerLev(k),rTrans,wVel,tracer,af,myThid)       &         bi,bj,k,dTtracerLev(k),rTrans,wVel,tracer,af,myThid)
385         ELSE          ELSE
386          STOP 'GAD_CALC_RHS: Bad vertAdvecScheme (R)'            STOP 'GAD_CALC_RHS: Bad vertAdvecScheme (R)'
387         ENDIF          ENDIF
388  C-     add the advective flux to fVerT  C-     add the advective flux to fVerT
389         DO j=1-Oly,sNy+Oly          DO j=1-Oly,sNy+Oly
390          DO i=1-Olx,sNx+Olx           DO i=1-Olx,sNx+Olx
391           fVerT(i,j,kUp) = fVerT(i,j,kUp) + af(i,j)            fVerT(i,j,kUp) = fVerT(i,j,kUp) + af(i,j)
392             ENDDO
393          ENDDO          ENDDO
394         ENDDO  #ifdef ALLOW_DIAGNOSTICS
395            IF ( useDiagnostics ) THEN
396              diagName = 'ADVr'//diagSufx
397              kk = -k
398              CALL DIAGNOSTICS_FILL(af,diagName, kk,1, 2,bi,bj, myThid)
399            ENDIF
400    #endif
401        ENDIF        ENDIF
402    
403  C-    Diffusive flux in R  C-    Diffusive flux in R
# Line 385  C *note* should update GMREDI_RTRANSPORT Line 431  C *note* should update GMREDI_RTRANSPORT
431         ENDDO         ENDDO
432        ENDDO        ENDDO
433    
434    #ifdef ALLOW_DIAGNOSTICS
435    C-    Diagnostics of Tracer flux in R dir (mainly Diffusive terms),
436    C       Explicit terms only & excluding advective terms:
437          IF ( useDiagnostics .AND.
438         &    (.NOT.implicitDiffusion .OR. useGMRedi) ) THEN
439              diagName = 'DFrE'//diagSufx
440              kk = -k
441              CALL DIAGNOSTICS_FILL(df,diagName, kk,1, 2,bi,bj, myThid)
442          ENDIF
443    #endif
444    
445  #ifdef ALLOW_KPP  #ifdef ALLOW_KPP
446  C-    Set non local KPP transport term (ghat):  C-    Set non local KPP transport term (ghat):
447        IF ( useKPP .AND. k.GE.2 ) THEN        IF ( useKPP .AND. k.GE.2 ) THEN
Legend:
Removed from v.1.31  
changed lines
  Added in v.1.32
  ViewVC Help
Powered by ViewVC 1.1.22