/[MITgcm]/MITgcm/model/src/diags_phi_rlow.F
ViewVC logotype

Diff of /MITgcm/model/src/diags_phi_rlow.F

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

revision 1.1 by jmc, Sun Feb 9 02:58:39 2003 UTC revision 1.6 by jmc, Fri Jan 15 22:39:54 2010 UTC
# Line 6  C $Name$ Line 6  C $Name$
6  CBOP  CBOP
7  C     !ROUTINE: DIAGS_PHI_RLOW  C     !ROUTINE: DIAGS_PHI_RLOW
8  C     !INTERFACE:  C     !INTERFACE:
9        SUBROUTINE DIAGS_PHI_RLOW(        SUBROUTINE DIAGS_PHI_RLOW(
10       I                       k, bi, bj, iMin,iMax, jMin,jMax,       I                       k, bi, bj, iMin,iMax, jMin,jMax,
11       I                       phiHyd, alphRho, tFld, sFld,       I                       phiHydF, phiHydC, alphRho, tFld, sFld,
12       I                       myTime, myIter, myThid)       I                       myTime, myIter, myThid)
13  C     !DESCRIPTION: \bv  C     !DESCRIPTION: \bv
14  C     *==========================================================*  C     *==========================================================*
15  C     | S/R DIAGS_PHI_RLOW                                      C     | S/R DIAGS_PHI_RLOW
16  C     | o Diagnose Phi-Hydrostatic at r-lower boundary  C     | o Diagnose Phi-Hydrostatic at r-lower boundary
17  C     |   = bottom pressure (ocean in z-coord) ;  C     |   = bottom pressure (ocean in z-coord) ;
18  C     |   = sea surface elevation (ocean in p-coord) ;  C     |   = sea surface elevation (ocean in p-coord) ;
19  C     |   = height at the top of atmosphere (in p-coord) ;  C     |   = height at the top of atmosphere (in p-coord) ;
20  C     *==========================================================*  C     *==========================================================*
# Line 32  C     == Global variables == Line 32  C     == Global variables ==
32    
33  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
34  C     == Routine Arguments ==  C     == Routine Arguments ==
35  C     bi,bj      :: tile index              C     bi,bj      :: tile index
36  C     iMin,iMax,jMin,jMax :: Loop counters  C     iMin,iMax,jMin,jMax :: Loop counters
37  C     phiHyd     :: Hydrostatic Potential anomaly  C     phiHydF    :: hydrostatic potential anomaly at middle between
38    C                   2 centers k & k+1 (interface k+1)
39    C     phiHydC    :: hydrostatic potential anomaly at cell center
40  C                  (atmos: =Geopotential ; ocean-z: =Pressure/rho)  C                  (atmos: =Geopotential ; ocean-z: =Pressure/rho)
41  C     alphRho    :: Density (z-coord) or specific volume (p-coord)  C     alphRho    :: Density (z-coord) or specific volume (p-coord)
42  C     tFld       :: Potential temp.  C     tFld       :: Potential temp.
43  C     sFld       :: Salinity  C     sFld       :: Salinity
44  C     myTime :: Current time  C     myTime     :: Current time
45  C     myIter :: Current iteration number  C     myIter     :: Current iteration number
46  C     myThid :: Instance number for this call of the routine.  C     myThid     :: my Thread Id number
47        INTEGER k, bi,bj, iMin,iMax, jMin,jMax        INTEGER k, bi,bj, iMin,iMax, jMin,jMax
48        _RL phiHyd(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL phiHydF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
49          _RL phiHydC(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
50        _RL alphRho(1-OLx:sNx+OLx,1-OLy:sNy+OLy)        _RL alphRho(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
51        _RL tFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)        _RL tFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
52        _RL sFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)        _RL sFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
# Line 57  C     == Local variables == Line 60  C     == Local variables ==
60  C     i,j :: Loop counters  C     i,j :: Loop counters
61        INTEGER i,j        INTEGER i,j
62        _RL zero, one, half        _RL zero, one, half
63        _RL dRloc        _RL ddRloc, ratioRm, ratioRp
64        PARAMETER ( zero= 0. _d 0 , one= 1. _d 0 , half= .5 _d 0 )        PARAMETER ( zero= 0. _d 0 , one= 1. _d 0 , half= .5 _d 0 )
65  CEOP  CEOP
66    
67          dRloc=drC(k)        IF ( usingZCoords ) THEN
         IF (k.EQ.1) dRloc=drF(1)  
68    
69        IF ( buoyancyRelation .eq. 'OCEANIC' ) THEN  C----- Compute bottom pressure deviation from gravity*rho0*H
70    C      Start from phiHyd at the (bottom) tracer point and add Del_h*g*rho'
71    C      with Del_h = distance from the bottom up to tracer point
72    
73    C--    Initialise to zero (otherwise phi0surf accumulate over land)
74           IF ( k.EQ.1 ) THEN
75             DO j=1-Oly,sNy+Oly
76              DO i=1-Olx,sNx+Olx
77                phiHydLow(i,j,bi,bj) = 0. _d 0
78              ENDDO
79             ENDDO
80           ENDIF
81    
82         IF (integr_GeoPot.EQ.1) THEN         IF (integr_GeoPot.EQ.1) THEN
83  C  --  Finite Volume Form  C  --  Finite Volume Form
# Line 72  C  --  Finite Volume Form Line 85  C  --  Finite Volume Form
85           DO j=jMin,jMax           DO j=jMin,jMax
86            DO i=iMin,iMax            DO i=iMin,iMax
87             IF ( k .EQ. kLowC(i,j,bi,bj) ) THEN             IF ( k .EQ. kLowC(i,j,bi,bj) ) THEN
88              phiHydLow(i,j,bi,bj) = phiHyd(i,j,k)               ddRloc = rC(k)-R_low(i,j,bi,bj)
89       &          + hFacC(i,j,k,bi,bj)               phiHydLow(i,j,bi,bj) = phiHydC(i,j)
90       &             *drF(K)*gravity*alphRho(i,j)*recip_rhoConst       &            + ddRloc*gravity*alphRho(i,j)*recip_rhoConst
      &          + Bo_surf(i,j,bi,bj)*etaN(i,j,bi,bj)  
      &          + phi0surf(i,j,bi,bj)  
91             ENDIF             ENDIF
92            ENDDO            ENDDO
93           ENDDO           ENDDO
# Line 84  C  --  Finite Volume Form Line 95  C  --  Finite Volume Form
95         ELSE         ELSE
96  C  --  Finite Difference Form  C  --  Finite Difference Form
97    
98  C---------- Compute bottom pressure deviation from gravity*rho0*H           ratioRm = one
99  C           This has to be done starting from phiHyd at the current           ratioRp = one
100  C           tracer point and .5 of the cell's thickness has to be           IF (k.GT.1 ) ratioRm = half*drC(k)/(rF(k)-rC(k))
101  C           substracted from hFacC           IF (k.LT.Nr) ratioRp = half*drC(k+1)/(rC(k)-rF(k+1))
102    
103           DO j=jMin,jMax           DO j=jMin,jMax
104            DO i=iMin,iMax            DO i=iMin,iMax
105             IF ( K .EQ. kLowC(i,j,bi,bj) ) THEN             IF ( k .EQ. kLowC(i,j,bi,bj) ) THEN
106              phiHydLow(i,j,bi,bj) = phiHyd(i,j,k)               ddRloc = rC(k)-R_low(i,j,bi,bj)
107       &          + (half*dRloc+(hFacC(i,j,k,bi,bj)-half)*drF(k))               phiHydLow(i,j,bi,bj) = phiHydC(i,j)
108       &                    *gravity*alphRho(i,j)*recip_rhoConst       &                  +( MIN(zero,ddRloc)*ratioRm
109       &          + Bo_surf(i,j,bi,bj)*etaN(i,j,bi,bj)       &                    +MAX(zero,ddRloc)*ratioRp
110       &          + phi0surf(i,j,bi,bj)       &                   )*gravity*alphRho(i,j)*recip_rhoConst
111             ENDIF             ENDIF
112            ENDDO            ENDDO
113           ENDDO           ENDDO
114    
115  C  --  end if integr_GeoPot = ...  C  --  end if integr_GeoPot = ...
116         ENDIF         ENDIF
           
 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  
       ELSEIF ( buoyancyRelation .eq. 'OCEANICP' ) THEN  
   
        IF (integr_GeoPot.EQ.1) THEN  
 C  --  Finite Volume Form  
   
          DO j=jMin,jMax  
           DO i=iMin,iMax  
            IF ( K .EQ. kLowC(i,j,bi,bj) ) THEN  
             phiHydLow(i,j,bi,bj) = phiHyd(i,j,k)  
      &          + hFacC(i,j,k,bi,bj)*drF(K)*alphRho(i,j)  
      &          + Bo_surf(i,j,bi,bj)*etaN(i,j,bi,bj)  
      &          + phi0surf(i,j,bi,bj)  
            ENDIF  
           ENDDO  
          ENDDO  
117    
118         ELSE  C  -- end usingZCoords
119  C  --  Finite Difference Form        ENDIF
120    
121  C---------- Compute gravity*(sea surface elevation) first        IF ( k.EQ.Nr ) THEN
122  C           This has to be done starting from phiHyd at the current  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
123  C           tracer point and .5 of the cell's thickness has to be  C  --  last level (bottom): rescale (r*) and add surface contribution
 C           substracted from hFacC  
124    
125           DO j=jMin,jMax         IF ( usingPCoords ) THEN
126            DO i=iMin,iMax  C  -- P coordinate : Phi(R_low) is simply at the top :
127             IF ( K .EQ. kLowC(i,j,bi,bj) ) THEN          DO j=jMin,jMax
128              phiHydLow(i,j,bi,bj) = phiHyd(i,j,k)           DO i=iMin,iMax
129       &          + ( half*dRloc+(hFacC(i,j,k,bi,bj)-half)*drF(k)             phiHydLow(i,j,bi,bj) = phiHydF(i,j)
      &            )*alphRho(i,j)  
      &          + Bo_surf(i,j,bi,bj)*etaN(i,j,bi,bj)  
      &          + phi0surf(i,j,bi,bj)  
            ENDIF  
           ENDDO  
130           ENDDO           ENDDO
131            ENDDO
132           ENDIF
133    
134  C  --  end if integr_GeoPot = ...  #ifdef NONLIN_FRSURF
135    c      IF ( select_rStar.GE.2 .AND. nonlinFreeSurf.GE.4 ) THEN
136           IF ( select_rStar.GT.0 .AND. nonlinFreeSurf.GE.4 ) THEN
137            DO j=jMin,jMax
138             DO i=iMin,iMax
139               phiHydLow(i,j,bi,bj) = phiHydLow(i,j,bi,bj)
140         &                          *rStarFacC(i,j,bi,bj)
141             ENDDO
142            ENDDO
143         ENDIF         ENDIF
144    #endif /* NONLIN_FRSURF */
145    
146            DO j=jMin,jMax
147             DO i=iMin,iMax
148               phiHydLow(i,j,bi,bj) = phiHydLow(i,j,bi,bj)
149         &            + Bo_surf(i,j,bi,bj)*etaN(i,j,bi,bj)
150         &            + phi0surf(i,j,bi,bj)
151             ENDDO
152            ENDDO
153    
 c     ELSEIF ( buoyancyRelation .eq. 'ATMOSPHERIC' ) THEN  
154  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
155    C  -- end if k=Nr.
156        ENDIF        ENDIF
157    
158  #endif /* INCLUDE_PHIHYD_CALCULATION_CODE */  #endif /* INCLUDE_PHIHYD_CALCULATION_CODE */
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.6
  ViewVC Help
Powered by ViewVC 1.1.22