/[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.12.8.1 by heimbach, Sun Mar 24 04:11:34 2002 UTC revision 1.17 by jmc, Sun Jan 26 21:08:36 2003 UTC
# Line 11  C !INTERFACE: ========================== Line 11  C !INTERFACE: ==========================
11       I           bi,bj,iMin,iMax,jMin,jMax,k,kM1,kUp,kDown,       I           bi,bj,iMin,iMax,jMin,jMax,k,kM1,kUp,kDown,
12       I           xA,yA,uTrans,vTrans,rTrans,maskUp,       I           xA,yA,uTrans,vTrans,rTrans,maskUp,
13       I           diffKh, diffK4, KappaRT, Tracer,       I           diffKh, diffK4, KappaRT, Tracer,
14       I           tracerIdentity, advectionScheme,       I           tracerIdentity, advectionScheme, calcAdvection,
15       U           fVerT, gTracer,       U           fVerT, gTracer,
16       I           myThid )       I           myThid )
17    
# 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    
47  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 63  C  KappaRT              :: 3-D array for Line 64  C  KappaRT              :: 3-D array for
64  C  Tracer               :: tracer field  C  Tracer               :: tracer field
65  C  tracerIdentity       :: identifier for the tracer (required only for KPP)  C  tracerIdentity       :: identifier for the tracer (required only for KPP)
66  C  advectionScheme      :: advection scheme to use  C  advectionScheme      :: advection scheme to use
67    C  calcAdvection        :: =False if Advec terms computed with multiDim scheme
68  C  myThid               :: thread number  C  myThid               :: thread number
69        INTEGER bi,bj,iMin,iMax,jMin,jMax        INTEGER bi,bj,iMin,iMax,jMin,jMax
70        INTEGER k,kUp,kDown,kM1        INTEGER k,kUp,kDown,kM1
# Line 77  C  myThid               :: thread number Line 79  C  myThid               :: thread number
79        _RL Tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)        _RL Tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
80        INTEGER tracerIdentity        INTEGER tracerIdentity
81        INTEGER advectionScheme        INTEGER advectionScheme
82          LOGICAL calcAdvection
83        INTEGER myThid        INTEGER myThid
84    
85  C !OUTPUT PARAMETERS: ==================================================  C !OUTPUT PARAMETERS: ==================================================
# Line 128  C--   Make local copy of tracer array Line 131  C--   Make local copy of tracer array
131    
132  C--   Unless we have already calculated the advection terms we initialize  C--   Unless we have already calculated the advection terms we initialize
133  C     the tendency to zero.  C     the tendency to zero.
134        IF (.NOT. multiDimAdvection .OR.        IF (calcAdvection) THEN
      &    advectionScheme.EQ.ENUM_CENTERED_2ND .OR.  
      &    advectionScheme.EQ.ENUM_UPWIND_3RD .OR.  
      &    advectionScheme.EQ.ENUM_CENTERED_4TH ) THEN  
135         DO j=1-Oly,sNy+Oly         DO j=1-Oly,sNy+Oly
136          DO i=1-Olx,sNx+Olx          DO i=1-Olx,sNx+Olx
137           gTracer(i,j,k,bi,bj)=0. _d 0           gTracer(i,j,k,bi,bj)=0. _d 0
# Line 154  C--   Initialize net flux in X direction Line 154  C--   Initialize net flux in X direction
154        ENDDO        ENDDO
155    
156  C-    Advective flux in X  C-    Advective flux in X
157        IF (.NOT. multiDimAdvection .OR.        IF (calcAdvection) THEN
      &    advectionScheme.EQ.ENUM_CENTERED_2ND .OR.  
      &    advectionScheme.EQ.ENUM_UPWIND_3RD .OR.  
      &    advectionScheme.EQ.ENUM_CENTERED_4TH ) THEN  
158        IF (advectionScheme.EQ.ENUM_CENTERED_2ND) THEN        IF (advectionScheme.EQ.ENUM_CENTERED_2ND) THEN
159         CALL GAD_C2_ADV_X(bi,bj,k,uTrans,localT,af,myThid)         CALL GAD_C2_ADV_X(bi,bj,k,uTrans,localT,af,myThid)
160        ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN        ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN
# Line 200  C-    GM/Redi flux in X Line 197  C-    GM/Redi flux in X
197  C *note* should update GMREDI_XTRANSPORT to use localT and set df  *aja*  C *note* should update GMREDI_XTRANSPORT to use localT and set df  *aja*
198          CALL GMREDI_XTRANSPORT(          CALL GMREDI_XTRANSPORT(
199       I     iMin,iMax,jMin,jMax,bi,bj,K,       I     iMin,iMax,jMin,jMax,bi,bj,K,
200       I     xA,Tracer,       I     xA,Tracer,tracerIdentity,
201       U     df,       U     df,
202       I     myThid)       I     myThid)
203        ENDIF        ENDIF
# Line 229  C--   Initialize net flux in Y direction Line 226  C--   Initialize net flux in Y direction
226        ENDDO        ENDDO
227    
228  C-    Advective flux in Y  C-    Advective flux in Y
229        IF (.NOT. multiDimAdvection .OR.        IF (calcAdvection) THEN
      &    advectionScheme.EQ.ENUM_CENTERED_2ND .OR.  
      &    advectionScheme.EQ.ENUM_UPWIND_3RD .OR.  
      &    advectionScheme.EQ.ENUM_CENTERED_4TH ) THEN  
230        IF (advectionScheme.EQ.ENUM_CENTERED_2ND) THEN        IF (advectionScheme.EQ.ENUM_CENTERED_2ND) THEN
231         CALL GAD_C2_ADV_Y(bi,bj,k,vTrans,localT,af,myThid)         CALL GAD_C2_ADV_Y(bi,bj,k,vTrans,localT,af,myThid)
232        ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN        ELSEIF (advectionScheme.EQ.ENUM_FLUX_LIMIT) THEN
# Line 275  C-    GM/Redi flux in Y Line 269  C-    GM/Redi flux in Y
269  C *note* should update GMREDI_YTRANSPORT to use localT and set df  *aja*  C *note* should update GMREDI_YTRANSPORT to use localT and set df  *aja*
270         CALL GMREDI_YTRANSPORT(         CALL GMREDI_YTRANSPORT(
271       I     iMin,iMax,jMin,jMax,bi,bj,K,       I     iMin,iMax,jMin,jMax,bi,bj,K,
272       I     yA,Tracer,       I     yA,Tracer,tracerIdentity,
273       U     df,       U     df,
274       I     myThid)       I     myThid)
275        ENDIF        ENDIF
# Line 296  C-    Bi-harmonic flux in Y Line 290  C-    Bi-harmonic flux in Y
290         ENDDO         ENDDO
291        ENDIF        ENDIF
292    
293    #ifdef NONLIN_FRSURF
294    C--   Compute vertical flux fVerT(kDown) at interface k+1 (between k & k+1):
295          IF ( calcAdvection .AND. K.EQ.Nr .AND.
296         &     useRealFreshWaterFlux .AND.
297         &     buoyancyRelation .EQ. 'OCEANICP' ) THEN  
298           DO j=1-Oly,sNy+Oly
299            DO i=1-Olx,sNx+Olx
300             fVerT(i,j,kDown) = convertEmP2rUnit*PmEpR(i,j,bi,bj)
301         &     *rA(i,j,bi,bj)*maskC(i,j,k,bi,bj)*Tracer(i,j,k,bi,bj)
302            ENDDO
303           ENDDO
304          ENDIF
305    #endif /* NONLIN_FRSURF */
306    
307    C--   Compute vertical flux fVerT(kUp) at interface k (between k-1 & k):
308  C-    Advective flux in R  C-    Advective flux in R
309        IF (.NOT. multiDimAdvection .OR.        IF (calcAdvection) THEN
      &    advectionScheme.EQ.ENUM_CENTERED_2ND .OR.  
      &    advectionScheme.EQ.ENUM_UPWIND_3RD .OR.  
      &    advectionScheme.EQ.ENUM_CENTERED_4TH ) THEN  
310  C     Note: wVel needs to be masked  C     Note: wVel needs to be masked
311        IF (K.GE.2) THEN        IF (K.GE.2) THEN
312  C-    Compute vertical advective flux in the interior:  C-    Compute vertical advective flux in the interior:
# Line 358  C           boundary condition. Line 364  C           boundary condition.
364        ELSE        ELSE
365         CALL GAD_DIFF_R(bi,bj,k,KappaRT,tracer,df,myThid)         CALL GAD_DIFF_R(bi,bj,k,KappaRT,tracer,df,myThid)
366        ENDIF        ENDIF
 c     DO j=1-Oly,sNy+Oly  
 c      DO i=1-Olx,sNx+Olx  
 c       fVerT(i,j,kUp) = fVerT(i,j,kUp) + df(i,j)*maskUp(i,j)  
 c      ENDDO  
 c     ENDDO  
367    
368  #ifdef ALLOW_GMREDI  #ifdef ALLOW_GMREDI
369  C-    GM/Redi flux in R  C-    GM/Redi flux in R
# Line 370  C-    GM/Redi flux in R Line 371  C-    GM/Redi flux in R
371  C *note* should update GMREDI_RTRANSPORT to set df  *aja*  C *note* should update GMREDI_RTRANSPORT to set df  *aja*
372         CALL GMREDI_RTRANSPORT(         CALL GMREDI_RTRANSPORT(
373       I     iMin,iMax,jMin,jMax,bi,bj,K,       I     iMin,iMax,jMin,jMax,bi,bj,K,
374       I     Tracer,       I     Tracer,tracerIdentity,
375       U     df,       U     df,
376       I     myThid)       I     myThid)
 c      DO j=1-Oly,sNy+Oly  
 c       DO i=1-Olx,sNx+Olx  
 c        fVerT(i,j,kUp) = fVerT(i,j,kUp) + df(i,j)*maskUp(i,j)  
 c       ENDDO  
 c      ENDDO  
377        ENDIF        ENDIF
378  #endif  #endif
379    
# Line 431  C--   Divergence of fluxes Line 427  C--   Divergence of fluxes
427         ENDDO         ENDDO
428        ENDDO        ENDDO
429    
430    #ifdef NONLIN_FRSURF
431    C-- account for 3.D divergence of the flow in rStar coordinate:
432          IF (calcAdvection .AND. select_rStar.GT.0) THEN
433           DO j=1-Oly,sNy+Oly-1
434            DO i=1-Olx,sNx+Olx-1
435             gTracer(i,j,k,bi,bj) = gTracer(i,j,k,bi,bj)
436         &     - (rStarExpC(i,j,bi,bj) - 1. _d 0)/deltaTfreesurf
437         &       *tracer(i,j,k,bi,bj)*maskC(i,j,k,bi,bj)
438            ENDDO
439           ENDDO
440          ENDIF
441          IF (calcAdvection .AND. select_rStar.LT.0) THEN
442           DO j=1-Oly,sNy+Oly-1
443            DO i=1-Olx,sNx+Olx-1
444             gTracer(i,j,k,bi,bj) = gTracer(i,j,k,bi,bj)
445         &     - rStarDhCDt(i,j,bi,bj)
446         &       *tracer(i,j,k,bi,bj)*maskC(i,j,k,bi,bj)
447            ENDDO
448           ENDDO
449          ENDIF
450    #endif /* NONLIN_FRSURF */
451          
452    
453        RETURN        RETURN
454        END        END
Legend:
Removed from v.1.12.8.1  
changed lines
  Added in v.1.17
  ViewVC Help
Powered by ViewVC 1.1.22