/[MITgcm]/MITgcm/pkg/generic_advdiff/gad_som_adv_y.F

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

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

-revision 1.1 by jmc,
Tue Jan 16 04:38:34 2007 UTC
+revision 1.5 by jmc,
Fri May  9 21:43:16 2008 UTC
 Line 9 
 C !ROUTINE: GAD_SOM_ADV_Y
  C !INTERFACE: ==========================================================
        SUBROUTINE GAD_SOM_ADV_Y(
       I           bi,bj,k, limiter,
+      I           overlapOnly, interiorOnly,
+      I           N_edge, S_edge, E_edge, W_edge,
       I           deltaTloc, vTrans,
       U           sm_v, sm_o, sm_x, sm_y, sm_z,
       U           sm_xx, sm_yy, sm_zz, sm_xy, sm_xz, sm_yz,
-Line 34 
 C---+----1----+----2----+----3----+----4
+Line 36 
 C---+----1----+----2----+----3----+----4
  C !USES: ===============================================================
        IMPLICIT NONE
  #include "SIZE.h"
- c #include "GRID.h"
  #include "GAD.h"
  C !INPUT PARAMETERS: ===================================================
- C  bi,bj        :: tile indices
+ C  bi,bj         :: tile indices
- C  k            :: vertical level
+ C  k             :: vertical level
- C  limiter      :: 0: no limiter ; 1: Prather, 1986 limiter
+ C  limiter       :: 0: no limiter ; 1: Prather, 1986 limiter
- C  vTrans       :: zonal volume transport
+ C  overlapOnly   :: only update the edges of myTile, but not the interior
- C  myThid       :: my Thread Id. number
+ C  interiorOnly  :: only update the interior of myTile, but not the edges
+ C [N,S,E,W]_edge :: true if N,S,E,W edge of myTile is an Edge of the cube
+ C  vTrans        :: zonal volume transport
+ C  myThid        :: my Thread Id. number
        INTEGER bi,bj,k
        INTEGER limiter
+       LOGICAL overlapOnly, interiorOnly
+       LOGICAL N_edge, S_edge, E_edge, W_edge
        _RL deltaTloc
        _RL vTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
- c     _RL tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
        INTEGER myThid
  C !OUTPUT PARAMETERS: ==================================================
-Line 70 
 C  vT           :: meridional advective
+Line 75 
 C  vT           :: meridional advective
        _RL sm_yz (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
        _RL vT    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
- #ifdef GAD_ALLOW_SOM_ADVECT
  C !LOCAL VARIABLES: ====================================================
- C  i,j          :: loop indices
+ C  i,j           :: loop indices
- C  vLoc         :: volume transported (per time step)
+ C  vLoc          :: volume transported (per time step)
+ C [iMin,iMax]Upd :: loop range to update tracer field
+ C [jMin,jMax]Upd :: loop range to update tracer field
+ C  nbStrips      :: number of strips (if region to update is splitted)
        _RL two, three
        PARAMETER( two   = 2. _d 0 )
        PARAMETER( three = 3. _d 0 )
        INTEGER i,j
+       INTEGER ns, nbStrips
+       INTEGER iMinUpd(2), iMaxUpd(2), jMinUpd(2), jMaxUpd(2)
+       _RL  recip_dT
        _RL  slpmax, s1max, s1new, s2new
        _RL  vLoc, alf1, alf1q, alpmn
        _RL  alfp, alpq, alp1, locTp
-Line 108 
 C  vLoc         :: volume transported (p
+Line 118 
 C  vLoc         :: volume transported (p
        _RL  fn_yz(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
  CEOP
+       recip_dT = 0.
+       IF ( deltaTloc.GT.0. _d 0 ) recip_dT = 1.0 _d 0 / deltaTloc
+ C-    Set loop ranges for updating tracer field (splitted in 2 strips)
+       nbStrips   = 1
+       iMinUpd(1) = 1-Olx
+       iMaxUpd(1) = sNx+Olx
+       jMinUpd(1) = 1-Oly+1
+       jMaxUpd(1) = sNy+Oly-1
+       IF ( overlapOnly ) THEN
+ C     update in overlap-Only
+         IF ( S_edge ) jMinUpd(1) = 1
+         IF ( N_edge ) jMaxUpd(1) = sNy
+         IF ( W_edge ) THEN
+           iMinUpd(1) = 1-Olx
+           iMaxUpd(1) = 0
+         ENDIF
+         IF ( E_edge ) THEN
+           IF ( W_edge ) nbStrips = 2
+           iMinUpd(nbStrips) = sNx+1
+           iMaxUpd(nbStrips) = sNx+Olx
+         ENDIF
+       ELSE
+ C     do not only update the overlap
+         IF ( interiorOnly .AND. W_edge ) iMinUpd(1) = 1
+         IF ( interiorOnly .AND. E_edge ) iMaxUpd(1) = sNx
+       ENDIF
+ C-    Internal exchange for calculations in Y
+ c     IF ( overlapOnly ) THEN
+ c         CALL GAD_SOM_FILL_CS_CORNER( .FALSE.,
+ c    U           sm_v,  sm_o,  sm_x,  sm_y,  sm_z,
+ c    U           sm_xx, sm_yy, sm_zz, sm_xy, sm_xz, sm_yz,
+ c    I           bi, bj, myThid )
+ c     ENDIF
+ C--   start 1rst loop on strip number "ns"
+       DO ns=1,nbStrips
        IF ( limiter.EQ.1 ) THEN
-        DO j=1-OLy,sNy+OLy
+        DO j=jMinUpd(1)-1,jMaxUpd(1)+1
-         DO i=1-OLx,sNx+OLx
+         DO i=iMinUpd(ns),iMaxUpd(ns)
  C     If flux-limiting transport is to be applied, place limits on
  C     appropriate moments before transport.
           slpmax = 0.
-Line 121 
 C     appropriate moments before transpo
+Line 170 
 C     appropriate moments before transpo
       &                MAX(ABS(s1new)-slpmax,sm_yy(i,j))  )
           sm_xy(i,j) = MIN( slpmax, MAX(-slpmax,sm_xy(i,j)) )
           sm_yz(i,j) = MIN( slpmax, MAX(-slpmax,sm_yz(i,j)) )
-          sm_y (i,j) = s1new ;
+          sm_y (i,j) = s1new
-          sm_yy(i,j) = s2new ;
+          sm_yy(i,j) = s2new
          ENDDO
         ENDDO
        ENDIF
  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
  C---  part.1 : calculate flux for all moments
-       DO i=1-OLx,sNx+OLx
+       DO j=jMinUpd(1),jMaxUpd(1)+1
-        vT(i,1-OLy)=0.
+        DO i=iMinUpd(ns),iMaxUpd(ns)
-       ENDDO
-       DO j=1-OLy+1,sNy+OLy
-        DO i=1-OLx,sNx+OLx
          vLoc = vTrans(i,j)*deltaTloc
  C--    Flux from (j-1) to (j) when V>0 (i.e., take right side of box j-1)
          fp_v (i,j) = MAX( 0. _d 0,  vLoc )
-Line 175 
 C       use same indexing as velocity, "
+Line 221 
 C       use same indexing as velocity, "
          fn_zz(i,j) = aln(i,j)*sm_zz(i, j )
          fn_xz(i,j) = aln(i,j)*sm_xz(i, j )
  C--    Save zero-order flux:
-         vT(i,j)    = fp_o(i,j) - fn_o(i,j)
+         vT(i,j) = ( fp_o(i,j) - fn_o(i,j) )*recip_dT
         ENDDO
        ENDDO
+ C--   end 1rst loop on strip number "ns"
+ c     ENDDO
+ C-    Internal exchange for next calculations in X
+ c     IF ( overlapOnly ) THEN
+ c         CALL GAD_SOM_FILL_CS_CORNER( .TRUE.,
+ c    U           sm_v,  sm_o,  sm_x,  sm_y,  sm_z,
+ c    U           sm_xx, sm_yy, sm_zz, sm_xy, sm_xz, sm_yz,
+ c    I           bi, bj, myThid )
+ c     ENDIF
  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
+ C--   start 2nd loop on strip number "ns"
+ c     DO ns=1,nbStrips
  C---  part.2 : re-adjust moments remaining in the box
  C      take off from grid box (j): negative V(j) and positive V(j+1)
-       DO j=1-OLy+1,sNy+OLy-1
+       DO j=jMinUpd(1),jMaxUpd(1)
-        DO i=1-OLx,sNx+OLx
+        DO i=iMinUpd(ns),iMaxUpd(ns)
          alf1  = 1. _d 0 - aln(i,j) - alp(i,j+1)
          alf1q = alf1*alf1
          alpmn = alp(i,j+1) - aln(i,j)
-Line 204 
 C      take off from grid box (j): negat
+Line 264 
 C      take off from grid box (j): negat
  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
  C---  part.3 : Put the temporary moments into appropriate neighboring boxes
  C      add into grid box (j): positive V(j) and negative V(j+1)
-       DO j=1-OLy+1,sNy+OLy-1
+       DO j=jMinUpd(1),jMaxUpd(1)
-        DO i=1-OLx,sNx+OLx
+        DO i=iMinUpd(ns),iMaxUpd(ns)
          sm_v (i,j) = sm_v (i,j) + fp_v (i,j) + fn_v (i,j+1)
          alfp = fp_v(i, j )/sm_v(i,j)
          alfn = fn_v(i,j+1)/sm_v(i,j)
-Line 243 
 C      add into grid box (j): positive V
+Line 303 
 C      add into grid box (j): positive V
         ENDDO
        ENDDO
- #endif /* GAD_ALLOW_SOM_ADVECT */
+ C--   end 2nd loop on strip number "ns"
+       ENDDO
        RETURN
        END

 Legend:



Removed from v.1.1
 


changed lines


 
Added in v.1.5
 Legend:



Removed from v.1.1
 


changed lines


 
Added in v.1.5
-Removed from v.1.1
+Added in v.1.5

	ViewVC Help
Powered by ViewVC 1.1.22