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

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

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


Revision 1.25 - (hide annotations) (download)
Mon Jul 9 16:46:29 2001 UTC (22 years, 11 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint40pre3, checkpoint40pre7, checkpoint40pre6, checkpoint40pre9, checkpoint40pre8, checkpoint40pre2, checkpoint40pre4, checkpoint40pre5, checkpoint40
Changes since 1.24: +30 -48 lines
exchange R_low and Ro_surf to get the same mask as before (cube corners).

1 jmc 1.25 C $Header: /u/gcmpack/models/MITgcmUV/model/src/ini_depths.F,v 1.24 2001/07/06 21:48:31 jmc Exp $
2 adcroft 1.23 C $Name: $
3 cnh 1.1
4 cnh 1.16 #include "CPP_OPTIONS.h"
5 cnh 1.1
6     CStartOfInterface
7     SUBROUTINE INI_DEPTHS( myThid )
8     C /==========================================================\
9     C | SUBROUTINE INI_DEPTHS |
10 jmc 1.24 C | o define R_position of Lower and Surface Boundaries |
11 cnh 1.1 C |==========================================================|
12 jmc 1.24 C |atmosphere orography: |
13     C | define either in term of P_topo or converted from Z_topo |
14     C |ocean bathymetry: |
15 cnh 1.1 C | The depths of the bottom of the model is specified in |
16     C | terms of an XY map with one depth for each column of |
17     C | grid cells. Depths do not have to coincide with the |
18 cnh 1.15 C | model levels. The model lopping algorithm makes it |
19 cnh 1.1 C | possible to represent arbitrary depths. |
20     C | The mode depths map also influences the models topology |
21     C | By default the model domain wraps around in X and Y. |
22     C | This default doubly periodic topology is "supressed" |
23     C | if a depth map is defined which closes off all wrap |
24     C | around flow. |
25     C \==========================================================/
26 adcroft 1.17 IMPLICIT NONE
27 cnh 1.1
28     C === Global variables ===
29     #include "SIZE.h"
30     #include "EEPARAMS.h"
31     #include "PARAMS.h"
32     #include "GRID.h"
33    
34     C == Routine arguments ==
35     C myThid - Number of this instance of INI_DEPTHS
36     INTEGER myThid
37     CEndOfInterface
38    
39 jmc 1.24 C == Local variables in common ==
40     C Hloc - Temporary array used to read surface topography
41     C has to be in common for multi threading
42     COMMON / LOCAL_INI_DEPTHS / Hloc
43     _RS Hloc(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
44    
45 cnh 1.1 C == Local variables ==
46     C iG, jG - Global coordinate index
47     C bi,bj - Loop counters
48 adcroft 1.21 C I,J,K
49 cnh 1.13 C oldPrec - Temporary used in controlling binary input dataset precision
50 jmc 1.24 C msgBuf - Informational/error meesage buffer
51 cnh 1.1 INTEGER iG, jG
52     INTEGER bi, bj
53 adcroft 1.23 INTEGER I, J
54 jmc 1.24 CHARACTER*(MAX_LEN_MBUF) msgBuf
55    
56     IF (groundAtK1 .AND. bathyFile .NE. ' '
57     & .AND. topoFile .NE. ' ' ) THEN
58     WRITE(msgBuf,'(A,A)')
59     & 'S/R INI_DEPTHS: both bathyFile & topoFile are specified:',
60     & ' select the right one !'
61     CALL PRINT_ERROR( msgBuf , myThid)
62     STOP 'ABNORMAL END: S/R INI_DEPTHS'
63     ENDIF
64    
65     C------
66     C 0) Initialize R_low and Ro_surf (define an empty domain)
67     C------
68     DO bj = myByLo(myThid), myByHi(myThid)
69     DO bi = myBxLo(myThid), myBxHi(myThid)
70     DO j=1-Oly,sNy+Oly
71     DO i=1-Olx,sNx+Olx
72     R_low(i,j,bi,bj) = 0.
73     Ro_surf(i,j,bi,bj) = 0.
74     ENDDO
75     ENDDO
76     ENDDO
77     ENDDO
78 cnh 1.1
79 jmc 1.24 C------
80     C 1) Set R_low = the Lower (in r sense) boundary of the fluid column :
81     C------
82 jmc 1.25 IF (groundAtK1 .OR. bathyFile .EQ. ' ') THEN
83 jmc 1.24 C- e.g., atmosphere : R_low = Top of atmosphere
84 jmc 1.25 C- ocean : R_low = Bottom
85 cnh 1.4 DO bj = myByLo(myThid), myByHi(myThid)
86     DO bi = myBxLo(myThid), myBxHi(myThid)
87 jmc 1.25 DO j=1,sNy
88     DO i=1,sNx
89 jmc 1.24 R_low(i,j,bi,bj) = rF(Nr+1)
90 cnh 1.4 ENDDO
91 cnh 1.1 ENDDO
92     ENDDO
93     ENDDO
94 cnh 1.4 ELSE
95 jmc 1.24 _BEGIN_MASTER( myThid )
96 adcroft 1.19 C Read the bathymetry using the mid-level I/O pacakage read_write_rec
97     C The 0 is the "iteration" argument. The 1 is the record number.
98 jmc 1.24 CALL READ_REC_XY_RS( bathyFile, R_low, 1, 0, myThid )
99 adcroft 1.19 C Read the bathymetry using the mid-level I/O pacakage read_write_fld
100     C The 0 is the "iteration" argument. The ' ' is an empty suffix
101 jmc 1.24 c CALL READ_FLD_XY_RS( bathyFile, ' ', R_low, 0, myThid )
102 adcroft 1.19 C Read the bathymetry using the low-level I/O package
103 jmc 1.24 c CALL MDSREADFIELD( bathyFile, readBinaryPrec,
104     c & 'RS', 1, R_low, 1, myThid )
105     _END_MASTER(myThid)
106    
107 jmc 1.25 ENDIF
108     C- end setup R_low in the interior
109 jmc 1.24
110 jmc 1.25 C- fill in the overlap :
111     _EXCH_XY_R4(R_low, myThid )
112 adcroft 1.21
113 jmc 1.24 c CALL PLOT_FIELD_XYRS(R_low,'Bottom depths (ini_depths)',1,myThid)
114     c _BEGIN_MASTER( myThid )
115     c CALL WRITE_FLD_XY_RS( 'R_low' ,' ', R_low, 0,myThid)
116     c _END_MASTER(myThid)
117    
118     c---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
119    
120     C------
121     C 2) Set R_surf = Surface boundary: ocean surface / ground for the atmosphere
122     C------
123    
124 jmc 1.25 IF ( groundAtK1 .AND. bathyFile.NE.' ' ) THEN
125     C------ read directly Po_surf from bathyFile (only for backward compatibility)
126    
127     _BEGIN_MASTER( myThid )
128     CALL READ_REC_XY_RS( bathyFile, Ro_surf, 1, 0, myThid )
129     _END_MASTER(myThid)
130     _BARRIER
131    
132     ELSEIF ( topoFile.EQ.' ' ) THEN
133     C------ set default value:
134 jmc 1.24
135     DO bj = myByLo(myThid), myByHi(myThid)
136     DO bi = myBxLo(myThid), myBxHi(myThid)
137 jmc 1.25 DO j=1,sNy
138     DO i=1,sNx
139 jmc 1.24 Ro_surf(i,j,bi,bj) = Ro_SeaLevel
140     ENDDO
141     ENDDO
142     ENDDO
143     ENDDO
144    
145     ELSE
146 jmc 1.25 C------ read from file:
147 jmc 1.24
148     C- read surface topography (in m) from topoFile (case topoFile.NE.' '):
149     _BEGIN_MASTER( myThid )
150     CALL READ_REC_XY_RS( topoFile, Hloc, 1, 0, myThid )
151     _END_MASTER(myThid)
152     _BARRIER
153    
154     IF (buoyancyRelation .EQ. 'ATMOSPHERIC') THEN
155 jmc 1.25 C----
156 jmc 1.24 C Convert Surface Geopotential to (reference) Surface Pressure
157     C according to Tref profile, using same discretisation as in calc_phi_hyd
158 jmc 1.25 C----
159 jmc 1.24 c _BEGIN_MASTER( myThid )
160     c CALL WRITE_FLD_XY_RS( 'topo_Z',' ',Hloc,0,myThid)
161     c _END_MASTER(myThid)
162    
163     CALL INI_P_GROUND( Hloc, Ro_surf, myThid )
164    
165     _BARRIER
166     _BEGIN_MASTER( myThid )
167     CALL WRITE_FLD_XY_RS( 'topo_P',' ',Ro_surf,0,myThid)
168     _END_MASTER(myThid)
169    
170     ELSE
171 jmc 1.25 C----
172 jmc 1.24 C Direct Transfer to Ro_surf :
173     DO bj = myByLo(myThid), myByHi(myThid)
174     DO bi = myBxLo(myThid), myBxHi(myThid)
175     DO j=1,sNy
176     DO i=1,sNx
177     Ro_surf(i,j,bi,bj) = Hloc(i,j,bi,bj)
178     ENDDO
179     ENDDO
180     ENDDO
181     ENDDO
182 jmc 1.25
183 jmc 1.24 ENDIF
184    
185 jmc 1.25 C------ end case "read topoFile"
186     ENDIF
187 jmc 1.24
188 jmc 1.25 C----- fill in the overlap :
189     _EXCH_XY_R4(Ro_surf, myThid )
190 jmc 1.24
191     c---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
192 adcroft 1.21
193 jmc 1.24 C------
194     C 3) Close the Domain (special configuration).
195     C------
196 adcroft 1.23 IF (groundAtK1) THEN
197 adcroft 1.21 DO bj = myByLo(myThid), myByHi(myThid)
198     DO bi = myBxLo(myThid), myBxHi(myThid)
199     DO j=1-Oly,sNy+Oly
200     DO i=1-Olx,sNx+Olx
201 jmc 1.24 iG = myXGlobalLo-1+(bi-1)*sNx+I
202     jG = myYGlobalLo-1+(bj-1)*sNy+J
203     C Test for eastern edge
204     c IF ( iG .EQ. Nx ) Ro_surf(i,j,bi,bj) = 0.
205     C Test for northern edge
206     c IF ( jG .EQ. Ny ) Ro_surf(i,j,bi,bj) = 0.
207 adcroft 1.23 IF (usingSphericalPolarGrid .AND. abs(yC(I,J,bi,bj)).GE.90. )
208 jmc 1.24 & Ro_surf(I,J,bi,bj) = rF(Nr+1)
209 adcroft 1.23 ENDDO
210     ENDDO
211     ENDDO
212     ENDDO
213     ELSE
214     DO bj = myByLo(myThid), myByHi(myThid)
215     DO bi = myBxLo(myThid), myBxHi(myThid)
216     DO j=1-Oly,sNy+Oly
217     DO i=1-Olx,sNx+Olx
218 jmc 1.24 iG = myXGlobalLo-1+(bi-1)*sNx+I
219     jG = myYGlobalLo-1+(bj-1)*sNy+J
220     C Test for eastern edge
221     c IF ( iG .EQ. Nx ) R_low(i,j,bi,bj) = 0.
222     C Test for northern edge
223     c IF ( jG .EQ. Ny ) R_low(i,j,bi,bj) = 0.
224 adcroft 1.23 IF (usingSphericalPolarGrid .AND. abs(yC(I,J,bi,bj)).GE.90. )
225 jmc 1.24 & R_low(I,J,bi,bj) = Ro_SeaLevel
226 adcroft 1.21 ENDDO
227     ENDDO
228     ENDDO
229     ENDDO
230     ENDIF
231 jmc 1.24
232     c _BEGIN_MASTER( myThid )
233     c CALL WRITE_FLD_XY_RS('Ro_surf',' ',Ro_surf,0,myThid)
234     c _END_MASTER(myThid)
235 adcroft 1.23
236 cnh 1.1 RETURN
237     END

  ViewVC Help
Powered by ViewVC 1.1.22