/[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.15 by heimbach, Tue Nov 12 20:42:24 2002 UTC revision 1.20 by jmc, Mon Sep 22 22:13:11 2003 UTC
# Line 41  C !USES: =============================== Line 41  C !USES: ===============================
41  #include "PARAMS.h"  #include "PARAMS.h"
42  #include "GRID.h"  #include "GRID.h"
43  #include "DYNVARS.h"  #include "DYNVARS.h"
44    #include "SURFACE.h"
45  #include "GAD.h"  #include "GAD.h"
46    #ifdef ALLOW_PTRACERS
47    #include "PTRACERS_OPTIONS.h"
48    #include "PTRACERS.h"
49    #endif
50    
51  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
52  #include "tamc.h"  #include "tamc.h"
# Line 61  C  diffKh               :: horizontal di Line 66  C  diffKh               :: horizontal di
66  C  diffK4               :: bi-harmonic diffusion coefficient  C  diffK4               :: bi-harmonic diffusion coefficient
67  C  KappaRT              :: 3-D array for vertical diffusion coefficient  C  KappaRT              :: 3-D array for vertical diffusion coefficient
68  C  Tracer               :: tracer field  C  Tracer               :: tracer field
69  C  tracerIdentity       :: identifier for the tracer (required only for KPP)  C  tracerIdentity       :: identifier for the tracer (required for KPP and GM)
70  C  advectionScheme      :: advection scheme to use  C  advectionScheme      :: advection scheme to use
71  C  calcAdvection        :: =False if Advec terms computed with multiDim scheme  C  calcAdvection        :: =False if Advec terms computed with multiDim scheme
72  C  myThid               :: thread number  C  myThid               :: thread number
# Line 130  C--   Make local copy of tracer array Line 135  C--   Make local copy of tracer array
135    
136  C--   Unless we have already calculated the advection terms we initialize  C--   Unless we have already calculated the advection terms we initialize
137  C     the tendency to zero.  C     the tendency to zero.
138        IF (calcAdvection) THEN  C     <== now done earlier at the beginning of thermodynamics.
139         DO j=1-Oly,sNy+Oly  c     IF (calcAdvection) THEN
140          DO i=1-Olx,sNx+Olx  c      DO j=1-Oly,sNy+Oly
141           gTracer(i,j,k,bi,bj)=0. _d 0  c       DO i=1-Olx,sNx+Olx
142          ENDDO  c        gTracer(i,j,k,bi,bj)=0. _d 0
143         ENDDO  c       ENDDO
144        ENDIF  c      ENDDO
145    c     ENDIF
146    
147  C--   Pre-calculate del^2 T if bi-harmonic coefficient is non-zero  C--   Pre-calculate del^2 T if bi-harmonic coefficient is non-zero
148        IF (diffK4 .NE. 0.) THEN        IF (diffK4 .NE. 0.) THEN
# Line 289  C-    Bi-harmonic flux in Y Line 295  C-    Bi-harmonic flux in Y
295         ENDDO         ENDDO
296        ENDIF        ENDIF
297    
298    #ifdef NONLIN_FRSURF
299    C--   Compute vertical flux fVerT(kDown) at interface k+1 (between k & k+1):
300          IF ( calcAdvection .AND. K.EQ.Nr .AND.
301         &     useRealFreshWaterFlux .AND.
302         &     buoyancyRelation .EQ. 'OCEANICP' ) THEN  
303           DO j=1-Oly,sNy+Oly
304            DO i=1-Olx,sNx+Olx
305             fVerT(i,j,kDown) = convertEmP2rUnit*PmEpR(i,j,bi,bj)
306         &     *rA(i,j,bi,bj)*maskC(i,j,k,bi,bj)*Tracer(i,j,k,bi,bj)
307            ENDDO
308           ENDDO
309          ENDIF
310    #endif /* NONLIN_FRSURF */
311    
312    C--   Compute vertical flux fVerT(kUp) at interface k (between k-1 & k):
313  C-    Advective flux in R  C-    Advective flux in R
314        IF (calcAdvection) THEN        IF (calcAdvection) THEN
315  C     Note: wVel needs to be masked  C     Note: wVel needs to be masked
# Line 386  C *note* should update KPP_TRANSPORT_T t Line 407  C *note* should update KPP_TRANSPORT_T t
407       I     iMin,iMax,jMin,jMax,bi,bj,k,km1,       I     iMin,iMax,jMin,jMax,bi,bj,k,km1,
408       I     KappaRT,       I     KappaRT,
409       U     df )       U     df )
410    #ifdef ALLOW_PTRACERS
411           ELSEIF (tracerIdentity .GE. GAD_TR1 .AND.
412         &         tracerIdentity .LE. (GAD_TR1+PTRACERS_numInUse-1)) THEN
413            CALL KPP_TRANSPORT_PTR(
414         I     iMin,iMax,jMin,jMax,bi,bj,k,km1,
415         I     tracerIdentity,KappaRT,
416         U     df )
417    #endif
418         ELSE         ELSE
419            PRINT*,'invalid tracer indentity: ', tracerIdentity
420          STOP 'GAD_CALC_RHS: Ooops'          STOP 'GAD_CALC_RHS: Ooops'
421         ENDIF         ENDIF
422         DO j=1-Oly,sNy+Oly         DO j=1-Oly,sNy+Oly
# Line 411  C--   Divergence of fluxes Line 441  C--   Divergence of fluxes
441         ENDDO         ENDDO
442        ENDDO        ENDDO
443    
444    #ifdef NONLIN_FRSURF
445    C-- account for 3.D divergence of the flow in rStar coordinate:
446          IF (calcAdvection .AND. select_rStar.GT.0) THEN
447           DO j=1-Oly,sNy+Oly-1
448            DO i=1-Olx,sNx+Olx-1
449             gTracer(i,j,k,bi,bj) = gTracer(i,j,k,bi,bj)
450         &     - (rStarExpC(i,j,bi,bj) - 1. _d 0)/deltaTfreesurf
451         &       *tracer(i,j,k,bi,bj)*maskC(i,j,k,bi,bj)
452            ENDDO
453           ENDDO
454          ENDIF
455          IF (calcAdvection .AND. select_rStar.LT.0) THEN
456           DO j=1-Oly,sNy+Oly-1
457            DO i=1-Olx,sNx+Olx-1
458             gTracer(i,j,k,bi,bj) = gTracer(i,j,k,bi,bj)
459         &     - rStarDhCDt(i,j,bi,bj)
460         &       *tracer(i,j,k,bi,bj)*maskC(i,j,k,bi,bj)
461            ENDDO
462           ENDDO
463          ENDIF
464    #endif /* NONLIN_FRSURF */
465          
466    
467        RETURN        RETURN
468        END        END
Legend:
Removed from v.1.15  
changed lines
  Added in v.1.20
  ViewVC Help
Powered by ViewVC 1.1.22