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

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

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


Revision 1.9 - (hide annotations) (download)
Tue Dec 5 05:25:08 2006 UTC (17 years, 5 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint58t_post, mitgcm_mapl_00, checkpoint58u_post, checkpoint58s_post
Changes since 1.8: +21 -22 lines
start to implement deep-atmosphere and/or anelastic formulation

1 jmc 1.9 C $Header: /u/gcmpack/MITgcm/model/src/calc_grad_phi_hyd.F,v 1.8 2005/12/08 15:44:33 heimbach Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "CPP_OPTIONS.h"
5    
6     CBOP
7     C !ROUTINE: CALC_GRAD_PHI_HYD
8     C !INTERFACE:
9 jmc 1.9 SUBROUTINE CALC_GRAD_PHI_HYD(
10 jmc 1.1 I k, bi, bj, iMin,iMax, jMin,jMax,
11 jmc 1.3 I phiHydC, alphRho, tFld, sFld,
12 jmc 1.1 O dPhiHydX, dPhiHydY,
13     I myTime, myIter, myThid)
14     C !DESCRIPTION: \bv
15     C *==========================================================*
16 jmc 1.9 C | S/R CALC_GRAD_PHI_HYD
17     C | o Calculate the gradient of Hydrostatic potential anomaly
18 jmc 1.1 C *==========================================================*
19     C \ev
20    
21     C !USES:
22     IMPLICIT NONE
23     C == Global variables ==
24     #include "SIZE.h"
25     #include "EEPARAMS.h"
26     #include "PARAMS.h"
27     #include "GRID.h"
28     #include "SURFACE.h"
29 jmc 1.3 #include "DYNVARS.h"
30 jmc 1.1
31     C !INPUT/OUTPUT PARAMETERS:
32     C == Routine Arguments ==
33 jmc 1.9 C bi,bj :: tile index
34 jmc 1.1 C iMin,iMax,jMin,jMax :: Loop counters
35 jmc 1.9 C phiHydC :: Hydrostatic Potential anomaly
36 jmc 1.1 C (atmos: =Geopotential ; ocean-z: =Pressure/rho)
37     C alphRho :: Density (z-coord) or specific volume (p-coord)
38     C tFld :: Potential temp.
39 jmc 1.9 C sFld :: Salinity
40 jmc 1.1 C dPhiHydX,Y :: Gradient (X & Y directions) of Hyd. Potential
41     C myTime :: Current time
42     C myIter :: Current iteration number
43     C myThid :: Instance number for this call of the routine.
44     INTEGER k, bi,bj, iMin,iMax, jMin,jMax
45 jmc 1.3 _RL phiHydC(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
46 jmc 1.1 _RL alphRho(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
47     _RL tFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
48     _RL sFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
49     _RL dPhiHydX(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
50     _RL dPhiHydY(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
51     _RL myTime
52     INTEGER myIter, myThid
53    
54     #ifdef INCLUDE_PHIHYD_CALCULATION_CODE
55    
56     C !LOCAL VARIABLES:
57     C == Local variables ==
58     C i,j :: Loop counters
59     INTEGER i,j
60     _RL varLoc(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
61 jmc 1.7 #ifdef NONLIN_FRSURF
62 jmc 1.3 _RL factorZ, factorP, conv_theta2T
63 jmc 1.5 _RL factPI
64     CHARACTER*(MAX_LEN_MBUF) msgBuf
65 jmc 1.7 #endif
66 jmc 1.1 CEOP
67    
68 jmc 1.3 #ifdef NONLIN_FRSURF
69 jmc 1.5 IF (select_rStar.GE.2 .AND. nonlinFreeSurf.GE.4 ) THEN
70 heimbach 1.8 # ifndef DISABLE_RSTAR_CODE
71 jmc 1.5 C- Integral of b.dr = rStarFac * Integral of b.dr* :
72 jmc 1.9 C and will add later (select_rStar=2) the contribution of
73 jmc 1.5 C the slope of the r* coordinate.
74     IF ( buoyancyRelation .EQ. 'ATMOSPHERIC' ) THEN
75     C- Consistent with Phi'= Integr[ theta'.dPi ] :
76     DO j=jMin-1,jMax
77     DO i=iMin-1,iMax
78     varLoc(i,j) = phiHydC(i,j)*rStarFacC(i,j,bi,bj)**atm_kappa
79     & + phi0surf(i,j,bi,bj)
80     ENDDO
81     ENDDO
82     ELSE
83     DO j=jMin-1,jMax
84     DO i=iMin-1,iMax
85     varLoc(i,j) = phiHydC(i,j)*rStarFacC(i,j,bi,bj)
86     & + phi0surf(i,j,bi,bj)
87     ENDDO
88     ENDDO
89     ENDIF
90     ELSEIF (select_rStar.GE.1 .AND. nonlinFreeSurf.GE.4 ) THEN
91     C- Integral of b.dr but scaled to correspond to a fixed r-level (=r*)
92     C no contribution of the slope of the r* coordinate (select_rStar=1)
93     IF ( buoyancyRelation .EQ. 'ATMOSPHERIC' ) THEN
94     C- Consistent with Phi'= Integr[ theta'.dPi ] :
95     DO j=jMin-1,jMax
96     DO i=iMin-1,iMax
97     IF (Ro_surf(i,j,bi,bj).EQ.rC(k)) THEN
98 jmc 1.9 factPI=atm_Cp*( ((etaH(i,j,bi,bj)+rC(k))/atm_Po)**atm_kappa
99     & -( rC(k) /atm_Po)**atm_kappa
100 jmc 1.5 & )
101     varLoc(i,j) = factPI*alphRho(i,j)
102     ELSEIF (Ro_surf(i,j,bi,bj).NE.0. _d 0) THEN
103     factPI = (rC(k)/Ro_surf(i,j,bi,bj))**atm_kappa
104     varLoc(i,j) = phiHydC(i,j)
105     & *(rStarFacC(i,j,bi,bj)**atm_kappa - factPI)
106     & /(1. _d 0 -factPI)
107     & + phi0surf(i,j,bi,bj)
108     ENDIF
109     ENDDO
110     ENDDO
111     ELSE
112     DO j=jMin-1,jMax
113     DO i=iMin-1,iMax
114     IF (Ro_surf(i,j,bi,bj).EQ.rC(k)) THEN
115 jmc 1.6 WRITE(msgBuf,'(3A)') 'CALC_GRAD_PHI_HYD: ',
116     & 'Problem when Ro_surf=rC',
117     & ' with select_rStar,integr_GeoPot=1,4'
118 jmc 1.5 CALL PRINT_ERROR( msgBuf , myThid)
119     STOP 'CALC_GRAD_PHI_HYD: Pb in r* options implementation'
120     ELSE
121     varLoc(i,j) = phiHydC(i,j)
122     & *(etaH(i,j,bi,bj)+Ro_surf(i,j,bi,bj)-rC(k))
123     & / (Ro_surf(i,j,bi,bj)-rC(k))
124     & + phi0surf(i,j,bi,bj)
125     ENDIF
126     ENDDO
127 jmc 1.3 ENDDO
128 jmc 1.5 ENDIF
129 heimbach 1.8 # endif /* DISABLE_RSTAR_CODE */
130 jmc 1.3 ELSE
131     #else /* NONLIN_FRSURF */
132     IF (.TRUE.) THEN
133     #endif /* NONLIN_FRSURF */
134     DO j=jMin-1,jMax
135     DO i=iMin-1,iMax
136     varLoc(i,j) = phiHydC(i,j)+phi0surf(i,j,bi,bj)
137     ENDDO
138     ENDDO
139     ENDIF
140 jmc 1.1
141 jmc 1.3 C-- Zonal & Meridional gradient of potential anomaly
142 jmc 1.1 DO j=jMin,jMax
143     DO i=iMin,iMax
144 jmc 1.9 dPhiHydX(i,j) = _recip_dxC(i,j,bi,bj)*recip_deepFacC(k)
145     & *( varLoc(i,j)-varLoc(i-1,j) )*recip_rhoFacC(k)
146     dPhiHydY(i,j) = _recip_dyC(i,j,bi,bj)*recip_deepFacC(k)
147     & *( varLoc(i,j)-varLoc(i,j-1) )*recip_rhoFacC(k)
148 jmc 1.1 ENDDO
149     ENDDO
150 jmc 1.3
151     #ifdef NONLIN_FRSURF
152     IF (select_rStar.GE.2 .AND. nonlinFreeSurf.GE.1 ) THEN
153     IF ( buoyancyRelation .EQ. 'OCEANIC' ) THEN
154     C-- z* coordinate slope term: rho'/rho0 * Grad_r(g.z)
155 jmc 1.9 factorZ = gravity*recip_rhoConst*recip_rhoFacC(k)*0.5 _d 0
156 jmc 1.3 DO j=jMin-1,jMax
157     DO i=iMin-1,iMax
158     varLoc(i,j) = etaH(i,j,bi,bj)
159     & *(1. _d 0 + rC(k)*recip_Rcol(i,j,bi,bj))
160     ENDDO
161     ENDDO
162     DO j=jMin,jMax
163     DO i=iMin,iMax
164     dPhiHydX(i,j) = dPhiHydX(i,j)
165     & +factorZ*(alphRho(i-1,j)+alphRho(i,j))
166     & *(varLoc(i,j)-varLoc(i-1,j))
167 jmc 1.9 & *recip_dxC(i,j,bi,bj)*recip_deepFacC(k)
168 jmc 1.3 dPhiHydY(i,j) = dPhiHydY(i,j)
169     & +factorZ*(alphRho(i,j-1)+alphRho(i,j))
170     & *(varLoc(i,j)-varLoc(i,j-1))
171 jmc 1.9 & *recip_dyC(i,j,bi,bj)*recip_deepFacC(k)
172 jmc 1.3 ENDDO
173     ENDDO
174     ELSEIF (buoyancyRelation .EQ. 'OCEANICP' ) THEN
175     C-- p* coordinate slope term: alpha' * Grad_r( p )
176     factorP = 0.5 _d 0
177     DO j=jMin,jMax
178     DO i=iMin,iMax
179     dPhiHydX(i,j) = dPhiHydX(i,j)
180     & +factorP*(alphRho(i-1,j)+alphRho(i,j))
181     & *(rStarFacC(i,j,bi,bj)-rStarFacC(i-1,j,bi,bj))
182 jmc 1.9 & *rC(k)*recip_dxC(i,j,bi,bj)*recip_deepFacC(k)
183 jmc 1.3 dPhiHydY(i,j) = dPhiHydY(i,j)
184     & +factorP*(alphRho(i,j-1)+alphRho(i,j))
185     & *(rStarFacC(i,j,bi,bj)-rStarFacC(i,j-1,bi,bj))
186 jmc 1.9 & *rC(k)*recip_dyC(i,j,bi,bj)*recip_deepFacC(k)
187 jmc 1.3 ENDDO
188     ENDDO
189     ELSEIF ( buoyancyRelation .EQ. 'ATMOSPHERIC' ) THEN
190     C-- p* coordinate slope term: alpha' * Grad_r( p )
191     conv_theta2T = (rC(k)/atm_Po)**atm_kappa
192     factorP = (atm_Rd/rC(k))*conv_theta2T*0.5 _d 0
193     DO j=jMin,jMax
194     DO i=iMin,iMax
195     dPhiHydX(i,j) = dPhiHydX(i,j)
196 jmc 1.5 & +factorP*(alphRho(i-1,j)+alphRho(i,j))
197 jmc 1.3 & *(rStarFacC(i,j,bi,bj)-rStarFacC(i-1,j,bi,bj))
198 jmc 1.9 & *rC(k)*recip_dxC(i,j,bi,bj)*recip_deepFacC(k)
199 jmc 1.3 dPhiHydY(i,j) = dPhiHydY(i,j)
200 jmc 1.5 & +factorP*(alphRho(i,j-1)+alphRho(i,j))
201 jmc 1.3 & *(rStarFacC(i,j,bi,bj)-rStarFacC(i,j-1,bi,bj))
202 jmc 1.9 & *rC(k)*recip_dyC(i,j,bi,bj)*recip_deepFacC(k)
203 jmc 1.3 ENDDO
204     ENDDO
205     ENDIF
206     ENDIF
207     #endif /* NONLIN_FRSURF */
208 jmc 1.1
209     #endif /* INCLUDE_PHIHYD_CALCULATION_CODE */
210    
211     RETURN
212     END

  ViewVC Help
Powered by ViewVC 1.1.22