/[MITgcm]/MITgcm/pkg/ptracers/ptracers_integrate.F
ViewVC logotype

Diff of /MITgcm/pkg/ptracers/ptracers_integrate.F

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

revision 1.51 by jmc, Fri Dec 6 01:56:52 2013 UTC revision 1.54 by jmc, Mon Jul 14 22:47:09 2014 UTC
# Line 8  C !ROUTINE: PTRACERS_INTEGRATE Line 8  C !ROUTINE: PTRACERS_INTEGRATE
8    
9  C !INTERFACE: ==========================================================  C !INTERFACE: ==========================================================
10        SUBROUTINE PTRACERS_INTEGRATE(        SUBROUTINE PTRACERS_INTEGRATE(
11       I                    bi, bj,       I                    bi, bj, recip_hFac,
12       I                    uFld, vFld, wFld,       I                    uFld, vFld, wFld,
13       U                    KappaRk,       U                    KappaRk,
14       I                    myTime, myIter, myThid )       I                    myTime, myIter, myThid )
# Line 23  C !USES: =============================== Line 23  C !USES: ===============================
23  #include "SIZE.h"  #include "SIZE.h"
24  #include "EEPARAMS.h"  #include "EEPARAMS.h"
25  #include "PARAMS.h"  #include "PARAMS.h"
26    #include "GRID.h"
27  #ifdef ALLOW_LONGSTEP  #ifdef ALLOW_LONGSTEP
28  #include "LONGSTEP_PARAMS.h"  #include "LONGSTEP_PARAMS.h"
29  #endif  #endif
# Line 37  C !USES: =============================== Line 38  C !USES: ===============================
38  #endif  #endif
39    
40  C !INPUT PARAMETERS: ===================================================  C !INPUT PARAMETERS: ===================================================
41  C  bi, bj               :: tile indices  C  bi, bj           :: tile indices
42  C  uFld, vFld, wFld     :: Local copy of velocity field (3 components)  C  recip_hFac       :: reciprocal of cell open-depth factor (@ next iter)
43  C  KappaRk              :: vertical diffusion used for one passive tracer  C  uFld, vFld, wFld :: Local copy of velocity field (3 components)
44  C  myTime               :: model time  C  KappaRk          :: vertical diffusion used for one passive tracer
45  C  myIter               :: time-step number  C  myTime           :: model time
46  C  myThid               :: thread number  C  myIter           :: time-step number
47    C  myThid           :: thread number
48        INTEGER bi, bj        INTEGER bi, bj
49        _RL uFld    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RS recip_hFac(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
50        _RL vFld    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL uFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
51        _RL wFld    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL vFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
52        _RL KappaRk (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL wFld      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
53          _RL KappaRk   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
54        _RL myTime        _RL myTime
55        INTEGER myIter        INTEGER myIter
56        INTEGER myThid        INTEGER myThid
# Line 58  C  none Line 61  C  none
61  #ifdef ALLOW_PTRACERS  #ifdef ALLOW_PTRACERS
62    
63  C !LOCAL VARIABLES: ====================================================  C !LOCAL VARIABLES: ====================================================
64  C  iTracer              :: tracer index  C  iTracer          :: tracer index
65  C  iMin,iMax,jMin,jMax  :: loop ranges  C  iMin, iMax       :: 1rst index loop range
66  C  k                    :: vertical level number  C  jMin, jMax       :: 2nd  index loop range
67  C  kUp,kDown            :: toggle indices for even/odd level fluxes  C  k                :: vertical level number
68  C  kM1                  :: =min(1,k-1)  C  kUp,kDown        :: toggle indices for even/odd level fluxes
69  C  GAD_TR               :: passive tracer id (GAD_TR1+iTracer-1)  C  kM1              :: =min(1,k-1)
70  C  xA                   :: face area at U points in level k  C  GAD_TR           :: passive tracer id (GAD_TR1+iTracer-1)
71  C  yA                   :: face area at V points in level k  C  xA               :: face area at U points in level k
72  C  maskUp               :: mask for vertical transport  C  yA               :: face area at V points in level k
73  C  uTrans               :: zonal transport in level k  C  maskUp           :: mask for vertical transport
74  C  vTrans               :: meridional transport in level k  C  uTrans           :: zonal transport in level k
75  C  rTrans               :: vertical volume transport at interface k  C  vTrans           :: meridional transport in level k
76  C  rTransKp             :: vertical volume transport at interface k+1  C  rTrans           :: vertical volume transport at interface k
77  C  fVerT                :: passive tracer vertical flux  C  rTransKp         :: vertical volume transport at interface k+1
78    C  fZon             :: passive tracer zonal flux
79    C  fMer             :: passive tracer meridional flux
80    C  fVer             :: passive tracer vertical flux
81    C  gTr_AB           :: Adams-Bashforth tracer tendency increment
82        INTEGER iTracer        INTEGER iTracer
83        INTEGER iMin,iMax,jMin,jMax        INTEGER iMin,iMax,jMin,jMax
84        INTEGER i, j, k        INTEGER i, j, k
# Line 84  C  fVerT                :: passive trace Line 91  C  fVerT                :: passive trace
91        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL vTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
92        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rTrans  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
93        _RL rTransKp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL rTransKp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
94        _RL fVerT   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)        _RL fZon    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
95          _RL fMer    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
96          _RL fVer    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,2)
97        _RL gTr_AB  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL gTr_AB  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
98        LOGICAL calcAdvection        LOGICAL calcAdvection
99        INTEGER iterNb        INTEGER iterNb
# Line 100  CEOP Line 109  CEOP
109    
110  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
111    
 C-    Loop ranges for daughter routines  
       iMin = 1-OLx+2  
       iMax = sNx+OLx-1  
       jMin = 1-OLy+2  
       jMax = sNy+OLy-1  
   
112  C-    Compute iter at beginning of ptracer time step  C-    Compute iter at beginning of ptracer time step
113  #ifdef ALLOW_LONGSTEP  #ifdef ALLOW_LONGSTEP
114        iterNb = myIter - LS_nIter + 1        iterNb = myIter - LS_nIter + 1
# Line 129  C-    Initialise tracer tendency to zero Line 132  C-    Initialise tracer tendency to zero
132         IF ( PTRACERS_StepFwd(iTracer) ) THEN         IF ( PTRACERS_StepFwd(iTracer) ) THEN
133          GAD_TR = GAD_TR1 + iTracer - 1          GAD_TR = GAD_TR1 + iTracer - 1
134    
135    C-    Loop ranges for daughter routines
136            iMin = 1-OLx+2
137            iMax = sNx+OLx-1
138            jMin = 1-OLy+2
139            jMax = sNy+OLy-1
140    
141  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
142            act0 = iTracer - 1            act0 = iTracer - 1
143            max0 = PTRACERS_num            max0 = PTRACERS_num
# Line 148  C-    Initialise tracer tendency to zero Line 157  C-    Initialise tracer tendency to zero
157    
158          DO j=1-OLy,sNy+OLy          DO j=1-OLy,sNy+OLy
159           DO i=1-OLx,sNx+OLx           DO i=1-OLx,sNx+OLx
160             fVerT(i,j,1) = 0. _d 0             fVer(i,j,1) = 0. _d 0
161             fVerT(i,j,2) = 0. _d 0             fVer(i,j,2) = 0. _d 0
162           ENDDO           ENDDO
163          ENDDO          ENDDO
164  #ifdef ALLOW_AUTODIFF  #ifdef ALLOW_AUTODIFF
# Line 162  C-    Initialise tracer tendency to zero Line 171  C-    Initialise tracer tendency to zero
171          ENDDO          ENDDO
172  #endif /* ALLOW_AUTODIFF */  #endif /* ALLOW_AUTODIFF */
173    
174          IF ( .NOT.implicitDiffusion ) THEN          CALL CALC_3D_DIFFUSIVITY(
           CALL CALC_3D_DIFFUSIVITY(  
175       I         bi, bj, iMin,iMax,jMin,jMax,       I         bi, bj, iMin,iMax,jMin,jMax,
176       I         GAD_TR,       I         GAD_TR,
177       I         PTRACERS_useGMRedi(iTracer), PTRACERS_useKPP(iTracer),       I         PTRACERS_useGMRedi(iTracer), PTRACERS_useKPP(iTracer),
178       O         kappaRk,       O         kappaRk,
179       I         myThid)       I         myThid)
         ENDIF  
180    
181  #ifndef DISABLE_MULTIDIM_ADVECTION  #ifndef DISABLE_MULTIDIM_ADVECTION
182  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 186  CADJ &      = comlev1_bibj_ptracers, key Line 193  CADJ &      = comlev1_bibj_ptracers, key
193       I                       PTRACERS_ImplVertAdv(iTracer),       I                       PTRACERS_ImplVertAdv(iTracer),
194       I                       PTRACERS_advScheme(iTracer),       I                       PTRACERS_advScheme(iTracer),
195       I                       PTRACERS_advScheme(iTracer),       I                       PTRACERS_advScheme(iTracer),
196       I                       GAD_TR1+iTracer-1,       I                       GAD_TR,
197       I                       PTRACERS_dTLev,       I                       PTRACERS_dTLev, uFld, vFld, wFld,
      I                       uFld, vFld, wFld,  
198       I                       pTracer(1-OLx,1-OLy,1,1,1,iTracer),       I                       pTracer(1-OLx,1-OLy,1,1,1,iTracer),
199       U                       _Ptracers_som(:,:,:,:,:,:,iTracer),       U                       _Ptracers_som(:,:,:,:,:,:,iTracer),
200       O                       gPtr(1-OLx,1-OLy,1,1,1,iTracer),       O                       gPtr(1-OLx,1-OLy,1,1,1,iTracer),
# Line 204  CADJ &      = comlev1_bibj_ptracers, key Line 210  CADJ &      = comlev1_bibj_ptracers, key
210       I                        PTRACERS_ImplVertAdv(iTracer),       I                        PTRACERS_ImplVertAdv(iTracer),
211       I                        PTRACERS_advScheme(iTracer),       I                        PTRACERS_advScheme(iTracer),
212       I                        PTRACERS_advScheme(iTracer),       I                        PTRACERS_advScheme(iTracer),
213       I                        GAD_TR1+iTracer-1,       I                        GAD_TR,
214       I                        PTRACERS_dTLev,       I                        PTRACERS_dTLev, uFld, vFld, wFld,
      I                        uFld, vFld, wFld,  
215       I                        pTracer(1-OLx,1-OLy,1,1,1,iTracer),       I                        pTracer(1-OLx,1-OLy,1,1,1,iTracer),
216       O                        gPtr(1-OLx,1-OLy,1,1,1,iTracer),       O                        gPtr(1-OLx,1-OLy,1,1,1,iTracer),
217       I                        bi, bj, myTime, myIter, myThid )       I                        bi, bj, myTime, myIter, myThid )
# Line 225  C-    Start vertical index (k) loop (Nr: Line 230  C-    Start vertical index (k) loop (Nr:
230            kDown= 1+MOD(k,2)            kDown= 1+MOD(k,2)
231    
232  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
233  CADJ STORE fVerT(:,:,:) = comlev1_bibj_k_ptracers,  CADJ STORE fVer(:,:,:) = comlev1_bibj_k_ptracers,
234  CADJ &     key = kkey, byte = isbyte, kind = isbyte  CADJ &     key = kkey, byte = isbyte, kind = isbyte
235  CADJ STORE gPtr(:,:,k,bi,bj,iTracer) = comlev1_bibj_k_ptracers,  CADJ STORE gPtr(:,:,k,bi,bj,iTracer) = comlev1_bibj_k_ptracers,
236  CADJ &     key = kkey, byte = isbyte, kind = isbyte  CADJ &     key = kkey, byte = isbyte, kind = isbyte
# Line 258  C      (advection, [explicit] diffusion, Line 263  C      (advection, [explicit] diffusion,
263       I             .FALSE., .FALSE.,       I             .FALSE., .FALSE.,
264       I             PTRACERS_useGMRedi(iTracer),       I             PTRACERS_useGMRedi(iTracer),
265       I             PTRACERS_useKPP(iTracer),       I             PTRACERS_useKPP(iTracer),
266       U             fVerT, gPtr(1-OLx,1-OLy,1,1,1,iTracer),       O             fZon, fMer,
267         U             fVer, gPtr(1-OLx,1-OLy,1,1,1,iTracer),
268       I             myTime, myIter, myThid )       I             myTime, myIter, myThid )
269    
270  C-    External forcing term(s)  C-    External forcing term(s)
271            IF ( tracForcingOutAB.NE.1 )            IF ( tracForcingOutAB.NE.1 )
272       &      CALL PTRACERS_FORCING(       &      CALL PTRACERS_APPLY_FORCING(
273       I                    bi,bj,iMin,iMax,jMin,jMax,k,iTracer,       U                    gPtr(1-OLx,1-OLy,k,bi,bj,iTracer),
274       U                    gPtr(1-OLx,1-OLy,1,1,1,iTracer),       I                    surfaceForcingPTr(1-OLx,1-OLy,bi,bj,iTracer),
275       I                    surfaceForcingPTr(1-OLx,1-OLy,1,1,iTracer),       I                    iMin,iMax,jMin,jMax, k, bi, bj,
276       I                    myIter, myTime, myThid )       I                    iTracer, myTime, myIter, myThid )
277    
278  C-    If using Adams-Bashforth II, then extrapolate tendencies  C-    If using Adams-Bashforth II, then extrapolate tendencies
279  C     gPtr is now the tracer tendency for explicit advection/diffusion  C     gPtr is now the tracer tendency for explicit advection/diffusion
# Line 294  C     prevent gPtr from being replaced b Line 300  C     prevent gPtr from being replaced b
300    
301  C-    External forcing term(s)  C-    External forcing term(s)
302            IF ( tracForcingOutAB.EQ.1 )            IF ( tracForcingOutAB.EQ.1 )
303       &      CALL PTRACERS_FORCING(       &      CALL PTRACERS_APPLY_FORCING(
304       I                    bi,bj,iMin,iMax,jMin,jMax,k,iTracer,       U                    gPtr(1-OLx,1-OLy,k,bi,bj,iTracer),
305       U                    gPtr(1-OLx,1-OLy,1,1,1,iTracer),       I                    surfaceForcingPTr(1-OLx,1-OLy,bi,bj,iTracer),
306       I                    surfaceForcingPTr(1-OLx,1-OLy,1,1,iTracer),       I                    iMin,iMax,jMin,jMax, k, bi, bj,
307       I                    myIter, myTime, myThid )       I                    iTracer, myTime, myIter, myThid )
308    
309  #ifdef NONLIN_FRSURF  #ifdef NONLIN_FRSURF
310  C-    Account for change in level thickness  C-    Account for change in level thickness
# Line 325  C-    Integrate forward in time, storing Line 331  C-    Integrate forward in time, storing
331  C-    end of vertical index (k) loop (Nr:1)  C-    end of vertical index (k) loop (Nr:1)
332          ENDDO          ENDDO
333    
334    #ifdef ALLOW_DOWN_SLOPE
335            IF ( PTRACERS_useDWNSLP(iTracer) ) THEN
336              IF ( usingPCoords ) THEN
337                CALL DWNSLP_APPLY(
338         I                  GAD_TR, bi, bj, kSurfC,
339         I                  recip_drF, recip_hFac, recip_rA,
340         I                  PTRACERS_dTLev,
341         I                  pTracer(1-OLx,1-OLy,1,1,1,iTracer),
342         O                  gPtr(1-OLx,1-OLy,1,1,1,iTracer),
343         I                  myTime, myIter, myThid )
344              ELSE
345                CALL DWNSLP_APPLY(
346         I                  GAD_TR, bi, bj, kLowC,
347         I                  recip_drF, recip_hFac, recip_rA,
348         I                  PTRACERS_dTLev,
349         I                  pTracer(1-OLx,1-OLy,1,1,1,iTracer),
350         O                  gPtr(1-OLx,1-OLy,1,1,1,iTracer),
351         I                  myTime, myIter, myThid )
352              ENDIF
353            ENDIF
354    #endif /* ALLOW_DOWN_SLOPE */
355    
356    C       All explicit advection/diffusion/sources should now be
357    C       done. The updated tracer field is in gPtr.
358    #ifdef ALLOW_MATRIX
359    C       Accumalate explicit tendency and also reset gPtr to initial
360    C       tracer field for implicit matrix calculation
361            IF ( useMATRIX ) CALL MATRIX_STORE_TENDENCY_EXP(
362         I                               iTracer, bi, bj,
363         I                               myTime, myIter, myThid )
364    #endif /* ALLOW_MATRIX */
365    
366    C--     Vertical advection & diffusion (implicit) for passive tracers
367            iMin = 0
368            iMax = sNx+1
369            jMin = 0
370            jMax = sNy+1
371    
372    #ifdef ALLOW_AUTODIFF_TAMC
373    CADJ STORE gPtr(:,:,:,bi,bj,iTracer) = comlev1_bibj_ptracers,
374    CADJ &     key=iptrkey, byte=isbyte
375    #endif /* ALLOW_AUTODIFF_TAMC */
376    
377    #ifdef INCLUDE_IMPLVERTADV_CODE
378            IF ( PTRACERS_ImplVertAdv(iTracer) ) THEN
379              CALL GAD_IMPLICIT_R(
380         I         PTRACERS_ImplVertAdv(iTracer),
381         I         PTRACERS_advScheme(iTracer), GAD_TR,
382         I         PTRACERS_dTLev, kappaRk, recip_hFac, wFld,
383         I         pTracer(1-OLx,1-OLy,1,1,1,iTracer),
384         U         gPtr(1-OLx,1-OLy,1,1,1,iTracer),
385         I         bi, bj, myTime, myIter, myThid )
386    
387            ELSEIF ( implicitDiffusion ) THEN
388    #else /* INCLUDE_IMPLVERTADV_CODE */
389            IF ( implicitDiffusion ) THEN
390    #endif /* INCLUDE_IMPLVERTADV_CODE */
391              CALL IMPLDIFF(
392         I         bi, bj, iMin, iMax, jMin, jMax,
393         I         GAD_TR, kappaRk, recip_hFac,
394         U         gPtr(1-OLx,1-OLy,1,1,1,iTracer),
395         I         myThid )
396            ENDIF
397    
398    #ifdef ALLOW_OBCS
399    C--     Apply open boundary conditions for each passive tracer
400            IF ( useOBCS ) THEN
401              CALL OBCS_APPLY_PTRACER(
402         I         bi, bj, 0, iTracer,
403         U         gPtr(1-OLx,1-OLy,1,bi,bj,iTracer),
404         I         myThid )
405            ENDIF
406    #endif /* ALLOW_OBCS */
407    
408  C--   end of tracer loop  C--   end of tracer loop
409         ENDIF         ENDIF
410        ENDDO        ENDDO
Legend:
Removed from v.1.51  
changed lines
  Added in v.1.54
  ViewVC Help
Powered by ViewVC 1.1.22