/[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.42 by jmc, Sun Nov 19 22:04:09 2006 UTC revision 1.48 by jmc, Tue Sep 18 15:27:56 2007 UTC
# Line 2  C $Header$ Line 2  C $Header$
2  C $Name$  C $Name$
3    
4  #include "GAD_OPTIONS.h"  #include "GAD_OPTIONS.h"
5    #ifdef ALLOW_PTRACERS
6    # include "PTRACERS_OPTIONS.h"
7    #endif
8    
9  CBOP  CBOP
10  C !ROUTINE: GAD_CALC_RHS  C !ROUTINE: GAD_CALC_RHS
# Line 44  C !USES: =============================== Line 47  C !USES: ===============================
47  #include "GRID.h"  #include "GRID.h"
48  #include "SURFACE.h"  #include "SURFACE.h"
49  #include "GAD.h"  #include "GAD.h"
50    #ifdef ALLOW_PTRACERS
51    # include "PTRACERS_SIZE.h"
52    # include "PTRACERS.h"
53    #endif
54    
55  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
56  #include "tamc.h"  #include "tamc.h"
# Line 75  C tracerIdentity   :: tracer identifier Line 82  C tracerIdentity   :: tracer identifier
82  C advectionScheme  :: advection scheme to use (Horizontal plane)  C advectionScheme  :: advection scheme to use (Horizontal plane)
83  C vertAdvecScheme  :: advection scheme to use (Vertical direction)  C vertAdvecScheme  :: advection scheme to use (Vertical direction)
84  C calcAdvection    :: =False if Advec computed with multiDim scheme  C calcAdvection    :: =False if Advec computed with multiDim scheme
85    C--- needs to pass those 2 flags as argument
86    C trUseGMRedi      :: true if this tracer uses GM-Redi
87    C trUseKPP         :: true if this tracer uses KPP
88    C---
89  C implicitAdvection:: =True if vertical Advec computed implicitly  C implicitAdvection:: =True if vertical Advec computed implicitly
90  C applyAB_onTracer :: apply Adams-Bashforth on Tracer (rather than on gTr)  C applyAB_onTracer :: apply Adams-Bashforth on Tracer (rather than on gTr)
91  C myTime           :: current time  C myTime           :: current time
# Line 98  C myThid           :: thread number Line 109  C myThid           :: thread number
109        _RL TracAB (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)        _RL TracAB (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
110        INTEGER tracerIdentity        INTEGER tracerIdentity
111        INTEGER advectionScheme, vertAdvecScheme        INTEGER advectionScheme, vertAdvecScheme
112        LOGICAL calcAdvection        LOGICAL calcAdvection, trUseGMRedi, trUseKPP
113        LOGICAL implicitAdvection, applyAB_onTracer        LOGICAL implicitAdvection, applyAB_onTracer
114        _RL     myTime        _RL     myTime
115        INTEGER myIter, myThid        INTEGER myIter, myThid
# Line 140  C--   the kDown is still required Line 151  C--   the kDown is still required
151        fVerT(1,1,kDown) = fVerT(1,1,kDown)        fVerT(1,1,kDown) = fVerT(1,1,kDown)
152  #endif  #endif
153    
154          trUseGMRedi = useGMRedi
155          trUseKPP    = useKPP
156    #ifdef ALLOW_PTRACERS
157          IF ( usePTRACERS ) THEN
158            i = MAX( 1, tracerIdentity - GAD_TR1 + 1 )
159            trUseGMRedi = useGMRedi .AND. ( tracerIdentity.LT.GAD_TR1
160         &                                  .OR. PTRACERS_useGMRedi(i) )
161            trUseKPP    = useKPP    .AND. ( tracerIdentity.LT.GAD_TR1
162         &                                  .OR. PTRACERS_useKPP(i)    )
163          ENDIF
164    #endif
165    
166  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
167  C--   Set diagnostic suffix for the current tracer  C--   Set diagnostic suffix for the current tracer
168        IF ( useDiagnostics ) THEN        IF ( useDiagnostics ) THEN
# Line 210  C-    Advective flux in X Line 233  C-    Advective flux in X
233            CALL GAD_C2_ADV_X(bi,bj,k,uTrans,locABT,af,myThid)            CALL GAD_C2_ADV_X(bi,bj,k,uTrans,locABT,af,myThid)
234          ELSEIF ( advectionScheme.EQ.ENUM_UPWIND_1RST          ELSEIF ( advectionScheme.EQ.ENUM_UPWIND_1RST
235       &          .OR. advectionScheme.EQ.ENUM_DST2 ) THEN       &          .OR. advectionScheme.EQ.ENUM_DST2 ) THEN
236            CALL GAD_DST2U1_ADV_X( bi,bj,k, advectionScheme,            CALL GAD_DST2U1_ADV_X( bi,bj,k, advectionScheme, .TRUE.,
237       I            dTtracerLev(k), uTrans, uFld, locABT,       I            dTtracerLev(k), uTrans, uFld, locABT,
238       O            af, myThid )       O            af, myThid )
239          ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN          ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN
240            CALL GAD_FLUXLIMIT_ADV_X( bi,bj,k, dTtracerLev(k),            CALL GAD_FLUXLIMIT_ADV_X( bi,bj,k, .TRUE., dTtracerLev(k),
241       I            uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,       I            uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,
242       O            af, myThid )       O            af, myThid )
243          ELSEIF (advectionScheme.EQ.ENUM_UPWIND_3RD ) THEN          ELSEIF (advectionScheme.EQ.ENUM_UPWIND_3RD ) THEN
# Line 222  C-    Advective flux in X Line 245  C-    Advective flux in X
245          ELSEIF (advectionScheme.EQ.ENUM_CENTERED_4TH) THEN          ELSEIF (advectionScheme.EQ.ENUM_CENTERED_4TH) THEN
246            CALL GAD_C4_ADV_X(bi,bj,k,uTrans,locABT,af,myThid)            CALL GAD_C4_ADV_X(bi,bj,k,uTrans,locABT,af,myThid)
247          ELSEIF (advectionScheme.EQ.ENUM_DST3 ) THEN          ELSEIF (advectionScheme.EQ.ENUM_DST3 ) THEN
248            CALL GAD_DST3_ADV_X( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3_ADV_X( bi,bj,k, .TRUE., dTtracerLev(k),
249       I            uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,       I            uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,
250       O            af, myThid )       O            af, myThid )
251          ELSEIF (advectionScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN          ELSEIF (advectionScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN
# Line 230  C-    Advective flux in X Line 253  C-    Advective flux in X
253  cph This block is to trick the adjoint:  cph This block is to trick the adjoint:
254  cph IF inAdExact=.FALSE., we want to use DST3  cph IF inAdExact=.FALSE., we want to use DST3
255  cph with limiters in forward, but without limiters in reverse.  cph with limiters in forward, but without limiters in reverse.
256            CALL GAD_DST3_ADV_X( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3_ADV_X( bi,bj,k, .TRUE., dTtracerLev(k),
257       I           uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,       I           uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,
258       O           af, myThid )       O           af, myThid )
259           ELSE           ELSE
260            CALL GAD_DST3FL_ADV_X( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3FL_ADV_X( bi,bj,k, .TRUE., dTtracerLev(k),
261       I           uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,       I           uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,
262       O           af, myThid )       O           af, myThid )
263           ENDIF           ENDIF
264    #ifndef ALLOW_AUTODIFF_TAMC
265            ELSEIF (advectionScheme.EQ.ENUM_OS7MP ) THEN
266              CALL GAD_OS7MP_ADV_X( bi,bj,k, .TRUE., dTtracerLev(k),
267         I            uTrans, uFld, maskW(1-Olx,1-Oly,k,bi,bj), locABT,
268         O            af, myThid )
269    #endif
270          ELSE          ELSE
271           STOP 'GAD_CALC_RHS: Bad advectionScheme (X)'           STOP 'GAD_CALC_RHS: Bad advectionScheme (X)'
272          ENDIF          ENDIF
# Line 272  C-    Add bi-harmonic diffusive flux in Line 301  C-    Add bi-harmonic diffusive flux in
301    
302  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
303  C-    GM/Redi flux in X  C-    GM/Redi flux in X
304        IF (useGMRedi) THEN        IF ( trUseGMRedi ) THEN
305  C *note* should update GMREDI_XTRANSPORT to set df  *aja*  C *note* should update GMREDI_XTRANSPORT to set df  *aja*
306          IF ( applyAB_onTracer ) THEN          IF ( applyAB_onTracer ) THEN
307            CALL GMREDI_XTRANSPORT(            CALL GMREDI_XTRANSPORT(
# Line 289  C *note* should update GMREDI_XTRANSPORT Line 318  C *note* should update GMREDI_XTRANSPORT
318          ENDIF          ENDIF
319        ENDIF        ENDIF
320  #endif  #endif
321    C     anelastic: advect.fluxes are scaled by rhoFac but hor.diff. flx are not
322        DO j=1-Oly,sNy+Oly        DO j=1-Oly,sNy+Oly
323         DO i=1-Olx,sNx+Olx         DO i=1-Olx,sNx+Olx
324          fZon(i,j) = fZon(i,j) + df(i,j)          fZon(i,j) = fZon(i,j) + df(i,j)*rhoFacC(k)
325         ENDDO         ENDDO
326        ENDDO        ENDDO
327    
# Line 299  C *note* should update GMREDI_XTRANSPORT Line 329  C *note* should update GMREDI_XTRANSPORT
329  C-    Diagnostics of Tracer flux in X dir (mainly Diffusive term),  C-    Diagnostics of Tracer flux in X dir (mainly Diffusive term),
330  C       excluding advective terms:  C       excluding advective terms:
331        IF ( useDiagnostics .AND.        IF ( useDiagnostics .AND.
332       &    (diffKh.NE.0. .OR. diffK4 .NE.0. .OR. useGMRedi) ) THEN       &    (diffKh.NE.0. .OR. diffK4 .NE.0. .OR. trUseGMRedi) ) THEN
333            diagName = 'DFxE'//diagSufx            diagName = 'DFxE'//diagSufx
334            CALL DIAGNOSTICS_FILL(df,diagName, k,1, 2,bi,bj, myThid)            CALL DIAGNOSTICS_FILL(df,diagName, k,1, 2,bi,bj, myThid)
335        ENDIF        ENDIF
# Line 318  C-    Advective flux in Y Line 348  C-    Advective flux in Y
348            CALL GAD_C2_ADV_Y(bi,bj,k,vTrans,locABT,af,myThid)            CALL GAD_C2_ADV_Y(bi,bj,k,vTrans,locABT,af,myThid)
349          ELSEIF ( advectionScheme.EQ.ENUM_UPWIND_1RST          ELSEIF ( advectionScheme.EQ.ENUM_UPWIND_1RST
350       &          .OR. advectionScheme.EQ.ENUM_DST2 ) THEN       &          .OR. advectionScheme.EQ.ENUM_DST2 ) THEN
351            CALL GAD_DST2U1_ADV_Y( bi,bj,k, advectionScheme,            CALL GAD_DST2U1_ADV_Y( bi,bj,k, advectionScheme, .TRUE.,
352       I            dTtracerLev(k), vTrans, vFld, locABT,       I            dTtracerLev(k), vTrans, vFld, locABT,
353       O            af, myThid )       O            af, myThid )
354          ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN          ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN
355            CALL GAD_FLUXLIMIT_ADV_Y( bi,bj,k, dTtracerLev(k),            CALL GAD_FLUXLIMIT_ADV_Y( bi,bj,k, .TRUE., dTtracerLev(k),
356       I            vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,       I            vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,
357       O            af, myThid )       O            af, myThid )
358          ELSEIF (advectionScheme.EQ.ENUM_UPWIND_3RD ) THEN          ELSEIF (advectionScheme.EQ.ENUM_UPWIND_3RD ) THEN
# Line 330  C-    Advective flux in Y Line 360  C-    Advective flux in Y
360          ELSEIF (advectionScheme.EQ.ENUM_CENTERED_4TH) THEN          ELSEIF (advectionScheme.EQ.ENUM_CENTERED_4TH) THEN
361            CALL GAD_C4_ADV_Y(bi,bj,k,vTrans,locABT,af,myThid)            CALL GAD_C4_ADV_Y(bi,bj,k,vTrans,locABT,af,myThid)
362          ELSEIF (advectionScheme.EQ.ENUM_DST3 ) THEN          ELSEIF (advectionScheme.EQ.ENUM_DST3 ) THEN
363            CALL GAD_DST3_ADV_Y( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3_ADV_Y( bi,bj,k, .TRUE., dTtracerLev(k),
364       I            vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,       I            vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,
365       O            af, myThid )       O            af, myThid )
366          ELSEIF (advectionScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN          ELSEIF (advectionScheme.EQ.ENUM_DST3_FLUX_LIMIT ) THEN
# Line 338  C-    Advective flux in Y Line 368  C-    Advective flux in Y
368  cph This block is to trick the adjoint:  cph This block is to trick the adjoint:
369  cph IF inAdExact=.FALSE., we want to use DST3  cph IF inAdExact=.FALSE., we want to use DST3
370  cph with limiters in forward, but without limiters in reverse.  cph with limiters in forward, but without limiters in reverse.
371            CALL GAD_DST3_ADV_Y( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3_ADV_Y( bi,bj,k, .TRUE., dTtracerLev(k),
372       I           vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,       I           vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,
373       O           af, myThid )       O           af, myThid )
374           ELSE           ELSE
375            CALL GAD_DST3FL_ADV_Y( bi,bj,k, dTtracerLev(k),            CALL GAD_DST3FL_ADV_Y( bi,bj,k, .TRUE., dTtracerLev(k),
376       I           vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,       I           vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,
377       O           af, myThid )       O           af, myThid )
378           ENDIF           ENDIF
379    #ifndef ALLOW_AUTODIFF_TAMC
380            ELSEIF (advectionScheme.EQ.ENUM_OS7MP ) THEN
381              CALL GAD_OS7MP_ADV_Y( bi,bj,k, .TRUE., dTtracerLev(k),
382         I            vTrans, vFld, maskS(1-Olx,1-Oly,k,bi,bj), locABT,
383         O            af, myThid )
384    #endif
385          ELSE          ELSE
386            STOP 'GAD_CALC_RHS: Bad advectionScheme (Y)'            STOP 'GAD_CALC_RHS: Bad advectionScheme (Y)'
387          ENDIF          ENDIF
# Line 380  C-    Add bi-harmonic flux in Y Line 416  C-    Add bi-harmonic flux in Y
416    
417  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
418  C-    GM/Redi flux in Y  C-    GM/Redi flux in Y
419        IF (useGMRedi) THEN        IF ( trUseGMRedi ) THEN
420  C *note* should update GMREDI_YTRANSPORT to set df  *aja*  C *note* should update GMREDI_YTRANSPORT to set df  *aja*
421          IF ( applyAB_onTracer ) THEN          IF ( applyAB_onTracer ) THEN
422            CALL GMREDI_YTRANSPORT(            CALL GMREDI_YTRANSPORT(
# Line 397  C *note* should update GMREDI_YTRANSPORT Line 433  C *note* should update GMREDI_YTRANSPORT
433          ENDIF          ENDIF
434        ENDIF        ENDIF
435  #endif  #endif
436    C     anelastic: advect.fluxes are scaled by rhoFac but hor.diff. flx are not
437        DO j=1-Oly,sNy+Oly        DO j=1-Oly,sNy+Oly
438         DO i=1-Olx,sNx+Olx         DO i=1-Olx,sNx+Olx
439          fMer(i,j) = fMer(i,j) + df(i,j)          fMer(i,j) = fMer(i,j) + df(i,j)*rhoFacC(k)
440         ENDDO         ENDDO
441        ENDDO        ENDDO
442    
# Line 407  C *note* should update GMREDI_YTRANSPORT Line 444  C *note* should update GMREDI_YTRANSPORT
444  C-    Diagnostics of Tracer flux in Y dir (mainly Diffusive terms),  C-    Diagnostics of Tracer flux in Y dir (mainly Diffusive terms),
445  C       excluding advective terms:  C       excluding advective terms:
446        IF ( useDiagnostics .AND.        IF ( useDiagnostics .AND.
447       &    (diffKh.NE.0. .OR. diffK4 .NE.0. .OR. useGMRedi) ) THEN       &    (diffKh.NE.0. .OR. diffK4 .NE.0. .OR. trUseGMRedi) ) THEN
448            diagName = 'DFyE'//diagSufx            diagName = 'DFyE'//diagSufx
449            CALL DIAGNOSTICS_FILL(df,diagName, k,1, 2,bi,bj, myThid)            CALL DIAGNOSTICS_FILL(df,diagName, k,1, 2,bi,bj, myThid)
450        ENDIF        ENDIF
# Line 456  cph with limiters in forward, but withou Line 493  cph with limiters in forward, but withou
493       I         dTtracerLev(k),rTrans,wFld,TracAB(1-Olx,1-Oly,1,bi,bj),       I         dTtracerLev(k),rTrans,wFld,TracAB(1-Olx,1-Oly,1,bi,bj),
494       O         af, myThid )       O         af, myThid )
495            ENDIF            ENDIF
496    #ifndef ALLOW_AUTODIFF_TAMC
497            ELSEIF (vertAdvecScheme.EQ.ENUM_OS7MP ) THEN
498               CALL GAD_OS7MP_ADV_R( bi,bj,k,
499         I         dTtracerLev(k),rTrans,wFld,TracAB(1-Olx,1-Oly,1,bi,bj),
500         O         af, myThid )
501    #endif
502          ELSE          ELSE
503            STOP 'GAD_CALC_RHS: Bad vertAdvecScheme (R)'            STOP 'GAD_CALC_RHS: Bad vertAdvecScheme (R)'
504          ENDIF          ENDIF
# Line 495  C           boundary condition. Line 538  C           boundary condition.
538    
539  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
540  C-    GM/Redi flux in R  C-    GM/Redi flux in R
541        IF (useGMRedi) THEN        IF ( trUseGMRedi ) THEN
542  C *note* should update GMREDI_RTRANSPORT to set df  *aja*  C *note* should update GMREDI_RTRANSPORT to set df  *aja*
543          IF ( applyAB_onTracer ) THEN          IF ( applyAB_onTracer ) THEN
544            CALL GMREDI_RTRANSPORT(            CALL GMREDI_RTRANSPORT(
# Line 523  C *note* should update GMREDI_RTRANSPORT Line 566  C *note* should update GMREDI_RTRANSPORT
566  C-    Diagnostics of Tracer flux in R dir (mainly Diffusive terms),  C-    Diagnostics of Tracer flux in R dir (mainly Diffusive terms),
567  C       Explicit terms only & excluding advective terms:  C       Explicit terms only & excluding advective terms:
568        IF ( useDiagnostics .AND.        IF ( useDiagnostics .AND.
569       &    (.NOT.implicitDiffusion .OR. useGMRedi) ) THEN       &    (.NOT.implicitDiffusion .OR. trUseGMRedi) ) THEN
570            diagName = 'DFrE'//diagSufx            diagName = 'DFrE'//diagSufx
571            CALL DIAGNOSTICS_FILL(df,diagName, k,1, 2,bi,bj, myThid)            CALL DIAGNOSTICS_FILL(df,diagName, k,1, 2,bi,bj, myThid)
572        ENDIF        ENDIF
# Line 531  C       Explicit terms only & excluding Line 574  C       Explicit terms only & excluding
574    
575  #ifdef ALLOW_KPP  #ifdef ALLOW_KPP
576  C-    Set non local KPP transport term (ghat):  C-    Set non local KPP transport term (ghat):
577        IF ( useKPP .AND. k.GE.2 ) THEN        IF ( trUseKPP .AND. k.GE.2 ) THEN
578         DO j=1-Oly,sNy+Oly         DO j=1-Oly,sNy+Oly
579          DO i=1-Olx,sNx+Olx          DO i=1-Olx,sNx+Olx
580           df(i,j) = 0. _d 0           df(i,j) = 0. _d 0
# Line 539  C-    Set non local KPP transport term ( Line 582  C-    Set non local KPP transport term (
582         ENDDO         ENDDO
583         IF (tracerIdentity.EQ.GAD_TEMPERATURE) THEN         IF (tracerIdentity.EQ.GAD_TEMPERATURE) THEN
584          CALL KPP_TRANSPORT_T(          CALL KPP_TRANSPORT_T(
585       I     iMin,iMax,jMin,jMax,bi,bj,k,km1,       I           iMin,iMax,jMin,jMax,bi,bj,k,km1,
586       O     df )       O           df,
587         I           myTime, myIter, myThid )
588         ELSEIF (tracerIdentity.EQ.GAD_SALINITY) THEN         ELSEIF (tracerIdentity.EQ.GAD_SALINITY) THEN
589          CALL KPP_TRANSPORT_S(          CALL KPP_TRANSPORT_S(
590       I     iMin,iMax,jMin,jMax,bi,bj,k,km1,       I           iMin,iMax,jMin,jMax,bi,bj,k,km1,
591       O     df )       O           df,
592         I           myTime, myIter, myThid )
593  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
594         ELSEIF (tracerIdentity .GE. GAD_TR1) THEN         ELSEIF (tracerIdentity .GE. GAD_TR1) THEN
595          CALL KPP_TRANSPORT_PTR(          CALL KPP_TRANSPORT_PTR(
596       I     iMin,iMax,jMin,jMax,bi,bj,k,km1,       I           iMin,iMax,jMin,jMax,bi,bj,k,km1,
597       I     tracerIdentity-GAD_TR1+1,       I           tracerIdentity-GAD_TR1+1,
598       O     df )       O           df,
599         I           myTime, myIter, myThid )
600  #endif  #endif
601         ELSE         ELSE
602          PRINT*,'invalid tracer indentity: ', tracerIdentity          PRINT*,'invalid tracer indentity: ', tracerIdentity
# Line 558  C-    Set non local KPP transport term ( Line 604  C-    Set non local KPP transport term (
604         ENDIF         ENDIF
605         DO j=1-Oly,sNy+Oly         DO j=1-Oly,sNy+Oly
606          DO i=1-Olx,sNx+Olx          DO i=1-Olx,sNx+Olx
607           fVerT(i,j,kUp) = fVerT(i,j,kUp) + df(i,j)*maskUp(i,j)           fVerT(i,j,kUp) = fVerT(i,j,kUp)
608         &                  + df(i,j)*maskUp(i,j)*rhoFacF(k)
609          ENDDO          ENDDO
610         ENDDO         ENDDO
611        ENDIF        ENDIF
612  #endif  #endif
613    
614  C--   Divergence of fluxes  C--   Divergence of fluxes
615    C     Anelastic: scale vertical fluxes by rhoFac and leave Horizontal fluxes unchanged
616        DO j=1-Oly,sNy+Oly-1        DO j=1-Oly,sNy+Oly-1
617         DO i=1-Olx,sNx+Olx-1         DO i=1-Olx,sNx+Olx-1
618          gTracer(i,j,k,bi,bj)=gTracer(i,j,k,bi,bj)          gTracer(i,j,k,bi,bj)=gTracer(i,j,k,bi,bj)
619       &   -_recip_hFacC(i,j,k,bi,bj)*recip_drF(k)*recip_rA(i,j,bi,bj)       &   -_recip_hFacC(i,j,k,bi,bj)*recip_drF(k)
620         &   *recip_rA(i,j,bi,bj)*recip_deepFac2C(k)*recip_rhoFacC(k)
621       &   *( (fZon(i+1,j)-fZon(i,j))       &   *( (fZon(i+1,j)-fZon(i,j))
622       &     +(fMer(i,j+1)-fMer(i,j))       &     +(fMer(i,j+1)-fMer(i,j))
623       &     +(fVerT(i,j,kDown)-fVerT(i,j,kUp))*rkSign       &     +(fVerT(i,j,kDown)-fVerT(i,j,kUp))*rkSign
Legend:
Removed from v.1.42  
changed lines
  Added in v.1.48
  ViewVC Help
Powered by ViewVC 1.1.22