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

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

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


Revision 1.10 - (hide annotations) (download)
Fri Nov 4 01:19:24 2005 UTC (18 years, 7 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint57y_post, checkpoint57y_pre, checkpoint57x_post
Changes since 1.9: +1 -2 lines
remove unused variables (reduces number of compiler warning)

1 jmc 1.10 C $Header: /u/gcmpack/MITgcm/model/src/ini_linear_phisurf.F,v 1.9 2003/02/26 03:11:32 jmc Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "CPP_OPTIONS.h"
5    
6 cnh 1.3 CBOP
7     C !ROUTINE: INI_LINEAR_PHISURF
8     C !INTERFACE:
9 jmc 1.1 SUBROUTINE INI_LINEAR_PHISURF( myThid )
10 cnh 1.3
11     C !DESCRIPTION: \bv
12     C *==========================================================*
13     C | SUBROUTINE INI_LINEAR_PHISURF
14     C | o Initialise the Linear Relation Phi_surf(eta)
15     C *==========================================================*
16 jmc 1.7 C | Initialise -Boyancy at surface level (Bo_surf)
17 cnh 1.3 C | to setup the Linear relation: Phi_surf(eta)=Bo_surf*eta
18 jmc 1.7 C | Initialise phi0surf = starting point for integrating
19     C | phiHyd (= phiHyd at r=RoSurf)
20 cnh 1.3 C *==========================================================*
21     C \ev
22    
23     C !USES:
24 jmc 1.1 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    
32 cnh 1.3 C !INPUT/OUTPUT PARAMETERS:
33 jmc 1.1 C === Routine arguments ===
34     C myThid - Thread no. that called this routine.
35     INTEGER myThid
36    
37 jmc 1.9 C == Local variables in common ==
38     C Hloc - Temporary array used to write surface topography
39     C has to be in common for multi threading
40     COMMON / LOCAL_INI_PHISURF / topoHloc
41     _RS topoHloc(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
42    
43 cnh 1.3 C !LOCAL VARIABLES:
44 jmc 1.1 C === Local variables ===
45     C bi,bj - Loop counters
46     C I,J,K
47     INTEGER bi, bj
48     INTEGER I, J, K
49 mlosch 1.4 _RL rhoLoc
50 jmc 1.1 _RL dPIdp
51 cnh 1.3 CEOP
52 heimbach 1.6
53     #ifdef ALLOW_AUTODIFF_TAMC
54     DO bj=myByLo(myThid),myByHi(myThid)
55     DO bi=myBxLo(myThid),myBxHi(myThid)
56     DO J=1-Oly,sNy+Oly
57     DO I=1-Olx,sNx+Olx
58     Bo_surf(I,J,bi,bj) = 0. _d 0
59     recip_Bo(I,J,bi,bj) = 0. _d 0
60     ENDDO
61     ENDDO
62     ENDDO
63     ENDDO
64     #endif
65 jmc 1.1
66     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
67    
68     C-- Initialise -Boyancy at surface level : Bo_surf
69     C Bo_surf is defined as d/dr(Phi_surf) and set to g/rtoz (linear free surface)
70     C with rtoz = conversion factor from r-unit to z-unit (=horiVertRatio)
71     C an accurate formulation includes P_surf and T,S_surf effects on rho_surf:
72     C (setting uniformLin_PhiSurf=.FALSE.):
73     C z-ocean (rtoz=1) : Bo_surf = - Boyancy = gravity * rho_surf/rho_0
74     C p-atmos (rtoz=rho_c*g) : Bo_surf = (1/rho)_surf
75     C Note on Phi_surf splitting : Non-linear Time-dependent effects on b_surf
76     C [through eta & (T-tRef)_surf] are included in PhiHyd rather than in Bo_surf
77     C--
78     IF ( buoyancyRelation .eq. 'OCEANIC' ) THEN
79     C- gBaro = gravity (except for External mode test with reduced gravity)
80     DO bj=myByLo(myThid),myByHi(myThid)
81     DO bi=myBxLo(myThid),myBxHi(myThid)
82     DO J=1-Oly,sNy+Oly
83     DO I=1-Olx,sNx+Olx
84     Bo_surf(I,J,bi,bj) = gBaro
85     recip_Bo(I,J,bi,bj) = 1. _d 0 / gBaro
86     ENDDO
87     ENDDO
88     ENDDO
89     ENDDO
90     ELSEIF ( uniformLin_PhiSurf ) THEN
91     C- use a linear (in ps) uniform relation : Phi'_surf = 1/rhoConst * ps'_surf
92     DO bj=myByLo(myThid),myByHi(myThid)
93     DO bi=myBxLo(myThid),myBxHi(myThid)
94     DO J=1-Oly,sNy+Oly
95     DO I=1-Olx,sNx+Olx
96     Bo_surf(I,J,bi,bj) = recip_rhoConst
97     recip_Bo(I,J,bi,bj) = rhoConst
98     ENDDO
99     ENDDO
100     ENDDO
101     ENDDO
102 mlosch 1.4 ELSEIF ( buoyancyRelation .eq. 'OCEANICP' ) THEN
103     DO bj=myByLo(myThid),myByHi(myThid)
104     DO bi=myBxLo(myThid),myBxHi(myThid)
105     DO J=1-Oly,sNy+Oly
106     DO I=1-Olx,sNx+Olx
107     IF ( Ro_surf(I,J,bi,bj).GT.0. _d 0
108     & .AND. ksurfC(I,J,bi,bj).LE.Nr ) THEN
109     k = ksurfC(I,J,bi,bj)
110     CALL FIND_RHO_SCALAR(
111     & tRef(k), sRef(k), Ro_surf(I,J,bi,bj),
112     & rhoLoc, myThid )
113 mlosch 1.5 rhoLoc = rhoLoc + rhoConst
114     if ( rhoLoc .eq. 0. _d 0 ) then
115 mlosch 1.4 Bo_surf(I,J,bi,bj) = 0. _d 0
116     else
117 mlosch 1.5 Bo_surf(I,J,bi,bj) = 1./rhoLoc
118 mlosch 1.4 endif
119 mlosch 1.5 recip_Bo(I,J,bi,bj) = rhoLoc
120 mlosch 1.4 ELSE
121     Bo_surf(I,J,bi,bj) = 0. _d 0
122     recip_Bo(I,J,bi,bj) = 0. _d 0
123     ENDIF
124     ENDDO
125     ENDDO
126     ENDDO
127     ENDDO
128     ELSEIF ( buoyancyRelation .eq. 'ATMOSPHERIC' ) THEN
129 jmc 1.1 C- use a linearized (in ps) Non-uniform relation : Bo_surf(Po_surf,tRef_surf)
130     C--- Bo = d/d_p(Phi_surf) = tRef_surf*d/d_p(PI) ; PI = Cp*(p/Po)^kappa
131     DO bj=myByLo(myThid),myByHi(myThid)
132     DO bi=myBxLo(myThid),myBxHi(myThid)
133     DO J=1-Oly,sNy+Oly
134     DO I=1-Olx,sNx+Olx
135 jmc 1.2 IF ( Ro_surf(I,J,bi,bj).GT.0. _d 0
136     & .AND. ksurfC(I,J,bi,bj).LE.Nr ) THEN
137 jmc 1.7 dPIdp = (atm_Cp*atm_kappa/atm_Po)*
138     & (Ro_surf(I,J,bi,bj)/atm_Po)**(atm_kappa-1. _d 0)
139 jmc 1.2 Bo_surf(I,J,bi,bj) = dPIdp*tRef(ksurfC(I,J,bi,bj))
140 jmc 1.1 recip_Bo(I,J,bi,bj) = 1. _d 0 / Bo_surf(I,J,bi,bj)
141     ELSE
142     Bo_surf(I,J,bi,bj) = 0.
143     recip_Bo(I,J,bi,bj) = 0.
144     ENDIF
145     ENDDO
146     ENDDO
147     ENDDO
148     ENDDO
149 mlosch 1.4 ELSE
150     STOP 'INI_LINEAR_PHISURF: We should never reach this point!'
151 jmc 1.1 ENDIF
152    
153     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
154    
155     C-- Update overlap regions
156     _EXCH_XY_R8(Bo_surf, myThid)
157     _EXCH_XY_R8(recip_Bo, myThid)
158    
159 mlosch 1.4 IF ( ( buoyancyRelation .eq. 'ATMOSPHERIC' .OR.
160     & buoyancyRelation .eq. 'OCEANICP' )
161     & .AND. .NOT.uniformLin_PhiSurf ) THEN
162 jmc 1.7
163     _BEGIN_MASTER( myThid )
164 jmc 1.1 CALL WRITE_FLD_XY_RL( 'Bo_surf',' ',Bo_surf,0,myThid)
165 jmc 1.7 _END_MASTER( myThid )
166    
167 jmc 1.1 ENDIF
168    
169 jmc 1.7 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
170    
171     C-- Initialise phi0surf: used for atmos. surf. P-loading (ocean, z-coord)
172     C or topographic geopotential anom. (p-coord)
173    
174     DO bj=myByLo(myThid),myByHi(myThid)
175     DO bi=myBxLo(myThid),myBxHi(myThid)
176     DO J=1-Oly,sNy+Oly
177     DO I=1-Olx,sNx+Olx
178     phi0surf(I,J,bi,bj) = 0.
179     ENDDO
180     ENDDO
181     ENDDO
182     ENDDO
183    
184     IF ( buoyancyRelation .eq. 'ATMOSPHERIC'
185 jmc 1.9 & .AND. topoFile.NE.' ' ) THEN
186 jmc 1.7
187 heimbach 1.8 #ifdef ALLOW_AUTODIFF_TAMC
188     STOP 'CANNOT PRESENTLY USE THIS OPTION WITH ADJOINT'
189     #else
190    
191 jmc 1.9 C-- Compute topoH = PhiRef(Po_surf)/g ; is different from original
192     C topoZ(read from file) because of truncation of Po_surf.
193     C NOTE: not clear for now which topoZ needs to be saved in common block
194     C-- AND set phi0surf = starting point for integrating Geopotential;
195    
196     CALL INI_P_GROUND( -2,
197     O topoHloc,
198 jmc 1.7 I Ro_surf, myThid )
199    
200 jmc 1.9 _BEGIN_MASTER( myThid )
201     CALL WRITE_FLD_XY_RS( 'topo_H',' ',topoHloc,0,myThid)
202     _END_MASTER( myThid )
203    
204     IF (selectFindRoSurf.NE.0) THEN
205 jmc 1.7 _EXCH_XY_RS(phi0surf, myThid)
206    
207     _BEGIN_MASTER( myThid )
208     CALL WRITE_FLD_XY_RS( 'phi0surf',' ',phi0surf,0,myThid)
209     _END_MASTER( myThid )
210 jmc 1.9 ENDIF
211 heimbach 1.8
212     #endif /* ALLOW_AUTODIFF_TAMC */
213 jmc 1.7
214     ENDIF
215    
216     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
217 jmc 1.1 RETURN
218     END

  ViewVC Help
Powered by ViewVC 1.1.22