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

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

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

Revision 1.9 - (hide annotations) (download)
Tue Jul 19 23:49:34 2011 UTC (12 years, 10 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint63g, checkpoint64, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint64a, checkpoint64b
Changes since 1.8: +41 -11 lines 
- include r* effect on reference-state (PhiRef) in totPhiHyd
  (so that it's really the potential anomaly at the cell center);
  affect solution only if z* with full pressure in EOS.
- keep the previous diagnostics as "PHIHYDcR" (closer to potential anomaly
  @ constant r) until interpolation at constant r is implemented.
1 jmc 1.9 C $Header: /u/gcmpack/MITgcm/model/src/diags_phi_rlow.F,v 1.8 2011/06/10 20:55:51 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "CPP_OPTIONS.h"
5    
6     CBOP
7     C !ROUTINE: DIAGS_PHI_RLOW
8     C !INTERFACE:
9 jmc 1.5 SUBROUTINE DIAGS_PHI_RLOW(
10 jmc 1.1 I k, bi, bj, iMin,iMax, jMin,jMax,
11 jmc 1.2 I phiHydF, phiHydC, alphRho, tFld, sFld,
12 jmc 1.1 I myTime, myIter, myThid)
13     C !DESCRIPTION: \bv
14     C *==========================================================*
15 jmc 1.5 C | S/R DIAGS_PHI_RLOW
16 jmc 1.1 C | o Diagnose Phi-Hydrostatic at r-lower boundary
17 jmc 1.5 C | = bottom pressure (ocean in z-coord) ;
18 jmc 1.1 C | = sea surface elevation (ocean in p-coord) ;
19     C | = height at the top of atmosphere (in p-coord) ;
20     C *==========================================================*
21     C \ev
22    
23     C !USES:
24     IMPLICIT NONE
25     C == Global variables ==
26     #include "SIZE.h"
27     #include "EEPARAMS.h"
28     #include "PARAMS.h"
29     #include "GRID.h"
30     #include "SURFACE.h"
31     #include "DYNVARS.h"
32    
33     C !INPUT/OUTPUT PARAMETERS:
34     C == Routine Arguments ==
35 jmc 1.5 C bi,bj :: tile index
36 jmc 1.1 C iMin,iMax,jMin,jMax :: Loop counters
37 jmc 1.5 C phiHydF :: hydrostatic potential anomaly at middle between
38 jmc 1.2 C 2 centers k & k+1 (interface k+1)
39     C phiHydC :: hydrostatic potential anomaly at cell center
40 jmc 1.1 C (atmos: =Geopotential ; ocean-z: =Pressure/rho)
41     C alphRho :: Density (z-coord) or specific volume (p-coord)
42     C tFld :: Potential temp.
43 jmc 1.5 C sFld :: Salinity
44     C myTime :: Current time
45     C myIter :: Current iteration number
46     C myThid :: my Thread Id number
47 jmc 1.1 INTEGER k, bi,bj, iMin,iMax, jMin,jMax
48 jmc 1.2 _RL phiHydF(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
49     _RL phiHydC(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
50 jmc 1.1 _RL alphRho(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
51     _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)
53     _RL myTime
54     INTEGER myIter, myThid
55    
56     #ifdef INCLUDE_PHIHYD_CALCULATION_CODE
57    
58     C !LOCAL VARIABLES:
59     C == Local variables ==
60     C i,j :: Loop counters
61     INTEGER i,j
62     _RL zero, one, half
63 jmc 1.2 _RL ddRloc, ratioRm, ratioRp
64 jmc 1.1 PARAMETER ( zero= 0. _d 0 , one= 1. _d 0 , half= .5 _d 0 )
65 jmc 1.9 #ifdef NONLIN_FRSURF
66     _RL facP, dPhiRef
67     #endif /* NONLIN_FRSURF */
68 jmc 1.1 CEOP
69    
70 jmc 1.6 IF ( usingZCoords ) THEN
71 jmc 1.1
72 jmc 1.2 C----- Compute bottom pressure deviation from gravity*rho0*H
73 jmc 1.7 C Start from phiHyd at the (bottom) tracer point and add Del_h*g*rho_prime
74 jmc 1.2 C with Del_h = distance from the bottom up to tracer point
75 jmc 1.1
76 jmc 1.6 C-- Initialise to zero (otherwise phi0surf accumulate over land)
77     IF ( k.EQ.1 ) THEN
78     DO j=1-Oly,sNy+Oly
79     DO i=1-Olx,sNx+Olx
80     phiHydLow(i,j,bi,bj) = 0. _d 0
81     ENDDO
82     ENDDO
83     ENDIF
84    
85 jmc 1.1 IF (integr_GeoPot.EQ.1) THEN
86     C -- Finite Volume Form
87    
88     DO j=jMin,jMax
89     DO i=iMin,iMax
90     IF ( k .EQ. kLowC(i,j,bi,bj) ) THEN
91 jmc 1.2 ddRloc = rC(k)-R_low(i,j,bi,bj)
92     phiHydLow(i,j,bi,bj) = phiHydC(i,j)
93     & + ddRloc*gravity*alphRho(i,j)*recip_rhoConst
94 jmc 1.1 ENDIF
95     ENDDO
96     ENDDO
97    
98     ELSE
99     C -- Finite Difference Form
100    
101 jmc 1.2 ratioRm = one
102     ratioRp = one
103     IF (k.GT.1 ) ratioRm = half*drC(k)/(rF(k)-rC(k))
104     IF (k.LT.Nr) ratioRp = half*drC(k+1)/(rC(k)-rF(k+1))
105 jmc 1.1
106     DO j=jMin,jMax
107     DO i=iMin,iMax
108 jmc 1.2 IF ( k .EQ. kLowC(i,j,bi,bj) ) THEN
109     ddRloc = rC(k)-R_low(i,j,bi,bj)
110     phiHydLow(i,j,bi,bj) = phiHydC(i,j)
111     & +( MIN(zero,ddRloc)*ratioRm
112     & +MAX(zero,ddRloc)*ratioRp
113     & )*gravity*alphRho(i,j)*recip_rhoConst
114 jmc 1.1 ENDIF
115     ENDDO
116     ENDDO
117    
118     C -- end if integr_GeoPot = ...
119     ENDIF
120 jmc 1.5
121 jmc 1.6 C -- end usingZCoords
122 jmc 1.2 ENDIF
123    
124 jmc 1.6 IF ( k.EQ.Nr ) THEN
125 jmc 1.1 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
126 jmc 1.2 C -- last level (bottom): rescale (r*) and add surface contribution
127 jmc 1.1
128 jmc 1.6 IF ( usingPCoords ) THEN
129 jmc 1.2 C -- P coordinate : Phi(R_low) is simply at the top :
130     DO j=jMin,jMax
131     DO i=iMin,iMax
132     phiHydLow(i,j,bi,bj) = phiHydF(i,j)
133 jmc 1.1 ENDDO
134 jmc 1.2 ENDDO
135     ENDIF
136 jmc 1.1
137 jmc 1.5 #ifdef NONLIN_FRSURF
138     c IF ( select_rStar.GE.2 .AND. nonlinFreeSurf.GE.4 ) THEN
139 jmc 1.8 IF ( select_rStar.GE.1 .AND. nonlinFreeSurf.GE.4 ) THEN
140     C- Integral of b.dr = rStarFac * Integral of b.dr* :
141 jmc 1.9 IF ( fluidIsAir ) THEN
142 jmc 1.8 C- Consistent with Phi'= Integr[ theta'.dPi ] :
143     DO j=jMin,jMax
144     DO i=iMin,iMax
145 jmc 1.9 facP = rStarFacC(i,j,bi,bj)**atm_kappa
146     dPhiRef = phiRef(2*k+1) - gravity*topoZ(i,j,bi,bj)
147     & - phi0surf(i,j,bi,bj)
148     phiHydLow(i,j,bi,bj) =
149     & phiHydLow(i,j,bi,bj)*facP
150     & + MAX( dPhiRef, zero )*( facP - 1. _d 0 )
151     & + phi0surf(i,j,bi,bj)
152     ENDDO
153     ENDDO
154     ELSEIF ( usingPCoords ) THEN
155     C- Note: dPhiRef*(rStarFacC -1) = 1/rho*PSo*((eta+PSo)/PSo -1) = Bo_surf*etaN
156     C so that this expression is the same as before (same as in "ELSE" block below)
157     DO j=jMin,jMax
158     DO i=iMin,iMax
159     dPhiRef = ( Ro_surf(i,j,bi,bj) - R_low(i,j,bi,bj) )
160     & *recip_rhoConst
161     phiHydLow(i,j,bi,bj) =
162     & phiHydLow(i,j,bi,bj)*rStarFacC(i,j,bi,bj)
163     & + dPhiRef*( rStarFacC(i,j,bi,bj) - 1. _d 0 )
164     & + phi0surf(i,j,bi,bj)
165 jmc 1.8 ENDDO
166     ENDDO
167     ELSE
168 jmc 1.9 C- Note: dPhiRef*(rStarFacC -1) = g*H*((eta+H)/H -1) = Bo_surf*etaN
169     C so that this expression is the same as before (same as in "ELSE" block below)
170 jmc 1.8 DO j=jMin,jMax
171     DO i=iMin,iMax
172 jmc 1.9 dPhiRef = ( Ro_surf(i,j,bi,bj)-R_low(i,j,bi,bj) )
173     & *gravity
174     phiHydLow(i,j,bi,bj) =
175     & phiHydLow(i,j,bi,bj)*rStarFacC(i,j,bi,bj)
176     & + dPhiRef*( rStarFacC(i,j,bi,bj) - 1. _d 0 )
177     & + phi0surf(i,j,bi,bj)
178 jmc 1.8 ENDDO
179 jmc 1.5 ENDDO
180 jmc 1.8 ENDIF
181 jmc 1.9 ELSE
182     #else /* NONLIN_FRSURF */
183     IF ( .TRUE. ) THEN
184 jmc 1.5 #endif /* NONLIN_FRSURF */
185    
186 jmc 1.9 DO j=jMin,jMax
187     DO i=iMin,iMax
188 jmc 1.2 phiHydLow(i,j,bi,bj) = phiHydLow(i,j,bi,bj)
189     & + Bo_surf(i,j,bi,bj)*etaN(i,j,bi,bj)
190 jmc 1.5 & + phi0surf(i,j,bi,bj)
191 jmc 1.9 ENDDO
192     ENDDO
193     ENDIF
194 jmc 1.1
195     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
196 jmc 1.2 C -- end if k=Nr.
197 jmc 1.1 ENDIF
198    
199     #endif /* INCLUDE_PHIHYD_CALCULATION_CODE */
200    
201     RETURN
202     END
  ViewVC Help
Powered by ViewVC 1.1.22